C
C Program SRFWRL to convert triangulated or polygonated surfaces into
C the Virtual Reality Modeling Language or GOCAD representation
C
C Version: 5.60
C Date: 2002, May 17
C
C Coded by: Ludek Klimes & Vaclav Bucha
C     Department of Geophysics, Charles University Prague,
C     Ke Karlovu 3, 121 16 Praha 2, Czech Republic,
C     E-mails: klimes@seis.karlov.mff.cuni.cz
C              bucha@seis.karlov.mff.cuni.cz
C
C References:
C     
C     VRML (Virtual Reality Modeling Language) version 1.0C
C     
C     VRML97 (Virtual Reality Modeling Language ISO/IEC 14772)
C     
C     GOCAD
C     
C     Persistence of Vision scene description language, version 3.1
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 Data specifying input files:
C     VRTX='string'... Name of the file with vertices of the polygons.
C             Description of file VRTX
C             Default: VRTX='vrtx.out'
C     TRGL='string'... Name of the file describing the triangles or
C             polygons.  Triangles are recommended (and obligatory if
C             VRML='GOCAD').
C             Description of file TRGL
C             Default: TRGL='trgl.out'
C     COLORS='string'... Name of the file containing the data describing
C             the colour map.
C             Description of file COLORS
C             Default: COLORS='hsv.dat'
C Input/output file:
C     WRL='string'... Name of the file to be supplemented with surfaces
C             or to be copied to the beginning of the output file.
C             If the filename is blank, output file starts from a
C             scratch (mostly not reasonable).
C             The default name of the output file is equal to WRL.
C             It is recommended to specify WRL rather than to use
C             the default name.
C             Default: WRL='out.wrl'
C     WRLOUT='string'... Name of the output file if different from WRL.
C             Default: WRLOUT=WRL
C Data specifying the form of the output file:
C     VRML='string'... Virtual reality scene description language.
C             VRML='VRML1': VRML (Virtual Reality Modeling Language)
C                           version 1.0.
C             VRML='VRML2': VRML97 according to ISO/IEC 14772 standard.
C             VRML='GOCAD': GOCAD description of surfaces (TSurf).
C             VRML='POV':   POV (Persistence Of Vision) scene
C                           description language, version 3.1.
C                           This option has not been used and is thus
C                           poorly debugged.
C             Default: VRML='VRML2' (recommended if not using GOCAD)
C     NAME='string'... String containing the GOCAD name of the surface.
C             Be sure to select different names for all objects within
C             the GOCAD file.
C             The same name is used for the corresponding colour scale,
C             written if KOLSRF is positive.
C             Used only if VRML='GOCAD'.  Obligatory parameter, must be
C             specified and cannot be blank if VRML='GOCAD'.
C Data specifying the values to be scaled in colours:
C     KOLSRF=integer... If zero, all surfaces will have the same colour
C             given by parameters R, G, B.  If positive, the values in
C             KOLSRF-th column of input file VRTX will be colour-coded
C             at each vertex of each triangle or polygon of the surface.
C             If VRML.NE.'GOCAD', this setting may be modified by
C             parameters KOLPOS and KOLNEG.
C             If both KOLPOS and KOLNEG are specified, KOLSRF is used
C             only if VRML='GOCAD'.
C             Default: KOLSRF=7
C     KOLPOS=integer... Analogous to KOLSRF, but applies just to the
C             positive side of the the surface.
C             Not used if VRML='GOCAD'.
C             Default: KOLPOS=KOLSRF
C     KOLNEG=integer... Analogous to KOLSRF, but applies just to the
C             negative side of the the surface.
C             Not used if VRML='GOCAD'.
C             Default: KOLNEG=KOLSRF
C     PROPERTIES='string'... String containing names of properties
C             corresponding to values Z1,Z2,Z3,V1,...,VN (see file
C             VRTX) which may be used to control the
C             colour of the surface.  The names are separated by blanks.
C             If the number of names is smaller than the number of
C             values, the leftmost values are considered.  PROPERTIES
C             must be specified if VRML='GOCAD' and KOLSRF is positive.
C             If KOLSRF is 1, 2 or 3, the last name is assumed to denote
C             the KOLSRFth coordinate rather than the quantity in the
C             corresponding column, and the value of the coordinate
C             copied into that column.
C             If PROPERTIES=' ', no values are considered and GOCAD atom
C             VRTX is used for the vertices (otherwise, GOCAD atom PVRTX
C             is used).
C             Used only if VRML='GOCAD'.
C             Default: PROPERTIES=' '
C Data specifying the colour scale:
C     VADD=real, VMUL=real, VPER=real, VREF=real, CREF=real, CREF1=real,
C     CREF2=real, CREF3=real, etc... Refer to file
C             colors.for.
C     R=real, G=real, B=real... Float numbers between 0 and 1 specifying
C             the colour of the surfaces if KOLSRF=0 or KOLPOS=0 or
C             KOLNEG=0.
C             Defaults: R=1, G=1, B=1 (white)
C     TRANSP=real... Transparency of the surfaces (sometimes called
C             transmit).  Values from 0 to 1.
C             Default: TRANSP=0.
C     AMBIENT=real... Float number between 0 and 1 specifying the
C             intensity of the ambient light.  The colour of the ambient
C             light is assumed white.  Applied to the surfaces only if
C             VRML='vrml1'.  Otherwise, the ambient light source of
C             intensity AMBIENT is prescribed by program
C             iniwrl.for.
C             Not used if VRML='GOCAD'.
C             Default:  AMBIENT=0.20 (default for VRML materials)
C     SPECULAR=real... Intensity of the specular reflections from
C             glossy surfaces.  Values from 0 to 1.
C             Not used if VRML='GOCAD'.
C             Default: SPECULAR=0 (default for VRML materials)
C     SHININESS=real... Shininess of the surfaces (sometimes called
C             transmit).  Values from 0 to 1.
C             Not used if VRML='GOCAD'.
C             Default: SHININESS=0.20 (default for VRML materials)
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,Z1,Z2,Z3,V1,...,VN/
C     'NAME'... Name of the vertex.  Not considered.  May be blank.
C     X1,X2,X3... Coordinates of the vertex.
C     Z1,Z2,Z3... Normal to the surface at the vertex.  The normals are
C             used for shading of the surface if VRML='VRML1' or
C             VRML='VRML2'.  If at least one normal is zero, shading
C             corresponds  to flat triangles or polygons.
C             Normals are transmitted to the GOCAD file if VRML='GOCAD'
C             and parameter PROPERTIES is specified, but do not
C             influence the surface appearance.
C     V1,...,VN...Optional values which may be used to control the
C             colour of the surface.
C     /...    None to several values terminated by a slash.
C (3) / or end of file.
C
C                                                    
C Input file TRGL with the triangles or polygons:
C (1) For each triangle data (1.1), or for each polygon data (1.2):
C (1.1) I1,I2,I3,/
C     I1,I2,I3... Indices of 3 vertices of the triangle, right-handed
C             with respect to the given surface normals.
C             The vertices in file VRTX are indexed by positive integers
C             according to their order.
C             For polygon, three indices I1,I2,I3 are replaced with more
C             ones.
C     /...    List of vertices is terminated by a slash.
C (1.2) I1,I2,...,IN,/
C     I1,I2,...,IN... Indices of N vertices of the polygon.
C             The vertices in file VRTX are indexed by positive integers
C             according to their order.
C     /...    List of vertices must be terminated by a slash.
C (2) / or end of file.
C
C=======================================================================
C
C Common block /RAMC/:
      INCLUDE 'ram.inc'
