C
C P R O G R A M B P L O T C C C PROGRAM BPLOT IS DESIGNED FOR THE PLOTTING OF SYNTHETIC C SEISMOGRAMS AND GENERATING A FILE FOR POLARPLOT C C ************************************************************ C CHARACTER*80 TEXT,PSTEXT,MPRINT,IPRINT,TITLE CHARACTER*80 FILEIN,FILEOU,FILE1,FILE2 DIMENSION SMAX(100),SEIS(2048) DIMENSION IS(2048),IEP(100) C C C************************************************** C LIN=5 LOU=6 LU8=1 LU4=2 FILEIN='blot.dat' FILEOU='bplot.out' FILE1='lu8.dat' FILE2='lu4.dat' WRITE(*,'(2A)') ' (BPLOT) SPECIFY NAMES OF INPUT AND OUTPUT', 1' FILES LIN, LOU, LU8, LU4: ' READ(*,*) FILEIN,FILEOU,FILE1,FILE2 IF(FILE1.EQ.' ') GO TO 99 IF(FILE2.EQ.' ') LU4=0 OPEN(LIN,FILE=FILEIN,FORM='FORMATTED',STATUS='OLD') OPEN(LOU,FILE=FILEOU,FORM='FORMATTED') OPEN(LU8,FILE=FILE1,FORM='FORMATTED',STATUS='OLD') IF(LU4.NE.0)OPEN(LU4,FILE=FILE2,FORM='FORMATTED') C C************************************************** C TEXT='BPLOT' PSTEXT=' ' IPR=0 READ(LIN,*)TEXT READ(LIN,*)IPR,PSTEXT WRITE(LOU,108)TEXT WRITE(LOU,104)IPR,PSTEXT IF(IPR.LT.0)THEN CALL PLOTN(PSTEXT,0) IPR=-IPR END IF CALL PLOTS(LDUM1,LDUM2,7) CALL PLOT(4.,6.,-3) C C READING FROM LU8 C IF(LU8.NE.0)REWIND LU8 IF(LU4.NE.0)REWIND LU4 READ(LU8,106)MPRINT READ(LU8,106)IPRINT READ(LU8,106)TITLE READ(LU8,152)XSOUR,YSOUR,ZSOUR,TSOUR,RSTEP,DT,DF READ(LU8,101)NDST,NT,MCOMP,ILOC WRITE(LOU,108)MPRINT WRITE(LOU,108)IPRINT WRITE(LOU,108)TITLE WRITE(LOU,152)XSOUR,YSOUR,ZSOUR,TSOUR,RSTEP,DT,DF WRITE(LOU,101)NDST,NT,MCOMP C C INPUT DATA C MRED=0 MEPIC=0 NTICX=1 NTICY=1 NDX=0 NDY=0 1 READ(LIN,*)MRED,MEPIC,NTICX,NTICY,NDX,NDY WRITE(LOU,101)MRED,MEPIC,NTICX,NTICY,NDX,NDY IF(MEPIC.EQ.0)GO TO 3 READ(LIN,*)NEPIC,(IEP(I),I=1,NEPIC) WRITE(LOU,101)NEPIC,(IEP(I),I=1,NEPIC) 3 CONTINUE READ(LIN,*)XMIN,XMAX,XLEN,DTICX,YMIN,YMAX,YLEN,DTICY WRITE(LOU,102)XMIN,XMAX,XLEN,DTICX,YMIN,YMAX,YLEN,DTICY NP=IFIX((YMAX-YMIN)/DT) NPMIN=MIN0(NP,NT) VRED=6. AMP=0. B1=1. EPICS=10. EPS=0. SC=1. READ(LIN,*)VRED,AMP,B1,EPICS,EPS,SC WRITE(LOU,102)VRED,AMP,B1,EPICS,EPS,SC SMAXIM=0. XMX=0. C C COMPUTATION OF SMAXIM C DO 40 I=1,NDST READ(LU8,162)XX,TMIN,AREDUC,NUM IF(NUM.EQ.0)GO TO 40 READ(LU8,109)(IS(M),M=1,NUM) IF(XX.LE.XMIN.OR.XX.GE.XMAX)GO TO 40 IF (MEPIC.EQ.0)GO TO 45 DO 46 J=1,NEPIC IF(I.EQ.IEP(J))GO TO 45 46 CONTINUE GO TO 40 45 CONTINUE TSTART=YMIN-TMIN IF(MRED.EQ.1)TSTART=YMIN+ABS(XX-XSOUR)/VRED-TMIN IPOM=IFIX(TSTART/DT)+1 47 IF(IPOM.GT.0)GO TO 41 IPOM=IPOM+NT GO TO 47 41 IF(IPOM.LE.NT)GO TO 42 IPOM=IPOM-NT GO TO 41 42 CONTINUE IIAUX=NT-IPOM DO 43 J=1,NT IF(J.GE.IPOM)JAUX=J-IPOM+1 IF(J.LT.IPOM)JAUX=J+IIAUX+1 43 SEIS(JAUX)=FLOAT(IS(J)) SMAX(I)=0. DO 44 J=1,NPMIN IF(ABS(SEIS(J)).GT.SMAX(I))SMAX(I)=ABS(SEIS(J)) 44 CONTINUE SMAX(I)=SMAX(I)*AREDUC/999.1 IF(SMAX(I).GT.SMAXIM)XMX=XX IF(SMAX(I).GT.SMAXIM)SMAXIM=SMAX(I) 40 CONTINUE WRITE(LOU,103)XMX,SMAXIM WRITE(LOU,102)(SMAX(I),I=1,NDST) REWIND LU8 READ(LU8,106)MPRINT READ(LU8,106)IPRINT READ(LU8,106)TITLE READ(LU8,152)XSOUR,YSOUR,ZSOUR,TSOUR,RSTEP,DT,DF READ(LU8,101)NDST,NT,MCOMP,ILOC C C END OF COMPUTATION OF SMAXIM C 2 CONTINUE IF(LU4.EQ.0)GO TO 12 WRITE(LU4,106)TEXT WRITE(LU4,100)NDST,MRED,MCOMP,ILOC,VRED,RSTEP,XSOUR,YSOUR,DT WRITE(LU4,110)XMX,SMAXIM 12 CONTINUE XMER=XLEN/(XMAX-XMIN) YMER=YLEN/(YMAX-YMIN) DDX=RSTEP*XMER ELM=.45*SC CALL BORDER(XLEN,DTICX,YLEN,DTICY,SC,TEXT,0,XMIN,XMAX, 1YMIN,YMAX,NTICX,NTICY,NDX,NDY) T=.5*(XLEN-6.3*SC) CALL SYMBOL(T,-1.6*SC,ELM,'DISTANCE IN KM',0.,14) T=.5*(YLEN-8.1*SC) U=-(1.6+.4*NDX)*SC IF(MRED.EQ.0) 1CALL SYMBOL(U,T,ELM,'TRAVEL TIME IN SEC',90.,18) IF(MRED.EQ.0)GO TO 9 CALL SYMBOL(U,T,ELM,'T-D/ ',90.,5) T=T+1.8*SC CALL NUMBER(U,T,ELM,VRED,90.,2) T=T+2.7*SC CALL SYMBOL(U,T,ELM,'(IN SEC)',90.,8) 9 CONTINUE IF(MCOMP.EQ.0)CALL SYMBOL(ELM,YLEN+SC,ELM,'VERTICAL',0.,8) IF(MCOMP.EQ.1)CALL SYMBOL(ELM,YLEN+SC,ELM,'X-COMPONENT',0.,11) IF(MCOMP.EQ.2)CALL SYMBOL(ELM,YLEN+SC,ELM,'Y-COMPONENT',0.,11) T=XLEN-7.5*SC CALL NUMBER(T,YLEN+.5*SC,.3*SC,AMP,0.,0) T=T+1.5*SC CALL NUMBER(T,YLEN+.5*SC,.3*SC,B1,0.,2) T=T+1.5*SC CALL NUMBER(T,YLEN+.5*SC,.3*SC,EPS,0.,1) T=T+1.5*SC CALL NUMBER(T,YLEN+.5*SC,.3*SC,SMAXIM,0.,5) C C LOOP FOR THE RECEIVER POSITIONS C 4 DO 10 I=1,NDST READ(LU8,162)XX,TMIN,AREDUC,NUM IF(IPR.GT.1)WRITE(LOU,162)XX,TMIN,AREDUC,NUM IF(NUM.EQ.0)GO TO 10 READ(LU8,109)(IS(M),M=1,NUM) IF(LU4.NE.0)WRITE(LU4,107)XX,AREDUC,TMIN,NUM IF(LU4.NE.0)WRITE(LU4,109)(IS(M),M=1,NUM) IF(XX.LE.(XMIN+0.001).OR.XX.GE.(XMAX-0.001))GO TO 10 IF(MEPIC.EQ.0)GO TO 5 DO 6 J=1,NEPIC IF(I.EQ.IEP(J))GO TO 5 6 CONTINUE GO TO 10 C C SHIFTING SEISMOGRAM INTO A REQUIRED TIME POSITION C 5 TSTART=YMIN-TMIN IF(MRED.EQ.1)TSTART=YMIN+ABS(XX-XSOUR)/VRED-TMIN IPOM=IFIX(TSTART/DT)+1 TL=TMIN+DT*FLOAT(IPOM-1) 37 IF(IPOM.GT.0)GO TO 31 IPOM=IPOM+NT GO TO 37 31 IF(IPOM.LE.NT)GO TO 32 IPOM=IPOM-NT GO TO 31 32 CONTINUE IIAUX=NT-IPOM DO 33 J=1,NT IF(J.GE.IPOM)JAUX=J-IPOM+1 IF(J.LT.IPOM)JAUX=J+IIAUX+1 33 SEIS(JAUX)=AREDUC*FLOAT(IS(J))/999.1 SMAX(I)=0. DO 39 J=1,NPMIN IF(ABS(SEIS(J)).GT.SMAX(I))SMAX(I)=ABS(SEIS(J)) 39 CONTINUE SMAXI=SMAX(I) IF(I.EQ.1)SMAX1=SMAXI C C SYNTHETIC SEISMOGRAM SEIS(J),J=1,NPMIN, CORRESPONDS TO TRAVEL C TIMES (OR REDUCED TRAVEL TIMES) FROM YMIN TO YMAX C C C COMPUTATION OF SCALING FACTORS C IF(SMAXI.LT.0.000001)GO TO 7 IF(ABS(AMP).LT.0.00001)FACTOR=B1*DDX/SMAXI IF(ABS(AMP).LT.0.00001)GO TO 8 7 FACTOR=0. IF(ABS(EPS).GT.0.00001)GO TO 20 IF(AMP.LT.(-0.00001))FACTOR=B1*DDX/SMAXIM IF(AMP.GT.0.00001.AND.AMP.LT.5.)FACTOR=B1 IF(AMP.GT.5.)FACTOR=B1*DDX/SMAX1 GO TO 8 20 IF(AMP.LT.0.00001)FACTOR=B1*DDX*((ABS(XX-XSOUR)/EPICS)**EPS) 1/SMAXIM IF(AMP.GT.0.00001)FACTOR=B1*(ABS(XX-XSOUR)/EPICS)**EPS 8 CONTINUE SFMAX=FACTOR*SMAXI IF(IPR.EQ.1)WRITE(LOU,120)XX,SMAXI,FACTOR,SFMAX IF(IPR.LE.1)GO TO 25 WRITE(LOU,121)XX,TSTART,SMAXI,FACTOR,SFMAX DO 26 J=1,NPMIN 26 IS(J)=IFIX(999.1*SEIS(J)/SMAXI) WRITE(LOU,109)(IS(J),J=1,NPMIN) 25 CONTINUE C C PLOTTING C XM=XMIN X0=(XX-XM)*XMER XNEW=X0 YNEW=0. CALL PLOT(XNEW,YNEW,3) IOLD=1 AMPL=SEIS(1)+(SEIS(2)-SEIS(1))*(TL-YMIN)/DT XNEW=X0-FACTOR*AMPL IF(XNEW.LT.0..OR.XNEW.GT.XLEN)GO TO 15 CALL PLOT(XNEW,YNEW,2) 15 CONTINUE DO 11 J=2,NPMIN INEW=0 U1=FACTOR*SEIS(J) IF(J.EQ.NPMIN)GO TO 50 IF(ABS(U1).LT.0.003*SFMAX)INEW=1 IFUT=0 U11=FACTOR*SEIS(J+1) IF(ABS(U11).LT.0.003*SFMAX)IFUT=1 IF(IOLD.EQ.1.AND.INEW.EQ.1.AND.IFUT.EQ.1)GO TO 11 50 XNEW=X0-U1 YNEW=(TL-YMIN+DT*FLOAT(J-1))*YMER CALL PLOT(XNEW,YNEW,2) 11 IOLD=INEW AMPL=SEIS(NPMIN)+(SEIS(NPMIN+1)-SEIS(NPMIN))*(TL-YMIN)/DT XNEW=X0-FACTOR*AMPL IF(XNEW.LT.0..OR.XNEW.GT.XLEN)GO TO 10 CALL PLOT(XNEW,YLEN,2) 10 CONTINUE C C END OF THE LOOP FOR RECEIVER POSITIONS C 100 FORMAT(4I5,5F10.5) 101 FORMAT(16I5) 102 FORMAT(8F10.5) 103 FORMAT(/,'SMAXIM=',F15.5,' AT XMX=' ,F10.5/) 104 FORMAT(I5,1X,A) 106 FORMAT(A) 107 FORMAT(F10.5,E15.8,F10.5,I5) 108 FORMAT(1X,A) 109 FORMAT(20I4) 110 FORMAT(22X,F10.5,9X,F15.9) 111 FORMAT(20X,F15.9) 112 FORMAT(6F15.9) 120 FORMAT(F10.5,3E15.4) 121 FORMAT(/1X,'SYNTHETIC SEISMOGRAM',2X,'X=',F10.5,2X,'T0=',F10.5, 1/1X,'SMAX=',1E15.6,2X,'FACTOR=',1E15.6,2X,'SFMAX=',F10.5) 152 FORMAT(5F10.5,2E15.7) 154 FORMAT(3E16.8) 162 FORMAT(2F10.3,1E12.5,I5) 99 CALL PLOT(0.,0.,999) STOP END C C======================================================================= C INCLUDE 'border.for' C border.for INCLUDE 'error.for' C error.for INCLUDE 'length.for' C length.for INCLUDE 'sep.for' C sep.for INCLUDE 'calcops.for' C calcops.for C C======================================================================= C