C
C Program TRGLSORT to sort triangles according to values at its vertices
C
C Version: 5.80
C Date: 2004, June 11
C
C Coded by: Petr Bulant
C     Department of Geophysics, Charles University Prague,
C     Ke Karlovu 3, 121 16 Praha 2, Czech Republic,
C     E-mail: bulant@seis.karlov.mff.cuni.cz
C
C.......................................................................
C This program reads file TRGL with triangles and the corresponding
C file VRTX with coordinates of the vertices of the triangles.
C The triangles are sorted into three output files TRGLN, TRGLN1 and
C TRGLN2 according one of the following criteria:
C
C If KOLUMN and VALUE are given:
C  The program compares the values in KOLUMN-th column of file VRTX
C with the value VALUE.  The triangles with the values in all their
C vertices equal to VALUE are written to files TRGLN and VRTXN.
C The triangles with the values in all their vertices less than VALUE
C are written to files TRGLN1 and VRTXN1.
C Remaining triangles are written to files TRGLN2 and VRTXN2.
C  This option may be used, e.g., to separate triangles corresponding
C to given interface, or to separate triangles according to the
C coordinates of their vertices.
C
C If KOLUM1 and KOLUM2 are given:
C  The program computes for each triangle the absolute value of
C the sum of values in column KOLUM1 and absolute value of the sum
C of values in column KOLUM2 of file VRTX.  Triangles with
C sum corresponding to column KOLUM1 equal the sum corresponding
C to column KOLUM2 are written to file TRGLN and their vertices to
C file VRTXN.  Triangles with sum corresponding to column KOLUM1 less
C then sum corresponding to column KOLUM2 are written to file TRGLN1
C and their vertices to file VRTXN1.  Remaining triangles and their
C vertices are written to files TRGLN2 and VRTXN2.
C  This option may be used, e.g., to sort the triangles according
C to their distances from two interfaces.  The distances must be
C calculated in advance by program 'intf.for', and written to columns
C KOLUM1 and KOLUM2 of the corresponding files.
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 Data specifying input files:
C     VRTX='string'... Name of the file with vertices of the triangles.
C             Description of file VRTX
C             Default: VRTX='vrtx.out'
C     TRGL='string'... Name of the file with the triangles.
C             Description of file TRGL
C             Default: TRGL='trgl.out'
C Data for sorting of the triangles:
C     KOLUMN=integer, KOLUM1=integer, KOLUM2=integer... Indices
C             of columns in input file VRTX, which contain the values,
C             according which the triangles are to be sorted.
C             Default: KOLUMN=0, KOLUM1=0, KOLUM2=0
C     VALUE=real ... Value according which the triangles are
C             to be sorted.
C             Default: VALUE=UNDEF
C                      For the value of UNDEF see below.
C     Either KOLUMN and VALUE, or KOLUM1 and KOLUM2 must be specified,
C     (i.e. one pair of parameters must be specified, specification of
C     both pairs of parameters results in error).
C Data specifying output files:
C     VRTXN='string',VRTXN1='string',VRTXN2='string' ...
C             Names of the output files with vertices of triangles.
C             If blank, files are not generated.
C             Description of files VRTXNi
C             Default: VRTXN=' ', VRTXN1=' ', VRTXN2=' '
C     TRGLN='string',TRGLN1='string',TRGLN2='string' ...
C             Names of the output files with the triangles.
C             If blank, files are not generated.
C             Description of files TRGLNi
C             Default: TRGLN=' ', TRGLN1=' ', TRGLN2=' '
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 Input file VRTX with the vertices:
C (1) None to several strings terminated by / (a slash)
C (2) For each vertex data (2.1):
C (2.1) 'NAME',X1,X2,X3,R1,R2,/
C     'NAME'... Name of the vertex.  Not considered.  May be blank.
C     X1,X2,X3... Coordinates of the vertex.
C     R1,R2,/ ... None to several values terminated by a slash.
C             Number of values should be the same for all the vertices.
C (3) / or end of file.
C
C                                                    
C Input file TRGL with the triangles:
C (1) For each triangle data (1.1):
C (1.1) I1,I2,I3,/
C     I1,I2,I3... Indices of 3 vertices of the triangle.
C             The vertices in file VRTX are indexed by positive integers
C             according to their order.
C     /...    List of vertices of the triangle is terminated by a slash.
C
C                                                   
C Output files VRTXN, VRTXN1 and VRTXN2 with the vertices:
C (1) / (a slash)
C (2) For each vertex data (2.1):
C (2.1) 'NAME',X1,X2,X3,R1,R2,/
C     'NAME'..Name of the vertex.  String in apostrophes containing
C             the index of the vertex corresponding to file TRGLN,
C             TRGLN1 or TRGLN2.
C     X1,X2,X3,R1,R2,/ ... Unchanged values from the file VRTX.
C (3) / (a slash)
C
C                                                   
C Output files TRGLN, TRGLN1 and TRGLN2 with the triangles:
C (1) For each triangle data (1.1):
C (1.1) I1,I2,I3,/
C     I1,I2,I3... Indices of 3 vertices of the triangle.
C             The vertices in corresponding file VRTXN, VRTXN1 or VRTXN2
C             are indexed by positive integers according to their order.
C     /...    List of vertices of the triangle is terminated by a slash.
C
C=======================================================================
C
C Common block /RAMC/:
      INCLUDE 'ram.inc'
