C
C Program SP (Seismogram Plotting) to plot seismograms previously stored
C in the GSE data exchange format.
C
C Version: 5.50
C Date: 2000, December 11
C
C Coded by: Ludek Klimes
C     Department of Geophysics, Charles University Prague,
C     Ke Karlovu 3, 121 16 Praha 2, Czech Republic,
C     E-mail: klimes@seis.karlov.mff.cuni.cz
C
C.......................................................................
C                                                    
C Description of data files:
C
C Input data read from the standard input device (*):
C     The data are read by the list directed input (free format) and
C     consist of a single string 'SEP':
C     'SEP'...String in apostrophes containing the name of the input
C             SEP parameter or history file with the input data.
C     No default, 'SEP' must be specified and cannot be blank.
C
C                                                     
C Input data file 'SEP':
C     File 'SEP' has the form of the SEP
C     parameter file.  The parameters, which do not differ from their
C     defaults, need not be specified in file 'SEP'.
C Names of the input files:
C     SS='string'... String with the name of the input data file in the
C             GSE data exchange format, containing the seismograms to be
C             plotted.
C             If the source names are specified in the comment lines
C               of the GSE file, the hypocentral coordinates of the
C               sources are taken, in the order of preference, from file
C               PTS (if the filename is specified), file SRC (if the
C               filename is specified), corresponding GSE file SS* (if
C               the coordinates are present), default to zeros.
C             If the source name is not specified in the comment lines
C               of the GSE file, the hypocentral coordinates of the
C               sources are taken, in the order of preference, from the
C               the first point in file SRC (if the filename is
C               specified), GSE file SS (if the coordinates are
C               present), default to zeros.
C             The coordinates of the receivers are taken, in the order
C               of preference, from file PTS (if the filename is
C               specified), file REC (if the filename is specified),
C               GSE file SS.
C             Only the first 6 characters of source and receiver names
C               are significant.
C             Description of GSE file SS
C             Default: SS='ss.gse'
C     SS1='string', SS2='string', ..., SS9='string'... Strings with the
C             additional, optional names of the input data files in the
C             GSE data exchange format, containing the seismograms to be
C             plotted simultaneously with seismograms given by parameter
C             SS into the same figures.  The order of plottting is SS,
C             SS1, SS2, ..., SS9, considering just nonblank filenames.
C             The seismograms are plotted in colours given by parameters
C             KOLOR, KOLOR1, KOLOR2, ..., KOLOR9, respectively.  Refer
C             to file calcops.rgb for the
C             definitions of the colours.  The frame and description
C             remain in colour number 1, which is probably black.
C             Defaults: SS1=' ', SS2=' ', SS3=' ', ..., SS9=' '
C     SRC='string'... String with the name of the input data file
C             with the name(s) and coordinates of the hypocentre(s).
C             If specified, this file serves (a) to select the sources
C             for plotting seismograms, (b) to update the hypocentre
C             coordinates if PTS=' ' or if the source name is not
C             specified in the comment lines of the GSE file and file
C             PTS thus cannot be used.
C             If SRC is not blank and the source names are specified
C               comment lines of the GSE file, only seismograms
C               corresponding to the sources listed in file SRC will be
C               plotted, otherwise no selection according to the source
C               names is performed.
C             If the source names are specified in the comment lines
C               of the GSE file, the hypocentral coordinates of the
C               sources are taken, in the order of preference, from file
C               PTS (if the filename is specified), file SRC (if the
C               filename is specified), corresponding GSE file SS* (if
C               the coordinates are present), default to zeros.
C             If the source name is not specified in the comment lines
C               of the GSE file, the hypocentral coordinates are taken,
C               in the order of preference, from the the first point in
C               file SRC (if the filename is specified), corresponding
C               GSE file SS* (if the coordinates are present), default
C               to zeros.
C             Only the first 6 characters of source names are
C               significant.
C             The hypocentral coordinates are used only for optional
C             travel-time reduction on a given velocity, or for
C             amplitude power scaling.  File SRC thus usually need not
C             be specified if the seismograms are generated by programs
C             'greenss.for' or 'ss.for', because those programs enter
C             the hypocentral coordinates directly to the comments of
C             the data section in GSE data exchange file.
C             Description of file SRC
C             Default: SRC=' '
C     REC='string'... String with the name of the input data file
C             containing the list of receivers.
C             If specified, this file serves (a) to select the receivers
C             for plotting seismograms, (b) to update the receiver
C             coordinates if PTS=' '.
C             If REC=' ' seismograms at all receivers will be plotted,
C               otherwise only seismograms at the receivers listed in
C               file REC will be plotted.
C             The coordinates of the receivers are taken, in the order
C               of preference, from file PTS (if the filename is
C               specified), file REC (if the filename is specified),
C               corresponding GSE file SS*.
C             File REC also enables to determine the number NREC of all
C               receivers, in which the seimograms may be plotted, for
C               scaling purposes.  If REC=' ', NREC=999 for scaling.
C             Only the first 6 characters of receiver names are
C               significant.
C             In most cases, file REC will be the same as for the
C             calculation of synthetic seismograms.
C             Parameter REC has to be specified and cannot be blank if
C             KODESP=0, see below.  If KODESP=0, the horizontal position
C             of a plotted seismogram corresponds to the position of the
C             corresponding receiver in file REC.
C             Description of file REC
C             Default: REC=' '
C     FTT='string'... String with the name of the input data file
C             with the list of travel times to be highlighted.
C             The travel times are specified in dependence on the
C             source and receiver names.  The coordinates of the source
C             and receivers are taken from file PTS.  If PTS=' ', the
C             coordinates of the source are taken from file SRC and the
C             coordinates of receivers are taken from file REC.
C             If PTS=' ' and SRC=' ', the source coordinates are assumed
C             zero.  If FTT.NE.' ' and PTS=' ', REC cannot be blank.
C             If FTT=' ', no travel times are highlighted.
C             Description of file FTT
C             Default: FTT=' '
C     PTS='string'... String with the name of the input data file
C             with the coordinates corresponding to the source and
C             receiver names.  Since this file is given just to specify
C             the coordinates, the coordinates from this file have the
C             highest priority.  This feature is especially useful when
C             changing the coordinate system.  The seismogram files
C             need not be changed and may preserve the old coordinates.
C             The same advantage applies to files SRC, REC and SPHILI.
C             Only the first 6 characters of point names are
C             significant.
C             Description of file PTS
C             Default: PTS=' '
C     SPHILI='string'... String with the name of the input data file
C             with the list of times to be highlighted.
C             The times are specified in dependence on the names and
C             coordinates of receivers.
C             If REC is not blank and the receiver names are specified
C               in file SPHILI, only times corresponding to the
C               receivers listed in file REC will be plotted, otherwise
C               no selection according to the receiver names is
C               performed.
C             The coordinates of the receivers are taken, in the order
C               of preference, from file PTS (if filename PTS is
C               specified and the point names in file SPHILI are not
C               blank), from file SPHILI.
C             If SPHILI=' ', no times are highlighted.
C             Description of file SPHILI
C             Default: SPHILI=' '
C Names of output files:
C     SP1='string'... String with the name of the first output
C             PostScript file, usually contaning the plot of the first
C             component of the seismograms.
C             If blank, the file is not created.
C             Default: SP1='ss1.gse'
C     SP2='string'... String with the name of the second output
C             PostScript file, usually contaning the plot of the second
C             component of the seismograms.
C             If blank, the file is not created.
C             Default: SP2='ss2.gse'
C     SP3='string'... String with the name of the third output
C             PostScript file, usually contaning the plot of the third
C             component of the seismograms.
C             If blank, the file is not created.
C             Default: SP3='ss3.gse'
C Component selection (mostly not needed):
C     KOMP1=integer... Component of the seismograms of file given by
C             parameter SS, plotted into the output file given by
C             parameter SP1.
C             Default: KOMP1=1
C     KOMP11=integer, KOMP12=integer, KOMP13=integer, KOMP14=integer,
C     KOMP15=integer, KOMP16=integer, KOMP17=integer, KOMP18=integer,
C     KOMP19=integer... Components of the seismograms of files given by
C             parameters SS1, SS2, SS3, SS4, SS5, SS6, SS7, SS8, SS9,
C             respectively, plotted into the output file given by
C             parameter SP1.
C             Defaults: KOMP11=KOMP1, KOMP12=KOMP1, KOMP13=KOMP1,
C                       KOMP14=KOMP1, KOMP15=KOMP1, KOMP16=KOMP1,
C                       KOMP17=KOMP1, KOMP18=KOMP1, KOMP19=KOMP1
C     KOMP2=integer... Component of the seismograms of file given by
C             parameter SS, plotted into the output file given by
C             parameter SP2.  Analogous to KOMP1.
C             Default: KOMP2=2
C     KOMP21=integer, KOMP22=integer, KOMP23=integer, KOMP24=integer,
C     KOMP25=integer, KOMP26=integer, KOMP27=integer, KOMP28=integer,
C     KOMP29=integer... Analogous to KOMP11 to KOMP19, but for file SP2.
C             Defaults equal the value of KOMP2.
C     KOMP3=integer... Component of the seismograms of file given by
C             parameter SS, plotted into the output file given by
C             parameter SP3.  Analogous to KOMP1.
C             Default: KOMP3=3
C     KOMP31=integer, KOMP32=integer, KOMP33=integer, KOMP34=integer,
C     KOMP35=integer, KOMP36=integer, KOMP37=integer, KOMP38=integer,
C     KOMP39=integer... Analogous to KOMP11 to KOMP19, but for file SP2.
C             Defaults equal the value of KOMP3.
C Data to control plotting:
C Colours:
C     KOLOR=positive integer, KOLOR1=positive integer,
C     KOLOR2=positive integer, ..., KOLOR9=positive integer... Colours
C             to plot seismograms of files SS,SS1, SS2, ..., SS9,
C             respectively.  Colour indices correspond to dummy
C             parameter INP of CalComp subroutine
C             NEWPEN.
C             The colours corresponding to the indices may be defined
C             or changed in file calcops.rgb.
C             Default: KOLOR=1, KOLOR1=2, KOLOR2=3, ..., KOLOR9=10
C     KOLORTT=positive integer... Colour to plot the travel times of
C             optional file FTT.  It is also used as the default colour
C             for optional file SPHILI.
C             Travel times are not plotted if KOLORTT is not positive.
C             Default: KOLORTT=1
C     KOLORTD=integer... Colour to plot the error bar of optional file
C             FTT.  It is also used as the default colour for optional
C             file SPHILI.
C             The error bar is not plotted if KOLORTD is negative.
C             If KOLORTD=0 (white), the contour of the rectangle is
C             plotted in colour KOLORTT.
C             Default: KOLORTD=0
C Data for the time axis (vertical):
C     SPTMIN=real... Minimum time (or reduced time), i.e.
C             the time corresponding to the bottom of the
C             seismogram plot.
C             Default: SPTMIN=0.
C     SPTMAX=real... Maximum time (or reduced time), i.e.
C             the time corresponding to the top of the
C             seismogram plot.
C             Default: SPTMAX=1.
C     SPTLEN=real... Length of the vertical time axis in cm.
C             Default: SPTLEN=10.
C     SPTDIV=real... ABS(SPTDIV) is the number of intervals along the
C             time axis, starting at the bottom.  Must be SPTDIV.NE.0.
C             SPTDIV.GT.0.: the marks at the endpoints of intervals will
C               be supplemented with corresponding times (or reduced
C               times).
C             SPTDIV.LT.0.: the time axis will have no description.
C             Default: SPTDIV=1.
C     SPTSUB=real... Number of subintervals to be marked in each
C             interval.  Must be SPTSUB.GT.0.
C             Default: SPTSUB=1.
C     SPVRED=real... Reduction velocity.  If non-zero, the time at each
C             receiver is reduced by the value of RR/SPVRED, where RR is
C             the hypocentral distance.
C             No time reduction is applied if SPVRED=0.
C             Default: SPVRED=0.
C Data for the distance axis (horizontal):
C     KODESP=integer... Specifies the distribution and description of
C             seismograms in the plot.  See the description of SPXMIN,
C             SPXMAX,SPYMIN,SPYMAX below.
C             Default: KODESP=0
C     SPXMIN=real, SPXMAX=real, SPYMIN=real, SPYMAX=real:
C             For KODESP=0: Horizontal axis represents the index of the
C               receiver, corresponding to the receiver position in
C               the file given by parameter REC.
C               SPXMIN and SPXMAX are the receiver indices corresponding
C               to the lef-hand and right-hand sides of the frame.
C               Example: For default values SPXMIN=0 and SPXMAX=NREC+1,
C                 where NREC is the number of receivers in file REC, the
C                 horizontal axis is divided into NREC+1 intervals.
C                 The seismogram at the Ith receiver is then plotted at
C                 the endpoint of the Ith interval.
C               SPYMIN and SPYMAX have no meaning.
C               If SPXDIV.GT.0., the horizontal axis is denoted by the
C               names of the receivers corresponding to the plotted
C               seismograms, otherwise it has no description.
C               Parameter REC has to be specified and cannot be blank in
C               this case.
C             For KODESP=1: Horizontal axis represents the profile with
C               endpoints (X1,X2)=(SPXMIN,SPYMIN) and
C               (X1,X2)=(SPXMAX,SPYMAX),
C               situated in a horizontal plane.  The seismograms are
C               located at the orthogonal projections of the receivers
C               onto the profile.  If SPXDIV.GT.0., the horizontal axis
C               is supplemented with the values of the X1 coordinates.
C             For KODESP=2: Horizontal axis represents the profile with
C               endpoints (X1,X2)=(SPXMIN,SPYMIN) and
C               (X1,X2)=(SPXMAX,SPYMAX),
C               situated in a horizontal plane.  The seismograms are
C               located at the orthogonal projections of the receivers
C               onto the profile.  If SPXDIV.GT.0., the horizontal axis
C               is supplemented with the values of both X1 and X2
C               coordinates.
C             For KODESP=3: Horizontal axis represents the vertical
C               profile with endpoints X3=SPXMIN and X3=SPXMAX.
C               The seismograms are located at the horizontal
C               projections of the receivers onto the profile.
C               If SPXDIV.GT.0., the horizontal axis is supplemented
C               with the values of X3 coordinate.
C               SPYMIN and SPYMAX have no meaning.
C             For KODESP=4: Horizontal axis represents the hypocentral
C               distance with endpoints RR=SPXMIN and RR=SPXMAX.
C               The seismograms are located according to the hypocentral
C               distances the receivers.
C               If SPXDIV.GT.0., the horizontal axis is supplemented
C               with the values of the hypocentral distance.
C               SPYMIN and SPYMAX have no meaning.
C             Default: SPXMIN=0., SPXMAX=NREC+1, SPYMIN=0., SPYMAX=0.,
C               where NREC is the number of receivers in file REC.
C     SPXMIN1=real, SPXMIN2=real, SPXMIN3=real, SPXMIN4=real,
C     SPXMIN5=real, SPXMIN6=real, SPXMIN7=real, SPXMIN8=real,
C     SPXMIN9=real... Analogous to SPXMIN, but for files SS1 to SS9.
C             Do not influence the description of the horizontal axis
C             nor highlighting of times.
C             Usually need not be specified.
C             Defaults equal the value of SPXMIN.
C     SPXMAX1=real, SPXMAX2=real, SPXMAX3=real, SPXMAX4=real,
C     SPXMAX5=real, SPXMAX6=real, SPXMAX7=real, SPXMAX8=real,
C     SPXMAX9=real... Analogous to SPXMAX, but for files SS1 to SS9.
C             Do not influence the description of the horizontal axis
C             nor highlighting of times.
C             Usually need not be specified.
C             Defaults equal the value of SPXMAX.
C     SPYMIN1=real, SPYMIN2=real, SPYMIN3=real, SPYMIN4=real,
C     SPYMIN5=real, SPYMIN6=real, SPYMIN7=real, SPYMIN8=real,
C     SPYMIN9=real... Analogous to SPYMIN, but for files SS1 to SS9.
C             Do not influence the description of the horizontal axis
C             nor highlighting of times.
C             Usually need not be specified.
C             Defaults equal the value of SPYMIN.
C     SPYMAX1=real, SPYMAX2=real, SPYMAX3=real, SPYMAX4=real,
C     SPYMAX5=real, SPYMAX6=real, SPYMAX7=real, SPYMAX8=real,
C     SPYMAX9=real... Analogous to SPYMAX, but for files SS1 to SS9.
C             Do not influence the description of the horizontal axis
C             nor highlighting of times.
C             Usually need not be specified.
C             Defaults equal the value of SPYMAX.
C     SPXLEN=real... Length of the horizontal axis in cm.
C             Default: SPXLEN=FLOAT(NREC+1), where NREC is the number of
C                      receivers in file REC.
C     SPXDIV=real... ABS(SPXDIV) is the number of intervals along the
C               horizontal axis, starting at the left.
C               Must be SPXDIV.NE.0.
C             SPXDIV.GT.0.: The marks at the endpoints of intervals will
C               be supplemented with the corresponding values.
C             SPXDIV.LT.0.: Horizontal axis will have no description.
C             Default: SPXDIV=1.
C     SPXSUB=real... Number of subintervals to be marked in each
C             interval.  Must be SPXSUB.GT.0.
C             Default: SPXSUB=1.
C Amplitude scaling:
C     NORMSP=integer... Type of amplitude scaling:
C             NORMSP.LT.0: Like NORMSP=0, but the maximum amplitude is
C               calculated only over the plotted part of seismogram.
C             NORMSP.EQ.0: Maximum amplitude at each trace normalized to
C               the given value.
C               Simple and universal option.  No other amplitude scaling
C               parameter usually needs to be specified for this option,
C               although the input parameters SPAMP and SPAMPi enable
C               for possible amplitude scaling.
C               Amplitude scale AA is
C                 AA=SPAMPi*XD/AMAX
C               where AMAX is the maximum amplitude in each seismogram
C               and XD is the average distance between seismograms, i.e.
C               XD=SPXLEN/(NREC+1), where NREC is the number of
C               receivers in file REC.  The receiver file given by input
C               parameter REC is thus required to determine NREC.
C             NORMSP.GT.0: All seismograms have the same scaling or
C               power scaling.
C               Amplitude scale AA is
C                 AA=SPAMPi*(RR/SPDIST)**SPOWER
C               where SPDIST is the hypocentral distance of the receiver
C               under consideration.
C             Default: NORMSP=0
C     SPAMP=real... Amplitude scale for all 3 components.
C             Default: SPAMP=1.
C     SPAMP1=real, SPAMP2=real, ..., SPAMP9=real... Amplitude scales
C             SPAMP individually set for optional input GSE files
C             given by parameters SS1, SS2, ..., SS9, respectively.
C             Defaults: SPAMP1=SPAMP, SPAMP2=SPAMP, ..., SPAMP9=SPAMP
C     SPAMPX1=real, SPAMPX2=real, SPAMPX3=real... Amplitude scale
C             multiplicative factors for individual output files SP1,
C             SP2 and SP3, usually corresponding to different
C             components.
C             SPAMP will be multiplied by SPAMPX1 for file SP1 (usually
C             component 1), by SPAMPX2 for file SP1 (usually component
C             2), by SPAMPX3 for file SP1 (usually component 3).
C             Default: SPAMPX1=1, SPAMPX2=1, SPAMPX3=1
C     SPOWER=real... Exponent of the amplitude power scaling.
C             Need not be specified if power scaling is not required.
C             Has no meaning if NORMSP.LE.0.
C             Default: SPOWER=0.
C     SPOWER1=real, SPOWER2=real, ..., SPOWER9=real... Exponents of the
C             amplitude power scaling for optional input GSE files
C             given by parameters SS1, SS2, ..., SS9, respectively.
C             Defaults: SPOWER1=SPOWER, SPOWER2=SPOWER, ...,
C                       SPOWER9=SPOWER
C     SPDIST=real... Reference hypocentral distance for the amplitude
C             power scaling.
C             Need not be specified if power scaling is not required.
C             Has no meaning if NORMSP.LE.0 or SPOWER=0.
C             Default: SPDIST=1.
C     SPEXP=real, SPEXPT=real... Exponential scaling of seismograms with
C             respect to time.  The amplitude scale is multiplied by the
C             factor of EXP(SPEXP*(t-SPEXPT)).
C             Example 1:  Exponential scaling with SPEXP=pi*F/Q, where
C               F is the dominant frequency and Q is the quality factor,
C               may roughly compensate for attenuation when plotting
C               late arrivals.
C             Example 2:  Exponential scaling with SPEXP equal to half
C               the sum of the positive Lyapunov exponents may
C               compensate for large geometrical spreading when plotting
C               late arrivals.
C             Default: SPEXP=0., SPEXPT=0.
C     SPEXP1=real, SPEXP2=real, ..., SPEXP9=real... Exponential scaling
C             set individually for optional input GSE files given by
C             parameters SS1, SS2, ..., SS9, respectively.
C             Defaults: SPEXP1=SPEXP, SPEXP2=SPEXP, ..., SPEXP9=SPEXP
C Other data to control plotting:
C     SPTEXT1='string'... Text to be written above the left-hand top
C             corner of the frame.
C             Default: SPTEXT1=' '
C     SPTEXT2='string'... Text to be written above the right-hand top
C             corner of the frame.
C             Default: SPTEXT2=' '
C     SPTEXT3='string'... Text to be written below the frame.
C             Default: SPTEXT3=' '
C     SPTEXT4='string'... Text to be written below text SPTEXT3.
C             Default: SPTEXT4=' '
C     SPCHRH=real... Character height in cm.
C             Default: SPCHRH=0.4
C     SPHIWI=real... Width of the highlighted area when plotting travel
C             times.
C             Default: SPHIWI=XD,
C               where XD is the average distance between seismograms,
C               i.e. XD=SPXLEN/(NREC+1), where NREC is the number of
C               receivers in file REC, see the amplitude scaling.
C
C                                                     
C Input GSE files SS* with the seismograms to plot:
C     File in the GSE data exchange format, see the description in file
C     'gse.for'.
C     The 'sp.for' program is looking in the comment lines of the
C     waveform identification section for the source name identified by
C     string 'NAMESRC', and for the hypocentral coordinates identified
C     by strings 'X1SRC', 'X2SRC' and 'X3SRC'.
C     Description of format GSE
C
C                                                     
C Input file SRC with the hypocentral coordinates:
C (1) /
C     None to 20 character strings terminated by a slash.
C (2) 'SRCNAME',X1SRC,X2SRC,X3SRC,/
C     'SRCNAME'... Name of the source.  Not used.
C     X1SRC,X2SRC,X3SRC... Coordinates of the hypocentre.
C     Default: X1SRC=0., X2SRC=0., X3SRC=0.
C
C                                                     
C Input file REC containing receiver coordinates:
C (1) Several strings terminated by / (a slash).
C             The simplest way is to submit just the /.
C (2) Several times (2.1):
C (2.1) 'NAMER(IR)',X1REC(IR),X2REC(IR),X3REC(IR),/
C     'NAMER(IR)'... String reserved for the name of the receiver.
C             No meaning here, anything in apostrophes may be submitted.
C     X1REC(IR),X2REC(IR),X3REC(IR)... Coordinates of the receiver.
C             The coordinates need not be present in the file.  It may
C             thus be comfortable to omit them if preparing the list for
C             the selection of particular receivers for seismogram
C             plotting.
C     Default: X1REC(IR)=0, X2REC(IR)=0, X3REC(IR)=0.
C (3) / (a slash).
C
C                                                  
C Input file SPHILI containing travel times to highlight:
C (1) Several strings terminated by / (a slash).
C             The simplest way is to submit just the /.
C (2) Several times (2.1):
C (2.1) 'NAMER(IR)',X1REC(IR),X2REC(IR),X3REC(IR),TT,TTERR,K1,K2,/
C     'NAMER(IR)',X1REC(IR),X2REC(IR),X3REC(IR)... Same meaning as in
C             the receiver file above.
C     TT...   Travel time.  If given, it will be plotted as horizontal
C             line of width SPHIWI and colour K1.  It is plotted after
C             the frame and the corresponding error bar, but before
C             seismograms.
C     TTERR...Travel time error.  If given, it will be plotted as
C             a solid rectangle of width SPHIWI and colour K2. It is
C             plotted after the frame but before the corresponding
C             travel time and before seismograms.
C     K1...   Colour to plot the travel time.  The travel time is not
C             plotted if K1 is not positive.
C     K2...   Colour to plot the error bar.  The error bar is not
C             plotted if K2 is negative.  If K2=0 (white), the contour
C             of the rectangle is plotted in colour K1.
C     Default: X1REC(IR)=0, X2REC(IR)=0, X3REC(IR)=0., K1=KOLORTD,
C             K2=KOLORTD
C (3) / (a slash).
C
C.......................................................................
C                                                
C This Fortran77 file consists of the following external procedures:
C     SP...   Main program to read and plot the seismograms.
C             SP
C     FRAME...Subroutine called by the main to plot the rectangular
C             frame around the seismograms and supplement it with simple
C             descriptions.
C             FRAME
C
C Other external procedures required:
C     RGSE1,RGSE2... Subroutines of the Fortran 77 file 'gse.for'
C             (package MODEL), designed to read seismograms in the GSE
C             data exchange format.
C     PLOTS,PLOT,SYMBOL,NUMBER... CALCOMP plotting subroutines.  For
C             example, Fortran 77 routines of file 'calcops.for'
C             (package MODEL) may be used to generate seismogram plots
C             in the PostScript files.
C
C=======================================================================
C
C     
C
C Common block /RAMC/:
      INCLUDE 'ram.inc'
