C SUBROUTINE FILE 'CODE.FOR' - CODES FOR ELEMENTARY WAVES. C C BY VLASTISLAV CERVENY, LUDEK KLIMES, IVAN PSENCIK C C THIS FILE CONSISTS OF: C GENERAL DESCRIPTION OF THE CODES FOR ELEMENTARY WAVES. C CODEB...BLOCK DATA SUBROUTINE DEFINING COMMON BLOCK /COD/ TO STORE C THE CODE OF THE COMPUTED ELEMENTARY WAVE. C CODE1...SUBROUTINE DESIGNED TO READ THE INPUT DATA FOR THE CODES C OF ELEMENTARY WAVES AND TO STORE THEM IN THE COMMON BLOCK C /COD/. C LCODE...INTEGER FUNCTION RETURNING THE LENGTH OF THE CODE OF THE C CURRENT ELEMENTARY WAVE. C NELEM...INTEGER FUNCTION RETURNING THE NUMBER OF RAY ELEMENTS C CORRESPONDING TO THE GIVEN POSITION IN THE CODE OF THE C ELEMENTARY WAVE. C CODE... SUBROUTINE DESIGNED TO TRANSFORM THE USED NUMERICAL CODE C OF THE ELEMENTARY WAVE UNDER CONSIDERATION INTO C INSTRUCTIONS SPECIFYING THE BEHAVIOUR OF THE WAVE AT THE C INITIAL POINTS OF RAYS AND AT ALL POINTS OF INCIDENCE OF C THE RAYS AT INTERFACES (BOUNDARIES OF COMPLEX BLOCKS). C C INPUT DATA FOR THE CODES OF ELEMENTARY WAVES: C THE DATA ARE READ IN BY THE LIST DIRECTED INPUT (FREE FORMAT). IN C THE LIST OF INPUT DATA BELOW, EACH NUMBERED PARAGRAPH INDICATES C THE BEGINNING OF A NEW INPUT OPERATION (NEW READ STATEMENT). IF C THE FIRST LETTER OF THE SYMBOLIC NAME OF THE INPUT VARIABLE IS C I-N, THE CORRESPONDING VALUE IN INPUT DATA MUST BE OF THE TYPE C INTEGER. OTHERWISE (EXCEPT TEXTC), THE INPUT PARAMETER IS OF THE C TYPE REAL. C (1) TEXTC C STRING DESCRIBING THE DATA. ONLY THE FIRST 80 CHARACTERS OF THE C STRING ARE SIGNIFICANT. C (2) NSKIP C THE NUMBER OF ELEMENTARY WAVES TO BE SKIPPED OVER. THE C (NSKIP+1)-TH ELEMENTARY WAVE IN THESE INPUT DATA WILL BE THE FIRST C COMPUTED WAVE. DEFAULT VALUE: 0. C FOR EACH ELEMENTARY WAVE IWAVE THE FOLLOWING DATA (3): C (3) KODTYP,KODE C KODTYP..DETERMINES THE TYPE OF THE CODE C 1... BLOCK-SURFACE CODE, C 2... BASIC BLOCK CODE, C 3... COMPOUND-ELEMENT BLOCK CODE, C 4... GENERALIZED LAYER CODE, C 5... TRANSMISSION-REFLECTION CODE. C KODE... THE CODE OF THE ELEMENTARY WAVE IN THE FORM OF A FINITE C SEQUENCE OF INTEGERS TERMINATED BY A SLASH. C ATTENTION: C GENERALLY, DIFFERENT TYPES OF THE CODE MAY BE COMBINED C WITHIN THE SAME INPUT DATA, BUT THE USER HAS THEN BE VERY C CAREFUL. IF COMBINING THE DIFFERENT TYPES OF THE CODE, C OUTPUT PARAMETERS IWAVE0 AND IKODE OF SUBROUTINE CODE1 MAY C HAVE NO SENSE. THEN: C IF THE PARAMETER MODCRT OF THE DATA SET 'DCRT' (SEE SOURCE C CODE FILE 'RAY.FOR') IS ZERO, THE OUTPUT PARAMETERS C IWAVE0 AND IKODE OF SUBROUTINE CODE1 ARE HOPEFULLY C UNUSED. C IF THE PARAMETER MODCRT OF THE DATA SET 'DCRT' IS NOT C ZERO, MIXING THE TYPES OF THE WAVE CODE SHOULD BE C AVOIDED UNLESS THE USER FINDS A WAY HOW TO DO IT C REASONABLY. C (4) A SLASH. C C STORAGE IN THE MEMORY: C THE INPUT DATA (2) FOR THE COMPUTED ELEMENTARY WAVE ARE STORED IN C THE COMMON BLOCK /COD/ DEFINED IN THE FOLLOWING SUBROUTINE: C ------------------------------------------------------------------ BLOCK DATA CODEB INTEGER MKODE PARAMETER (MKODE=128) INTEGER KODTYP,KODE(MKODE) COMMON /COD/ KODTYP,KODE SAVE /COD/ END C ------------------------------------------------------------------ C KODTYP..DETERMINES THE TYPE OF THE CODE, SEE THE INPUT DATA (2). C KODE... ARRAY CONTAINING THE CODE OF THE ELEMENTARY WAVE. C THE CODE IS A FINITE SEQUENCE OF NONZERO INTEGERS. C ZERO INDICATES THE END OF A CODE. C COMMON BLOCK /COD/ IS INCLUDED IN EXTERNAL PROCEDURES COD1, CODE, C WRIT1, AND MAY BE INCLUDED IN ANY OTHER SUBROUTINE. IF MKODE IS C CHANGED, IT MUST BE ADJUSTED IN ALL SUBROUTINES WHICH INCLUDE C COMMON BLOCK /COD/. C C DATE: 1993, AUGUST 3 C CODED BY LUDEK KLIMES C C======================================================================= C C CODES FOR ELEMENTARY WAVES - GENERAL DESCRIPTION C C WE CONSIDER ORDINARY SEISMIC BODY WAVES, SUCH AS REFRACTED, C PRIMARILY OR MULTIPLY REFLECTED, POSSIBLY CONVERTED WAVES. IN C GENERAL, INCIDENCE OF A WAVE AT AN INTERFACE (BOUNDARY OF A C COMPLEX BLOCK) PRODUCES FOUR WAVES, REFLECTED P AND S, AND C TRANSMITTED P AND S WAVES. WHEN PERFORMING COMPLETE RAY TRACING, C WE MUST KNOW A PRIORI WHICH OF THE FOUR GENERATED WAVES TO FOLLOW. C THIS DECISION MUST BE MADE AT EACH INTERFACE. THE ALPHANUMERIC C STRING SPECIFYING THE BEHAVIOUR OF A RAY FROM ITS INITIAL POINT TO C ITS ENDPOINT IS A 'CODE'. IN THE FOLLOWING, WE SHALL CONSIDER THE C CODE TO BE A SEQUENCE OF NONZERO INTEGERS. C C THE TERM 'ELEMENTARY WAVE' DOES NOT HAVE UNIQUE MEANING IN THE C LITERATURE. HERE WE APPLY THE TERM ELEMENTARY WAVE TO THAT PART C OF THE WAVEFIELD THAT IS DESCRIBED BY ONE SPECIFIC CODE. SINCE C THERE MAY BE VARIOUS TYPES OF CODES, THERE IS ALSO A VARIETY OF C DIVISIONS OF THE WAVE FIELD INTO ELEMENTARY WAVES. C C WE INTRODUCE THE TERM 'ELEMENT OF A RAY', WHICH HAS AN IMPORTANT C MEANING IN THE CONSTRUCTION OF CODES. BY AN ELEMENT OF A RAY, WE C DENOTE THAT PART OF THE RAY THAT IS SITUATED IN ONE COMPLEX BLOCK C BETWEEN TWO SUCCESSIVE POINTS OF REFLECTION/TRANSMISSION, OR C BETWEEN THE INITIAL POINT OR ENDPOINT OF THE RAY AND THE CLOSEST C POINT OF REFLECTION/TRANSMISSION, OR BETWEEN THE INITIAL POINT AND C THE ENDPOINT OF THE RAY, IF THE RAY IS ENTIRELY SITUATED IN ONE C COMPLEX BLOCK. C C IN THE FOLLOWING, WE PRESENT SEVERAL POSSIBLE TYPES OF NUMERICAL C CODES OF ELEMENTARY WAVES. C C EXAMPLES OF CODE TYPES C (1) BLOCK-SURFACE CODE C TWO INTEGERS ARE USED FOR ANY ELEMENT OF THE RAY. THE C ABSOLUTE VALUE OF THE FIRST INTEGER SPECIFIES THE INDEX OF C THE COMPLEX BLOCK, IN WHICH THE ELEMENT OF THE RAY IS C SITUATED, THE SIGN SPECIFIES THE TYPE OF THE WAVE (PLUS C SIGN... P WAVE, MINUS SIGN... S WAVE). THE SECOND C NUMBER IS THE INDEX OF THE SMOOTH SURFACE, AT WHICH THE C ELEMENT TERMINATES. THE CODE OF THE RAY IS A CHAIN OF THE C ABOVE DOUBLETS DESCRIBING THE RAY FROM ITS INITIAL POINT C TO ITS ENDPOINT. THE INDEX OF A SURFACE, AT WHICH THE C LAST ELEMENT OF A RAY TERMINATES, NEED NOT BE SPECIFIED. C IN THIS CASE THE RAY IS ALLOWED TO TERMINATE ON ANY C SURFACE BOUNDING THE COMPLEX BLOCK IN WHICH THE LAST C ELEMENT OF THE RAY IS SITUATED. APPLICATION OF THIS CODE C IS STRAIGHTFORWARD. C C IF ANY COMPLEX BLOCK IS BOUNDED BY SEVERAL SMOOTH C SURFACES, THE BLOCK-SURFACE CODE MAY DIVIDE ARTIFICIALLY C THE WAVE FIELD INTO SEVERAL ELEMENTARY WAVES. THIS WOULD C RESULT IN REPEATED APPLICATION OF THE COMPLETE RAY TRACING C TO ALL OF THE ELEMENTARY WAVES, WHAT MAY BE TIME C CONSUMING. IN SUCH SITUATIONS, IT IS THEREFORE DESIRABLE C TO USE OTHER CODES, FOR WHICH THE ELEMENTARY WAVE HAS A C MORE GENERAL MEANING. C C (2) BASIC BLOCK CODE C THIS CODE MAY BE OBTAINED FROM THE BLOCK-SURFACE CODE IF C THE INDICES OF THE SURFACES, AT WHICH INDIVIDUAL ELEMENTS C OF A RAY TERMINATE, ARE OMITTED. THEN, ANY ELEMENT OF THE C RAY IS DESCRIBED BY A SINGLE INTEGER SPECIFYING THE C COMPLEX BLOCK AND WAVE TYPE, AND IT MAY TERMINATE AT ANY C SURFACE BOUNDING THE COMPLEX BLOCK. THE CODE OF THE RAY C IS AGAIN A CHAIN OF THE NUMBERS DESCRIBING SUCCESIVELY ITS C INDIVIDUAL ELEMENTS. C C ELEMENTARY WAVES SPECIFIED BY THE BLOCK-SURFACE CODE, C PASSING THROUGH THE SAME COMPLEX BLOCKS, AND REFLECTED C FROM DIFFERENT SURFACES BOUNDING A COMPLEX BLOCK, ARE C UNITED INTO ONE ELEMENTARY WAVE IN THE BASIC BLOCK CODE. C SIMILARLY, ELEMENTARY WAVES SPECIFIED BY THE BLOCK-SURFACE C CODE AND TRANSMITTED FROM A BLOCK INTO THE NEIGHBOURING C BLOCK ACROSS DIFFERENT SURFACES SEPARATING THESE BLOCKS, C ARE UNITED INTO ONE ELEMENTARY WAVE BY THE BASIC BLOCK C CODE. C C (3) COMPOUND-ELEMENT BLOCK CODE C LET US INTRODUCE THE TERMS SIMPLE AND COMPOUND ELEMENT OF C A RAY. FOR THIS PURPOSE, WE CALL THE 'LOWER-INDEX C BOUNDARY' OF A COMPLEX BLOCK THAT PART OF ITS BOUNDARY, C WHICH SEPARATES THE COMPLEX BLOCK EITHER FROM COMPLEX C BLOCKS WITH LOWER INDICES OR FROM A FREE SPACE. THE C REMAINING PART OF THE BOUNDARY IS CALLED THE C 'HIGHER-INDEX BOUNDARY' OF THE COMPLEX BLOCK. THE INITIAL C POINT OF A RAY IS TREATED IN THE SAME WAY AS THE POINTS C SITUATED AT A LOWER-INDEX BOUNDARY. C C AN ELEMENT OF A RAY IS CALLED 'SIMPLE ELEMENT' IF ITS C INITIAL POINT (E.G. THE POINT WHERE THE WAVE ENTERS THE C CORRESPONDING COMPLEX BLOCK) AND ITS ENDPOINT (E.G. THE C POINT AT WHICH THE WAVE LEAVES THE COMPLEX BLOCK) ARE C SITUATED ONE ON THE LOWER-INDEX AND THE OTHER ON THE C HIGHER-INDEX BOUNDARY OF THE COMPLEX BLOCK. THE C 'COMPOUND ELEMENT' IS SUCH AN ELEMENT FOR WHICH ITS C INITIAL POINT AND ITS ENDPOINT ARE BOTH SITUATED EITHER C ON THE LOWER-INDEX OR ON THE HIGHER-INDEX BOUNDARY OF THE C COMPLEX BLOCK. THE COMPOUND ELEMENT IS FORMALLY C CONSIDERED AS TWO SIMPLE ELEMENTS. C C OFTEN IT IS CONVENIENT TO WORK WITH WAVES, THE RAYS OF C WHICH HAVE EITHER A COMPOUND ELEMENT IN A COMPLEX BLOCK OR C TWO SIMPLE ELEMENTS OF THE SAME WAVE TYPE (UNCONVERTED C REFLECTION) IN THE SAME BLOCK, AS WITH ONE ELEMENTARY C WAVE. THE COMPOUND-ELEMENT BLOCK CODE MAKES A DIVISION OF C THE WAVE FIELD INTO SUCH ELEMENTARY WAVES POSSIBLE. C C WE INTRODUCE THE COMPOUND-ELEMENT BLOCK CODE AS FOLLOWS. C FOR ANY SIMPLE ELEMENT OF A RAY ONE NUMBER IS USED. THE C ABSOLUTE VALUE OF THIS NUMBER SPECIFIES THE INDEX OF THE C COMPLEX BLOCK, IN WHICH THE SIMPLE ELEMENT IS SITUATED. C THE SIGN OF THIS NUMBER SPECIFIES THE TYPE OF THE WAVE C (PLUS SIGN... P WAVE, MINUS SIGN... S WAVE). FOR ANY C COMPOUND ELEMENT OF A RAY, TWO IDENTICAL NUMBERS ARE USED. C A COMPOUND ELEMENT AND A DOUBLET OF SIMPLE ELEMENTS, C DESCRIBED BY THE SAME NUMBERS, ARE NOT DISTINGUISHED. THE C CODE OF A RAY IS A CHAIN OF THESE NUMBERS DESCRIBING THE C RAY FROM ITS INITIAL POINT TO ITS ENDPOINT. THE USE OF C THE COMPOUND-ELEMENT BLOCK CODE LEADS TO MORE GENERAL C ELEMENTARY WAVES. C C (4) GENERALIZED LAYER CODE C ANOTHER CODE MAY BE OBTAINED FROM THE COMPOUND-ELEMENT C BLOCK CODE IF THE MODIFIED INTERPRETATION OF THE MODEL C STRUCTURE IS USED. THE MODIFIED INTERPRETATION CONSISTS C IN ASSUMING AN EXISTENCE OF FICTITIOUS PARTS OF BLOCKS OF C A ZERO THICKNESS SITUATED BETWEEN EVERY NEIGHBOURING C COMPLEX BLOCKS, THE INDICES OF WHICH ARE NOT SEQUENTIALLY C ORDERED. THE NUMBER OF FICTITIOUS BLOCKS IS CONSIDERED TO C BE JUST THE NUMBER, WHICH IS NECESSARY TO FILL THE GAP C BETWEEN INDICES OF THE TWO NEIGHBOURING COMPLEX BLOCKS. C THE RAY CAN ONLY PASS THROUGH THE FICTITIOUS BLOCKS, NO C REFLECTION IN FICTITIOUS BLOCKS BEING ALLOWED. SIMILARLY, C NO CONVERSION IS ALLOWED NEITHER AT INTERFACES BETWEEN C FICTITIOUS BLOCKS NOR AT THE INTERFACE AT WHICH THE RAY C LEAVES A FICTITIOUS BLOCK (CONVERSION MAY TAKE PLACE ONLY C AT THE INTERFACE AT WHICH THE RAY ENTERS THE FICTITIOUS C BLOCK). IF ANY OF THESE PROHIBITED SITUATIONS IS C SPECIFIED BY THE CODE, IT LEADS TO THE TERMINATION OF C COMPUTATIONS OF THE CORRESPONDING RAY. C C THIS CODE IS A GENERALIZATION OF THE CODE USED IN THE 2-D C PROGRAM PACKAGE 'SEIS83' AND IS DESCRIBED IN CERVENY, C MOLOTKOV AND PSENCIK: RAY METHOD IN SEISMOLOGY, CHARLES C UNIVERSITY, PRAGUE 1977. IT ENABLES ANY BLOCK STRUCTURE C TO BE INTERPRETED AS LOCALLY LAYERED STRUCTURE WITH C FICTITIOUS LAYERS. THUS, THE ROUTINES FOR AUTOMATIC OR C SEMIAUTOMATIC GENERATION OF RAY CODES FOR LAYERED C STRUCTURES MAY BE USED FOR THIS CODE EVEN IN GENERAL BLOCK C STRUCTURES. ITS EFFECTIVE USE IS CONDITIONED BY PROPER C INDEXING OF COMPLEX BLOCKS, WHICH MINIMIZES THE NUMBER OF C FICTITIOUS LAYERS. FOR EXAMPLE, IN A LAYERED STRUCTURE, C THE BLOCKS (LAYERS) SHOULD BE INDEXED SEQUENTIALLY FROM C THE TOP TO THE BOTTOM. C C (5) REFLECTION-TRANSMISSION CODE C FOR ANY ELEMENT OF A RAY ONE INTEGER IS USED. THE FIRST C ELEMENT OF THE RAY, I.E. THE ELEMENT CONTAINING THE C INITIAL POINT OF THE RAY, IS DENOTED BY 1 FOR P AND BY -1 C FOR S WAVE. ANY OTHER ELEMENT IS DENOTED BY 1 OR -1 (P OR C S WAVE) IF THE RAY IS TRANSMITTED AT THE INITIAL POINT OF C THE ELEMENT, AND BY 2 OR -2 (P OR S WAVE) IF THE RAY IS C REFLECTED AT THE INITIAL POINT OF THE ELEMENT. THE CODE C OF A RAY IS A CHAIN OF THE ABOVE NUMBERS CORRESPONDING TO C INDIVIDUAL ELEMENTS OF THE RAY FROM ITS INITIAL POINT TO C ITS ENDPOINT. IN A LAYERED MODEL, THIS CODE IS VERY C CONVENIENT FOR THE DESCRIPTION OF REFRACTED AND PRIMARY C REFLECTED WAVES. THE APPLICATION OF THIS CODE IS C STRAIGHTFORWARD. C C SPECIFICATION OF AN ELEMENTARY WAVE C AN ELEMENTARY WAVE IS SPECIFIED BY THE FOLLOWING DATA: C (1) INTEGER KODTYP, WHICH DETERMINES THE TYPE OF THE CODE C 1... BLOCK-SURFACE CODE, C 2... BASIC BLOCK CODE, C 3... COMPOUND-ELEMENT BLOCK CODE, C 4... GENERALIZED LAYER CODE, C 5... TRANSMISSION-REFLECTION CODE. C (2) INTEGER ARRAY KODE CONTAINING THE CODE OF THE ELEMENTARY C WAVE. THE CODE IS A FINITE SEQUENCE OF NONZERO INTEGERS. C A ZERO INDICATES THE END OF A CODE IN THE INPUT DATA. C THE DATA (1) AND (2) SHOULD BE STORED IN COMMON BLOCK C -------------------------------- C COMMON /COD/ KODTYP,KODE(MKODE). C -------------------------------- C THE DIMENSION MKODE OF THE ARRAY KODE MAY BE ADJUSTED BY THE USER. C C DATE: 1988, JUNE 3 C WRITTEN BY VLASTISLAV CERVENY, LUDEK KLIMES, IVAN PSENCIK C C======================================================================= C SUBROUTINE CODE1(LUN,IWAVE,IWAVE0,IKODE) INTEGER LUN,IWAVE,IWAVE0,IKODE C C THIS SUBROUTINE IS CALLED WHEN STARTING THE COMPUTATION OF A NEW C ELEMENTARY WAVE. IT STORES THE CODE OF THE ELEMENTARY WAVE INTO THE C COMMON BLOCK /COD/. C C INPUT: C LUN... LOGICAL UNIT NUMBER OF THE EXTERNAL INPUT DEVICE C CONTAINING THE INPUT DATA. C IWAVE...ZERO WHEN STARTING THE COMPLETE RAY TRACING PROGRAM, C OTHERWISE THE INDEX OF THE LAST ELEMENTARY WAVE WHICH HAS C BEEN COMPUTED (I.E. THE OUTPUT FROM THE PREVIOUS C INVOCATION OF THE SUBROUTINE CODE1). C C OUTPUT: C IWAVE...ZERO IF ALL REQUIRED ELEMENTARY WAVES ARE COMPUTED AND THE C COMPLETE RAY TRACING PROGRAM WILL BE TERMINATED. C OTHERWISE, THE INDEX OF THE ELEMENTARY WAVE WHICH WILL BE C COMPUTED (I.E. NSKIP+1 FOR THE FIRST COMPUTED ELEMENTARY C WAVE, WHERE NSKIP IS THE NUMBER OF ELEMENTARY WAVES HAVING C BEEN SKIPPED OVER, OTHERWISE THE INPUT VALUE INCREASED BY C ONE). C IWAVE0..INDEX OF THE ALREADY COMPUTED ELEMENTARY WAVE HAVING THE C MOST NUMEROUS COMMON ELEMENTS WITH THE CURRENT ELEMENTARY C WAVE. IN THE CASE OF SEVERAL POSSIBILITIES, THE FIRST C COMPUTED WAVE OF THEM IS TAKEN. C IKODE...THE LENGTH OF THE COMMON PART OF THE CODES OF THE IWAVE-TH C AND IWAVE0-TH ELEMENTARY WAVES. C C COMMON BLOCK /COD/: INTEGER MKODE PARAMETER (MKODE=128) INTEGER KODTYP,KODE(MKODE) COMMON /COD/ KODTYP,KODE C ALL THE STORAGE LOCATIONS OF THE COMMON BLOCK ARE DEFINED IN THIS C SUBROUTINE. C C SUBROUTINES AND EXTERNAL FUNCTIONS REQUIRED: EXTERNAL CODEB,LCODE C CODEB.. BLOCK DATA SUBROUTINE OF THIS FILE. C LCODE.. THIS FILE. C C ERROR MESSAGES: C 401... INSUFFICIENT MEMORY FOR CODES: C THE DIMENSION MKODES OF THE ARRAY KODES(0:MKODES) IN C THIS SUBROUTINE SHOULD BE INCREASED. THE DIMENSION MKODES C SHOULD AT LEAST EQUAL TO THE NUMBER OF ELEMENTARY WAVES C PLUS THE TOTAL NUMBER OF ALL INDECES FORMING THE CODES OF C THE WAVES. C NOTE: C IF THE PARAMETER MODCRT OF THE DATA SET 'DCRT' (SEE SOURCE C CODE FILE 'RAY.FOR') IS ZERO, THE OUTPUT PARAMETERS IWAVE0 C AND IKODE OF THIS SUBROUTINE ARE LIKELY UNUSED. THEN, C INSTEAD OF INCREASING MKODES, THE PAUSE AND STOP C STATEMENTS GENERATING THIS ERROR MESSAGE MAY BE DISABLED C TO SAVE THE COMPUTER MEMORY. THEN THIS SUBROUTINE WILL C WORK, BUT CONSIDERING ONLY THE FIRST SEVERAL ELEMENTARY C WAVES WHEN CALCULATING OUTPUT VALUES OF IWAVE0 AND IKODE. C C DATE: 1993, DECEMBER 18 C CODED BY LUDEK KLIMES C C----------------------------------------------------------------------- C C STORAGE FOR PREVIOUS CODES: INTEGER MKODES PARAMETER (MKODES=1024) INTEGER KODES(0:MKODES) SAVE KODES C C AUXILIARY STORAGE LOCATIONS: CHARACTER*80 TEXTC INTEGER NSKIP,I,J,K,L SAVE NSKIP C C TEXTC...THE NAME OF THE DATA. STRING OF 80 CHARACTERS. C NSKIP...THE NUMBER OF ELEMENTARY WAVES TO BE SKIPPED OVER. C I,J... AUXILIARY LOOP VARIABLES. C K... AUXILIARY VARIABLE - INDEX IN ARRAY KODES. C L... AUXILIARY VARIABLE - LENGTH OF THE CURRENT CODE. C C....................................................................... C IF(IWAVE.EQ.0) THEN C READING THE NAME OF THE INPUT DATA READ(LUN,*) TEXTC C C READING THE NUMBER OF ELEMENTARY WAVES TO BE SKIPPED OVER NSKIP=0 READ(LUN,*) NSKIP C NO CODES OF ELEMENTARY WAVES ARE READ AND STORED NOW KODES(0)=0 END IF C C READING THE CODE OF THE ELEMENTARY WAVE 10 CONTINUE KODTYP=0 DO 11 I=1,MKODE KODE(I)=0 11 CONTINUE READ(LUN,*) KODTYP,KODE IF(KODTYP.EQ.0) THEN IWAVE=0 RETURN ELSE K=KODES(IWAVE)+1 L=LCODE() IWAVE=IWAVE+1 IF(K+L.GT.MKODES) THEN PAUSE 'ERROR 401 IN CODE1: INSUFFICIENT MEMORY FOR CODES' STOP ELSE C STORING THE CODES DO 12 I=L,1,-1 KODES(K+I)=KODE(I) 12 CONTINUE DO 13 I=K-1,IWAVE,-1 KODES(I+1)=KODES(I) 13 CONTINUE KODES(IWAVE)=K+L DO 14 I=IWAVE-1,0,-1 KODES(I)=KODES(I)+1 14 CONTINUE END IF END IF IF(IWAVE.LE.NSKIP) GO TO 10 C C DETERMINING ELEMENTARY WAVE HAVING THE MOST NUMEROUS COMMON C ELEMENTS WITH THE CURRENT ELEMENTARY WAVE IKODE=0 IWAVE0=0 L=LCODE() DO 39 J=1,MIN0(KODES(0),IWAVE-1) K=KODES(J-1) DO 31 I=1,MIN0(KODES(J)-K,L) IF(KODE(I).NE.KODES(K+I)) THEN GO TO 32 END IF 31 CONTINUE 32 CONTINUE I=I-1 IF(I.GT.IKODE) THEN IKODE=I IWAVE0=J END IF 39 CONTINUE RETURN END C C======================================================================= C INTEGER FUNCTION LCODE() C C INTEGER FUNCTION KOOR IS DESIGNED TO RETURN THE LENGTH OF THE CODE OF C THE CURRENT ELEMENTARY WAVE. C C NO INPUT. C C OUTPUT: C LCODE...LENGTH OF THE CODE OF THE CURRENT ELEMENTARY WAVE, I.E. C THE COUNT OF THE NUMERIC ITEMS IN THE ARRAY KODE WHICH C DESCRIBES THE BEHAVIOUR OF A RAY AT INTERFACES. THE ARRAY C IS CALLED HERE THE CODE OF ELEMENTARY WAVE. C C COMMON BLOCK /COD/: INTEGER MKODE PARAMETER (MKODE=128) INTEGER KODTYP,KODE(MKODE) COMMON /COD/ KODTYP,KODE C NONE OF THE STORAGE LOCATIONS OF THE COMMON BLOCK ARE ALTERED. C C NO SUBROUTINES AND EXTERNAL FUNCTIONS REQUIRED. C C DATE: 1989, DECEMBER 19 C CODED BY LUDEK KLIMES C C----------------------------------------------------------------------- C C AUXILIARY STORAGE LOCATIONS: INTEGER I C DO 1 I=1,MKODE IF(KODE(I).EQ.0) THEN LCODE=I-1 GO TO 2 END IF 1 CONTINUE LCODE=MKODE 2 CONTINUE RETURN END C C======================================================================= C INTEGER FUNCTION NELEM(KODIND) INTEGER KODIND C C THIS FUNCTION IS DESIGNED TO RETURN THE MAXIMUM POSSIBLE NUMBER OF RAY C ELEMENTS CORRESPONDING TO THE GIVEN POSITION IN THE CODE OF THE C ELEMENTARY WAVE. C C INPUT: C KODIND..POSITION IN THE CODE (INDEX IN ARRAY KODE) C NONE OF THE INPUT PARAMETERS ARE ALTERED. C C OUTPUT: C NELEM...NUMBER OF RAY ELEMENTS BETWEEN THE INITIAL POINT OF THE C RAY AND THE END OF THE ELEMENT CORRESPONDING TO THE GIVEN C POSITION IN THE CODE. NOTE THAT, FOR A PARTICULAR RAY, C SOME OF THESE POSSIBLE RAY ELEMENTS MAY BE LEFT OUT OR C COUPLED INTO ONE VIRTUAL ELEMENT. C C COMMON BLOCK /COD/: INTEGER MKODE PARAMETER (MKODE=128) INTEGER KODTYP,KODE(MKODE) COMMON /COD/ KODTYP,KODE C NONE OF THE STORAGE LOCATIONS OF THE COMMON BLOCK ARE ALTERED. C C NO SUBROUTINES AND EXTERNAL FUNCTIONS REQUIRED. C C DATE: 1989, DECEMBER 19 C CODED BY LUDEK KLIMES C C----------------------------------------------------------------------- C C NO AUXILIARY STORAGE LOCATIONS. C IF(KODTYP.EQ.1) THEN NELEM=(KODIND+1)/2 ELSE NELEM=KODIND END IF RETURN END C C======================================================================= C SUBROUTINE CODE(IY,KODIND,ICBNEW,IEND) INTEGER IY(12),KODIND,ICBNEW,IEND C C THIS SUBROUTINE TRANSFORMS THE USED NUMERICAL CODE OF THE ELEMENTARY C WAVE UNDER CONSIDERATION INTO INSTRUCTIONS SPECIFYING THE BEHAVIOUR C OF THE WAVE AT THE INITIAL POINTS OF RAYS AND AT ALL POINTS OF C INCIDENCE OF THE RAYS AT INTERFACES (BOUNDARIES OF COMPLEX BLOCKS). C THUS, IT IS FIRST CALLED AT THE INITIAL POINT OF A RAY AND THEN C SUCCESSIVELY AT ALL POINTS OF INCIDENCE. C C INPUT: C IY... INTEGER ARRAY. ITS NUMERICAL STORAGE UNITS 2, 3, 5, 6, 8 C MUST BE DEFINED AS FOLLOWS: C IY(2)... POSITION IN THE CODE (INDEX IN ARRAY KODE) C CORRESPONDING TO THE LAST COMPUTED ELEMENT OF A RAY. C IY(2)=0 AT THE INITIAL POINT OF THE RAY, C IY(2)=KODIND FROM THE LAST INVOCATION OF SUBROUTINE CODE C AT OTHER POINTS OF THE RAY. C IY(3)=ICB0... INDEX OF THE NEIGHBOURING COMPLEX BLOCK, C FROM WHICH THE RAY ENTERED THE COMPLEX BLOCK IN WHICH C THE LAST COMPUTED ELEMENT OF THE RAY IS SITUATED. C IY(3)=0 BEFORE LEAVING THE COMPLEX BLOCK, IN WHICH THE C INITIAL POINT OF THE RAY IS SITUATED. C AT THE INITIAL POINT OF THE RAY (IY(2)=0), IY(3) IS C IGNORED. C IY(5)=ICB1... INDEX OF THE COMPLEX BLOCK CONTAINING THE C COMPUTED ELEMENT OF THE RAY, SUPPLEMENTED BY THE SIGN C '+' FOR P WAVE AND SIGN '-' FOR S WAVE. C ICB1 IS IGNORED AT THE INITIAL POINT OF THE RAY C (IY(2)=0). C IY(6)=ISRF... INDEX OF THE SURFACE AT WHICH THE ENDPOINT C OF THE LAST COMPUTED ELEMENT OF THE RAY IS SITUATED. C THE SIGN OF IY(6) IS IGNORED. C IY(6)=0 AT THE INITIAL POINT OF THE RAY. C IY(8)=ICB2... INDEX OF THE COMPLEX BLOCK TOUCHING THE C COMPLEX BLOCK ICB1 FROM THE OTHER SIDE OF THE SURFACE C ISRF AT THE ENDPOINT OF THE LAST COMPUTED ELEMENT OF THE C RAY. IY(8)=0 FOR A FREE SPACE. C AT THE INITIAL POINT OF THE RAY, ICB2 IS THE INDEX OF C THE COMPLEX BLOCK CONTAINING THE INITIAL POINT. C THE INPUT PARAMETER IS NOT ALTERED. C C OUTPUT: C KODIND..POSITION IN THE CODE (INDEX IN THE ARRAY KODE) CORRES- C PONDING TO THE NEXT ELEMENT OF THE RAY. C ICBNEW..THE INDEX OF THE COMPLEX BLOCK IN WHICH THE NEXT C ELEMENT OF THE RAY IS TO BE SITUATED, SUPPLEMENTED C BY THE SIGN "+" FOR P WAVE OR "-" FOR S WAVE. C IEND... INFORMATION ON THE PROCESS OF THE INTERPRETATION C OF THE CODE: C IEND.EQ.0... COMPUTATION OF THE RAY CONTINUES. C IEND.NE.0... THE COMPUTATION OF THE RAY TERMINATES. C DIFFERENT VALUES OF IEND SPECIFY THE REASON FOR THE C TERMINATION: C 10... RAY SATISFIES THE WHOLE CODE (REGULAR TERMINATION C OF THE RAY). C 21... THE POINT OF INCIDENCE IS SITUATED AT A DIFFERENT C SURFACE THAN THAT REQUIRED BY THE CODE. C 22... THE NEXT ELEMENT OF THE RAY IS REQUIRED BY THE C CODE TO BE SITUATED IN A COMPLEX BLOCK THAT DOES NOT C TOUCH THE POINT OF INCIDENCE. C 23... TRANSMISSION IS REQUIRED BY THE CODE AT A FREE C SURFACE. C 30... REFLECTION OR WAVE CONVERSION AT THE FICTITIOUS C PART OF THE INTERFACE. C C COMMON BLOCK /COD/: C ------------------------------------------------------------------ INTEGER MKODE PARAMETER (MKODE=128) INTEGER KODTYP,KODE(MKODE) COMMON /COD/ KODTYP,KODE C ------------------------------------------------------------------ C KODTYP..DETERMINES THE TYPE OF THE CODE: C 1... BLOCK-SURFACE CODE, C 2... BASIC BLOCK CODE, C 3... COMPOUND-ELEMENT BLOCK CODE, C 4... GENERALIZED LAYER CODE, C 5... TRANSMISSION-REFLECTION CODE. C KODE... ARRAY CONTAINING THE CODE OF THE ELEMENTARY WAVE. C THE CODE IS A FINITE SEQUENCE OF NONZERO INTEGERS. C ZERO INDICATES THE END OF A CODE. C NONE OF THE STORAGE LOCATIONS OF THE COMMON BLOCK ARE ALTERED. C C NO SUBROUTINES AND EXTERNAL FUNCTIONS REQUIRED. C C DATE: 1990, NOVEMBER 18 C CODED BY LUDEK KLIMES C C----------------------------------------------------------------------- C C AUXILIARY STORAGE LOCATIONS: INTEGER I,J C KODIND=IY(2) GO TO (10,20,30,30,50),KODTYP C C BLOCK-SURFACE CODE: 10 CONTINUE IF(KODIND.GT.0) THEN C CHECK ON THE SURFACE: KODIND=KODIND+1 IF(KODE(KODIND).EQ.0) THEN C END OF CODE IEND=10 RETURN ELSE IF(KODE(KODIND).NE.IABS(IY(6))) THEN C WRONG SURFACE IEND=21 RETURN END IF END IF C THE REST OF INTERPRETATION IS THE SAME AS IN THE BASIC BLOCK CODE C C BASIC BLOCK CODE: 20 CONTINUE KODIND=KODIND+1 ICBNEW=KODE(KODIND) IF(ICBNEW.EQ.0) THEN C END OF CODE IEND=10 ELSE IF(IABS(ICBNEW).EQ.IY(8)) THEN C TRANSMISSION IEND=0 ELSE IF(KODIND.GT.1.AND.IABS(ICBNEW).EQ.IABS(IY(5))) THEN C REFLECTION IEND=0 ELSE C REQUIRED COMPLEX BLOCK IS NOT ATTAINABLE IEND=22 END IF RETURN C C COMPOUND-ELEMENT BLOCK CODE: 30 CONTINUE IF(KODIND.GT.0) THEN I=IABS(IY(5)) J=(I-IABS(IY(3)))*(I-IY(8)) DO 31 I=KODIND-1,1,-1 IF(KODE(I).NE.KODE(KODIND)) GO TO 32 J=-J 31 CONTINUE 32 CONTINUE IF(J.GT.0) THEN C COMPOUND ELEMENT: KODIND=KODIND+1 IF(KODE(KODIND).NE.IY(5)) THEN C WRONG SURFACE IEND=21 RETURN END IF END IF END IF C REST OF INTERPRETATION IS THE SAME AS IN BASIC BLOCK CODE IF(KODTYP.EQ.3) GO TO 20 C C GENERALIZED LAYER CODE: C INTERPRETATION OF COMPOUND ELEMENTS IS THE SAME AS IN C COMPOUND-ELEMENT BLOCK CODE, THEN: KODIND=KODIND+1 ICBNEW=KODE(KODIND) IF(ICBNEW.EQ.0) THEN C END OF CODE IEND=10 ELSE IF(KODIND.EQ.1) THEN C INITIAL POINT OF A RAY: IF(IABS(ICBNEW).EQ.IY(8)) THEN C INITIAL POINT OF A RAY IS SITUATED IN THE REQUIRED C.B. IEND=0 ELSE C INITIAL POINT OF A RAY IS NOT SITUATED IN THE REQUIRED C.B. IEND=22 END IF ELSE IF(IABS(ICBNEW).EQ.IABS(IY(5))) THEN C REFLECTION IEND=0 ELSE C POSSIBLE TRANSMISSION: J=ISIGN(1,IY(8)-IABS(IY(5))) IF(IABS(ICBNEW).EQ.IABS(IY(5))+J) THEN C LOOP FOR FICTITIOUS PARTS OF BLOCKS: DO 41 I=IABS(IY(5))+J+J,IY(8),J KODIND=KODIND+1 ICBNEW=KODE(KODIND) IF(ICBNEW.EQ.0) THEN C END OF CODE IN THE FICTITIOUS PART OF BLOCK IEND=10 RETURN ELSE IF(ICBNEW.EQ.ISIGN(I,KODE(KODIND-1))) THEN C TRANSMISSION FROM THE FICTITIOUS PART OF BLOCK IEND=0 ELSE C TERMINATION OF THE RAY COMPUTATION IN THE FICTITIOUS PART C OF THE BLOCK IEND=30 RETURN END IF 41 CONTINUE ELSE C REQUIRED COMPLEX BLOCK IS NOT ATTAINABLE IEND=22 END IF END IF RETURN C C TRAMSMISSION-REFLECTION CODE: 50 CONTINUE KODIND=KODIND+1 I=KODE(KODIND) IF(IABS(I).EQ.1) THEN C TRANSMISSION: IF(IY(8).GT.0) THEN C TRANSMISSION INTO MATERIAL BLOCK ICBNEW=ISIGN(IY(8),I) IEND=0 ELSE C TRANSMISSION INTO FREE SPACE IEND=23 END IF ELSE IF(IABS(I).EQ.2.AND.KODIND.GT.1) THEN C REFLECTION ICBNEW=ISIGN(IABS(IY(5)),I) IEND=0 ELSE C END OF CODE IEND=10 END IF RETURN END C C======================================================================= C