C PROGRAM S Y N T P L C ********************* C C PROGRAM SYNTPL IS DESIGNED FOR THE COMPUTATION OF RAY-SYNTHETIC C SEISMOGRAMS FROM THE DATA STORED IN THE FILE LU2 GENERATED IN C THE PROGRAM SEIS88. THE FILE LU2 CONTAINS TRAVEL TIMES, AMPLITU- C DES AND PHASE SHIFTS FOR ALL SPECIFIED RECEIVERS AND ALL ELEMEN- C TARY WAVES. THE FILE ALSO CONTAINS THE INITIAL ANGLES OF CORRES- C PONDING RAYS AND CORRESPONDING GLOBAL ABSORPTION FACTORS. PROGRAM C SYNTPL GENERATES A NEW FILE LU3 AND AN AUXILIARY FILE LU4 WITH C SYNTHETIC SEISMOGRAMS. THESE SYNTHETIC SEISMOGRAMS CAN BE PLOTTED C BY STANDARD PLOTTING ROUTINES, AVAILABLE AT MANY INSTITUTIONS. C FOR PLOTTING OF SEISMOGRAMS, IT IS ALSO POSSIBLE TO USE THE PRO- C GRAM SEISPLOT INCLUDED IN THIS PACKAGE. C THE SOURCE-TIME FUNCTION USED IN SYNTPL IS A HARMONIC CARRIER C MODULATED BY A GAUSSIAN ENVELOPE, C F(T)=EXP(-(OMEGA*T/GAMA)**2)COS(OMEGA*T+PSI), (1) C WHERE T IS TIME, OMEGA=2.*PI*FREQ, AND FREQ, GAMA AND PSI ARE C THREE FREE PARAMETERS OF THE SOURCE-TIME FUNCTION. BY A PROPER C CHOICE OF THESE PARAMETERS, THE USER CAN SIMULATE A BROAD VARIETY C OF WAVELETS WHICH ARE CLOSE TO THOSE OBSERVED IN SEISMOLOGY AND C SEISMIC PROSPECTING (SEE DETAILS IN CERVENY,MOLOTKOV AND PSENCIK, C RAY METHOD IN SEISMOLOGY, UNIVERZITA KARLOVA, PRAGUE 1977). C NOTE THAT THE ZERO TIME OF THE SIGNAL (1) CORRESPONDS TO THE C MAXIMUM OF ITS ENVELOPE, NOT TO ITS FIRST ARRIVAL. THIS GIVES C A TIME SHIFT OF THIS SIGNAL IN SYNTHETIC SEISMOGRAMS IN COMPARI- C SON WITH OTHER SIMILAR SIGNALS (MUELLER SIGNAL,ETC.). THE USER C WHO IS NOT SATISFIED WITH THIS SHIFT CAN USE THE PARAMETERS C MSHIFT AND SHIFT READ FROM INPUT DATA TO REMOVE IT. C WHEN COMPUTING THE SYNTHETIC SEISMOGRAMS FROM THE DATA STORED C IN THE FILE LU2, IT IS POSSIBLE TO SELECT ONLY CERTAIN DATA C FROM LU2. C FIRST SELECTION IS PERFORMED USING THE INTERNAL WAVE CODE OF C ELEMENTARY WAVES, WHICH CORRESPONDS TO THE SUCCESSIVE NUMBER OF C THE ELEMENTARY WAVE AS IT IS GENERATED BY THE ROUTINE GENER IN C THE PROGRAM SEIS88, SEE DESCRIPTION OF THIS ROUTINE FOR DETAILS. C IN THIS WAY, IT IS POSSIBLE TO SELECT, E.G., ONLY PP REFLECTIONS, C OR PS REFLECTIONS, OR TO REMOVE THE DIRECT WAVE, ETC. C THE SECOND SELECTION APPLIES TO EPICENTRAL DISTANCES. IT IS C POSSIBLE TO CONSTRUCT THE SYNTHETIC SEISMOGRAMS ONLY FOR SELEC- C TED RECEIVERS,AND NOT FOR THE WHOLE SYSTEM OF RECEIVERS FOR C WHICH THE DATA IN SEIS88 WERE COMPUTED. C C INPUT DATA C ********** C C INPUT DATA CONSIST PARTIALLY OF THE DATA OBTAINED IN THE PROGRAM C SEIS88, WHICH ARE STORED IN FORMATTED OR BINARY FORM IN THE FILE C LU2, AND PARTIALLY OF THE ADDITIONAL DATA CONTROLLING THE COM- C PUTATIONS, INTRODUCED DIRECTLY (CARD READER, TERMINAL). C C THE DATA STORED AT LU2 C ********************** C C THE DATA MAY BE STORED IN LU2 IN A FORMATTED FORM (THEN SPE- C CIFY LTAPE=0 IN CARD NO.1) OR IN A BINARY FORM (SPECIFY C LTAPE=1). FOR DETAIL SEE THE DESCRIPTION OF THE CONTENT OF C THE FILE LU2 IN PROGRAM SEIS88. C C 1) MTEXT FORMAT(17A4) C 2) NDST,KSH,ITPR FORMAT(26I3) C 3) XSOUR,ZSOUR,TSOUR,RSTEP,ROS,VPS,VSS FORMAT(8F10.5) C 4) DST(I),I=1,NDST FORMAT(8F10.5) C 5) NCODE,II,T,AX,AY,AZ,PHX,PHY,PHZ,ANGLE,TAST C FORMAT(2I3,F10.5,3E12.6,5F10.6) C C THE DATA IN 5 CAN BE REPEATED ON LU2 MANY TIMES. IN THE PROGRAM C SYNTPL, THIS NUMBER MAY NOT EXCEED 2000. OTHERWISE IT WOULD BE C NECESSARY TO CHANGE THE DIMENSIONS IN THE PROGRAM. THE DATA 1-5 C CAN BE READ IN FROM DIFFERENT FILES (SEE INPUT DATA CARD NO.7) C GENERATED IN THE PROGRAM SEIS88 FOR A SPECIFIC MODEL, SOURCE AND C RECEIVER POSITIONS. NOTE THAT THIS POSSIBILITY CAN BE USED ONLY C IN OPERATING SYSTEMS ALLOWING SPECIFICATION OF THE NUMBER OF A C FILE OUTSIDE THE SOURCE PROGRAM (E.G. IBM OPERATING SYSTEM). C NOTE THAT THE MAXIMUM NUMBER OF THE RECEIVER POSITIONS IS 100. C C C THE ADDITIONAL INPUT DATA C ************************* C C THESE DATA SHOULD BE INTRODUCED DIRECTLY (CARDS,TERMINAL). THEY C CONTROL THE COMPUTATION OF SYNTHETIC SEISMOGRAMS. IN THE DESCRIP- C TION OF ADDITIONAL INPUT DATA, WE SHALL DENOTE PLACES WHERE THE C DATA FROM LU2 ARE READ IN BY **LU2/1, **LU2/2, ETC. C C 1) ONE CARD C LU2,LU3,MPRINT,LU4,LTAPE FORMAT(26I3) C LU2... NUMBER OF THE FILE IN WHICH THE IMPULSE SYNTHETIC C SEISMOGRAMS COMPUTED IN THE PROGRAM SEIS88 ARE C STORED. IF LU2 IS NOT SPECIFIED, THEN THE COMPUTA- C TION TERMINATES. C LU3... NUMBER OF THE FILE IN WHICH THE SYNTHETIC SEISMO- C GRAMS COMPUTED IN THIS PROGRAM ARE STORED. C MPRINT... THE PARAMETER WHICH CONTROLS THE PRINTOUT ON THE C LINE PRINTER, SEE DETAILS IN THE SECTION "OUTPUT C TABLES". C LU4... NUMBER OF AN AUXILIARY FILE USED IN SYNTPL FOR C STORAGE OF SYNTHETIC SEISMOGRAMS IN UNFORMATTED C FORM FOR PLOTTING OF SEISMOGRAMS ON SPECIAL SEIS- C MIC PLOTTERS. IF LU4 IS NOT SPECIFIED, THEN LU4=10. C C LTAPE... LTAPE=0 FOR DATA STORED IN LU2 IN A FORMATTED C FORM. C LTAPE=1 FOR DATA STORED IN LU2 IN A BINARY FORM. C C **LU2/1 C **LU2/2 C **LU2/3 C **LU2/4 C C 2) ONE CARD, VARIOUS SWITCHES. C MCOMP,MRED,MSELEC,MEPIC,MSHIFT,KABS,MTYPE FORMAT(26I3) C MCOMP... CONTROLS THE CHOICE OF THE COMPONENT OF THE C DISPLACEMENT VECTOR. C MCOMP=0: VERTICAL COMPONENT. C MCOMP=1: RADIAL COMPONENT. C MCOMP=2: TRANSVERSE COMPONENT. C MRED... CONTROLS REDUCTION OF TRAVEL TIMES. IF VSP C SYNTHETICS ARE COMPUTED, MRED SET AUTOMATI- C CALLY, MRED=0. C MRED=0: TRAVEL TIMES NOT REDUCED. C MRED=1: TRAVEL TIMES REDUCED. THE REDUCTION C VELOCITY IS GIVEN IN THE CARD NO.3. C MSELEC... CONTROLS SELECTION OF ELEMENTARY WAVES. C MSELEC=0: ALL WAVES STORED AT LU2 ARE USED TO C CONSTRUCT SYNTHETIC SEISMOGRAMS. C MSELEC=1: ONLY SELECTED WAVES ARE USED TO CON- C STRUCT SYNTHETIC SEISMOGRAMS, SEE CARD NO.5. C MEPIC... CONTROLS SELECTION OF RECEIVERS. C MEPIC=0: THE SYNTHETIC SEISMOGRAMS ARE CONSTRUC- C TED AT ALL RECEIVER POSITIONS. C MEPIC=1: THE SYNTHETIC SEISMOGRAMS ARE CONSTRUC- C TED ONLY AT SELECTED RECEIVER POSITIONS, SEE C CARD NO.4. C MSHIFT... CONTROLS THE TIME SHIFT OF THE SIGNAL. C MSHIFT=0: NO TIME SHIFT. ZERO TIME CORRESPONDS C TO THE MAXIMUM OF THE ENVELOPE OF THE SIGNAL (1). C MSHIFT=1: THE TIME SHIFT IS INTRODUCED AUTOMATI- C CALLY. IT CORRESPONDS TO THE TIME INTERVAL BET- C WEEN THE MAXIMUM OF THE ENVELOPE AND THE TIME C WHERE THE ENVELOPE IS 0.1 FROM ITS MAXIMUM VALUE. C MSHIFT=2: THE TIME SHIFT IS INTRODUCED MANUALLY, C BY THE PARAMETER 'SHIFT' IN CARD NO.3. C KABS... SPECIFIES THE DISSIPATION TYPE. C KABS=0: PERFECTLY ELASTIC MEDIUM. C KABS=1: SLIGHTLY ABSORBING MEDIUM WITH NON-CAUSAL C ABSORPTION. C KABS=2: SLIGHTLY ABSORBING MEDIUM WITH CAUSAL C ABSORPTION (FUTTERMAN'S MODEL). C MTYPE... CONTROLS THE TYPE OF THE POINT SOURCE. C MTYPE=0: UNIT ISOTROPIC RADIATION PATTERN. C MTYPE=1: SINGLE FORCE. C MTYPE=2: DOUBLE COUPLE. C MTYPE=3: EXPLOSIVE (IMPLOSIVE) SOURCE. C NOTE: SINGLE FORCE AND DOUBLE COUPLE POINT SOUR- C CES CAN BE ARBITRARILY ORIENTED IN THE 3-D SPACE C (THEY DO NOT HAVE TO BE CONFINED TO THE VERTICAL C (X,Z) PLANE). C C 3) ONE CARD, SPECIFICATION OF SYNTHETIC SEISMOGRAMS. C TMIN,DT,TMAX,FREQ,GAMA,PSI,VRED,SHIFT FORMAT(8F10.5) C TMIN,DT,TMAX... MINIMUM TIME, TIME STEP AND THE MAXIMUM C TIME OF COMPUTED SYNTHETIC SEISMOGRAMS. IN CASE C OF MRED=1, THESE QUANTITIES CORRESPOND TO THE C REDUCED TRAVEL TIME. C FREQ,GAMA,PSI... PARAMETERS WHICH SPECIFY THE SOURCE-TIME C FUNCTION (1). C WHEN FREQ NOT SPECIFIED, FREQ=4. C WHEN GAMA NOT SPECIFIED, GAMA=4. C VRED... REDUCTION VELOCITY. C WHEN VRED NOT SPECIFIED, VRED=8. C SHIFT... WHEN MSHIFT=2, SHIFT GIVES THE TIME SHIFT OF THE C SIGNAL (1). OTHERWISE IT MAY BE ARBITRARY. C C 4) ONE CARD, SELECTION OF RECEIVER POSITIONS. C INCLUDED ONLY WHEN MEPIC.NE.0. C NEPIC,(IEP(I),I=1,NEPIC) FORMAT(26I3) C NEPIC...NUMBER OF SELECTED RECEIVER POSITIONS AT WHICH C THE SYNTHETIC SEISMOGRAMS SHOULD BE CONSTRUCTED. C IEP(1),IEP(2),...,IEP(NEPIC)... SEQUENTIAL NUMBERS OF C SELECTED RECEIVER POSITIONS. THE RECEIVER POSI- C TIONS ARE NUMBERED FROM THE LEFT TO THE RIGHT. C C 5) ONE CARD, SELECTION OF WAVES C INCLUDED ONLY WHEN MSELEC.NE.0. C NSELEC,(ISEL(I),I=1,NSELEC) FORMAT(26I3) C NSELEC...NUMBER OF SELECTED WAVES USED FOR THE CONSTRUC- C TION OF SYNTHETIC SEISMOGRAMS. C ISEL(1),ISEL(2),...,ISEL(NSELEC)... INTERNAL CODES OF C SELECTED WAVES. THEY CORRESPOND TO THE SEQUENTIAL C NUMBER OF THE WAVE AS IT IS GENERATED BY THE C ROUTINE GENER IN SEIS88. FOR DETAILS, SEE C DESCRIPTION OF ROUTINE GENER. C C 6) ONE CARD, REFERENCE FREQUENCY. C INCLUDED ONLY IF KABS.NE.0. C FREF,QRED FORMAT(8F10.5) C FREF... REFERENCE FREQUENCY FOR ABSORPTION COMPUTATIONS. C THE VELOCITY VALUES AND QUALITY FACTOR USED IN C THE MODEL ARE SPECIFIED FOR THE FREQUENCY FREF. C DEFAULT VALUE, FREF=1. C QRED... FACTOR BY WHICH GLOBAL ABSORPTION FACTORS "TAST" C ALONG ALL RAYS ARE MULTIPLIED. IT CORRESPONDS TO C DIVISION OF THE CORRESPONDING Q-FACTOR BY QRED. C DEFAULT VALUE, QRED=1. C C 7) ONE CARD, SPECIFICATION OF PARAMETERS OF THE SOURCE. C INCLUDED ONLY IF MTYPE.NE.0. C IPAR(1),...,IPAR(4),PAR(1),...,PAR(6) FORMAT(4I5,6F10.5) C C MTYPE=1: SINGLE FORCE. C C IPAR(1)-IPAR(4)... FREE PARAMETERS. C PAR(1)... FORCE AZIMUTH - THE ANGLE, IN RADIANS, MADE BY C THE PROJECTION OF THE FORCE VECTOR INTO A HORI- C ZONTAL PLANE AND THE POSITIVE X-AXIS. POSITIVE C CLOCKWISE. C PAR(2)... MAGNITUDE OF THE FORCE. C PAR(3)... FORCE DECLINATION - THE ANGLE, IN RADIANS, MADE C BY THE FORCE VECTOR AND A HORIZONTAL PLANE. C PAR(4)-PAR(6)... FREE PARAMETERS. C C MTYPE=2: DOUBLE COUPLE. C C IPAR(1)-IPAR(4)... FREE PARAMETERS. C PAR(1)... DIP ANGLE, IN RADIANS. C PAR(2)... SOURCE MOMENT. C PAR(3)... STRIKE ANGLE, IN RADIANS. C PAR(4)... RAKE ANGLE, IN RADIANS. C PAR(5)-PAR(6)... FREE PARAMETERS. C NOTE: THE FORMULAE FOR THE DOUBLE COUPLE FOLLOW THE DE- C FINITION IN AKI AND RICHARDS (1980). C C MTYPE=3: EXPLOSIVE (IMPLOSIVE) SOURCE. C C IPAR(1)=1 (-1)... FOR EXPLOSIVE (IMPLOSIVE) SOURCE. C IPAR(2)-IPAR(4)... FREE PARAMETERS. C PAR(1)... MAGNITUDE OF THE SOURCE. C PAR(2)-PAR(6)... FREE PARAMETERS. C C **LU2/5 C C 8) ONE CARD C LU FORMAT(26I3) C LU... NUMBER OF AN ADDITIONAL FILE IN WHICH QUANTITIES, C COMPUTED IN PROGRAM SEIS88 FOR OTHER ELEMENTARY C WAVES THAN THOSE IN LU2, ARE STORED. NOTE THAT C THIS POSSIBILITY CAN BE USED ONLY IN OPERATING C SYSTEMS ALLOWING SPECIFICATION OF A FILE OUTSIDE C THE CORRESPONDING SOURCE PROGRAM (E.G. IN IBM C OPERATING SYSTEM). IF LU.NE.0, THE WHOLE CONTENT C OF THE FILE LU IS READ IN. THE INPUT DATA CARD C NO.7 CAN BE READ IN ANY NUMBER TIMES, UNLESS LU=0. C IF LU=0, COMPUTATION OF SYNTHETIC SEISMOGRAMS C STARTS. C C ***************************************************************** C C TERMINATION OF THE COMPUTATIONS C ******************************* C C THE SYSTEM OF CARDS 1-7 CAN BE REPEATED ANY NUMBER TIMES. C EACH SYSTEM OF CARDS 1-7 PRODUCES A NEW FILE LU3 (IF LU3.NE.0). C THIS POSSIBILITY CAN BE USED AGAIN ONLY IN IBM TYPE OPERATING C SYSTEMS, SEE THE DESCRIPTION OF THE PARAMETER LU. C COMPUTATIONS TERMINATE IF LU2=0 IN INPUT DATA CARD NO.1. C C *************************************************************** C C OUTPUT TABLES C ************* C C OUTPUT ON THE LINE PRINTER C ************************** C ALL THE ADDITIONAL INPUT DATA ARE REPRODUCED ON THE LINE PRINTER. C THE OTHER PRINTOUT IS CONTROLLED BY THE PARAMETER MPRINT, SEE C CARD NO.1. C FOR MPRINT=0, THE DATA LU2/1, LU2/2, LU2/3 AND LU2/4 ARE REPRO- C DUCED. FOR EACH RECEIVER POSITION THE FOLLOWING QUANTITIES ARE C ALSO PRINTED (IN ONE LINE): THE COORDINATE OF THE RECEIVER (X- C COORDINATE IN CASE OF RECEIVERS SITUATED ON THE EARTH'S SURFACE C OR AN INTERFACE, Z-COORDINATE IN CASE OF RECEIVERS SITUATED ON C A VERTICAL PROFILE), THE MAXIMUM AMPLITUDE SMAX OF THE SEISMO- C GRAM, TIME TM CORRESPONDING TO THE FIRST NON-ZERO POINT OF THE C SEISMOGRAM AND THE NUMBER NPS OF POINTS IN THE SEISMOGRAM. C AT THE END, THE QUANTITY SMAXIM (THE MAXIMUM AMPLITUDE IN THE C WHOLE SYSTEM OF SYNTHETIC SEISMOGRAMS) IS PRINTED. C FOR MPRINT=1, ALSO THE DATA LU2/5 ARE PRINTED. C FOR MPRINT=2, THE SEISMOGRAMS ARE ALSO PRINTED. C C OUTPUT TO THE FILE LU3 C ********************** C THE FILE LU3 CONTAINS THE COMPUTED RAY SYNTHETIC SEISMOGRAMS C AT INDIVIDUAL RECEIVER POSITIONS TOGETHER WITH SOME OTHER IM- C PORTANT INFORMATION. THE ZERO VALUES IN THE INITIAL AND FINAL C PARTS OF EACH SEISMOGRAM ARE NOT STORED. THE DATA ARE STORED C IN THE FORMATTED FORM TO ALLOW FOR AN INSPECTION OF COMPUTATIONS. C TO PERFORM THE PLOTTING OF SYNTHETIC SEISMOGRAMS FROM THE FILE C LU3, THE PROGRAM SEISPLOT, INCLUDED IN THIS PACKAGE, CAN BE USED. C THE DATA ON LU3 ARE STORED IN THE FOLLOWING ORDER: C C 1) MTEXT FORMAT(17A4) C ARBITRARY ALPHANUMERIC TEXT DESCRIBING THE COMPUTATIONS. THIS C TEXT WILL APPEAR UNDER THE PLOTS. THE TEXT WAS SPECIFIED IN C THE PROGRAM SEIS88. C 2) MMD,MRED,MCOMP,ITPR,VRED,RSTEP,XSOUR,DT FORMAT(4I5,4F10.5) C MMD... THE NUMBER OF SELECTED RECEIVER POSITIONS C MRED... MRED=0: NON REDUCED TRAVEL TIMES C MRED=0: REDUCED TRAVEL TIMES C MCOMP... MCOMP=0: VERTICAL COMPONENT C MCOMP=1: RADIAL COMPONENT C MCOMP=2: TRANSVERSE COMPONENT C ITPR... ITPR=0: RECEIVERS ALONG THE EARTH'S SURFACE C ITPR=1: RECEIVERS ALONG A VERTICAL PROFILE C ITPR=2: RECEIVERS ALONG AN INTERFACE C VRED... REDUCTION VELOCITY C RSTEP... AVERAGE DIFFERENCE BETWEEN X-COORDINATES OF C NEIGHBOURING RECEIVER POSITIONS. C XSOUR... X-COORDINATE OF THE SOURCE C DT... TIME STEP IN SYNTHETIC SEISMOGRAM. C 3) XMX,SMAXIM FORMAT(22X,F10.5,E15.9) C XMX... COORDINATE OF THE RECEIVER AT WHICH MAXIMUM C AMPLITUDE SMAXIM WAS RECORDED. C SMAXIM... THE MAXIMUM AMPLITUDE OF ALL STORED SYNTHETIC C SEISMOGRAMS. C THE FOLLOWING CARDS 4,5 ARE SUCCESSIVELY REPEATED FOR ALL RECEI- C VER POSITIONS. THEY CONTAIN INFORMATION ABOUT INDIVIDUAL SEIS- C MOGRAMS: C 4) XX,SMAX,TM,NPS FORMAT(F10.5,E15.9,F10.5,I5) C XX... COORDINATE OF THE RECEIVER C SMAX... MAXIMUM AMPLITUDE IN THE SYNTHETIC SEISMOGRAM C TM... TM IS THE TIME CORRESPONDING TO THE FIRST NON- C ZERO POINT IN THE SEISMOGRAM. IF MRED.NE.0, TM C IS A REDUCED TIME. C NPS... NUMBER OF POINTS IN SYNTHETIC SEISMOGRAMS. C 5) SYNTHETIC SEISMOGRAM FORMAT(20I4) C C NOTE: THE MAXIMUM LENGTH OF THE SYNTHETIC SEISMOGRAMS IS 3001 C POINTS. C C C OUTPUT TO THE FILE LU4 C ********************** C THE FILE LU4 IS AN AUXILIARY FILE. DATA STORED IN IT ARE STORED C IN AN UNFORMATTED FORM. THE FILE CAN BE USED FOR PLOTTING ON C SPECIAL SEISMIC PLOTTERS USING STANDARD ROUTINES SUCH AS THE C ROUTINE TRPLOT IN SEP STANFORD. IN CONTRAST TO LU3, COMPLETE C AND NON-NORMALIZED TRACES (ALL EQUALLY LONG) ARE STORED IN LU4. C THE STRUCTURE OF THE FILE IS AS FOLLOWS: C C 1) NPTS,TMIN,DT,NDST,DIST(1),RSTEP C NPTS,TMIN,DT... NUMBER OF POINTS, MINIMUM TIME AND TIME STEP C IN SYNTHETIC SEISMOGRAMS (EQUAL FOR ALL TRACES). C NDST,DIST(1),RSTEP... NUMBER OF TRACES, COORDINATE OF THE C FIRST TRACE AND DISTANCE BETWEEN TRACES. C 2) SAMPLED TRACES ONE AFTER OTHER, FROM THE LEFT TO THE RIGHT. C C C **************************************************************** C