PROGRAM TRIPLT INTEGER PLTFIL, RAYFIL PARAMETER ( PLTFIL = 1, RAYFIL = 21 ) COMMON /BOTTOM/ RB(8000), ZB(8000) C DIMENSION PLTITL(20), CODES(5000), NRPLOT(65), ZR(1000), TR(1000), 1 SS(1000), TIME(1000), PHASE(1000), NBR(1000), NTU(1000), 2 NSR(1000), NTO(1000), R(8000) COMMON /TRAJCT/ ZSTOR(36000), TSTOR(36000) DATA ZSTOR/36000*-1./ COMMON /HITS/ NBHITS, NSHITS INTEGER*2 NBHIT(36000), NSHIT(36000) COMMON /GRAPH/ R1KM, DRKM, R2KM, Z1, DZ, Z2, XLG, YHGT DATA YHGT /6./, TENSE/.1/ REAL RBKM(8000) C 900 FORMAT(20A4) 901 FORMAT(16I5) 902 FORMAT(10F8.0) 903 FORMAT(' PLOT AXIS RANGES(KM), DEPTH(M), LGTH/HGT(IN)=', 5F8.1) 904 FORMAT(' DB, NSR, NBR LIMITS, TITLE FLAG=', F8.1, 8I5) 905 FORMAT(5H0NRAY,20I5) 906 FORMAT(1X, I6, ' RAYS;', I10, ' BLOCKS OF LENGTH=', I6) C C***** THIS PROGRAM READS THE RAY-HISTORY FILE CREATED BY TRIMAIN, C FOR001.DAT, AND PLOTS SELECTED RAY PATHS IN THE RANGE-DEPTH PLANE. C***** PLOT IS FROM RANGES R1KM TO R2KM AND DEPTHS 0 TO ZMAX. C HEIGHT IS 8 INCHES, LENGTH IN INCHES=XLG. C ALIM=TRANSMISSION LOSS AT WHICH A RAY IS DROPPED FROM PLOT. C READ( PLTFIL, * ) R1KM, DRKM, R2KM, ZMAX, XLG, ALIM C C***** NSRS AND NBRS ARE SURFACE AND BOTTOM-HIT LIMITS AT WHICH C TO DROP A RAY FROM THE PLOT. NRAYS=NUMBER OF RAYS TO PLOT. C READ( PLTFIL, * ) NSRS, NBRS, NRAYS C R1 = 1000.*R1KM DR = 1000.*DRKM R2 = 1000.*R2KM IF(ALIM.EQ.0.) ALIM = 200. WRITE(6,903) R1KM, R2KM, ZMAX, XLG, YHGT WRITE(6,904) ALIM, NSRS, NBRS SSMIN = 10.**(-.1*ALIM) IF(XLG.LE.0..OR.XLG.GT.120.) STOP IF(ZMAX.LE.0.) STOP Z1 = 0. DZ = 1000. IF(ZMAX.LT.4000.) DZ = 500. IF(ZMAX.LT.2000.) DZ = 250. NDZ = (ZMAX + DZ - 5.)/DZ Z2 = NDZ*DZ C IF(NRAYS.EQ.0) GO TO 5 READ( PLTFIL, * ) (NRPLOT(I),I=1,NRAYS) C GO TO 26 C5 IL = 1 C6 IH = IL + 15 C READ( PLTFIL, 901, END=10 ) (NRPLOT(I),I=IL,IH) C GO TO 15 C10 IH=IL-1 C GO TO 25 C 15 WRITE(6,905) (NRPLOT(I),I=IL,IH) C IF(NRPLOT(IH).EQ.0) GO TO 20 C IL=IL+15 C GO TO 6 C20 IH=IH-1 C IF(NRPLOT(IH).NE.0) GO TO 25 C IF(IH.GT.0) GO TO 20 C25 NRAYS=IH 26 IF(NRAYS.LE.0 .OR. NRAYS.GT.64) STOP LBLK = 36000/NRAYS NRS = LBLK MBLK=(NRS-1)/LBLK+1 WRITE(6,906) NRAYS, MBLK, LBLK IF(MBLK.GT.1) STOP NR=NRS M=1 IREC=0 ZBMIN = 10000. DO 500 IR=1,LBLK IREC=IREC+1 IFEF=0 33 READ( RAYFIL, END=36) NRAY,R(IREC), ZB(IREC), (TR(I),I=1,NRAY), 1 (ZR(I), I=1, NRAY), (SS(I), I=1, NRAY), (TIME(I), I=1, NRAY), 2 (PHASE(I),I=1,NRAY), (NBR(I),I=1,NRAY), (NTU(I),I=1,NRAY), 3 (NSR(I),I=1,NRAY), (NTO(I),I=1,NRAY) RB(IREC) = R(IREC) IF(IREC.LT.5) TYPE *, R(IREC), ZB(IREC), ZR(3) GO TO 40 36 IFEF=1 IF(NRS.LT.IREC-1) GO TO 501 37 NRS = IREC-1 NR=NRS GO TO 501 40 IF(R(IREC).LT.R1) GO TO 33 IF(IREC.EQ.2 .AND. R(1).EQ.R(2)) GO TO 33 IF(R(IREC).GT.R2) GO TO 37 ZBMIN = AMIN1(ZBMIN,ZB(IREC)) DO 200 I=1,NRAYS K=NRPLOT(I) C ------ CHECK THAT THE RAY ACTUALLY EXISTS IF ( K .GT. NRAY ) THEN TYPE *, '*** WARNING: RAY ',K,' DOES NOT EXIST ***' NRAYS = I-1 GOTO 500 ENDIF J=LBLK*(I-1)+IREC ZSTOR(J) = -1. IF(SS(K).LE.SSMIN .OR. NSR(K).GT.NSRS .OR. & NBR(K).GT.NBRS) GO TO 200 ZSTOR(J) = ZR(K) NBHIT(J) = NBR(K) NSHIT(J) = NSR(K) TSTOR(J) = -TR(K) 200 CONTINUE 500 CONTINUE 501 RMAX=R(NR) C C***** THESE DISSPLA CALLS DRAW THE AXES AND THE BOTTOM PROFILE. C CALL COMPRS CALL COMPLX C CALL SHDCHR(0.,1,.005,1) C CALL BASALF('STANDARD') CALL NOBRDR CALL PAGE(XLG+4., YHGT+3.) CALL INTAXS CALL YAXANG(0.) CALL XTICKS(2) CALL YTICKS(2) CALL XNONUM CALL HEIGHT(.25) CALL TITLE( 'TRIMAIN RAY PATHS$', 100, 'Range (km)$', 100, & 'Depth (m)$', 100, XLG, YHGT ) CALL FRAME CALL GRAF(R1KM, DRKM, R2KM, Z2, -DZ, Z1) CALL GRACE(0.) CALL NOCHEK IF(ZBMIN.GT.ZMAX) GO TO 550 CALL THKCRV(.02) C ------ CONVERT TO KM DO 507 IR = 1, NR RBKM(IR) = RB(IR)/1000.0 507 CONTINUE CALL CURVE(RB, ZB, NR, 0) CALL RESET('THKCRV') C C***** PLOT ONE RAY PATH EACH PASS THROUGH LOOP. C 550 CONTINUE DO 1000 I=1,NRAYS L=(I-1)*LBLK+1 L1 = L + NRS DO 800 J=1,NRS JJ = L1 - J IF(ZSTOR(JJ).LT.0) GO TO 800 NPTS = JJ GO TO 801 800 CONTINUE 801 NBHITS = NBHIT(NPTS) NSHITS = NSHIT(NPTS) NPTS = NPTS - L + 1 CALL RAYPLT(R, NPTS, ZSTOR(L), TSTOR(L), NBHIT(L), NSHIT(L)) 1000 CONTINUE CALL RESET('XNONUM') CALL ENDPL(0) CALL DONEPL C STOP END SUBROUTINE RAYPLT(R, N, Z, T, NBHIT, NSHIT) COMMON /BOTTOM/ RB(8000), ZB(8000) COMMON /HITS/ NBHITS, NSHITS DIMENSION R(8000), Z(2000), T(2000) INTEGER *2 NBHIT(2000), NSHIT(2000) REAL RKM(8000) C C***** SUBROUTINE TO PLOT THE TRAJECTORY OF A SINGLE RAY. C THE N POINTS ON THE PATH ARE STORED IN ARRAYS R=RANGE(KM), C Z=DEPTH(M), T=TANGENT W/HORIZ(+DOWN). RANGES EQUALLY SPACED. C NBHIT,NSHIT=COUNTS OF BOTTOM AND SURFACE HITS AT EACH RANGE. C THE BOTTOM CONTOUR IS IN RB(KM) AND ZB(M); R(I)=RB(I) AT ENTRY. C THIS ROUTINE COMPUTES APPROXIMATE LOCATION OF EACH SURFACE OR C BOTTOM HIT IN ORDER TO PLOT A SHARP REFLECTION, WITH SMOOTH C SEGMENT BETWEEN REFLECTIONS. C N1 = 1 C C***** IF THE RAY NEVER HITS BOUNDARY, INTERPOLATE IT ALL IN ONE PIECE. C IF(NBHITS.GT.0 .OR. NSHITS.GT.0) GO TO 50 NPT = N-1 GO TO 1001 C C***** OTHERWISE GO THROUGH THE LOOP TO FIND THE END OF THE SEGMENT. C THE DIRECTION OF THE RAY CHANGES CONTINUOUSLY WITHIN A SEGMENT. C 50 N11 = N1 + 1 NPT = 0 DO 1000 I=N11,N NPT = NPT + 1 DR = R(I) - R(I-1) IF(NSHIT(I).EQ.NSHIT(I-1)) GO TO 200 C C***** SURFACE HIT C Z0 = (DR*T(I)*T(I-1) - Z(I)*T(I-1) + Z(I-1)*T(I))/ 1 (T(I) - T(I-1)) RHIT = R(I) + (Z0-Z(I))/T(I) ZHIT = 0. ZSAVE = Z(I) TSAVE = T(I) R(I) = RHIT Z(I) = ZHIT T(I) = 0. NSHIT(I-1) = NSHIT(I) GO TO 1001 200 IF(NBHIT(I).EQ.NBHIT(I-1)) GO TO 1000 C C***** BOTTOM HIT C TB = (ZB(I)-ZB(I-1))/DR IF(TB.NE.0.) GO TO 210 C C***** PROJECT INCOMING & OUTGOING RAYS TO THEIR INTERSECTION AND C USE THE INTERSECTION RANGE AS THE RANGE OF THE BOTTOM HIT. C Z0 = (DR*T(I)*T(I-1) - Z(I)*T(I-1) + Z(I-1)*T(I))/ 1 (T(I) - T(I-1)) DRHIT = (Z(I)-Z0)/T(I) GO TO 220 C C***** PROJECT INCOMING & OUTGOING RAYS ONTO THE BOTTOM AND USE THE C MEAN OF THE INTERSECTION RANGES AS THE RANGE OF THE BOTTOM HIT. C 210 RH1 = (ZB(I-1) - Z(I-1))/(T(I-1) - TB) RH2 = (ZB(I) - Z(I))/(TB - T(I)) DRHIT = .5*(DR + RH2 - RH1) 220 ZHIT = ZB(I) - TB*DRHIT RHIT = R(I) - DRHIT ZSAVE = Z(I) TSAVE = T(I) R(I) = RHIT Z(I) = ZHIT T(I) = 0. NBHIT(I-1) = NBHIT(I) GO TO 1001 1000 CONTINUE C C***** NPTS IS NUMBER OF POINTS ON THE SEGMENT TO BE PLOTTED. C NOTE THAT AN EXTRA POINT IS INSERTED AT EACH BOUNDARY HIT, AND C THE ORIGINAL COORDINATES ARE RESTORED TO R AND Z FOR NEXT SEGMENT. C 1001 NPTS = NPT + 1 C ------ CONVERT TO KM DO 901 IR = N1, N1+NPTS-1 RKM(IR) = R(IR)/1000.0 901 CONTINUE CALL CURVE( RKM(N1), Z(N1), NPTS, 0 ) R(N1) = RB(N1) N1 = N1 + NPT IF(N1.GE.N .AND. T(N1).NE.0.) GO TO 1010 C C***** IF LAST POINT OF RAY PATH HAS NOT BEEN PLOTTED, THERE MUST C HAVE BEEN A SURFACE OR BOTTOM HIT CAUSING EXIT FROM LOOP. C N2 = N1 N1 = N2 - 1 R(N1) = RHIT Z(N1) = ZHIT T(N1) = 0. R(N2) = RB(N2) Z(N2) = ZSAVE T(N2) = TSAVE GO TO 50 1010 N1 = MIN0(N1,N) R(N1) = RB(N1) C RETURN END