PROGRAM   S Y N T P L
*********************

  PROGRAM SYNTPL IS DESIGNED FOR THE COMPUTATION OF RAY-SYNTHETIC
SEISMOGRAMS FROM THE  DATA STORED IN THE FILE LU2 GENERATED IN
THE PROGRAM SEIS. THE FILE LU2 CONTAINS TRAVEL TIMES, AMPLITUDES 
AND PHASE SHIFTS FOR ALL SPECIFIED RECEIVERS AND ALL ELEMENTARY 
WAVES. THE FILE ALSO CONTAINS THE INITIAL ANGLES OF CORRESPONDING 
RAYS AND CORRESPONDING GLOBAL ABSORPTION FACTORS. PROGRAM
SYNTPL GENERATES A NEW FILE LU3 AND AN AUXILIARY FILE LU4 WITH
SYNTHETIC SEISMOGRAMS. THESE SYNTHETIC SEISMOGRAMS CAN BE PLOTTED
BY STANDARD PLOTTING ROUTINES, AVAILABLE AT MANY INSTITUTIONS.
FOR PLOTTING OF SEISMOGRAMS, IT IS ALSO POSSIBLE TO USE THE 
PROGRAM SEISPLOT INCLUDED IN THIS PACKAGE.
  THE SOURCE-TIME FUNCTION USED IN SYNTPL IS A HARMONIC CARRIER
MODULATED BY A GAUSSIAN ENVELOPE,
  F(T)=EXP(-(OMEGA*T/GAMA)**2)COS(OMEGA*T+PSI),          (1)
WHERE T IS TIME, OMEGA=2.*PI*FREQ, AND FREQ, GAMA AND PSI ARE
THREE FREE PARAMETERS OF THE SOURCE-TIME FUNCTION. BY A PROPER
CHOICE OF THESE PARAMETERS, THE USER CAN SIMULATE A BROAD VARIETY
OF WAVELETS WHICH ARE CLOSE TO THOSE OBSERVED IN SEISMOLOGY AND
SEISMIC PROSPECTING (SEE DETAILS IN CERVENY,MOLOTKOV AND PSENCIK,
RAY METHOD IN SEISMOLOGY, UNIVERZITA KARLOVA, PRAGUE 1977).
  NOTE THAT THE ZERO TIME OF THE SIGNAL (1) CORRESPONDS TO THE
MAXIMUM OF ITS ENVELOPE, NOT TO ITS FIRST ARRIVAL. THIS GIVES
A TIME SHIFT OF THIS SIGNAL IN SYNTHETIC SEISMOGRAMS IN COMPARISON 
WITH OTHER SIMILAR SIGNALS (MUELLER SIGNAL,ETC.). THE USER
WHO IS NOT SATISFIED WITH THIS SHIFT CAN USE THE PARAMETERS
MSHIFT AND SHIFT READ FROM INPUT DATA TO REMOVE IT.
  WHEN COMPUTING THE SYNTHETIC SEISMOGRAMS FROM THE DATA STORED
IN THE FILE LU2, IT IS POSSIBLE TO SELECT ONLY CERTAIN DATA
FROM LU2.
  FIRST SELECTION IS PERFORMED  USING THE INTERNAL WAVE CODE OF
ELEMENTARY WAVES, WHICH CORRESPONDS TO THE SUCCESSIVE NUMBER OF
THE ELEMENTARY WAVE AS IT IS GENERATED BY THE ROUTINE GENER IN
THE PROGRAM SEIS, SEE DESCRIPTION OF THIS ROUTINE FOR DETAILS.
IN THIS WAY, IT IS POSSIBLE TO SELECT, E.G., ONLY PP REFLECTIONS,
OR PS REFLECTIONS, OR TO REMOVE THE DIRECT WAVE, ETC.
  THE SECOND SELECTION APPLIES TO EPICENTRAL DISTANCES. IT IS
POSSIBLE TO CONSTRUCT THE SYNTHETIC SEISMOGRAMS ONLY FOR SELECTED 
RECEIVERS,AND NOT FOR THE WHOLE SYSTEM OF RECEIVERS FOR WHICH 
THE DATA IN SEIS WERE COMPUTED.

