C 00000010 C 00000020 C 00000030 C THIS PROGRAM IS A COMPILATION OF PROGRAMS PRESENTED IN 00000040 C USCGS TECHNICAL BULLETINS NO 21 23 29 34 AND 00000050 C A PROGRAM BY SCHUT OF NRC COMPILATION BY CARLIN 1969 00000060 C 00000070 C S T R I P A E R O T R I A N G U L A T I O N CA00000080 C 00000090 IMPLICIT REAL*8(A-H,O-Z) 00000100 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00000110 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00000120 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00000130 2JCON,INSTRU 00000140 DEFINE FILE 13(800,25,U,KK3) 00000150 JS1=1 00000160 JS2=2 00000170 KD3=13 00000180 C DATA SET REFERENCE NUMBERS 00000190 NREAD=5 00000200 NPRNT=6 00000210 NPUNC=7 00000220 300 READ (NREAD,70)IJOB 00000230 C JOB CONTROL CARD 00000240 IF(IJOB.EQ.9)GO TO 71 00000250 IF(IJOB.EQ.8)CALL EXIT 00000260 GO TO 300 00000270 71 I50=0 00000280 REWIND JS1 00000290 READ (NREAD,122) 00000310 C PROJECT NAME 00000320 WRITE (NPRNT,520) 00000330 ISTP=0 00000340 WRITE (NPRNT,122) 00000350 READ (NREAD,152)IS1A0,IPOLY,IEARTH,ICOORD,IFROM,ITO,ISHAR,ITURN, 00000360 1INSTRU 00000370 C INSTRUCTION PARAMETERS 00000380 WRITE(NPRNT,152)IS1A0,IPOLY,IEARTH,ICOORD,IFROM,ITO,ISHAR,ITURN, 00000390 1INSTRU 00000400 IF(IEARTH.EQ.1)GO TO 402 00000410 IF(ISHAR.EQ.1)GO TO 80 00000420 C IF ISHAR EQUALS 1 RESET FOR POLYNOMIAL ADJUSTMENT 00000430 204 IF(IS1A0.EQ.0)GO TO 350 00000440 IF(IS1A0.EQ.1)GO TO 385 00000450 IF(IS1A0.EQ.2)GO TO 60 00000460 IF(ICOORD.EQ.2.AND.IS1A0.EQ.3)READ(NREAD,151)N5 00000470 C IF ICOORD EQUALS 2 AND IS1A0 EQUALS 3 RESET FOR COORD. TRANSFORMA00000480 IF(ICOORD.EQ.2)CALL TCOORD 00000490 GO TO 300 00000500 350 IF(ITURN.EQ.0)GO TO 102 00000510 IGO=1 00000520 GO TO 104 00000530 C TRIPLET ASSEMBLY OF STRIP 00000540 102 CALL PART1 00000550 201 WRITE (NPRNT,200)(ISTRP(I),IRCSP(I),I=1,ISTP) 00000560 END FILE JS1 00000570 IF(I50.EQ.0)GO TO 401 00000580 WRITE (NPRNT,151)(ICAN(I),I=1,I50) 00000590 GO TO 401 00000600 385 IF(ITURN.EQ.0)GO TO 384 00000610 IGO=2 00000620 GO TO 104 00000630 C INPUT STRIP MODEL COORDINATES 00000640 384 IREC=0 00000650 READ (NREAD,151)ISTRIP 00000660 IF(ISTRIP.EQ.0)GO TO 201 00000670 ISTP=ISTP+1 00000680 383 READ (NREAD,381)IFK,X,Y,D,IT 00000690 IREC=IREC+1 00000700 WRITE (JS1)IFK,X,Y,D 00000710 IF (IT.GT.0)GO TO 382 00000720 GO TO 383 00000730 382 ISTRP(ISTP)=ISTRIP 00000740 IRCSP(ISTP)=IREC 00000750 IF(IREC.GT.500)GO TO 399 00000760 GO TO 384 00000770 80 READ (NREAD,151)ISTP,I50 00000780 C RERUN STORED DATA 00000790 READ (NREAD,200)(ISTRP(I),IRCSP(I),I=1,ISTP) 00000800 IF(I50.EQ.0)GO TO 203 00000810 READ (NREAD,151)(ICAN(I),I=1,I50) 00000820 203 IF(ITURN.EQ.0)GO TO 401 00000830 C IF ITURN EQUAL 1 CALCULATE ADDITIONAL MODEL COORDINATES 00000840 GO TO 204 00000850 104 ITOT=0 00000860 DO 100 I=1,ISTP 00000870 100 ITOT=ITOT+IRCSP(I) 00000880 C SET DATA SET NUMBER JS1 FOR ADDITIONAL STORAGE 00000890 DO 101 I=1,ITOT 00000900 101 READ (JS1,ERR=300) 00000910 GO TO (102,384,103),IGO 00000920 60 IF(ITURN.EQ.0)GO TO 103 00000930 IGO=3 00000940 GO TO 104 00000950 C INDEPENDENT MODEL ASSEMBLY BY SCHUT METHOD 00000960 103 CALL INMODL 00000970 GO TO 201 00000980 C INPUT PROJECT GROUND CONTROL 00000990 401 I=0 00001000 403 I=I+1 00001010 READ (NREAD,404)IGRND(I),X,Y,D,IT 00001020 WRITE (KD3'I)X,Y,D 00001030 IF(IT.GT.0)GO TO 405 00001040 GO TO 403 00001050 405 JCON=I 00001060 IF(JCON.GT.300)CALL EXIT 00001070 IN=0 00001080 C INPUT STRIP TIE POINT STATION NUMBERS 00001090 450 READ (NREAD,151)(ITIE(I),I=101,107) 00001100 DO 452 I=101,107 00001110 IF(ITIE(I).LE.0)GO TO 453 00001120 IN=IN+1 00001130 452 ITIE(IN)=ITIE(I) 00001140 GO TO 450 00001150 453 IF(IN.GT.100)GO TO 301 00001160 ICON=IN 00001170 402 IF(IPOLY.EQ.1)CALL ASTRIP 00001180 C POLYNOMIAL ADJUSTMENT OF STRIP BY USC&GS METHOD 00001190 WRITE (NPRNT,200)JS2,N7 00001200 105 IF(ICOORD.EQ.1)CALL TCOORD 00001210 C TRANSFORMATION OF GROUND COORDINATE SYSTEMS USC&GS; DNR 00001220 387 CONTINUE GO TO 300 00001260 399 WRITE (NPRNT,398) 00001270 I50=I50+1 00001280 IRCSP(ISTP)=IREC 00001290 C TOO MANY POINTS FOR STRIP ADJUSTMENT PROGRAM 00001300 ICAN(I50)=ISTRIP 00001310 GO TO 384 00001320 301 WRITE (NPRNT,302) 00001330 CALL EXIT 00001340 70 FORMAT (79X,I1) 00001350 88 FORMAT (/T40,'PROGRAM RUN TIME WAS ',F9.2,' SECONDS') 00001360 122 FORMAT (' JOB HEADING 00001370 1 ') 00001380 151 FORMAT (7I10) 00001390 152 FORMAT ( 9I5) 00001400 200 FORMAT (3X,2I10) 00001410 302 FORMAT (T35,'OVER WRITING ITIE ARRAY'/) 00001420 381 FORMAT (I10,3F16.2,21X,I1) 00001430 398 FORMAT (/T20,'STRIP CONTAINS MORE THAN 500 POINTS') 00001440 404 FORMAT (2X,I8,3F15.3,24X,I1) 00001450 520 FORMAT (1H1) 00001460 END 00001470 SUBROUTINE IMCORE (*,*,*) 00001480 C 00001490 C REFINEMENT OF IMAGE COORDINATES READ ON COMPARATOR USC&GS00001500 C 00001510 REAL*8 A(50,16),B(50,9),C(21,9),E(18,3),R(56,2),FK,D,SQR,X,Y,F, 00001520 1T(151) 00001530 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00001540 COMMON/AREA3/A,B,C,E,R,FK,D,SQR,X,Y,T,N,NM1,NP1,NP2,NCOL,NROW, 00001550 1ITEST,IGO,LINE,IPLATE,IBLOCK,NBAD,L,ITRANS,LROW,JP2,NDEV,M,IRR,ME,00001560 2IP1,NMI,I,J,K,IEND 00001570 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00001580 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00001590 2JCON,INSTRU 00001600 248 ITRANS=1 00001610 C MEANING COMPARATOR READINGS OF FIDUCIALS 00001620 C 1 MEANS FIDUCIALS; 2 MEANS PHOTO IMAGE 00001630 READ (NREAD,520)IPLATE 00001640 C NO. OF PHOTO IN PROCESS 00001650 IF(IPLATE)238,605,238 00001660 605 RETURN 3 00001670 C IF 0 PROCESS IMAGES ON RT. SIDE OF LAST TRI00001680 238 LINE=31 00001690 I=31 00001700 204 READ (NREAD,502)(A(I,J),J=12,14),ITEST 00001710 C OBSERVED FIDUCIAL COORDINATES 00001720 IF(INSTRU.EQ.1)A(I,14)=-A(I,14) 00001730 241 GO TO (222,239),ITRANS 00001740 222 IF(A(31,12)-A(I,12))151,150,203 00001750 C SEQUENCE CHECK AND ERROR MESSAGE 00001760 203 WRITE (NPRNT,504) 00001770 RETURN 2 00001780 150 IF(ITEST)202,202,199 00001790 C FINDING LAST OF MULTIPLE OBSERVATIONS 00001800 202 I=I+1 00001810 GO TO 204 00001820 199 I=I+1 00001830 151 LROW=I-1 00001840 K=LROW-30 00001850 C NUMBER OF MULTIPLE OBSERVATIONS 00001860 FK=K 00001870 NDEV=1 00001880 C COUNT OF LARGE DEVIATIONS 00001890 219 IGO=1 00001900 A(LINE,2)=0.D0 00001910 A(LINE,3)=0.D0 00001920 DO 205 I=31,LROW 00001930 A(LINE,2)=A(LINE,2) +A(I,13) 00001940 205 A(LINE,3)=A(LINE,3) +A(I,14) 00001950 A(LINE,2)=A(LINE,2)/FK 00001960 A(LINE,3)=A(LINE,3)/FK 00001970 C MEAN X AND Y 00001980 251 DO 208 I=31,LROW 00001990 208 A(I,11)=DABS(A(I,13)-A(LINE,2)) 00002000 C X DEVIATION FROM MEAN 00002010 IF(2-NDEV)206,206,209 00002020 206 A(L,11)=0.D0 00002030 GO TO 209 00002040 207 DO 152 I=31,LROW 00002050 152 A(I,11)=DABS(A(I,14)-A(LINE,3)) 00002060 C Y DEVIATION FROM MEAN 00002070 IGO=2 00002080 IF(2-NDEV)153,153,209 00002090 153 A(M,11)=0.D0 00002100 209 L=31 00002110 J=31 00002120 C FINDING LARGEST DEVIATION 00002130 214 J=J+1 00002140 IF(J-LROW)155,155,154 00002150 155 IF(A(L,11)-A(J,11))215,214,214 00002160 215 L=J 00002170 GO TO 214 00002180 154 IF(A(L,11)-0.000025D0)217,217,216 00002190 C TEST FOR EXCESSIVE DEVIATION 00002200 216 IF(2-NDEV)203,203,218 00002210 C ONLY ONE IS PERMITTED 00002220 218 IOUT3=A(31,12) 00002230 WRITE (NPRNT,505)IPLATE,IOUT3 00002240 C DELETE THE OBSERVATIONS FOR EXCESSIVE IMAGE00002250 A(L,13)=0.D0 00002260 A(L,14)=0.D0 00002270 FK=FK-1.D0 00002280 NDEV=2 00002290 M=L 00002300 GO TO 219 00002310 217 GO TO (207,220),IGO 00002320 220 GO TO (198,240),ITRANS 00002330 C FIDUCIALS OR IMAGES 00002340 198 LINE=LINE+1 00002350 IF(N+30-LINE)221,197,197 00002360 C ACCOUNTING FOR ALL FIDUCIALS 00002370 197 A(31,12)=A(LROW+1,12) 00002380 A(31,13)=A(LROW+1,13) 00002390 A(31,14)=A(LROW+1,14) 00002400 C IMAGE NUMBER, X AND Y OF NEXT POINT 00002410 I=31 00002420 GO TO 241 00002430 C CORRECT X,Y FOR COMPARATOR ERROR BY USING MATRIX FOR FILM 00002440 C SHRINKAGE AND ORIENTATION OF PLATE 00002450 221 LROW=N+30 00002460 DO 223 I=31,LROW 00002470 A(I,NP1)=A(I,NP1)-A(I,2) 00002480 A(I,NP2)=A(I,NP2)-A(I,3) 00002490 A(I,1)=1.D0 00002500 A(I,4)=A(I,2)*A(I,3) 00002510 IF(LROW-34)223,223,224 00002520 224 A(I,5)=A(I,2)*A(I,2) 00002530 A(I,6)=A(I,3)*A(I,3) 00002540 A(I,7)=A(I,4)*A(I,3) 00002550 A(I,8)=A(I,5)*A(I,3) 00002560 C COEFFICIENTS FOR OBSERVATION EQUATIONS 00002570 223 CONTINUE 00002580 C FORMATION AND SOLUTION OF NORMAL EQUATIONS 00002590 DO 225 I=1,N 00002600 DO 225 J=I,NP2 00002610 A(I,J)=0.D0 00002620 DO 225 K=1,N 00002630 225 A(I,J)=A(I,J) + A(K+30,I)*A(K+30,J) 00002640 C FORMATION AND FORWARD SOLUTION 00002650 DO 227 I=1,N 00002660 702 SQR=DSQRT(A(I,I)) 00002670 DO 228 J=I,NP2 00002680 228 A(I,J)=A(I,J)/SQR 00002690 IF(I-N)229,230,230 00002700 229 IP1=I+1 00002710 DO 227 L=IP1,N 00002720 DO 227 J=L,NP2 00002730 227 A(L,J)=A(L,J)-A(I,L)*A(I,J) 00002740 C BACK SOLUTION 00002750 230 A(N,NP1)=A(N,NP1)/A(N,N) 00002760 A(N,NP2)=A(N,NP2)/A(N,N) 00002770 DO 231 I=1,NM1 00002780 NMI=N-I 00002790 NMIP1=NMI+1 00002800 DO 232 J=NMIP1,N 00002810 A(NMI,NP1)=A(NMI,NP1)-A(J,NP1)*A(NMI,J) 00002820 232 A(NMI,NP2)=A(NMI,NP2)-A(J,NP2)*A(NMI,J) 00002830 A(NMI,NP1)=A(NMI,NP1)/A(NMI,NMI) 00002840 231 A(NMI,NP2)=A(NMI,NP2)/A(NMI,NMI) 00002850 C PROCESSING NONFIDUCIAL IMAGE READINGS FROM COMPARATOR 00002860 237 ITRANS=2 00002870 NROW=0 00002880 C ROW IN B WHERE MEAN X AND Y STORED 00002890 GO TO 238 00002900 C TO READ IMAGE COORDINATES 00002910 239 IF(A(31,12)-A(I,12))151,150,151 00002920 C SEQUENCE CHECK 00002930 240 NROW=NROW+1 00002940 B(NROW,NCOL)=A(31,12) 00002950 B(NROW,NCOL+1)=A(31,2) 00002960 B(NROW,NCOL+2)=A(31,3) 00002970 C DATA TRANSFER 00002980 IF(ITEST)197,197,242 00002990 C TEST FOR LAST IMAGE 00003000 C CORRECTION FOR COMPARATOR ERROR, FILM SHRINKAGE AND ORIENTATION00003010 242 DO 243 I=1,NROW 00003020 DO 243 K=NP1,NP2 00003030 C K DESIGNATES X OR Y 00003040 IF(N-4)244,244,245 00003050 244 A(31,1)=A(1,K) +A(2,K)*B(I,NCOL+1) +A(3,K)*B(I,NCOL+2) +A(4,K)* 00003060 1B(I,NCOL+1)*B(I,NCOL+2) 00003070 GO TO 246 00003080 245 A(31,1)=A(1,K) +A(2,K)*B(I,NCOL+1) +A(3,K)*B(I,NCOL+2) +A(4,K)* 00003090 1B(I,NCOL+1)*B(I,NCOL+2) +A(5,K)*B(I,NCOL+1)**2 +A(6,K)*B(I,NCOL+2)00003100 2**2 +A(7,K)*B(I,NCOL+1)*B(I,NCOL+2)**2 +A(8,K)*B(I,NCOL+1)**2*B(I,00003110 3NCOL+2) 00003120 246 IF(K-NP1)261,261,262 00003130 261 X=A(31,1)+B(I,NCOL+1) 00003140 GO TO 243 00003150 262 B(I,NCOL+2)=A(31,1)+B(I,NCOL+2) 00003160 B(I,NCOL+1)= X 00003170 243 CONTINUE 00003180 C CORRECTION FOR FOCAL PLANE TILT 00003190 257 DO 247 I=1,NROW 00003200 X=A(42,3)*B(I,NCOL+1) +A(42,4)*B(I,NCOL+2) 00003210 Y=-A(42,4)*B(I,NCOL+1) +A(42,3)*B(I,NCOL+2) 00003220 B(I,NCOL+1)=X +A(42,5)*X*X 00003230 B(I,NCOL+2)=Y +A(42,5)*X*Y 00003240 X=A(42,3)* B(I,NCOL+1) - A(42,4)*B(I,NCOL+2) 00003250 Y=A(42,4)* B(I,NCOL+1) + A(42,3)*B(I,NCOL+2) 00003260 C CORRECTION FOR LENS DISTORTION AND ATMOSPHERIC REFRACTION 00003270 FK=DSQRT(X*X + Y*Y) 00003280 K=FK*1000.D0 00003290 F=K 00003300 D=1.D0+((1000.D0*FK-F)*(T(K+2)-T(K+1)) + T(K+1))*.000001D0 +A(42,100003310 1 ) + FK*FK*A(42,2) 00003320 B(I,NCOL+1) =X *D 00003330 247 B(I,NCOL+2) =Y *D 00003340 C FINAL REFINED COORDINATES 00003350 B(50,NCOL)=NROW 00003360 IF(A(41,5)-1.D0)259,146,259 00003370 146 K=NCOL+2 00003380 DO 148 I=1,NROW 00003390 148 WRITE (NPUNC,509)(B(I,J),J=NCOL,K) 00003400 C PROCESSING PLATES 2 AND 3 00003410 259 DO 157 I=1,N 00003420 A(I+30,NP1)=A(43,I) 00003430 157 A(I+30,NP1+1)=A(46,I) 00003440 C FIDUCIAL COORDINATES FOR NEW PHOTO 00003450 NCOL=NCOL+3 00003460 163 IF(NCOL-7)248,248,600 00003470 600 RETURN 1 00003480 502 FORMAT (8X,F7.0,F10.6,F10.6,44X,I1) 00003490 504 FORMAT (' MEASURED FIDUCIAL MARK CARDS NOT IN S5QUENCE OR 2 BAD 00003500 1READINGS ON POINT -EXIT') 00003510 505 FORMAT(' DEV FROM MEAN OVER 25 MICRONS. PLATE',I9,' POINT',I9) 509 FORMAT (2X,F9.0,1X,D16.10,2X,D16.10) 00003530 520 FORMAT ( I10) 00003540 END 00003550 C 00003560 SUBROUTINE PART1 00003570 C 00003580 C A N A L Y T I C A E R O T R I A N G U L A T I O N USC&GS00003590 C 00003600 REAL*8 A(50,16),B(50,9),C(21,9),E(18,3),R(56,2),FK,D,SQR,X,Y,F, 00003610 1T(151) 00003620 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00003630 COMMON/AREA3/A,B,C,E,R,FK,D,SQR,X,Y,T,N,NM1,NP1,NP2,NCOL,NROW, 00003640 1ITEST,IGO,LINE,IPLATE,IBLOCK,NBAD,L,ITRANS,LROW,JP2,NDEV,M,IRR,ME,00003650 2IP1,NMI,I,J,K,IEND 00003660 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00003670 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRN4,IEARTH,ICON, 00003680 2JCON,INSTRU 00003690 350 READ (NREAD,151)ISTRIP,ICAM 00003700 IF(ISTRIP.EQ.0)RETURN 00003710 400 IEXIT=0 00003720 J90=KLOCK(J) 00003730 IREC=0 00003740 ISTP=ISTP+1 00003750 ISTRP(ISTP)=ISTRIP 00003760 IF(ICAM.NE.1)GO TO 200 00003770 READ (NREAD,508)(A(42,J),J=1,2),A(41,1),A(41,5) 00003780 C REFRACTION; REJECTION LIMIT 00003790 READ (NREAD,501)N,(A(42,J),J=3,6) 00003800 C NO. OF MARKS; ASYM. DIST.; FOCAL LENGTH 00003810 200 NM1=N-1 00003820 NP1=N+1 00003830 NP2=N+2 00003840 IF(ICAM.NE.1)GO TO 201 00003850 READ (NREAD,503)(T(I),I=1,151) 00003860 C SYMMETRIC RADIAL LENS DISTORTION TABLE 00003870 READ (NREAD,502)(A(I+30,1),(A(I+30,J+1),J=N,NP1),ITEST,I=1,N) 00003880 C CORRECT FIDUCIAL COORDINATES 00003890 201 IBLOCK=1 00003900 C TEST NUMBER FOR INITIAL TRIPLET 00003910 NCOL=1 00003920 NBAD=0 00003930 C NUMBER OF RESIDUALS THAT EXCEED LIMIT 00003940 C STORE CORRECT COORDINATES OF FIDUCIALS 00003950 DO 156 I=1,N 00003960 A(43,I)=A(I+30,NP1) 00003970 156 A(46,I)=A(I+30,NP2) 00003980 603 CALL IMCORE (&260,&352,&355) 00003990 C BEGIN PHOTO ORIENTATION 00004000 602 DO 157 I=1,N 00004010 A(I+30,NP1)=A(43,I) 00004020 157 A(I+30,NP1+1)=A(46,I) 00004030 NCOL=NCOL+3 00004040 163 IF(NCOL-7)603,603,260 00004050 260 CALL PART2 (&352,&168) 00004060 C SOLUTION OF NORMAL EQUATIONS GAUSS 00004070 298 LINE=1 00004080 147 NROW=0 00004090 DO 310 M=1,12 00004100 LROW=LINE+2 00004110 DO 309 I=1,3 00004120 NROW=NROW+1 00004130 731 SQR=DSQRT(A(NROW,I)) 00004140 DO 308 J=I,16 00004150 308 A(NROW,J)=A(NROW,J)/SQR 00004160 IF(I-3)307,310,307 00004170 307 K=NROW+1 00004180 ME=I 00004190 DO 309 L=K,LROW 00004200 ME =ME+1 00004210 DO 309 J=ME,16 00004220 309 A(L,J)=A(L,J)-A(NROW,ME)*A(NROW,J) 00004230 310 LINE=LINE+3 00004240 K=3 00004250 DO 301 L=37,48 00004260 K=K+1 00004270 DO 301 I=1,36 00004280 DO 301 J=K,16 00004290 301 A(L,J)=A(L,J)-A(I,K)*A(I,J) 00004300 149 DO 302 I=37,48 00004310 K=I-33 00004320 IF(A(I,K))299,302,299 00004330 299 SQR=DSQRT(A(I,K)) 00004340 DO 303 J=K,16 00004350 303 A(I,J)=A(I,J)/SQR 00004360 IF(I-48)304,306,306 00004370 304 IP1=I+1 00004380 DO 335 L=IP1,48 00004390 K=K+1 00004400 DO 335 J=K,16 00004410 335 A(L,J)=A(L,J)-A(I,K)*A(I,J) 00004420 302 CONTINUE 00004430 C BACK SOLUTION 00004440 306 A(48,16)=A(48,16)/A(48,15) 00004450 DO 311 I=1,11 00004460 NMI=48-I 00004470 K=NMI-32 00004480 DO 305 J=K,15 00004490 IF(A(NMI,K-1))305,311,305 00004500 305 A(NMI,16)=A(NMI,16) -A(J+33,16)*A(NMI,J) 00004510 A(NMI,16)=A(NMI,16)/A(NMI,K-1) 00004520 311 CONTINUE 00004530 LROW=11 00004540 313 DO 312 I=1,3 00004550 L=I+LROW 00004560 K=4-I 00004570 NROW=48-L 00004580 DO 315 J=4,15 00004590 315 A(NROW,16)=A(NROW,16)-A(J+33,16)*A(NROW,J) 00004600 IF(I-2)312,316,317 00004610 316 A(NROW,16)=A(NROW,16)-A(NROW+1,16)*A(NROW,3) 00004620 GO TO 312 00004630 317 A(NROW,16)=A(NROW,16)-A(NROW+2,16)*A(NROW,3) 00004640 A(NROW,16)=A(NROW,16)-A(NROW+1,16)*A(NROW,2) 00004650 312 A(NROW,16)=A(NROW,16)/A(NROW,K) 00004660 C OBJECT COORDINATES 00004670 LROW=L 00004680 IF(LROW-44)313,313,314 00004690 314 DO 144 I=37,39 00004700 144 A(I,16)=A(I,16)*10.D0 00004710 C RESCALING ANGULAR PARAMETERS 00004720 DO 145 I=43,45 00004730 145 A(I,16)=A(I,16)*10.D0 00004740 C ADD LEAST SQUARES RESULTS TO C AND E ARRAYS 00004750 DO 321 J=4,6 00004760 C APPLICATION OF 48 COMPUTED ADJUSTMENTS 00004770 C(2,J)=C(2,J)+A(J+33,16) 00004780 321 C(1,J)=C(1,J)+A(J+36,16) 00004790 DO 322 J=7,9 00004800 C(2,J)=C(2,J)+A(J+36,16) 00004810 322 C(1,J)=C(1,J)+A(J+39,16) 00004820 M=1 00004830 DO 325 I=7,18 00004840 DO 325 J=1,3 00004850 E(I,J)=E(I,J)+A(M,16) 00004860 325 M=M+1 00004870 C TEST MAGNITUDE OF CORRECTIONS FOR ORIENTATION PARAMETERS 00004880 DO 318 I=37,39 00004890 IF(DABS(A(I,16))-.00001D0)318,318,319 00004900 318 CONTINUE 00004910 DO 320 I=43,45 00004920 IF(DABS(A(I,16))-.00001D0)320,320,319 00004930 320 CONTINUE 00004940 IEND=1 00004950 IRR=0 00004960 319 ITRANS=2 00004970 CALL PART3 (&352,&168) 00004980 GO TO 298 00004990 C X Y Z OUTPUT FOR PASS POINTS OBJECT COORDINATES 00005000 168 WRITE (NPRNT,512) 00005010 L=0 00005020 C SIGNAL FOR FINDING RELIABILITES 00005030 GO TO (323,332),IBLOCK 00005040 323 K=1 00005050 GO TO 160 00005060 332 K=7 00005070 160 M=12 00005080 DO 161 I=K,M 00005090 J=B(I,4)/100000.D0 00005100 FK=J*100000 00005110 FK=B(I,4)-FK 00005120 IFK=FK 00005130 IF(I-6)179,179,180 00005140 180 NCOL=I+L 00005150 X=DSQRT((R(NCOL,2)**2+R(NCOL+1,2)**2 +R(NCOL+12,2)**2 +R(NCOL+13, 00005160 12)**2 +R(NCOL+24,2)**2 +R(NCOL+25,2)**2)/6.D0) 00005170 C OBJECT RELIABILITY 00005180 IREC=IREC+1 00005190 WRITE (JS1)IFK,(E(I,J),J=1,3) 00005200 WRITE (NPRNT,529)IFK,(E(I,J),J=1,3),X 00005210 IF(IEARTH.EQ.9)WRITE (NPUNC,529)IFK,(E(I,J),J=1,3) 00005220 L=L+1 00005230 GO TO 161 00005240 179 WRITE (NPRNT,529)IFK,(E(I,J),J=1,3),R(I,2) 00005250 IREC=IREC+1 00005260 WRITE (JS1)IFK,(E(I,J),J=1,3) 00005270 IF(IEARTH.EQ.9)WRITE (NPUNC,529)IFK,(E(I,J),J=1,3) 00005280 C OBJECT NO. X Y Z AND SIGMA RELIABILITY 00005290 161 CONTINUE 00005300 R(55,2)=0.D0 00005310 DO 181 I=1,54 00005320 181 R(55,2)=R(55,2)+R(I,2)**2 00005330 R(55,2)=DSQRT(R(55,2)/42.D0) 00005340 WRITE (NPRNT,510)R(55,2) 00005350 C TRIPLET RELIABILITY RMS RESIDUAL PARALLAXES00005360 C COMPUTE X Y Z OF OTHER OBJECTS, IMAGES NOT USED IN ADJUSTMENT, 00005370 357 WRITE (NPRNT,515) 00005380 L=B(50,4) 00005390 M=B(50,1) 00005400 IF(L-18)333,360,333 00005410 333 IF(M-18)226,360,226 00005420 226 DO 347 K=19,M 00005430 J=B(K,1)/100000.D0 00005440 FK=J*100000 00005450 FK=B(K,1)-FK 00005460 DO 349 NCOL=19,L 00005470 J=B(NCOL,4)/100000.D0 00005480 D=J*100000 00005490 D=B(NCOL,4)-D 00005500 IF(FK-D)349,351,349 00005510 349 CONTINUE 00005520 C BEGIN INTERSECTION ROUTINE 00005530 GO TO 347 00005540 351 DO 348 I=1,3 00005550 348 C(16,I)=C(4,I)*B(K,2) + C(5,I)*B(K,3) + C(6,I)*(-A(42,6)) 00005560 DO 255 I=4,6 00005570 255 C(16,I)=C(4,I)*B(NCOL,5)+C(5,I)*B(NCOL,6)+C(6,I)*(-A(42,6)) 00005580 A(40,3)=C(16,1)/C(16,3) 00005590 A(40,4)=C(16,2)/C(16,3) 00005600 A(40,2)=C(16,4)-A(40,3)*C(16,6) 00005610 A(40,1)=C(1,6)*C(16,4)-C(1,4)*C(16,6)+C(1,1)*C(16,6)-C(1,3)*A(40,300005620 1)*C(16,6) 00005630 E(1,3)=A(40,1)/A(40,2) 00005640 E(1,2)=C(1,2)+A(40,4)*(E(1,3)-C(1,3)) 00005650 E(1,1)=C(1,1)+A(40,3)*(E(1,3)-C(1,3)) 00005660 E(2,2)=C(1,5)+(C(16,5)/C(16,6))*(E(1,3)-C(1,6)) 00005670 Y=E(2,2)-E(1,2) 00005680 C Y DIFFERENCE 00005690 E(1,2)=(E(1,2)+E(2,2))/2.D0 00005700 IFK=FK 00005710 IREC=IREC+1 00005720 WRITE (JS1)IFK,(E(1,J),J=1,3) 00005730 WRITE (NPRNT,529)IFK,(E(1,J),J=1,3),Y 00005740 IF(IEARTH.EQ.9)WRITE (NPUNC,529)IFK,(E(1,J),J=1,3) 00005750 C OBJECT NO. X Y Z AND Y^DIFF 00005760 347 CONTINUE 00005770 360 IF(IPLATE)358,359,358 00005780 C OUTPUT OF RESIDUAL PARALLAXES FOR PASS POINTS IN TRIPLET 00005790 358 WRITE (NPRNT,517) 00005800 WRITE (NPRNT,530)(R(I,1),R(I,2),I=1,54) 00005810 C IMAGE NO. AND RESIDUAL PARALLAX 00005820 C PROCEED TO NEW TRIPLET 00005830 WRITE (NPRNT,516)IPLATE 00005840 IBLOCK=2 00005850 DO 169 I=1,6 00005860 DO 169 J=1,3 00005870 C SHIFT DATA IN E ARRAY FOR NEW PHOTO 00005880 169 E(I,J)=E(I+6,J) 00005890 DO 162 I=1,50 00005900 DO 162 J=1,6 00005910 162 B(I,J)=B(I,J+3) 00005920 DO 166 I=1,15 00005930 DO 166 J=1,6 00005940 166 C(I,J)=C(I,J+3) 00005950 NCOL=4 00005960 GO TO 602 00005970 C TO REFINE COORDS. ON NEW PHOTO 00005980 352 I50=I50+1 00005990 IRCSP(ISTP)=IREC 00006000 ICAN(I50)=ISTRIP 00006010 IF(IEXIT.EQ.2)GO TO 350 00006020 405 READ (NREAD,518)I10,I80 00006030 IF(I10.NE.0)GO TO 405 00006040 353 IF(I80.NE.0)GO TO 405 00006050 GO TO 350 00006060 C OUTPUT DATA FOR RIGHT SIDE OF LAST TRIPLET IN STRIP 00006070 355 WRITE (NPRNT,512) 00006080 C THESE REQUIRE SPECIAL OUTPUT INSTRUCTIONS 00006090 DO 356 I=13,18 00006100 J=B(I,1)/100000.D0 00006110 FK=J*100000 00006120 FK=B(I,1)-FK 00006130 IFK=FK 00006140 IREC=IREC+1 00006150 WRITE (JS1)IFK,(E(I,J),J=1,3) 00006160 IF(IEARTH.EQ.9)WRITE (NPUNC,529)IFK,(E(I,J),J=1,3) 00006170 356 WRITE (NPRNT,529)IFK,(E(I,J),J=1,3),R(I+30,2) 00006180 GO TO 357 00006190 359 J91=KLOCK(J) 00006200 TIME=(J91-J90)/100. 00006210 WRITE (NPRNT,519)TIME 00006220 IF(IREC.GT.500)GO TO 399 00006230 IRCSP(ISTP)=IREC 00006240 GO TO 350 00006250 399 WRITE (NPRNT,398) 00006260 IEXIT=2 00006270 GO TO 352 00006280 151 FORMAT (7I10) 00006290 398 FORMAT (/T20,'STRIP CONTAINS MORE THAN 500 POINTS') 00006300 501 FORMAT (1X,I3,4(1X,D13.7)) 00006310 502 FORMAT (3X,F7.0,2F10.7,49X,I1) 00006320 503 FORMAT (15F5.1) 00006330 508 FORMAT (2F10.7,F10.6,1X,F1.0) 00006340 509 FORMAT (3X,F7.0,4(2X,D14.8)) 00006350 510 FORMAT (' RMS PARALLAX FOR TRIPLET IS',D16.8/) 512 FORMAT (/T6,'PASS PT',T20,'X',T36,'Y',T52,'Z',T64,'PRECISION'/) 515 FORMAT (/T7,'OBJECT',T20,'X',T36,'Y',T52,'Z',T64,'PRECISION'/) 516 FORMAT (/7X,'N E W T R I P L E T',T40,'LAST PLATE USED',I10/) 00006390 517 FORMAT(/T6,'IMAGE RESIDUAL X,Y PARALLAX'/) 00006400 518 FORMAT (I10,69X,I1) 00006410 519 FORMAT (/' ** STRIP ANALYTIC AEROTRIANGULATION COMPLETED ** 00006420 1RUN TIME WAS',F9.2,' SEC.'/) 00006430 521 FORMAT (40X,2I10,D14.7) 00006440 529 FORMAT (2X,I8,4(2X,D14.8)) 00006450 530 FORMAT (4(1X,F10.0,2X,D12.6)) 00006460 END 00006470 SUBROUTINE PART2(*,*) 00006480 C 00006490 C P R O C E S S I N G T R I P L E T S T R I P USC&GS 00006500 C 00006510 REAL*8 A(50,16),B(50,9),C(21,9),E(18,3),R(56,2),FK,D,SQR,X,Y,F, 00006520 1T(151) 00006530 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00006540 COMMON/AREA3/A,B,C,E,R,FK,D,SQR,X,Y,T,N,NM1,NP1,NP2,NCOL,NROW, 00006550 1ITEST,IGO,LINE,IPLATE,IBLOCK,NBAD,L,ITRANS,LROW,JP2,NDEV,M,IRR,ME,00006560 2IP1,NMI,I,J,K,IEND 00006570 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00006580 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00006590 2JCON,INSTRU 00006600 260 IF(NBAD)211,210,211 00006610 211 DO 329 I=1,6 00006620 L=I/2 00006630 J=B(I,7)/100000.D0 00006640 FK=J*100000 00006650 FK=B(I,7)-FK 00006660 IF(R(55,1)-FK)213,212,213 00006670 213 IF(R(56,1)-FK)329,212,329 00006680 212 IF(I-2*L)121,121,122 00006690 121 DO 123 J=7,9 00006700 123 B(I,J)=B(I-1,J) 00006710 GO TO 329 00006720 122 DO 124 J=7,9 00006730 124 B(I,J)=B(I+1,J) 00006740 C PUTING COMPANION DATA IN PLACE OF UNACCEPTA00006750 329 CONTINUE 00006760 210 DO 130 J=4,7,3 00006770 DO 130 K=1,12 00006780 I=B(K,J)/100000.D0 00006790 C SORTING TEST 00006800 FK=I*100000 00006810 FK=B(K,J)-FK 00006820 I=B(K+6,J-3)/100000.D0 00006830 D=I*100000 00006840 D=B(K+6,J-3)-D 00006850 IF(FK-D)131,130,131 00006860 131 IOUT1=B(K,J) 00006870 IOUT2=B(K+6,J-3) 00006880 WRITE (NPRNT,506)IOUT1,IOUT2 00006890 400 RETURN 1 00006900 130 CONTINUE 00006910 NBAD=0 00006920 72 IEND=0 00006930 ITRANS=1 00006940 GO TO (164,324),IBLOCK 00006950 C THREE PHOTO ORIENTATION 00006960 C C ARRAY ROWS 1 THROUGH 15 00006970 164 DO 265 J=1,9 00006980 C ONLY FOR FIRST TRIPLET 00006990 C(1,J)=0.D0 00007000 265 C(2,J)=0.D0 00007010 ENTRY PART3 (*,*) 00007020 324 DO 132 I=1,33 00007030 DO 132 J=1,16 00007040 132 A(I,J)=0.D0 00007050 DO 133 I=34,48 00007060 K=I-33 00007070 DO 133 J=K,16 00007080 133 A(I,J)=0.D0 00007090 LROW=0 00007100 IF(IBLOCK-1)196,196,269 00007110 C SKIP ON FIRST SOLUTION 00007120 196 GO TO (195,269),ITRANS 00007130 195 J=1 00007140 IGO=1 00007150 C FOR FIRST PHOTO OF FIRST TRIPLET 00007160 264 JP2=J+2 00007170 C ELEMENTS OF ROTATION MATRIX AND THEIR PART.00007180 DO 266 K=J,JP2 00007190 266 C(3,K)=DSQRT(1.D0-C(2,K)*C(2,K)) 00007200 C(4,J)=C(3,J+1)*C(3,JP2) 00007210 C(5,J)=-C(3,J+1)*C(2,JP2) 00007220 C(6,J)=C(2,J+1) 00007230 C(10,J)=-C(2,J+1)*C(3,JP2) 00007240 C(11,J)=C(2,J+1)*C(2,JP2) 00007250 C(12,J)=C(3,J+1) 00007260 C(10,J+1)=C(4,J)*C(2,J) 00007270 C(11,J+1)=C(5,J)*C(2,J) 00007280 C(12,J+1)=C(2,J)*C(2,J+1) 00007290 C(10,JP2)=-C(4,J)*C(3,J) 00007300 C(11,JP2)=-C(5,J)*C(3,J) 00007310 C(12,JP2)=-C(3,J)*C(2,J+1) 00007320 C(4,J+1)=C(3,J)*C(2,JP2) + C(12,J+1)*C(3,JP2) 00007330 C(5,J+1)=C(3,J)*C(3,JP2) - C(12,J+1)*C(2,JP2) 00007340 C(6,J+1)=-C(2,J)*C(3,J+1) 00007350 C(4,JP2)=C(2,J)*C(2,JP2) + C(10,J)*C(3,J) 00007360 C(5,JP2)=C(2,J)*C(3,JP2) + C(11,J)*C(3,J) 00007370 C(6,JP2)=C(3,J)*C(3,J+1) 00007380 DO 267 I=7,9 00007390 C(I,J)=0.D0 00007400 C(I,J+1)=-C(I-3,JP2) 00007410 267 C(I,JP2)=C(I-3,J+1) 00007420 DO 268 I=J,JP2 00007430 C(13,I)=C(5,I) 00007440 C(14,I)=-C(4,I) 00007450 268 C(15,I)=0.D0 00007460 271 GO TO (269,263,272),IGO 00007470 269 J=4 00007480 IGO=2 00007490 C FOR SECOND PHOTO OF TRIPLET 00007500 IF(IBLOCK-1)170,170,264 00007510 170 GO TO (331,264),ITRANS 00007520 C RESET DURING INITIAL TRIPLET 00007530 331 C(1,J)=B(17,2) 00007540 GO TO 264 00007550 263 J=7 00007560 IGO=3 00007570 C FOR THIRD PHOTO 00007580 GO TO (330,264),ITRANS 00007590 330 C(1,J)=B(17,5)+C(1,J-3) 00007600 GO TO 264 00007610 C C ARRAY ROWS 16 THROUGH 21 AND E ARRAY 00007620 C THE AB PRODUCTS AND MODEL COORD. LIST 00007630 272 IGO=1 00007640 NDEV=37 00007650 LINE=42 00007660 J=4 00007670 JP2=J+2 00007680 K=1 00007690 C COUNT OF IMAGE PROCESSED 00007700 328 GO TO (94,289),IBLOCK 00007710 94 M=4 00007720 DO 334 L=18,20 00007730 M=M+3 00007740 DO 334 I=J,JP2 00007750 334 C(L,I-3)=C(M,I)*B(K,5) +C(M+1,I)*B(K,6) +C(M+2,I)*(-A(42,6)) 00007760 275 DO 274 I=4,6 00007770 C BEGIN INTERSECTION ROUTINE 00007780 274 C(16,I)=C(4,I)*B(K,5) + C(5,I)*B(K,6) + C(6,I)*(-A(42,6)) 00007790 IF(IGO-4)129,129,142 00007800 129 DO 277 I=1,3 00007810 277 C(16,I)=C(4,I)*B(K+6,2) + C(5,I)*B(K+6,3) + C(6,I)*(-A(42,6)) 00007820 GO TO 278 00007830 142 DO 143 I=7,9 00007840 143 C(16,I)=C(4,I)*B(K-6,8) + C(5,I)*B(K-6,9) +C(6,I)*(-A(42,6)) 00007850 J=7 00007860 JP2=J+2 00007870 278 A(40,3)=C(16,J-3)/C(16,J-1) 00007880 A(40,4)=C(16,J-2)/C(16,J-1) 00007890 A(40,2)=C(16,J)-A(40,3)*C(16,JP2) 00007900 A(40,1)=C(1,JP2)*C(16,J)-C(1,J)*C(16,JP2)+C(1,J-3)*C(16,JP2)- 00007910 1C(1,J-1)*A(40,3)*C(16,JP2) 00007920 E(K,3)=A(40,1)/A(40,2) 00007930 C Z Y AND X 00007940 E(K,2)=C(1,J-2)+A(40,4)*(E(K,3)-C(1,J-1)) 00007950 E(K,1)=C(1,J-3)+A(40,3)*(E(K,3)-C(1,J-1)) 00007960 233 J=4 00007970 300 JP2=J+2 00007980 289 C(17,J)=E(K,1)-C(1,J) 00007990 C(17,J+1)=E(K,2)-C(1,J+1) 00008000 C(17,JP2)=E(K,3)-C(1,JP2) 00008010 C B VECTOR 00008020 M=4 00008030 DO 279 L=18,21 00008040 DO 279 I=J,JP2 00008050 C(L,I)=C(M,J)*C(17,J) + C(M,J+1)*C(17,J+1) + C(M,JP2)*C(17,JP2) 00008060 279 M=M+1 00008070 C AB VECTOR PRODUCT 00008080 254 GO TO (273,273,290,292,273,292),IGO 00008090 C COMPUTE P CONSTANTS COEFF. OF OBS. EQUATIONS 00008100 273 NCOL=1 00008110 DO 280 L=18,21 00008120 A(44,NCOL)=(B(K,J+1)*C(L,JP2)-(-A(42,6))*C(L,J))/C(18,JP2) 00008130 A(45,NCOL)=(B(K,JP2)*C(L,JP2)-(-A(42,6))*C(L,J+1))/C(18,JP2) 00008140 280 NCOL=NCOL+1 00008150 DO 281 I=J,JP2 00008160 A(44,NCOL)=(B(K,J+1)*C(6,I) - (-A(42,6))*C(4,I))/C(18,JP2) 00008170 A(45,NCOL)=(B(K,JP2)*C(6,I) - (-A(42,6))*C(5,I))/C(18,JP2) 00008180 281 NCOL=NCOL+1 00008190 IF(IEND)344,344,171 00008200 C IF 1 COMPUTE RESIDUAL PARALLAX 00008210 171 IRR=IRR+1 00008220 R(IRR,1)=B(K,J) 00008230 R(IRR,2)=A(45,1) 00008240 C RESIDUAL PARALLAXES 00008250 GO TO (137,172,172,172,295,297),IGO 00008260 172 IRR=IRR+1 00008270 R(IRR,1)=B(K,J) 00008280 R(IRR,2)=A(44,1) 00008290 GO TO (137,137,138,293),IGO 00008300 344 GO TO (270,92),IBLOCK 00008310 270 IF(IGO-1)89,89,92 00008320 C COMPUTE CONSTANTS FOR FIRST TRIPLET 00008330 89 DO 326 I=44,45 00008340 326 A(I,8)=(1.D0/A(40,2))*(A(I,5)*A(40,3)+A(I,6)*A(40,4)+A(I,7)) 00008350 M=2 00008360 DO 327 L=18,20 00008370 A(41,M)=C(17,J)*C(L,J-1) - C(17,JP2)*C(L,J-3) 00008380 327 M=M+1 00008390 92 A(49,16)=-A(44,1) 00008400 A(50,16)=-A(45,1) 00008410 C CLEAR AREA TO STORE 2 OBS. EQUATIONS 00008420 DO 282 I=49,50 00008430 DO 282 M=1,15 00008440 282 A(I,M)=0.D0 00008450 GO TO (134,134,287,139,134,139),IGO 00008460 C COMUPTE CONTRIBUTION TO NORMAL EQUATIONS 00008470 134 DO 283 M=4,6 00008480 A(49,M)=A(44,M-2)*10.D0 00008490 283 A(50,M)=A(45,M-2)*10.D0 00008500 C SCALE FACTOR OF ANGULAR PARAMETERS COEFF 00008510 C PHOTO 2 00008520 DO 284 M=7,9 00008530 A(49,M)=-A(44,M-2) 00008540 284 A(50,M)=-A(45,M-2) 00008550 GO TO (285,286),IBLOCK 00008560 285 A(49,7)=0.D0 00008570 A(50,7)=0.D0 00008580 C FOR FIRST TRIPLET ONLY 00008590 286 IF(IGO-1)88,88,287 00008600 287 DO 288 M=1,3 00008610 A(49,M)=A(44,M+4) 00008620 288 A(50,M)=A(45,M+4) 00008630 C TRANSFER DX DY DZ COEFF. FOR COMPRESSED FOR00008640 GO TO 127 00008650 139 DO 140 M=10,12 00008660 C PHOTO 3 00008670 A(49,M)=A(44,M-8)*10.D0 00008680 140 A(50,M)=A(45,M-8)*10.D0 00008690 DO 141 M=13,15 00008700 A(49,M)=-A(44,M-8) 00008710 141 A(50,M)=-A(45,M-8) 00008720 GO TO 287 00008730 127 DO 128 I=1,3 00008740 LROW=LROW+1 00008750 DO 128 M=I,16 00008760 A(LROW,M)=A(LROW,M) + A(49,I)*A(49,M) 00008770 128 A(LROW,M)=A(LROW,M) + A(50,I)*A(50,M) 00008780 C CONTRIBUTION TO NORMAL EQUATIONS 00008790 256 IF(IGO-3)125,138,125 00008800 88 GO TO (165,125),IBLOCK 00008810 165 DO 93 I=49,50 00008820 A(I,9)=A(I,9)+A(I-5,8)*C(16,J) 00008830 DO 93 M=4,6 00008840 93 A(I,M)=A(I,M) +A(I-5,8)*A(41,M-2)*10.D0 00008850 125 DO 126 I=NDEV,LINE 00008860 L=I-33 00008870 DO 126 M=L,16 00008880 A(I,M)=A(I,M) + A(49,L)*A(49,M) 00008890 126 A(I,M)=A(I,M) + A(50,L)*A(50,M) 00008900 C CONTRIBUTION TO NORMAL EQUATIONS 00008910 GO TO (137,137,137,293,295,297),IGO 00008920 C INDEX VALUES FOR DIFFERENT CASES 00008930 137 K=K+1 00008940 IF(K-6)328,328,276 00008950 276 IF(K-12)135,135,136 00008960 135 IGO=2 00008970 GO TO (275,289),ITRANS 00008980 136 LROW=0 00008990 IGO=3 00009000 K=7 00009010 J=1 00009020 GO TO 300 00009030 290 K=K+6 00009040 GO TO 273 00009050 138 K=K-5 00009060 IF(K-12)289,289,291 00009070 291 IGO=4 00009080 LROW=0 00009090 K=7 00009100 J=7 00009110 NDEV=43 00009120 LINE=48 00009130 GO TO 300 00009140 292 K=K-6 00009150 GO TO 273 00009160 293 K=K+7 00009170 IF(K-12)289,289,294 00009180 294 IGO=5 00009190 NDEV=37 00009200 LINE=42 00009210 GO TO (275,233),ITRANS 00009220 295 K=K+1 00009230 GO TO (345,346),ITRANS 00009240 345 IF(K-18)275,275,296 00009250 346 IF(K-18)289,289,296 00009260 296 IGO=6 00009270 LROW=18 00009280 J=7 00009290 K=13 00009300 NDEV=43 00009310 LINE=48 00009320 GO TO 300 00009330 297 K=K+7 00009340 IF(K-18)289,289,194 00009350 194 IF(IEND)298,298,173 00009360 C TO LEAST SQUARES OR RESIDUALS 00009370 173 IF(IBLOCK-1)192,192,234 00009380 C ANALYZE FOR RESIDUALS IN PROPER CASES 00009390 192 I=1 00009400 GO TO 235 00009410 234 I=7 00009420 C FINDING LARGEST ABSOLUTE RESIDUAL 00009430 235 M=I 00009440 177 M=M+1 00009450 IF(M-54)174,174,175 00009460 174 IF(DABS(R(I,2)) - DABS(R(M,2)))176,177,177 00009470 176 I=M 00009480 GO TO 177 00009490 175 IF(DABS(R(I,2))-A(41,1))168,168,178 00009500 C 168 IF OK 00009510 168 RETURN 2 00009520 178 IF(I-6)339,339,340 00009530 C IDENTIFICATION OF BAD RESIDUAL 00009540 340 IF(I-42)341,341,339 00009550 341 L=I/2 00009560 IF(I-2*L)342,342,339 00009570 342 IOUT1=R(I,1) 00009580 WRITE (NPRNT,514)IOUT1,R(I,2) 00009590 GO TO 343 00009600 339 IOUT2=R(I,1) 00009610 WRITE (NPRNT,507)IOUT2,R(I,2) 00009620 343 NBAD=NBAD+1 00009630 GO TO (182,183,184),NBAD 00009640 183 J=R(I,1)/100000.D0 00009650 FK=J*100000 00009660 FK=R(I,1)-FK 00009670 IF(DABS(R(55,1) - FK) - 1.D0)185,185,182 00009680 C BOTH IMAGES BAD ? 00009690 185 WRITE (NPRNT,511) 00009700 GO TO 400 00009710 184 WRITE (NPRNT,513) 00009720 GO TO 400 00009730 182 J=R(I,1)/100000.D0 00009740 FK=J*100000 00009750 R(55,1)=R(I,1)-FK 00009760 GO TO (253,236),NBAD 00009770 253 R(56,1)=R(55,1) 00009780 236 DO 187 I=1,18 00009790 C SEARCH FOR BAD DATA AND REPLACE WITH COMPAN00009800 L=I/2 00009810 M=1 00009820 186 NCOL=M+2 00009830 J=B(I,M)/100000.D0 00009840 FK=J*100000 00009850 FK=B(I,M)-FK 00009860 IF(M-4)338,188,189 00009870 338 IF(R(55,1)-FK)190,193,190 00009880 190 M=4 00009890 GO TO 186 00009900 188 IF(R(55,1)-FK)249,193,249 00009910 249 M=7 00009920 GO TO 186 00009930 189 IF(R(55,1)-FK)187,193,187 00009940 193 IF(I-2*L)158,158,159 00009950 158 DO 336 J=M,NCOL 00009960 336 B(I,J)=B(I-1,J) 00009970 GO TO 187 00009980 159 DO 337 J=M,NCOL 00009990 337 B(I,J)=B(I+1,J) 00010000 187 CONTINUE 00010010 GO TO 72 00010020 298 RETURN 00010030 506 FORMAT (' SORTING ERROR PREVENTS MATCHING',I9,' EXIT') 507 FORMAT (/' YPARALLAX FOR'I9,' IS',D16.8,' AND EXCEEDS LIMIT') 511 FORMAT (' TWO BAD RESIDUAL PARALLAXES IN SAME LOCATION. EXIT') 00010060 513 FORMAT (' THREE BAD RESIDUAL PARALLAXES IN SAME TRIPLET. EXIT') 00010070 514 FORMAT (/' XPARALLAX FOR',I9,' IS',D16.8,' AND EXEEDS LIMIT') 521 FORMAT (40X,2I10,D14.7) 00010090 END 00010100 C 00010110 C 00010120 SUBROUTINE ASTRIP 00010130 C 00010140 C P O L Y N O M I A L A D J U S T M E N T O F S T R I P 00010150 C 00010160 IMPLICIT REAL*8(A-H,O-Z) 00010170 REAL*8 C(40,14),A(40,8),EN(7,8),VAC(7),HAC(7),SUM1,SUM2,SQR,SCALE,00010180 1ADX(100),ADY(100),ADZ(100),XM(500),YM(500),ZM(500) 00010190 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00010200 1,IFKA(500) 00010210 COMMON/AREA3/A,C,EN,VAC,HAC,ADX,ADY,ADZ,XM,YM,ZM,SUM1,SUM2,SQR 00010220 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00010230 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00010240 2JCON,INSTRU 00010250 DEFINE FILE 13(800,25,U,KK3) 00010260 ISTOP=0 00010270 IST8=0 00010280 DO 412 I=1,100 00010290 412 ADX(I)=0.D0 00010300 N7=0 00010310 REWIND JS2 00010320 C INSTRUCTION PARAMETERS 00010330 1 READ (NREAD,100)IST9,NHC,NVC,IDHA,IDVA,IREP,PCONST,XAX,YAY 00010340 IF(IST9.EQ.0)GO TO 171 00010350 IF(I50.EQ.0)GO TO 337 00010360 DO 336 I51=1,I50 00010370 IF(IST9.EQ.ICAN(I51))GO TO 170 00010380 336 CONTINUE 00010390 337 I30=KLOCK(J) 00010400 ITNC=NHC+NVC 00010410 IF(IST8.EQ.IST9)GO TO 411 00010420 DO 410 K=1,ITNC 00010430 DO 410 J=5,7 00010440 410 C(K,J)=0.D0 00010450 411 IFIND=0 00010460 I40=0 00010470 REWIND JS1 00010480 DO 306 I=1,ISTP 00010490 IF(IST9.EQ.ISTRP(I))GO TO 312 00010500 306 CONTINUE 00010510 WRITE (NPRNT,307)IST9 00010520 GO TO 170 00010530 312 IREC=IRCSP(I) 00010540 I=I-1 00010550 IF(I.GT.0)GO TO 310 00010560 316 DO 311 I2=1,IREC 00010570 311 READ (JS1,ERR=313,END=313)IFKA(I2),XM(I2),YM(I2),ZM(I2) 00010580 GO TO 313 00010590 310 ITOT=0 00010600 DO 314 I3=1,I 00010610 314 ITOT=ITOT+IRCSP(I3) 00010620 IF(ITOT.EQ.0)GO TO 316 00010630 DO 315 I4=1,ITOT 00010640 315 READ (JS1,ERR=170,END=170) 00010650 GO TO 316 00010660 313 IF(IS1A0.EQ.1)GO TO 348 00010670 IF(IS1A0.GT.2)GO TO 348 00010680 IF(IEARTH.EQ.1)GO TO 348 00010690 WRITE (NPRNT,103) 00010700 348 WRITE (NPRNT,130)IST9 00010710 WRITE (NPRNT,341)NHC,NVC,IS1A0,IDHA,IDVA,IREP,PCONST,XAX,YAY 00010720 L22=ITNC 00010730 IF(IREP.EQ.0)GO TO 377 00010740 DO 378 I=1,ITNC 00010750 378 C(I,8)=C(I,1) 00010760 377 IF(ITNC.GT.40)GO TO 349 00010770 IF(IS1A0.EQ.1)GO TO 2 00010780 IF(IS1A0.GT.2)GO TO 2 00010790 C 1 OR 2 INSTRUMENTAL, 0 ANALYTIC 00010800 DO 322 I=1,IREC 00010810 C MODEL COORDINATES IN MILLIMETERS 00010820 XM(I)=XM(I)*1000.D0 00010830 YM(I)=YM(I)*1000.D0 00010840 C REVERSE Y DIRECTION OF INDEPENDENT MODEL STRIP COORDINATES 00010850 IF(IS1A0.EQ.2)YM(I)=5000.D0-YM(I) 00010860 322 ZM(I)=ZM(I)*1000.D0 00010870 2 READ (NREAD,151)LINE1,IPN,LINE2,IPN2 00010880 C TERMINAL PHOTO CENTERS 00010890 WRITE (NPRNT,151)LINE1,IPN,LINE2,IPN2 00010900 ISTEP=1 00010910 C LOCATE MODEL COORDINATES 00010920 321 DO 317 I=1,IREC 00010930 IF(IPN.EQ.IFKA(I))GO TO 318 00010940 317 CONTINUE 00010950 IF(IFIND.EQ.1)GO TO 329 00010960 WRITE (NPRNT,319)IPN 00010970 GO TO 170 00010980 318 GO TO (323,320,326,331),ISTEP 00010990 323 X1=XM(I) 00011000 Y1=YM(I) 00011010 ISTEP=2 00011020 IPN=IPN2 00011030 GO TO 321 00011040 320 X2=XM(I) 00011050 Y2=YM(I) 00011060 4 IF(LINE2-LINE1)5,5,6 00011070 C CENTERS IN PROPER ORDER 00011080 5 WRITE (NPRNT,102) 00011090 GO TO 170 00011100 C CONSTANTS FOR TRANSFORMATION OF MODEL COORDINATES TO AXIS OF 00011110 C FLIGHT SYSTEM 00011120 6 XDIF=X1-X2 00011130 YDIF=Y1-Y2 00011140 DIST=DSQRT(XDIF*XDIF+YDIF*YDIF) 00011150 A1=-XDIF/DIST 00011160 B1=YDIF/DIST 00011170 C1=(-DIST/2.D0)-(A1*X1) +(B1*Y1) 00011180 D1=(-B1*X1)-(A1*Y1) 00011190 C READING CONTROL STATIONS 00011200 7 IF(IREP.EQ.1)GO TO 370 00011210 READ (NREAD,324)(C(I,1),I=1,ITNC) 00011220 370 ISTEP=3 00011230 II=0 00011240 C LOCATE MODEL COORDINATES FOR STRIP CONTROL POINTS 00011250 327 II=II+1 00011260 IF(II.GT.ITNC)GO TO 9 00011270 IPN=C(II,1) 00011280 GO TO 321 00011290 326 C(II,2)=XM(I) 00011300 C(II,3)=YM(I) 00011310 C(II,4)=ZM(I) 00011320 GO TO 327 00011330 C INPUT ALL CONTROL VALUES USED IN THE ADJUST00011340 9 WRITE (NPRNT,105)(C(I,1),(C(I,J),J=2,4),I=1,ITNC) 00011350 IF(IREP.EQ.1)GO TO 371 00011360 C LOCATE GROUND COORDINATES FOR STRIP CONTROL POINTS 00011370 DO 10 I=1,ITNC 00011380 C(I,8)=C(I,1) 00011390 IT8=C(I,1) 00011400 DO 11 K=1,JCON 00011410 IF(IT8.EQ.IGRND(K))GO TO 400 00011420 11 CONTINUE 00011430 GO TO 10 00011440 400 READ (KD3'K)(C(I,J),J=5,7) 00011450 10 CONTINUE 00011460 371 IFIND =1 00011470 C ARE COORDINATES DEFINED FOR CONTROL 00011480 DO 160 I=1,NHC 00011490 IF(C(I,5).EQ.0.D0)GO TO 351 00011500 GO TO 160 00011510 351 IT8=C(I,1) 00011520 DO 158 N8=1,ICON 00011530 IF(IT8.EQ.ITIE(N8))GO TO 159 00011540 158 CONTINUE 00011550 GO TO 162 00011560 159 C(I,5)=ADX(N8) 00011570 C(I,6)=ADY(N8) 00011580 160 CONTINUE 00011590 II=NHC+1 00011600 DO 157 I=II,ITNC 00011610 IF(C(I,7).EQ.0.D0)GO TO 352 00011620 GO TO 157 00011630 352 IT8=C(I,1) 00011640 DO 161 N8=1,ICON 00011650 IF(IT8.EQ.ITIE(N8))GO TO 353 00011660 161 CONTINUE 00011670 GO TO 162 00011680 353 C(I,7)=ADZ(N8) 00011690 157 CONTINUE 00011700 WRITE (NPRNT,105)(C(I,8),(C(I,J),J=5,7),I=1,ITNC) 00011710 GO TO 76 00011720 162 WRITE (NPRNT,163)C(I,8) 00011730 GO TO 170 00011740 C TRANSFORMATION OF MODEL COORDINATES TO AXIS OF FLIGHT SYSTEM 00011750 76 IF(XAX.EQ.0.D0.AND.YAY.EQ.0.D0)GO TO 380 00011760 DO 379 I=1,L22 00011770 C(I,5)=C(I,5)-XAX 00011780 379 C(I,6)=C(I,6)-YAY 00011790 380 DO 12 I=1,ITNC 00011800 SUM1=A1*C(I,2)-B1*C(I,3)+C1 00011810 SUM2=B1*C(I,2)+A1*C(I,3)+D1 00011820 C(I,2)=SUM1 00011830 12 C(I,3)=SUM2 00011840 C VERTICAL ADJUSTMENT FOR SLOPE FACTORS 00011850 XDIF=C(1,2)-C(NHC,2) 00011860 YDIF=C(1,3)-C(NHC,3) 00011870 DGX=C(1,5)-C(NHC,5) 00011880 DGY=C(1,6)-C(NHC,6) 00011890 SCALE=DSQRT((DGX*DGX+DGY*DGY)/(XDIF*XDIF+YDIF*YDIF)) 00011900 WRITE (NPRNT,123)SCALE 00011910 SUM2=0.D0 00011920 FTNC=ITNC 00011930 DO 13 I=1,ITNC 00011940 13 SUM2=SUM2+C(I,4) 00011950 AVIZ=SUM2/FTNC 00011960 C AVERAGE VERT. INSTRUMENT (MODEL) ELEV. 00011970 NHCP1=NHC+1 00011980 FNVC=NVC 00011990 SUM1=0.D0 00012000 DO 14 I=NHCP1,ITNC 00012010 14 SUM1=SUM1+C(I,7) 00012020 SUM1=SUM1/FNVC 00012030 EINDEX=AVIZ-(SUM1/SCALE) 00012040 IGO=1 00012050 31 DO 15 I=NHCP1,ITNC 00012060 C(I,8)=C(I,7)/SCALE+EINDEX 00012070 15 C(I,9)=C(I,8)-C(I,4) 00012080 LINE1=NHCP1 00012090 LINE2=ITNC 00012100 ITRANS=1 00012110 45 I=1 00012120 C 1ST, 2ND, 3RD, DEG. VERT. ADJUSTMENT? 00012130 IF(IDVA-2)16,17,18 00012140 C 3RD DEGREE VERTICAL OBSERVATION EQUATIONS 00012150 18 DO 19 J=LINE1,LINE2 00012160 A(I,3)=C(J,2) 00012170 A(I,6)=C(J,3) 00012180 A(I,5)=A(I,3)*A(I,6) 00012190 A(I,2)=A(I,3)*A(I,3) 00012200 A(I,4)=A(I,2)*A(I,6) 00012210 A(I,1)=A(I,3)*A(I,2) 00012220 A(I,7)=1.D0 00012230 A(I,8)=C(J,9) 00012240 19 I=I+1 00012250 N=7 00012260 GO TO (77,74),ITRANS 00012270 C 2ND DEGREE VERTICAL OBSERVATION EQUATIONS 00012280 17 DO 20 J=LINE1,LINE2 00012290 A(I,2)=C(J,2) 00012300 A(I,4)=C(J,3) 00012310 A(I,3)=A(I,2)*A(I,4) 00012320 A(I,1)=A(I,2)*A(I,2) 00012330 A(I,5)=1.D0 00012340 A(I,6)=C(J,9) 00012350 20 I=I+1 00012360 N=5 00012370 GO TO (77,74),ITRANS 00012380 C 1ST DEGREE VERTICAL OBSERVATION EQUATIONS 00012390 16 DO 21 J=LINE1,LINE2 00012400 A(I,1)=C(J,2) 00012410 A(I,3)=C(J,3) 00012420 A(I,2)=A(I,1)*A(I,3) 00012430 A(I,4)=1.D0 00012440 A(I,5)=C(J,9) 00012450 21 I=I+1 00012460 N=4 00012470 GO TO (77,74),ITRANS 00012480 77 CONTINUE 00012490 C FORMATION AND SOLUTION OF NORMAL EQUATIONS 00012500 22 LINE=NVC 00012510 40 NP1=N+1 00012520 NM1=N-1 00012530 DO 23 I=1,N 00012540 DO 23 J=I,NP1 00012550 EN(I,J)=0.D0 00012560 DO 23 K=1,LINE 00012570 23 EN(I,J)=EN(I,J)+A(K,I)*A(K,J) 00012580 C FORWARD SOLUTION 00012590 DO 99 I=1,N 00012600 SQR=DSQRT(EN(I,I)) 00012610 DO 98 J=I,NP1 00012620 98 EN(I,J)=EN(I,J)/SQR 00012630 IF(I-N)97,96,96 00012640 97 IP1=I+1 00012650 DO 99 L=IP1,N 00012660 DO 99 J=L,NP1 00012670 99 EN(L,J)=EN(L,J)-EN(I,L)*EN(I,J) 00012680 96 EN(N,NP1)=EN(N,NP1)/EN(N,N) 00012690 C BACK SOLUTION 00012700 DO 91 I=1,NM1 00012710 NMI=N-I 00012720 NMIP1=NMI+1 00012730 DO 90 J=NMIP1,N 00012740 90 EN(NMI,NP1)=EN(NMI,NP1)-EN(J,NP1)*EN(NMI,J) 00012750 91 EN(NMI,NP1)=EN(NMI,NP1)/EN(NMI,NMI) 00012760 GO TO (24,32,41),IGO 00012770 C SLOPE FACTORS APPLIED TO LIST OF HOTIZONTAL & VERT. CONTROL POINTS00012780 24 IF(IDVA-2)25,26,27 00012790 27 DO 28 I=1,ITNC 00012800 XS=3.D0*EN(1,8)*C(I,2)*C(I,2)+2.D0*EN(2,8)*C(I,2)+EN(3,8) 00012810 YS=EN(4,8)*C(I,2)*C(I,2)+EN(5,8)*C(I,2)+EN(6,8) 00012820 ZS=DSQRT(1.D0+XS*XS+YS*YS) 00012830 C(I,2)=C(I,2)-(C(I,4)-AVIZ)*XS 00012840 C(I,3)=C(I,3)-(C(I,4)-AVIZ)*YS 00012850 28 C(I,4)=C(I,4)*ZS 00012860 GO TO 201 00012870 26 DO 29 I=1,ITNC 00012880 XS=2.D0*EN(1,6)*C(I,2) + EN(2,6) 00012890 YS=EN(3,6)*C(I,2) + EN(4,6) 00012900 ZS=DSQRT(1.D0+XS*XS+YS*YS) 00012910 C(I,2)=C(I,2)-(C(I,4)-AVIZ)*XS 00012920 C(I,3)=C(I,3)-(C(I,4)-AVIZ)*YS 00012930 29 C(I,4)=C(I,4)*ZS 00012940 GO TO 201 00012950 25 DO 200 I=1,ITNC 00012960 XS=EN(1,5) 00012970 YS=EN(2,5)*C(I,2) + EN(3,5) 00012980 ZS=DSQRT(1.D0+XS*XS+YS*YS) 00012990 C(I,2)=C(I,2)-(C(I,4)-AVIZ)*XS 00013000 C(I,3)=C(I,3)-(C(I,4)-AVIZ)*YS 00013010 200 C(I,4)=C(I,4)*ZS 00013020 201 CONTINUE 00013030 C CONSTANTS FOR TRANSFORMATION OF GROUND PLANE COORDINATES TO AXIS 00013040 C OF FLIGHT SYSTEM 00013050 80 XDIF=C(1,2)-C(NHC,2) 00013060 YDIF=C(1,3)-C(NHC,3) 00013070 DIST=XDIF*XDIF+YDIF*YDIF 00013080 A2=(DGX * XDIF + DGY * YDIF)/DIST 00013090 B2=(XDIF * DGY - YDIF * DGX)/DIST 00013100 C2=C(1,5)-A2*C(1,2)+B2*C(1,3) 00013110 D2=C(1,6)-B2*C(1,2)-A2*C(1,3) 00013120 DIST=1.D0/(A2*A2 + B2*B2) 00013130 A21=A2*DIST 00013140 B21=B2*DIST 00013150 C21=(A2*C2 + B2*D2)*DIST 00013160 D21=(A2*D2 - B2*C2)*DIST 00013170 C INVERSE 00013180 C CONTROL TRANSFORMED INTO AXIS OF FLIGHT SYSTEM 00013190 DO 30 I=1,NHC 00013200 C(I,8)=A21*C(I,5) + B21*C(I,6) - C21 00013210 C(I,9)=-B21*C(I,5) + A21*C(I,6) - D21 00013220 C(I,10)=C(I,8)-C(I,2) 00013230 30 C(I,11)=C(I,9)-C(I,3) 00013240 SCALE=DSQRT((DGX*DGX + DGY*DGY)/(XDIF*XDIF + YDIF*YDIF)) 00013250 81 CONTINUE 00013260 82 IGO=2 00013270 GO TO 31 00013280 C VERTICAL ADJUSTMENT COEFFICIENTS Z POLYNOMIAL 00013290 32 DO 33 I=1,N 00013300 33 VAC(I)=EN(I,NP1) 00013310 DO 47 I=1,NVC 00013320 SUM1=0.D0 00013330 DO 48 J=1,N 00013340 48 SUM1=SUM1+VAC(J)*A(I,J) 00013350 K=I+NHC 00013360 47 C(K,10)=C(K,9)-SUM1 00013370 C HORIZONTAL ADJUSTMENT COEFFICIENTS X AND Y POLYNOMIALS 00013380 LINE=NHC*2 00013390 LINE1=1 00013400 LINE2=NHC 00013410 ITRANS=1 00013420 IGO =3 00013430 52 I=1 00013440 IF(IDHA-2)34,35,36 00013450 C 3RD DEGREE HORIZONTAL OBSERVATION EQUATIONS 00013460 36 N=7 00013470 DO 37 J=LINE1,LINE2 00013480 A(I,3)=C(J,2) 00013490 A(I,5)=-C(J,3) 00013500 A(I,2)=A(I,3)*A(I,3) 00013510 A(I,1)=A(I,3)*A(I,2) 00013520 A(I,4)=2.D0*A(I,3)*A(I,5) 00013530 A(I,6)=1.D0 00013540 A(I,7)=0.D0 00013550 A(I,8)=C(J,10) 00013560 I=I+1 00013570 A(I,3)=-A(I-1,5) 00013580 A(I,5)=A(I-1,3) 00013590 A(I,4)=A(I-1,2) 00013600 A(I,2)=-A(I-1,4) 00013610 A(I,6)=0.D0 00013620 A(I,7)=1.D0 00013630 A(I,1)=3.D0*A(I,4)*A(I,3) 00013640 A(I,8)=C(J,11) 00013650 37 I=I+1 00013660 GO TO (83,54),ITRANS 00013670 C 2ND DEGREE HORIZONTAL OBSERVATION EQUATIONS 00013680 35 N=6 00013690 DO 38 J=LINE1,LINE2 00013700 A(I,2)=C(J,2) 00013710 A(I,4)=-C(J,3) 00013720 A(I,1)=A(I,2)*A(I,2) 00013730 A(I,3)=2.D0*A(I,2)*A(I,4) 00013740 A(I,5)=1.D0 00013750 A(I,6)=0.D0 00013760 A(I,7)=C(J,10) 00013770 I=I+1 00013780 A(I,2)=-A(I-1,4) 00013790 A(I,4)=A(I-1,2) 00013800 A(I,3)=A(I,4)*A(I,4) 00013810 A(I,1)=-A(I-1,3) 00013820 A(I,5)=0.D0 00013830 A(I,6)=1.D0 00013840 A(I,7)=C(J,11) 00013850 38 I=I+1 00013860 GO TO (83,54),ITRANS 00013870 C 1ST DEGREE HORIZONTAL OBSERVATION EQUATIONS 00013880 34 N=4 00013890 DO 39 J=LINE1,LINE2 00013900 A(I,1)=C(J,2) 00013910 A(I,2)=-C(J,3) 00013920 A(I,3)=1.D0 00013930 A(I,4)=0.D0 00013940 A(I,5)=C(J,10) 00013950 I=I+1 00013960 A(I,1)=-A(I-1,2) 00013970 A(I,2)=A(I-1,1) 00013980 A(I,3)=0.D0 00013990 A(I,4)=1.D0 00014000 A(I,5)=C(J,11) 00014010 39 I=I+1 00014020 GO TO (83,54),ITRANS 00014030 83 CONTINUE 00014040 GO TO 40 00014050 C EQUATION SOLUTION 00014060 41 DO 42 I=1,N 00014070 42 HAC(I)=EN(I,NP1) 00014080 CXBOW=HAC(N-1) 00014090 CYBOW=HAC(N) 00014100 C DATA FOR PRINT OUT OF HORIZONTAL CONTROL USED FOR ADJUSTMENT 00014110 WRITE (NPRNT,107) 00014120 WRITE (NPRNT,108) 00014130 C JOB TITLE COLUMN HEADINGS 00014140 53 L=1 00014150 DO 43 I=1,NHC 00014160 SUM1=0.D0 00014170 SUM2=0.D0 00014180 DO 44 J=1,N 00014190 SUM1=SUM1+HAC(J)*A(L,J) 00014200 44 SUM2=SUM2+HAC(J)*A(L+1,J) 00014210 L=L+2 00014220 C(I,12)=C(I,10)-SUM1 00014230 43 C(I,13)=C(I,11)-SUM2 00014240 ITRANS=2 00014250 GO TO 45 00014260 74 DO 49 I=1,NHC 00014270 SUM1=0.D0 00014280 DO 46 J=1,N 00014290 46 SUM1=SUM1+VAC(J)*A(I,J) 00014300 C(I,14)=C(I,4)+SUM1 00014310 C(I,14)=SCALE*(C(I,14)-EINDEX) 00014320 C COMPUTED GROUND ELEVATION 00014330 IT8=C(I,1) 00014340 DO 401 K=1,JCON 00014350 IF(IT8.EQ.IGRND(K))GO TO 402 00014360 401 CONTINUE 00014370 GO TO 49 00014380 C IF VERTICAL CONTROL MISSING USE CALCULATED VALUE 00014390 402 READ (KD3'K)PP,PP2,PP3 00014400 IF(PP3.EQ.0.D0)PP3=C(I,14) 00014410 WRITE (KD3'K)PP,PP2,PP3 00014420 49 CONTINUE 00014430 WRITE (NPRNT,109)(C(I,1),(C(I,J),J=10,14),I=1,NHC) 00014440 SUM1=0.D0 00014450 SUM2=0.D0 00014460 DO 50 I=1,NHC 00014470 SUM1=SUM1+C(I,12)*C(I,12) 00014480 50 SUM2=SUM2+C(I,13)*C(I,13) 00014490 DIST=NHC-1 00014500 STDX=DSQRT(SUM1/DIST) 00014510 STDY=DSQRT(SUM2/DIST) 00014520 STDXY=DSQRT(STDX*STDX + STDY*STDY) 00014530 WRITE (NPRNT,110)STDX,STDY,STDXY 00014540 WRITE (NPRNT,111)CXBOW,CYBOW 00014550 85 CONTINUE 00014560 C DATA FOR PRINT OUT OF VERTICAL CONTROL USED IN ADJUSTMENT 00014570 86 SUM1=0.D0 00014580 DO 51 I=NHCP1,ITNC 00014590 C(I,11)=0.D0 00014600 51 SUM1=SUM1+C(I,10)*C(I,10) 00014610 DIST=NVC-1 00014620 STDZ=DSQRT(SUM1/DIST) 00014630 LINE1=NHCP1 00014640 LINE2=ITNC 00014650 GO TO 52 00014660 54 L=1 00014670 DO 55 I=NHCP1,ITNC 00014680 SUM1=0.D0 00014690 SUM2=0.D0 00014700 DO 73 J=1,N 00014710 SUM1=SUM1+HAC(J)*A(L,J) 00014720 73 SUM2=SUM2+HAC(J)*A(L+1,J) 00014730 L=L+2 00014740 C(I,11)=SUM1+C(I,2) 00014750 C(I,12)=SUM2+C(I,3) 00014760 SUM1=A2*C(I,11)-B2*C(I,12)+C2 00014770 SUM2=B2*C(I,11)+A2*C(I,12)+D2 00014780 C(I,11)=SUM1 00014790 C(I,12)=SUM2 00014800 IT8=C(I,1) 00014810 DO 403 K=1,JCON 00014820 IF(IT8.EQ.IGRND(K))GO TO 404 00014830 403 CONTINUE 00014840 GO TO 55 00014850 C IF HORIZONTAL CONTROL MISSING USE CALCULATED VALUE 00014860 404 READ (KD3'K)PP,PP2,PP3 00014870 IF(PP.EQ.0.D0)PP=C(I,11) 00014880 IF(PP2.EQ.0.D0)PP2=C(I,12) 00014890 WRITE (KD3'K)PP,PP2,PP3 00014900 55 CONTINUE 00014910 WRITE (NPRNT,112) 00014920 WRITE (NPRNT,113) 00014930 WRITE (NPRNT,123)(C(I,1),(C(I,J),J=9,12),I=NHCP1,ITNC) 00014940 WRITE (NPRNT,114)STDZ 00014950 IF(IEARTH.EQ.0.AND.ITO.EQ.13)GO TO 71 00014960 WRITE (NPRNT,117) 00014970 WRITE (NPRNT,118) 00014980 C VERT. AND OTHER HORIZONTAL CONTROL USED 00014990 C DATA FOR PUNCH OUT OF OTHER CONTROL AND BRIDGE POINTS 00015000 56 READ (NREAD,325)(C(I,2),I=1,7),ITYPE,ITEST 00015010 II=0 00015020 ISTEP=4 00015030 ITNC=0 00015040 DO 328 J=1,7 00015050 IF(C(J,2).EQ.0.D0)GO TO 329 00015060 328 ITNC=ITNC+1 00015070 329 IF(ITNC.EQ.0.AND.ITYPE.EQ.2)GO TO 59 00015080 IF(ITNC.EQ.0)GO TO 330 00015090 ITNC=ITNC-1 00015100 II=II+1 00015110 IPN=C(II,2) 00015120 GO TO 321 00015130 331 X=XM(I) 00015140 Y=YM(I) 00015150 Z=ZM(I) 00015160 GO TO 57 00015170 57 IF(ITYPE-1)62,58,59 00015180 58 WRITE (NPRNT,119) 00015190 WRITE (NPRNT,118) 00015200 ITYPE=0 00015210 GO TO 62 00015220 59 WRITE (NPRNT,120) 00015230 WRITE (NPRNT,118) 00015240 ITYPE=0 00015250 I40=1 00015260 I41=0 00015270 335 I41=I41+1 00015280 IF(I41.GT.IREC)GO TO 71 00015290 IPN=IFKA(I41) 00015300 X=XM(I41) 00015310 Y=YM(I41) 00015320 Z=ZM(I41) 00015330 62 SUM1=A1*X - B1*Y + C1 00015340 SUM2=B1*X + A1*Y + D1 00015350 C AXIS OF FLIGHT TRANSFORMATION 00015360 IF(IDVA-2)210,211,212 00015370 212 XS=3.D0*VAC(1)*SUM1*SUM1 + 2.D0*VAC(2)*SUM1 + VAC(3) 00015380 YS=VAC(4)*SUM1*SUM1 + VAC(5)*SUM1 + VAC(6) 00015390 GO TO 213 00015400 211 XS=2.D0*VAC(1)*SUM1 + VAC(2) 00015410 YS=VAC(3)*SUM1 + VAC(4) 00015420 GO TO 213 00015430 210 XS=VAC(1) 00015440 YS=VAC(2)*SUM1 + VAC(3) 00015450 213 ZS=DSQRT(1.D0+XS*XS+YS*YS) 00015460 X=SUM1-(Z-AVIZ)*XS 00015470 Y=SUM2-(Z-AVIZ)*YS 00015480 Z=Z*ZS 00015490 IF(IDHA-2)63,64,65 00015500 63 SUM1=HAC(1)*X + HAC(2)*(-Y) + HAC(3) 00015510 SUM2=HAC(1)*Y + HAC(2)*X + HAC(4) 00015520 GO TO 66 00015530 C POLYNOMIAL TRANSFORMATION 00015540 64 SUM1=HAC(1)*X*X + HAC(2)*X + HAC(3)*(-2.D0)*X*Y + HAC(4)*(-Y) + 00015550 1HAC(5) 00015560 SUM2=HAC(1)*2.D0*X*Y + HAC(2)*Y + HAC(3)*X*X + HAC(4)*X + HAC(6) 00015570 C 2ND DEGREE 00015580 GO TO 66 00015590 65 SUM1=HAC(1)*X*X*X + HAC(2)*X*X + HAC(3)*X + HAC(4)*(-2.D0)*X*Y + 00015600 1HAC(5)*(-Y) + HAC(6) 00015610 SUM2=HAC(1)*3.D0*X*X*Y + HAC(2)*2.D0*X*Y + HAC(3)*Y + HAC(4)*X*X +00015620 1HAC(5)*X + HAC(7) 00015630 66 TRUEX=SUM1+X 00015640 TRUEY=SUM2+Y 00015650 IF(IDVA-2)67,68,69 00015660 69 SUM1=VAC(1)*X*X*X +VAC(2)*X*X +VAC(3)*X + VAC(4)*X*X*Y + VAC(5) 00015670 1*X*Y + VAC(6)*Y + VAC(7) 00015680 GO TO 70 00015690 68 SUM1=VAC(1)*X*X + VAC(2)*X + VAC(3)*X*Y + VAC(4)*Y + VAC(5) 00015700 GO TO 70 00015710 67 SUM1=VAC(1)*X + VAC(2)*X*Y + VAC(3)*Y +VAC(4) 00015720 70 TRUEZ=SUM1+Z 00015730 X=A2*TRUEX-B2*TRUEY+C2 00015740 Y=B2*TRUEX+A2*TRUEY+D2 00015750 Z=SCALE*(TRUEZ-EINDEX) 00015760 PLOTX=PCONST * X 00015770 PLOTY=PCONST * Y 00015780 C FORM DIFFERENCE BETWEEN CONTROL & COMPUTED VALUES 00015790 DO 381 M=1,JCON 00015800 IF(IPN.EQ.IGRND(M))GO TO 382 00015810 381 CONTINUE 00015820 GO TO 383 00015830 382 READ (KD3'M)PP,PP2,PP3 00015840 IF(PP.NE.0.D0)PP=PP-X 00015850 IF(PP2.NE.0.D0)PP2=PP2-Y 00015860 IF(PP3.NE.0.D0)PP3=PP3-Z 00015870 WRITE (NPRNT,384)IPN,PP,PP2,PP3 00015880 383 IF(ICON.EQ.0)GO TO 155 00015890 DO 154 I=1,ICON 00015900 IF(IPN.EQ.ITIE(I))GO TO 156 00015910 154 CONTINUE 00015920 GO TO 155 00015930 C GENERATE CONTROL VALUES FOR STRIP TIE POINTS 00015940 156 IF(ADX(I).NE.0.D0)GO TO 155 00015950 ADX(I)=X 00015960 ADY(I)=Y 00015970 ADZ(I)=Z 00015980 155 IF(ICOORD.EQ.1.AND.I40.EQ.1)GO TO 342 00015990 343 WRITE (NPRNT,139)IPN,X,Y,Z,PLOTX,PLOTY 00016000 IF(IEARTH.EQ.1)WRITE (NPUNC,407)IPN,X,Y,Z 00016010 IF(I40.EQ.1)GO TO 335 00016020 IF(ITEST.EQ.0)GO TO 329 00016030 330 IF(ITEST)56,56,71 00016040 342 N7=N7+1 00016050 WRITE (JS2)IPN,X,Y,Z 00016060 GO TO 343 00016070 71 WRITE (NPRNT,121) 00016080 IST8=IST9 00016090 I31=KLOCK(J) 00016100 TIME= (I31-I30)/100. 00016110 WRITE (NPRNT,88)TIME 00016120 GO TO 1 00016130 C REJECT STRIP DATA 00016140 170 IF(IREP.EQ.1)GO TO 180 00016150 ISTOP=ISTOP+1 00016160 ISTOP2=ISTP-ISTOP 00016170 IF(ISTOP2.LE.0)ICOORD=0 00016180 180 READ (NREAD,169)ITEST 00016190 IF(ITEST.EQ.1)GO TO 1 00016200 GO TO 180 00016210 171 IF(ICON.EQ.0)GO TO 304 00016220 IF(IEARTH.EQ.0.AND.ITO.EQ.13)GO TO 304 00016230 WRITE (NPRNT,89) 00016240 DO 164 I=1,ICON 00016250 164 WRITE (NPRNT,139)ITIE(I),ADX(I),ADY(I),ADZ(I) 00016260 304 END FILE JS2 00016270 RETURN 00016280 349 WRITE (NPRNT,350)ITNC 00016290 GO TO 170 00016300 88 FORMAT (T40,'STRIP RUN TIME WAS',F9.2,' SECONDS') 00016310 89 FORMAT (/T7,'TIE POINTS X,Y,Z IN GROUND COORDINATE SYSTEM'/) 00016320 100 FORMAT (6I5,F10.9,2F10.3) 00016330 101 FORMAT (I2,1X,I7,3(2X,D14.8)) 00016340 102 FORMAT (/' INITIAL - TERMINAL CARDS REVERSED EXIT') 00016350 103 FORMAT (1H1) 00016360 105 FORMAT (3X,F7.0,3F16.3) 00016370 106 FORMAT (I2,1X,I7,3F16.2) 00016380 107 FORMAT (/5X,'HORIZONTAL CONTROL USED FOR ADJUSTMENT'/) 00016390 108 FORMAT (/5X,'POINT',T19,'CX',T30,'CY',T41,'RX',T52,'RY',T60,'COMP 00016400 1Z GROUND'/) 00016410 109 FORMAT (1X,F11.0,2X,D9.3,2X,D9.3,2X,D9.3,2X,D9.3,2X,D16.10) 00016420 110 FORMAT (/T5,'STDX= ',G14.8,' STDY= ',G14.8,' STDXY= ',G14.8,/) 00016430 111 FORMAT (T5,'BOW ERROR AT AXIS OF FLIGHT ORIGIN CX= ',G14.8,' CY= 00016440 1',G14.8/) 00016450 112 FORMAT (//T5,'VERTICAL CONTROL USED FOR ADJUSTMENT'/) 00016460 113 FORMAT (/3X,' POINT',T19,'CZ',T30,'RZ',T38,'COMP X GROUND COM00016470 1P Y GROUND'/) 00016480 114 FORMAT (/T5,'STDZ= ',G14.8) 00016490 115 FORMAT (3X,I7,3(2X,D14.8),20X,2I1) 00016500 116 FORMAT (3X,I7,3F16.2,20X,2I1) 00016510 117 FORMAT (//T5,'OTHER HORIZONTAL CONTROL'/) 00016520 118 FORMAT (/T6,'POINT COMP X GROUND COMP Y GROUND COMP Z G00016530 1ROUND PLOT X PLOT Y'/) 00016540 119 FORMAT (//T5,'OTHER VERTICAL CONTROL'/) 00016550 120 FORMAT (//T5,'ALL POINTS IN STRIP'/) 00016560 121 FORMAT(/' ** STRIP ADJUSTMENT COMPLETED **') 00016570 123 FORMAT (1X,F11.0,2X,D9.3,2X,D9.3,2X,D16.10,2X,D16.10) 00016580 125 FORMAT (F2.0,F8.0,3D16.8 ) 00016590 126 FORMAT (F2.0,F8.0,3F16.2 ) 00016600 130 FORMAT (///T20,'COMMENCE ADJUSTMENT OF STRIP NO. ',I10//) 00016610 139 FORMAT (1X,I 9,5(2X,D16.10)) 00016620 151 FORMAT (7I10) 00016630 163 FORMAT (3X,'CONTROL NOT DEFINED FOR ',F10.0) 00016640 169 FORMAT (79X,I1) 00016650 195 FORMAT (1X,F4.0,F12.0,3D16.8) 00016660 196 FORMAT (1X,F4.0,F12.0,3F16.2) 00016670 307 FORMAT (T20,'STRIP NO.',I12,' NOT IN LIST,REJECT STRIP') 00016680 319 FORMAT (T20,'POINT NO.',I10,' NOT IN LIST, REJECT STRIP') 00016690 324 FORMAT (7F10.0) 00016700 325 FORMAT (7F10.0,8X,2I1) 00016710 341 FORMAT (6I5,F15.8,2F15.3) 00016720 350 FORMAT (T10,'CONTROL FOR',I4,' POINTS EXCEEDS LIMIT') 00016730 384 FORMAT (T10,I10,3F15.3) 00016740 407 FORMAT (1X,I9,3(3X,D18.12)) 00016750 END 00016760 C 00016770 C 00016780 SUBROUTINE TCOORD 00016790 C 00016800 C USCGS T R A N S F O R M A T I O N O F C O O R D I N A T E00016810 C 00016820 IMPLICIT REAL*8(A-H,O-Z) 00016830 REAL*8G(4,500),CONV(500),T(11) 00016840 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00016850 COMMON/AREA3/G,T,PHDG,PHMIN,PHSEC,PLGDG,PLGMN,PLSEC,AMTFT,FTMT, 00016860 1SECRD,ESQ,CONV,A,SNONS 00016870 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00016880 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00016890 2JCON,INSTRU 00016900 C SET ELLIPSOID PARAMETERS AND CONVERSION CONSTANTS CLARKE 1866 00016910 DEFINE FILE 13(800,25,U,KK3) 00016920 DATA PI/3.1415926535898D0/ 00016930 REWIND JS2 00016940 ESQ=.6768658D-2 00016950 FTMT=3.28083989501312D0 00016960 AMTFT=.3048D0 00016970 A=6378206.4D0 00016980 SECRD=648000.D0/PI 00016990 SNONS=1.D0/SECRD 00017000 IRUN2=0 00017010 IRUN=0 00017020 IF(ICOORD.NE.2)WRITE (NPRNT,998) 00017030 I9=KLOCK(J) 00017040 IF(IEARTH.EQ.0.AND.ITO.EQ.13)GO TO 220 00017050 IF(ICOORD.EQ.1)N5=N7 00017060 IF(ICOORD.EQ.2)N7=N5 00017070 C INPUT VALUES FROM SEQUENTIAL DATA SET 2 00017080 1 IF(N7-500)30,30,31 00017090 30 NN=N7 00017100 IRUN=1 00017110 33 DO 32 I=1,NN 00017120 READ (JS2,ERR=34,END=34)IPN,(G(J,I),J=2,4) 00017130 32 G(1,I)=IPN 00017140 GO TO 34 00017150 31 NN=500 00017160 N7=N7-NN 00017170 GO TO 33 00017180 220 IRUN=1 00017190 DO 221 I=1,JCON 00017200 221 READ (KD3'I)(G(J,I),J=2,4) 00017210 NN=JCON 00017220 34 IF(IRUN2.EQ.1)GO TO 53 00017230 C WRITE ORIGIN AND DESTINATION COORDINATE SYSTEMS 00017240 WRITE (NPRNT,1021)IFROM,ITO 00017250 WRITE (NPRNT,1020)ESQ,A 00017260 53 L=0 00017270 C IF SECANT PLANE SYSTEM IS INVOLVED, READ AND WRITE THE LATITUDE AN00017280 C LONGITUDE OF THE SYSTEM ORIGIN 00017290 IF(IFROM.NE.13.AND.ITO.NE.13) GO TO 2 00017300 IF(IRUN2.EQ.1)GO TO 2 00017310 READ (NREAD,1000)PHDG,PHMIN,PHSEC,PLGDG,PLGMN,PLSEC 00017320 IOUT1=PHDG 00017330 IOUT2=PHMIN 00017340 IOUT3=PLGDG 00017350 IOUT4=PLGMN 00017360 WRITE (NPRNT,1022)IOUT1,IOUT2,PHSEC,IOUT3,IOUT4,PLSEC 00017370 2 IF(IFROM.EQ.1) GO TO 5 00017380 IF(IFROM.EQ.2) GO TO 15 00017390 IF(ITO.EQ.1) GO TO 5 00017400 GO TO 15 00017410 5 IF(IRUN2.EQ.1)GO TO 54 00017420 READ (NREAD,1003)(T(M),M=1,6) 00017430 WRITE (NPRNT,1023) 00017440 WRITE (NPRNT,1024)(T(M),M=1,6) 00017450 54 IF(L.EQ.1) GO TO 20 00017460 GO TO 15 00017470 10 IF(IRUN2.EQ.1)GO TO 36 00017480 READ (NREAD,1003)(T(M),M=1,11) 00017490 WRITE (NPRNT,1023) 00017500 WRITE (NPRNT,1024)(T(M),M=1,11) 00017510 36 IF(L.EQ.1)GO TO 20 00017520 15 IF(IFROM.NE.12)GO TO 130 00017530 19 I=ITO 00017540 IDIR=1 00017550 C CALL APPROPRIATE TRANSFORMATION SUBROUTINE 00017560 20 IF(IABS(I).GT.100)GO TO 80 00017570 GO TO (50,60,70,80,80,80,80,80,80,80,80,1,90,100),I 00017580 50 CONTINUE 00017590 GO TO 110 00017600 60 CALL STERGP 00017610 GO TO 110 00017620 70 CONTINUE 00017630 GO TO 110 00017640 80 CALL UTMGP 00017650 GO TO 110 00017660 90 CALL GPSP 00017670 GO TO 110 00017680 100 CALL GEOCGP 00017690 110 IF(IDIR.EQ.0)GO TO 160 00017700 IF(ITO.EQ.13)GO TO 120 00017710 IF(ITO.EQ.2)GO TO 200 00017720 C PRINT PLANE COORDINATE OUTPUT 00017730 IF(IRUN2.EQ.1)GO TO 37 00017740 WRITE (NPRNT,1027) 00017750 37 DO 141 J=1,NN 00017760 IOUT=G(1,J) 00017770 C WRITE (NPUNC,1036)IOUT,(G(M,J),M=2,4) 00017780 141 WRITE (NPRNT,1036)IOUT,(G(M,J),M=2,4) 00017790 GO TO 900 00017800 C PRINT SECANT PLANE OUTPUT 00017810 120 IF(IRUN2.EQ.1)GO TO 38 00017820 IF(IEARTH.EQ.0.AND.ITO.EQ.13)GO TO 222 00017830 WRITE (NPRNT,1028) 00017840 38 DO 142 J=1,NN 00017850 IOUT=G(1,J) 00017860 WRITE (NPUNC,1036)IOUT,(G(M,J),M=2,4) 00017870 142 WRITE (NPRNT,1036)IOUT,(G(M,J),M=2,4) 00017880 GO TO 900 00017890 200 IF(IRUN2.EQ.1)GO TO 39 00017900 WRITE (NPRNT,1040) 00017910 39 DO 201 J=1,NN 00017920 IOUT=G(1,J) 00017930 CC=CONV(J) 00017940 CALL DEG(CC,IC,MC,SC) 00017950 C WRITE (NPUNC,1041)IOUT,(G(M,J),M=2,4) 00017960 201 WRITE (NPRNT,1041)IOUT,(G(M,J),M=2,4),IC,MC,SC 00017970 GO TO 900 00017980 130 IDIR=0 00017990 IF(IFROM.NE.13)GO TO 150 00018000 C TRANSFORM SECANT PLANE COORDINATES TO GEOGRAPHIC POSITIONS AND ELE00018010 CALL GPSP 00018020 IF(ITO.NE.12)GO TO 19 00018030 135 IF(IRUN2.EQ.1)GO TO 136 00018040 WRITE (NPRNT,1026) 00018050 C CONVERT SECONDS OF GEOGRAPHIC POSITIONS TO DEGREES,MINUTES,SECONDS00018060 136 DO 140 J=1,NN 00018070 LAD=G(2,J)/3600.D0 00018080 XS=LAD 00018090 G(2,J)=G(2,J)-3600.D0*XS 00018100 LAM=G(2,J)/60.D0 00018110 XS=LAM 00018120 XS=G(2,J)-60.D0*XS 00018130 LOD=G(3,J)/3600.D0 00018140 YS=LOD 00018150 G(3,J)=G(3,J)-3600.D0*YS 00018160 LOM=G(3,J)/60.D0 00018170 YS=LOM 00018180 YS=G(3,J)-60.D0*YS 00018190 LOM=IABS(LOM) 00018200 YS=DABS(YS) 00018210 LAM=IABS(LAM) 00018220 XS=DABS(XS) 00018230 IOUT=G(1,J) 00018240 C PRINT GEOGRAPHIC POSITION AND ELEVATION 00018250 IF(IFROM.EQ.2)GO TO 205 00018260 C WRITE (NPUNC,1018)IOUT,LAD,LAM,XS,LOD,LOM,YS,G(4,J) 00018270 WRITE (NPRNT,1018)IOUT,LAD,LAM,XS,LOD,LOM,YS,G(4,J) 00018280 GO TO 140 00018290 205 CC=DABS(CONV(J)) 00018300 CALL DEG(CC,IC,MC,SC) 00018310 WRITE (NPRNT,1018)IOUT,LAD,LAM,XS,LOD,LOM,YS,G(4,J),IC,MC,SC 00018320 140 CONTINUE 00018330 GO TO 900 00018340 150 I=IFROM 00018350 GO TO 20 00018360 160 IF(ITO.EQ.12)GO TO 135 00018370 IF(ITO.NE.13)GO TO 170 00018380 IDIR=1 00018390 CALL GPSP 00018400 GO TO 120 00018410 170 IDIR=1 00018420 I=ITO 00018430 L=1 00018440 K=IFROM+ITO 00018450 IF(K.GT.3)GO TO 20 00018460 GO TO (5,10),ITO 00018470 222 DO 223 I=1,JCON 00018480 223 WRITE (KD3'I)(G(J,I),J=2,4) 00018490 900 IRUN2=1 00018500 IF(IRUN.EQ.0)GO TO 1 00018510 I=KLOCK(J) 00018520 TIME=(I-I9)/100. 00018530 WRITE (NPRNT,1009)TIME 00018540 RETURN 00018550 998 FORMAT (1H1) 00018560 C SECANT PLANE SYSTEM ORIGIN LATITUDE AND LONGITUDE FORMAT 00018570 1000 FORMAT (2(F4.0,F2.0,F8.5,6X)) 00018580 C ORIGIN AND DESTINATION SYSTEMS FORMAT 00018590 1001 FORMAT (3I5) 00018600 C PLANE CONSTANTS INPUT FORMAT 00018610 1003 FORMAT (5D16.10) 00018620 1009 FORMAT (///' TRANSFORMATION RUN TIME WAS ',F7.2,' SECONDS'/) 00018630 C GEOGRAPHIC POSITION AND ELEVATION OUTPUT FORMAT 00018640 1018 FORMAT (1X,I10,1X,I4,I3,1X,F8.5,5X,I4,I3,1X,F8.5,4X,F11.5,I8,I3,F500018650 1.1) 00018660 1020 FORMAT (/' ECCENTRECITY =',D14.8,' MAJOR SEMI-DIAMETER = ',D16.100018670 10/) 00018680 1021 FORMAT (' ORIGINAL SYSTEM = ',I4,' FINAL SYSTEM = ',I4/) 00018690 1022 FORMAT (' SECANT ORIGIN LAT ',I4,I3,F9.5,' LONG ',I4,I3,F9.5/) 00018700 1023 FORMAT (/' PLANE CONSTANTS'/) 00018710 1024 FORMAT (1X,D16.10) 00018720 1026 FORMAT (///T5,'STATION',T18,'LATITUDE',T38,'LONGITUDE',T54,'ELEVAT00018730 1ION FT.',T74,'CONV.'/) 00018740 1027 FORMAT (///' PLANE, GEOCENTRIC OR UTM OUTPUT COORDINATES'/) 00018750 1028 FORMAT (///' SECANT PLANE COORDINATE OUTPUT X,Y,Z IN METERS'/) 00018760 C SECANT PLANE, UTM, GEOCENTRIC OR PLANE COORDINATE OUTPUT FORMAT 00018770 1036 FORMAT (1X,I9,3(3X,D18.12)) 00018780 1040 FORMAT (///T5,'STATION',T18,'X COORD',T33,'Y COORD',T48,'ELEVATION00018790 1',T66,'CONV.'/) 00018800 1041 FORMAT (1X,I10,3F15.3,I8,I3,F5.1) 00018810 END 00018820 C 00018830 C 00018840 SUBROUTINE STERGP 00018850 C 00018860 IMPLICIT REAL*8(A-H,O-Z) 00018870 REAL*8G(4,500),CONV(500),T(11) 00018880 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00018890 COMMON/AREA3/G,T,PHDG,PHMIN,PHSEC,PLGDG,PLGMN,PLSEC,AMTFT,FTMT, 00018900 1SECRD,ESQ,CONV,A,SNONS 00018910 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00018920 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00018930 2JCON,INSTRU 00018940 DATA RAD/57.2957795D0/,E1/676.8658D-5/,E2/681.4785D-5/,FU/.434294500018950 1D0/ 00018960 C TRANSFORMATION AND INVERSE FOR STEREOGRAPHIC AND GEOGRAPHIC POSITION 00018970 IF(IDIR.EQ.1)GO TO 400 00018980 C STEREOGRAPHIC TO GEOGRAPHIC POSITION 00018990 DO 350 N=1,NN 00019000 Z=(G(2,N) -1.D6)/1.D6 00019010 W=(G(3,N) -1.D6)/1.D6 00019020 G(3,N)=.665D2+(-.142952945665D5*Z-.7186600649D3*Z*W-.7206D-3*Z*Z 00019030 1-.442990316D2*Z*W*W+.148300635D2*Z**3-.25941871D1*Z*W**3+.2630264800019040 2D1*Z**3*W-.1245733D0*Z**5)/3600.D0 00019050 G(2,N)=.465D2+(.987200337240D4*W-.2481452114D3*Z*Z-.2400224D1*W*W-00019060 1.179973728D2*Z*Z*W-.18776066D1*W**3+.2988811D0*Z**4-.923441D0*Z*Z 00019070 2*W*W+.573046D-1*Z**4*W)/3600.D0 00019080 CONV(N)=((G(3,N)-.665D2)*3600.D0*DSIN(((G(2,N)+.465D2)/2.D0)/RAD))00019090 1/3600.D0 00019100 C WEST LONGITUDE IS NEGATIVE 00019110 G(3,N)=-G(3,N)*3600.D0 00019120 350 G(2,N)=G(2,N)*3600.D0 00019130 RETURN 00019140 C GEOGRAPHIC POSITION TO STEROGRAPHIC 00019150 400 T3=FU/(24.D0*RAD*RAD) 00019160 T31=T3*2.D0 00019170 T32=T3*3.D0 00019180 DO 450 N=1,NN 00019190 DT=G(2,N)/3600.D0 00019200 DG=-G(3,N)/3600.D0 00019210 DT1=DT-.465D2 00019220 DG1=.665D2-DG 00019230 IF(DT1)42,41,42 00019240 41 IF(DG1)42,43,42 00019250 43 CONV(N)=0.D0 00019260 G(2,N)=0.D0 00019270 G(3,N)=0.D0 00019280 GO TO 450 00019290 42 DT12=DT1*DT1 00019300 DG12=DG1*DG1 00019310 SL=DSIN(((DT+.465D2)/2.D0)/RAD) 00019320 SL2=SL*SL 00019330 CL=DCOS(((DT+.465D2)/2.D0)/RAD) 00019340 CL2=CL*CL 00019350 TL=SL/CL 00019360 TL2=TL*TL 00019370 FN2=E2*CL2 00019380 V2=1.D0+FN2 00019390 V4=V2*V2 00019400 BN=A/DSQRT(1.D0-E1*SL2) 00019410 BM=A*(1.D0-E1)/DSQRT((1.D0-E1*SL2)**3) 00019420 T1=RAD/BM 00019430 T2=RAD/BN 00019440 T4=(1.D0+FN2-9.D0*FN2*TL2)/V4*T3 00019450 T5=T3*(1.D0-2.D0*FN2) 00019460 T6=T32*(FN2*(TL2-1.D0-FN2-4.D0*FN2*TL2)/V4) 00019470 T7=T31*V2 00019480 T8=T3*(3.D0+8.D0*FN2+5.D0*FN2*FN2)/V4 00019490 SSA=DG1*CL/T2*(1.D0-T3*DG12*SL2/FU+T4*DT12/FU) 00019500 SCA=DT1*DCOS((DG1/2.D0)/RAD)/T1*(1.D0+T5*DG12*CL2/FU-T6*DT12/FU) 00019510 DA=DG1*SL*(1.D0+T7*DG12*CL2/FU+T8*DT12/FU) 00019520 IF(DT1)30,31,30 00019530 31 IF(DG1)32,33,33 00019540 32 BG=90.D0 00019550 GO TO 4 00019560 33 BG=270.D0 00019570 GO TO 4 00019580 30 BG=DATAN(SSA/SCA)*RAD 00019590 IF(DT1)1,2,2 00019600 1 BG=BG+180.D0 00019610 2 IF(DG1)3,4,4 00019620 3 BG=BG+360.D0 00019630 4 SMA=DSIN(BG/RAD) 00019640 CMA=DCOS(BG/RAD) 00019650 IF(DABS(SMA)-DABS(CMA))5,6,6 00019660 5 S=SCA/CMA 00019670 GO TO 7 00019680 6 S=SSA/SMA 00019690 7 GA=BG-DA/2.D0 00019700 IF(GA)8,9,9 00019710 8 GA=GA+360.D0 00019720 9 GD=S/AMTFT 00019730 GD=GD+1.902D0*1.D-16*GD**3-8.846*1.D-5*GD+4.343D0*1.D-32*GD**5 00019740 G(2,N)=GD*DSIN(GA/RAD)+1.D6 00019750 G(3,N)=GD*DCOS(GA/RAD)+1.D6 00019760 CONV(N)=DABS(DA) 00019770 450 CONTINUE 00019780 RETURN 00019790 END 00019800 C 00019810 C 00019820 SUBROUTINE DEG(ANG,ID,IM,SS) 00019830 C 00019840 IMPLICIT REAL*8(A-H,O-Z) 00019850 ID=ANG 00019860 AID=ID 00019870 FID=(ANG-AID)*60.D0 00019880 IM=FID 00019890 FIS=IM 00019900 SS=(FID-FIS)*60.D0 00019910 RETURN 00019920 END 00019930 SUBROUTINE GPSP 00019940 C 00019950 IMPLICIT REAL*8(A-H,O-Z) 00019960 REAL*8G(4,500),CONV(500),T(11) 00019970 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00019980 COMMON/AREA3/G,T,PHDG,PHMIN,PHSEC,PLGDG,PLGMN,PLSEC,AMTFT,FTMT, 00019990 1SECRD,ESQ,CONV,A,SNONS 00020000 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00020010 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00020020 2JCON,INSTRU 00020030 DATA PI/3.1415926535898D0/ 00020040 C TRANSFORMATION AND INVERSE FOR GEOGRAPHIC POSITIONS AND SECANT PLA00020050 C SOLVE FOR THE TRANSFORMATION CONSTANTS DEPENDENT ON ELLIPSOID 00020060 C PARAMETERS AND SECANT PLANE SYSTEM ORIGIN 00020070 PANME=1.D0-ESQ 00020080 IF(PHDG)3,4,4 00020090 3 PSGLAT=-1.D0 00020100 PHMIN=-PHMIN 00020110 PHSEC=-PHSEC 00020120 GO TO 11 00020130 4 PSGLAT=1.D0 00020140 11 PHRAD=(PHDG*3600.D0+PHMIN*60.D0+PHSEC)*SNONS 00020150 IF(PLGDG)5,6,6 00020160 5 PSGLNG=-1.D0 00020170 PLGMN=-PLGMN 00020180 PLSEC=-PLSEC 00020190 PNT1=0.D0 00020200 PNT2=-PI 00020210 PNT3=-PI/2.D0 00020220 GO TO 13 00020230 6 PSGLNG=1.D0 00020240 PNT1=PI 00020250 PNT2=0.D0 00020260 PNT3=PI/2.D0 00020270 13 PLRAD=(PLGDG*3600.D0+PLGMN*60.D0+PLSEC)*SNONS 00020280 PSNPHO=DSIN(PHRAD) 00020290 PCSPHO=DCOS(PHRAD) 00020300 PSNLGO=DSIN(PLRAD) 00020310 PCSLGO=DCOS(PLRAD) 00020320 PEPHM=DSQRT(1.D0-ESQ*PSNPHO*PSNPHO) 00020330 C COMPUTE LENGTH OF NORMAL THROUGH SYSTEM ORIGIN 00020340 PENO=A/PEPHM 00020350 PENSYT=PENO*ESQ*PSNPHO 00020360 PSNPL=PSNPHO*PSNLGO 00020370 PCSPL=PCSPHO*PSNLGO 00020380 PSCPLA=PSNPHO*PCSLGO 00020390 PCSPLA=PCSPHO*PCSLGO 00020400 C COMPUTE MAGNITUDE OF TRANSLATION FOR X-Y PLANE 00020410 IBASE=PENO/1000.D0 00020420 PBASE=IBASE*1000 00020430 IF(IDIR.EQ.1)GO TO 99 00020440 C SECANT PLANE TO GEOGRAPHIC POSITION 00020450 50 DO 98 N=1,NN 00020460 C TRANSFORM SECANT PLANE COORDINATES TO MODIFIED GEOCENTRIC COORDS, 00020470 PZPB=G(4,N)+PBASE 00020480 PX=G(2,N)*PCSLGO-G(3,N)*PSNPL+PZPB*PCSPL 00020490 PY=-(G(2,N)*PSNLGO+G(3,N)*PSCPLA)+PZPB*PCSPLA 00020500 PZ=G(3,N)*PCSPHO+PZPB*PSNPHO 00020510 PZAPP=PZ 00020520 PDEM=DSQRT(PX*PX+PY*PY) 00020530 DO 80 IX=1,5 00020540 PTNPHP=PZ/PDEM 00020550 PHR=DATAN(PTNPHP) 00020560 PSNPHP=DSIN(PHR) 00020570 PCSPHP=DCOS(PHR) 00020580 PENPH=DSQRT(1.D0-ESQ*PSNPHP*PSNPHP) 00020590 PENPT=A/PENPH 00020600 80 PZ=-PENSYT+(PENPT*ESQ*PSNPHP)+PZAPP 00020610 81 PALAM=DATAN(PX/PY) 00020620 IF(PALAM)74,75,76 00020630 74 PANGLE=PALAM+PNT1 00020640 GO TO 77 00020650 75 PANGLE=PNT3 00020660 GO TO 77 00020670 76 PANGLE=PALAM+PNT2 00020680 77 G(3,N)=PANGLE*SECRD 00020690 PSINL=DSIN(PANGLE) 00020700 PH=(PX/(PCSPHP*PSINL))-PENPT 00020710 G(4,N)=PH*FTMT 00020720 G(2,N)=PHR*SECRD 00020730 98 CONTINUE 00020740 IF (PHDG)30,31,31 00020750 30 PHMIN=-PHMIN 00020760 PHSEC=-PHSEC 00020770 31 IF(PLGDG)32,33,33 00020780 32 PLGMN=-PLGMN 00020790 PLSEC=-PLSEC 00020800 33 RETURN 00020810 C GEOGRAPHIC POSITION TO SECANT PLANE 00020820 99 DO 199 N=1,NN 00020830 PELEV=G(4,N) 00020840 C CONVERT LATITUDE AND LONGITUDE TO RADIANS 00020850 PPRAD=G(2,N)*SNONS 00020860 PLGRD=G(3,N)*SNONS 00020870 PSNPHI=DSIN(PPRAD) 00020880 PCSPHI=DCOS(PPRAD) 00020890 PSNLAM=DSIN(PLGRD) 00020900 PCSLAM=DCOS(PLGRD) 00020910 PH=G(4,N)*AMTFT 00020920 PAMSQ=DSQRT(1.D0-ESQ*PSNPHI*PSNPHI) 00020930 C COMPUTE NORMAL OF STATION 00020940 PRN=A/PAMSQ 00020950 C COMPUTE CLASSICAL GEOCENTRIC COORDINATES 00020960 PSNH=PRN+PH 00020970 PX=PSNLAM*PCSPHI*PSNH 00020980 PY=PCSLAM*PCSPHI*PSNH 00020990 PZ=((PRN*PANME+PH)*PSNPHI)+PENSYT 00021000 G(2,N)=PX*PCSLGO-PY*PSNLGO 00021010 G(3,N)=PZ*PCSPHO-(PX*PSNPL+PY*PSCPLA) 00021020 G(4,N)=PX*PCSPL+PY*PCSPLA+PZ*PSNPHO-PBASE 00021030 199 CONTINUE 00021040 IF(PHDG)34,35,35 00021050 34 PHMIN=-PHMIN 00021060 PHSEC=-PHSEC 00021070 35 IF(PLGDG)36,200,200 00021080 36 PLGMN=-PLGMN 00021090 PLSEC=-PLSEC 00021100 200 RETURN 00021110 END 00021120 C 00021130 SUBROUTINE GEOCGP 00021140 IMPLICIT REAL*8(A-H,O-Z) 00021150 REAL*8G(4,500),CONV(500),T(11) 00021160 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00021170 COMMON/AREA3/G,T,PHDG,PHMIN,PHSEC,PLGDG,PLGMN,PLSEC,AMTFT,FTMT, 00021180 1SECRD,ESQ,CONV,A,SNONS 00021190 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00021200 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00021210 2JCON,INSTRU 00021220 DATA PI/3.1415926535898D0/ 00021230 IF(IDIR.EQ.1) GO TO 800 00021240 DO 750 N=1,NN 00021250 X=G(3,N)/G(2,N) 00021260 IF(X.LT.10.D0)GO TO 710 00021270 FLONG=DATAN(G(2,N)/G(3,N)) 00021280 X=0.5D0*PI-DABS(FLONG) 00021290 IF(FLONG.LT.0.D0)FLONG=-X 00021300 IF(FLONG.GE.0.D0)FLONG=X 00021310 GO TO 715 00021320 710 FLONG=DATAN(X) 00021330 IF(G(2,N).LT.0.D0.AND.G(3,N).GT.0.D0)FLONG=PI+FLONG 00021340 IF(G(2,N).LT.0.D0.AND.G(3,N).LT.0.D0)FLONG=-PI+FLONG 00021350 715 PHIEST=DATAN(G(4,N)/DSQRT(G(2,N)*G(2,N)+G(3,N)*G(3,N))) 00021360 DO 720 I=1,10 00021370 IF(I.EQ.10)WRITE (NPRNT,1000) 00021380 S=DSIN(PHIEST) 00021390 S=S*S 00021400 TEMP=A*DSQRT(1.D0-((ESQ*S*(1.D0-ESQ))/(1.D0-ESQ*S))) 00021410 H=DSQRT(G(2,N)*G(2,N)+G(3,N)*G(3,N)+G(4,N)*G(4,N)) -TEMP 00021420 BSQ=A*A*(1.D0-ESQ) 00021430 TEMP=BSQ+A*H*DSQRT(1.D0-ESQ*S) 00021440 TEMP=(A*A+A*H*DSQRT(1.D0-ESQ*S))/TEMP 00021450 PHI=DATAN(G(4,N)*TEMP/DSQRT(G(2,N)*G(2,N)+G(3,N)*G(3,N))) 00021460 IF(DABS(PHIEST-PHI).LT.1.D-15)GO TO 730 00021470 PHIEST=PHI 00021480 720 CONTINUE 00021490 730 G(2,N)=PHI*SECRD 00021500 G(3,N)=FLONG*SECRD 00021510 750 G(4,N)=H*FTMT 00021520 RETURN 00021530 C GEOGRAPHIC POSITION TO GEOCENTRIC 00021540 800 DO 850 N=1,NN 00021550 H=G(4,N)*AMTFT 00021560 SP=DSIN(G(2,N)*SNONS) 00021570 CP=DCOS(G(2,N)*SNONS) 00021580 SL=DSIN(G(3,N)*SNONS) 00021590 CL=DCOS(G(3,N)*SNONS) 00021600 D=DSQRT(1.D0-ESQ*SP*SP) 00021610 G(2,N)=(A*CP*CL/D)+H*CP*CL 00021620 G(3,N)=(A*CP*SL/D)+H*CP*SL 00021630 850 G(4,N)=(A*(1.D0-ESQ)*SP/D)+H*SP 00021640 RETURN 00021650 1000 FORMAT(//' 10 ITERATIONS ON LATITUDE NECESSARY'//) 00021660 END 00021670 C 00021680 C 00021690 SUBROUTINE UTMGP 00021700 C 00021710 IMPLICIT REAL*8(A-H,O-Z) 00021720 REAL*8G(4,500),CONV(500),T(11) 00021730 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00021740 COMMON/AREA3/G,T,PHDG,PHMIN,PHSEC,PLGDG,PLGMN,PLSEC,AMTFT,FTMT, 00021750 1SECRD,ESQ,CONV,A,SNONS 00021760 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00021770 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00021780 2JCON,INSTRU 00021790 DATA CONS1/.24682D0/,CONS2/30.02335D0/,CONS3/5078.64977D0/, 00021800 1CONS4/.1570499810D0/,SEPQ/.0408887094D0/,Y1/.11422D0/,Y2/21.73607 00021810 2D0/,Y3/5104.57388D0/,Y4/6367399.689D0/,CNII/.2424068406D-1/, 00021820 3CNIII/.0195870255D-21/,CNV/.0391740509D-2/, 00021830 4CNA6/.0153460626D-6/,CNB5/.0460381879D-6/,CNI3/7.83481018D0/, 00021840 5CNC5/.36830550D-2/ 00021850 C TRANSFORMATION AND INVERSE FOR UNIVERSAL TRANSVERSE MERCATOR 00021860 C TO GEOGRAPHIC POSITION ABOVE CONSTANTS FOR CLARKE 1866 00021870 IF(IABS(IFROM).LE.14.AND.IDIR.EQ.0)GO TO 703 00021880 IF(IABS(IFROM).GT.14.AND.IDIR.EQ.0)GO TO 700 00021890 IF(IABS(ITO).LE.14.AND.IDIR.EQ.1)GO TO 704 00021900 I=ITO 00021910 GO TO 701 00021920 700 I=IFROM 00021930 C SET SCALE AND ORIGIN PARAMETERS FOR UTM SYSTEM 00021940 701 ALORG=0.D0 00021950 SCALE=.9996D0 00021960 FE=500000.D0*FTMT 00021970 IX=IABS(I)-100 00021980 X=IX 00021990 CM=(6.D0*X-183.D0)*3600.D0 00022000 GO TO 706 00022010 703 I=IFROM 00022020 GO TO 705 00022030 704 I=ITO 00022040 C ALASKA ZONES 00022050 705 ALORG=5985100.923D0 00022060 SCALE=.9999 00022070 IX=I-2 00022080 FE=500000.D0 00022090 IF(IX.EQ.7)FE=700000.D0 00022100 IF(IX.EQ.9)FE=600000.D0 00022110 X=IX 00022120 CM=482400.D0+X*14400.D0 00022130 706 EPSQ=ESQ/(1.D0-ESQ) 00022140 IF(IDIR.EQ.1)GO TO 800 00022150 C UTM TO GEOGRAPHIC 00022160 DO 750 N=1,NN 00022170 IF(I.LE.0)G(3,N)=10000000.D0-G(3,N) 00022180 C CONVERT FROM METERS TO FEET IF UTM 00022190 IF(IABS(I).LE.14)GO TO 702 00022200 G(2,N)=FTMT*G(2,N) 00022210 G(3,N)=FTMT*G(3,N) 00022220 702 PPRD=((G(3,N)*AMTFT+ALORG)*CONS4*10.D-7)/SCALE 00022230 SNLT=DSIN(PPRD) 00022240 CSLT=DCOS(PPRD) 00022250 CSSQ=CSLT*CSLT 00022260 PHRD=((CONS1*CSSQ+CONS2)*CSSQ+CONS3)*(SNLT*CSLT)*10.D-7+PPRD 00022270 PHIS=PHRD*SECRD 00022280 Q=((G(2,N)-FE)*AMTFT)*10.D-7 00022290 IF(Q)720,710,720 00022300 710 DLAM=0.D0 00022310 G(2,N)=PHIS 00022320 GO TO 740 00022330 720 SNLAT=DSIN(PHRD) 00022340 CSLAT=DCOS(PHRD) 00022350 SNSQ=SNLAT*SNLAT 00022360 CCSQ=CSLAT*CSLAT 00022370 TNLAT=SNLAT/CSLAT 00022380 TNSQ=TNLAT*TNLAT 00022390 ENU=A/DSQRT(1.D0-ESQ*SNSQ) 00022400 ENSNS=ENU*SNONS 00022410 EPCS=EPSQ*CCSQ 00022420 EPCSQ=1.D0+EPCS 00022430 QSQ=Q*Q 00022440 QCU=QSQ*Q 00022450 QFR=QCU*Q 00022460 QFV=QFR*Q 00022470 QSX=QFV*Q 00022480 SCLAT=1.D0/CSLAT 00022490 ENSN5=ENU**4*ENSNS 00022500 SVN=(((TNLAT/(2.D0*ENU*ENSNS))*EPCSQ)/SCALE**2)*10.D11 00022510 EGHT=((TNLAT/(24.D0*ENU**3*ENSNS))*(5.D0+3.D0*TNSQ+SEPQ*(CCSQ-SNSQ00022520 1)-(3.D0*EPSQ**2*CCSQ)*(CCSQ+3.D0*SNSQ))/SCALE**4)*10.D23 00022530 D6=(QSX*(TNLAT/(720.D0*ENU**5*ENSNS))*(61.D0+(45.D0*TNSQ)*(2.D0+ 00022540 1TNSQ-EPSQ*SNSQ)+EPSQ*(107.D0*CCSQ-162.D0*SNSQ))/SCALE**6)*10.D35 00022550 ANINE=(SCLAT/ENSNS)/SCALE*10.D5 00022560 TEN=(SCLAT/(6.D0*ENU**2*ENSNS))*(1.D0+2.D0*TNSQ+EPCS)/SCALE**3*10.00022570 1D17 00022580 E5=QFV*(SCLAT/(120.D0*ENSN5))*(5.D0+(4.0D0*TNSQ)*(7.D0+6.D0*TNSQ) 00022590 1+(2.D0*EPSQ)*(3.D0*CCSQ+4.D0*SNSQ))/SCALE**5*10.D29 00022600 G(2,N)=PHIS-SVN*QSQ+EGHT*QFR-D6 00022610 IF(IFROM.LT.0)G(2,N)=-G(2,N) 00022620 DLAM=ANINE*Q-TEN*QCU+E5 00022630 740 G(3,N)=CM-DLAM 00022640 IF(IABS(IFROM).GE.15)G(3,N)=CM+DLAM 00022650 IF(IABS(IFROM).LE.11)G(3,N)=-G(3,N) 00022660 750 CONTINUE 00022670 RETURN 00022680 C GEOGRAPHIC POSITION TO UTM 00022690 800 DO 850 N=1,NN 00022700 IF(IABS(ITO).LE.11)G(3,N)=-G(3,N) 00022710 PHRD=G(2,N)*SNONS 00022720 SNLAT=DSIN(PHRD) 00022730 CSLAT=DCOS(PHRD) 00022740 SNSQ=SNLAT*SNLAT 00022750 CCSQ=CSLAT*CSLAT 00022760 TNLAT=SNLAT/CSLAT 00022770 TNSQ=TNLAT*TNLAT 00022780 ENU=A/DSQRT(1.D0-ESQ*SNSQ) 00022790 ENSNS=ENU*SNONS 00022800 EPCS=EPSQ*CCSQ 00022810 EPCSQ=1.D0+EPCS 00022820 RWON=(Y4*(PHRD-(SNLAT*CSLAT*10.D-7)*(Y3-CCSQ*(Y2-Y1*CCSQ))))*SCALE00022830 P=(CM-G(3,N))*10.D-5 00022840 IF (IABS(ITO).GE.15) P=-P 00022850 IF(P)819,820,819 00022860 820 G(2,N)=FE 00022870 G(3,N)=(RWON-ALORG)*FTMT 00022880 GO TO 825 00022890 819 PSQ=P*P 00022900 PCU=PSQ*P 00022910 PFR=PCU*P 00022920 PFV=PFR*P 00022930 PSX=PFV*P 00022940 FOUR=ENSNS*CSLAT*SCALE*10.D3 00022950 TWO=FOUR*SNLAT*CNII 00022960 EPCTN=EPCSQ-TNSQ 00022970 THRE=TWO*CNIII*CCSQ*(4.D0*EPCSQ*EPCSQ+EPCTN) 00022980 FIVE=FOUR*CNV*CCSQ*EPCTN 00022990 CSFR=CCSQ*CCSQ 00023000 CSFRR=1.D0/CSFR 00023010 A6=PSX*TWO*CNA6*CSFR*(60.D0*EPCTN+CSFRR+540.D0*EPCS-330.D0*EPSQ) 00023020 B5=PFV*FOUR*CNB5*CSFR*(20.D0*EPCTN+CSFRR+52.D0*EPCS-58.D0*EPSQ-16.00023030 1D0) 00023040 G(2,N)=(P*FOUR+PCU*FIVE+B5)*FTMT+FE 00023050 G(3,N)=((PSQ*TWO+PFR*THRE+RWON+A6)-ALORG)*FTMT 00023060 C CONVERT FEET TO METERS IF UTM CONVERSION 00023070 825 IF(IABS(ITO).LE.14)GO TO 850 00023080 G(2,N)=G(2,N)*AMTFT 00023090 G(3,N)=G(3,N)*AMTFT 00023100 850 IF(ITO.LE.0)G(3,N)=10000000.D0-G(3,N) 00023110 RETURN 00023120 END 00023130 C 00023140 C 00023150 SUBROUTINE INMODL 00023160 C 00023170 C SCHUT INDEPENDENT MODEL STRIP FORMATION NRC 00023180 C 00023190 IMPLICIT REAL*8(A-H,O-Z) 00023200 REAL*8 AX(22),AY(22),AZ(22),AF(11),S(13) 00023210 INTEGER ISTRP(20),IRCSP(20),ICAN(20),ITIE(107),IGRND(300) 00023220 2,LS(22),LP(22) 00023230 COMMON/AREA2/ISTRP,IRCSP,ISTP,IS1A0,ICAN,I50,ICOORD,IFROM,ITO,N7, 00023240 1IDIR,NN,NPRNT,NREAD,NPUNC,JS1,JS2,KD3,N5,ITIE,IGRND,IEARTH,ICON, 00023250 2JCON,INSTRU 00023260 COMMON/AREA3/AX,AY,AZ,AF 00023270 EQUIVALENCE (S11,S(1)),(S12,S(4)),(S13,S(7)),(S14,S(10)), 00023280 1 (G1,S21,S(2)),(S22,S(5)),(S23,S(8)),(S24,S(11)), 00023290 2 (G2,S31,S(3)),(G3,S32,S(6)),(S33,S(9)),(S34,S(12)), 00023300 3 (G4,S(13)) 00023310 1 FORMAT (I4, I5, 3F9.5) 00023320 2 FORMAT (I4, I5, 3I9) 00023330 3 FORMAT (1H 2I7, 3P6F12.2) 00023340 4 FORMAT (1H1) 00023350 5 FORMAT (14H0ERROR AT CARD 2I7) 00023360 140 FORMAT (8I10) 00023370 153 FORMAT (/T20,'STRIP CONTAINS MORE THAN 500 POINTS') 00023380 154 FORMAT (T20,'INDEPENDENT MODEL RUN TIME WAS',F9.2,' SECONDS') 00023390 360 FORMAT (/T4,'MODEL',T11,'POINT',T24,'X',T36,'Y',T48,'Z',T61,'DX', 00023400 1T73,'DY',T85,'DZ'//) 00023410 WRITE (NPRNT,360) 00023420 IRUN=0 00023430 IEXIT=1 00023440 I20=KLOCK(J) 00023450 100 READ (NREAD,140)ISTRIP 00023460 IF(IRUN.EQ.0)GO TO 151 00023470 IRCSP(ISTP)=IREC 00023480 WRITE (NPRNT,4) 00023490 WRITE (NPRNT,360) 00023500 IF(IREC.GT.500)GO TO 152 00023510 151 IREC=0 00023520 ISTP=ISTP+1 00023530 ISTRP(ISTP)=ISTRIP 00023540 IRUN=2 00023550 J = 1 00023560 MA = 1 00023570 MB = 1 00023580 MC = 1 00023590 GO TO 801 00023600 101 MBB = 3 00023610 IF (J-1) 999, 110, 201 00023620 110 MB = 2 00023630 MC = 2 00023640 GO TO 801 00023650 200 MA = 1 00023660 MBB = 4 00023670 201 IF (J-12) 202, 202, 990 00023680 202 K = J - 1 00023690 MC = 3 00023700 XP = AX(1) 00023710 YP = AY(1) 00023720 ZP = AZ(1) 00023730 DO 203 L = 1,13 00023740 203 S(L) = 0.D0 00023750 GO TO 801 00023760 204 IF (J-K-1) 999, 205, 210 00023770 205 XQ = AX(K+1) 00023780 YQ = AY(K+1) 00023790 ZQ = AZ(K+1) 00023800 GO TO 800 00023810 210 JK = J - K 00023820 X1 = AX(JK) - XP 00023830 Y1 = AY(JK) - YP 00023840 Z1 = AZ(JK) - ZP 00023850 X2 = AX(J) - XQ 00023860 Y2 = AY(J) - YQ 00023870 Z2 = AZ(J) - ZQ 00023880 D1 = DSQRT(X1*X1 + Y1*Y1 + Z1*Z1) 00023890 D2 = DSQRT(X2*X2 + Y2*Y2 + Z2*Z2) 00023900 X1 = X1 / D1 00023910 Y1 = Y1 / D1 00023920 Z1 = Z1 / D1 00023930 X2 = X2 / D2 00023940 Y2 = Y2 / D2 00023950 Z2 = Z2 / D2 00023960 AF(JK) = D1 / D2 00023970 G4 = G4 + (X1+X2)**2 + (Y1+Y2)**2 + (Z1+Z2)**2 00023980 S11 = S11 + X1 * X2 00023990 S22 = S22 + Y1 * Y2 00024000 S33 = S33 + Z1 * Z2 00024010 S14 = S14 + Z1 * Y2 00024020 S24 = S24 + X1 * Z2 00024030 S34 = S34 + Y1 * X2 00024040 G1 = G1 + Z2 * Y1 00024050 G2 = G2 + X2 * Z1 00024060 G3 = G3 + Y2 * X1 00024070 IF (JK - K) 800, 300, 999 00024080 300 S44 = G4 -4.D0* (S11 + S22 + S33) 00024090 S11 = G4 -4.D0* S11 00024100 S22 = G4 -4.D0* S22 00024110 S33 = G4 -4.D0* S33 00024120 S12 =-2.D0* (G3 + S34) 00024130 S13 =-2.D0* (G2 + S24) 00024140 S23 =-2.D0* (G1 + S14) 00024150 S14 =2.D0 * (G1 - S14) 00024160 S24 =2.D0 * (G2 - S24) 00024170 S34 =2.D0 * (G3 - S34) 00024180 S22 = S22 - S12 * S12 / S11 00024190 S23 = S23 - S13 * S12 / S11 00024200 S24 = S24 - S14 * S12 / S11 00024210 S33 = S33 - S13 * S13 / S11 - S23 * S23 / S22 00024220 S34 = S34 - S14 * S13 / S11 - S24 * S23 / S22 00024230 S44 = S44 - S14 * S14 / S11 - S24 * S24 / S22 00024240 IF (S44 - S33) 301, 301, 302 00024250 301 D = 1.D0 00024260 C = -S34 / S33 00024270 GO TO 303 00024280 302 C = 1.D0 00024290 D = -S34 / S44 00024300 303 B = -(S24 * D + S23 * C) / S22 00024310 A = -(S14 * D + S13 * C + S12 * B) /S11 00024320 S11 = D*D + A*A - B*B - C*C 00024330 S22 = D*D - A*A + B*B - C*C 00024340 S33 = D*D - A*A - B*B + C*C 00024350 S21 =2.D0* (A*B + C*D) 00024360 S12 =2.D0* (A*B - C*D) 00024370 S13 =2.D0* (A*C + B*D) 00024380 S31 =2.D0* (A*C - B*D) 00024390 S32 =2.D0* (B*C + A*D) 00024400 S23 =2.D0* (B*C - A*D) 00024410 F = 0.D0 00024420 G4 = K - 1 00024430 DO 304 L = 2,K 00024440 304 F = F + AF(L) / G4 00024450 G4 = F / (D*D + A*A + B*B + C*C) 00024460 DO 305 L = 1,9 00024470 305 S(L) = S(L) * G4 00024480 J = 1 00024490 MC = 4 00024500 LS(1) = LS(K+1) 00024510 GO TO 820 00024520 401 MB = MBB 00024530 J = K + 1 00024540 402 J = J + 1 00024550 IF (J - 2 * K) 810, 810, 410 00024560 410 J = 1 00024570 MA = 2 00024580 MB = 2 00024590 GO TO 801 00024600 411 GO TO (401, 801, 402, 402), MB 00024610 500 MB = 1 00024620 MC = 5 00024630 GO TO 801 00024640 800 J = J + 1 00024650 801 READ (NREAD,1) LS(J), LP(J), AX(J), AY(J), AZ(J) 00024660 IF (LS(J)) 902, 901, 802 00024670 802 GO TO (811, 810), MA 00024680 810 X2 = AX(J) - XQ 00024690 Y2 = AY(J) - YQ 00024700 Z2 = AZ(J) - ZQ 00024710 AX(J) = S11 * X2 + S12 * Y2 + S13 * Z2 + XP 00024720 AY(J) = S21 * X2 + S22 * Y2 + S23 * Z2 + YP 00024730 AZ(J) = S31 * X2 + S32 * Y2 + S33 * Z2 + ZP 00024740 JK = J - K 00024750 811 GO TO (900, 820, 831, 830), MB 00024760 820 WRITE (NPRNT,3) LS(J), LP(J), AX(J), AY(J), AZ(J) 00024770 GO TO 840 00024780 830 AX(J) =.5D0* (AX(J) + AX(JK)) 00024790 AY(J)=.5D0*(AY(J)+AY(JK)) 00024800 AZ(J) =.5D0* (AZ(J) + AZ(JK)) 00024810 831 DX = AX(J) - AX(JK) 00024820 DY = AY(J) - AY(JK) 00024830 DZ = AZ(J) - AZ(JK) 00024840 IF (LP(J) - LP(JK)) 832, 833, 832 00024850 832 LP(J) = 0 00024860 833 WRITE (NPRNT,3) LS(J),LP(J),AX(J),AY(J),AZ(J),DX,DY,DZ 00024870 840 IREC=IREC+1 00024880 WRITE (JS1)LP(J),AX(J),AY(J),AZ(J) 00024890 900 GO TO (800, 801, 204, 411, 800), MC 00024900 901 GO TO (800, 990, 990, 999, 990), MC 00024910 902 IF (LS(J) + 2) 905, 904, 903 00024920 903 GO TO (101, 990, 990, 990, 200), MC 00024930 904 GO TO (990, 500, 990, 500, 990), MC 00024940 905 GO TO (990, 990, 990, 100, 990), MC 00024950 990 WRITE (NPRNT,5) LS(J), LP(J) 00024960 999 IRCSP(ISTP)=IREC 00024970 I21=KLOCK(J) 00024980 TIME=(I21-I20)/100. 00024990 WRITE (NPRNT,154)TIME 00025000 IF(IREC.LE.500)RETURN 00025010 IEXIT=2 00025020 152 WRITE (NPRNT,153) 00025030 I50=I50+1 00025040 ICAN(I50)=ISTRIP 00025050 IF(IEXIT.EQ.2)RETURN 00025060 GO TO 151 00025070 END 00025080