V3FIT
sortstuff.f
1  SUBROUTINE sortr(xp,zp,gsq,npc,xmx,zmx,xset,zset,j,jsep)
2  USE precision
3  IMPLICIT NONE
4  INTEGER :: j, jsep
5  REAL(rprec) :: xp(*),zp(*),gsq(*), xmx, zmx, xset, zset
6  REAL(rprec) :: dr1, dr2, dz1, dz2, DOt, cross
7  INTEGER :: i, npc
8 c ---------------------------------------
9 c order points counterclockwise IF needed
10 c ---------------------------------------
11  i=0
12  10 i=i+1
13  dr1=xp(i)-xmx
14  dr2=xp(i+1)-xmx
15  dz1=zp(i)-zmx
16  dz2=zp(i+1)-zmx
17  dot=dr1*dr2+dz1*dz2
18  cross=dr1*dz2-dr2*dz1
19  IF(cross.eq.0.) go to 10
20  IF(dot.eq.0.) go to 10
21  IF(cross/dot.gt.0.) go to 20
22 c ------------------------------------------------
23 c reorder the points to make them counterclockwise
24 c ------------------------------------------------
25  CALL sorter(xp,zp,gsq,npc)
26  20 CONTINUE
27 c ---------------------
28 c delete the last point
29 c ---------------------
30  npc=npc-1
31 c --------------------------------------------------------------------
32 c sort so first poINT is one with smallest positive angle about center
33 c xmx, zmx for j .lt. jsep
34 c pick poINT CLOSEst to xset, zset for j .ge. jsep
35 c --------------------------------------------------------------------
36  CALL sortog(xp,zp,gsq,npc,xmx,zmx,xset,zset,j,jsep)
37 c -----------------------------------------
38 c reset new last poINT to equal first point
39 c -----------------------------------------
40  npc=npc+1
41  xp(npc)=xp(1)
42  zp(npc)=zp(1)
43  gsq(npc)=gsq(1)
44  END SUBROUTINE sortr
45  SUBROUTINE sorter(x,y,d,n)
46  USE precision
47  IMPLICIT NONE
48 cray lcm (x),(y),(d)
49  REAL(rprec) :: x(*),y(*),d(*)
50  REAL(rprec) :: t
51  INTEGER :: n, nsort, nl, i
52  nsort=n/2
53  nl=n+1
54  DO 100 i=1,nsort
55  nl=nl-1
56  t=x(i)
57  x(i)=x(nl)
58  x(nl)=t
59  t=y(i)
60  y(i)=y(nl)
61  y(nl)=t
62  t=d(i)
63  d(i)=d(nl)
64  d(nl)=t
65  100 CONTINUE
66  END SUBROUTINE sorter
67  SUBROUTINE sortog(x,y,d,n,xs,ys,xset,zset,j,jsep)
68  USE precision
69  IMPLICIT NONE
70  INTEGER :: j, jsep,n
71  REAL(rprec) :: xs, ys, xset, zset
72  REAL(rprec) :: pi, ang, angp, tx, ty, td, dist, dt
73  INTEGER :: i, npc, isp, im1, nm1, jj
74  REAL(rprec) :: x(*),y(*),d(*)
75  pi=4*atan(1._dbl)
76  IF(j.ge.jsep) go to 110
77  ang=1.e10
78  DO 50 i=1,n
79  angp=atan2(y(i)-ys,x(i)-xs)
80  angp=abs(angp)
81  IF(angp.gt.ang) go to 50
82  30 isp=i
83  ang=angp
84  50 CONTINUE
85  60 IF(isp.eq.1) RETURN
86  im1=isp-1
87  nm1=n-1
88  DO 100 i=1,im1
89  tx=x(1)
90  ty=y(1)
91  td=d(1)
92  DO 90 jj=1,nm1
93  x(jj)=x(jj+1)
94  y(jj)=y(jj+1)
95  d(jj)=d(jj+1)
96  90 CONTINUE
97  x(n)=tx
98  y(n)=ty
99  d(n)=td
100  100 CONTINUE
101  RETURN
102 c ----------------------------------------------
103 c sorting PROCEDURE for points inside separatrix
104 c find points CLOSEst to xset, zset
105 c ----------------------------------------------
106  110 dist=1.e10
107  DO 130 i=1,n
108  dt=(y(i)-zset)**2+(x(i)-xset)**2
109  IF(dt.gt.dist) go to 130
110  isp=i
111  dist=dt
112  130 CONTINUE
113  go to 60
114  END SUBROUTINE sortog
115