C
C Program GREEN converting the unformatted output of program CRT into a C formatted file containing the ray-theory elastodynamic Green function. C Program GREEN can consider coupling ray-theory in weakly anisotropic C models. For S waves, second-order perturbations of anisotropic travel C times may be optionally calculated. C C Version: 5.60 C Date: 2002, May 20 C C Coded by: Petr Bulant & Ludek Klimes C Department of Geophysics, Charles University Prague, C Ke Karlovu 3, 121 16 Praha 2, Czech Republic, C E-mails: bulant@seis.karlov.mff.cuni.cz C klimes@seis.karlov.mff.cuni.cz C C....................................................................... C C C Description of data files: C C Main input data read from external interactive device (*): C The data consist of character strings read by list directed (free C format) input. The strings have thus to be enclosed in C apostrophes. The interactive * external unit may be redirected to C the file containing the data. C (1) 'SEP'/ C 'SEP'...Name of the file with input parameters. C Description of file SEP C If blank, default values of the corresponding data are C considered. C Default: 'SEP'=' '. C C C Input file 'SEP' in the SEP format: C The file has the form of a SEP parameter file. C For the description of the SEP format refer to file C 'sep.for'. C Names of optional input files (related to ray tracing): C CRTOUT='string'...File with the names of the output files of C program CRT. C If blank, default names are considered. C Description of file CRTOUT C Default: CRTOUT=' ' (usually sufficient for the first C elementary wave) C REC='string'... If non-blank, the name of the file with the names C and coordinates of the receiver points. The names are C used within the strings describing the points of two-point C rays and the coordinates for the linear interpolation of C travel times from ray endpoints to the receivers. C If blank, the two-point rays are denoted by the receiver C index and the interpolation is not performed. C Description of file REC C Default: REC=' ' C SRC='string'... If non-blank, the name of the file with the name C and coordinates of the source point. The name is used C within the strings describing the two-point rays and the C coordinates for the linear interpolation of travel times C from ray initial points to the source. C If blank, the interpolation is not performed. C Description of file SRC C Default: SRC=' ' (mostly sufficient) C MODEL='string'... Name of the input formatted file with the input C data for the anisotropic model. If blank, the model is C assumed isotropic and the data for the model are not read. C Description of MODEL C Default: MODEL=' ' (sufficient in isotropic models) C The name of the output file: C GREEN='string'... Name of the output formatted file with the Green C tensor. C Description of file GREEN C Default: GREEN='green.out' (mostly sufficient) C Data describing the frequency domain common with program 'ss.for': C These data are used only if MODEL is not blank and the model is C anisotropic, i.e., when the green function is frequency C dependent. These data common with 'ss.for' are used to C generate the optimum default values of input parameters C DF, OF and NF compatible with data for 'ss.for'. C DT=real... Time step to digitize the source time function and C seismograms. C Default: DT=1. C NFFT=integer... Number of the time samples for the fast Fourier C transform. Will be used to convert the time step to the C frequency step. C NFFT must be an integer power of 2. C Default: NFFT=1 C FMIN=real... The lowest frequency of the cosine band-pass filter. C Default: FMIN=0 (mostly sufficient) C FMAX=real... The highest frequency of the cosine band-pass filter. C Default: FMAX=0.5/DT (Nyquist frequency, mostly C sufficient) C Data describing the frequency domain specific to this program: C These data are used only if MODEL is not blank and the model is C anisotropic. C DF=real... Frequency step. C Default: DF=1/(NFFT*DT) (recommended) C OF=real... The elementary Green functions are calculated for NF C frequencies OF,OF+DF,OF+2*DF,...,OF+(NF-1)*DF. C Default: OF=DF*NINT(FMIN/DF) (i.e. FMIN rounded to the C nearest multiple of the frequency step, recommended) C NF=integer... Number of frequencies. C Default: NF=NINT((FMAX-OF)/DF)+1 (recommended) C Hint: Do not specify DF, OF and NF. Use the data for program C 'ss.for' instead and leave the default values of DF, OF and NF C if you are going to generate synthetic seismograms by 'ss.for'. C Specify explicitly DF, OF and NF only if you generate the Green C function for another purpose than the generation of synthetic C seismograms by programs 'greenss.for' and 'ss.for'. C Auxiliary parameters for anisotropic models, need not be specified: C ERRWAN=real... Maximum error in propagator matrix along whole ray. C Default: ERRWAN=0.001 C ERRQI=real... Maximum error of the quasi-isotropic projection of C the polarization vectors. If ERRQI is exceeded, the C warning message is generated. C This parameter has no meaning for QIPV=0. C Default: ERRQI=ERRWAN C QIPV=integer... Optional quasi-isotropic projection of C the polarization vectors. C QIPV=0: No modification of the polarization vectors of the C coupling ray theory. Eigenvectors of the Christoffel C matrix are used as the polarization vectors in C anisotropic models (default). C QIPV=1: Quasi-isotropic projection of the polarization C vectors onto the basis vectors of the ray-centred C coordinate system. This option corresponds to versions C 5.20 to 5.50 of this program, and may also be useful for C comparison of the results with package ANRAY. C Both options should yield equivalent results if QICHM=1. C Default: QIPV=0 C QICHM=integer... Optional quasi-isotropic modification of the C Christoffel matrix. C QICHM=0: No modification of the Christoffel matrix C (default). C QICHM=1: Quasi-isotropic modification of the Christoffel C matrix. Mixed components Gamma_12, Gamma_21, Gamma_13 C and Gamma_31 of the Christoffel matrix in the C ray-centred coordinates are put to zero, in order to C keep the reference isotropic polarization vectors C unperturbed. This modification is designed to compare C the results with package ANRAY. C Default: QICHM=0 C QITT=integer... Optional quasi-isotropic approximation of the C average anisotropic travel time. C QITT=0: No modification of the average anisotropic travel C time calculated along the reference ray. Anisotropic C travel times are linearized with respect to the C Christoffel matrix powered to -1/2, which represents C the slowness (default). C QITT=1: Quasi-isotropic modification of the average C anisotropic travel time calculated along the reference C ray. The anisotropic travel times are linearized with C respect to the Christoffel matrix, which represents C the squared velocity. This modification is designed to C compare the results with package ANRAY, for the default C value of parameter INEWB of program 'weakan.for'. C QITT=1 corresponds to INEWB=0, QITT=0 to INEWB=1. C Default: QITT=0 C QIDT=integer... Optional quasi-isotropic approximation of the C difference between the anisotropic travel times. C QIDT=0: No modification of the difference between the C anisotropic travel times (default). C QIDT=1: Quasi-isotropic modification of the difference C between the anisotropic travel times. The anisotropic C travel times are linearized with respect to the C Christoffel matrix instead of its inverse square root. C This modifications is designed to compare the results C with package ANRAY, for the default value of parameter C INEWB of program 'weakan.for'. C Default: QIDT=0 C Parameters for optional second-order perturbations of anisotropic C travel times of S waves. The perturbations may be calculated only C along the rays which do not cross any interface. The perturbations C may not be calculated for P waves. Furthermore, in more complex C models in the vicinity of caustics, the second-order perturbations C may be infinitely large and should be avoided. C QIRAY=integer ... Optional second-order perturbations C of anisotropic travel times of S waves. C QIRAY=0: No second-order perturbations (default). C QIRAY=1: Second-order perturbations of anisotropic travel C times of S waves are calculated and are taken into C account when calculating the output Green function. C Default: QIRAY=0 C QILST='string' ... Name of the optional output file containing C second-order perturbations of anisotropic travel times C of S waves. C Specification of QILST does not influence the output C file with the Green function. C Description of file QILST. C Default: QILST=' ' (the file is not generated) C C Example of file 'SEP' C C C Input formatted file CRTOUT: C For general description of file CRTOUT refer to the file C writ.for. C Description specific to this program: C 'CRT-R'.. The file is not read if MODEL=' ' or if the model C is isotropic. If the file is read, only two-point rays C stored in both files 'CRT-R' and 'CRT-S' are selected. C 'CRT-S'.. Only two-point rays (rays incident to the vicinities of C the receivers situated at the reference surface) are read C from file 'CRT-S'. The output elementary Green functions C are then evaluated at the receivers. C 'CRT-I'.. If file 'CRT-R' is not read, file 'CRT-I' must contain C at least initial points of all rays contained in file C 'CRT-S'. If file 'CRT-R' is read, file 'CRT-I' must C contain at least the initial points of all rays common C to files 'CRT-R' and 'CRT-S'. C 'CRT-T'.. Not used. C Only the first data line is read in the current version. C Defaults for the first data line: C 'CRT-R'='r01.out', 'CRT-S'='s01.out', 'CRT-I'='s01i.out' C C C Output formatted file GREEN: C (1) / (a slash). C (2) For each two-point ray (2.1): C (2.1) 'R','S',(GREEN(I),I=1,NRAYPT),/ C Here NRAYPT=32 except for S waves in anisotropic model, C where NRAYPT=14+18*NF C 'R'... String in apostrophes describing the receiver. C 'S'... String in apostrophes describing the source. C GREEN(1)... Travel time between the receiver and the source. C GREEN(2)... Imaginary part of the complex-valued travel time C between receiver and source due to attenuation. C GREEN(3:8)... Coordinates of the receiver and coordinates of the C source. C GREEN(9:14)... Derivatives of the travel time with respect to the C coordinates of the receiver and coordinates of the source. C GREEN(15:32)... 1000000 times enlarged amplitude of the Green C function: contravariant components of the complex-valued C 3*3 matrix Gij in model coordinates, where the first C subscript corresponds to the receiver and the second C subscript corresponds to the source. The components are C ordered as C Re(G11),Im(G11),Re(G21),Im(G21),Re(G31),Im(G31), C Re(G12),Im(G12),Re(G22),Im(G22),Re(G32),Im(G32), C Re(G13),Im(G13),Re(G23),Im(G23),Re(G33),Im(G33). C GREEN(33:*)... Analogy of GREEN(15:32) for other frequencies. C Written just for frequency-dependent elementary Green C functions. C /... An obligatory slash at the end of line. It may also stand C for zeros at the end of line. C (3) / (a slash). C File form GREEN is, to some extent, an extension of file form C Travel Times. C C C Output formatted file QILST with S-wave travel-time perturbations: C (1) / (a slash). C (2) For each two-point ray (2.1): C (2.1) IREC,T1,T2,DTL1,DTL2,DTQ1,DTQ2,DTQC,/ C IREC... Number of the receiver corresponding to the two-point ray. C T1,T2 ... Anisotropic travel times of QS1 and QS2 waves. C DTL1,DTL2 ... First-order perturbations of travel time, C from the isotropic S-wave reference ray to the C anisotropic ray theory S1 and S2 rays. C DTQ1,DTQ2 ... Second-order perturbations of travel time, C from the isotropic S-wave reference ray to the C anisotropic ray theory S1 and S2 rays. C DTQC ... Estimated second-order perturbations of travel time, C from the anisotropic S-wave common reference ray to the C anisotropic ray theory S1 and S2 rays. C /... A slash at the end of line. C (3) / (a slash). C C----------------------------------------------------------------------- C C Subroutines and external functions required: EXTERNAL KOOR,NSRFC,AP00,AP21,TXT1,TXT2,REC,SRC,TTSORT, * FORM2,LENGTH,WAN INTEGER KOOR,NSRFC,LENGTH C AP00,AP21... File 'ap.for'. C AP03,AP03A... File 'ap.for' (called by AP21,AP03). C TXT1,TXT2,REC,SRC... File 'crtout.for'. C FORM2...File 'forms.for'. C FORM1...File 'forms.for' (called by FORM2). C LENGTH..File 'length.for'. C TTSORT..File 'ttsort.for'. C INDEXX..File 'indexx.for' (called by TTSORT). C WAN ... File 'wan.for'. C WAPTS,WAREAD,WACHRI,WAMAT,WACHAN,WAPROJ,WAPERT,WASUM,WASUM4,WASUM5 C ... File 'wan.for' (called by WAN). C C----------------------------------------------------------------------- C C Common block /RAMC/: INCLUDE 'ram.inc' C ram.inc C C Common block /POINTC/ to store the results of complete ray tracing: INCLUDE 'pointc.inc' C pointc.inc C None of the storage locations of the common block are altered. C C----------------------------------------------------------------------- C C Allocation of working arrays: INTEGER IRAM(MRAM) EQUIVALENCE (IRAM,RAM) C C Auxiliary storage locations: INTEGER LU0,LU1,LU2,LU3,LU4,LUQI PARAMETER (LU0=1,LU1=2,LU2=3,LU3=4,LU4=5) CHARACTER*80 FILSEP,FILMOD,FILCRT,FILREC,FILSRC,FGREEN,FQILST CHARACTER*80 FILE1,FILE2,FILE3,FILO1,FILO2,FILO3 CHARACTER*260 FORMAT CHARACTER*80 RAYTXT,TEXT INTEGER MQ,MTPR,NTPR,IW,IGREEN,I,K INTEGER NFFT,NF,NQ REAL DT,FMIN,FMAX,DF,OF,ERRWAN,ERRQI,ERR(2) C MQ... Maximum number of quantities stored for a two-point ray: C IRAM(MQ*NTPR+1)... Number of quantities stored for the C (NTPR+1)-th two-point ray. C IRAM(MQ*NTPR+2)... Index of the receiver. C RAM(MQ*NTPR+3) to RAM(MQ*NTPR+16) ... Travel time, C coordinates, gradient of travel time. C RAM(MQ*NTPR+17) to RAM(MQ*NTPR+IRAM(MQ*NTPR+1))... C Amplitudes. C MTPR... Maximum number of two-point rays. C Note: MTPR=1 indicates no sorting according to receivers. C NTPR... Number of two-point rays. C IW... Loop variable: index of the elementary wave. C C----------------------------------------------------------------------- C C Reading name of SEP file with input data: WRITE(*,'(A)') '+GREEN: Enter input filename: ' FILSEP=' ' READ(*,*) FILSEP C C Reading all data from the SEP file into the memory: IF (FILSEP.NE.' ') THEN CALL RSEP1(LU0,FILSEP) ELSE C GREEN-06 CALL ERROR('GREEN-06: 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)') '+GREEN: Working ... ' C C Name of the model file CALL RSEP3T('MODEL',FILMOD,' ') C Name of the receiver file CALL RSEP3T('REC',FILREC,' ') C Name of the source file CALL RSEP3T('SRC',FILSRC,' ') C Name of the output file CALL RSEP3T('GREEN',FGREEN,'green.out') C C Reading the data for a possibly anisotropic model: IF(FILMOD.NE.' ') THEN OPEN(LU0,FILE=FILMOD,STATUS='OLD') CALL MODEL1(LU0) CLOSE(LU0) C Checking whether the model is anisotropic CALL PARM4(I) IF(I.GT.0) THEN C Isotropic model FILMOD=' ' END IF END IF IF(FILMOD.NE.' ') THEN IF(KOOR().NE.0) THEN C GREEN-07 CALL ERROR * ('GREEN-07: Anisotropic model in curvilinear coordinates') C Current version of subroutine file 'wan.for' for weakly C anisotropic models is coded only for Cartesian coordinates. END IF IF(NSRFC().NE.0) THEN C GREEN-51 CALL WARN('GREEN-51: Anisotropic model with interfaces') C This version of program GREEN is not debugged for anisotropic C models with structural interfaces. END IF END IF C C Number of frequencies, maximum number of values per two-point ray: IF(FILMOD.EQ.' ') THEN NF=1 MQ=34 ELSE CALL RSEP3I('NFFT',NFFT,1) CALL RSEP3R('DT ',DT ,1.) CALL RSEP3R('FMIN',FMIN,0.) CALL RSEP3R('FMAX',FMAX,0.5/DT) CALL RSEP3R('DF ',DF ,1./(FLOAT(NFFT)*DT)) CALL RSEP3R('OF ',OF ,DF*NINT(FMIN/DF)) CALL RSEP3I('NF ',NF ,NINT((FMAX-OF)/DF)+1) MQ=16+18*NF CALL RSEP3R('ERRWAN',ERRWAN,0.001) CALL RSEP3R('ERRQI' ,ERRQI ,ERRWAN) ERR(1)=0. ERR(2)=0. CALL RSEP3T('QILST' ,FQILST ,' ') IF (FQILST.EQ.' ') THEN LUQI=0 ELSE LUQI=6 OPEN(LUQI,FILE=FQILST) WRITE(LUQI,'(A)') '/' ENDIF END IF C C Reading the source and receiver names and preparing them for TXT2: CALL TXT1(LU0,FILSRC,FILREC) C C Initializing the output file: OPEN(LU4,FILE=FGREEN) WRITE(LU4,'(A)') '/' C C Maximum number of two-point rays: MTPR=MRAM/(MQ+2) C Note: MQ is the maximum number of quantities per a two-point ray C (including two integers), and the space for 2 working arrays C (required for sorting) is reserved at the end of array RAM. C Note: If MTPR=1, the rays are not sorted according to receivers. C C Number of stored two-point rays: NTPR=0 C C....................................................................... C C Loop for elementary waves: FILO1=' ' FILO2=' ' FILO3=' ' DO 70 IW=1,1 C C....................................................................... C C Reading output filenames of the CRT program: CALL RSEP3T('CRTOUT',FILCRT,' ') IF(IW.EQ.1) THEN FILE1='r01.out' FILE2='s01.out' FILE3='s01i.out' ELSE FILE1=' ' FILE2=' ' FILE3=' ' END IF IF(FILCRT.NE.' ') THEN IF(IW.EQ.1) THEN OPEN(LU0,FILE=FILCRT,STATUS='OLD') END IF READ(LU0,*,END=90) FILE1,FILE2,FILE3 END IF IF(FILE1.EQ.' '.AND.FILE2.EQ.' '.AND.FILE3.EQ.' ') THEN C No new elementary wave GO TO 90 END IF C C Opening output files of the CRT program for input: IF(FILMOD.NE.' ') THEN IF(FILE1.EQ.' ') THEN C GREEN-01 CALL ERROR('GREEN-01: No input file with rays') C Calculation of the coupling ray theory Green function C in the anisotropic model, specified by input parameter MODEL, C requires the two-point rays stored by program 'crt.for' into C the file which name is taken from the file given by parameter C CRTOUT. END IF IF(FILE1.NE.FILO1) THEN IF(IW.GT.1) THEN CLOSE(LU1) END IF OPEN(LU1,FILE=FILE1,FORM='UNFORMATTED',STATUS='OLD',ERR=2) END IF END IF IF(FILE2.NE.FILO2) THEN IF(IW.GT.1) THEN CLOSE(LU2) END IF OPEN(LU2,FILE=FILE2,FORM='UNFORMATTED',STATUS='OLD',ERR=3) END IF IF(FILE3.NE.FILO3) THEN IF(IW.GT.1) THEN CLOSE(LU3) END IF OPEN(LU3,FILE=FILE3,FORM='UNFORMATTED',STATUS='OLD',ERR=4) END IF GO TO 9 C Error messages: 2 CONTINUE C GREEN-02 CALL ERROR('GREEN-02: Cannot open the input file with rays') C The default name of the file is 'r01.out'. 3 CONTINUE C GREEN-03 CALL ERROR('GREEN-03: Cannot open the input file with points') C The default name of the file is 's01.out'. 4 CONTINUE C GREEN-04 CALL ERROR * ('GREEN-04: Cannot open the input file with initial points') C The default name of the file is 's01i.out'. 9 CONTINUE C C....................................................................... C C Loop for the points of rays: 10 CONTINUE C Checking free memory for the next two-point ray: IF (NTPR.GE.MTPR) THEN C GREEN-05 CALL ERROR('GREEN-05: Too many Elementary Green functions') C Elementary Green functions require 36 storage locations for C each two-point ray in an isotropic model and 18+18*NF storage C locations for each two-point ray in an anisotropic model. C The Green functions are stored in array RAM declared in C include file ram.inc. END IF C C Reading the results of the complete ray tracing IGREEN=MQ*NTPR IF(FILMOD.EQ.' ') THEN C Isotropic model 20 CONTINUE CALL AP00(0,LU2,LU3) IF(IWAVE.LT.1) THEN C End of rays GO TO 60 END IF IF (IREC.LE.0) GO TO 20 C Determination of the Green function CALL AP21(1,RAM(1),RAM(1),RAM(IGREEN+3)) NQ=32 ELSE C Anisotropic model CALL WAN(LU1,LU2,LU3,RAM(IGREEN+3),MQ-2,RAM(IGREEN+MQ+1), * IRAM(IGREEN+MQ+1),MRAM-IGREEN-MQ,NQ,ERR,LUQI) IF(NQ.LE.0) THEN C End of rays GO TO 60 END IF END IF C IRAM(IGREEN+1)=NQ+2 IRAM(IGREEN+2)=IREC C Linear Taylor expansion of travel time CALL RECSRC(FILREC,FILSRC,RAM(IGREEN+3)) C Rescaling the amplitudes of the Green function DO 30 I=IGREEN+17,IGREEN+2+NQ RAM(I)=1000000.*RAM(I) 30 CONTINUE C C Writing unsorted two-point rays if MTPR=1: NTPR=NTPR+1 IF(MTPR.EQ.1) THEN CALL WGREEN(LU4,IRAM(1)-2,IRAM(2),RAM(3)) NTPR=0 END IF GO TO 10 60 CONTINUE C C....................................................................... C FILO1=FILE1 FILO2=FILE2 FILO3=FILE3 70 CONTINUE C End of loop for elementary waves. C C....................................................................... C C Sorting and writing two-point rays: 90 CONTINUE IF(MTPR.GT.1) THEN C Collecting the indices of the receivers for sorting: DO 91 K=1,NTPR IRAM(MQ*NTPR+K)=IRAM(MQ*(K-1)+2) 91 CONTINUE C C Sorting elementary Green functions according to the receivers: CALL TTSORT(MQ,NTPR,3,RAM,IRAM(MQ*NTPR+1), * RAM(MQ*NTPR+1),IRAM((MQ+1)*NTPR+1)) C C Writing elementary Green functions: DO 92 K=(MQ+1)*NTPR+1,(MQ+2)*NTPR IGREEN=MQ*(IRAM(K)-1) CALL WGREEN(LU4,IRAM(IGREEN+1)-2,IRAM(IGREEN+2),RAM(IGREEN+3)) 92 CONTINUE END IF C C....................................................................... C WRITE(LU4,'(A)') '/' CLOSE(LU4) CLOSE(LU3) CLOSE(LU2) IF(FILMOD.NE.' ') THEN CLOSE(LU1) IF (FQILST.NE.' ') THEN WRITE(LUQI,'(A)') '/' CLOSE(LUQI) ENDIF END IF IF(FILCRT.NE.' ') THEN CLOSE(LU0) END IF IF(FILMOD.NE.' ') THEN IF(ERR(1).GT.ERRQI.OR.ERR(2).GT.ERRQI) THEN C GREEN-52 WRITE(TEXT,'(2(A,F8.6))') * 'GREEN-52: QI projection error: source ',ERR(1), * ', receiver ',ERR(2) CALL WARN(TEXT(1:LENGTH(TEXT))) END IF END IF WRITE(*,'(A)') '+GREEN: Done. ' STOP END C C======================================================================= C SUBROUTINE RECSRC(FILREC,FILSRC,GREEN) CHARACTER*(*) FILREC,FILSRC REAL GREEN(14) C C----------------------------------------------------------------------- C C Common block /POINTC/ to store the results of complete ray tracing: INCLUDE 'pointc.inc' C pointc.inc C None of the storage locations of the common block are altered. C C----------------------------------------------------------------------- C REAL COOR(3) C C....................................................................... C C Linear Taylor expansion of travel time: C Receiver: IF(FILREC.NE.' ') THEN C Receiver: COOR(1)=GREEN(3) COOR(2)=GREEN(4) COOR(3)=GREEN(5) CALL REC(IREC,COOR(1),COOR(2),COOR(3)) GREEN(1)=GREEN(1)+(COOR(1)-GREEN(3))*GREEN( 9) * +(COOR(2)-GREEN(4))*GREEN(10) * +(COOR(3)-GREEN(5))*GREEN(11) GREEN(3)=COOR(1) GREEN(4)=COOR(2) GREEN(5)=COOR(3) END IF IF(FILSRC.NE.' ') THEN C Source: COOR(1)=GREEN(6) COOR(2)=GREEN(7) COOR(3)=GREEN(8) CALL SRC(1,COOR(1),COOR(2),COOR(3)) GREEN(1)=GREEN(1)+(COOR(1)-GREEN(6))*GREEN(12) * +(COOR(2)-GREEN(7))*GREEN(13) * +(COOR(3)-GREEN(8))*GREEN(14) GREEN(6)=COOR(1) GREEN(7)=COOR(2) GREEN(8)=COOR(3) END IF RETURN END C C======================================================================= C SUBROUTINE WGREEN(LU4,NGREEN,IREC,GREEN) INTEGER LU4,NGREEN,IREC REAL GREEN(NGREEN) C C----------------------------------------------------------------------- C CHARACTER*260 FORMAT CHARACTER*80 FILREC,RAYTXT INTEGER LENTXT,I,L C C....................................................................... C C Text strings: CALL TXT2(0,1,1,1,IREC,LENTXT,RAYTXT) L=INDEX(RAYTXT(1:LENTXT),'''') L=INDEX(RAYTXT(L+1:LENTXT),'''')+L CALL RSEP3T('REC',FILREC,' ') IF(FILREC.EQ.' ') THEN C Shortening the receiver string part from 8 to 6 characters LENTXT=LENTXT-2 DO 71 I=L+4,LENTXT RAYTXT(I:I)=RAYTXT(I+2:I+2) 71 CONTINUE END IF C C Writing: FORMAT(1:4)='(4A,' IF(NGREEN.LE.32) THEN CALL FORM2(32,GREEN,GREEN,FORMAT(5:260)) WRITE(LU4,FORMAT) RAYTXT(L+2:LENTXT),' ',RAYTXT(1:L), * (' ',GREEN(I),I=1,NGREEN),' /' ELSE CALL FORM2(14,GREEN,GREEN,FORMAT(5:260)) WRITE(LU4,FORMAT) RAYTXT(L+2:LENTXT),' ',RAYTXT(1:L), * (' ',GREEN(I),I=1,14) L=15 DO 80 L=L,NGREEN-18,18 FORMAT(1:4)='(1A,' CALL FORM2(18,GREEN(L),GREEN(L),FORMAT(5:260)) WRITE(LU4,FORMAT) (' ',GREEN(I),I=L,L+17) 80 CONTINUE FORMAT(1:4)='(1A,' CALL FORM2(NGREEN-L+1,GREEN(L),GREEN(L),FORMAT(5:260)) WRITE(LU4,FORMAT) (' ',GREEN(I),I=L,NGREEN),' /' END IF C RETURN 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 INCLUDE 'indexx.for' C indexx.for of Numerical Recipes INCLUDE 'ap.for' C ap.for INCLUDE 'crtout.for' C crtout.for INCLUDE 'ttsort.for' C ttsort.for INCLUDE 'wan.for' C wan.for INCLUDE 'model.for' C model.for INCLUDE 'metric.for' C metric.for INCLUDE 'srfc.for' C srfc.for INCLUDE 'parm.for' C parm.for INCLUDE 'val.for' C val.for INCLUDE 'fit.for' C fit.for C C======================================================================= C