INPUT DATA
**********

INPUT DATA CONSIST PARTIALLY OF THE DATA OBTAINED IN THE PROGRAM
SEIS, WHICH ARE STORED IN THE FILE LU2, AND PARTIALLY OF THE 
ADDITIONAL DATA CONTROLLING THE COMPUTATIONS, INTRODUCED DIRECTLY. 
THE INPUT DATA ARE READ IN FROM STANDARD INPUT BY LIST-DIRECTED 
INPUT (FREE FORMAT):


THE DATA STORED AT LU2
**********************

FOR DETAILS SEE THE DESCRIPTION OF THE CONTENT OF THE FILE LU2 
IN PROGRAM SEIS.

1) MTEXT                                  FORMAT(17A4)
2) NDST,KSH,ITPR                          FORMAT(26I3)
3) XSOUR,ZSOUR,TSOUR,RSTEP,ROS,VPS,VSS    FORMAT(8F10.5)
4) DST(I),I=1,NDST                        FORMAT(8F10.5)
5) NCODE,II,T,AX,AY,AZ,PHX,PHY,PHZ,ANGLE,TAST
                            FORMAT(2I3,F10.5,3E12.6,5F10.6)

  THE DATA IN 5 CAN BE REPEATED ON LU2 MANY TIMES. IN THE PROGRAM
SYNTPL, THIS NUMBER MAY NOT EXCEED 2000. OTHERWISE IT WOULD BE
NECESSARY TO CHANGE THE DIMENSIONS IN THE PROGRAM. THE DATA 1-5
CAN BE READ IN FROM DIFFERENT FILES (SEE INPUT DATA CARD NO.7)
GENERATED IN THE PROGRAM SEIS FOR A SPECIFIC MODEL, SOURCE AND
RECEIVER POSITIONS. NOTE THAT THIS POSSIBILITY CAN BE USED ONLY
IN OPERATING SYSTEMS ALLOWING SPECIFICATION OF THE NUMBER OF A
FILE OUTSIDE THE SOURCE PROGRAM (E.G. IBM OPERATING SYSTEM).
  NOTE THAT THE MAXIMUM NUMBER OF THE RECEIVER POSITIONS IS 100.

ADDITIONAL INPUT DATA
*********************

THESE DATA SHOULD BE INTRODUCED DIRECTLY (CARDS,TERMINAL). THEY
CONTROL THE COMPUTATION OF SYNTHETIC SEISMOGRAMS. IN THE 
DESCRIPTION OF ADDITIONAL INPUT DATA, WE SHALL DENOTE PLACES WHERE 
THE DATA FROM LU2 ARE READ IN BY **LU2/1, **LU2/2, ETC.

1) ONE CARD
   NEXT,MPRINT
      NEXT ...  NEXT=1: CONTINUATION OF COMPUTATIONS
                NEXT=0: TERMINATION OF COMPUTATIONS.
      MPRINT... THE PARAMETER WHICH CONTROLS THE PRINTOUT ON THE
                LINE PRINTER, SEE DETAILS IN THE SECTION "OUTPUT
                TABLES".

**LU2/1
**LU2/2
**LU2/3
**LU2/4