C     ram.inc
      INTEGER IRAM(MRAM)
      EQUIVALENCE (RAM,IRAM)
C
C     Allocation of working arrays:
      INTEGER MPTS,MSS
      PARAMETER (MPTS=3000,MSS=MRAM-3*MPTS)
      CHARACTER*6 PTS(MPTS)
      COMMON/PTSC/PTS
C
C-----------------------------------------------------------------------
C
C     External functions and subroutines:
      EXTERNAL LENGTH
      INTEGER  LENGTH
C
C     Input and output data filenames and logical unit numbers:
      INTEGER LU,LUPAR
      PARAMETER (LU=1,LUPAR=2)
      CHARACTER*80 FILSEP,FILPAR,FILPTS,FILSRC,FILREC,FILFTT,FILHIL
      CHARACTER*80 FILOLD(0:9),FILESS(0:9),FILEPS(3)
C
C     Storing seismograms in memory
      INTEGER IFILO(0:9),IFILE(0:9),ISSRAM(0:10)
C     IFILO(I:I)... Index of the data corresponding to file FILOLD(I:I).
C     IFILE(I:I)... Index of the data corresponding to file FILESS(I:I).
C     ISSRAM(IFILE:IFILE)... End of IFILEth data stored in RAM.
C
C     Parameters and small working arrays:
      REAL UNDEF
      PARAMETER (UNDEF=-999999.)
      INTEGER KOLOR(0:9),KOMP(0:9,3)
      INTEGER ISS,ISP
      REAL SPXMIN(0:9),SPXMAX(0:9),SPYMIN(0:9),SPYMAX(0:9)
      REAL SPAMP(0:9,3),SPOWER(0:9),SPEXP(0:9)
      REAL SPAMPX(3),XPTS(4),YPTS(4)