C     ram.inc
C
      INTEGER IRAM(MRAM)
      EQUIVALENCE (IRAM,RAM)
C
C.......................................................................
C
C     External functions and subroutines:
      EXTERNAL LENGTH,RSEP1,RSEP3T,RSEP3I,ERROR,FORM2,COLOR1,COLOR2
      INTEGER  LENGTH
C
C     Filenames and parameters:
      CHARACTER*80 FSEP,FVRTX,FTRGL,FCOLS,FIN,FOUT
      INTEGER LU1,LU2,IUNDEF,MVRTX,MQ
      PARAMETER (LU1=1,LU2=2,IUNDEF=-999999,MVRTX=99,MQ=30)
C     MVRTX...  Maximum number of vertices of a single polygon.
C
C     Other variables:
      CHARACTER*(8+8*MQ) FORMAT
      CHARACTER*5   VRML
      CHARACTER*255 NAME,PROP,TEXT
      LOGICAL LNORM
      INTEGER KOLSRF,KOLPOS,KOLNEG,KQ,NQ
      INTEGER NVRTX,NPLGN,IREF,IRGB,I0,I1,I2,I,N
      REAL AMBI,TRANSP,SPEC,SHIN,RED,GREEN,BLUE
      REAL OUTMIN(MQ),OUTMAX(MQ),R,G,B,AUX,AUXA(1)
C     LNORM.. Says whether the surface normals are specified.
C
C.......................................................................
C
C     Reading main input data:
      WRITE(*,'(A)') '+SRFWRL: Enter input filename: '
      FSEP=' '
      READ (*,*) FSEP
      IF(FSEP.EQ.' ') THEN
C       SRFWRL-08
        CALL ERROR('SRFWRL-08: No input file specified')
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.
      END IF
      CALL RSEP1(LU1,FSEP)
      WRITE(*,'(A)') '+SRFWRL: Working...            '
C
C     Reading input and output filenames:
      CALL RSEP3T('VRTX'  ,FVRTX,'vrtx.out')
      CALL RSEP3T('TRGL'  ,FTRGL,'trgl.out')
      CALL RSEP3T('COLORS',FCOLS,'hsv.dat' )
      CALL RSEP3T('WRL'   ,FIN  ,'out.wrl' )
      CALL RSEP3T('WRLOUT',FOUT ,FIN       )
      CALL RSEP3T('VRML'  ,VRML ,'VRML2'   )
      CALL LOWER(VRML)
C
C     Reading input parameters for surface appearance:
      CALL RSEP3I('KOLSRF',KOLSRF,7)
      CALL RSEP3I('KOLPOS',KOLPOS,KOLSRF)
      CALL RSEP3I('KOLNEG',KOLNEG,KOLSRF)
      CALL RSEP3R('AMBIENT'  ,AMBI  ,0.20)
      CALL RSEP3R('TRANSP'   ,TRANSP,0.00)
      CALL RSEP3R('SPECULAR' ,SPEC  ,0.00)
      CALL RSEP3R('SHININESS',SHIN  ,0.20)
      CALL RSEP3R('R'        ,RED   ,1.)
      CALL RSEP3R('G'        ,GREEN ,1.)
      CALL RSEP3R('B'        ,BLUE  ,1.)
C
C     Opening the output file and writing its beginning:
      CALL WRL1(LU1,LU2,FIN,FOUT,VRML,1)