2) ONE CARD, VARIOUS SWITCHES.
   MCOMP,MRED,MSELEC,MEPIC,MSHIFT,KABS,MSOUR    
      MCOMP...  CONTROLS THE CHOICE OF THE COMPONENT OF THE
                DISPLACEMENT VECTOR.
             MCOMP=0: VERTICAL COMPONENT.
             MCOMP=1: RADIAL COMPONENT.
             MCOMP=2: TRANSVERSE COMPONENT.
      MRED...   CONTROLS REDUCTION OF TRAVEL TIMES. IF VSP
                SYNTHETICS ARE COMPUTED, MRED SET AUTOMATICALLY, 
                MRED=0.
             MRED=0: TRAVEL TIMES NOT REDUCED.
             MRED=1: TRAVEL TIMES REDUCED. THE REDUCTION
                VELOCITY IS GIVEN IN THE CARD NO.3.
      MSELEC... CONTROLS SELECTION OF ELEMENTARY WAVES.
             MSELEC=0: ALL WAVES STORED AT LU2 ARE USED TO
                CONSTRUCT SYNTHETIC SEISMOGRAMS.
             MSELEC=1: ONLY SELECTED WAVES ARE USED TO CONSTRUCT 
                SYNTHETIC SEISMOGRAMS, SEE CARD NO.5.
      MEPIC...  CONTROLS SELECTION OF RECEIVERS.
             MEPIC=0: THE SYNTHETIC SEISMOGRAMS ARE CONSTRUCTED AT 
                ALL RECEIVER POSITIONS.
             MEPIC=1: THE SYNTHETIC SEISMOGRAMS ARE CONSTRUCTED ONLY 
                AT SELECTED RECEIVER POSITIONS, SEE CARD NO.4.
      MSHIFT... CONTROLS THE TIME SHIFT OF THE SIGNAL.
             MSHIFT=0: NO TIME SHIFT. ZERO TIME CORRESPONDS
                TO THE MAXIMUM OF THE ENVELOPE OF THE SIGNAL (1).
             MSHIFT=1: THE TIME SHIFT IS INTRODUCED AUTOMATICALLY. 
                IT CORRESPONDS TO THE TIME INTERVAL BETWEEN 
                THE MAXIMUM OF THE ENVELOPE AND THE TIME
                WHERE THE ENVELOPE IS 0.1 FROM ITS MAXIMUM VALUE.
             MSHIFT=2: THE TIME SHIFT IS INTRODUCED MANUALLY,
                BY THE PARAMETER 'SHIFT' IN CARD NO.3.
      KABS...   SPECIFIES THE DISSIPATION TYPE.
             KABS=0: PERFECTLY ELASTIC MEDIUM.
             KABS=1: SLIGHTLY ABSORBING MEDIUM WITH NON-CAUSAL
                ABSORPTION.
             KABS=2: SLIGHTLY ABSORBING MEDIUM WITH CAUSAL
                ABSORPTION (FUTTERMAN'S MODEL).
      MSOUR...  CONTROLS THE TYPE OF THE POINT SOURCE.
             MSOUR=0: UNIT ISOTROPIC RADIATION PATTERN.
             MSOUR=1: SINGLE FORCE.
             MSOUR=2: DOUBLE COUPLE.
             MSOUR=3: EXPLOSIVE (IMPLOSIVE) SOURCE.

NOTE: SINGLE FORCE AND DOUBLE COUPLE POINT SOURCES CAN BE ARBITRARILY 
ORIENTED IN THE 3-D SPACE (THEY DO NOT HAVE TO BE CONFINED TO 
THE VERTICAL (X,Z) PLANE).

3) ONE CARD, SPECIFICATION OF SYNTHETIC SEISMOGRAMS.
   TMIN,DT,TMAX,FREQ,GAMA,PSI,VRED,SHIFT      
      TMIN,DT,TMAX... MINIMUM TIME, TIME STEP AND THE MAXIMUM
                TIME OF COMPUTED SYNTHETIC SEISMOGRAMS. IN CASE
                OF MRED=1, THESE QUANTITIES CORRESPOND TO THE
                REDUCED TRAVEL TIME.
      FREQ,GAMA,PSI...  PARAMETERS WHICH SPECIFY THE SOURCE-TIME
                FUNCTION (1). WHEN FREQ NOT SPECIFIED, FREQ=4.
                WHEN GAMA NOT SPECIFIED, GAMA=4.
      VRED...   REDUCTION VELOCITY. WHEN VRED NOT SPECIFIED, VRED=8.
      SHIFT...  WHEN MSHIFT=2, SHIFT GIVES THE TIME SHIFT OF THE
                SIGNAL (1). OTHERWISE IT MAY BE ARBITRARY.