C
      CHARACTER*1  CHAN,TEXT
      CHARACTER*6  NAMSRC,NAMREC
C
C     Line of optional SEP parameter file or comments in the GSE file
      CHARACTER*80 LINE
C
C     Data specifying labels in the plot:
      CHARACTER*80 TEXT1,TEXT2,TEXT3,TEXT4
      CHARACTER*20 KXTEXT(5)
C
C     Lists of point coordinates, sources and receivers:
C     INTEGER NPTS,NSRC,NREC
C
      DATA FILESS/10*' '/
      DATA KXTEXT/' ','X1','X1','X3','HYPOCENTRAL DISTANCE'/
C
C.......................................................................
C
C     Reading name of SEP file with input data:
      FILSEP=' '
      WRITE(*,'(A)') '+SP: Enter input filename: '
      READ(*,*) FILSEP
      IF (FILSEP.EQ.' ') THEN
C       SP-01
        CALL ERROR('SP-01: SEP file not given')
C       Input file in the form of the SEP (Stanford Exploration Project)
C       parameter or history file must be specified.
C       There is no default filename.
      ENDIF
C
C     Reading all data from the SEP file into the memory:
      CALL RSEP1(LU,FILSEP)
      WRITE(*,'(A)') '+SP: Working...            '
C
C     Loop over SP executions in optional SP parameter file:
      CALL RSEP3T('SPPAR',FILPAR,' ')
      IF(FILPAR.NE.' ') THEN
        OPEN(LUPAR,FILE=FILPAR,STATUS='OLD')
        NSSRAM=0
        ISSRAM(0)=0
      END IF
  100 CONTINUE
      IF(FILPAR.NE.' ') THEN
C       Loop over lines the SP parameter file
  110   CONTINUE
          READ(LUPAR,'(A)',END=999) LINE
          CALL RSEP2(LINE)
          I=INDEX(LINE,'#')
          IF(I.EQ.0) THEN
            I=LENGTH(LINE)
          END IF
          I=INDEX(LINE(1:I),'sp:')
          IF(I.GT.1) THEN
            IF(LINE(I-1:I-1).NE.' ') THEN
              I=0
            END IF
          END IF
        IF(I.EQ.0) GO TO 110
C       SP execution is prescribed at the current positions in the SP
C       parameter file.
C       Copying filenames corresponding to the data stored in RAM:
        DO 111 I2=0,9
          FILOLD(I2)=FILESS(I2)
          IFILO (I2)=IFILE (I2)
  111   CONTINUE
      END IF
C
C     Input and output filenames:
      CALL RSEP3T('SS' ,FILESS(0),'ss.gse')
      CALL RSEP3T('SS1',FILESS(1),' ')
      CALL RSEP3T('SS2',FILESS(2),' ')
      CALL RSEP3T('SS3',FILESS(3),' ')
      CALL RSEP3T('SS4',FILESS(4),' ')
      CALL RSEP3T('SS5',FILESS(5),' ')
      CALL RSEP3T('SS6',FILESS(6),' ')
      CALL RSEP3T('SS7',FILESS(7),' ')
      CALL RSEP3T('SS8',FILESS(8),' ')
      CALL RSEP3T('SS9',FILESS(9),' ')
      CALL RSEP3T('SP1',FILEPS(1),'ss1.ps')
      CALL RSEP3T('SP2',FILEPS(2),'ss2.ps')
      CALL RSEP3T('SP3',FILEPS(3),'ss3.ps')
C###
      IF(FILPAR.EQ.' ') THEN
        ISS0=0
      ELSE
        DO 129 I2=0,9
          IF(FILOLD(I2).NE.' ') THEN
            DO 121 I1=0,9
              IF(FILESS(I1).EQ.FILOLD(I2)) THEN
                GO TO 128
              END IF
  121       CONTINUE
C             Removing seismograms from the memory
              I=ISSRAM(IFILE(I2))-ISSRAM(IFILE(I2)-1)
              DO 122 I1=ISSRAM(IFILE(I2))+1,ISSRAM(NSSRAM)
                RAM(I1-I)=RAM(I1)
  122         CONTINUE
              DO 123 I1=IFILE(I2),NSSRAM
                ISSRAM(I1-1)=ISSRAM(I1)-I
  123         CONTINUE
              NSSRAM=NSSRAM-1
              DO 124 I1=I2+1,9
                IF(FILOLD(I1).EQ.FILOLD(I2)) THEN
                  FILOLD(I1)=' '
                END IF
  124         CONTINUE
              FILOLD(I2)=' '
  128       CONTINUE
          END IF
  129   CONTINUE
        DO 139 I2=0,9
          IF(FILESS(I2).NE.' ') THEN
            DO 131 I1=0,9
              IF(FILESS(I2).EQ.FILOLD(I1)) THEN
                IFILE(I2)=IFILO(I1)
                GO TO 138
              END IF
  131       CONTINUE
            DO 132 I1=0,I2-1
              IF(FILESS(I2).EQ.FILESS(I1)) THEN
                IFILE(I2)=IFILE(I1)
                GO TO 138
              END IF
  132       CONTINUE