C
C     Writing the prolog for the surface (part 1):
      IF (VRML.EQ.'vrml1') THEN
        CONTINUE
      ELSE IF (VRML.EQ.'vrml2') THEN
        CONTINUE
      ELSE IF (VRML.EQ.'gocad') THEN
        KOLPOS=KOLSRF
        KOLNEG=KOLSRF
        CALL RSEP3T('NAME',NAME,' ')
C       Subroutine WRL has already checked that NAME is not blank.
        WRITE(LU2,'(A)')
     *   'GOCAD TSurf 1.0'
        WRITE(LU2,'(2A)')
     *   'HDR name:',NAME(1:LENGTH(NAME))
        WRITE(LU2,'(A)')
     *   'HDR *visible:true'
        CALL RSEP3T('PROPERTIES',PROP,' ')
        I0=1
        KQ=3
        DO 11 I=1,LEN(PROP)-1
          IF (PROP(I:I).EQ.' '.AND.PROP(I+1:I+1).NE.' ') THEN
            I0=I+1
          END IF
          IF (PROP(I:I).NE.' '.AND.PROP(I+1:I+1).EQ.' ') THEN
            KQ=KQ+1
            IF (KQ.EQ.KOLSRF.OR.(1.LE.KOLSRF.AND.KOLSRF.LE.3)) THEN
              I1=I0
              I2=I
            END IF
          END IF
   11   CONTINUE
        IF (KOLSRF.LE.0) THEN
          WRITE(LU2,'(3(A,F4.2))')
     *     'HDR color: ',RED,' ',GREEN,' ',BLUE
          WRITE(LU2,'(A,F4.2)')
     *     'HDR *solid*transparency:',TRANSP
        ELSE
          IF (KQ.LT.KOLSRF.OR.KQ.LT.4) THEN
C           SRFWRL-09
            CALL ERROR('SRFWRL-09: GOCAD property name not specified')
C           If KOLSRF is not zero, list PROPERTIES of property names
C           must contain MAX(1,KOLSRF-3) names at the least, see the
C           description of the input data.
          END IF
          WRITE(LU2,'(A)')
     *     'HDR *painted:true'
     *    ,'HDR *shaded_painted:true'
     *    ,'HDR *precise_painted:true'
          WRITE(LU2,'(2A)')
     *     'HDR *painted*variable:',PROP(I1:I2)
        END IF
        IF (KQ.GT.3) THEN
          WRITE(LU2,'(2A)')
     *     'PROPERTIES ',PROP(1:LENGTH(PROP))
        END IF
        IF (KOLSRF.NE.0) THEN
          IF (LENGTH(PROP)+(KQ-3)*LENGTH(NAME).GT.LEN(TEXT)) THEN
C           SRFWRL-10
            CALL ERROR('SRFWRL-10: Too long property class names')
C           Each property class name is composed of the object name
C           given by input parameter NAME and the property name.
C           The property names are given by input parameter PROPERTIES.
C           All property class names should fit into character variable
C           TEXT.  The length of TEXT thus should not be smaller than
C           the length of the value of PROPERTIES, plus the number of
C           properties times the length of the value of NAME.
          END IF
          I0=0
          DO 12 I=1,LENGTH(PROP)
            IF (I.EQ.1.AND.PROP(1:1).NE.' ') THEN
              TEXT(I0+1:I0+LENGTH(NAME))=NAME(1:LENGTH(NAME))
              I0=I0+LENGTH(NAME)
            END IF
            I0=I0+1
            TEXT(I0:I0)=PROP(I:I)
            IF (PROP(I:I).EQ.' '.AND.PROP(I+1:I+1).NE.' ') THEN
              TEXT(I0+1:I0+LENGTH(NAME))=NAME(1:LENGTH(NAME))
              I0=I0+LENGTH(NAME)
            END IF
   12     CONTINUE
          WRITE(LU2,'(2A)')
     *     'PROPERTY_CLASSES ',TEXT(1:I0)
          WRITE(LU2,'(4A)')
     *    'PROPERTY_CLASS_HEADER ',NAME(1:LENGTH(NAME)),PROP(I1:I2),' {'
C         The output file now waits for the colour scale.
        END IF
C       KQ is the number of coordinates and properties at each point.
      ELSE
C       SRFWRL-11
        CALL ERROR('SRFWRL-11: No valid string in VRML')
C       Valid string specifying the form of the output file is:
C       VRML='VRML1' or 'VRML2' or 'GOCAD' or 'POV'.  Default and
C       recommended value is 'VRML2'.
      END IF
C
C     Reading vertices:
      LNORM=.TRUE.
      IF(VRML.EQ.'gocad') THEN
        NQ=KQ
      ELSE
        KQ=MAX0(6,KOLPOS,KOLNEG)
        IF(KOLPOS.EQ.0.AND.KOLNEG.EQ.0) THEN
          NQ=6
        ELSE IF(KOLPOS.EQ.KOLNEG) THEN
          NQ=7
        ELSE
          NQ=8
        END IF
C       Values to be displayed will be shifted to the 7th and 8th column
      END IF
      IF(NQ.GT.MQ) THEN
C       SRFWRL-12
        CALL ERROR('SRFWRL-12: Too small arrays OUTMIN and OUTMAX')
      END IF
      OPEN(LU1,FILE=FVRTX,STATUS='OLD')
      READ(LU1,*) (TEXT,I=1,20)
      NVRTX=0
   20 CONTINUE
        IF(NVRTX+KQ.GT.MRAM) THEN
