C 
C Program GREEN converting the unformatted output of program CRT into a
C formatted file containing the ray-theory elastodynamic Green function.
C
C Version: 5.10
C Date: 1997, September 27
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                                                    
C Description of data files:
C
C Main input data read from external interactive device (*):
C     The data consist of character strings read by list directed (free
C     format) input.  The strings have thus to be enclosed in
C     apostrophes.  The interactive * external unit may be redirected to
C     the file containing the data.
C (1) 'REC','SRC','CRT-S','CRT-I','GREEN'
C     'REC'...If non-blank, the name of the file with the names of the
C             receiver points.  The names are then used within the
C             strings describing the points of two-point rays.
C             Otherwise, the two-point rays are denoted by the receiver
C             index.
C             Description of file REC
C     'SRC'...If non-blank, the name of the file with the name of the
C             source point.  The name is then used within the strings
C             describing the rays.
C             Description of file SRC
C     'CRT-S'... File with the quantities stored at the points of
C             intersections of rays with the specified reference surface
C             (see C.R.T.5.5.2).  Only two-point rays incident to the
C             vicinities of the receivers situated at the reference
C             surface are considered.  The output elementary Green
C             functions are then evaluated at the receivers.
C             Description of file CRT-S
C     'CRT-I'... File with the quantities at the initial points of rays,
C             corresponding to file 'CRT-S' (see C.R.T.6.1).
C             Description of file CRT-I
C     'GREEN'... Name of the output formatted file with the Green
C             tensor.
C             Description of file GREEN
C Default: 'REC'=' ', 'SRC'=' ', 'CRT-S'='s01.out', 'CRT-I'='s01i.out',
C     'GREEN'='green.out'.
C
C Input formatted files 'REC' and 'SRC', if specified, must correspond
C     to the receiver and source files used during Complete Ray Tracing.
C
C Input unformatted file 'CRT-S':
C     See the description within source code file 'writ.for'.
C     Description of file CRT-S
C
C Input unformatted file 'CRT-I':
C     See the description within source code file 'writ.for'.
C     Description of file CRT-I
C
C                                                   
C Output formatted file 'GREEN':
C (1) / (a slash).
C (2) For each two-point ray (2.1):
C (2.1) 'R','S',(GREEN(I),I=1,32),/
C     'R'...  String in apostrophes describing the receiver.
C     'S'...  String in apostrophes describing the source.
C     GREEN(1)... Travel time between receiver and source.
C     GREEN(2)... Imaginary part of the complex-valued travel time
C             between receiver and source due to attenuation.
C     GREEN(3:8)... Coordinates of the receiver and coordinates of the
C             source.
C     GREEN(9:14)... Derivatives of the travel time with respect to the
C             coordinates of the receiver and coordinates of the source.
C     GREEN(15:32)... 1000000 times enlarged amplitude of the Green
C             function: contravariant components of the complex-valued
C             3*3 matrix Gij in model coordinates, where the first
C             subscript corresponds to the receiver and the second
C             subscript corresponds to the source.  The components are
C             ordered as
C             Re(G11),Im(G11),Re(G21),Im(G21),Re(G31),Im(G31),
C             Re(G12),Im(G12),Re(G22),Im(G22),Re(G32),Im(G32),
C             Re(G13),Im(G13),Re(G23),Im(G23),Re(G33),Im(G33).
C     /...    An obligatory slash after at the end of line, in place
C             where the slowness vector components could be written.
C (3) / (a slash).
C File form GREEN is, to some extent, an extension of file form
C Travel Times.
C
C-----------------------------------------------------------------------
C
C Common block /POINTC/ to store the results of complete ray tracing:
      INCLUDE 'pointc.inc'
C     pointc.inc
C None of the storage locations of the common block are altered.
C
C Subroutines and external functions required:
      EXTERNAL AP00,AP21,TXT1,TXT2,REC,SRC,TTSORT,FORM2