C             Reading seismograms into the memory
              WRITE(*,'(2A)') '+SP: Reading ',FILESS(I2)(1:66)
              OPEN(LU,FILE=FILESS(I2),STATUS='OLD')
              CALL RGSE1(LU,TEXT)
              ISS0=ISSRAM(NSSRAM)+1
  133         CONTINUE
                CALL RGSE2(LU,NAMREC,CHAN,I,X1R,X2R,X3R,T0,TD,
     *                                      NSS,MSS-ISS0-22,RAM(ISS0+1))
                IF(NSS.LE.-1) THEN
C                 End of the GSE file
                  GO TO 137
                END IF
  134           CONTINUE
                  CALL RGSE2C(LINE,*135)
                  CALL RSEP2(LINE)
                GO TO 134
  135           CONTINUE
                CALL RSEP3T('NAMESRC',NAMSRC,' ')
                CALL RSEP3R('X1SRC',X1S,0.)
                CALL RSEP3R('X2SRC',X2S,0.)
                CALL RSEP3R('X3SRC',X3S,0.)
                CALL WRAM2(ISS0,NAMSRC,X1S,X2S,X3S,
     *                          NAMREC,X1R,X2R,X3R,I,T0,TD,NSS,IRAM,RAM)
              GO TO 133
  137         CONTINUE
              NSSRAM=NSSRAM+1
              ISSRAM(NSSRAM)=ISS0
              IFILE(I2)=NSSRAM
              CLOSE(LU)
  138       CONTINUE
          END IF
  139   CONTINUE
      END IF
C^^^
C     Reading lists of point coordinates, sources and receivers:
C     Point coordinates
      NPTS=0
      CALL RSEP3T('PTS',FILPTS,' ')
      IF(FILPTS.NE.' ') THEN
        OPEN(LU,FILE=FILPTS,STATUS='OLD')
        READ(LU,*) (TEXT,I=1,20)
        DO 11 I=1,MPTS
          I0=MSS+3*I
          PTS(I)='$$$$$$'
          RAM(I0-2)=0.
          RAM(I0-1)=0.
          RAM(I0  )=0.
          READ(LU,*,END=12) PTS(I),RAM(I0-2),RAM(I0-1),RAM(I0)
          IF(PTS(I).EQ.'$$$$$$') THEN
            GO TO 12
          END IF
   11   CONTINUE
C         SP-02
          CALL ERROR('SP-02: Array dimension MPTS small for points')
   12   CONTINUE
        NPTS=I-1
        CLOSE(LU)
      END IF
C     Sources
      NSRC=0
      CALL RSEP3T('SRC',FILSRC,' ')
      IF(FILSRC.NE.' ') THEN
        OPEN(LU,FILE=FILSRC,STATUS='OLD')
        READ(LU,*) (TEXT,I=1,20)
        DO 13 I=NPTS+1,MPTS
          I0=MSS+3*I
          PTS(I)='$$$$$$'
          RAM(I0-2)=0.
          RAM(I0-1)=0.
          RAM(I0  )=0.
          READ(LU,*,END=14) PTS(I),RAM(I0-2),RAM(I0-1),RAM(I0)
          IF(PTS(I).EQ.'$$$$$$') THEN
            GO TO 14
          END IF
   13   CONTINUE
C         SP-03
          CALL ERROR('SP-03: Array dimension MPTS small for sources')
   14   CONTINUE
        NSRC=I-1
        CLOSE(LU)
      END IF
C     Receivers
      NREC=0
      RECNUM=999.
      CALL RSEP3T('REC',FILREC,' ')
      IF(FILREC.NE.' ') THEN
        OPEN(LU,FILE=FILREC,STATUS='OLD')
        READ(LU,*) (TEXT,I=1,20)
        DO 15 I=NPTS+NSRC+1,MPTS
          I0=MSS+3*I
          PTS(I)='$$$$$$'
          RAM(I0-2)=0.
          RAM(I0-1)=0.
          RAM(I0  )=0.
          READ(LU,*,END=16) PTS(I),RAM(I0-2),RAM(I0-1),RAM(I0)
          IF(PTS(I).EQ.'$$$$$$') THEN
            GO TO 16
          END IF
   15   CONTINUE
C         SP-04
          CALL ERROR('SP-04: Array dimension MPTS small for points')
   16   CONTINUE
        NREC=I-NPTS-NSRC-1
        RECNUM=FLOAT(NREC)
        CLOSE(LU)
      END IF
C
C     Components:
      CALL RSEP3I('KOMP1 ',KOMP(0,1),1)
      CALL RSEP3I('KOMP11',KOMP(1,1),KOMP(0,1))
      CALL RSEP3I('KOMP12',KOMP(2,1),KOMP(0,1))
      CALL RSEP3I('KOMP13',KOMP(3,1),KOMP(0,1))
      CALL RSEP3I('KOMP14',KOMP(4,1),KOMP(0,1))
      CALL RSEP3I('KOMP15',KOMP(5,1),KOMP(0,1))
      CALL RSEP3I('KOMP16',KOMP(6,1),KOMP(0,1))
      CALL RSEP3I('KOMP17',KOMP(7,1),KOMP(0,1))
      CALL RSEP3I('KOMP18',KOMP(8,1),KOMP(0,1))
      CALL RSEP3I('KOMP19',KOMP(9,1),KOMP(0,1))
      CALL RSEP3I('KOMP2 ',KOMP(0,2),2)
      CALL RSEP3I('KOMP21',KOMP(1,2),KOMP(0,2))
      CALL RSEP3I('KOMP22',KOMP(2,2),KOMP(0,2))
      CALL RSEP3I('KOMP23',KOMP(3,2),KOMP(0,2))
      CALL RSEP3I('KOMP24',KOMP(4,2),KOMP(0,2))
      CALL RSEP3I('KOMP25',KOMP(5,2),KOMP(0,2))
      CALL RSEP3I('KOMP26',KOMP(6,2),KOMP(0,2))
      CALL RSEP3I('KOMP27',KOMP(7,2),KOMP(0,2))
      CALL RSEP3I('KOMP28',KOMP(8,2),KOMP(0,2))
      CALL RSEP3I('KOMP29',KOMP(9,2),KOMP(0,2))
      CALL RSEP3I('KOMP3 ',KOMP(0,3),3)
      CALL RSEP3I('KOMP31',KOMP(1,3),KOMP(0,3))
      CALL RSEP3I('KOMP32',KOMP(2,3),KOMP(0,3))
      CALL RSEP3I('KOMP33',KOMP(3,3),KOMP(0,3))
      CALL RSEP3I('KOMP34',KOMP(4,3),KOMP(0,3))
      CALL RSEP3I('KOMP35',KOMP(5,3),KOMP(0,3))
      CALL RSEP3I('KOMP36',KOMP(6,3),KOMP(0,3))
      CALL RSEP3I('KOMP37',KOMP(7,3),KOMP(0,3))
      CALL RSEP3I('KOMP38',KOMP(8,3),KOMP(0,3))
      CALL RSEP3I('KOMP39',KOMP(9,3),KOMP(0,3))
C
C     Colours:
      CALL RSEP3I('KOLOR ',KOLOR(0),1)
      CALL RSEP3I('KOLOR1',KOLOR(1),2)
      CALL RSEP3I('KOLOR2',KOLOR(2),3)
      CALL RSEP3I('KOLOR3',KOLOR(3),4)
      CALL RSEP3I('KOLOR4',KOLOR(4),5)
      CALL RSEP3I('KOLOR5',KOLOR(5),6)
      CALL RSEP3I('KOLOR6',KOLOR(6),7)
      CALL RSEP3I('KOLOR7',KOLOR(7),8)
      CALL RSEP3I('KOLOR8',KOLOR(8),9)
      CALL RSEP3I('KOLOR9',KOLOR(9),10)
      CALL RSEP3I('KOLORTT',KOLORT,1)
      CALL RSEP3I('KOLORTD',KOLORD,0)
C
C     Initial values for plotting frame:
C     Time axis:
      CALL RSEP3R('SPTMIN',SPTMIN, 0.)
      CALL RSEP3R('SPTMAX',SPTMAX, 1.)
      CALL RSEP3R('SPTLEN',SPTLEN,10.)
      CALL RSEP3R('SPTDIV',SPTDIV, 1.)
      CALL RSEP3R('SPTSUB',SPTSUB, 1.)
      CALL RSEP3R('SPVRED',SPVRED, 0.)