C         SRFWRL-01
          CALL ERROR('SRFWRL-01: Too small array RAM')
        END IF
        TEXT='$'
        DO 21 I=NVRTX+2,NVRTX+KQ
          RAM(I)=0.
   21   CONTINUE
        READ(LU1,*,END=29) TEXT,(RAM(I),I=NVRTX+1,NVRTX+KQ)
        IF(TEXT.EQ.'$') THEN
          GO TO 29
        END IF
C       Relocating the values to be displayed to the 7th and 8th columns
        IF(VRML.EQ.'gocad') THEN
          IF(1.LE.KOLSRF.AND.KOLSRF.LE.3) THEN
            RAM(NVRTX+KQ)=RAM(NVRTX+KOLSRF)
          END IF
        ELSE
          IF(KOLNEG.GT.0) THEN
            AUX=RAM(NVRTX+KOLNEG)
          END IF
          IF(KOLPOS.GT.0) THEN
            RAM(NVRTX+7)=RAM(NVRTX+KOLPOS)
          END IF
          IF(KOLNEG.GT.0.AND.KOLPOS.NE.KOLNEG) THEN
            RAM(NVRTX+8)=AUX
          END IF
        END IF
C       Normalizing the normal
        AUX=SQRT(RAM(NVRTX+4)**2+RAM(NVRTX+5)**2+RAM(NVRTX+6)**2)
        IF(AUX.GT.0.) THEN
          AUX=0.999/AUX
          RAM(NVRTX+4)=RAM(NVRTX+4)*AUX
          RAM(NVRTX+5)=RAM(NVRTX+5)*AUX
          RAM(NVRTX+6)=RAM(NVRTX+6)*AUX
        ELSE
          LNORM=.FALSE.
        END IF
C       Determining the minimum and maximum values
        IF(NVRTX.EQ.0) THEN
          DO 22 I=1,NQ
            OUTMIN(I)=RAM(NVRTX+I)
            OUTMAX(I)=RAM(NVRTX+I)
   22     CONTINUE
        ELSE
          DO 23 I=1,NQ
            OUTMIN(I)=AMIN1(OUTMIN(I),RAM(NVRTX+I))
            OUTMAX(I)=AMAX1(OUTMAX(I),RAM(NVRTX+I))
   23     CONTINUE
        END IF
C       Number of storage locations in RAM used for the vertices
        NVRTX=NVRTX+NQ
      GO TO 20
   29 CONTINUE
      CLOSE(LU1)
C     NVRTX is the number of storage locations in RAM used for vertices
      IF(VRML.NE.'gocad') THEN
        IF(KOLNEG.NE.0) THEN
          IF(KOLPOS.EQ.KOLNEG) THEN
            KOLNEG=7
          ELSE
            KOLNEG=8
          END IF
        END IF
        IF(KOLPOS.NE.0) THEN
          KOLPOS=7
        END IF
        IF(NQ.GE.8) THEN
          OUTMIN(7)=AMIN1(OUTMIN(7),OUTMIN(8))
          OUTMAX(7)=AMAX1(OUTMAX(7),OUTMAX(8))
        END IF
      END IF
C     Values to be displayed have been shifted to the 7th or 8th columns
C
C     Determining the colour map:
      IF(KOLPOS.GT.0.OR.KOLNEG.GT.0) THEN
        IF(VRML.EQ.'gocad') THEN
          CALL COLOR1(LU1,MRAM-NVRTX,IRAM(NVRTX+1),RAM(NVRTX+1),
     *                                  1,OUTMIN(KOLSRF),OUTMAX(KOLSRF))
          WRITE(LU2,'(2A)')
     *     ' *colormap:',NAME(1:LENGTH(NAME))
          FORMAT='(A,'
          CALL FORM2(1,OUTMIN(KOLSRF),OUTMAX(KOLSRF),FORMAT(4:11))
          FORMAT(9:11)=')  '
          IF(OUTMAX(KOLSRF).GT.OUTMIN(KOLSRF)) THEN
            WRITE(LU2,FORMAT)
     *       ' *low_clip: ',OUTMIN(KOLSRF)
     *      ,' *high_clip:',OUTMAX(KOLSRF)
          ELSE
            WRITE(LU2,FORMAT)
     *       ' *low_clip: ',OUTMIN(KOLSRF)
     *      ,' *high_clip:',OUTMIN(KOLSRF)+1.
          END IF
          WRITE(LU2,'(4A)')
     *     ' *colormap*',NAME(1:LENGTH(NAME)),'*colors: ',CHAR(92)
          AUX=(OUTMAX(KOLSRF)-OUTMIN(KOLSRF))/255.
          DO 31 I=0,255
            AUXA(1)=OUTMIN(KOLSRF)+FLOAT(I)*AUX
            CALL COLOR2(MRAM-NVRTX,IRAM(NVRTX+1),RAM(NVRTX+1),
     *                                                     1,AUXA,R,G,B)
            IF (I.LT.255) THEN
              WRITE(LU2,'(I5,3(1X,F4.2),2A)')
     *         I,R,G,B,' ',CHAR(92)
            ELSE
              WRITE(LU2,'(I5,3(1X,F4.2),2A)')
     *         I,R,G,B
            END IF
   31     CONTINUE
          IF(TRANSP.GT.0.) THEN
            WRITE(LU2,'(2A)')
     *       ' *colormap*alphas: ',CHAR(92)
            DO 32 I=0,255
              IF (I.LT.255) THEN
                WRITE(LU2,'(I5,1X,F4.2,2A)')
     *           I,TRANSP,' ',CHAR(92)
              ELSE
                WRITE(LU2,'(I5,1X,F4.2,2A)')
     *           I,TRANSP
              END IF
   32       CONTINUE
          END IF
          WRITE(LU2,'(A)')
     *     '}'
        ELSE IF (VRML.EQ.'pov') THEN
          AUX=0.01/SHIN
          WRITE(LU2,'(A)')
     *     '#default {'
          WRITE(LU2,'(A,2(F4.2,A))')
     *     '  finish { ambient 1.00 specular ',SPEC,
     *                                            ' roughness ',AUX,' }'
          WRITE(LU2,'(A)')
     *     '  pigment {'
     *    ,'    color_map {'
          CALL COLOR3(MRAM-NVRTX,IRAM(NVRTX+1),RAM(NVRTX+1),1,IREF,IRGB)
          I=NVRTX+1+IRAM(NVRTX+1)
          IREF=NVRTX+IREF
          IRGB=NVRTX+IRGB
          DO 57 I2=1,IRAM(NVRTX+2)-IRAM(NVRTX+1)
            WRITE(LU2,'(A,F8.6,A,4(F4.2,A))')
     *       '      [',RAM(I+I2),' rgbt <',
     *                     (RAM(IRGB+I1),',',I1=3*I2-2,3*I2),TRANSP,'>]'
   57     CONTINUE
          WRITE(LU2,'(A)')
     *     '    }'
     *    ,'  }'
     *    ,'}'
          WRITE(LU2,'(A,G13.6,A)')
     *     '#declare CREF = ',RAM(IREF+1),';'
     *    ,'#declare VREF = ',RAM(IREF+2),';'
     *    ,'#declare VPER = ',RAM(IREF+3),';'
        ELSE
          CALL COLOR1(LU1,MRAM-NVRTX,IRAM(NVRTX+1),RAM(NVRTX+1),
     *                                            1,OUTMIN(7),OUTMAX(7))
        END IF
      END IF
