C
C PROGRAM S Y N T P L C ********************** C C PROGRAM SYNTPL IS DESIGNED FOR THE COMPUTATION OF RAY C SYNTHETIC SEISMOGRAMS FROM THE IMPULSE SYNTHETIC SEISMOGRAMS C STORED IN THE FILE LU2 GENERATED BY THE PROGRAM SEIS C C ************************************************************** C CHARACTER*80 MTEXT,FILEIN,FILEOU,FILE1,FILE2,FILE3 DIMENSION IDIST(2000),TIME(2000),A(2000),PH(2000), 1SEIS(3001),TAS(2000),DIST(100),SMAX(100),ISEL(100),XX(100), 2IEP(100),JS(20) COMMON/SOUR/ROS,VPS,VSS C C *** C LIN=5 LOU=6 LU2=1 LU3=2 LU4=3 FILEIN='synt.dat' FILEOU='synt.out' FILE1='lu2.dat' FILE2='lu3.dat' FILE3='lu4.dat' WRITE(*,'(2A)') ' (SYNTPL) SPECIFY NAMES OF INPUT AND OUTPUT', 1' FILES LIN, LOU, LU2, LU3, LU4: ' READ(*,*) FILEIN,FILEOU,FILE1,FILE2,FILE3 IF(FILE1.EQ.' ') GO TO 99 IF(FILE2.EQ.' ') LU3=0 IF(FILE3.EQ.' ') LU4=0 OPEN(LIN,FILE=FILEIN,FORM='FORMATTED',STATUS='OLD') OPEN(LOU,FILE=FILEOU,FORM='FORMATTED') OPEN(LU2,FILE=FILE1,FORM='FORMATTED',STATUS='OLD') IF(LU3.NE.0)OPEN(LU3,FILE=FILE2,FORM='FORMATTED') IF(LU4.NE.0)OPEN(LU4,FILE=FILE3,FORM='UNFORMATTED') C C *** C PI=3.14159265 SHF=0. LTAPE=0 NEXT=0 1 READ(lin,*)NEXT,MPRINT IF(NEXT.EQ.0)GO TO 99 WRITE(lou,100)NEXT,MPRINT C C *** C REWIND LU2 REWIND LU4 IF(LU3.NE.0)REWIND LU3 2 IF(LTAPE.EQ.0)READ(LU2,101)MTEXT IF(LTAPE.EQ.1)READ(LU2)MTEXT IF(LTAPE.EQ.0)READ(LU2,100)NDST,KSH,itpr IF(LTAPE.EQ.1)READ(LU2)NDST,KSH,itpr IF(LTAPE.EQ.0)READ(LU2,102)XSOUR,ZSOUR,TSOUR,RSTEP,ROS,VPS,VSS IF(LTAPE.EQ.1)READ(LU2)XSOUR,ZSOUR,TSOUR,RSTEP,ROS,VPS,VSS IF(LTAPE.EQ.0)READ(LU2,102)(DIST(K),K=1,NDST) IF(LTAPE.EQ.1)READ(LU2)(DIST(K),K=1,NDST) C WRITE(lou,101)MTEXT WRITE(lou,100)NDST,ksh,itpr WRITE(lou,102)XSOUR,ZSOUR,TSOUR,RSTEP,ROS,VPS,VSS WRITE(lou,102)(DIST(K),K=1,NDST) C READ(lin,*)MCOMP,MRED,MSELEC,MEPIC,MSHIFT,KABS,MSOUR WRITE(lou,100)MCOMP,MRED,MSELEC,MEPIC,MSHIFT,KABS,MSOUR if(itpr.eq.1)mred=0 IF(MCOMP.EQ.2.AND.KSH.EQ.0)WRITE(lou,115) IF(MCOMP.EQ.2.AND.KSH.EQ.0)GO TO 99 IF(MCOMP.NE.2.AND.KSH.EQ.1)WRITE(lou,117) IF(MCOMP.NE.2.AND.KSH.EQ.1)GO TO 99 READ(lin,*)TMIN,DT,TMAX,FREQ,GAMA,PSI,VRED,SHIFT WRITE(lou,102)TMIN,DT,TMAX,FREQ,GAMA,PSI,VRED,SHIFT IF(MEPIC.EQ.0)GO TO 11 READ(lin,*)NEPIC,(IEP(I),I=1,NEPIC) WRITE(lou,100)NEPIC,(IEP(I),I=1,NEPIC) 11 IF(MSELEC.EQ.0)GO TO 12 READ(lin,*)NSELEC,(ISEL(I),I=1,NSELEC) WRITE(lou,100)NSELEC,(ISEL(I),I=1,NSELEC) 12 CONTINUE IF(KABS.EQ.0)GO TO 5 READ(lin,*) FREF,QRED WRITE(lou,102)FREF,QRED IF(FREF.LT..00001)FREF=1. IF(QRED.LT..00001)QRED=1. RF=FREQ/FREF C 5 CALL SOURCE(lin,lou,0,0,MSOUR,0.,ANGLE,AMSOUR,PHSOUR) C NCD=0 I=0 3 I=I+1 IF(I.GT.2000)WRITE(lou,104) IF(I.GT.2000)GO TO 99 17 IF(LTAPE.EQ.0)READ(LU2,116,END=7) 1NC,IDIST(I),TIME(I),AX,AY,AZ,PHX,PHY,PHZ,ANGLE,TAST IF(LTAPE.EQ.1)READ(LU2,END=7) 1NC,IDIST(I),TIME(I),AX,AY,AZ,PHX,PHY,PHZ,ANGLE,TAST 18 NC1=IABS(NC) NC1=NC1+NCD IF(MSELEC.EQ.0)GO TO 14 DO 13 J=1,NSELEC IF(ISEL(J).EQ.NC1)GO TO 14 13 CONTINUE GO TO 3 14 IF(MEPIC.EQ.0)GO TO 16 DO 15 J=1,NEPIC IF(IEP(J).EQ.IDIST(I))GO TO 16 15 CONTINUE GO TO 3 16 MWAVE=3 IF(NC.GT.0.AND.MCOMP.EQ.2)GO TO 17 IF(NC.LT.0)MWAVE=1 IF(NC.LT.0.AND.MCOMP.EQ.2)MWAVE=2 CALL SOURCE(lin,lou,1,MWAVE,MSOUR,0.,ANGLE,AMSOUR,PHSOUR) IF(KSH.EQ.1)GO TO 10 IF(MCOMP.NE.0)GO TO 8 A(I)=AZ*AMSOUR PH(I)=PHZ+PHSOUR GO TO 9 8 IF(MCOMP.NE.1)GO TO 10 A(I)=AX*AMSOUR PH(I)=PHX+PHSOUR GO TO 9 10 A(I)=AY*AMSOUR PH(I)=PHY+PHSOUR 9 IF(MPRINT.GE.1)WRITE(lou,114) 1NC,IDIST(I),TIME(I),A(I),PH(I),ANGLE,TAST IF(KABS.NE.0)TAS(I)=TAST IF(KABS.EQ.2)PH(I)=PH(I)-2.*FREQ*TAST*ALOG(RF) GO TO 3 C C 7 CONTINUE NEXT=0 READ(lin,*)NEXT,LTAPE LTAPE=0 IF(NEXT.EQ.0)GO TO 4 WRITE(lou,100)NEXT,LTAPE REWIND LU2 LU2=LU REWIND LU2 NCD=NC1 GO TO 2 C C 4 ISUM=I-1 IF(VRED.LT.0.0001)VRED=8. IF(FREQ.LT.0.0001)FREQ=4. IF(GAMA.LT.0.0001)GAMA=4.0 IF(MSHIFT.EQ.0)SHF=0. IF(MSHIFT.EQ.1)SHF=.241506*GAMA/FREQ IF(MSHIFT.EQ.2)SHF=SHIFT OMEGA=2.*PI*FREQ OMRF=2.*PI*FREF OG=OMEGA/GAMA OGG=OG*OG TSH=0.45*GAMA/FREQ C C CONSTRUCTION OF SYNTHETIC SEISMOGRAMS C IF(NDST.GT.100)WRITE(lou,105) IF(NDST.GT.100)GO TO 99 NPTS=(TMAX-TMIN)/DT+1.5 IF(NPTS.GT.3001)WRITE(lou,106) IF(NPTS.GT.3001)GO TO 99 C C MMD=NDST IF(MEPIC.NE.0)MMD=NEPIC WRITE(LU4)NPTS,TMIN,DT,NDST,DIST(1),RSTEP IF(LU3.EQ.0)GO TO 6 WRITE(LU3,101)MTEXT WRITE(LU3,108)MMD,MRED,MCOMP,itpr,VRED,RSTEP,XSOUR,DT C C LOOP FOR RECEIVER POSITIONS C 6 CONTINUE DO 30 I=1,MMD JJ=I IF(MEPIC.NE.0)JJ=IEP(I) XX(I)=DIST(JJ) XXD=ABS(XX(I)-XSOUR) XXV=0. IF(MRED.NE.0)XXV=XXD/VRED AMAX=0. DO 20 J1=1,NPTS 20 SEIS(J1)=0. DO 21 K=1,ISUM IF(IDIST(K).NE.JJ)GO TO 21 TR=TIME(K)-XXV+SHF IF((TR+TSH).LT.TMIN)GO TO 21 IF((TR-TSH).GT.TMAX)GO TO 21 AMPL=A(K) PHASE=PH(K) IF(AMPL.LT.0.005*AMAX)GO TO 21 IF(AMPL.GT.AMAX)AMAX=AMPL IF1=IFIX((TR-TSH-TMIN)/DT+0.1) IF2=IFIX((TR+TSH-TMIN)/DT+0.1) IFF1=1 IF(IF1.GT.0)IFF1=IF1 IF1=IFF1 IFF2=NPTS IF(IF2.LT.NPTS)IFF2=IF2 IF2=IFF2 IF(KABS.EQ.0)GO TO 27 TAST=TAS(K)*QRED OMAS=OMEGA*(1.-OMEGA*TAST/(GAMA*GAMA)) IF(OMAS.LE.0.)WRITE(lou,112) IF(OMAS.LE.0.)GO TO 32 AABS1=OMEGA*TAST/GAMA AABS1=.25*(AABS1*AABS1) AABS2=.5*OMAS*TAST IF(KABS.NE.2)GO TO 27 AABS4=OMAS/OMRF AABS5=ALOG(AABS4) AABS3=(TAST/PI)*(1.+AABS5) 27 DO 24 KK=IF1,IF2 TT=TMIN+DT*FLOAT(KK-1) TD=TT-TR IF(KABS.EQ.2)TD=TD+AABS3 AA=PSI-PHASE BB=-OGG*TD*TD IF(KABS.EQ.0)AA=AA+OMEGA*TD IF(KABS.NE.0)AA=AA+OMAS*TD IF(KABS.EQ.2)AA=AA-OMAS*TAST/PI IF(KABS.NE.0)BB=BB-AABS1-AABS2 AA=AMPL*COS(AA)*EXP(BB) 24 SEIS(KK)=SEIS(KK)+AA 21 CONTINUE SMAX(I)=0. DO 25 KK=1,NPTS IF(SMAX(I).GT.ABS(SEIS(KK)))GO TO 25 SMAX(I)=ABS(SEIS(KK)) 25 CONTINUE C C WRITE(LU4)(SEIS(J),J=1,NPTS) 30 CONTINUE REWIND LU2 REWIND LU4 C C END OF THE LOOP FOR RECEIVERS C SMAXIM=0. DO 26 I=1,MMD IF(SMAXIM.GT.SMAX(I))GO TO 26 SMAXIM=SMAX(I) XXX=XX(I) 26 CONTINUE IF(LU3.NE.0)WRITE(LU3,113)XXX,SMAXIM WRITE(lou,113)XXX,SMAXIM IF(LU3.EQ.0)GO TO 1 C READ(LU4)NPTS,TMIN,DT,NDST,DIST(1),RSTEP DO 31 I=1,MMD READ(LU4)(SEIS(J),J=1,NPTS) C C NORMALIZATION OF SEISMOGRAMS C IF1=0 IF2=0 SMAXI=SMAX(I) IF(SMAXI.LT..000001)GO TO 35 DO 34 J=1,NPTS SEIS(J)=999.1*SEIS(J)/SMAXI IF(ABS(SEIS(J)).LT.1.)GO TO 33 IF2=0 IF(IF1.EQ.0.AND.J.EQ.1)IF1=1 IF(IF1.EQ.0)IF1=J-1 GO TO 34 33 IF(IF1.EQ.0)GO TO 34 IF(IF2.EQ.0)IF2=J 34 CONTINUE TM=TMIN+FLOAT(IF1-1)*DT IF(IF2.EQ.0)IF2=NPTS NPS=IF2-IF1+1 35 IF(IF1.EQ.0)NPS=0 IF(IF1.EQ.0)TM=TMIN IF(IF1.EQ.0)IF1=1 C WRITE(LU3,109)XX(I),SMAX(I),TM,NPS WRITE(lou,107)XX(I),SMAX(I),TM,NPS ISS=20 IS=IF1-1 41 IS1=IF2-IS IF(IS1.LT.20)ISS=IS1 DO 40 K=1,ISS 40 JS(K)=SEIS(IS+K) WRITE(LU3,111)(JS(K),K=1,ISS) IF(MPRINT.EQ.2.AND.NPS.NE.0) 1WRITE(lou,110)(JS(K),K=1,ISS) IF(IS1.LE.20)GO TO 31 IS=IS+20 GO TO 41 31 CONTINUE REWIND LU3 32 CONTINUE C C *** C *** C GO TO 1 C C 100 FORMAT(26I3) 101 FORMAT(A) 102 FORMAT(8F10.5) 103 FORMAT(2I3,F10.5,2E12.6,4F10.6) 104 FORMAT(1X,'NUMBER OF READINGS FROM LU2 GREATER THAN 2000.'/, 11X,'COMPUTATION TERMINATED'/) 105 FORMAT(1X,'NUMBER OF RECEIVER POSITIONS GREATER THAN 100.',/, 11X,'COMPUTATION TERMINATED'/) 106 FORMAT(1X,'NUMBER OF POINTS IN THE SYNTHETIC SEISMOGRAM 1GREATER THAN 3001.'/1X,'COMPUTATION TERMINATED'/) 107 FORMAT(1X,'EPICENTRAL DISTANCE =',F10.5,E15.9,1X,F10.5,I5) 108 FORMAT(4I5,4F10.5) 109 FORMAT(F10.5,E15.8,F10.5,I5) 110 FORMAT(1X,20I4) 111 FORMAT(20I4) 112 FORMAT(1X,'FREQUENCY TOO HIGH OR GAMA TOO SMALL'/ 11X,'COMPUTATION TERMINATED'/) 113 FORMAT(1X,'EPICENTRAL DISTANCE =',F10.5,1X,'SMAXIM =',E15.9) 114 FORMAT(2I3,F10.5,E12.6,3F10.6) 115 FORMAT(1X,'THERE ARE NO DATA IN THE FILE FOR SH SEISMOGRAMS' 1/1X,'COMPUTATION TERMINATED - USE MCOMP.NE.2'/) 116 FORMAT(2I3,F10.5,3E12.6,5F10.6) 117 FORMAT(1X,'THERE ARE DATA IN THE FILE ONLY FOR SH SEISMOGRAMS' 1/1X,'COMPUTATION TERMINATED - USE MCOMP.EQ.2'/) C 99 CONTINUE C C *** C *** C STOP END C C======================================================================= C INCLUDE 'source.for' C source.for C C======================================================================= C