C     AP00,AP21... File 'ap.for'.
C     AP03,AP03A... File 'ap.for' (called by AP21,AP03).
C     TXT1,TXT2,REC,SRC... File 'crtout.for'.
C     LENGTH..File 'length.for' (called by TXT2).
C     TTSORT..File 'ttsort.for'.
C     INDEXX..File 'indexx.for' (called by TTSORT).
C     FORM2...File 'forms.for'.
C     FORM1...File 'forms.for' (called by FORM2).
C
C-----------------------------------------------------------------------
C
C Common block /RAMC/:
      INCLUDE 'ram.inc'
C     ram.inc
C
C     Allocation of working arrays:
      INTEGER MPTS
      PARAMETER (MPTS=MRAM/34)
      INTEGER IRECS(MPTS),INDX(MPTS)
      REAL RECS(MPTS),GREEN(32,MPTS)
      EQUIVALENCE (IRECS,RAM          )
      EQUIVALENCE (RECS ,RAM          )
      EQUIVALENCE (INDX ,RAM(  MPTS+1))
      EQUIVALENCE (GREEN,RAM(2*MPTS+1))
C
C-----------------------------------------------------------------------
C
C     Auxiliary storage locations:
      INTEGER LU1,LU2,LU3,NPTS
      PARAMETER (LU1=1,LU2=2,LU3=3)
C     MPTS... maximum number of two-point rays.
C     Note: MPTS=1 indicates no sorting according to receivers.
      CHARACTER*80 FILREC,FILSRC,FILE1,FILE2,FILE3
      CHARACTER*260 FORMAT
      CHARACTER*80 RAYTXT
      INTEGER LENTXT,I,J,K,L
      REAL COOR(3)
C
C.......................................................................
C
C     Opening input and output files:
      FILREC=' '
      FILSRC=' '
      FILE1='s01.out'
      FILE2='s01i.out'
      FILE3='green.out'
      WRITE(*,'(A)') ' Enter 5 filenames (REC,SRC,S01,S01I,GREEN): '
      READ(*,*) FILREC,FILSRC,FILE1,FILE2,FILE3
C
      CALL TXT1(LU1,FILSRC,FILREC)
      OPEN(LU1,FILE=FILE1,FORM='UNFORMATTED',STATUS='OLD')
      OPEN(LU2,FILE=FILE2,FORM='UNFORMATTED',STATUS='OLD')
      OPEN(LU3,FILE=FILE3)
      WRITE(LU3,'(A)') '/'
C
C.......................................................................
C
C     Loop for the points of rays:
      NPTS=0
   10 CONTINUE
C       Reading the results of the complete ray tracing
        CALL AP00(0,LU1,LU2)
        IF(IWAVE.LT.1)THEN
C         End of rays
          GO TO 60
        ELSE IF (IREC.GT.0) THEN
          NPTS=NPTS+1
          CALL AP21(GREEN(1,NPTS))
          IRECS(NPTS)=IREC
C
C         Linear Taylor expansion of travel time:
C         Receiver:
          IF(FILREC.NE.' ') THEN
C           Receiver:
            COOR(1)=GREEN(3,NPTS)
            COOR(2)=GREEN(4,NPTS)
            COOR(3)=GREEN(5,NPTS)
            CALL REC(IREC,COOR(1),COOR(2),COOR(3))
            GREEN(1,NPTS)=GREEN(1,NPTS)
     *                   +(COOR(1)-GREEN(3,NPTS))*GREEN( 9,NPTS)
     *                   +(COOR(2)-GREEN(4,NPTS))*GREEN(10,NPTS)
     *                   +(COOR(3)-GREEN(5,NPTS))*GREEN(11,NPTS)
            GREEN(3,NPTS)=COOR(1)
            GREEN(4,NPTS)=COOR(2)
            GREEN(5,NPTS)=COOR(3)
          END IF
          IF(FILSRC.NE.' ') THEN