C
C     Writing the prolog for the surface (part 2):
      IF (VRML.EQ.'vrml1') THEN
        WRITE(LU2,'(A)')
     *   'DEF SurfaceMaterial Material {'
        WRITE(LU2,'(3(A,F4.2))')
     *   '  diffuseColor    ',RED,' ',GREEN,' ',BLUE
     *  ,'  ambientColor    ',RED*AMBI,' ',GREEN*AMBI,' ',BLUE*AMBI
     *  ,'  specularColor   ',SPEC,' ',SPEC,' ',SPEC
        WRITE(LU2,'(A,F4.2)')
     *   '  shininess       ',SHIN
     *  ,'  transparency    ',TRANSP
        WRITE(LU2,'(A)')
     *   '  emissiveColor    0.00 0.00 0.00'
     *  ,'}'
        WRITE(LU2,'(A)')
     *   'Separator {'
     *  ,'USE SurfaceMaterial'
        IF(LNORM) THEN
          WRITE(LU2,'(A)') 'NormalBinding { value PER_VERTEX }'
        ELSE
          WRITE(LU2,'(A)') 'NormalBinding { value PER_FACE }'
        END IF
      ELSE IF (VRML.EQ.'vrml2') THEN
        WRITE(LU2,'(A)')
     *   'Shape {'
     *  ,'  appearance DEF SurfaceAppearance Appearance {'
     *  ,'    material Material {'
        WRITE(LU2,'(3(A,F4.2))')
     *   '      diffuseColor     ',RED,' ',GREEN,' ',BLUE
     *  ,'      specularColor    ',SPEC,' ',SPEC,' ',SPEC
        WRITE(LU2,'(A,F4.2)')
     *   '      shininess        ',SHIN
     *  ,'      transparency     ',TRANSP
        WRITE(LU2,'(A)')
     *   '      ambientIntensity 1.00'
     *  ,'      emissiveColor    0.00 0.00 0.00'
     *  ,'    }'
     *  ,'  }'
     *  ,'}'
     *  ,'Surface {'
     *  ,'appearance USE SurfaceAppearance'
      ELSE IF (VRML.EQ.'pov') THEN
        WRITE(LU2,'(A,I6,A)')
     *   '#declare NVRTX =',NVRTX/NQ,';'
        WRITE(LU2,'(A)')
     *   '#declare PTS = array[NVRTX][7]'
     *  ,'#declare IVRTX = 0;'
     *  ,'#macro VRTX(X1,X2,X3,Z1,Z2,Z3,V1)'
     *  ,'  #declare PTS[IVRTX][0] = X1;'
     *  ,'  #declare PTS[IVRTX][1] = X2;'
     *  ,'  #declare PTS[IVRTX][2] = X3;'
     *  ,'  #declare PTS[IVRTX][3] = Z1;'
     *  ,'  #declare PTS[IVRTX][4] = Z2;'
     *  ,'  #declare PTS[IVRTX][5] = Z3;'
     *  ,'  #declare PTS[IVRTX][6] = V1;'
     *  ,'  #declare IVRTX = IVRTX + 1;'
     *  ,'#end'
     *  ,'#macro TRGL(I1,I2,I3)'
     *  ,'  #local X1=;'
     *  ,'  #local X2=;'
     *  ,'  #local X3=;'
     *  ,'  #local Z1=;'
     *  ,'  #local Z2=;'
     *  ,'  #local Z3=;'
     *  ,'  #local V1=PTS[I1][6]-PTS[I3][6];'
     *  ,'  #local V2=PTS[I2][6]-PTS[I3][6];'
     *  ,'  #local V3=           PTS[I3][6];'
     *  ,'  #if (V1=0 & V2=0)'
     *  ,'    #local V1=VPER/999999;'
     *  ,'  #end'
     *  ,'  #local D1=X1-X3;'
     *  ,'  #local D2=X2-X3;'
     *  ,'  #local D11=vdot(D1,D1);'
     *  ,'  #local D12=vdot(D1,D2);'
     *  ,'  #local D22=vdot(D2,D2);'
     *  ,'  #local D  =D11*D22-D12*D12;'
     *  ,'  #local G =(D1*(D22*V1-D12*V2)+D2*(-D12*V1+D11*V2))/D;'
     *  ,'  #local GN= vlength(G);'
     *  ,'  #local G0= G*VPER/GN/GN;'
     *  ,'  #local G1= V2*D1-V1*D2;'
     *  ,'  #local G2= vcross(G0,G1);'
     *  ,'  smooth_triangle {'
     *  ,'    X1,Z1,X2,Z2,X3,Z3'
     *  ,'    texture {'
     *  ,'      pigment {'
     *  ,'        gradient x'
     *  ,'        translate ((VREF-V3)/VPER-CREF-100)*x'
     *  ,'        matrix '
     *  ,'        translate X3'
     *  ,'      }'
     *  ,'    }'
     *  ,'  }'
     *  ,'#end'
      END IF