C     ram.inc
C
      INTEGER IRAM(MRAM)
      EQUIVALENCE (IRAM,RAM)
C
C.......................................................................
C
      EXTERNAL LENGTH,ERROR,FORM1,WARN,RSEP1,RSEP3T,RSEP3R,RSEP3I
      INTEGER  LENGTH
C
C.......................................................................
C
C     Auxiliary storage locations:
      INTEGER LU,IUNDEF,NOUT
      PARAMETER (LU=1,IUNDEF=-999999,NOUT=3)
      REAL UNDEF
      PARAMETER (UNDEF=9.9E9)
      CHARACTER*80 FSEP,FVRTXO,FTRGLO,FVRTXN(NOUT),FTRGLN(NOUT)
      CHARACTER*10 FORMA1
      CHARACTER*26 FORMA2
      CHARACTER*80 TEXT
      INTEGER NVRTX,NTRGL,I,I1,I2,I3,J1,J2,J3,NQ
      INTEGER KOLUMN,KOLUM1,KOLUM2
      REAL A1,A2,A3,B1,B2,B3,W1,W2,OUTMIN,OUTMAX,VALUE
C
C     NVRTX...Last storage location with the vertices,
C             i.e. NQ+NOUT times the number of vertices.
C     NTRGL...Last storage location with triangles,
C             i.e. NVRTX + 3+NOUT times the number of triangles.
C
C-----------------------------------------------------------------------
C
C     Reading name of SEP file with input data:
      WRITE(*,'(A)') '+TRGLSORT: Enter input filename: '
      FSEP=' '
      READ (*,*) FSEP
      WRITE(*,'(A)') '+TRGLSORT: Working...            '
C
C     Reading all data from the SEP file into the memory:
      IF (FSEP.NE.' ') THEN
        CALL RSEP1(LU,FSEP)
      ELSE