C     Distance axis:
      CALL RSEP3I('KODESP',KODESP, 0 )
      CALL RSEP3R('SPXLEN',SPXLEN,RECNUM+1.)
      CALL RSEP3R('SPXDIV',SPXDIV, 1.)
      CALL RSEP3R('SPXSUB',SPXSUB, 1.)
      CALL RSEP3R('SPXMIN ',SPXMIN(0), 0.)
      CALL RSEP3R('SPXMIN1',SPXMIN(1),SPXMIN(0))
      CALL RSEP3R('SPXMIN2',SPXMIN(2),SPXMIN(0))
      CALL RSEP3R('SPXMIN3',SPXMIN(3),SPXMIN(0))
      CALL RSEP3R('SPXMIN4',SPXMIN(4),SPXMIN(0))
      CALL RSEP3R('SPXMIN5',SPXMIN(5),SPXMIN(0))
      CALL RSEP3R('SPXMIN6',SPXMIN(6),SPXMIN(0))
      CALL RSEP3R('SPXMIN7',SPXMIN(7),SPXMIN(0))
      CALL RSEP3R('SPXMIN8',SPXMIN(8),SPXMIN(0))
      CALL RSEP3R('SPXMIN9',SPXMIN(9),SPXMIN(0))
      CALL RSEP3R('SPXMAX ',SPXMAX(0),RECNUM+1.)
      CALL RSEP3R('SPXMAX1',SPXMAX(1),SPXMAX(0))
      CALL RSEP3R('SPXMAX2',SPXMAX(2),SPXMAX(0))
      CALL RSEP3R('SPXMAX3',SPXMAX(3),SPXMAX(0))
      CALL RSEP3R('SPXMAX4',SPXMAX(4),SPXMAX(0))
      CALL RSEP3R('SPXMAX5',SPXMAX(5),SPXMAX(0))
      CALL RSEP3R('SPXMAX6',SPXMAX(6),SPXMAX(0))
      CALL RSEP3R('SPXMAX7',SPXMAX(7),SPXMAX(0))
      CALL RSEP3R('SPXMAX8',SPXMAX(8),SPXMAX(0))
      CALL RSEP3R('SPXMAX9',SPXMAX(9),SPXMAX(0))
      CALL RSEP3R('SPYMIN ',SPYMIN(0), 0.)
      CALL RSEP3R('SPYMIN1',SPYMIN(1),SPYMIN(0))
      CALL RSEP3R('SPYMIN2',SPYMIN(2),SPYMIN(0))
      CALL RSEP3R('SPYMIN3',SPYMIN(3),SPYMIN(0))
      CALL RSEP3R('SPYMIN4',SPYMIN(4),SPYMIN(0))
      CALL RSEP3R('SPYMIN5',SPYMIN(5),SPYMIN(0))
      CALL RSEP3R('SPYMIN6',SPYMIN(6),SPYMIN(0))
      CALL RSEP3R('SPYMIN7',SPYMIN(7),SPYMIN(0))
      CALL RSEP3R('SPYMIN8',SPYMIN(8),SPYMIN(0))
      CALL RSEP3R('SPYMIN9',SPYMIN(9),SPYMIN(0))
      CALL RSEP3R('SPYMAX ',SPYMAX(0), 0.)
      CALL RSEP3R('SPYMAX1',SPYMAX(1),SPYMAX(0))
      CALL RSEP3R('SPYMAX2',SPYMAX(2),SPYMAX(0))
      CALL RSEP3R('SPYMAX3',SPYMAX(3),SPYMAX(0))
      CALL RSEP3R('SPYMAX4',SPYMAX(4),SPYMAX(0))
      CALL RSEP3R('SPYMAX5',SPYMAX(5),SPYMAX(0))
      CALL RSEP3R('SPYMAX6',SPYMAX(6),SPYMAX(0))
      CALL RSEP3R('SPYMAX7',SPYMAX(7),SPYMAX(0))
      CALL RSEP3R('SPYMAX8',SPYMAX(8),SPYMAX(0))
      CALL RSEP3R('SPYMAX9',SPYMAX(9),SPYMAX(0))
C     Characters:
      CALL RSEP3T('SPTEXT1',TEXT1,' ')
      CALL RSEP3T('SPTEXT2',TEXT2,' ')
      CALL RSEP3T('SPTEXT3',TEXT3,' ')
      CALL RSEP3T('SPTEXT4',TEXT4,' ')
      CALL RSEP3R('SPCHRH',SPCHRH, 0.4)