C
C     Writing the vertices:
      IF (VRML.EQ.'vrml1') THEN
        WRITE(LU2,'(A)') 'Coordinate3 { point ['
      ELSE IF (VRML.EQ.'vrml2') THEN
        WRITE(LU2,'(A)') 'point ['
      END IF
C     ------
      IF (VRML.EQ.'vrml1'.OR.VRML.EQ.'vrml2') THEN
        FORMAT='('
        CALL FORM2(3,OUTMIN(1),OUTMAX(1),FORMAT(2:25))
        DO 60 I=1,NVRTX,NQ
          WRITE(LU2,FORMAT) RAM(I),' ',RAM(I+1),' ',RAM(I+2),','
   60   CONTINUE
      ELSE IF (VRML.EQ.'gocad') THEN
C       Writing the vertices with normals and values:
        FORMAT='(A,I0,A,'
        FORMAT(5:5)=CHAR(ICHAR('1')+INT(ALOG10(FLOAT(NVRTX/NQ)+0.5)))
        IF (NQ.EQ.3) THEN
          CALL FORM2(3,OUTMIN(1),OUTMAX(1),FORMAT(9:8+8*3))
          DO 62 I0=1,NVRTX,NQ
            WRITE(LU2,FORMAT) 'VRTX ',I0/NQ+1,(' ',RAM(I),I=I0,I0+2)
   62     CONTINUE
        ELSE
          CALL FORM2(NQ,OUTMIN(1),OUTMAX(1),FORMAT(9:8+8*NQ))
          DO 63 I0=1,NVRTX,NQ
            WRITE(LU2,FORMAT) 'PVRTX ',I0/NQ+1,(' ',RAM(I),I=I0,I0+NQ-1)
   63     CONTINUE
        END IF
      ELSE IF (VRML.EQ.'pov') THEN
C       Writing the vertices with normals and values:
        IF(KOLNEG.NE.KOLPOS) THEN
C         SRFWRL-51
          CALL WARN('SRFWRL-51: POV surface sides differently coloured')
C         POV scene description language does not allow for different
C         colours at the positive and negative side of a surface.
        END IF
        FORMAT='(A,'
        CALL FORM2(3,OUTMIN(1),OUTMAX(1),FORMAT(4:27))
        FORMAT(27:38)=',3(F5.3,A),'
        CALL FORM2(1,OUTMIN(7),OUTMAX(7),FORMAT(39:46))
        DO 65 I=1,NVRTX,NQ
          WRITE(LU2,FORMAT) 'VRTX(',(RAM(I1),',',I1=I,I+5),RAM(I+6),')'
   65   CONTINUE
      END IF
C     ------
      IF (VRML.EQ.'vrml1') THEN
        WRITE(LU2,'(A)') '] }'
      ELSE IF (VRML.EQ.'vrml2') THEN
        WRITE(LU2,'(A)') ']'
      END IF
C
C     Writing the right-handed normals (positive surface side):
      IF(LNORM) THEN
        IF (VRML.EQ.'vrml1') THEN
          WRITE(LU2,'(A)') 'DEF SurfaceNormal Normal { vector ['
        ELSE IF (VRML.EQ.'vrml2') THEN
          WRITE(LU2,'(A)') 'normalPos Normal { vector ['
        END IF
C       ------
        IF (VRML.EQ.'vrml1'.OR.VRML.EQ.'vrml2') THEN
          FORMAT='(3(F5.3,A))'
          DO 66 I=4,NVRTX,NQ
            WRITE(LU2,FORMAT) RAM(I),' ',RAM(I+1),' ',RAM(I+2),','
   66     CONTINUE
        END IF
C       ------
        IF (VRML.EQ.'vrml1') THEN
          WRITE(LU2,'(A)') '] }'
        ELSE IF (VRML.EQ.'vrml2') THEN
          WRITE(LU2,'(A)') '] }'
        END IF
      END IF