C       TRGLSORT-01
        CALL ERROR('TRGLSORT-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 input and output filenames:
      CALL RSEP3T('VRTX',FVRTXO,'vrtx.out')
      CALL RSEP3T('TRGL',FTRGLO,'trgl.out')
      CALL RSEP3T('VRTXN', FVRTXN(1),' ')
      CALL RSEP3T('VRTXN1',FVRTXN(2),' ')
      CALL RSEP3T('VRTXN2',FVRTXN(3),' ')
      CALL RSEP3T('TRGLN', FTRGLN(1),' ')
      CALL RSEP3T('TRGLN1',FTRGLN(2),' ')
      CALL RSEP3T('TRGLN2',FTRGLN(3),' ')
C
C     Reading the columns with the values and the limits of coordinates:
      CALL RSEP3I('KOLUMN',KOLUMN,0)
      CALL RSEP3R('VALUE',VALUE,UNDEF)
      CALL RSEP3I('KOLUM1',KOLUM1,0)
      CALL RSEP3I('KOLUM2',KOLUM2,0)
      IF ((.NOT.(((KOLUMN.NE.0).AND.(VALUE.NE.UNDEF)).OR.
     *           ((KOLUM1.NE.0).AND.(KOLUM2.NE.0   ))))     .OR.
     *    (     (((KOLUMN.NE.0).AND.(VALUE.NE.UNDEF)).AND.
     *           ((KOLUM1.NE.0).AND.(KOLUM2.NE.0   ))))) THEN
C       TRGLSORT-02
        CALL ERROR('TRGLSORT-02: Wrong sorting criterion.')
C       Either KOLUMN and VALUE, or KOLUM1 and KOLUM2 must be specified,
C       see input data.
      ENDIF
C
C     Reading vertices:
      OPEN(LU,FILE=FVRTXO,STATUS='OLD')
      READ(LU,*) (TEXT,I=1,20)
      DO 2, I=1,MRAM
        RAM(I)=UNDEF
   2  CONTINUE
      TEXT='$'
      READ(LU,*,END=18) TEXT,(RAM(I),I=1,MRAM)
      IF(TEXT.EQ.'$') GO TO 18
      NQ=0
      DO 4, I=MRAM,1,-1
        IF (RAM(I).NE.UNDEF) THEN
          NQ=I
          GOTO 5
        ENDIF
   4  CONTINUE
   5  CONTINUE
      IF (NQ.LT.MAX0(KOLUM1,KOLUM2)) THEN
C       TRGLSORT-03
        CALL ERROR('TRGLSORT-03: Missing values in file VRTX')
C       Each line of file VRTX must contain at least MAX(KOLUM1,KOLUM2)
C       quantities.
      ENDIF
      CLOSE(LU)
      OPEN(LU,FILE=FVRTXO,STATUS='OLD')
      READ(LU,*) (TEXT,I=1,20)
      NVRTX=0
  10  CONTINUE
        IF(NVRTX+NQ+NOUT.GT.MRAM) THEN
C         TRGLSORT-04
          CALL ERROR('TRGLSORT-04: Too small array RAM')
        END IF
        TEXT='$'
        READ(LU,*,END=18) TEXT,(RAM(I),I=NVRTX+1,NVRTX+NQ)
        IF(TEXT.EQ.'$') THEN
          GO TO 18
        END IF
        DO 16, I=1,NOUT
          IRAM(NVRTX+NQ+I)=IUNDEF
  16    CONTINUE
        NVRTX=NVRTX+NQ+NOUT
      GO TO 10
  18  CONTINUE
      CLOSE(LU)
C
C     Reading triangles:
      DO 19 I=NVRTX+1,MRAM
        IRAM(I)=IUNDEF
  19  CONTINUE
      OPEN(LU,FILE=FTRGLO,STATUS='OLD')
      NTRGL=NVRTX
  20  CONTINUE
        IF (NTRGL+3+NOUT.GT.MRAM) THEN
C         TRGLSORT-05
          CALL ERROR('TRGLSORT-05: Too small array RAM')
        ENDIF
        READ(LU,*,END=22) (IRAM(I),I=NTRGL+1,NTRGL+3)
        NTRGL=NTRGL+3+NOUT
      GOTO 20
  22  CONTINUE
      CLOSE(LU)
C
C     Sorting the triangles:
      DO 30, I1=NVRTX,NTRGL-3-NOUT,3+NOUT
        J1=(NQ+NOUT)*(IRAM(I1+1)-1)
        J2=(NQ+NOUT)*(IRAM(I1+2)-1)
        J3=(NQ+NOUT)*(IRAM(I1+3)-1)
        W1=0.
        W2=0.
        IF ((KOLUM1.NE.0).AND.(KOLUM2.NE.0)) THEN
          A1=RAM(J1+KOLUM1)
          A2=RAM(J2+KOLUM1)
          A3=RAM(J3+KOLUM1)
          B1=RAM(J1+KOLUM2)
          B2=RAM(J2+KOLUM2)
          B3=RAM(J3+KOLUM2)
          W1=ABS(A1+A2+A3)
          W2=ABS(B1+B2+B3)
        ELSE
          A1=RAM(J1+KOLUMN)
          A2=RAM(J2+KOLUMN)
          A3=RAM(J3+KOLUMN)
          IF ((A1.EQ.VALUE).AND.(A2.EQ.VALUE).AND.(A3.EQ.VALUE)) THEN
            CONTINUE
          ELSEIF((A1.LT.VALUE).AND.(A2.LT.VALUE).AND.(A3.LT.VALUE)) THEN
            W2=1.
          ELSE
            W1=1.
          ENDIF
        ENDIF
        IF     (W1.EQ.W2) THEN
          IRAM(I1+3+1)=1
          IRAM(J1+NQ+1)=1
          IRAM(J2+NQ+1)=1
          IRAM(J3+NQ+1)=1
        ELSEIF (W1.LT.W2) THEN
          IRAM(I1+3+2)=1
          IRAM(J1+NQ+2)=1
          IRAM(J2+NQ+2)=1
          IRAM(J3+NQ+2)=1
        ELSE
          IRAM(I1+3+3)=1
          IRAM(J1+NQ+3)=1
          IRAM(J2+NQ+3)=1
          IRAM(J3+NQ+3)=1
        ENDIF
  30  CONTINUE
C
C     Output format for the output files:
      FORMA1='(99(I0,A))'
      FORMA2='(A,1I0.0,A,00(F00.0,1X),A)'
      I=INT(ALOG10(FLOAT(NVRTX/(NQ+NOUT))))+1
      IF (I.GT.9) THEN
C       TRGLSORT-06
        CALL ERROR('TRGLSORT-06: Too many vertices in file VRTX')
C       This format specification allows for maximum of 100 000 000
C       of vertices in file VRTX
      ENDIF
      FORMA1(6:6)=CHAR(ICHAR('0')+I)
      FORMA2(6:6)=FORMA1(6:6)
      FORMA2(8:8)=FORMA1(6:6)
      FORMA2(13:13)=CHAR(ICHAR('0')+MOD(NQ/1,10))
      FORMA2(12:12)=CHAR(ICHAR('0')+MOD(NQ/10,10))
C
C     Indices of the vertices in the output files:
      DO 40, I1=1,NOUT
        I3=0
        DO 39, I2=NQ,NVRTX-NOUT,NQ+NOUT
          IF (IRAM(I2+I1).EQ.1) THEN
            I3=I3+1
            IRAM(I2+I1)=I3
          ENDIF
  39    CONTINUE
  40  CONTINUE
C
C     Writing the output files:
      DO 60, I1=1,3
C       Writing the vertices:
        IF (FVRTXN(I1).NE.' ') THEN
          OPEN(LU,FILE=FVRTXN(I1))
          WRITE(LU,'(A)') '/'
          DO 44, I2=0,NVRTX-NQ-NOUT,NQ+NOUT
            I3=IRAM(I2+NQ+I1)
            IF (I3.NE.IUNDEF) THEN
              OUTMIN=0.
              OUTMAX=0.
              DO 42, I=I2+1,I2+NQ
                IF(RAM(I).LT.OUTMIN) OUTMIN=RAM(I)
                IF(RAM(I).GT.OUTMAX) OUTMAX=RAM(I)
  42          CONTINUE
              CALL FORM1(OUTMIN,OUTMAX,FORMA2(15:22))
              FORMA2(21:24)= '1X),'
              WRITE(LU,FORMA2)
     *              ' ''',I3,''' ',(RAM(I2+I),I=1,NQ),'/'
            ENDIF
  44      CONTINUE
          WRITE(LU,'(A)') '/'
          CLOSE(LU)
        ENDIF
C       Writing the triangles:
        IF (FTRGLN(I1).NE.' ') THEN
          OPEN(LU,FILE=FTRGLN(I1))
          DO 52, I2=NVRTX,NTRGL-3-NOUT,3+NOUT
            I3=IRAM(I2+3+I1)
            IF (I3.NE.IUNDEF) THEN
              WRITE(LU,FORMA1) IRAM((NQ+NOUT)*(IRAM(I2+1)-1)+NQ+I1),' ',
     *                         IRAM((NQ+NOUT)*(IRAM(I2+2)-1)+NQ+I1),' ',
     *                         IRAM((NQ+NOUT)*(IRAM(I2+3)-1)+NQ+I1),' /'
            ENDIF
  52      CONTINUE
          CLOSE(LU)
        ENDIF
  60  CONTINUE
C
      WRITE(*,'(A)') '+TRGLSORT: Done.                 '
C
      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