C
C     Amplitude scaling:
      CALL RSEP3I('NORMSP',NORMSP,0)
      CALL RSEP3R('SPAMP ',SPAMP(0,1),1.)
      CALL RSEP3R('SPAMP1',SPAMP(1,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP2',SPAMP(2,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP3',SPAMP(3,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP4',SPAMP(4,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP5',SPAMP(5,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP6',SPAMP(6,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP7',SPAMP(7,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP8',SPAMP(8,1),SPAMP(0,1))
      CALL RSEP3R('SPAMP9',SPAMP(9,1),SPAMP(0,1))
      CALL RSEP3R('SPAMPX1',SPAMPX(1),1.)
      CALL RSEP3R('SPAMPX2',SPAMPX(2),1.)
      CALL RSEP3R('SPAMPX3',SPAMPX(3),1.)
      DO 18 I2=3,1,-1
        DO 17 I1=0,9
          SPAMP(I1,I2)=SPAMP(I1,1)*SPAMPX(I2)
   17   CONTINUE
   18 CONTINUE
      CALL RSEP3R('SPDIST' ,SPDIST,1.)
      CALL RSEP3R('SPOWER' ,SPOWER(0),0.)
      CALL RSEP3R('SPOWER1',SPOWER(1),SPOWER(0))
      CALL RSEP3R('SPOWER2',SPOWER(2),SPOWER(0))
      CALL RSEP3R('SPOWER3',SPOWER(3),SPOWER(0))
      CALL RSEP3R('SPOWER4',SPOWER(4),SPOWER(0))
      CALL RSEP3R('SPOWER5',SPOWER(5),SPOWER(0))
      CALL RSEP3R('SPOWER6',SPOWER(6),SPOWER(0))
      CALL RSEP3R('SPOWER7',SPOWER(7),SPOWER(0))
      CALL RSEP3R('SPOWER8',SPOWER(8),SPOWER(0))
      CALL RSEP3R('SPOWER9',SPOWER(9),SPOWER(0))
      CALL RSEP3R('SPEXP ',SPEXP(0),0.)
      CALL RSEP3R('SPEXP1',SPEXP(1),SPEXP(0))
      CALL RSEP3R('SPEXP2',SPEXP(2),SPEXP(0))
      CALL RSEP3R('SPEXP3',SPEXP(3),SPEXP(0))
      CALL RSEP3R('SPEXP4',SPEXP(4),SPEXP(0))
      CALL RSEP3R('SPEXP5',SPEXP(5),SPEXP(0))
      CALL RSEP3R('SPEXP6',SPEXP(6),SPEXP(0))
      CALL RSEP3R('SPEXP7',SPEXP(7),SPEXP(0))
      CALL RSEP3R('SPEXP8',SPEXP(8),SPEXP(0))
      CALL RSEP3R('SPEXP9',SPEXP(9),SPEXP(0))
      CALL RSEP3R('SPEXPT',SPEXPT,0.)
C
      SCY= SPTLEN/(SPTMAX-SPTMIN)
      XA = SPXMAX(0)-SPXMIN(0)
      YA = SPYMAX(0)-SPYMIN(0)
C
C     Higlighting given areas (e.g., travel times with error bars):
      CALL RSEP3T('FTT'   ,FILFTT,' ')
      CALL RSEP3T('SPHILI',FILHIL,' ')
      CALL RSEP3R('SPHIWI',SPHIWI,SPXLEN/(RECNUM+1.))
C
C.......................................................................
C
C     Loop over 3 seismogram files:
      DO 99 ISP=1,3
        IF(FILEPS(ISP).NE.' ') THEN
          WRITE(*,'(2A)') '+SP: Plotting ',FILEPS(ISP)(1:65)
C
C         Initialization of CALCOMP:
          CALL PLOTN(FILEPS(ISP),0)
          CALL PLOTS(0,0,0)
C
C         Plotting frame:
          CALL NEWPEN(1)
          IX = KODESP
          IT = -1
          IF(KODESP.GE.3)  IX=1
          IF(SPXDIV.LT.0.) IX=0
          IF(SPTDIV.LT.0.) IT=0
          CALL FRAME(SPXLEN,SPTLEN,ABS(SPXDIV),SPXSUB,
     *                             ABS(SPTDIV),SPTSUB,0,IX,IT,SPCHRH,
     *               SPXMIN(0),SPXMAX(0),KXTEXT(KODESP+1),
     *               SPYMIN(0),SPYMAX(0),'X2',
     *               SPTMIN,SPTMAX,'TIME',
     *               TEXT1,0,0.,0,TEXT2,0,0.,0,TEXT3,TEXT4)
C
C         Higlighting travel times of file FTT:
          IF(FILFTT.NE.' ') THEN
            OPEN(LU,FILE=FILFTT,STATUS='OLD')
            READ(LU,*) (TEXT,I=1,20)
C           Loop for areas to highlight:
   20       CONTINUE
              NAMSRC='$$$$$$'
              T0=0.
              TD=0.
              READ(LU,*,END=29) NAMSRC,NAMREC,T0,TD
              IF(NAMSRC.EQ.'$$$$$$') THEN
C               End of travel-time file
                GO TO 29
              END IF
C             Selecting the receiver:
              IF(FILREC.NE.' ') THEN
C               Loop for receivers
                DO 21 I=NPTS+NSRC+1,NPTS+NSRC+NREC
                  IF(PTS(I).EQ.NAMREC) THEN
                    IREC=I-NPTS-NSRC
                    GO TO 22
                  END IF
   21           CONTINUE
                GO TO 20
              END IF
   22         CONTINUE
C             Selecting the source:
              IF(FILSRC.NE.' ') THEN
C               Loop for sources
                DO 23 I=NPTS+1,NPTS+NSRC
                  IF(PTS(I).EQ.NAMSRC) GO TO 24
   23           CONTINUE
                GO TO 20
              END IF
   24         CONTINUE
C             Finding the receiver coordinates:
              DO 25 I=1,NPTS
                IF(PTS(I).EQ.NAMREC) THEN
                  I0=MSS+3*I
                  X1R=RAM(I0-2)
                  X2R=RAM(I0-1)
                  X3R=RAM(I0)
                  GO TO 26
                END IF
   25         CONTINUE
C               SP-05
                LINE='SP-05: Receiver '//NAMREC(1:LENGTH(NAMREC))
     *                                 //' not found in file PTS'
                CALL ERROR(LINE)
C               If file FTT with travel times is given, file PTS has
C               to contain all receiver names of file REC (if REC is
C               specified) or all receiver names of file FTT (if
C               REC=' ').
   26         CONTINUE
C             Finding the source coordinates:
              DO 27 I=1,NPTS
                IF(PTS(I).EQ.NAMSRC) THEN
                  I0=MSS+3*I
                  X1S=RAM(I0-2)
                  X2S=RAM(I0-1)
                  X3S=RAM(I0)
                  GO TO 28
                END IF
   27         CONTINUE
C               SP-06
                LINE='SP-06: Source '//NAMSRC(1:LENGTH(NAMSRC))
     *                               //' not found in file PTS'
                CALL ERROR(LINE)
C               If file FTT with travel times is given, file PTS has
C               to contain all receiver names of file REC (if REC is
C               specified) or all receiver names of file FTT (if
C               REC=' ').
   28         CONTINUE
C             Reduction of the travel time
              IF(SPVRED.NE.0.) THEN
                RR=SQRT((X1R-X1S)**2+(X2R-X2S)**2+(X3R-X3S)**2)
                T0=T0-RR/SPVRED
              END IF
              IF(KODESP.EQ.0) THEN
                IF(FILREC.EQ.' ') THEN
C                 SP-07
                  CALL ERROR('SP-07: No receiver file specified')
C                 For KODESP=0, filename REC must be specified in the
C                 input data.
                END IF
                X=SPXLEN*(FLOAT(IREC)-SPXMIN(0))
     *                  /(SPXMAX(0)-SPXMIN(0))
              ELSE IF(KODESP.EQ.1.OR.KODESP.EQ.2) THEN
                X=SPXLEN*((X1R-SPXMIN(0))*XA+(X2R-SPYMIN(0))*YA)
     *                  /(XA*XA+YA*YA)
              ELSE IF(KODESP.EQ.3) THEN
                X=SPXLEN*(X3R-SPXMIN(0))/XA
              ELSE
                X=SPXLEN*(RR-SPXMIN(0))/XA
              END IF
C             Plotting the highlight:
              IF(KOLORD.GE.0) THEN
                XPTS(1)=X-0.5*SPHIWI
                XPTS(2)=X+0.5*SPHIWI
                XPTS(3)=X+0.5*SPHIWI
                XPTS(4)=X-0.5*SPHIWI
                YPTS(1)=SCY*(T0-TD-SPTMIN)
                YPTS(2)=SCY*(T0-TD-SPTMIN)
                YPTS(3)=SCY*(T0+TD-SPTMIN)
                YPTS(4)=SCY*(T0+TD-SPTMIN)
                IF(KOLORD.GT.0) THEN
                  CALL NEWPEN(KOLORD)
                  CALL FILL(XPTS,YPTS,4)
                ELSE
                  CALL NEWPEN(KOLORT)
                  CALL PLOT(XPTS(1),YPTS(1),3)
                  CALL PLOT(XPTS(2),YPTS(2),2)
                  CALL PLOT(XPTS(3),YPTS(3),2)
                  CALL PLOT(XPTS(4),YPTS(4),2)
                  CALL PLOT(XPTS(1),YPTS(1),2)
                END IF
              END IF
              IF(KOLORT.GT.0) THEN
                CALL NEWPEN(KOLORT)
                CALL PLOT(X-0.5*SPHIWI,SCY*(T0-SPTMIN),3)
                CALL PLOT(X+0.5*SPHIWI,SCY*(T0-SPTMIN),2)
              END IF
            GO TO 20
   29       CONTINUE
C           Closing highlighting file
            CLOSE(LU)
          END IF
C
C         Higlighting times given in file SPHILI:
          IF(FILHIL.NE.' ') THEN
            OPEN(LU,FILE=FILHIL,STATUS='OLD')
            READ(LU,*) (TEXT,I=1,20)
C           Loop for areas to highlight:
   30       CONTINUE
              NAMREC='$$$$$$'
              T0=UNDEF
              TD=UNDEF
              K1=KOLORT
              K2=KOLORD
              READ(LU,*,END=39) NAMREC,X1R,X2R,X3R,T0,TD,K1,K2
              IF(NAMREC.EQ.'$$$$$$') THEN
                GO TO 39
              END IF
C             Selecting the receiver:
              IF(FILREC.NE.' ') THEN
C               Loop for receivers
                DO 31 I=NPTS+NSRC+1,NPTS+NSRC+NREC
                  IF(PTS(I).EQ.NAMREC) THEN
                    IREC=I-NPTS-NSRC
                    GO TO 32
                  END IF
   31           CONTINUE
                GO TO 30
              END IF
   32         CONTINUE
C             Updating the coordinates:
              IF(FILPTS.NE.' '.AND.NAMREC.NE.' ') THEN
C               Receiver coordinates
                DO 33 I=1,NPTS
                  IF(PTS(I).EQ.NAMREC) THEN
                    I0=MSS+3*I
                    X1R=RAM(I0-2)
                    X2R=RAM(I0-1)
                    X3R=RAM(I0)
                    GO TO 34
                  END IF
   33           CONTINUE
C                 SP-08
                  LINE='SP-08: Receiver '//NAMREC(1:LENGTH(NAMREC))
     *                                   //' not found in file PTS'
                  CALL ERROR(LINE)
C                 If file PTS with the coordinates of points is given
C                 and file SPHILI contains receiver names, file PTS has
C                 to contain all receiver names of file REC (if REC is
C                 specified) or all receiver names of file SPHILI (if
C                 REC=' ').
   34           CONTINUE
              END IF
              IF(SPVRED.NE.0.) THEN
                IF(FILSRC.EQ.' ') THEN
                  X1S=0.
                  X2S=0.
                  X3S=0.
                ELSE
                  I0=MSS+3*NPTS
                  X1S=RAM(I0+1)
                  X2S=RAM(I0+2)
                  X3S=RAM(I0+3)
                END IF
                RR=SQRT((X1R-X1S)**2+(X2R-X2S)**2+(X3R-X3S)**2)
                T0=T0-RR/SPVRED
              END IF
              IF(KODESP.EQ.0) THEN
                IF(FILREC.EQ.' ') THEN
C                 SP-09
                  CALL ERROR('SP-09: No receiver file specified')
C                 For KODESP=0, filename REC must be specified in the
C                 input data.
                END IF
                X=SPXLEN*(FLOAT(IREC)-SPXMIN(0))
     *                  /(SPXMAX(0)-SPXMIN(0))
              ELSE IF(KODESP.EQ.1.OR.KODESP.EQ.2) THEN
                X=SPXLEN*((X1R-SPXMIN(0))*XA+(X2R-SPYMIN(0))*YA)
     *                  /(XA*XA+YA*YA)
              ELSE IF(KODESP.EQ.3) THEN
                X=SPXLEN*(X3R-SPXMIN(0))/XA
              ELSE
                X=SPXLEN*(RR-SPXMIN(0))/XA
              END IF
C             Plotting the highlight:
              IF(T0.EQ.UNDEF) THEN
                K1=-1
                K2=-1
              END IF
              IF(TD.EQ.UNDEF) THEN
                K2=-1
              END IF
              IF(K2.GE.0) THEN
                XPTS(1)=X-0.5*SPHIWI
                XPTS(2)=X+0.5*SPHIWI
                XPTS(3)=X+0.5*SPHIWI
                XPTS(4)=X-0.5*SPHIWI
                YPTS(1)=SCY*(T0-TD-SPTMIN)
                YPTS(2)=SCY*(T0-TD-SPTMIN)
                YPTS(3)=SCY*(T0+TD-SPTMIN)
                YPTS(4)=SCY*(T0+TD-SPTMIN)
                IF(K2.GT.0) THEN
                  CALL NEWPEN(K2)
                  CALL FILL(XPTS,YPTS,4)
                ELSE
                  CALL NEWPEN(K1)
                  CALL PLOT(XPTS(1),YPTS(1),3)
                  CALL PLOT(XPTS(2),YPTS(2),2)
                  CALL PLOT(XPTS(3),YPTS(3),2)
                  CALL PLOT(XPTS(4),YPTS(4),2)
                END IF
              END IF
              IF(K1.GT.0) THEN
                CALL NEWPEN(K1)
                CALL PLOT(X-0.5*SPHIWI,SCY*(T0-SPTMIN),3)
                CALL PLOT(X+0.5*SPHIWI,SCY*(T0-SPTMIN),2)
              END IF
            GO TO 30
   39       CONTINUE
C           Closing highlighting file
            CLOSE(LU)
          END IF
C
C         Plotting seismograms:
C         Loop for GSE files
          DO 90 ISS=0,9
            IF(FILESS(ISS).NE.' ') THEN
              CALL NEWPEN(KOLOR(ISS))
              XA=SPXMAX(ISS)-SPXMIN(ISS)
              YA=SPYMAX(ISS)-SPYMIN(ISS)
C###
C             Opening input GSE file with the seismograms:
              IF(FILPAR.EQ.' ') THEN
                OPEN(LU,FILE=FILESS(ISS),STATUS='OLD')
                CALL RGSE1(LU,TEXT)
              ELSE
                ISS0=ISSRAM(IFILE(ISS)-1)+1
              END IF
C^^^
C             Loop for seismograms
              IREC=0
   40         CONTINUE
   41           CONTINUE
C###
C               Selecting the component:
   42           CONTINUE
                  IF(FILPAR.EQ.' ') THEN
                    CALL RGSE2
     *                  (LU,NAMREC,CHAN,I,X1R,X2R,X3R,T0,TD,NSS,MSS,RAM)
                  ELSE
                    ISS1=ISS0
                    CALL RRAM2(ISS0,NAMSRC,X1S,X2S,X3S,
     *                          NAMREC,X1R,X2R,X3R,I,T0,TD,NSS,IRAM,RAM)
                  END IF
                  IF(NSS.LE.-1) THEN
C                   End of the GSE file
                    GO TO 80
                  END IF
                IF(I.NE.KOMP(ISS,ISP)) GO TO 42
C^^^
C               Selecting the receiver:
                IF(FILREC.EQ.' ') THEN
                  IREC=IREC+1
                ELSE
C                 Loop for receivers
                  DO 51 I=NPTS+NSRC+1,NPTS+NSRC+NREC
                    IF(PTS(I).EQ.NAMREC) THEN
                      I0=MSS+3*I
                      X1R=RAM(I0-2)
                      X2R=RAM(I0-1)
                      X3R=RAM(I0)
                      IREC=I-NPTS-NSRC
                      GO TO 52
                    END IF
   51             CONTINUE
                  GO TO 41
                END IF
   52           CONTINUE
C###
C               Reading the source information:
                IF(FILPAR.EQ.' ') THEN
   53             CONTINUE
                    CALL RGSE2C(LINE,*54)
                    CALL RSEP2(LINE)
                  GO TO 53
   54             CONTINUE
                  CALL RSEP3T('NAMESRC',NAMSRC,' ')
                  CALL RSEP3R('X1SRC',X1S,0.)
                  CALL RSEP3R('X2SRC',X2S,0.)
                  CALL RSEP3R('X3SRC',X3S,0.)
                END IF
C^^^
C               Selecting the source:
                IF(FILSRC.NE.' ') THEN
                  IF(NAMSRC.EQ.' ') THEN
                    I0=MSS+3*NPTS
                    X1S=RAM(I0+1)
                    X2S=RAM(I0+2)
                    X3S=RAM(I0+3)
                  ELSE
C                   Loop for sources
                    DO 55 I=NPTS+1,NPTS+NSRC
                      IF(PTS(I).EQ.NAMSRC) THEN
                        I0=MSS+3*I
                        X1S=RAM(I0-2)
                        X2S=RAM(I0-1)
                        X3S=RAM(I0)
                        GO TO 56
                      END IF
   55               CONTINUE
                    GO TO 41
                  END IF
                END IF
   56           CONTINUE
C
C               Updating the coordinates:
                IF(FILPTS.NE.' ') THEN
C                 Receiver coordinates
                  DO 61 I=1,NPTS
                    IF(PTS(I).EQ.NAMREC) THEN
                      I0=MSS+3*I
                      X1R=RAM(I0-2)
                      X2R=RAM(I0-1)
                      X3R=RAM(I0)
                      GO TO 62
                    END IF
   61             CONTINUE
C                   SP-10
                    LINE='SP-10: Receiver '//NAMREC(1:LENGTH(NAMREC))
     *                                     //' not found in file PTS'
                    CALL ERROR(LINE)
C                   If file PTS with the coordinates of points is given,
C                   it has to contain all receiver names of file REC
C                   (if REC is specified) or all receiver names of the
C                   GSE file (if REC=' ').
   62             CONTINUE
C                 Source coordinates
                  IF(NAMSRC.NE.' ') THEN
                    DO 63 I=1,NPTS
                      IF(PTS(I).EQ.NAMSRC) THEN
                        I0=MSS+3*I
                        X1S=RAM(I0-2)
                        X2S=RAM(I0-1)
                        X3S=RAM(I0)
                        GO TO 64
                      END IF
   63               CONTINUE
C                   SP-11
                    CALL ERROR('SP-11: Source not found in file PTS')
C                   If file PTS with the coordinates of points is given
C                   and the GSE file contains source names, file PTS has
C                   to contain all source names of file SRC (if SRC is
C                   specified) or all source names of the GSE file (if
C                   SRC=' ').
   64               CONTINUE
                  END IF
                END IF
C
C               Plotting the seismogram:
                RR=SQRT((X1R-X1S)**2+(X2R-X2S)**2+(X3R-X3S)**2)
                IF(SPEXP(ISS).NE.0.) THEN
                  DO 76 I=1,NSS
                   RAM(ISS1+I)=RAM(ISS1+I)*EXP(SPEXP(ISS)
     *                                    *(T0+TD*FLOAT(I-1)-SPEXPT))
   76             CONTINUE
                END IF
                IF(SPVRED.NE.0.) THEN
                  T0=T0-RR/SPVRED
                END IF
C
                IF(NORMSP.LE.0) THEN
                  IF(NORMSP.LT.0) THEN
                    I1=MAX0(INT((SPTMIN-T0)/TD+2.),1)
                    I2=MIN0(INT((SPTMAX-T0)/TD+1.),NSS)
                  ELSE
                    I1=1
                    I2=NSS
                  END IF
                  AA=0.
                  DO 77 I=I1,I2
                    AA=AMAX1(ABS(RAM(ISS1+I)),AA)
   77             CONTINUE
                  IF(AA.NE.0.) THEN
                    AA=SPAMP(ISS,ISP)*SPXLEN/(RECNUM+1.)/AA
                  END IF
                ELSE
                  IF(SPOWER(ISS).EQ.0.) THEN
                    AA=SPAMP(ISS,ISP)
                  ELSE
                    AA=SPAMP(ISS,ISP)*((RR/SPDIST)**SPOWER(ISS))
                  END IF
                END IF
C
                IF(KODESP.LE.0) THEN
                  IF(FILREC.EQ.' ') THEN
C                   SP-12
                    CALL ERROR('SP-12: No receiver file specified')
C                   For KODESP=0, filename REC must be specified in the
C                   input data.
                  END IF
                  X=SPXLEN*(FLOAT(IREC)-SPXMIN(ISS))
     *                    /(SPXMAX(ISS)-SPXMIN(ISS))
                ELSE IF(KODESP.EQ.1.OR.KODESP.EQ.2) THEN
                  X=SPXLEN*((X1R-SPXMIN(ISS))*XA+(X2R-SPYMIN(ISS))*YA)
     *                    /(XA*XA+YA*YA)
                ELSE IF(KODESP.EQ.3) THEN
                  X=SPXLEN*(X3R-SPXMIN(ISS))/XA
                ELSE
                  X=SPXLEN*(RR-SPXMIN(ISS))/XA
                END IF
                Y=0.
C
C               Denoting the horizontal axis by the receiver name:
C               IF(ISS.EQ.0) THEN
                  IF(KODESP.LE.0.AND.SPXDIV.GT.0.) THEN
                    CALL SYMBOL(X-0.5*(LENGTH(NAMREC)-0.43)*SPCHRH,
     *                     Y-1.8*SPCHRH,SPCHRH,NAMREC,0.,LENGTH(NAMREC))
                  END IF
C               END IF
C
                CALL PLOT(X,Y,3)
                DO 78 I=1,NSS
                  T = T0+TD*FLOAT(I-1)
                  IF (SPTMIN.LE.T.AND.T.LE.SPTMAX) THEN
                    A = X-AA*RAM(ISS1+I)
                    Y = SCY*(T-SPTMIN)
                    IF (T.LT.SPTMIN+TD) THEN
C                     Straight line before the seismogram
                      CALL PLOT(X,Y,2)
                    END IF
                    CALL PLOT(A,Y,2)
                  END IF
   78           CONTINUE
C
C               Straight line after the seismogram
                CALL PLOT(X,Y,2)
                Y = SPTLEN
                CALL PLOT(X,Y,2)
              GO TO 40
   80         CONTINUE
C             End of loop for receivers
C###
C             Closing GSE file
              IF(FILPAR.EQ.' ') THEN
                CLOSE(LU)
              END IF
            END IF
   90     CONTINUE
C         End of loop for GSE files
C^^^
          CALL PLOT(0.,0.,999)
        END IF
   99 CONTINUE
C
C     End of loop over SP executions in optional SP parameter file:
      IF(FILPAR.NE.' ') THEN
        GO TO 100
      END IF
  999 CONTINUE
      IF(FILPAR.NE.' ') THEN
        CLOSE(LUPAR)
      END IF
C
      WRITE(*,'(A)') '+SP: Done.                                       '
      STOP
      END
C
C=======================================================================
C
C     
C
      SUBROUTINE FRAME(XX,YY,XM,XN,YM,YN,MM,IX,IY,HEIGHT,
     *                 XL1,XR1,KX1,XL2,XR2,KX2,YL,YR,KY,
     *                 K1,M1,A1,N1,K2,M2,A2,N2,K3,K4)
C
      CHARACTER*(*) KX1,KX2,KY,K1,K2,K3,K4
      REAL XX,YY,XM,XN,YM,YN,HEIGHT,XL1,XR1,XL2,XR2,YL,YR,A1,A2
      INTEGER MM,IX,IY,M1,N1,M2,N2
C
C Input:
C     XX...   Length of the horizontal axis in cm.
C     YY...   Length of the vertical axis in cm.
C     XM...   The number of intervals along the horizontal axis,
C             starting at the left.  The intervals are marked with
C             long ticks and are supplemented with the numerical
C             values if IX is positive.  XM must be positive.
C     XN...   Number of subintervals to be marked in each interval
C             with short ticks.  XN must be positive.
C     YM...   The number of intervals along the vertical axis,
C             starting at the bottom.  The intervals are marked with
C             long ticks and are supplemented with the numerical
C             values if IY is positive.  YM must be positive.
C     YN...   Number of subintervals to be marked in each interval
C             with short ticks.  YN must be positive.
C     MM...   If negative, the top line of the frame is not plotted.
C     IX...   IX=0: No labeling of the horizontal axis.
C             IX=1: Labeling of the horizontal axis with the first
C                   variable only.
C             IX=2: Labeling of the horizontal axis with both variables.
C     IY...   IY=0: No labeling of the vertical axis.
C             IY=1: Labeling of the vertical axis.
C     HEIGHT... Height of the characters in cm.
C     XL1,XR1... Values of the first variable along the horizontal axis,
C             corresponding to left-hand and right-hand corners.
C     KX1...  First label of the horizontal axis.  String to be written
C             below the horizontal axis.
C     XL2,XR2... Values of the second variable along the horizontal
C             axis, corresponding to left-hand and right-hand corners.
C     KX2...  Second label of the horizontal axis.  String to be written
C             below the horizontal axis.
C     YL,YR...Values of the variable along the vertical axis,
C             corresponding to left-hand and right-hand corners.
C     KY...   Label of the vertical axis.  String to be written to the
C             left of the vertical axis.
C     K1,A1...String and the number to be written above the frame,
C             starting 0.5*HEIGHT above the left top corner.
C     M1...   Number A1 is written only if M2 is positive.
C     N1...   Number of decimal places of A1.
C     K2,A2...String and the number to be written above the frame,
C             ending 0.5*HEIGHT above the right top corner.
C     M2...   Width of A2 in characters.  A2 is written only if M2 is
C             positive.
C     N2...   Number of decimal places of A2.
C     K3...   String to be written below the frame, starting 5.3*HEIGHT
C             below the left bottom corner.
C     K4...   String to be written below the frame, starting 7.0*HEIGHT
C             below the left bottom corner.
C
C-----------------------------------------------------------------------
C
      EXTERNAL LENGTH
      INTEGER  LENGTH
C
      INTEGER LX1,LX2,LY,L1,L2,L3,L4
C
C     Lengths of input strings
      LX1=LENGTH(KX1)
      LX2=LENGTH(KX2)
      LY =LENGTH(KY)
      L1 =LENGTH(K1)
      L2 =LENGTH(K2)
      L3 =LENGTH(K3)
      L4 =LENGTH(K4)
C
      HEIGH0=.215*HEIGHT
      HEIGH2=.500*HEIGHT
C
C     Initial values for plotting frame:
      I0= 0
      JX= IABS(IX)-1
      MX= INT(XM*XN+0.001)
      NX= INT(XN+0.5)
      MY= INT(YM*YN+0.001)
      NY= INT(YN+0.5)
      XD= XX/XM/XN
      YD= YY/YM/YN
      XH1= (XR1-XL1)/XM
      XH2= (XR2-XL2)/XM
      YH = (YR -YL )/YM
C
C     Plotting border 29.7*21.0 cm (landscape):
C     CALL PLOT( 0. , 0. ,3)
C     CALL PLOT(29.7, 0. ,2)
C     CALL PLOT(29.7,21.0,2)
C     CALL PLOT( 0. ,21.0,2)
C     CALL PLOT( 0. , 0. ,2)
C
C     Plotting frame XX*YY cm centred on the A4 sheet:
C     Landscape
*     X = (29.7-XX)/2.
*     Y = (21.0-YY)/2.
C     Portrait
      X = (21.0-XX)/2.
      Y = (29.7-YY)/2.
C     Leaving 2 cm belts for description of axes
      IF(IX.GE.0) THEN
        X=X+2.5*HEIGHT
      END IF
      IF (IY.NE.0) THEN
        Y=Y+2.5*HEIGHT
      END IF
C     Shifting the origin and plotting the frame
      CALL PLOT(X,Y,-3)
      CALL PLOT(XX,0.,2)
      CALL PLOT(XX,YY,2)
      IF(MM.GE.0) THEN
        I = 2
      ELSE
        I = 3
      END IF
      CALL PLOT(0.,YY,I)
      CALL PLOT(0.,0.,2)
C
C     Description of axes:
C
C     Bottom horizontal axis:
      X = 0.
      X1= XL1-XH1
      X2= XL2-XH2
      DO 16 I=I0,MX
        IF (MOD(I,NX).NE.0) THEN
          A = 0.1
        ELSE
          A = 0.2
          IF (JX.GE.0) THEN
            X1= X1+XH1
            X2= X2+XH2
            IF(IX.GE.0.OR.MOD(I,MX).NE.0) THEN
C             Determination of the number of decimal places:
ccc           J = INT(.99-ALOG10(AMAX1(ABS(XH1),ABS(XH2),0.001)))
              DO 11 J=0,5
                IF(AMOD(ABS(XH1)+0.000001,0.1**J).LT.0.000002) THEN
                  GO TO 12
                END IF
   11         CONTINUE
   12         CONTINUE
C             J is the preferable number of decimal places.
              JMAX=MAX0(INT(ALOG10(ABS(X1)+0.5*0.1**J))+1,1)
              IF(X1.LT.0.) JMAX=JMAX+1
C             JMAX is the number of digits left to the decimal point.
              IF (JX.GT.0) THEN
                DO 13 J=J,5
                  IF(AMOD(ABS(XH2)+0.000001,0.1**J).LT.0.000002) THEN
                    GO TO 14
                  END IF
   13           CONTINUE
   14           CONTINUE
C               J is the preferable number of decimal places.
                JMAX=MAX0(INT(ALOG10(ABS(X2)+0.5*0.1**J))+1,JMAX)
                IF(X1.GE.0..AND.X2.LT.0.) JMAX=JMAX+1
C               JMAX is the number of digits left to the decimal point.
              END IF
              JMAX=INT(XX/XM/HEIGHT)-JMAX-1
C             JMAX is the maximum number of decimal places.
              J=MIN0(J,JMAX)
              IF(J.LE.0) THEN
                J=-1
              END IF
C             J is the number of decimal places.
C
              B = X-( .5*FLOAT(1+JX)*AINT(ALOG10(ABS(X1)+.5))+0.285
     *                  -FLOAT(JX)+.5*FLOAT(J+1) )*HEIGHT
              IF(X1.LT.0.) THEN
                B=B-HEIGHT
              END IF
              CALL NUMBER(B,-1.8*HEIGHT,HEIGHT,X1,0.,J)
              IF (JX.GT.0) THEN
                B = X-HEIGHT*AINT(ALOG10(ABS(X2)+.5))+.715*HEIGHT
                IF(X2.LT.0.) B=B-HEIGHT
                CALL NUMBER(B,-3.3*HEIGHT,HEIGHT,X2,0.,J)
              END IF
            END IF
          END IF
        END IF
        IF (MOD(I,MX).NE.0) THEN
          CALL PLOT(X,0.,3)
          CALL PLOT(X,A ,2)
        END IF
        X = X+XD
   16 CONTINUE
      IF (JX.EQ.0) THEN
        B = XX-XX/XM/2.-HEIGH2*FLOAT(LX1)+HEIGH0
        CALL SYMBOL(B,-3.3*HEIGHT,HEIGHT,KX1,0.,LX1)
      ELSE IF (JX.GT.0) THEN
        CALL SYMBOL(XX+2.715*HEIGHT,-1.8*HEIGHT,HEIGHT,KX1,0.,LX1)
        CALL SYMBOL(XX+2.715*HEIGHT,-3.3*HEIGHT,HEIGHT,KX2,0.,LX2)
      END IF
C
C     Right-hand vertical axis:
      Y = 0.
      M = MY
      DO 23 I=1,M
        Y = Y+YD
        A = 0.1
        IF (MOD(I,NY).EQ.0) THEN
          A = 0.2
        END IF
        CALL PLOT(XX,Y,3)
        CALL PLOT(XX-A,Y,2)
   23 CONTINUE
C
C     Top horizontal axis:
      IF (MM.GE.0) THEN
        X = XX
        M = MX-1
        DO 33 I=1,M
          X = X-XD
          IF (MOD(I,NX).NE.0) THEN
            A = 0.1
          ELSE
            A = 0.2
          END IF
          CALL PLOT(X,YY,3)
          CALL PLOT(X,YY-A,2)
   33   CONTINUE
      END IF
C
C     Left-hand vertical axis:
      Y = 0.
      Y1= YL-YH
      DO 45 I=I0,MY
        IF (MOD(I,NY).NE.0) THEN
          A = 0.1
        ELSE
          A = 0.2
          IF (IY.NE.0) THEN
            Y1= Y1+YH
C
C           Determination of the number of decimal places:
ccc         J = INT(.99-ALOG10(AMAX1(ABS(YH),0.001)))
            DO 41 J=0,5
              IF(AMOD(ABS(YH)+0.000001,0.1**J).LT.0.000002) THEN
                GO TO 42
              END IF
   41       CONTINUE
   42       CONTINUE
            IF(J.LE.0) THEN
              J=-1
            END IF
ccc         IF(J.LE.0.OR.IY.GT.0) THEN
ccc           J=IY
ccc         END IF
C           J is the number of decimal places.
C
            B = -( AINT(ALOG10(ABS(Y1)+.5))+2.785+FLOAT(J) )*HEIGHT
            IF(Y1.LT.0.) THEN
              B=B-HEIGHT
            END IF
            CALL NUMBER(B,Y-HEIGH2,HEIGHT,Y1,0.,J)
          END IF
        END IF
        IF (I.NE.0) THEN
          CALL PLOT(0.,Y,3)
          CALL PLOT(A,Y ,2)
        END IF
        Y = Y+YD
   45 CONTINUE
      IF (IY.NE.0) THEN
        A = -HEIGHT*FLOAT(LY)-1.285*HEIGHT
        IF (YM-FLOAT(MY/NY).GE.0.25) THEN
          B = YY-HEIGHT
        ELSE
          B = YY-YY/YM/2.-HEIGH2
        END IF
        CALL SYMBOL(A,B,HEIGHT,KY,0.,LY)
      END IF
C
C     Writing texts:
      CALL SYMBOL(HEIGH0,YY+HEIGH2,HEIGHT,K1,0.,L1)
      B = HEIGHT*FLOAT(L1)+1.215*HEIGHT
      IF(M1.GT.0) THEN
        CALL NUMBER(B,YY+HEIGH2,HEIGHT,A1,0.,N1)
      END IF
      B = XX-HEIGHT*FLOAT(L2+M2)-.785*HEIGHT
      CALL SYMBOL( B,YY+HEIGH2,HEIGHT,K2,0.,L2)
      B = XX-HEIGHT*FLOAT(M2)   +.215*HEIGHT
      IF(M2.GT.0) THEN
        CALL NUMBER(B,YY+HEIGH2,HEIGHT,A2,0.,N2)
      END IF
      CALL SYMBOL(HEIGH0,-5.3*HEIGHT,HEIGHT,K3,0.,L3)
      CALL SYMBOL(HEIGH0,-7.0*HEIGHT,HEIGHT,K4,0.,L4)
C
      RETURN
      END
C
C=======================================================================
C###
      SUBROUTINE WRAM2(ISS0,NAMSRC,X1S,X2S,X3S,
     *                      NAMREC,X1R,X2R,X3R,KOMP,T0,TD,NSS,IRAM,RAM)
      CHARACTER*(*) NAMSRC,NAMREC
      INTEGER IRAM(*)
      REAL RAM(*)
C
      CHARACTER*6 NAME
      IRAM(ISS0)=22+NSS
      I=ISS0+IRAM(ISS0)
      NAME=NAMSRC
      IRAM(I-21)=ICHAR(NAME(1:1))
      IRAM(I-20)=ICHAR(NAME(2:2))
      IRAM(I-19)=ICHAR(NAME(3:3))
      IRAM(I-18)=ICHAR(NAME(4:4))
      IRAM(I-17)=ICHAR(NAME(5:5))
      IRAM(I-16)=ICHAR(NAME(6:6))
      NAME=NAMREC
      IRAM(I-15)=ICHAR(NAME(1:1))
      IRAM(I-14)=ICHAR(NAME(2:2))
      IRAM(I-13)=ICHAR(NAME(3:3))
      IRAM(I-12)=ICHAR(NAME(4:4))
      IRAM(I-11)=ICHAR(NAME(5:5))
      IRAM(I-10)=ICHAR(NAME(6:6))
      RAM(I-9)=X1SRC
      RAM(I-8)=X2SRC
      RAM(I-7)=X3SRC
      RAM(I-6)=X1REC
      RAM(I-5)=X2REC
      RAM(I-4)=X3REC
      IRAM(I-3)=KOMP
      RAM(I-2)=T0
      RAM(I-1)=TD
      ISS0=I
      IRAM(ISS0)=0
      RETURN
      END
C
C=======================================================================
C
      SUBROUTINE RRAM2(ISS0,NAMSRC,X1S,X2S,X3S,
     *                      NAMREC,X1R,X2R,X3R,KOMP,T0,TD,NSS,IRAM,RAM)
      CHARACTER*(*) NAMSRC,NAMREC
      INTEGER IRAM(*)
      REAL RAM(*)
C
      CHARACTER*6 NAME
      NSS=IRAM(ISS0)-22
      I=ISS0+IRAM(ISS0)
      IF(NSS.GT.-1) THEN
        NAME(1:1)=CHAR(IRAM(I-21))
        NAME(2:2)=CHAR(IRAM(I-20))
        NAME(3:3)=CHAR(IRAM(I-19))
        NAME(4:4)=CHAR(IRAM(I-18))
        NAME(5:5)=CHAR(IRAM(I-17))
        NAME(6:6)=CHAR(IRAM(I-16))
        NAMSRC=NAME
        NAME(1:1)=CHAR(IRAM(I-15))
        NAME(2:2)=CHAR(IRAM(I-14))
        NAME(3:3)=CHAR(IRAM(I-13))
        NAME(4:4)=CHAR(IRAM(I-12))
        NAME(5:5)=CHAR(IRAM(I-11))
        NAME(6:6)=CHAR(IRAM(I-10))
        NAMREC=NAME
        X1SRC=RAM(I-9)
        X2SRC=RAM(I-8)
        X3SRC=RAM(I-7)
        X1REC=RAM(I-6)
        X2REC=RAM(I-5)
        X3REC=RAM(I-4)
        KOMP=IRAM(I-3)
        T0=RAM(I-2)
        TD=RAM(I-1)
      END IF
      ISS0=I
      RETURN
      END
C^^^
C=======================================================================
C
      INCLUDE 'error.for'
C     error.for
      INCLUDE 'sep.for'
C     sep.for
      INCLUDE 'gse.for'
C     gse.for
      INCLUDE 'length.for'
C     length.for
      INCLUDE 'calcops.for'
C     calcops.for
C
C=======================================================================
C