C
C     Writing the left-handed normals (negative surface side):
      IF(LNORM) THEN
        IF (VRML.EQ.'vrml1') THEN
          WRITE(LU2,'(A)') 'Normal { vector ['
        ELSE IF (VRML.EQ.'vrml2') THEN
          WRITE(LU2,'(A)') 'normalNeg Normal { vector ['
        END IF
C       ------
        IF (VRML.EQ.'vrml1'.OR.VRML.EQ.'vrml2') THEN
          DO 67 I=4,NVRTX,NQ
            WRITE(LU2,FORMAT) -RAM(I),' ',-RAM(I+1),' ',-RAM(I+2),','
   67     CONTINUE
        END IF
C       ------
        IF (VRML.EQ.'vrml1') THEN
          WRITE(LU2,'(A)') '] }'
        ELSE IF (VRML.EQ.'vrml2') THEN
          WRITE(LU2,'(A)') '] }'
        END IF
      END IF
C
C     Writing the colours of the positive surface side:
      IF(KOLPOS.GT.0) THEN
        IF (VRML.EQ.'vrml1') THEN
          WRITE(LU2,'(A)') 'DEF SurfaceColor Material { diffuseColor ['
        ELSE IF (VRML.EQ.'vrml2') THEN
          WRITE(LU2,'(A)') 'colorPos DEF SurfaceColor Color { color ['
        END IF
C       ------
        IF (VRML.EQ.'vrml1'.OR.VRML.EQ.'vrml2') THEN
          DO 71 I=KOLPOS,NVRTX,NQ
            CALL COLOR2(MRAM-NVRTX,IRAM(NVRTX+1),RAM(NVRTX+1),
     *                                                   1,RAM(I),R,G,B)
            WRITE(LU2,'(3(F4.2,A))') R,' ',G,' ',B,','
   71     CONTINUE
        END IF
C       ------
        IF (VRML.EQ.'vrml1') THEN
          WRITE(LU2,'(A)') '] }'
        ELSE IF (VRML.EQ.'vrml2') THEN
          WRITE(LU2,'(A)') '] }'
        END IF
      END IF
C
C     Writing the colours of the negative surface side:
      IF(KOLNEG.GT.0) THEN
        IF(KOLNEG.EQ.KOLPOS) THEN
          IF (VRML.EQ.'vrml1') THEN
            CONTINUE
          ELSE IF (VRML.EQ.'vrml2') THEN
            WRITE(LU2,'(A)') 'colorNeg USE SurfaceColor'
          END IF
        ELSE
          IF (VRML.EQ.'vrml1') THEN
            WRITE(LU2,'(A)') 'Material { diffuseColor ['
          ELSE IF (VRML.EQ.'vrml2') THEN
            WRITE(LU2,'(A)') 'colorNeg Color { color ['
          ELSE IF (VRML.EQ.'pov') THEN
          END IF
C         ------
          IF (VRML.EQ.'vrml1'.OR.VRML.EQ.'vrml2') THEN
            DO 72 I=KOLNEG,NVRTX,NQ
              CALL COLOR2(MRAM-NVRTX,IRAM(NVRTX+1),RAM(NVRTX+1),
     *                                                   1,RAM(I),R,G,B)
              WRITE(LU2,'(3(F4.2,A))') R,' ',G,' ',B,','
   72       CONTINUE
          END IF
C         ------
          IF (VRML.EQ.'vrml1') THEN
            WRITE(LU2,'(A)') '] }'
          ELSE IF (VRML.EQ.'vrml2') THEN
            WRITE(LU2,'(A)') '] }'
          END IF
        END IF
      END IF
C
C     Reading the polygons (usually triangles):
      DO 81 I=1,MRAM
        IRAM(I)=0
   81 CONTINUE
      OPEN(LU1,FILE=FTRGL)
      NPLGN=0
   82 CONTINUE
        IF(NPLGN+MVRTX+1.GT.MRAM) THEN
C         SRFWRL-02
          CALL ERROR('SRFWRL-02: Too small array RAM')
        END IF
        IRAM(NPLGN+1)=IUNDEF
        READ(LU1,*,END=89) (IRAM(I),I=NPLGN+1,NPLGN+MVRTX+1)
        IF(IRAM(NPLGN+1).EQ.IUNDEF) THEN
          GO TO 89
        END IF
        DO 83 I=NPLGN+1,NPLGN+MVRTX+1
          IF(IRAM(I).LE.0) THEN
C           Number of polygon vertices
            N=I-1-NPLGN
            GO TO 84
          ELSE IF(IRAM(I).GT.NVRTX/NQ) THEN
C           SRFWRL-03
            WRITE(TEXT,'(A,I6)')'SRFWRL-03: Wrong vertex index:',IRAM(I)
            CALL ERROR(TEXT(1:LENGTH(TEXT)))
          END IF
   83   CONTINUE
C         SRFWRL-04
          CALL ERROR('SRFWRL-04: Too many vertices in polygons')
   84   CONTINUE
        IF(N.LT.3) THEN
C         SRFWRL-52
          CALL WARN('SRFWRL-52: Polygon of less than 3 vertices')
        END IF
C       Checking vertex indices:
        DO 86 I2=NPLGN+1,NPLGN+N
          DO 85 I1=I2+1,NPLGN+N
          IF(IRAM(I2).EQ.IRAM(I1)) THEN
C           SRFWRL-05
            WRITE(TEXT,'(A,I6)')
     *        'SRFWRL-05: The same vertex twice in a polygon:',IRAM(I2)
            CALL ERROR(TEXT(1:LENGTH(TEXT)))
C           All vertices of a polygon must be different.
          END IF
   85     CONTINUE
   86   CONTINUE
