C PROGRAM INV3 TO UPDATE THE INPUT DATA FOR A FUNCTION DESCRIBING THE
C MODEL.
C
C JUST A PRELIMINARY DEMO VERSION, GENERATING A TABLE OF VALUES OF A
C GIVEN FUNCTION DESCRIBING THE MODEL.  THE GIVEN FUNCTION MAY BE A
C FUNCTION DESCRIBING A SOOTH SURFACE COVERING A STRUCTURAL INTERFACE,
C OR A FUNCTION DESCRIBING THE SPATIAL DISTRIBUTION OF A MATERIAL
C PARAMETER.  THE FUNCTION IS EVALUATED AT GRIDPOINTS OF A GIVEN
C RECTANGULAR GRID OF POINTS.  THE MODEL PARAMETERS (COEFFICIENTS OF
C FUNCTIONS) MAY BE UPDATED BY INCREMENTS RESULTING FROM AN INVERSION OF
C THE SYSTEM OF EQUATIONS GENERATED BY THE INV1 PROGRAM.  CONSEQUENTLY,
C THE TABLES GENERATED BY THE INV3 PROGRAM MAY BE USED TO UPDATE THE
C INPUT DATA FOR THE MODEL BY MEANS OF MANUAL EDITTING.  THE INVERSION
C PROGRAM (READING RESULTS 'DATA' AND 'SOFT' OF PROGRAM INV1, AND
C GENERATING INPUT 'ANSWER' OF PROGRAM INV3) SHOULD BE CONTRIBUTED BY A
C USER.
C
C PROGRAM INV3 ASSUMES ALL MODEL PARAMETERS (COEFFICIENTS) STORED IN THE
C COMMON BLOCK /VALC/ AS IN THE SUBMITTED VERSIONS OF USER-DEFINED MODEL
C SPECIFICATION FORTRAN77 SOURCE CODE FILES 'SRFC.FOR', 'PARM.FOR', AND
C 'VAL.FOR'.  THUS, UNLIKE THE OTHER PARTS OF THE COMPLETE RAY TRACING,
C THE INV3 PROGRAM CANNOT WORK WITH USER'S MODIFICATIONS OF THE
C SUBROUTINES SRFC1, SRFC2, PARM1, AND PARM2.
C
C MAIN INPUT DATA READ FROM THE INTERACTIVE DEVICE (*):
C (1) 'MODEL','ANSWER','GRID',/
C     'MODEL','ANSWER','GRID'... NAMES OF THE INPUT AND OUTPUT
C             FILES DESCRIBED BELOW.
C     'MODEL'... INPUT DATA FILE CONTAINING THE MODEL PARAMETERS.
C     'ANSWER'... UPDATES TO THE MODEL.  IF BLANK, MODEL IS NOT UPDATED.
C     'GRID'..INPUT-OUTPUT FILE:
C             INPUT: SPECIFICATION OF THE GRID IN WHICH A SELECTED
C               FUNCTION IS DISCRETIZED.
C             OUTPUT: FUNCTION VALUES AT GRIDPOINTS ARE ADDED.
C     /...    OBLIGATORY SLASH TO ENABLE FUTURE COMPATIBLE EXTENSIONS.
C     DEFAULT: 'ANSWER'='ANSWER.OUT', 'GRID'='GRID.OUT'.
C
C INPUT FILE 'MODEL':
C INPUT DATA FILE CONTAINING THE MODEL PARAMETERS.  SEE THE FILE
C 'MODEL.FOR' OF THE PACKAGE 'MODEL'.
C
C INPUT FILE 'ANSWER' (SEE PROGRAM INV2):
C (1) NM
C     NM...   NUMBER OF MODEL PARAMETERS.
C             NM=0: INITIAL MODEL IS DISCRETIZED INTO THE GIVEN GRID.
C             DEFAULT: NM=0.
C (2) NM-TIMES THE FOLLOWING DATA:
C     INDM,RM
C     INDM... INDEX OF A MODEL PARAMETER.
C     RM...   INCREMENT OF THE MODEL PARAMETER.
C (3) (CM(I,J),I=1,J),J=1,NM)
C     CM...   MODEL PARAMETER COVARIANCE MATRIX INCLUDING THE SUBJECTIVE
C             PRIOR INFORMATION.
C
C FILE 'GRID':
C INPUT:
C (1) TEXTG,IGROUP
C     IDENTIFICATION OF THE GROUP.
C     TEXTG...A STRING.  IF ITS FIRST CHARACTER IS 'I' OR 'S', THE
C             COMPUTED FUNCTION DESCRIBES A SURFACE.  OTHERWISE, IT
C             DESCRIBES A MATERIAL PARAMETER.
C     IGROUP..INDEX OF A SURFACE, OR OF A COMPLEX BLOCK.
C (2) TEXTF
C     THIS INPUT IS NOT PERFORMED FOR A SURFACE, I.E. IF THE FIRST
C     CHARACTER OF TEXTG (SEE ABOVE) IS 'I' OR 'S'.
C     TEXTF...STRING IDENTIFYING A MATERIAL PARAMETER.
C (3) K1,K2,K3
C     K1,K2,K3... INDICES OF COORDINATES.
C (4) N1,N2,N3
C     N1,N2,N3... NUMBERS OF GRID LINES.
C (5) X1(1),...,X1(N1)
C     THE GRID COORDINATES CORRESPONDING TO THE FIRST INDEPENDENT
C     VARIABLE.
C (6) X2(1),...,X2(N2)
C     THE GRID COORDINATES CORRESPONDING TO THE SECOND INDEPENDENT
C     VARIABLE.
C (7) X3(1),...,X3(N3)
C     THE GRID COORDINATES CORRESPONDING TO THE THIRD INDEPENDENT
C     VARIABLE.
C OUTPUT:
C (8) (((W(I1,I2,I3),I1=1,N1),I2=1,N2),I3=1,N3)
C     THE VALUES OF FUNCTION W AT GRID POINTS. FUNCTION VALUE
C     W(I1,I2,I3) CORRESPONDS TO POINT (X1(I1),X2(I2),X3(I3)).
C (9) 'STANDARD DEVIATIONS ....................' (A STRING)
C (10) (((E(I1,I2,I3),I1=1,N1),I2=1,N2),I3=1,N3)
C     STANDARD DEVIATIONS OF FUNCTION VALUES.
C
C DATE: 1992, DECEMBER 31
C CODED BY LUDEK KLIMES
C
C-----------------------------------------------------------------------
C
C COMMON BLOCK /VALC/:
      INTEGER NPAR
      PARAMETER (NPAR=8191)
      INTEGER IPAR(0:NPAR)
      REAL RPAR(0:NPAR)
      EQUIVALENCE (IPAR,RPAR)
      COMMON/VALC/IPAR
