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','INV3IN','INV3OUT',/ 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 'INV3IN'... Input file specifying the grid in which a selected C function is discretized. C 'INV3OUT'... Copy of the specification of the grid in which the C selected function is discretized from 'INV3IN', C supplemented with the function values at gridpoints. C /... Obligatory slash to enable future compatible extensions. C Default: 'MODEL'='model.dat', 'ANSWER'=' ', 'INV3IN'='inv3.dat', C 'INV3OUT'='inv3.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': 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 Input file 'INV3IN': 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 C Output file 'INV3OUT': C (1)-(7): Copy of (1)-(7) from input file 'INV3IN'. 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: 1996, September 30 C Coded by Ludek Klimes C C----------------------------------------------------------------------- C C Common block /VALC/: INCLUDE 'val.inc' C None of the storage locations of the common block are altered. C C Subroutines and external functions required: EXTERNAL LENGTH INTEGER LENGTH C C----------------------------------------------------------------------- C C Filenames: CHARACTER*80 FILE1,FILE2,FILE3,FILE4 C C Logical unit numbers: INTEGER LU1,LU2,LU3 PARAMETER (LU1=11) PARAMETER (LU2=12) PARAMETER (LU3=13) C C A line of the copied file: CHARACTER*160 LINE 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='model.dat' FILE2=' ' FILE3='inv3.dat' FILE4='inv3.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 1 I1=1,NM READ(LU2,*) INDM(I1),AUX0 RPAR(INDM(I1))=RPAR(INDM(I1))+AUX0 1 CONTINUE IF(NM.GT.0) THEN READ(LU2,*) (CM(I1),I1=1,NM*(NM+1)/2) END IF CLOSE(LU2) END IF C C Copying input file 'INV3IN' to output file 'INV3OUT': OPEN(LU2,FILE=FILE3,STATUS='OLD') OPEN(LU3,FILE=FILE4) 9 CONTINUE READ(LU2,'(A)',END=10) LINE WRITE(LU3,'(A)') LINE(1:LENGTH(LINE)) GO TO 9 10 CONTINUE CLOSE(LU2) REWIND(LU3) C CLOSE(LU3) C C Reading function and grid specifications: C OPEN(LU3,FILE=FILE4,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 INCLUDE 'modelv.for' INCLUDE 'metric.for' INCLUDE 'srfc.for' INCLUDE 'parmv.for' INCLUDE 'valv.for' INCLUDE 'fitv.for' INCLUDE 'var.for' INCLUDE 'length.for' C C======================================================================= C