C
C Program PTSGRD to generate grid file containing undefined values
C at the gridpoints closest to the given points and zeros at the other
C gridpoints
C
C Version: 5.40
C Date: 2000, January 26
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 Imput and output filenames:
C     PTS='string'... String in apostrophes containing the name of the
C             input file with the coordinates of the points.
C             If PTS=' ', zeros at all gridpoints are generated.
C             Description of the file PTS
C             Default: PTS=' '
C     PTSGRD='string'... String in apostrophes containing the name of
C             the output file with the grid values.
C             For general description of files with gridded data refer
C             to file forms.htm.
C             No default, 'PTSGRD' must be specified and cannot be blank
C Data specifying grid dimensions:
C     N1=positive integer... Number of gridpoints along the X1 axis.
C             Default: N1=1
C     N2=positive integer... Number of gridpoints along the X2 axis.
C             Default: N2=1
C     N3=positive integer... Number of gridpoints along the X3 axis.
C             Default: N3=1
C     O1=real... First coordinate of the grid origin (first point of the
C             grid).
C             Default: O1=0.
C     O2=real... Second coordinate of the grid origin.
C             Default: O2=0.
C     O3=real... Third coordinate of the grid origin.
C             Default: O3=0.
C     D1=real... Grid interval in the direction of the first coordinate
C             axis.
C             Default: D1=1.
C     D2=real... Grid interval in the direction of the second coordinate
C             axis.
C             Default: D2=1.
C     D3=real... Grid interval in the direction of the third coordinate
C             axis.
C             Default: D3=1.
C
C                                                     
C Input file PTS with given points:
C (1) None to several strings terminated by / (a slash)
C (2) For each gridpoint data (2.1):
C (2.1) 'NAME',X1,X2,X3,/
C     'NAME'... Name of the point.  Any string different from '$'.
C      X1,X2,X3... Coordinates of the point.
C (3) / or end of file
C
C-----------------------------------------------------------------------
C
C Common block /RAMC/:
      INCLUDE 'ram.inc'
C     ram.inc
C
C.......................................................................
C
C     Filenames and parameters:
      CHARACTER*80 FSEP,FPTS,FGRD
      INTEGER LU
      REAL UNDEF
      PARAMETER (LU=1,UNDEF=-999999.)
C
C     Data:
      CHARACTER*1 TEXT
      INTEGER I,I1,I2,I3,N1,N2,N3
      REAL D1,D2,D3,O1,O2,O3,X1,X2,X3
C
C-----------------------------------------------------------------------
C
C     Reading name of SEP file with input data:
      WRITE(*,'(A)') '+PTSGRD: Enter input filename: '
      FSEP=' '
      READ(*,*) FSEP
C
C     Reading all data from the SEP file into the memory:
      IF (FSEP.NE.' ') THEN
        CALL RSEP1(LU,FSEP)
      ELSE
C       PTSGRD-01
        CALL ERROR('PTSGRD-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
      WRITE(*,'(A)') '+PTSGRD: Working...            '
C
C     Recalling the data specifying grid dimensions
      CALL RSEP3I('N1',N1,1)
      CALL RSEP3I('N2',N2,1)
      CALL RSEP3I('N3',N3,1)
      CALL RSEP3R('O1',O1,0.)
      CALL RSEP3R('O2',O2,0.)
      CALL RSEP3R('O3',O3,0.)
      CALL RSEP3R('D1',D1,1.)
      CALL RSEP3R('D2',D2,1.)
      CALL RSEP3R('D3',D3,1.)
C
C     Input and output filenames:
      CALL RSEP3T('PTS'   ,FPTS,' ')
      CALL RSEP3T('PTSGRD',FGRD,' ')
      IF (FGRD.EQ.' ') THEN
C       PTSGRD-02
        CALL ERROR('PTSGRD-02: Output file PTSGRD not given')
C       Output file PTSGRD must be specified.
C       There is no default filename.
      ENDIF
C
      DO 10 I=1,N1*N2*N3
        RAM(I)=0.
   10 CONTINUE
C
      IF(FPTS.NE.' ') THEN
C
        OPEN(LU,FILE=FPTS,STATUS='OLD')
        READ(LU,*) (TEXT,I=1,20)
   20   CONTINUE
          TEXT='$'
          X1=0.
          X2=0.
          X3=0.
          READ(LU,*,END=29) TEXT,X1,X2,X3
          IF(TEXT.EQ.'$') THEN
            GO TO 29
          END IF
          I1=NINT((X1-O1)/D1)
          I2=NINT((X2-O2)/D2)
          I3=NINT((X3-O3)/D3)
          IF(0.LE.I1.AND.I1.LT.N1.AND.
     *       0.LE.I2.AND.I2.LT.N2.AND.
     *       0.LE.I3.AND.I3.LT.N3) THEN
            RAM(1+I1+N1*(I2+N2*I3))=UNDEF
          END IF
        GO TO 20
   29   CONTINUE
        CLOSE(LU)
      END IF
C
C     Writing output grid values:
      CALL WARRAY(LU,FGRD,'FORMATTED',.TRUE.,UNDEF,.FALSE.,0.,
     *                                                     N1*N2*N3,RAM)
      WRITE(*,'(A)') '+PTSGRD: Done.                 '
      STOP
      END
C
C=======================================================================
C
      INCLUDE 'error.for'
C     error.for
      INCLUDE 'sep.for'
C     sep.for
      INCLUDE 'forms.for'
C     forms.for
      INCLUDE 'length.for'
C     length.for
C
C=======================================================================
C