C           Source:
            COOR(1)=GREEN(6,NPTS)
            COOR(2)=GREEN(7,NPTS)
            COOR(3)=GREEN(8,NPTS)
            CALL SRC(1,COOR(1),COOR(2),COOR(3))
            GREEN(1,NPTS)=GREEN(1,NPTS)
     *                   +(COOR(1)-GREEN(6,NPTS))*GREEN(12,NPTS)
     *                   +(COOR(2)-GREEN(7,NPTS))*GREEN(13,NPTS)
     *                   +(COOR(3)-GREEN(8,NPTS))*GREEN(14,NPTS)
            GREEN(6,NPTS)=COOR(1)
            GREEN(7,NPTS)=COOR(2)
            GREEN(8,NPTS)=COOR(3)
          END IF
C
C         Storing or writing:
          DO 20 I=15,32
            GREEN(I,NPTS)=1000000.*GREEN(I,NPTS)
   20     CONTINUE
          IF(MPTS.EQ.1) THEN
C           Text strings:
            CALL TXT2(0,1,IWAVE,IRAY,IREC,LENTXT,RAYTXT)
            L=INDEX(RAYTXT(1:LENTXT),'''')
            L=INDEX(RAYTXT(L+1:LENTXT),'''')+L
            IF(FILREC.EQ.' ') THEN
C             Shortening the receiver string part from 8 to 6 characters
              LENTXT=LENTXT-2
              DO 21 I=L+4,LENTXT
                RAYTXT(I:I)=RAYTXT(I+2:I+2)
   21         CONTINUE
            END IF
C
C           Writing:
            FORMAT(1:4)='(4A,'
            CALL FORM2(32,GREEN,GREEN,FORMAT(5:260))
            WRITE(LU3,FORMAT) RAYTXT(L+2:LENTXT),
     *                      ' ',RAYTXT(1:L),(' ',GREEN(I,1),I=1,32),' /'
            NPTS=0
          END IF
        END IF
      GO TO 10
   60 CONTINUE
C
C.......................................................................
C
C     Sorting and writing two-point rays:
      IF(MPTS.GT.1) THEN
        CALL TTSORT(32,NPTS,1,GREEN,IRECS,RECS,INDX)
        DO 80 K=1,NPTS
          J=INDX(K)
C
C         Text strings:
          CALL TXT2(0,1,1,1,IRECS(J),LENTXT,RAYTXT)
          L=INDEX(RAYTXT(1:LENTXT),'''')
          L=INDEX(RAYTXT(L+1:LENTXT),'''')+L
          IF(FILREC.EQ.' ') THEN
C           Shortening the receiver string part from 8 to 6 characters
            LENTXT=LENTXT-2
            DO 71 I=L+4,LENTXT
              RAYTXT(I:I)=RAYTXT(I+2:I+2)
   71       CONTINUE
          END IF
C
C         Writing:
          FORMAT(1:4)='(4A,'
          CALL FORM2(32,GREEN(1,J),GREEN(1,J),FORMAT(5:260))
          WRITE(LU3,FORMAT) RAYTXT(L+2:LENTXT),
     *                      ' ',RAYTXT(1:L),(' ',GREEN(I,J),I=1,32),' /'
   80   CONTINUE
      END IF
C
C.......................................................................
C
      WRITE(LU3,'(A)') '/'
      CLOSE(LU3)
      CLOSE(LU2)
      CLOSE(LU1)
      STOP
      END
C
C=======================================================================
C
      INCLUDE 'ap.for'
C     ap.for
      INCLUDE 'crtout.for'
C     crtout.for
      INCLUDE 'ttsort.for'
C     ttsort.for
      INCLUDE 'forms.for'
C     forms.for
      INCLUDE 'length.for'
C     length.for
      INCLUDE 'indexx.for'
C     indexx.for
C
C=======================================================================
C