C       Terminating polygon by zero
        IF(N.GE.3) THEN
          NPLGN=NPLGN+N+1
          IRAM(NPLGN)=0
        END IF
      GO TO 82
   89 CONTINUE
      CLOSE(LU1)
C
C     Writing the polygons (usually triangles):
      IF(VRML.EQ.'vrml1') THEN
        IF(KOLNEG.GT.0) THEN
          WRITE(LU2,'(A)') 'MaterialBinding { value PER_VERTEX }'
        ELSE
          WRITE(LU2,'(A)') 'MaterialBinding { value OVERALL }'
        END IF
        WRITE(LU2,'(A)') 'ShapeHints {'
        WRITE(LU2,'(A)') '  vertexOrdering CLOCKWISE'
        WRITE(LU2,'(A)') '  shapeType SOLID'
        WRITE(LU2,'(A)') '}'
        WRITE(LU2,'(A)') 'DEF Surface IndexedFaceSet { coordIndex ['
      ELSE IF(VRML.EQ.'vrml2') THEN
        WRITE(LU2,'(A)') 'coordIndex ['
      END IF
C     ------
      N=0
      IF(VRML.EQ.'vrml1'.OR.VRML.EQ.'vrml2') THEN
        FORMAT='(99(I0,A))'
        I=INT(ALOG10(FLOAT(NVRTX/NQ)-0.5))+1
        FORMAT(6:6)=CHAR(ICHAR('0')+I)
        DO 91 I2=1,NPLGN
          IF(IRAM(I2).LE.0) THEN
            WRITE(LU2,FORMAT)
     *                (IRAM(I1)-1,', ',I1=N+1,I2-2),IRAM(I2-1)-1,', -1,'
            N=I2
          END IF
   91   CONTINUE
      ELSE IF (VRML.EQ.'gocad') THEN
        FORMAT='(A,3(A,I0))'
        I=INT(ALOG10(FLOAT(NVRTX/NQ)+0.5))+1
        FORMAT(9:9)=CHAR(ICHAR('0')+I)
        DO 92 I2=1,NPLGN
          IF(IRAM(I2).LE.0) THEN
            IF(I2-N.GT.4) THEN
C             SRFWRL-06
              CALL ERROR('SRFWRL-06: More than 3 vertices in polygon')
C             In this version of the SRFWRL program, only triangles are
C             allowed for GOCAD.  Polygons should be divided into
C             triangles using program 'trgl.for'.
            END IF
            WRITE(LU2,FORMAT) 'TRGL',(' ',IRAM(I1),I1=N+1,I2-1)
            N=I2
          END IF
   92   CONTINUE
      ELSE IF(VRML.EQ.'pov') THEN
        FORMAT='(99(A,I0))'
        I=INT(ALOG10(FLOAT(NVRTX/NQ)-0.5))+1
        FORMAT(8:8)=CHAR(ICHAR('0')+I)
        DO 93 I2=1,NPLGN
          IF(IRAM(I2).LE.0) THEN
            IF(I2-N.GT.4) THEN
C             SRFWRL-07
              CALL ERROR('SRFWRL-07: More than 3 vertices in polygon')
C             In this version of the SRFWRL program, only triangles are
C             allowed for the POV scene description language.  Polygons
C             should be divided into triangles using program 'trgl.for'.
            END IF
            WRITE(LU2,FORMAT)
     *             'TRGL(',(IRAM(I1)-1,',',I1=N+1,I2-2),IRAM(I2-1)-1,')'
            N=I2
          END IF
   93   CONTINUE
      END IF
C     ------
      IF(VRML.EQ.'vrml1') THEN
        WRITE(LU2,'(A)') '] }'
        IF(LNORM) THEN
          WRITE(LU2,'(A)') 'USE SurfaceNormal'
        END IF
        IF(KOLPOS.GT.0) THEN
          WRITE(LU2,'(A)') 'USE SurfaceColor'
          WRITE(LU2,'(A)') 'MaterialBinding { value PER_VERTEX }'
        ELSE
          WRITE(LU2,'(A)') 'MaterialBinding { value OVERALL }'
        END IF
        WRITE(LU2,'(A)') 'ShapeHints {'
        WRITE(LU2,'(A)') '  vertexOrdering COUNTERCLOCKWISE'
        WRITE(LU2,'(A)') '  shapeType SOLID'
        WRITE(LU2,'(A)') '}'
        WRITE(LU2,'(A)') 'USE Surface'
      ELSE IF(VRML.EQ.'vrml2') THEN
        WRITE(LU2,'(A)') ']'
      END IF
C
C     Writing the trailor for the surface:
      IF (VRML.EQ.'vrml1') THEN
        WRITE(LU2,'(A)') '}'
      ELSE IF (VRML.EQ.'vrml2') THEN
        WRITE(LU2,'(A)') '}'
      ELSE IF (VRML.EQ.'gocad') THEN
        WRITE(LU2,'(A)') 'END'
      END IF
      CLOSE(LU2)
      WRITE(*,'(A)') '+SRFWRL: Done.                 '
      STOP
      END
C
C=======================================================================
C
      INCLUDE 'error.for'
C     error.for
      INCLUDE 'sep.for'
C     sep.for
      INCLUDE 'length.for'
C     length.for
      INCLUDE 'forms.for'
C     forms.for
      INCLUDE 'colors.for'
C     colors.for
      INCLUDE 'wrl.for'
C     wrl.for
C
C=======================================================================
C