4) ONE CARD, SELECTION OF RECEIVER POSITIONS.
   INCLUDED ONLY WHEN MEPIC.NE.0.
   NEPIC,(IEP(I),I=1,NEPIC)                
      NEPIC...NUMBER OF SELECTED RECEIVER POSITIONS AT WHICH
                THE SYNTHETIC SEISMOGRAMS SHOULD BE CONSTRUCTED.
      IEP(1),IEP(2),...,IEP(NEPIC)... SEQUENTIAL NUMBERS OF
                SELECTED RECEIVER POSITIONS. THE RECEIVER POSITIONS 
                ARE NUMBERED FROM THE LEFT TO THE RIGHT.

5) ONE CARD, SELECTION OF WAVES
   INCLUDED ONLY WHEN MSELEC.NE.0.
   NSELEC,(ISEL(I),I=1,NSELEC)              
      NSELEC...NUMBER OF SELECTED WAVES USED FOR THE CONSTRUCTION 
                OF SYNTHETIC SEISMOGRAMS.
      ISEL(1),ISEL(2),...,ISEL(NSELEC)... INTERNAL CODES OF
                SELECTED WAVES. THEY CORRESPOND TO THE SEQUENTIAL
                NUMBER OF THE WAVE AS IT IS GENERATED BY THE
                ROUTINE GENER IN SEIS. FOR DETAILS, SEE THE
                DESCRIPTION OF ROUTINE GENER.

6) ONE CARD, REFERENCE FREQUENCY.
   INCLUDED ONLY IF KABS.NE.0.
   FREF,QRED                                
      FREF... REFERENCE FREQUENCY FOR ABSORPTION COMPUTATIONS.
                THE  VELOCITY VALUES AND QUALITY FACTOR USED IN
                THE MODEL ARE SPECIFIED FOR THE FREQUENCY FREF.
                DEFAULT VALUE, FREF=1.
      QRED... FACTOR BY WHICH GLOBAL ABSORPTION FACTORS "TAST"
                ALONG ALL RAYS ARE MULTIPLIED. IT CORRESPONDS TO
                DIVISION OF THE CORRESPONDING Q-FACTOR BY QRED.
                DEFAULT VALUE, QRED=1.

7) ONE CARD, SPECIFICATION OF PARAMETERS OF THE SOURCE.
   INCLUDED ONLY IF MSOUR.NE.0.
   IPAR(1),...,IPAR(4),PAR(1),...,PAR(6)     

   MSOUR=1: SINGLE FORCE.

      IPAR(1)-IPAR(4)... FREE PARAMETERS.
      PAR(1)... FORCE AZIMUTH - THE ANGLE, IN RADIANS, MADE BY
                THE PROJECTION OF THE FORCE VECTOR INTO A HORIZONTAL 
                PLANE AND THE POSITIVE X-AXIS. POSITIVE CLOCKWISE.
      PAR(2)... MAGNITUDE OF THE FORCE.
      PAR(3)... FORCE DECLINATION - THE ANGLE, IN RADIANS, MADE
                BY THE FORCE VECTOR AND A HORIZONTAL PLANE.
      PAR(4)-PAR(6)... FREE PARAMETERS.

   MSOUR=2: DOUBLE COUPLE.

      IPAR(1)-IPAR(4)... FREE PARAMETERS.
      PAR(1)... DIP ANGLE, IN RADIANS.
      PAR(2)... SOURCE MOMENT.
      PAR(3)... STRIKE ANGLE, IN RADIANS.
      PAR(4)... RAKE ANGLE, IN RADIANS.
      PAR(5)-PAR(6)... FREE PARAMETERS.
NOTE: THE FORMULAE FOR THE DOUBLE COUPLE FOLLOW THE DEFINITION IN 
AKI AND RICHARDS (1980).

   MSOUR=3: EXPLOSIVE (IMPLOSIVE) SOURCE.

      IPAR(1)=1 (-1)... FOR EXPLOSIVE (IMPLOSIVE) SOURCE.
      IPAR(2)-IPAR(4)... FREE PARAMETERS.
      PAR(1)... MAGNITUDE OF THE SOURCE.
      PAR(2)-PAR(6)... FREE PARAMETERS.

**LU2/5