C NONE OF THE STORAGE LOCATIONS OF THE COMMON BLOCK ARE ALTERED.
C
C-----------------------------------------------------------------------
C
C     FILENAMES:
      CHARACTER*80 FILE1,FILE2,FILE3
C
C     LOGICAL UNIT NUMBERS:
      INTEGER LU1,LU2,LU3
      PARAMETER (LU1=11)
      PARAMETER (LU2=12)
      PARAMETER (LU3=13)
C
C     INPUT DATA:
      INTEGER MM
      PARAMETER (MM=256)
      INTEGER NM,INDM(MM)
      REAL CM(MM*(MM+1)/2)
      INTEGER MX
      PARAMETER (MX=100)
      INTEGER IGROUP,K1,K2,K3,N1,N2,N3
      REAL X(MX,3)
C
C     OUTPUT DATA:
      INTEGER ME
      PARAMETER (ME=8000)
      REAL W(MX),E(ME)
C
C     AUXILIARY STORAGE LOCATIONS:
      CHARACTER*3 TEXT
      INTEGER MFUN
      PARAMETER (MFUN=64)
      INTEGER I1,I2,I3,IVAL,IFUN(MFUN),NFUN,IE,II,JJ,IJJ
      REAL COOR(3),UP(10),US(10),AUX0,AUX1,AUX2
      REAL FUN(MFUN),B0I,B1I,B2I,B3I
C
C.......................................................................
C
C     READING MAIN INPUT DATA:
      WRITE(*,'(A)') ' ENTER NAMES OF INPUT AND OUTPUT FILES: '
      FILE1=' '
      FILE2='ANSWER.OUT'
      FILE3='GRID.OUT'
      READ(*,*) FILE1,FILE2,FILE3
      WRITE(*,'(A)') '+                                       '
C
C     READING INPUT DATA FOR MODEL:
      OPEN(LU1,FILE=FILE1,STATUS='OLD')
      CALL MODEL1(LU1)
      CLOSE(LU1)
C
C     UPDATING THE MODEL:
      NM=0
      IF(FILE2.NE.' ') THEN
        OPEN(LU2,FILE=FILE2,STATUS='OLD')
        READ(LU2,*) NM
        IF(NM.GT.MM) THEN
          PAUSE 'ERROR: TOO MANY MODEL PARAMETERS'
          STOP
        END IF
        DO 10 I1=1,NM
          READ(LU2,*) INDM(I1),AUX0
          RPAR(INDM(I1))=RPAR(INDM(I1))+AUX0
   10   CONTINUE
        IF(NM.GT.0) THEN
          READ(LU2,*) (CM(I1),I1=1,NM*(NM+1)/2)
        END IF
        CLOSE(LU2)
      END IF
