C
C Program GRDMERGE to merge two nonoverlapping sets of values given on
C the same grid into a single set.
C
C Version: 5.40
C Date: 2000, January 22
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 and output files:
C     GRD1='string', GRD2='string'... Names of the input ASCII files
C             with the grid values.  No gridpoint may have the value
C             defined in both the files if 'GRD4' is blank (default).
C             Default: GRD1='grd1.out', GRD2='grd2.out'
C     GRD3='string'... Name of the output ASCII file containing the
C             grid values collected from both the input files.
C             Default: GRD3='grd.out'
C     GRD4='string'... Name of the auxiliary output file.
C             If specified and nonblank, the input grid values may
C             overlap.  Of each pair of overlaping values, the value
C             of file GRD2 will be collected in this file.  This file
C             is not created if there is no overlap.
C             Default: GRD4=' ' (auxiliary file not generated)
C     For general description of the files with gridded data refer
C     to file forms.htm.
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 Optional parameters specifying the form of the real quantities
C written in the output formatted files:
C     MINDIG,MAXDIG=positive integers ... See the description in file
C             forms.for.
C
C=======================================================================
C
C Common block /RAMC/:
      INCLUDE 'ram.inc'
C     ram.inc
C
C.......................................................................
C
      INTEGER LU,N1,N2,N3,N1N2N3,I
      REAL UNDEF
      PARAMETER (LU=1,UNDEF=-999999.)
      CHARACTER*80 FSEP,FGRD1,FGRD2,FGRD,FGRDA
      LOGICAL LGRDA
C
C-----------------------------------------------------------------------
C
C     Reading name of SEP file with input data:
      WRITE(*,'(A)') '+GRDMERGE: Enter input filename: '
      FSEP=' '
      READ(*,*) FSEP
      WRITE(*,'(A)') '+GRDMERGE: Working ...           '
C
C     Reading all data from the SEP file into the memory:
      IF (FSEP.NE.' ') THEN
        CALL RSEP1(LU,FSEP)
      ELSE
C       GRDMERGE-03
        CALL ERROR('GRDMERGE-03: 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 input parameters from the SEP file:
      CALL RSEP3T('GRD1',FGRD1,'grd1.out')
      CALL RSEP3T('GRD2',FGRD2,'grd2.out')
      CALL RSEP3T('GRD3',FGRD,'grd.out')
      CALL RSEP3T('GRD4',FGRDA,' ')
C
C     Recalling the data specifying grid dimensions
C     (arguments: Name of value in input data, Variable, Default):
      CALL RSEP3I('N1',N1,1)
      CALL RSEP3I('N2',N2,1)
      CALL RSEP3I('N3',N3,1)
      N1N2N3=N1*N2*N3
      IF(2*N1N2N3.GT.MRAM) THEN
C       GRDMERGE-01
        CALL ERROR('GRDMERGE-01: Too small array RAM(MRAM)')
C       Array RAM(MRAM) allocated in include file 'ram.inc' is too small
C       to contain two input grids (2*N1*N2*N3 values).  You may wish to
C       increse the dimension MRAM in file 'ram.inc'.
C       ram.inc
      END IF
C
C     Reading input grids:
      CALL RARRAY(LU,FGRD1,'FORMATTED',.TRUE.,UNDEF,N1N2N3,RAM)
      CALL RARRAY(LU,FGRD2,'FORMATTED',.TRUE.,UNDEF,N1N2N3,
     *                                                    RAM(N1N2N3+1))
C
C     Merging the grids:
      LGRDA=.FALSE.
      DO 10 I=1,N1N2N3
        IF(RAM(I).EQ.UNDEF) THEN
          RAM(I)=RAM(N1N2N3+I)
          RAM(N1N2N3+I)=UNDEF
        ELSE
          IF(RAM(N1N2N3+I).NE.UNDEF) THEN
            IF(FGRDA.EQ.' ') THEN
C             GRDMERGE-02
              CALL ERROR('GRDMERGE-02: Overlapping grid values')
C             Value at the same gridpoint is defined in both input files
C             and the auxiliary file is not specified.
            ELSE
              LGRDA=.TRUE.
            END IF
          END IF
        END IF
   10 CONTINUE
C
C     Writing output grid:
      CALL WARRAY(LU,FGRD,'FORMATTED',.TRUE.,UNDEF,.FALSE.,0.,N1N2N3,
     *                                                              RAM)
      IF(LGRDA) THEN
        CALL WARRAY(LU,FGRDA,'FORMATTED',.TRUE.,UNDEF,.FALSE.,0.,N1N2N3,
     *                                                    RAM(N1N2N3+1))
      END IF
      WRITE(*,'(A)') '+GRDMERGE: 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