Re: iso curve extraction

On Jul 22, 9:11 am, "Wade Ward" <zaxf...@xxxxxxxxxxx> wrote:
"veda" <kedar.hardi...@xxxxxxxxx> wrote in message

news:1185061741.160634.283180@xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx





On Jul 22, 4:02 am, fj <francois.j...@xxxxxxx> wrote:
On 22 juil, 01:33, glen herrmannsfeldt <g...@xxxxxxxxxxxxxxxx> wrote:

veda wrote:

SUBROUTINE TITQ4N (X ,Y ,F ,V )
C&F
C&F Drawing of isovalues in a Q4 (node field)
C&F
C&A X : abscissa
C&A Y : ordinates
C&A F : field at nodes
C&A V : wanted value
C&A
DIMENSION X(4),Y(4),F(4),FV(4),A(2,4),XTEMP(2),YTEMP(2)
CCC Normalisation of the field
DO 10 I=1,4
FV(I)=F(I)-V
10 CONTINUE
NA=0
DO 20 I=1,4
I1=MOD(I,4)+1
IF(FV(I)*FV(I1).LT.0) THEN
CCC The curve has a point in the segment I,I+1
NA=NA+1
COEF=FV(I1)/(FV(I1)-FV(I))
A(1,NA)=COEF*X(I)+(1-COEF)*X(I1)
A(2,NA)=COEF*Y(I)+(1-COEF)*Y(I1)
ELSE IF(FV(I).EQ.0.) THEN
CCC The node I is a point of the curve
NA=NA+1
A(1,NA)=X(I)
A(2,NA)=Y(I)
ENDIF
20 CONTINUE
IF(NA.EQ.2) THEN
CCC Two points of the boundary belong to the curve
XTEMP(1)=A(1,1)
XTEMP(2)=A(1,2)
YTEMP(1)=A(2,1)
YTEMP(2)=A(2,2)
CALL TGPL(2,2,XTEMP,YTEMP)
ELSE IF(NA.EQ.4) THEN
CCC Four points of the boundary belong to the curve
CCC => in fact two curves are going through the element
CCC The form depends on the field value at the center of the
element
VAL=(FV(1)+FV(2)+FV(3)+FV(4))*FV(1)
IF(VAL.GT.0.) THEN
CCC Positive value => 2 lines 1-2 and 3-4
XTEMP(1)=A(1,1)
XTEMP(2)=A(1,2)
YTEMP(1)=A(2,1)
YTEMP(2)=A(2,2)
CALL TGPL(2,2,XTEMP,YTEMP)
XTEMP(1)=A(1,3)
XTEMP(2)=A(1,4)
YTEMP(1)=A(2,3)
YTEMP(2)=A(2,4)
CALL TGPL(2,2,XTEMP,YTEMP)
ELSE IF(VAL.LT.0.) THEN
CCC Negative value => 2 lines 2-3 and 4-1
XTEMP(1)=A(1,2)
XTEMP(2)=A(1,3)
YTEMP(1)=A(2,2)
YTEMP(2)=A(2,3)
CALL TGPL(2,2,XTEMP,YTEMP)
XTEMP(1)=A(1,4)
XTEMP(2)=A(1,1)
YTEMP(1)=A(2,4)
YTEMP(2)=A(2,1)
CALL TGPL(2,2,XTEMP,YTEMP)
ELSE
CCC Null value => 2 lines 1-3 and 2-4
XTEMP(1)=A(1,1)
XTEMP(2)=A(1,3)
YTEMP(1)=A(2,1)
YTEMP(2)=A(2,3)
CALL TGPL(2,2,XTEMP,YTEMP)
XTEMP(1)=A(1,2)
XTEMP(2)=A(1,4)
YTEMP(1)=A(2,2)
YTEMP(2)=A(2,4)
CALL TGPL(2,2,XTEMP,YTEMP)
ENDIF
ENDIF
END

This is great! In fact doing it "element-by-element" is perfect for
what I need. I will take some time to see what this routine does. I
am not trying to plot the segments , but rather do some computations
on the zero set (linked list of segments). So, it will help if you
tell me about the variables in this routine -
I know it is an old routine but please see if you can provide some
comments on the following -

1. What XTMP, YTMP hold (of course the ends of the segments) - I am
assuming two points (x1,y1) (x2, y2) form the segment. Is this
correct? And then TGPL just plots the segment. Please confirm. If so,
I will replace call to TGPL by a call to add the segment to the list.
I will manage orientations separately.

2. There are some comments in the routine for the special case (two
segments). What is the numbering scheme referred to (2-3 and 4-1 etc)?

If you find it appropriate I can ask my stupid questions offline.
Understanding this routine will save me a lot of time.

I can't thank you enough :)

There's a happy customer.

I'll start off with the stupid questions. What does your caller look like?
--
Wade Ward- Hide quoted text -

- Show quoted text -

I have dx, dy - increments in x and y directions for the grid
(preferred uniform grid) or alternatively arrays x(:) and y(:) with
coordinate discretizations along x and y axes. Point (i,j) of teh grid
has coordinates (x(i), y(j)). I also have a supposedly preexisting
linked list that defines a "curve" as set of segmens (seg is two
points p1, p2 with direction p1 to p2 - consistent with some criterion
for normal. p2-p1 is my "tangent"). Also, I have a function f(i,j) at
all grid points and f(i,j) changes sign over domain. I want to call
the "isocontour" as follows -

Now, I recently learnt that I d o not need "tailfront" for my linked
list but whatever.....
=======================
prgram main
somewhere in my "main"
............
call isocontour(x, y, f, HeadFront, TailFront (?))

.....
end program
===============
Subroutine isocontour(.....)

first "free" the linked list defined by HeadFront
f should be processed to obtain the curve f=0
While returning HeadFront and TailFront should be changed to point to
a new linked list defining the front (list of segments) characterized
by f=0

end subroutine

Thanks and yes...Glenns reply has made me really happy and any help on
this will make me happier :)

Veda.

.

Loading