C >>>ANOTB 22:55:23 08/29/73 (73.241) ANOT 1 C ..................................................................ANOT 2 C ANOT 3 C SUBROUTINE ANOTB ANOT 4 C ANOT 5 C PURPOSE ANOT 6 C GIVEN 2 POLYGONS, TO FORM THEIR UNION, INTERSECTION, OR ANY ANOT 7 C ONE OF THE 16 POSSIBLE RESULTANT POLYGONS THAT MAY BE ANOT 8 C FORMED WHEN 2 POLYGONS ARE COMBINED. ANOT 9 C ANOT 10 C USAGE ANOT 11 C CALL ANOTB (AX,AY,NOA,BX,BY,NOB,ABX,ABY,NAB,ITT,DEBUG) ANOT 12 C ANOT 13 C DESCRIPTION OF THE PARAMETERS ANOT 14 C AX - NOA LONG VECTOR WITH X-COORDINATES OF VERTICES OF ANOT 15 C POLYGON #1. ANOT 16 C AY - NOA LONG VECTOR WITH Y-COORDINATES OF VERTICES OF ANOT 17 C POLYGON #1. ANOT 18 C NOA - NUMBER OF VERTICES IN POLYGON #1. ANOT 19 C BX - NOB LONG VECTOR WITH X-COORDINATES OF VERTICES OF ANOT 20 C POLYGON #2. ANOT 21 C BY - NOB LONG VECTOR WITH Y-COORDINATES OF VERTICES OF ANOT 22 C POLYGON #2. ANOT 23 C NOB - NUMBER OF VERTICES IN POLYGON #2. ANOT 24 C ABX - NAB LONG VECTOR WITH X-COORDINATES OF VERTICES OF ANOT 25 C RESULTANT POLYGON. ANOT 26 C ABY - NAB LONG VECTOR WITH Y-COORDINATES OF VERTICES OF ANOT 27 C RESULTANT POLYGON. ANOT 28 C ITT CODE FOR OPERATION WANTED ON POLYGONS. LET 'A' ANOT 29 C BE POLYGON #1, 'B' BE POLYGON #2, '.' STAND FOR ANOT 30 C 'AND', '+' STAND FOR 'OR', AND '¬' STAND FOR 'NOT' ANOT 31 C CODE OPERATION ANOT 32 C 1 A.B ANOT 33 C 2 A.¬B ANOT 34 C 3 A ANOT 35 C 4 B.¬A ANOT 36 C 5 B ANOT 37 C 6 (A+B).¬(A.B) THAT IS EXCLUSIVE OR ANOT 38 C 7 A+B ANOT 39 C 8 ¬(A+B) THAT IS NOR ANOT 40 C 9 (A.B)+¬(A+B) THAT IS EXCLUSIVE NOR ANOT 41 C 10 ¬B ANOT 42 C 11 A+¬B ANOT 43 C 12 ¬A ANOT 44 C 13 B+¬A ANOT 45 C 14 ¬(A.B) THAT IS NAND ANOT 46 C 15 THE WHOLE PLANE; THE UNIVERSAL POLYGON ANOT 47 C 16 NOTHING; THE NULL POLYGON ANOT 48 C ANOT 49 C DEBUG - CONTROLS PRINTING BY ANOTB ANOT 50 C 0 - SEVERE ERRORS (ANOTB FAILURE) PRINTED ANOT 51 C 1 - MINOR ERRORS (BAD DATA) ALSO PRINTED ANOT 52 C 2,3 - FURTHUR DIAGNOSTIC INFO PRINTED ANOT 53 C ANOT 54 C ANOT 55 C DEFINITIONS ANOT 56 C DEFINE A SIMPLE POLYGON TO BE A LIST OF VERTICES IN ANOT 57 C CARTESIAN COORDINATES IN 2 DIMENSIONS WITH THE FIRST VERTEX ANOT 58 C THE SAME AS THE LAST. THE EDGES ARE BETWEEN CONSECUTIVE ANOT 59 C VERTICES. NO 2 EDGES INTERSECT. IF THE VERTICES RUN ANOT 60 C CLOCKWISE, THE REGION COVERED BY THE POLYGON IS CONSIDERED ANOT 61 C TO BE THE PLANE INSIDE THE EDGES AND THE AREA IS POSITIVE. ANOT 62 C IF THE VERTICES RUN COUNTERCLOCKWISE, THE REGION IS THE ANOT 63 C WHOLE PLANE OUTSIDE THE EDGES, AND THE AREA OF THE POLYGON ANOT 64 C IS THE NEGATIVE OF THE AREA INSIDE THE EDGES. ANOT 65 C ANOT 66 C DEFINE A COMPOUND POLYGON TO BE A SET OF SIMPLE POLYGONS, ANOT 67 C CALLED SUBPOLYGONS, WITH NO 2 SUBPOLYGONS INTERSECTING. ANOT 68 C THE SUBPOLYGONS MAY BE NESTED TO ANY DEPTH. HOWEVER A ANOT 69 C SUBPOLYGON IMMEDIATELY INSIDE A CLOCKWISE SUBPOLYGON MUST ANOT 70 C RUN COUNTERCLOCKWISE AND VICE-VERSA. AN ANALOGY IN ANOT 71 C GEOGRAPHY OBTAINS IF YOU CONSIDER LANDMASSES TO BE POSITIVE ANOT 72 C SUBPOLYGONS AND LAKES TO BE NEGATIVE SUBPOLYGONS. YOU CAN ANOT 73 C HAVE ISLANDS IN LAKES IN ISLANDS IN LAKES IN ISLANDS...... ANOT 74 C BUT YOU CANNOT HAVE AN ISLAND IMMEDIATELY IN AN ISLAND. ANOT 75 C HERE THE RESULTING COMPOUND POLYGON IS EXACTLY THE TOTAL ANOT 76 C LANDMASS. THE COMPOUND POLYGON IS DEFINED BY THE VERTICES ANOT 77 C OF THE SUBPOLYGONS LISTED ONE AFTER ANOTHER. THE ANOT 78 C SUBPOLYGONS CAN BE LISTED IN ANY ORDER. SINCE THE FIRST AND ANOT 79 C LAST VERTICES OF EACH SUBPOLYGON ARE EQUAL, THE SUBPOLYGONS ANOT 80 C CAN BE SEPARATED FROM THE COMPOUND POLYGON. THE COMPLEMENT ANOT 81 C OF A POLYGON, SIMPLE OR COMPOUND, THAT IS THE REGION OF THE ANOT 82 C WHOLE PLANE EXCEPT THE POLYGON CAN BE FOUND BY LISTING THE ANOT 83 C VERTICES BACKWARDS. FOR COMPLEAT GENERALITY, 2 DEGENERATE ANOT 84 C POLYGONS ARE ALLOWED. THE NULL POLYGON WHICH OCCUPIES NO ANOT 85 C REGION IS DEFINED BY ONE VERTEX AT (0.,0.). THE UNIVERSAL ANOT 86 C POLYGON WHOSE REGION IS THE ENTIRE PLANE IS DEFINED BY ONE ANOT 87 C VERTEX AT (1.,1.). OTHER THAN THESE 2 POLYGONS, EVERY ANOT 88 C POLYGON HAS AT LEAST 4 VERTICES SINCE IT IS AT LEAST A ANOT 89 C TRIANGLE AND HAS ONE VERTEX REPEATED. ANOT 90 C ANOT 91 C MACHINE DEPENDENCE ANOT 92 C TO SAVE STORAGE, A BIT ARRAY IS STORED IN LCR. IT IS ANOT 93 C MANIPULATED BY ISON WHICH CALLS LOR AND LAND. FOR NON-IBM ANOT 94 C MACHINES, ISON WOULD PROBABLY HAVE TO BE MODIFIED AND LOR ANOT 95 C AND LAND REWRITTEN. ANOT 96 C ANOT 97 C METHOD ANOT 98 C FIRST ASSUME THAT THE INPUT POLYGONS, A AND B, ARE SIMPLE ANOT 99 C POLYGONS. ANOTB IMMEDIATELY CHECKS WHETHER THE OPERATION ANOT 100 C DESIRED IS TRIVIAL (ITT=3, 5, 10, 12, 15 OR 16) OR ONE OF ANOT 101 C THE INPUT POLYGONS ARE TRIVIAL (THE NULL OR UNIVERSAL ANOT 102 C POLYGON). IF SO, THE RESULT IS EASY TO CALCULATE AND ANOT 103 C RETURN. ELSE, THE EQUATIONS OF THE EDGES ARE CALCULATED AND ANOT 104 C THE INTERSECTIONS BETWEEN EDGES OF A AND B ARE STORED IN ANOT 105 C BIT ARRAY LCR. IF THERE ARE NO INTERSECTIONS, THEN EITHER ANOT 106 C ONE POLYGON IS INSIDE THE OTHER, OR THEY ARE DISJOINT. IN ANOT 107 C ANY CASE, THE RESULT IS EITHER A OR B, OR THE COMPLEMENT OF ANOT 108 C ONE OF THEM, OR THE NULL OR UNIVERSAL POLYGON. THE CORRECT ANOT 109 C RESULT IS RETURNED ACCORDING TO A MESSY RULE STORED IN ANOT 110 C DECISION TABLE TABNO. ANOT 111 C ANOT 112 C IN THE DIFFICULT CASE, A AND B INTERSECT. THE VERTICES OF ANOT 113 C THE RESULTING POLYGON ARE ALL THE INTERSECTION POINTS, PLUS ANOT 114 C SOME OF THE ORIGINAL VERTICES OF A AND B. ANOTB STARTS AT ANOT 115 C ANY INTERSECTION POINT AND TRAVELS ALONG ONE OF THE 4 EDGE ANOT 116 C SEGMENTS LEADING FROM THE INTERSECTION. WHICH ONE IS ANOT 117 C DETERMINED BY DECISION TABLE TABYES. TRAVELLING ALONG THE ANOT 118 C EDGE, IF IT COMES TO A VERTEX OF A OR B IT INCLUDES THAT ANOT 119 C VERTEX IN THE RESULTING POLYGON AND CONTINUES ALONG THE ANOT 120 C OTHER EDGE INCIDENT ON THAT VERTEX. IF IT COMES TO AN ANOT 121 C INTERSECTION POINT, IT MUST SWITCH FROM THE POLYGON IT IS ANOT 122 C ON TO THE OTHER POLYGON. HOWEVER THERE ARE 2 WAYS TO DO ANOT 123 C THIS. WHICH WAY IS DECIDED BY TABYES. THE LISTS OF ANOT 124 C INTERSECTION POINTS ALONG EVERY EDGE AND WHAT EDGES THESE ANOT 125 C INTERSECTION POINT ARE WITH HAVE BEEN STORED IN ARRAYS INO ANOT 126 C AND NOR. ANOTB USES THESE TO DETERMINE WHERE THE NEXT ANOT 127 C INTERSECTION POINT IS. NOTICE THAT EVEN IF A AND B ARE ANOT 128 C SIMPLE, THE RESULTING POLYGON MAY BE COMPOUND. IF ANOTB ANOT 129 C TRACES AROUND TO THE INTERSECTION POINT IT STARTED FROM, ANOT 130 C THEN ONE SUBPOLYGON OF THE RESULT HAS BEEN FOUND. NEXT ANOT 131 C ANOTB SEARCHES FOR ANOTHER UNUSED INTERSECTION POINT AND ANOT 132 C STARTS THE PROCESS AGAIN. WHEN NO MORE INTERSECTION POINTS ANOT 133 C REMAIN UNUSED, THE RESULT HAS BEEN COMPLEATLY FOUND. IN THE ANOT 134 C CASES OF ITT=6 OR 9 (EXCLUSIVE OR/NOR), ANOTB REPEATS ITS ANOT 135 C INTERSECTION ALGORITHM AND CATENATES THE RESULTS. IF ITT=6, ANOT 136 C THEN IT DOES ITT=2 AND 4. IF ITT=9, THEN IT DOES ITT= 11 ANOT 137 C AND 13 ANOT 138 C ANOT 139 C NOW ASSUME THAT THE INPUT POLYGONS ARE COMPOUND. HERE SOME ANOT 140 C SUBPOLYGONS OF EACH POLYGON HAVE INTERSECTION POINTS WITH ANOT 141 C THE OTHER POLYGON, AND SOME DO NOT. FIRST THE ONES THAT DO ANOT 142 C NOT ARE CONSIDERED USING THE SAME ALGORITHM AS FOR ANOT 143 C NONINTERSECTING SIMPLE POLYGONS DESCRIBED ABOVE. IF THE ANOT 144 C RESULT FROM THIS STAGE IS THE NULL OR UNIVERSAL POLYGON, ANOT 145 C THEN THIS IS RETURNED. ELSE THE ALGORITHM FOR INTERSECTING ANOT 146 C POLYGONS IS NOW APPLIED AND THE RESULT CATENATED TO THE ANOT 147 C PREVIOUS RESULT. SINCE THIS METHOD STARTS EACH NEW ANOT 148 C SUBPOLYGON AT AN INTERSECTION POINT, IT IGNORES THE ANOT 149 C NONINTERSECTING POLYGONS THAT WERE TREATED EARLIER. ANOT 150 C ANOT 151 C ANOTB PERFORMS SEVERAL NEATENING OPERATIONS ON THE DATA. ANOT 152 C FIRST IF AN INPUT POLYGON HAS CONSECUTIVE VERTICES TOO ANOT 153 C CLOSE TOGETHER, THEN THE SECOND VERTEX IS DELETED. ALSO IF ANOT 154 C THE LAST SUBPOLYGON IS NOT CLOSED, ANOTB CLOSES IT. ANY ANOT 155 C SUBPOLYGON LEFT WITH FEWER THAN 3 DIFFERENT VERTICES ANOT 156 C DELETED. NEXT AS EACH NEW VERTEX OF THE RESULTING POLYGON ANOT 157 C IS FOUND, IT IS NOT USED IF IT IS THE SAME AS THE PREVIOUS ANOT 158 C VERTEX OF THAT SUBPOLYGON OF THE RESULTANT POLYGON. LAST, ANOT 159 C AS EACH COMPLEAT SUBPOLYGON IS FOUND, ALL GROUPS OF THREE ANOT 160 C CONSECUTIVE VERTICES ARE CHECKED. IF THEY ARE COLLINEAR, ANOT 161 C THEN THE MIDDLE VERTEX IS DELETED. THIS IS BECAUSE IF ANOT 162 C THREE VERTICES ARE COLLINEAR, THEN THE MIDDLE ONE IS EITHER ANOT 163 C ON THE STRAIGHT LINE JOINING THE OTHER TWO, OR IT IS ON THE ANOT 164 C TIP OF AN ISOLATED SPIKE. IN EITHER CASE IT IS REDUNDANT. ANOT 165 C AFTER THIS IS DONE, THE FIRST VERTEX OF THE SUBPOLYGON IS ANOT 166 C COPIED AS THE LAST VERTEX TO MARK THE END OF THAT ANOT 167 C SUBPOLYGON. ANOT 168 C ANOT 169 C SIZE LIMITATIONS ANOT 170 C AFTER EACH LIMITATION IS GIVEN THE ARRAYS WHOSE DIMENSIONS ANOT 171 C MUST BE INCREASED AND SIMPLE VARIABLES THAT MUST BE CHANGED ANOT 172 C TO INCREASE THAT LIMITATION. ANOT 173 C ANOT 174 C 99 IS THE MAXIMUM NUMBER OF VERTICES IN EACH INPUT POLYGON. ANOT 175 C GROUP(99,2),A(2,99,2),EQ(3,99,2),SAME(99,2),NOLIN(99,2), ANOT 176 C NOR(10,99,2),INO(10,99,2) ANOT 177 C ANOT 178 C 40 IS THE MAXIMUM NUMBER OF SUBPOLYGONS IN EACH INPUT POLYGON. ANOT 179 C NEWGR(40,2),LMEET(40,2) MMEET=40 ANOT 180 C ANOT 181 C 313 IS (31+THE PRODUCT OF THE MAXIMUM NUMBER OF VERTICES IN ANOT 182 C EACH INPUT POLYGON) / 32 ANOT 183 C LCR(313) ANOT 184 C ANOT 185 C 8=10-2 IS THE MAXIMUM NUMBER OF INTERSECTIONS AN EDGE OF ONE ANOT 186 C POLYGON MAY HAVE WITH EDGES OF THE OTHER POLYGON. ANOT 187 C INO(10,99,2),NOR(10,99,2),PX(10) MCROSS=10 ANOT 188 C ANOT 189 C 300 IS THE MAXIMUM NUMBER OF VERTICES IN THE RESULT PLUS ANOT 190 C THE NUMBER OF VERTICES IN THE INPUT POLYGONS. ANOT 191 C IF(300),CROSS(2,300) LIMIT=300 ANOT 192 C ANOT 193 C REMARKS ANOT 194 C 1) IF THE INPUT POLYGONS ARE VERY LARGE, PASS ANOTB THE ANOT 195 C ARRAY A DIRECTLY INSTEAD OF PASSING AX, AY, BX AND BY SINCE ANOT 196 C ANOTB IMMEDIATELY COPIES THEM INTO A AND THEN DOESN'T USE ANOT 197 C THEM AGAIN. IF NOA=NOB=1000, THIS WILL SAVE 32K. ALSO, THE ANOT 198 C INFORMATION STORED IN INO AND NOR MIGHT BE STORED IN A LIST ANOT 199 C WITH POINTERS TO THE INFORMATION FOR EACH EDGE, INSTEAD OF ANOT 200 C IN A 3-D ARRAY. IF THE INPUT POLYGONS ARE LARGE THIS WILL ANOT 201 C SAVE QUITE A LOT OF STORAGE. ANOT 202 C ANOT 203 C 2) ANOTB CALCULATES MUCH SUPPLEMENTARY INFORMATION ANOT 204 C (STORED IN ARRAYS GROUP, EQ, SAME AND NEWGR) ABOUT EACH ANOT 205 C INPUT POLYGON BEFORE IT STARTS TO CALCULATE THE RESULT. IF ANOT 206 C THERE ARE MANY POLYGONS THAT ARE BEING INTERSECTED OR ANOT 207 C WHATEVER MANY TIMES, IT MIGHT BE WORTHWHILE TO STORE THIS ANOT 208 C INFORMATION INSTEAD OF CALCULATING IT AGAIN EACH TIME. ANOT 209 C ANOT 210 C AUTHOR ANOT 211 C WM. RANDOLPH FRANKLIN, ANOT 212 C DEPT. OF GEOGRAPHY, OR 458 THESSALY CIRCLE, ANOT 213 C SIMON FRASER UNIVERSITY, OTTAWA, CANADA, ANOT 214 C BURNABY 2, B. C. K1H 5W5 ANOT 215 C CANADA. ANOT 216 C ANOT 217 C VERSION ANOT 218 C THIS VERSION JUNE, 1973. THE DIFFERENCES BETWEEN THIS AND ANOT 219 C THE PREVIOUS VERSION ARE BETTER COMMENTS, THE CORRECTION OF ANOT 220 C AN ERROR WHICH OCCURED WHEN A VERTEX OF ONE POLYGON WAS ON ANOT 221 C AN EDGE OF THE OTHER POLYGON, THE STORAGE OF LCR AS A BIT ANOT 222 C ARRAY AND NEATENING THE RESULTANT POLYGON. ANOT 223 C ANOT 224 C KNOWN BUGS ANOT 225 C THE RESULTS MAY BE WRONG IF: ANOT 226 C 1) AN INPUT POLYGON INTERSECTS ITSELF OR HAS AN ILLEGAL ANOT 227 C COMBINATION OF SUBPOLYGONS (E.Q. A CLOCKWISE SUBPOLYGON ANOT 228 C IMMEDIATELY INSIDE A CLOCKWISE SUBPOLYGON). ANOT 229 C 2) AN EDGE OF ONE POLYGON COINCIDES WITH AN EDGE OF THE ANOT 230 C OTHER POLYGON, BUT THE 2 POLYGONS DON'T OTHERWISE INTERSECT. ANOT 231 C 3) A VERTEX OF ONE POLYGON IS ON AN EDGE OF THE OTHER POLYGON. ANOT 232 C ANOT 233 C ANOT 234 C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED ANOT 235 C AREA - FUNCTION TO CALCULATE SIGNED AREA OF A POLYGON. ANOT 236 C PNPOLY - SUBROUTINE TO TELL IF A POINT IS IN A POLYGON. ANOT 237 C ISON - SUBROUTINE TO MANIPULATE BIT ARRAYS ANOT 238 C LOR - ASSEMBLER SUBROUTINE FOR LOGICAL OR ANOT 239 C LAND - ASSEMBLER SUBROUTINE FOR LOGICAL AND ANOT 240 C ANOT 241 C ..................................................................ANOT 242 C ANOT 243 SUBROUTINE ANOTB (AX,AY,NOA,BX,BY,NOB,ABX,ABY,NAB,ITT,DEBUG) ANOT 244 INTEGER*2 GROUP(99,2),NEWGR(40,2),IF(300),NN(2),NGR(2), ANOT 245 1 SAME(99,2),NOLIN(99,2),NOR(10,99,2),INO(10,99,2) ANOT 246 LOGICAL*1 LMEET(40,2) ANOT 247 INTEGER LCR( 313),GO,GOB,O ANOT 248 LOGICAL LBUG,LBUG2,LSTART,LERR ANOT 249 REAL A(2,99,2),EQ(3,99,2),CROSS(2,300),PX(10) ANOT 250 1 , AX(NOA),AY(NOA),BX(NOB),BY(NOB),ABX(300),ABY(300),AREAS(2) ANOT 251 C NGR(I) IS THE NUMBER OF SUBPOLYGONS IN POLYGON #I. ANOT 252 C LMEET(J,I)=.TRUE. IFF SUBPOLYGON #J OF POLYGON #I INTERSECTS ANOT 253 C POLYGON #3-I. ANOT 254 C GROUP(J,I) TELLS WHICH SUBPOLYGON OF POLYGON #I, VERTEX #J OF ANOT 255 C POLYGON #I IS IN. ANOT 256 C NEWGR(J,I) GIVES THE LAST VERTEX IN POLYGON #I WHICH IS IN THAT ANOT 257 C POLYGON'S J-TH SUBPOLYGON. ANOT 258 C ANOT 259 INTEGER*2 ISTART(16),TABYES(2,16),TABNO(2,2,16) ANOT 260 DATA ISTART /1,-1,0,-1,0,0,-1,1,0,0,1,0,1,-1,0,0/ ANOT 261 #, TABYES /-1,-1,1,-1,2*30000,-1,1,4*30000,1,1,1,1, ANOT 262 14*30000,-1,1,2*30000,1,-1,-1,-1,4*30000/ ANOT 263 #, TABNO /1,3,1,3,3,1,2,3,4*5,2,3,3,1,4*5,2,1,2,1,4ANOT 264 1,1,4,1,3,2,3,2,1,2,1,2,4*5,1,4,4,2,4*5,4,2,1,4,2,4,2,4,8*5/ ANOT 265 #,O,KTEMP/6,4/ ANOT 266 C ANOT 267 C TABYES IS A DECISION TABLE USED WHEN 2 SUBPOLYGONS INTERSECT. ANOT 268 C THE ENTRIES WITH 30000 SHOULDN'T BE USED. ANOT 269 C TABNO IS THE DECISION TABLE USED FOR A DISJOINT SUBPOLYGON. ANOT 270 C LET I= NUMBER OF POLYGON TO WHICH SUBPOLYGON BELONGS. ANOT 271 C K= 1 IF THIS SUBPOLYGON IS INSIDE POLYGON#(3-I), ELSE K=2. ANOT 272 C THEN TABNO(K,I,IT)=1 IF THIS SUBPOLYGON IS TO BE USED, ANOT 273 C 2 IF ITS REVERSE IS TO BE USED, ANOT 274 C 3 IF IT IS NOT TO BE USED AT ALL, ANOT 275 C 4 IF THE RESULT IS THE WHOLE PLANE AS FAR ANOT 276 C AS THIS SUBPOLYGON IS CONCERNED, ANOT 277 C 5 IF THIS ENTRY IS NOT TO BE USED. ANOT 278 C ANOT 279 C ANOT 280 C IF THE RESULTANT POLYGON HAS NO VERTICES, NAB=1 IS RETURNED. ANOT 281 C IN THIS CASE, IF THE POLYGON IS NOTHING, THEN ABX(1)=ABY(1)=0.0 ANOT 282 C OR IF THE POLYGON IS THE WHOLE PLANE, ABX(1)=ABY(1)=1.0 ANOT 283 C NOLIN(I,J) : NO OF LINES THAT CROSS THE I-TH LINE OF THE J-TH ANOT 284 C POLYGON. ANOT 285 C INO(K, I, J) : K-TH LINE TO CROSS THE I-TH LINE OF THE J-TH ANOT 286 C POLYGON. ANOT 287 C CROSS(*, NOR (K,I,J) ): COORDINATES OF THE ABOVE INTERSECTION. ANOT 288 C LAST SUBSCRIPT OF AB IS NAB. ANOT 289 C O IS OUTPUT UNIT FOR PRINTED MESSAGES ANOT 290 C KTEMP IS A TEMPORARY ARRAY USED TO STORE LCR WHEN ITT=6 OR 9 ANOT 291 C FOR KTEMP, DCB=(LRECL=1254,BLKSIZE=1258,RECFM=VBS) ANOT 292 KSA=1 ANOT 293 TOLER=.00001 ANOT 294 MMEET=40 ANOT 295 LERR=DEBUG.GE.1. ANOT 296 LBUG2=DEBUG.GE.2. ANOT 297 LBUG=DEBUG.GE.3. ANOT 298 MCROSS=10 ANOT 299 IT=ITT ANOT 300 LIMIT=300 ANOT 301 C ANOT 302 C ON ERROR GO TO 990 ANOT 303 C ANOT 304 IF(.NOT.LBUG) GO TO 304 ANOT 305 WRITE(O,305) NOA,NOB,ITT ANOT 306 305 FORMAT('-SUBROUTINE ANOTB CALL: NOA=',I5,' NOB=',I5,' ITT=',I2) ANOT 307 K=1 ANOT 308 IF(NOA.GE.1) WRITE(O,306) (K,I,(AX(I),AY(I),J=1,2),I=1,NOA) ANOT 309 306 FORMAT(//(' POLYGON #',I1,' VERTEX #',I5,' : ',2G15.7,5X2Z10)) ANOT 310 K=2 ANOT 311 IF(NOB.GE.1) WRITE(O,306) (K,I,(BX(I),BY(I),J=1,2),I=1,NOB) ANOT 312 304 NANB=(NOA*NOB+31)/32 ANOT 313 NAB=0 ANOT 314 DO 1 I=1,2 ANOT 315 DO 1 J=1,MMEET ANOT 316 1 LMEET(J,I)=.FALSE. ANOT 317 AREAS(1)=AREA(AX,AY,NOA) ANOT 318 AREAS(2)=AREA(BX,BY,NOB) ANOT 319 C TOLER IS AN ABSOLUTE FUDGE FACTOR. ANOT 320 C FUDGE IS A RELATIVE FUDGE FACTOR. ANOT 321 FUDGE=TOLER*SQRT(AMAX1(ABS(AREAS(1)),ABS(AREAS(2)))) ANOT 322 NCR=0 ANOT 323 C <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<0 IF THE ANOT 612 C POINT IS IN THE POLYGON, AND <0 IF IT IS OUTSIDE. ANOT 613 C ANOT 614 DO 33 I=1,N ANOT 615 J=1+I ANOT 616 IF ((A(1,I,L)-A(2,I,L)+A(2,J,L))*EQ(1,I,L)+(A(2,I,L)+A(1,I,L)-A(1,ANOT 617 1J,L))*EQ(2,I,L).GT.EQ(3,I,L)) GO TO 33 ANOT 618 DO 32 K=1,3 ANOT 619 EQ(K,I,L)=-EQ(K,I,L) ANOT 620 32 CONTINUE ANOT 621 33 CONTINUE ANOT 622 34 CONTINUE ANOT 623 C ANOT 624 C THE FOLLOWING DETERMINES WHICH SIDES OF THE 2 POLYGONS CROSS. ANOT 625 C LCR(I,J)=.TRUE. IF THE I-TH SIDE OF THE 1-ST POLYGON CROSSES THE ANOT 626 C J-TH SIDE OF THE 2-ND POLYGON. ANOT 627 C ANOT 628 DO 43 I=1,NA ANOT 629 DO 43 J=1,NB ANOT 630 C THE YK CALCULATED FOR SOME I MUST BE THE SAME TO THE LAST ANOT 631 C BIT AS THE XK CALCULATED FOR THE NEXT HIGHER I. HENCE THE ANOT 632 C EXTRA PARENTHESES. ANOT 633 XK=(A(1,I,1)*EQ(1,J,2)+A(2,I,1)*EQ(2,J,2))-EQ(3,J,2) ANOT 634 YK=(A(1,I+1,1)*EQ(1,J,2)+A(2,I+1,1)*EQ(2,J,2))-EQ(3,J,2) ANOT 635 C N.B. SIGN(1.,0.) IS UNPREDICTABLE. ANOT 636 IK=SIGN(1.,XK) ANOT 637 IF (XK.EQ.0.) IK=1 ANOT 638 JK=SIGN(1.,YK) ANOT 639 IF (YK.EQ.0.) JK=1 ANOT 640 K=IK+JK ANOT 641 IF (K.NE.0) GO TO 42 ANOT 642 XL=(A(1,J,2)*EQ(1,I,1)+A(2,J,2)*EQ(2,I,1))-EQ(3,I,1) ANOT 643 YL=(A(1,J+1,2)*EQ(1,I,1)+A(2,J+1,2)*EQ(2,I,1))-EQ(3,I,1) ANOT 644 IL=SIGN(1.,XL) ANOT 645 IF (XL.EQ.0.) IL=-1 ANOT 646 JL=SIGN(1.,YL) ANOT 647 IF (YL.EQ.0.) JL=-1 ANOT 648 L=IL+JL ANOT 649 IF (L.NE.0) GO TO 42 ANOT 650 C ANOT 651 C EDGE(I,1) MEANS THE I-TH EDGE OF THE FIRST POLYGON, ETC. ANOT 652 C IF A POINT IS ON AN EDGE OF POLYGON #1,THEN CONSIDER THAT POINT ANOT 653 C TO BE OUTSIDE POLYGON #1. HOWEVER, IF A POINT IS ON AN EDGE OF ANOT 654 C POLYGON #2, CONSIDER THAT POINT TO BE INSIDE POLYGON #2. ANOT 655 C K=-2 MEANS THAT BOTH ENDS OF EDGE(I,1) ARE OUTSIDE POLYGON #2. ANOT 656 C 0 ONE END IS INSIDE AND ONE OUTSIDE. ANOT 657 C 2 BOTH ENDS ARE INSIDE. ANOT 658 C L=-2 MEANS THAT BOTH ENDS OF EDGE(J,2) ARE OUTSIDE POLYGON #1. ANOT 659 C 0 ONE END IS INSIDE AND ONE OUTSIDE. ANOT 660 C 2 BOTH ENDS ARE INSIDE. ANOT 661 C I AM USING 'INSIDE' AND 'OUTSIDE' TO MEAN LOCALLY WITH RESPECT TO ANOT 662 C THE EDGES THAT I AM CURRENTLY CONSIDERING. IF I SAY THAT POINT 'P'ANOT 663 C IS INSIDE A POLYGON WITH RESPECT TO EDGE 'E', I MEAN THAT IF A ANOT 664 C LINE WERE DRAWN FROM P TO ANY POINT 'Q' ON E, THEN POINTS ON PQ ANOT 665 C ARBITRARILY CLOSE TO E ARE INSIDE THE POLYGON. IF PQ IS EXTENDED ANOT 666 C ACROSS E, THEN POINTS ON PQ EXTENDED ARBITRARILY CLOSE TO E ARE ANOT 667 C OUTSIDE THE POLYGON. THAT INTERSECTIONS OF PQ WITH OTHER EDGES OF ANOT 668 C THE POLYGON WOULD PUT P REALLY OUTSIDE THE POLYGON DOES NOT ANOT 669 C MATTER HERE. ANOT 670 C ANOT 671 C EDGE(I,1) CUTS EDGE(J,2) IFF BOTH ENDS OF EDGE(I,1) ARE ON ANOT 672 C OPPOSITE SIDES OF EDGE(J,2) (THAT IS K=0) AND VICE VERSA, (THAT ISANOT 673 C L=0). THUS FOR AN INTERSECTION K=0 AND L=0. ANOT 674 C ANOT 675 JJ=ISON(LCR,I,J,0) ANOT 676 IF (LBUG) WRITE (O,35) I,J,XK,YK,IK,JK,K,XL,YL,IL,JL,L,JJ ANOT 677 35 FORMAT (' I,J,XK,YK,IK,JK,K,XL,YL,IL,JL,L,JJ=',2I4,2G13.5,3I3,1X2GANOT 678 113.5,4I5) ANOT 679 IF (JJ.EQ.0) GO TO 43 ANOT 680 NI=I ANOT 681 NJ=J ANOT 682 Y=EQ(1,I,1)*EQ(2,J,2)-EQ(2,I,1)*EQ(1,J,2) ANOT 683 C Y IS THE SINE OF THE ANGLE BETWEEN THE EDGES. ANOT 684 IF (ABS(Y).GT.TOLER) GO TO 38 ANOT 685 C EDGES ARE NEARLY PARALLEL. ANOT 686 C ARBITRARILY PUT THE INTERSECTION POINT HALFWAY BETWEEN THE MIDDLE ANOT 687 C 2 ENDPOINTS. ANOT 688 DO 36 K=1,2 ANOT 689 36 PX(K)=A(K,I,1) + A(K,I+1,1) + A(K,J,2) + A(K,J+1,2) - ANOT 690 #AMAX1 (A(K,I,1), A(K,I+1,1), A(K,J,2), A(K,J+1,2) )- ANOT 691 #AMIN1 (A(K,I,1), A(K,I+1,1), A(K,J,2), A(K,J+1,2) ) ANOT 692 X=PX(1) ANOT 693 Y=PX(2) ANOT 694 GO TO 39 ANOT 695 38 X=(EQ(3,I,1)*EQ(2,J,2)-EQ(2,I,1)*EQ(3,J,2))/Y ANOT 696 Y=(EQ(1,I,1)*EQ(3,J,2)-EQ(3,I,1)*EQ(1,J,2))/Y ANOT 697 39 K=NOLIN(I,1)+1 ANOT 698 IF (K.GT.MCROSS) GO TO 40 ANOT 699 NOLIN(I,1)=K ANOT 700 INO(K,I,1)=J ANOT 701 NCR=NCR+1 ANOT 702 NOR(K,I,1)=NCR ANOT 703 CROSS(1,NCR)=X ANOT 704 CROSS(2,NCR)=Y ANOT 705 K=NOLIN(J,2)+1 ANOT 706 IF (K.GT.MCROSS) GO TO 40 ANOT 707 NOLIN(J,2)=K ANOT 708 INO(K,J,2)=I ANOT 709 NOR(K,J,2)=NCR ANOT 710 LMEET(GROUP(I,1),1)=.TRUE. ANOT 711 LMEET(GROUP(J,2),2)=.TRUE. ANOT 712 GO TO 43 ANOT 713 40 WRITE (O,41) I,J ANOT 714 41 FORMAT ('-ERROR: IN SUBROUTINE ANOTB, MCROSS EXCEEDED. I=',I5,' J=ANOT 715 1',I5/'-SEE COMMENTS ON SIZE LIMITATIONS AT START OF ANOTB.') ANOT 716 RETURN ANOT 717 42 JJ=ISON(LCR,I,J,2) ANOT 718 43 CONTINUE ANOT 719 C ANOT 720 C NI=-1 IF THE 2 POLYGONS DON'T INTERSECT AT ALL. ANOT 721 C FILL ARRAY NEWGR. ANOT 722 C NEWGR(J,I)=LAST VERTEX IN SUBPOLYGON #J OF POLYGON #I. ANOT 723 C ANOT 724 DO 45 I=1,2 ANOT 725 N=NN(I)+1 ANOT 726 K=1 ANOT 727 DO 44 J=1,N ANOT 728 IF (GROUP(J,I).EQ.K) GO TO 44 ANOT 729 NEWGR(K,I)=J-1 ANOT 730 K=K+1 ANOT 731 44 CONTINUE ANOT 732 NEWGR(K,I)=N ANOT 733 45 CONTINUE ANOT 734 IF (.NOT.LBUG) GO TO 57 ANOT 735 WRITE(O,450) ISTART,TABYES,TABNO ANOT 736 450 FORMAT(' ISTART:',16I6,2(/' TABYES:',16I6),4(/' TABNO:',16I6)) ANOT 737 ASSIGN 57 TO KANT ANOT 738 46 WRITE (O,47) (LCR(I),I=1,NANB) ANOT 739 47 FORMAT (' LCR:',10Z10) ANOT 740 DO 51 L=1,2 ANOT 741 M=NN(L) ANOT 742 WRITE (O,48) L,(SAME(J,L),J=1,M) ANOT 743 48 FORMAT (' POLYGON #',I3/' SAME= ',(20I5)) ANOT 744 WRITE (O,49) (NOLIN(J,L),J=1,M) ANOT 745 49 FORMAT (' NOLIN=',(20I5)) ANOT 746 WRITE (O,50) ((INO(J,K,L),J=1,10),(NOR(J,K,L),J=1,10),K=1,M) ANOT 747 50 FORMAT (' INO:',60X,'NOR:'/(10I5,15X10I5)) ANOT 748 51 CONTINUE ANOT 749 WRITE (O,52) NCR,(J,(CROSS(I,J),I=1,2),J=1,NCR) ANOT 750 52 FORMAT (' NCR,CROSS:',I5/(3(I10,2G15.7))) ANOT 751 DO 53 L=1,2 ANOT 752 N=NN(L) ANOT 753 WRITE (O,54) L,(J,(EQ(I,J,L),I=1,3),((A(I,J,L),I=1,2),K=1,2), ANOT 754 1GROUP(J,L),J=1,N) ANOT 755 53 CONTINUE ANOT 756 54 FORMAT ('0EDGE EQUATIONS FOR POLYGON',I3,T70,'VERTEX',T125,'GROUP'ANOT 757 1/(' SIDE',I3,':',G15.5,'*X+',G15.5,'*Y=',G15.5,T65,2G15.7,2Z10, ANOT 758 2T125,I3)) ANOT 759 DO 55 I=1,2 ANOT 760 L=NGR(I) ANOT 761 WRITE (O,56) I,L,(J,LMEET(J,I),NEWGR(J,I),J=1,L) ANOT 762 55 CONTINUE ANOT 763 56 FORMAT (' POLYGON #',I1,' HAS',I3,' GROUPS.',8(T30,I2,':',L5,I5/))ANOT 764 GO TO KANT, (57,97) ANOT 765 C ANOT 766 C INO(*,I,J) CONTAINS THE LINE NUMBERS OF THE (3-J)-TH POLYGON ANOT 767 C THAT INTERSECT THE I-TH LINE OF THE J-TH POLYGON. ANOT 768 C ANOT 769 57 CONTINUE ANOT 770 C ANOT 771 C SOME SUBPOLYGONS MAY BE DISJOINT, THEY DO NOT INTERSECT WITH THE ANOT 772 C OTHER POLYGON. TREAT THEM FIRST. ANOT 773 C ANOT 774 ASSIGN 64 TO GO ANOT 775 ASSIGN 114 TO GOB ANOT 776 C *****THE FOLLOWING STATEMENTS SIMULATE 2 DO LOOPS SINCE WATFIV ANOT 777 C DOES NOT RECKOGNIZE THE EXTENDED RANGE OF A DO LOOP. ANOT 778 I=0 ANOT 779 58 I=I+1 ANOT 780 IF (I.GT.2) GO TO 65 ANOT 781 N=NGR(I) ANOT 782 J=0 ANOT 783 59 J=J+1 ANOT 784 IF (J.GT.N) GO TO 58 ANOT 785 IF (LMEET(J,I)) GO TO 64 ANOT 786 KSA=1 ANOT 787 IF (J.GT.1) KSA=NEWGR(J-1,I)+1 ANOT 788 M=NEWGR(J,I)+1-KSA ANOT 789 C KSA IS THE FIRST VERTEX OF SUBPOLYGON #J OF POLYGON #I. ANOT 790 IF (I.EQ.2) GO TO 60 ANOT 791 C ANOT 792 C I AM TRYING TO SELECT A POINT ON POLYGON A WHICH IS FAR ENOUGH ANOT 793 C FROM AN EDGE OF POLYGON B TO AVOID ROUND OFF ERROR. THIS RANDOM ANOT 794 C POINT SHOULD BE OK. ANOT 795 C ANOT 796 TEMP1=A(1,KSA,1)*0.1554304+A(1,KSA+1,1)*0.8445696 ANOT 797 TEMP2=A(2,KSA,1)*0.1554304+A(2,KSA+1,1)*0.8445696 ANOT 798 CALL PNPOLY (TEMP1,TEMP2,BX,BY,NOB,K) ANOT 799 GO TO 61 ANOT 800 60 TEMP1=A(1,KSA,2)*0.1554304+A(1,KSA+1,2)*0.8445696 ANOT 801 TEMP2=A(2,KSA,2)*0.1554304+A(2,KSA+1,2)*0.8445696 ANOT 802 CALL PNPOLY (TEMP1,TEMP2,AX,AY,NOA,K) ANOT 803 61 IF(LBUG2) WRITE(O,600)I,J,K,L,M,N,KSA,TEMP1,TEMP2 ANOT 804 600 FORMAT(' I=',I2,' J=',I3,' K=',I2,' L=',I6,' M=',I4,' N=', ANOT 805 * I3,' KSA=',I5,' TEMP1,TEMP2=',2G15.7) ANOT 806 IF(K.EQ.0) GO TO 994 ANOT 807 K=TABNO((3-SIGN(1.0,K*AREAS(3-I)))/2,I,IT) ANOT 808 IF(LBUG2)WRITE(O,601)K,IT ANOT 809 601 FORMAT(' K=',I5,' IT=',I5) ANOT 810 GO TO (66,68,63,62,993), K ANOT 811 62 ASSIGN 70 TO GOB ANOT 812 GO TO 64 ANOT 813 63 ASSIGN 72 TO GOB ANOT 814 64 CONTINUE ANOT 815 GO TO 59 ANOT 816 65 CONTINUE ANOT 817 IF (NI.GT.0) GO TO 74 ANOT 818 GO TO GOB, (114,70,72) ANOT 819 C INCLUDE THIS (SUB)POLYGON. ANOT 820 66 DO 67 L=1,M ANOT 821 NAB=NAB+1 ANOT 822 ABX(NAB)=A(1,KSA-1+L,I) ANOT 823 ABY(NAB)=A(2,KSA-1+L,I) ANOT 824 67 CONTINUE ANOT 825 GO TO GO, (64,114) ANOT 826 C INCLUDE THE COMPLEMENT OF THIS (SUB)POLYGON. ANOT 827 68 DO 69 L=1,M ANOT 828 NAB=NAB+1 ANOT 829 ABX(NAB)=A(1,M+KSA-L,I) ANOT 830 ABY(NAB)=A(2,M+KSA-L,I) ANOT 831 69 CONTINUE ANOT 832 GO TO GO, (64,114) ANOT 833 C THE RESULT IS THE WHOLE PLANE. ANOT 834 70 IF (NAB.GT.1) GO TO 114 ANOT 835 71 ABX(1)=1.0 ANOT 836 ABY(1)=1.0 ANOT 837 NAB=1 ANOT 838 GO TO 114 ANOT 839 C THE RESULT IS NOTHING. ANOT 840 72 IF (NAB.GT.1) GO TO 114 ANOT 841 73 ABX(1)=0.0 ANOT 842 ABY(1)=0.0 ANOT 843 NAB=1 ANOT 844 GO TO 114 ANOT 845 C ANOT 846 C THE 2 POLYGONS INTERSECT. $$$$$ ANOT 847 C ANOT 848 74 CONTINUE ANOT 849 ASSIGN 114 TO GO ANOT 850 DO 78 I=1,2 ANOT 851 K=NN(I) ANOT 852 DO 77 J=1,K ANOT 853 L=NOLIN(J,I) ANOT 854 IF (L.EQ.2) GO TO 77 ANOT 855 INO(2,J,I)=INO(L,J,I) ANOT 856 INO(L,J,I)=-2 ANOT 857 M=NOR(2,J,I) ANOT 858 NOR(2,J,I)=NOR(L,J,I) ANOT 859 NOR(L,J,I)=M ANOT 860 IF (L.EQ.3) GO TO 77 ANOT 861 C ANOT 862 C SORT INO, NOR ANOT 863 C ANOT 864 L=L-1 ANOT 865 N=NOR(1,J,I) ANOT 866 DO 75 M=2,L ANOT 867 II=NOR(M,J,I) ANOT 868 PX(M)=ABS(CROSS(1,II)-CROSS(1,N))+ABS(CROSS(2,II)-CROSS(2,N)) ANOT 869 75 CONTINUE ANOT 870 N=L-1 ANOT 871 DO 76 II=2,N ANOT 872 JJ=II+1 ANOT 873 DO 76 IJ=JJ,L ANOT 874 IF (PX(IJ).GT.PX(II)) GO TO 76 ANOT 875 DEN=PX(II) ANOT 876 PX(II)=PX(IJ) ANOT 877 PX(IJ)=DEN ANOT 878 M=INO(II,J,I) ANOT 879 INO(II,J,I)=INO(IJ,J,I) ANOT 880 INO(IJ,J,I)=M ANOT 881 M=NOR(II,J,I) ANOT 882 NOR(II,J,I)=NOR(IJ,J,I) ANOT 883 NOR(IJ,J,I)=M ANOT 884 76 CONTINUE ANOT 885 77 CONTINUE ANOT 886 78 CONTINUE ANOT 887 IF (.NOT.LBUG) GO TO 80 ANOT 888 DO 79 L=1,2 ANOT 889 M=NN(L) ANOT 890 WRITE (O,48) L,(SAME(J,L),J=1,M) ANOT 891 WRITE (O,49) (NOLIN(J,L),J=1,M) ANOT 892 WRITE (O,50) ((INO(J,K,L),J=1,10),(NOR(J,K,L),J=1,10),K=1,M) ANOT 893 79 CONTINUE ANOT 894 80 CONTINUE ANOT 895 LSTART=.TRUE. ANOT 896 KF=0 ANOT 897 IF (IT.EQ.6) GO TO 81 ANOT 898 IF (IT.NE.9) GO TO 83 ANOT 899 IT=1 ANOT 900 ASSIGN 118 TO GO ANOT 901 GO TO 82 ANOT 902 81 IT=2 ANOT 903 ASSIGN 119 TO GO ANOT 904 82 CONTINUE ANOT 905 REWIND KTEMP ANOT 906 WRITE (KTEMP) LCR ANOT 907 REWIND KTEMP ANOT 908 NIB=NI ANOT 909 NJB=NJ ANOT 910 IF (.NOT.LBUG2) GO TO 83 ANOT 911 WRITE(O,820) (LCR(I),I=1,NANB) ANOT 912 820 FORMAT(' LCR AFTER WRITE: ',10Z10) ANOT 913 83 L=1 ANOT 914 C ANOT 915 C L IS THE POLYGON WHOSE EDGE I AM ON. ANOT 916 C NI,NJ ARE THE 2 INTERSECTING SIDES. ANOT 917 C SINCE THE DIRECTION OF THE CYCLE IS IMPORTANT, THE FOLLOWING GETS ANOT 918 C ME STARTED IN THE RIGHT DIRECTION. ANOT 919 C ANOT 920 IF (INT(SIGN(1.0,(2.0*A(1,NI,1)-A(1,NI+1,1))*EQ(1,NJ,2)+(2.0*A(2,NANOT 921 1I,1)-A(2,NI+1,1))*EQ(2,NJ,2)-EQ(3,NJ,2)))*ISTART(IT)) 84,98,90 ANOT 922 84 L=2 ANOT 923 K=NJ ANOT 924 NJ=NI ANOT 925 NI=K ANOT 926 C ANOT 927 C THE FOLLOWING FINDS 1 CLOSED CYCLE OF VERTICES FOR THE NEW ANOT 928 C POLYGON. ANOT 929 C ANOT 930 85 IF (NJ.GT.0) GO TO 90 ANOT 931 C ANOT 932 C THIS NEW POINT IS THE CORNER OF AN OLD POLYGON. ANOT 933 C ANOT 934 IF (NJ.NE.-2) GO TO 86 ANOT 935 K=SAME(NI,L) ANOT 936 IF (K.GE.0) K=1 ANOT 937 NI=NI+K ANOT 938 86 CONTINUE ANOT 939 NAB=NAB+1 ANOT 940 IF (NAB.EQ.LIMIT) GO TO 112 ANOT 941 IF (NI.LE.0) GO TO 992 ANOT 942 ABX(NAB)=A(1,NI,L) ANOT 943 ABY(NAB)=A(2,NI,L) ANOT 944 C IF THIS VERTEX IS THE SAME AS THE LAST VERTEX, AND IT IS ANOT 945 C NOT THE FIRST VERTEX OF A SUBPOLYGON, DELETE IT. ANOT 946 IF (LSTART) GO TO 87 ANOT 947 IF (ABS(ABX(NAB)-ABX(NAB-1)).GT.FUDGE.OR. ABS(ABY(NAB)-ABY(NAB-1))ANOT 948 # .GT.FUDGE)GO TO 87 ANOT 949 IF(LBUG2) WRITE(O,870) K,L,NI,NJ,NAB ANOT 950 870 FORMAT(' NEW VERTEX TOO CLOSE TO PREV. ONE: K,L,NI,NJ,NAB=', ANOT 951 1 5I10) ANOT 952 NAB=NAB-1 ANOT 953 87 LSTART=.FALSE. ANOT 954 IF (LBUG2) WRITE (O,88) K,L,NI,NJ,I,NAB,ABX(NAB),ABX(NAB) ANOT 955 1 ,ABY(NAB),ABY(NAB) ANOT 956 88 FORMAT (' NEW VERTEX: K=',I3,' L=',I3,' NI=',I3,' NJ=',I3,' ANOT 957 * I=',I3,' NAB=',I3,' ABX=',G15.7,Z10,' ABY=',G15.7,Z10) ANOT 958 IF (NJ.EQ.-1) GO TO 89 ANOT 959 NJ=INO(2,NI,L) ANOT 960 GO TO 85 ANOT 961 89 CONTINUE ANOT 962 K=SAME(NI,L) ANOT 963 IF (K.LE.0) K=-1 ANOT 964 NI=NI+K ANOT 965 NJ=INO(NOLIN(NI,L)-1,NI,L) ANOT 966 GO TO 85 ANOT 967 C ANOT 968 C THIS NEW POINT IS THE INTERSECTION OF 2 OLD EDGES. ANOT 969 C ANOT 970 90 L=3-L ANOT 971 N=NOLIN(NJ,L) ANOT 972 DO 91 I=1,N ANOT 973 IF (INO(I,NJ,L).EQ.NI) GO TO 100 ANOT 974 91 CONTINUE ANOT 975 C ANOT 976 C SOMETHING IS WRONG. ANOT 977 C ANOT 978 II=7772 ANOT 979 92 WRITE (O,93) ANOT 980 93 FORMAT ('-*******************HELP? SUBROUTINE ANOTB IN TROUBLE. ANOT 981 1 ',/,20X,'PLEASE SEND YOUR LISTING TO R. FRANKLIN, ',/,20X,'458 THANOT 982 2ESSALY CIRCLE, OTTAWA, CANADA, K1H 5W5') ANOT 983 WRITE (O,94) I,J,K,L,M,N,II,IJ,JJ,NA,NB,NI,NJ,IT,ITT,KSA,NN,NGR,NOANOT 984 1A,NOB,NAB,NIB,NJB,LIMIT,NCR,MMEET,O,IG,LBUG,LBUG2,LSTART,KF ANOT 985 #,TEMP1,TEMP2 ANOT 986 94 FORMAT (' I=',I5,' J=',I5,' K=',I5,' L=',I5,' M=',I5,' N=',I5,' IIANOT 987 1=',I5,' IJ=',I5,' JJ=',I5,' NA=',I5,' NB=',I5,' NI=',I5/' NJ=',I5,ANOT 988 2' IT=',I5,' ITT=',I5,' KSA=',I5,' NN=',2I5,' NGR=',2I5,/' NOA=',I5ANOT 989 3,' NOB=',I5,' NAB=',I5/' NIB=',I5,' NJB=',I5,' LIMIT=',I5,' ANOT 990 4NCR=',I5,' MMEET=',I5,' O=',I5,' IG=',I5,' LBUG=',L1,' LBUG2=',LANOT 991 51,' LSTART=',L1,' KF=',I4/' TEMP1=',G15.7,' TEMP2=',G15.7) ANOT 992 IF(NAB.GE.1) WRITE (O,95) (ABX(I),ABY(I),I=1,NAB) ANOT 993 95 FORMAT ((' ABX,ABY:',4(5X,2F10.4))) ANOT 994 WRITE (O,96) IF ANOT 995 96 FORMAT (' IF: ',20I5) ANOT 996 ASSIGN 97 TO KANT ANOT 997 GO TO 46 ANOT 998 97 RETURN ANOT 999 98 II=7771 ANOT1000 GO TO 92 ANOT1001 99 II=7773 ANOT1002 GO TO 92 ANOT1003 990 II=7774 ANOT1004 GO TO 92 ANOT1005 991 II=7775 ANOT1006 GO TO 92 ANOT1007 992 II=7776 ANOT1008 GO TO 92 ANOT1009 993 II=7777 ANOT1010 GO TO 92 ANOT1011 994 II=7778 ANOT1012 GO TO 92 ANOT1013 100 M=NOR(I,NJ,L) ANOT1014 NAB=NAB+1 ANOT1015 IF (NAB.EQ.LIMIT) GO TO 112 ANOT1016 ABX(NAB)=CROSS(1,M) ANOT1017 ABY(NAB)=CROSS(2,M) ANOT1018 IF (LSTART) GO TO 101 ANOT1019 IF (ABS(ABX(NAB)-ABX(NAB-1)).GT.FUDGE.OR. ABS(ABY(NAB)-ABY(NAB-1))ANOT1020 # .GT.FUDGE)GO TO 101 ANOT1021 IF(LBUG)WRITE(O,870)K,L,NI,NJ,NAB ANOT1022 NAB=NAB-1 ANOT1023 101 LSTART=.FALSE. ANOT1024 IF (LBUG2) WRITE (O,88) K,L,NI,NJ,I,NAB,ABX(NAB),ABX(NAB) ANOT1025 * ,ABY(NAB),ABY(NAB) ANOT1026 IF (L.EQ.1) GO TO 102 ANOT1027 JJ=ISON(LCR,NI,NJ,2) ANOT1028 IF (JJ.EQ.0) GO TO 104 ANOT1029 GO TO 103 ANOT1030 102 JJ=ISON(LCR,NJ,NI,2) ANOT1031 IF (JJ.EQ.0) GO TO 104 ANOT1032 103 J=NOR(1,NJ,L) ANOT1033 C ANOT1034 C WHICH WAY TO TURN? ANOT1035 C ANOT1036 M=SIGN(1.0,CROSS(1,J)*EQ(1,NI,3-L)+CROSS(2,J)*EQ(2,NI,3-L)-EQ(3,NIANOT1037 1,3-L)) ANOT1038 IF (M.EQ.0) GO TO 99 ANOT1039 NI=NJ ANOT1040 NJ=INO(I+M*TABYES(L,IT),NJ,L) ANOT1041 GO TO 85 ANOT1042 C ANOT1043 C FIND THE START OF ANOTHER CYCLE OF VERTICES. ANOT1044 C ANOT1045 104 CONTINUE ANOT1046 IF(LBUG2) WRITE(O,1041) NAB,KF,NI,NJ,I,J,K,L ANOT1047 1041 FORMAT(' AT START OF NEW SUBPOLYGON, NAB,KF,NI,NJ=',4I10 ANOT1048 * /' I,J,K,L=',4I10) ANOT1049 C ANOT1050 C DELETE POINTS ON THE MIDDLE OF STRAIGHT LINES OR ON ISOLATED SPIKEANOT1051 C ANOT1052 N=NAB-KF ANOT1053 C IF N<=2 THEN SUBPOLYGON HAS DEGENERATED TO NOTHING. ANOT1054 IF(N.GT.3)GO TO 1040 ANOT1055 1049 NAB=KF ANOT1056 GO TO 110 ANOT1057 1040 M=N ANOT1058 IF(1)=1 ANOT1059 DO 106 J=2,N ANOT1060 L=KF+J ANOT1061 I=L ANOT1062 1042 I=I-1 ANOT1063 IF(IF(I).EQ.0) GO TO 1042 ANOT1064 K=I+2 ANOT1065 IF (J.EQ.N) K=KF+2 ANOT1066 IF(I+1)=1 ANOT1067 C IF IF(I)=0 THEN DELETE THIS POINT ANOT1068 X=ABY(I)-ABY(K) ANOT1069 Y=ABX(K)-ABX(I) ANOT1070 Z=SQRT(X**2+Y**2) ANOT1071 IF (Z.LT.FUDGE) GO TO 105 ANOT1072 X=X/Z ANOT1073 Y=Y/Z ANOT1074 Z=-ABX(I)*X-ABY(I)*Y ANOT1075 IF (ABS(ABX(L )*X+ABY(L )*Y+Z).GT.FUDGE) GO TO 106 ANOT1076 105 IF(L )=0 ANOT1077 M=M-1 ANOT1078 106 CONTINUE ANOT1079 IF(LBUG2) WRITE(O,1060) M,N,KF,NAB,I,J,K,L,(IF(I),I=1,N) ANOT1080 1060 FORMAT(' SPIKETEST: M=',I5,' N=',I5,' KF=',I5,' NAB=',I5, ANOT1081 * ' I,J,K,L=',4I5/ ANOT1082 #(' IF=',10I5)) ANOT1083 IF(LBUG2) WRITE(O,1061) (ABX(I),ABY(I),I=1,NAB) ANOT1084 1061 FORMAT(' ABX,ABY: ',10F10.3) ANOT1085 C HAVE ANY POINTS BEEN DELETED ? ANOT1086 IF (M.EQ.N) GO TO 110 ANOT1087 IF(M.LE.3) GO TO 1049 ANOT1088 C YES. ANOT1089 J=KF+2 ANOT1090 DO 109 K=2,M ANOT1091 I=KF+K ANOT1092 107 IF (IF(J).NE.0) GO TO 108 ANOT1093 J=J+1 ANOT1094 GO TO 107 ANOT1095 108 IF (J.EQ.I) GO TO 109 ANOT1096 ABX(I)=ABX(J) ANOT1097 ABY(I)=ABY(J) ANOT1098 109 J=J+1 ANOT1099 NAB=KF+M ANOT1100 ABX(KF+1)=ABX(NAB) ANOT1101 ABY(KF+1)=ABY(NAB) ANOT1102 110 CONTINUE ANOT1103 C LSTART=.TRUE. FOR THE FIRST VERTEX OF A SUBPOLYGON. ANOT1104 LSTART=.TRUE. ANOT1105 C KF+1 IS THE FIRST VERTEX OF THE NEXT SUBPOLYGON. ANOT1106 KF=NAB ANOT1107 DO 111 NI=1,NA ANOT1108 DO 111 NJ=1,NB ANOT1109 JJ=ISON(LCR,NI,NJ,0) ANOT1110 IF (JJ.EQ.1) GO TO 83 ANOT1111 111 CONTINUE ANOT1112 IF(NAB.GT.3) GO TO 1101 ANOT1113 IF(LBUG2) WRITE(O,1102) ANOT1114 1102 FORMAT(' ***POLYGON HAS DEGENERATED TO NOTHING DUE TO FUDGE FACTORANOT1115 #. SET LMEET=.FALSE. AND TRY AGAIN.') ANOT1116 C IF NAB=0, THEN DUE TO POINTS CLOSER THAN THE FUDGE FACTOR BEING ANOT1117 C COALESCED, THE RESULTANT POLYGON HAS DEGENERATED TO NOTHING. THUS ANOT1118 C SET LMEET=.FALSE. AND RUN THE ALGORITHM AGAIN, THIS TIME ASSUMING ANOT1119 C THE POLYGONS DON'T MEET. ANOT1120 DO 1100 I=1,2 ANOT1121 DO 1100 J=1,MMEET ANOT1122 1100 LMEET(J,I)=.FALSE. ANOT1123 NI=0 ANOT1124 NJ=0 ANOT1125 GO TO 57 ANOT1126 C ANOT1127 C NEW POLYGON HAS BEEN FOUND COMPLETELY ANOT1128 C ANOT1129 1101 GO TO GO, (114,118,119) ANOT1130 112 CONTINUE ANOT1131 C ANOT1132 C THE NEW POLYGON IS TOO LARGE FOR ITS ARRAY. ANOT1133 C ANOT1134 WRITE (O,113) NAB,I,J,K,L,M,N,NA,NB,NI,NJ,NN,NOA,NOB,MCROSS,LIMIT ANOT1135 113 FORMAT (' NEW POLYGON FOUND BY ANOTB TOO LARGE FOR ITS ARRAY'/(10IANOT1136 110)) ANOT1137 WRITE (O,47) LCR ANOT1138 114 CONTINUE ANOT1139 IF (NAB.LE.1.OR.NCR.LE.1.OR..NOT.LBUG2) RETURN ANOT1140 DO 116 I=1,NAB ANOT1141 DO 115 J=1,NCR ANOT1142 IF (CROSS(1,J).EQ.ABX(I).AND.CROSS(2,J).EQ.ABY(I)) GO TO 116 ANOT1143 115 CONTINUE ANOT1144 GO TO 92 ANOT1145 116 IF (I) =J ANOT1146 WRITE (O,117) (IF(I),I=1,NAB) ANOT1147 117 FORMAT (/(25I5)) ANOT1148 RETURN ANOT1149 118 ASSIGN 114 TO GO ANOT1150 IT=8 ANOT1151 GO TO 120 ANOT1152 119 ASSIGN 114 TO GO ANOT1153 IT=4 ANOT1154 120 READ (KTEMP) LCR ANOT1155 IF (LBUG2) WRITE(O,121) (LCR(I),I=1,NANB) ANOT1156 121 FORMAT(' LCR AFTER READ2: ',10Z10) ANOT1157 NI=NIB ANOT1158 NJ=NJB ANOT1159 GO TO 83 ANOT1160 END ANOT1161 ANOT1162 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1163 C ANOT1164 FUNCTION LOR(I,J) ANOT1165 LOGICAL LBUG/.FALSE./ ANOT1166 EQUIVALENCE (K,A),(L,B) ANOT1167 K=I ANOT1168 L=J ANOT1169 A=OR(A,B) ANOT1170 IF(LBUG) WRITE(6,1) I,J,K ANOT1171 1 FORMAT(' *****LOR CALL, I=',Z10,' J=',Z10,' LOR=',Z10) ANOT1172 LOR=K ANOT1173 RETURN ANOT1174 ENTRY LAND(I,J) ANOT1175 K=I ANOT1176 L=J ANOT1177 A=AND(A,B) ANOT1178 LAND=K ANOT1179 IF(LBUG) WRITE(6,2) I,J,LAND ANOT1180 2 FORMAT(' *****LAND CALL, I=',Z10,' J=',Z10,' LAND=',Z10) ANOT1181 RETURN ANOT1182 END ANOT1183 ANOT1184 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1185 C ANOT1186 C >>>ISON *ANOT1187 C ..................................................................ANOT1188 C ANOT1189 C FUNCTION ISON ANOT1190 C ANOT1191 C PURPOSE ANOT1192 C TO MANIPULATE A 2-D BIT ARRAY FOR FORTRAN ANOT1193 C ANOT1194 C USAGE ANOT1195 C I = ISON (IARRAY, NI, NJ, IND) ANOT1196 C ANOT1197 C DESCRIPTION OF THE PARAMETERS ANOT1198 C IARRAY - ARRAY TO PROVIDE THE STORAGE FOR THE BIT ARRAY ANOT1199 C NI, NJ - IF IND<0, THEN THE DIMENSIONS OF THE BIT ARRAY. IARRAYANOT1200 C MUST BE LARGE ENOUGH TO CONTAIN NI*NJ BITS. ANOT1201 C IF IND>0, THEN (NI,NJ) IS THE BIT UNDER CONSIDERATION ANOT1202 C IND - AN INDICATOR OF THE OPERATION TO BE PERFORMED ANOT1203 C -2: INITIALIZE THE ARRAY TO 1'S. ANOT1204 C -1: INITIALIZE THE ARRAY TO 0'S. ANOT1205 C 0: RETURN 0 OR 1 DEPENDING ON THE (NI,NJ) ELEMENT. ANOT1206 C 1: SET THE ELEMENT TO 1, AND RETURN ITS ORIGINAL VALUANOT1207 C 2: SET THE ELEMENT TO 0, AND RETURN ITS ORIGINAL VALUANOT1208 C ANOT1209 C SUBPROGRAMS REQUIRED ANOT1210 C LOR,LAND - ASSEMBLER FUNCTIONS THAT TAKE 2 1-BYTE ARGUMENTS ANOT1211 C AND RETURN THEIR LOGICAL 'OR' OR 'AND'. ANOT1212 C ANOT1213 C MACHINE DEPENDENCE ANOT1214 C ISON ASSUMES A WORD SIZE OF 32 BITS AND A BYTE SIZE OF 8 BITS. ANOT1215 C IN ANY NON-IBM SYSTEM, ISON WOULD HAVE TO BE CHANGED, AND LOR ANOT1216 C AND LAND REWRITTEN. ALSO THE SIZE OF IARRAY PASSED TO ISON ANOT1217 C WOULD HAVE TO BE CHANGED OF THE WORD SIZE WERE DIFFERENT. ANOT1218 C ANOT1219 C ..................................................................ANOT1220 C ANOT1221 FUNCTION ISON (IARRAY,NI,NJ,IND) ANOT1222 DIMENSION IARRAY( 313) ANOT1223 DIMENSION IMASK(32), JMASK(32) ANOT1224 LOGICAL LBUG ANOT1225 EXTERNAL LOR,LAND ANOT1226 DATA IMASK/Z1,Z2,Z4,Z8,Z10,Z20,Z40,Z80,Z100,Z200,Z400,Z800,Z1000,ZANOT1227 12000,Z4000,Z8000,Z10000,Z20000,Z40000,Z80000,Z100000,Z200000,Z4000ANOT1228 200,Z800000,Z1000000,Z2000000,Z4000000,Z8000000,Z10000000,Z20000000ANOT1229 3,Z40000000,Z80000000/,JMASK/ZFFFFFFFE,ZFFFFFFFD,ZFFFFFFFB,ZFFFFFFFANOT1230 47,ZFFFFFFEF,ZFFFFFFDF,ZFFFFFFBF,ZFFFFFF7F,ZFFFFFEFF,ZFFFFFDFF,ZFFFANOT1231 5FFBFF,ZFFFFF7FF,ZFFFFEFFF,ZFFFFDFFF,ZFFFFBFFF,ZFFFF7FFF,ZFFFEFFFF,ANOT1232 6ZFFFDFFFF,ZFFFBFFFF,ZFFF7FFFF,ZFFEFFFFF,ZFFDFFFFF,ZFFBFFFFF,ZFF7FFANOT1233 7FFF,ZFEFFFFFF,ZFDFFFFFF,ZFBFFFFFF,ZF7FFFFFF,ZEFFFFFFF,ZDFFFFFFF,ZBANOT1234 8FFFFFFF,Z7FFFFFFF/,LBUG/.FALSE./ ANOT1235 IF (LBUG) WRITE (6,1) NI,NJ,IND ANOT1236 1 FORMAT (' *****ISON CALL NI=',I4,' NJ=',I4,' IND=',I4) ANOT1237 ISON=0 ANOT1238 KIND=IND+3 ANOT1239 GO TO (2,3,7,7,7), KIND ANOT1240 C ILLEGAL IND ANOT1241 GO TO 10 ANOT1242 2 K=-1 ANOT1243 GO TO 4 ANOT1244 3 K=0 ANOT1245 4 NIJ=(NI*NJ+31)/32 ANOT1246 DO 5 I=1,NIJ ANOT1247 IARRAY(I)=K ANOT1248 5 CONTINUE ANOT1249 NN=NI ANOT1250 IF (LBUG) WRITE (6,6) (IARRAY(I),I=1,NIJ) ANOT1251 6 FORMAT (' *****ISON IARRAY:',10Z10) ANOT1252 GO TO 10 ANOT1253 7 IPOS=(NJ-1)*NN+NI-1 ANOT1254 IWORD=IPOS/32 ANOT1255 IBIT=IPOS-IWORD*32+1 ANOT1256 IWORD=IWORD+1 ANOT1257 ITEMP=LAND(IARRAY(IWORD),IMASK(IBIT)) ANOT1258 IF (ITEMP.NE.0) ISON=1 ANOT1259 GO TO (10,10,10,8,9), KIND ANOT1260 C SET ELEMENT = 1 ANOT1261 8 IARRAY(IWORD)=LOR(IARRAY(IWORD),IMASK(IBIT)) ANOT1262 GO TO 10 ANOT1263 C SET ELEMENT = 0 ANOT1264 9 IARRAY(IWORD)=LAND(IARRAY(IWORD),JMASK(IBIT)) ANOT1265 10 CONTINUE ANOT1266 IF (LBUG) WRITE (6,11) ISON,NN ANOT1267 11 FORMAT (' *****ISON RETURN ISON=',I4,' NN=',I4) ANOT1268 RETURN ANOT1269 END ANOT1270 ANOT1271 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1272 C ANOT1273 C>>>PNPOLX 17:25:49 07/26/73 (73.207) ANOT1274 C ..................................................................ANOT1275 C ANOT1276 C SUBROUTINE PNPOLY ANOT1277 C ANOT1278 C PURPOSE ANOT1279 C TO DETERMINE WHETHER A POINT IS INSIDE A POLYGON ANOT1280 C ANOT1281 C USAGE ANOT1282 C CALL PNPOLY (PX, PY, X, Y, N, INOUT ) ANOT1283 C ANOT1284 C DESCRIPTION OF THE PARAMETERS ANOT1285 C PX - X-COORDINATE OF POINT IN QUESTION. ANOT1286 C PY - Y-COORDINATE OF POINT IN QUESTION. ANOT1287 C X - N LONG VECTOR CONTAINING X-COORDINATES OF ANOT1288 C VERTICES OF POLYGON. ANOT1289 C Y - N LONG VECTOR CONTAINING Y-COORDINATES OF ANOT1290 C VERTICES OF POLYGON. ANOT1291 C N - NUMBER OF VERTICES IN THE POLYGON. ANOT1292 C INOUT - THE SIGNAL RETURNED: ANOT1293 C -1 IF THE POINT IS OUTSIDE OF THE POLYGON, ANOT1294 C 0 IF THE POINT IS ON AN EDGE OR AT A VERTEX, ANOT1295 C 1 IF THE POINT IS INSIDE OF THE POLYGON. ANOT1296 C ANOT1297 C REMARKS ANOT1298 C THE VERTICES MAY BE LISTED IN CLOCKWISE OR ANTICLOCKWISE ANOT1299 C ORDER. FOR THIS SUBROUTINE A POINT IS CONSIDERED INSIDE ANOT1300 C THE POLYGON IF IT IS LOCATED IN THE ENCLOSED AREA DEFINED ANOT1301 C BY THE LINE FORMING THE POLYGON. ANOT1302 C THE FIRST POINT MAY OPTIONALLY BE REPEATED, IF SO N MAY ANOT1303 C OPTIONALLY BE INCREASED BY 1. ANOT1304 C THE INPUT POLYGON MAY BE A COMPOUND POLYGON CONSISTING ANOT1305 C OF SEVERAL SEPARATE SUBPOLYGONS. IF SO, THE FIRST VERTEX ANOT1306 C OF EACH SUBPOLYGON MUST BE REPEATED, AND WHEN CALCULATING ANOT1307 C N, THESE FIRST VERTICES MUST BE COUNTED TWICE. ANOT1308 C INOUT IS THE ONLY PARAMETER WHOSE VALUE IS CHANGED. ANOT1309 C PNPOLY CAN HANDLE ANY NUMBER OF VERTICES IN THE POLYGON. ANOT1310 C ANOT1311 C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED ANOT1312 C NONE ANOT1313 C ANOT1314 C METHOD ANOT1315 C A VERTICAL SEMI-INFINITE LINE IS DRAWN UP FROM THE POINT ANOT1316 C IN QUESTION. IF IT CROSSES THE POLYGON AN ODD NUMBER OF ANOT1317 C TIMES, THEN THE POINT IS INSIDE THE POLYGON. ANOT1318 C WRITTEN BY W. RANDOLPH FRANKLIN ANOT1319 C THIS PNPOLY CONTAINS ADDITIONS FOR VIEWPLOT. ANOT1320 C ANOT1321 C ..................................................................ANOT1322 C ANOT1323 SUBROUTINE PNPOLX (PX,PY,X,Y,N,INOUT,P,ID,PIN) ANOT1324 C ANOT1325 COMMON/COMMAX/ MXV,MXE,MXF,MXLF,MXLIVI,MXRSUP,MXISUP,MXISP2,MXARTAANOT1326 C ANOT1327 REAL X(N),Y(N),P(3,MXV),PIN(2) ANOT1328 INTEGER*2 ID(N) ANOT1329 LOGICAL IX,IY,JX,JY,EOR ANOT1330 C EXCLUSIVE OR FUNCTION. ANOT1331 EOR(IX,IY)=(IX.OR.IY).AND..NOT.(IX.AND.IY) ANOT1332 6 INOUT=-1 ANOT1333 DO 4 I=1,N ANOT1334 GO TO IGO,(5,7) ANOT1335 5 XI=X(I)-PX ANOT1336 YI=Y(I)-PY ANOT1337 GO TO 10 ANOT1338 7 XI=P(1,ID(I))-PIN(1) ANOT1339 YI=P(2,ID(I))-PIN(2) ANOT1340 C CHECK WHETHER THE POINT IN QUESTION IS AT THIS VERTEX. ANOT1341 10 IF (XI.EQ.0.0.AND.YI.EQ.0.0) GO TO 2 ANOT1342 C J IS NEXT VERTEX NUMBER OF POLYGON. ANOT1343 J=1+MOD(I,N) ANOT1344 GO TO JGO,(8,9) ANOT1345 8 XJ=X(J)-PX ANOT1346 YJ=Y(J)-PY ANOT1347 GO TO 11 ANOT1348 9 XJ=P(1,ID(J))-PIN(1) ANOT1349 YJ=P(2,ID(J))-PIN(2) ANOT1350 C IS THIS LINE OF 0 LENGTH ? ANOT1351 11 IF (XI.EQ.XJ.AND.YI.EQ.YJ) GO TO 4 ANOT1352 IX=XI.GE.0.0 ANOT1353 IY=YI.GE.0.0 ANOT1354 JX=XJ.GE.0.0 ANOT1355 JY=YJ.GE.0.0 ANOT1356 C CHECK WHETHER (PX,PY) IS ON VERTICAL SIDE OF POLYGON. ANOT1357 IF (XI.EQ.0.0.AND.XJ.EQ.0.0.AND.EOR(IY,JY)) GO TO 2 ANOT1358 C CHECK WHETHER (PX,PY) IS ON HORIZONTAL SIDE OF POLYGON. ANOT1359 IF (YI.EQ.0.0.AND.YJ.EQ.0.0.AND.EOR(IX,JX)) GO TO 2 ANOT1360 C CHECK WHETHER BOTH ENDS OF THIS SIDE ARE COMPLETELY 1) TO RIGHT ANOT1361 C OF, 2) TO LEFT OF, OR 3) BELOW (PX,PY). ANOT1362 IF (.NOT.((IY.OR.JY).AND.EOR(IX,JX))) GO TO 4 ANOT1363 C DOES THIS SIDE OBVIOUSLY CROSS LINE RISING VERTICALLY FROM (PX,PY)ANOT1364 IF (.NOT.(IY.AND.JY.AND.EOR(IX,JX))) GO TO 1 ANOT1365 INOUT=-INOUT ANOT1366 GO TO 4 ANOT1367 1 IF ((YI*XJ-XI*YJ)/(XJ-XI)) 4,2,3 ANOT1368 2 INOUT=0 ANOT1369 RETURN ANOT1370 3 INOUT=-INOUT ANOT1371 4 CONTINUE ANOT1372 RETURN ANOT1373 C ANOT1374 C ANOT1375 ENTRY PNPOLY (PX,PY,X,Y,N,INOUT) ANOT1376 C ANOT1377 C ANOT1378 ASSIGN 5 TO IGO ANOT1379 ASSIGN 8 TO JGO ANOT1380 GO TO 6 ANOT1381 C ANOT1382 C ANOT1383 ENTRY PNP2 (PIN,P,ID,N,INOUT) ANOT1384 C ANOT1385 C ANOT1386 ASSIGN 7 TO IGO ANOT1387 ASSIGN 9 TO JGO ANOT1388 GO TO 6 ANOT1389 END ANOT1390 ANOT1391 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1392 C ANOT1393 C ..................................................................ANOT1394 C ANOT1395 C FUNCTION AREA ANOT1396 C ANOT1397 C PURPOSE ANOT1398 C TO FIND THE SIGNED AREA OF A POLYGON. ANOT1399 C ANOT1400 C USAGE ANOT1401 C A = AREA (X, Y, N) ANOT1402 C ANOT1403 C DESCRIPTION OF THE PARAMETERS ANOT1404 C X - N LONG VECTOR CONTAINING X-COORDINATES OF THE ANOT1405 C VERTICES OF THE POLYGON. ANOT1406 C Y - N LONG VECTOR CONTAINING Y-COORDINATES OF THEM. ANOT1407 C N - NUMBER OF VERTICES. ANOT1408 C ANOT1409 C REMARKS ANOT1410 C THE FIRST VERTEX MAY OPTIONALLY BE REPEATED. IF SO, N MAY ANOT1411 C OPTIONALLY BE INCREASED BY 1. ANOT1412 C IF THE VERTICES ARE LISTED CLOCKWISE, THE AREA RETURNED ANOT1413 C WILL BE POSITIVE; IF THEY ARE LISTED COUNTERCLOCKWISE, ANOT1414 C THE AREA RETURNED WILL BE NEGATIVE. ANOT1415 C THUS THIS FUNCTION MAY BE USED TO DETERMINE WHICH WAY ANOT1416 C THE VERTICES ARE LISTED. ANOT1417 C THE INPUT POLYGON MAY BE A COMPOUND POLYGON, CONSISTING OF ANOT1418 C SEVERAL DISJOINT SUBPOLYGONS. IN THIS CASE THE FIRST VERTEX ANOT1419 C OF EACH SUBPOLYGON MUST BE REPEATED, AND ALL THESE FIRST ANOT1420 C VERTICES MUST BE COUNTED TWICE WHEN CALCULATING N. ANOT1421 C IN THIS CASE THE AREA RETURNED WILL BE THE SUM OF THE AREAS ANOT1422 C OF THE SEPARATE SUBPOLYGONS. ANOT1423 C WRITTEN BY W. RANDOLPH FRANKLIN ANOT1424 C ANOT1425 C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED ANOT1426 C NONE ANOT1427 C ANOT1428 C METHOD ANOT1429 C THE POLYGON IS BROKEN INTO QUADRILATERALS WITH A POSSIBLE ANOT1430 C EXTRA TRIANGLE, AND THE AREAS OF THESE CALCULATED. ANOT1431 C ANOT1432 C THIS AREA CONTAINS ADDITIONS FOR VIEWPLOT. ANOT1433 C ..................................................................ANOT1434 C ANOT1435 FUNCTION AREAXX (X,Y,N,P,ID) ANOT1436 REAL X(N),Y(N),P(3,50) ANOT1437 INTEGER*2 ID(N) ANOT1438 C ANOT1439 C ANOT1440 ENTRY AREA (X,Y,N) ANOT1441 C ANOT1442 C ANOT1443 A=0.0 ANOT1444 W=X(1) ANOT1445 Z=Y(1) ANOT1446 IF (N-3) 4,3,1 ANOT1447 1 DO 2 I=4,N,2 ANOT1448 A=A+(X(I-1)-W)*(Y(I-2)-Y(I))-(X(I-2)-X(I))*(Y(I-1)-Z) ANOT1449 2 CONTINUE ANOT1450 3 IF (MOD(N,2).EQ.1) A=A+(X(N)-W)*(Y(N-1)-Z)-(X(N-1)-W)*(Y(N)-Z) ANOT1451 A=A/2.0 ANOT1452 AREA=A ANOT1453 RETURN ANOT1454 4 AREA=0.0 ANOT1455 RETURN ANOT1456 C ANOT1457 C ANOT1458 ENTRY AREA2 (P,ID,N) ANOT1459 C ANOT1460 C ANOT1461 A=0.0 ANOT1462 W=P(1,ID(1)) ANOT1463 Z=P(2,ID(1)) ANOT1464 IF(N-3) 4,13,11 ANOT1465 11 DO 12 I=4,N,2 ANOT1466 A=A+(P(1,ID(I-1))-W)*(P(2,ID(I-2))-P(2,ID(I)))-(P(1,ID(I-2))-P(1,IANOT1467 @D ANOT1468 @(I)))*(P(2,ID(I-1))-Z) ANOT1469 12 CONTINUE ANOT1470 13 IF (MOD(N,2).EQ.1) A=A+(P(1,ID(N))-W)*(P(2,ID(N-1))-Z)-(P(1,ID(N-1ANOT1471 @))-W)*(P(2,ID(N))-Z) ANOT1472 AREA=A/2.0 ANOT1473 RETURN ANOT1474 END ANOT1475 ANOT1476 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1477 C ANOT1478 LOR CSECT ANOT1479 USING LOR,15 ANOT1480 SAVE (14,10),T,* ANOT1481 L 2,0(0,1) LOAD ADDRESS OF ARGUMENT LIST IN R2 ANOT1482 L 0,0(2) ANOT1483 L 2,4(0,1) LOAD ADDRESS OF SECOND ARGUMENT INTO R2 ANOT1484 O 0,0(0,2) LOGICAL OR ANOT1485 LM 1,10,24(13) ANOT1486 RETURN (14,15),T,RC=0 ANOT1487 LAND CSECT ANOT1488 USING LAND,15 ANOT1489 SAVE (14,10),T,* ANOT1490 L 2,0(0,1) ANOT1491 L 0,0(0,2) ANOT1492 L 2,4(0,1) ANOT1493 N 0,0(0,2) LOGICAL AND ANOT1494 LM 1,10,24(13) ANOT1495 RETURN (14,15),T,RC=0 ANOT1496 END ANOT1497 ANOT1498 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1499 C ANOT1500 //A472ANOT JOB (1406,R3084),'A¬B A¬B A¬B A¬B ' ANOT1501 //* TEST ANOTB 16:58:47 07/18/73 (73.199) ANOT1502 //WATFIV PROC PRINT=?,PUNCH=D,NOPUNCH='DUMMY,',LIB='AS.A1195.WATLIB' ANOT1503 //* RANDOLPH FRANKLIN, GEOGRAPHY 291-4227 ANOT1504 //GO EXEC PGM=WATFIV,REGION=180K ANOT1505 //SYSPRINT DD SYSOUT=&PRINT ANOT1506 //WATLIB DD DSN=AS.A1195.FUNLIB,DISP=SHR,DCB=BLKSIZE=3520 ANOT1507 // DD DSN=&LIB,DISP=SHR ANOT1508 // DD DSN=SC.A1406.PDS,DISP=SHR ANOT1509 // DD DSN=SC.A1406.VPLT,DISP=SHR ANOT1510 //FT01F001 DD SPACE=(1024,(120,120)),UNIT=SYSDA,DSN=&SYSUT1 ANOT1511 //FT02F001 DD SPACE=(1024,(120,120)),UNIT=SYSDA,DSN=&SYSUT2 ANOT1512 //FT03F001 DD SPACE=(1024,(120,120)),UNIT=SYSDA,DSN=&SYSUT3 ANOT1513 //FT04F001 DD SPACE=(TRK,(1,1)),UNIT=SYSDA,DSN=&SYSUT4, ANOT1514 // DCB=(LRECL=1254,BLKSIZE=1258,RECFM=VBS) ANOT1515 //FT05F001 DD DDNAME=SYSIN ANOT1516 //FT06F001 DD SYSOUT=&PRINT ANOT1517 //FT07F001 DD &NOPUNCH.SYSOUT=&PUNCH ANOT1518 // PEND ANOT1519 //VPLT EXEC WATFIV ANOT1520 //FT05F001 DD * ANOT1521 $JOB 1406/FRANKLINRANDOLPH,TIME=1,PAGES=50,RUN=FREE,LIBLIST ANOT1522 C TEST ANOTB 17:00:36 07/18/73 (73.199) ANOT1523 REAL AX(99),AY(99),BX(99),BY(99),ABX(300),ABY(300) ANOT1524 INTEGER KIT(3)/1,4,7/ ANOT1525 LOGICAL LBUG/.TRUE./ ANOT1526 DO 1 INPUT=1,100 ANOT1527 WRITE(6,5) INPUT ANOT1528 5 FORMAT('1**********NEW SET OF INPUT POLYGONS #',I2,'**********') ANOT1529 READ,NOA,NOB ANOT1530 DO 7 I=1,NOA ANOT1531 7 READ,AX(I),AY(I) ANOT1532 DO8 I=1,NOB ANOT1533 8 READ,BX(I),BY(I) ANOT1534 WRITE(6,2) NOA,(AX(I),AY(I),I=1,NOA) ANOT1535 2 FORMAT('-NOA=',I3,' AX,AY='/(2G15.7)) ANOT1536 WRITE(6,3) NOB,(BX(I),BY(I),I=1,NOB) ANOT1537 3 FORMAT('-NOB=',I3,' BX,BY='/(2G15.7)) ANOT1538 DO 1 K=1,3 ANOT1539 ITT=KIT(K) ANOT1540 WRITE(6,6)ITT ANOT1541 6 FORMAT('-----------ITT=',I2,'----------') ANOT1542 CALL ANOTB(AX,AY,NOA,BX,BY,NOB,ABX,ABY,NAB,ITT,3.00) ANOT1543 IF(NAB.GE.1) WRITE(6,4) NAB,(ABX(I),ABY(I),I=1,NAB) ANOT1544 4 FORMAT('-NAB=',I3,' ABX,ABY='/(2G15.7)) ANOT1545 1 CONTINUE ANOT1546 RETURN ANOT1547 C DATA:15.36399 5.117999 ANOT1548 C DATA:14.63399 4.879999 ANOT1549 C DATA:13.65899 7.880999 ANOT1550 C DATA:13.50000 7.828999 ANOT1551 C DATA:13.76299 8.802999 ANOT1552 C DATA:14.54999 8.169999 ANOT1553 C DATA:14.38899 8.117999 ANOT1554 C DATA:15.36399 5.117999 ANOT1555 END ANOT1556 $ENTRY ANOT1557 4 9 ANOT1558 0 0 ANOT1559 0 1 ANOT1560 0 1 ANOT1561 1 0 ANOT1562 0 0 ANOT1563 0 1 ANOT1564 0 1 ANOT1565 1 0 ANOT1566 0 0 ANOT1567 0 0 ANOT1568 0 1 ANOT1569 1 1 ANOT1570 1 1 ANOT1571 8 10 ANOT1572 15.31999 5.183999 ANOT1573 14.67799 4.813999 ANOT1574 13.09199 7.560999 ANOT1575 12.94899 7.477999 ANOT1576 13.00000 8.462999 ANOT1577 13.87699 8.013999 ANOT1578 13.73399 7.930999 ANOT1579 15.31999 5.183999 ANOT1580 18.31898 6.078999 ANOT1581 18.91198 5.830999 ANOT1582 18.09028 4.957999 ANOT1583 18.05918 5.149999 ANOT1584 15.09999 4.522999 ANOT1585 14.99899 4.998999 ANOT1586 15.19199 4.565999 ANOT1587 17.96299 5.799999 ANOT1588 18.04199 5.625000 ANOT1589 18.31898 6.078999 ANOT1590 5 3 ANOT1591 0 0 ANOT1592 1 0 ANOT1593 1 1 ANOT1594 0 0 ANOT1595 0 0 ANOT1596 2 2 ANOT1597 2 3 ANOT1598 3 2 ANOT1599 6 1 ANOT1600 0 0 ANOT1601 1 0 ANOT1602 0 1 ANOT1603 0 0 ANOT1604 0 0 ANOT1605 0 0 ANOT1606 1 1 ANOT1607 9 10 ANOT1608 2 2 ANOT1609 2 3 ANOT1610 3 2 ANOT1611 2 2 ANOT1612 2 2 ANOT1613 2 2 ANOT1614 2 1 ANOT1615 1 1 ANOT1616 2 2 ANOT1617 5 6 ANOT1618 5 6 ANOT1619 6 5 ANOT1620 5 5 ANOT1621 7 7 ANOT1622 7 8 ANOT1623 7 7 ANOT1624 8 1 ANOT1625 1 2 ANOT1626 1 1 ANOT1627 12 1 ANOT1628 2 2 ANOT1629 2 3 ANOT1630 3 2 ANOT1631 2 2 ANOT1632 5 5 ANOT1633 5 6 ANOT1634 5 5 ANOT1635 6 6 ANOT1636 6 7 ANOT1637 6 6 ANOT1638 7 7 ANOT1639 7 8 ANOT1640 0 0 ANOT1641 $STOP ANOT1642 ANOT1643 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1644 C ANOT1645 //DOUANOTB JOB (1406,R3084),ANOTB.FORTG,MSGCLASS=R ANOT1646 //* TEST ANOTB 16:58:47 07/18/73 (73.199) ANOT1647 //A EXEC FORTGCLG,CPARM='ID,MAP', ANOT1648 // REGION=120K,LPARM='LET,MAP', ANOT1649 // LIB='SC.A1406.LOAD' ANOT1650 //FORT.SYSIN DD * ANOT1651 C TEST ANOTB 17:00:36 07/18/73 (73.199) ANOT1652 REAL AX(99),AY(99),BX(99),BY(99),ABX(300),ABY(300) ANOT1653 INTEGER KIT(3)/1,4,7/ ANOT1654 LOGICAL LBUG/.TRUE./ ANOT1655 DO 1 INPUT=1,100 ANOT1656 WRITE(6,5) INPUT ANOT1657 5 FORMAT('1**********NEW SET OF INPUT POLYGONS #',I2,'**********') ANOT1658 READ(5,10)NOA,NOB ANOT1659 10 FORMAT(2I5) ANOT1660 DO 7 I=1,NOA ANOT1661 7 READ(5,12)AX(I),AY(I) ANOT1662 12 FORMAT(2F10.4) ANOT1663 DO8 I=1,NOB ANOT1664 8 READ(5,12)BX(I),BY(I) ANOT1665 WRITE(6,2) NOA,(AX(I),AY(I),I=1,NOA) ANOT1666 2 FORMAT('-NOA=',I3,' AX,AY='/(2G15.7)) ANOT1667 WRITE(6,3) NOB,(BX(I),BY(I),I=1,NOB) ANOT1668 3 FORMAT('-NOB=',I3,' BX,BY='/(2G15.7)) ANOT1669 DO 1 K=1,3 ANOT1670 ITT=KIT(K) ANOT1671 WRITE(6,6)ITT ANOT1672 6 FORMAT('-----------ITT=',I2,'----------') ANOT1673 CALL ANOTB(AX,AY,NOA,BX,BY,NOB,ABX,ABY,NAB,ITT,3.) ANOT1674 IF(NAB.GE.1) WRITE(6,4) NAB,(ABX(I),ABY(I),I=1,NAB) ANOT1675 4 FORMAT('-NAB=',I3,' ABX,ABY='/(2G15.7)) ANOT1676 1 CONTINUE ANOT1677 RETURN ANOT1678 C DATA:15.36399 5.117999 ANOT1679 C DATA:14.63399 4.879999 ANOT1680 C DATA:13.65899 7.880999 ANOT1681 C DATA:13.50000 7.828999 ANOT1682 C DATA:13.76299 8.802999 ANOT1683 C DATA:14.54999 8.169999 ANOT1684 C DATA:14.38899 8.117999 ANOT1685 C DATA:15.36399 5.117999 ANOT1686 END ANOT1687 //GO.SYSIN DD * ANOT1688 8 10 ANOT1689 15.31999 5.183999 ANOT1690 14.67799 4.813999 ANOT1691 13.09199 7.560999 ANOT1692 12.94899 7.477999 ANOT1693 13.00000 8.462999 ANOT1694 13.87699 8.013999 ANOT1695 13.73399 7.930999 ANOT1696 15.31999 5.183999 ANOT1697 18.31898 6.078999 ANOT1698 18.91198 5.830999 ANOT1699 18.09028 4.957999 ANOT1700 18.05918 5.149999 ANOT1701 15.09999 4.522999 ANOT1702 14.99899 4.998999 ANOT1703 15.19199 4.565999 ANOT1704 17.96299 5.799999 ANOT1705 18.04199 5.625000 ANOT1706 18.31898 6.078999 ANOT1707 5 3 ANOT1708 0.0 0.0 ANOT1709 1.0 0.0 ANOT1710 1.0 1.0 ANOT1711 0.0 0.0 ANOT1712 0.0 0.0 ANOT1713 2.0 2.0 ANOT1714 2.0 3.0 ANOT1715 3.0 2.0 ANOT1716 6 1 ANOT1717 0.0 0.0 ANOT1718 1.0 0.0 ANOT1719 0.0 1.0 ANOT1720 0.0 0.0 ANOT1721 0.0 0.0 ANOT1722 0.0 0.0 ANOT1723 1.0 1.0 ANOT1724 9 10 ANOT1725 2.0 2.0 ANOT1726 2.0 3.0 ANOT1727 3.0 2.0 ANOT1728 2.0 2.0 ANOT1729 2.0 2.0 ANOT1730 2.0 2.0 ANOT1731 2.0 1.0 ANOT1732 1.0 1.0 ANOT1733 2.0 2.0 ANOT1734 5.0 6.0 ANOT1735 5.0 6.0 ANOT1736 6.0 5.0 ANOT1737 5.0 5.0 ANOT1738 7.0 7.0 ANOT1739 7.0 8.0 ANOT1740 7.0 7.0 ANOT1741 8.0 1.0 ANOT1742 1.0 2.0 ANOT1743 1.0 1.0 ANOT1744 12 1 ANOT1745 2.0 2.0 ANOT1746 2.0 3.0 ANOT1747 3.0 2.0 ANOT1748 2.0 2.0 ANOT1749 5.0 5.0 ANOT1750 5.0 6.0 ANOT1751 5.0 5.0 ANOT1752 6.0 6.0 ANOT1753 6.0 7.0 ANOT1754 6.0 6.0 ANOT1755 7.0 7.0 ANOT1756 7.0 8.0 ANOT1757 0.0 0.0 ANOT1758 //* ANOT1759 ANOT1760 C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ANOT1761 C ANOT1762