C
C     READING FUNCTION AND GRID SPECIFICATIONS:
      OPEN(LU3,FILE=FILE3,STATUS='OLD')
      READ(LU3,*) TEXT,IGROUP
      IF(TEXT(1:1).NE.'I'.AND.TEXT(1:1).NE.'S') THEN
        READ(LU3,*) TEXT
      END IF
      READ(LU3,*) K1,K2,K3
      READ(LU3,*) N1,N2,N3
      IF(MAX0(N1,N2,N3).GT.MX) THEN
        PAUSE 'ERROR: TOO MANY GRID LINES'
        STOP
      END IF
      IF(N1*N2*N3.GT.ME) THEN
        PAUSE 'ERROR: TOO LARGE GRID'
        STOP
      END IF
      READ(LU3,*) (X(I1,K1),I1=1,N1)
      READ(LU3,*) (X(I2,K2),I2=1,N2)
      READ(LU3,*) (X(I3,K3),I3=1,N3)
C
C     EVALUATING GRID VALUES OF THE GIVEN FUNCTION:
      IE=0
      DO 23 I3=1,N3
        COOR(K3)=X(I3,K3)
        DO 22 I2=1,N2
          COOR(K2)=X(I2,K2)
          DO 21 I1=1,N1
C
C           EVALUATING THE FUNCTIONAL VALUE:
            COOR(K1)=X(I1,K1)
            IF(TEXT(1:1).EQ.'I'.OR.TEXT(1:1).EQ.'S') THEN
              CALL SRFC2(IGROUP,COOR,UP)
              IVAL=1
              W(I1)=UP(1)
            ELSE
              CALL PARM2(IGROUP,COOR,UP,US,AUX0,AUX1,AUX2)
              IF(TEXT.EQ.'VP ') THEN
                IVAL=1
                W(I1)=UP(1)
              ELSE IF(TEXT.EQ.'VS ') THEN
                IVAL=2
                W(I1)=US(1)
              ELSE IF(TEXT.EQ.'DEN') THEN
                IVAL=3
                W(I1)=AUX0
              ELSE IF(TEXT.EQ.'QP ') THEN
                IVAL=4
                W(I1)=AUX1
              ELSE IF(TEXT.EQ.'QS ') THEN
                IVAL=5
                W(I1)=AUX2
              ELSE
                PAUSE 'ERROR: NAME OF MEDIUM PARAMETER NOT RECOGNIZED'
                STOP
              END IF
            END IF
C
C           EVALUATING THE STANDARD DEVIATION:
            IF(NM.GT.0) THEN
              II=0
   11         CONTINUE
                II=II+1
                CALL VAR6(IVAL,II,NFUN,IFUN(II),B0I,B1I,B2I,B3I)
                IF(II.LE.NFUN) THEN
                  IF(NFUN.GT.MFUN) THEN
                    PAUSE 'ERROR: ARRAY INDEX OUT OF RANGE'
                    STOP
                  END IF
                  FUN(II)=B0I
                END IF
              IF(II.LT.NFUN) GO TO 11
              DO 14 JJ=1,NFUN
                DO 12 II=1,NM
                  IF(INDM(II).EQ.IFUN(JJ)) THEN
                    IFUN(JJ)=II
                    GO TO 13
                  END IF
   12           CONTINUE
                PAUSE 'ERROR: MODEL PARAMETER INDEX NOT RECOGNISED'
                STOP
   13           CONTINUE
   14         CONTINUE
              AUX0=0.
              DO 17 JJ=1,NFUN
                IJJ=IFUN(JJ)*(IFUN(JJ)-1)/2
                AUX2=0.
                DO 16 II=1,JJ-1
                  AUX2=AUX2+     CM(IJJ+IFUN(II))*FUN(II)
   16           CONTINUE
                AUX0=AUX0+FUN(JJ)*(CM(IJJ+IFUN(JJ))*FUN(JJ)+2.*AUX2)
   17         CONTINUE
              IE=IE+1
              E(IE)=SQRT(AUX0)
            END IF
C
   21     CONTINUE
          WRITE(LU3,'(8F10.6)') (W(I1),I1=1,N1)
   22   CONTINUE
        IF(N1.NE.1.AND.N2.NE.1) WRITE(LU3,*)
   23 CONTINUE
C
      IF(NM.GT.0) THEN
        WRITE(LU3,'(A)') '''STANDARD DEVIATIONS ....................'''
        IE=0
        DO 33 I3=1,N3
          DO 32 I2=1,N2
            WRITE(LU3,'(8F10.6)') (E(I1),I1=IE+1,IE+N1)
            IE=IE+N1
   32     CONTINUE
          IF(N1.NE.1.AND.N2.NE.1) WRITE(LU3,*)
   33   CONTINUE
      END IF
C
      STOP
      END
C
C=======================================================================
C
