C
C *** PHASE AND GROUP VELOCITY AND 
C     AMPLITUDE SPECTRUM DUE TO VERTICAL POINT-FOORCE ON FREE SURFACE 
C     OF RAYLEIGH WAVE
C
      PROGRAM MRCUAV
C
C ** MARAMETER **
C MOM   : NOM IN MAIN PROGRAM DIMENSION
C       : IF MOM > 2**13(=8192), CAN NOT USE IN SUBROUTINE !!!
C MMD   : NEMD IN MAIN PROGRAM DIMENSION
C ML    : NL (LAYER NUMBER)
C ** INPUT DATA **
C DOM   : INCREMENT OF CIRCULAR FREQUENCY ( DELTA OMEGA )
C NOM   : TOTAL NUMBER OF CIRCULAR FREQUENCY
C DNS(I): DENSITY OF Ith LAYER 
C VS (I): VS OF Ith LAYER      
C DEP(I): DEPTH ( THICKNESS ) OF Ith LAYER   
C NL    : TOTAL NUMBER OF LAYERS             
C NDC   : DIVISION NUMBER FROM MIN VS TO MAX VS TO PHASE VELOCITY
C NSMD  : START OUTPUT MODE OF RAYLEIGH WAVE
C NEMD  : END OUTPUT MODE OF RAYLEIGH WAVE
C NDOM  : DIVISION NUMBER FROM 0 TO MAX OMEGA TO CUT-OFF OMEGA
C ** CALCULATE DATA **
C CFOM(I) : CUT-OFF OMEGA OF Ith MODE RAYLEIGH WAVE
C NCFOM   : NUMBER OF CUT-OFF OMEGA
C CC(I,J) : RAYLEIGH WAVE PHASE VELOCITY OF Jth MODE AT Ith OMEGA
C NMODE(I): OUTPUT HIGHEST RAYLEIGH MODE NUMBER AT Ith OMEGA
C UU(I,J) : GROUP VELOCITY OF Jth MODE RAYLEIGH WAVE AT Ith OMEGA
C ARK(I,J): AMPLITUDE SPECTRUM (AL/K) OF Jth MODE RAYLEIGH WAVE
C           AT Ith OMEGA
C UW0(I,J): AMPLITUDE RATIO ( IMAG(RK/RL) ) OF Jth MODE RAYLEIGH WAVE
C           AT Ith OMEGA
C
      PARAMETER(MOM=100,MMD=10,ML=30)
      IMPLICIT REAL*8(A-H,O-Z)
      DIMENSION CFOM(0:MMD),CC(MOM,0:MMD),NMODE(MOM)
     *         ,UU(MOM,0:MMD),ARK(MOM,0:MMD),UW0(MOM,0:MMD)
     *         ,DNS(ML),VP(ML),VS(ML),DEP(ML) 
C
C *** DATA INPUT ***
C
      CALL RAYINN(DOM,NOM,ML,NL,DNS,VP,VS,DEP,NDOM,NSMD,NEMD,NDC)
C
C *** PARAMETER CHECK ***
C
      IF(MOM.LT.NOM)THEN
        WRITE(*,*)'FROM MAIN MRCUAV :'
        WRITE(*,*)'MOM < NOM THUS STOP !!!'
        WRITE(*,*)'MOM,NOM =',MOM,NOM
        STOP
      END IF
C
      IF(ML.LT.NL)THEN
        WRITE(*,*)'FROM MAIN MRCUAV :'
        WRITE(*,*)'ML < NL THUS STOP !!!'
        WRITE(*,*)'ML,NL =',ML,NL
        STOP
      END IF
C
      IF(MMD.LT.NEMD)THEN
        WRITE(*,*)'FROM MAIN MRCUAV :'
        WRITE(*,*)'MMD < NEMD THUS STOP !!!'
        WRITE(*,*)'MMD,NEMD =',MMD,NEMD
        STOP
      END IF
C
C *** CUT-OFF OMEGA ***
C
       CALL RAYCTN(DOM,NOM,ML,NL,DNS,VP,VS,DEP,NDOM,CFOM,NEMD,MMD,NCFOM)
C
C *** PHASE VELOCITY ***
C
      CALL RAYCCN(DOM,NOM,MOM,ML,NL,DNS,VP,VS,DEP
     *           ,NDC,CFOM,NSMD,NEMD,MMD,CC,NMODE)
C
C *** GROUP VELOCITY AND AMPLITUDE SPECTRUM ***
C
      CALL RAYUAN(DOM,NOM,MOM,ML,NL,DNS,VP,VS,DEP
     *           ,CC,NMODE,CFOM,NSMD,MMD,UU,ARK,UW0)
C
C *** DATA OUTPUT ON SCREEN ***
C
      CALL RAYOTN(DOM,NOM,MOM,ML,NL,VS
     *           ,CFOM,NCFOM,CC,NSMD,MMD,UU,ARK,UW0)
C
C *** DATA OUTPUT ON PLOTTER ***
C
      CALL RAYPTN(DOM,NOM,MOM,ML,NL,VS
     *           ,CFOM,NCFOM,CC,NSMD,MMD,UU,ARK,UW0)
C
      STOP
      END