8) ONE CARD
   LU                                        
      LU...     NUMBER OF AN ADDITIONAL FILE IN WHICH QUANTITIES,
                COMPUTED IN PROGRAM SEIS FOR OTHER ELEMENTARY
                WAVES THAN THOSE IN LU2, ARE STORED. NOTE THAT
                THIS POSSIBILITY CAN BE USED ONLY IN OPERATING
                SYSTEMS ALLOWING SPECIFICATION OF A FILE OUTSIDE
                THE CORRESPONDING SOURCE PROGRAM (E.G. IN IBM
                OPERATING SYSTEM). IF LU.NE.0, THE WHOLE CONTENT
                OF THE FILE LU IS READ IN. THE INPUT DATA CARD
                NO.7 CAN BE READ IN ANY NUMBER TIMES, UNLESS LU=0.
                IF LU=0, COMPUTATION OF SYNTHETIC SEISMOGRAMS STARTS.

*****************************************************************

TERMINATION OF THE COMPUTATIONS
*******************************

THE SYSTEM OF CARDS 1-7 CAN BE REPEATED ANY NUMBER TIMES.
EACH SYSTEM OF CARDS 1-7 PRODUCES A NEW FILE LU3 (IF LU3.NE.0).
THIS POSSIBILITY CAN BE USED AGAIN ONLY IN IBM TYPE OPERATING
SYSTEMS, SEE THE DESCRIPTION OF THE PARAMETER LU.
COMPUTATIONS TERMINATE IF NEXT=0 IN INPUT DATA CARD NO.1.

***************************************************************

OUTPUT TABLES
*************

OUTPUT ON THE LINE PRINTER
**************************
  ALL THE ADDITIONAL INPUT DATA ARE REPRODUCED ON THE LINE PRINTER.
THE OTHER PRINTOUT IS CONTROLLED BY THE PARAMETER MPRINT, SEE
CARD NO.1.
  FOR MPRINT=0, THE DATA LU2/1, LU2/2, LU2/3 AND LU2/4 ARE REPRODUCED. 
FOR EACH RECEIVER POSITION THE FOLLOWING QUANTITIES ARE ALSO PRINTED 
(IN ONE LINE): THE COORDINATE OF THE RECEIVER (X-COORDINATE IN CASE 
OF RECEIVERS SITUATED ON THE EARTH'S SURFACE OR AN INTERFACE, 
Z-COORDINATE IN CASE OF RECEIVERS SITUATED ON A VERTICAL PROFILE), 
THE MAXIMUM AMPLITUDE SMAX OF THE SEISMOGRAM, TIME TM CORRESPONDING 
TO THE FIRST NON-ZERO POINT OF THE SEISMOGRAM AND THE NUMBER NPS OF 
POINTS IN THE SEISMOGRAM. AT THE END, THE QUANTITY SMAXIM (THE MAXIMUM 
AMPLITUDE IN THE WHOLE SYSTEM OF SYNTHETIC SEISMOGRAMS) IS PRINTED.
  FOR MPRINT=1, ALSO THE DATA LU2/5 ARE PRINTED.
  FOR MPRINT=2, THE SEISMOGRAMS ARE ALSO PRINTED.

OUTPUT TO THE FILE LU3
**********************
  THE FILE LU3 CONTAINS THE COMPUTED RAY SYNTHETIC SEISMOGRAMS
AT INDIVIDUAL RECEIVER POSITIONS TOGETHER WITH SOME OTHER IMPORTANT 
INFORMATION. THE ZERO VALUES IN THE INITIAL AND FINAL PARTS OF EACH 
SEISMOGRAM ARE NOT STORED. THE DATA ARE STORED IN THE FORMATTED FORM 
TO ALLOW FOR AN INSPECTION OF COMPUTATIONS. 
  TO PERFORM THE PLOTTING OF SYNTHETIC SEISMOGRAMS FROM THE FILE
LU3, THE PROGRAM SEISPLOT, INCLUDED IN THIS PACKAGE, CAN BE USED.
THE DATA ON LU3 ARE STORED IN THE FOLLOWING ORDER:

1) MTEXT                           FORMAT(17A4)
   ARBITRARY ALPHANUMERIC TEXT DESCRIBING THE COMPUTATIONS. THIS
   TEXT WILL APPEAR UNDER THE PLOTS. THE TEXT WAS SPECIFIED IN
   THE PROGRAM SEIS.
2) MMD,MRED,MCOMP,ITPR,VRED,RSTEP,XSOUR,DT  FORMAT(4I5,4F10.5)
      MMD...    THE NUMBER OF SELECTED RECEIVER POSITIONS
      MRED.. MRED=0: NON REDUCED TRAVEL TIMES
             MRED=0: REDUCED TRAVEL TIMES
      MCOMP..MCOMP=0: VERTICAL COMPONENT
             MCOMP=1: RADIAL COMPONENT
             MCOMP=2: TRANSVERSE COMPONENT
      ITPR.. ITPR=0: RECEIVERS ALONG THE EARTH'S SURFACE
             ITPR=1: RECEIVERS ALONG A VERTICAL PROFILE
             ITPR=2: RECEIVERS ALONG AN INTERFACE
      VRED...   REDUCTION VELOCITY
      RSTEP...  AVERAGE DIFFERENCE BETWEEN X-COORDINATES OF
                NEIGHBOURING RECEIVER POSITIONS.
      XSOUR...  X-COORDINATE OF THE SOURCE
      DT...     TIME STEP IN SYNTHETIC SEISMOGRAM.
3) XMX,SMAXIM                        FORMAT(22X,F10.5,E15.9)
      XMX...    COORDINATE OF THE RECEIVER AT WHICH MAXIMUM
                AMPLITUDE SMAXIM WAS RECORDED.
      SMAXIM... THE MAXIMUM AMPLITUDE OF ALL STORED SYNTHETIC
                SEISMOGRAMS.
THE FOLLOWING CARDS 4,5 ARE SUCCESSIVELY REPEATED FOR ALL RECEIVER 
POSITIONS. THEY CONTAIN INFORMATION ABOUT INDIVIDUAL SEISMOGRAMS:
4) XX,SMAX,TM,NPS                 FORMAT(F10.5,E15.9,F10.5,I5)
      XX...     COORDINATE OF THE RECEIVER
      SMAX...   MAXIMUM AMPLITUDE IN THE SYNTHETIC SEISMOGRAM
      TM...     TM IS THE TIME CORRESPONDING TO THE FIRST NON-
                ZERO POINT IN THE SEISMOGRAM. IF MRED.NE.0, TM
                IS A REDUCED TIME.
      NPS...    NUMBER OF POINTS IN SYNTHETIC SEISMOGRAMS.
5) SYNTHETIC SEISMOGRAM              FORMAT(20I4)

NOTE: THE MAXIMUM LENGTH OF THE SYNTHETIC SEISMOGRAMS IS 3001
POINTS.

OUTPUT TO THE FILE LU4
**********************

  THE FILE LU4 IS AN AUXILIARY FILE. DATA STORED IN IT ARE STORED
IN AN UNFORMATTED FORM. THE FILE CAN BE USED FOR PLOTTING ON
SPECIAL SEISMIC PLOTTERS USING STANDARD ROUTINES SUCH AS THE
ROUTINE TRPLOT IN SEP STANFORD. IN CONTRAST TO LU3, COMPLETE
AND NON-NORMALIZED TRACES (ALL EQUALLY LONG) ARE STORED IN LU4.
THE STRUCTURE OF THE FILE IS AS FOLLOWS:

1) NPTS,TMIN,DT,NDST,DIST(1),RSTEP
   NPTS,TMIN,DT... NUMBER OF POINTS, MINIMUM TIME AND TIME STEP
                IN SYNTHETIC SEISMOGRAMS (EQUAL FOR ALL TRACES).
   NDST,DIST(1),RSTEP... NUMBER OF TRACES, COORDINATE OF THE
                FIRST TRACE AND DISTANCE BETWEEN TRACES.
2) SAMPLED TRACES ONE AFTER OTHER, FROM THE LEFT TO THE RIGHT.

 ****************************************************************