C
C Program TCCOMP to read two files with propagator matrices U2 and U3 C written in the form of file C GREENTC, C to compute matrix DU=(U2-U3), and to compute its 'norm' as specified C by input parameter INORM. C C Version: 5.60 C Date: 2002, May 23 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 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 a SEP C parameter file. The parameters, which do not differ from their C defaults, need not be specified in file 'SEP'. C Name of the input file with the propagation matrices: C GREENTC2='string'... Name of the file with matrix U2. C Description of file C GREENTC2 C Default: GREENTC2=' ' C GREENTC='string'... Name of the file with matrix U3. C Description of file C GREENTC3 C Default: GREENTC=' ' C Name of the output file: C TCCOMP='string'... Name of the output formatted file with C the norm of matrix DU specified by INORM. C Default: TCCOMP='tccomp.out' C Specification of the norm to be used: C INORM=integer C INORM=1 ... output quantity ADU is the norm C of the matrix DU divided by square root of 2. C otherwise ... output quantity ADU is calculated as C ADU=SQRT(2.*|DU|/(|U2|+|U3|)), where |DU| stands for C the square of the norm of the matrix DU. C Default: INORM=1 C Output formatted file TCCOMP: C For each frequency following line: C F,ADU,TEXT C where F is the frequency, ADU is the 'norm' specified by INORM, C and text is a character enabling the output file to be included C into PostScript file. Thus the first line terminates with C character M, and other lines terminate with character L. C C----------------------------------------------------------------------- C C Subroutines and external functions required: EXTERNAL ERROR,RSEP1,RSEP3T C ERROR ... File C error.for. C RSEP1,RSEP3T ... File C sep.for. C C----------------------------------------------------------------------- C C Input and output data files: CHARACTER*80 FSEP,FIN2,FIN3,FOUT CHARACTER*2 TEXT INTEGER LU,LU2,LU3,LU4,INORM PARAMETER (LU=5,LU2=2,LU3=3,LU4=4) C Auxiliary storage locations: REAL RV11,RV21,RV31,RV12,RV22,RV32,RV13,RV23,RV33 REAL RW11,RW21,RW31,RW12,RW22,RW32,RW13,RW23,RW33 REAL RX11,RX21,RX31,RX12,RX22,RX32,RX13,RX23,RX33 REAL IV11,IV21,IV31,IV12,IV22,IV32,IV13,IV23,IV33 REAL IW11,IW21,IW31,IW12,IW22,IW32,IW13,IW23,IW33 REAL IX11,IX21,IX31,IX12,IX22,IX32,IX13,IX23,IX33 REAL F,ADU C C....................................................................... C C Main input data: FSEP=' ' WRITE(*,'(A)') '+TCCOMP: Enter input filename: ' READ(*,*) FSEP IF (FSEP.EQ.' ') THEN C TCCOMP-01 CALL ERROR('TCCOMP-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 WRITE(*,'(A)') '+TCCOMP: Working... ' C C Reading all the data from file FSEP to the memory C (SEP parameter file form): CALL RSEP1(LU,FSEP) C C Recalling the data: C (arguments: Name of value in input data, Variable, Default): CALL RSEP3T('GREENTC2',FIN2,' ') CALL RSEP3T('GREENTC' ,FIN3,' ') CALL RSEP3T('TCCOMP',FOUT,'tccomp.out') CALL RSEP3I('INORM',INORM,1) OPEN(LU2,FILE=FIN2,STATUS='OLD') OPEN(LU3,FILE=FIN3,STATUS='OLD') OPEN(LU4,FILE=FOUT) TEXT=' M' 10 CONTINUE READ(LU2,*,END=100) F,RV11,IV11,RV21,IV21,RV31,IV31,RV12,IV12, * RV22,IV22,RV32,IV32,RV13,IV13,RV23,IV23,RV33,IV33 READ(LU3,*,END=100) F,RW11,IW11,RW21,IW21,RW31,IW31,RW12,IW12, * RW22,IW22,RW32,IW32,RW13,IW13,RW23,IW23,RW33,IW33 RX11=RV11-RW11 IX11=IV11-IW11 RX21=RV21-RW21 IX21=IV21-IW21 RX31=RV31-RW31 IX31=IV31-IW31 RX12=RV12-RW12 IX12=IV12-IW12 RX22=RV22-RW22 IX22=IV22-IW22 RX32=RV32-RW32 IX32=IV32-IW32 RX13=RV13-RW13 IX13=IV13-IW13 RX23=RV23-RW23 IX23=IV23-IW23 RX33=RV33-RW33 IX33=IV33-IW33 ADU=RX11**2+IX11**2+RX21**2+IX21**2+RX31**2+IX31**2+ * RX12**2+IX12**2+RX22**2+IX22**2+RX32**2+IX32**2+ * RX13**2+IX13**2+RX23**2+IX23**2+RX33**2+IX33**2 IF (INORM.EQ.1) THEN ADU=SQRT(0.5*ADU) ELSE ADV=RV11**2+IV11**2+RV21**2+IV21**2+RV31**2+IV31**2+ * RV12**2+IV12**2+RV22**2+IV22**2+RV32**2+IV32**2+ * RV13**2+IV13**2+RV23**2+IV23**2+RV33**2+IV33**2 ADW=RW11**2+IW11**2+RW21**2+IW21**2+RW31**2+IW31**2+ * RW12**2+IW12**2+RW22**2+IW22**2+RW32**2+IW32**2+ * RW13**2+IW13**2+RW23**2+IW23**2+RW33**2+IW33**2 ADU=SQRT(2.*ADU/(ADV+ADW)) ENDIF WRITE(LU4,'(1(F9.6,1X),F9.6,A)') * F,ADU,TEXT TEXT=' L' GOTO 10 100 CONTINUE CLOSE(LU2) CLOSE(LU3) CLOSE(LU4) WRITE(*,'(A)') '+TCCOMP: Done. ' STOP END C C======================================================================= C INCLUDE 'error.for' C error.for INCLUDE 'sep.for' C sep.for INCLUDE 'length.for' C length.for C======================================================================= C