C
C *** OUTPUT OF RAYLEIGH WAVE PHASE AND GROUP VELOCITY AND 
C     AMPLITUDE SPECTRA ON PLOTTER DATA FILE
C
      SUBROUTINE RAYPT(DOM,NOM,MOM,CFOM,NCFOM,CC,NSMD,MMD,UU,ARK,UW0)
C
C ** INPUT DATA **
C DOM     : DELTA OMEGA
C NOM     : TOTAL NUMBER OF CIRCULAR FREQUENCY
C MOM     : NOM IN MAIN PROGRAM DIMENSION
C         : IF MOM > 2**13(=8192), CAN NOT USE IN SUBROUTINE !!!
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 NSMD    : OUTPUT START RAYLEIGH MODE
C MMD     : NSMD IN MAIN PROGRAM DIMENSION
C UU(I,J) : GROUP VELOCITY OF Jth MODE RAYLEIGH WAVE AT Ith OMEGA
C ARK(I,J): AMPLITUDE SPECTRA (AL/K) OF Jth MODE RAYLEIGH WAVE 
C           AT Ith OMEGA
C UW0(I,J): AMPLITUDE RATIO ( IMAG(U0/W0) ) OF Jth MODE RAYLEIGH WAVE 
C           AT Ith OMEGA
C
      IMPLICIT REAL*8 (A-H,O-Z)
      DIMENSION CFOM(0:MMD),CC(MOM,0:MMD)
     *         ,UU(MOM,0:MMD),ARK(MOM,0:MMD),UW0(MOM,0:MMD)
      COMMON / LAYER / NL,VP(50),VS(50),DEP(50),SM(50)
     *       / MINMAX / VSMIN,VSMAX,CCMIN,CCMAX,MINIL,MAXIL
C
C *** FINAL OMEGA ***
C
      EOM=DOM*DFLOAT(NOM)
C
C *** OUTPUT OF DATA INFORMATION ***
C
      WRITE(*,*)
      WRITE(*,1000)
   10 WRITE(*,1010)
      WRITE(*,1020)DOM,NOM
      WRITE(*,1030)CFOM(NSMD),EOM
      WRITE(*,1040)CCMIN,VSMAX
 1000 FORMAT(' From SUBROUTINE RAYPT : ')
 1010 FORMAT(/' DATA INFORMATION :')
 1020 FORMAT(' Delta Omega = ',F8.4,' ,Omega Number = ',I5)
 1030 FORMAT(' First Omega = ',F8.4,' ,Final Omega  = ',F8.4)
 1040 FORMAT(' Minimum C   = ',F8.3,' ,Maximum Vs   = ',F8.3)
C
C *** AXIS DATA INPUT ***
C
      WRITE(*,1200)
 1200 FORMAT(/' Axis Off ? ( Yes = 0 Enter ) :',\)
      READ(*,*)NAXY
C
      IF(NAXY.EQ.0)THEN
C
        OPEN(5,FILE='AX.DAT')
          READ(5,1210)XOR,YOR,XFAC,YFAC,AMPFAC,XVAL,IYAMP
        CLOSE(5)
 1210   FORMAT(6D15.7,I5)
C
      ELSE
C
        WRITE(*,1300)
        WRITE(*,1310)
        READ(*,*)XVAL
        WRITE(*,1320)
        READ(*,*)XLEN
        XLEN=XLEN*100.0D0
        WRITE(*,1330)
        READ(*,*)NUX
        XFAC=XLEN/XVAL
        WRITE(*,1340)
        READ(*,*)YVAL
        WRITE(*,1350)
        READ(*,*)YLEN
        YLEN=YLEN*100.D0
        WRITE(*,1360)
        READ(*,*)NUY
        YFAC=YLEN/YVAL
        WRITE(*,1370)
        WRITE(*,1380)
        READ(*,*)YAMP
        YAMP=YAMP*100.D0
        WRITE(*,1390)
        READ(*,*)IYAMP
 1300   FORMAT(/' Axis Data Input for Phase and Group Velocity :')
 1310   FORMAT(/' Enter X-Axis Max ( Max Omega : rad/sec )   :',\)
 1320   FORMAT( ' Enter X-Axis Length ( cm )                 :',\)
 1330   FORMAT( ' Enter X-Axis Division Number               :',\)
 1340   FORMAT(/' Enter Y-Axis Max Value ( Max Vs : m/sec )  :',\)
 1350   FORMAT( ' Enter Y-Axis Length ( cm )                 :',\)
 1360   FORMAT( ' Enter Y-Axis Division Number               :',\)
 1370   FORMAT(/' Axis Data Input for Amplitude Spectrum :')
 1380   FORMAT(/' Enter Y-Axis Length ( cm )                 :',\)
 1390   FORMAT(/' Superpose Amp W and U ? ( Yes = 0 Enter )  :',\)
C
        WRITE(*,1500)
        READ(*,*)XOR
        XOR=(2.0D0+XOR)*100.0D0
        WRITE(*,1510)
        READ(*,*)YOR
        YOR=(2.0D0+YOR)*100.0D0
        WRITE(*,1520)
        READ(*,*)IYES
        IF(IYES.NE.0)GOTO 10
 1500   FORMAT(/' Enter X Origin Position ( cm ) :',\)
 1510   FORMAT( ' Enter Y Origin Position ( cm ) :',\)
 1520   FORMAT(/' All Input Data OK ? ( Yes = Enter 0 ) :',\)
C
      END IF
C
C *** PLOTTER DATA FILE OPEN ***
C
      OPEN(6,FILE='PS.DAT')
C
C    ** OUTPUT OF PHASE AND GROUP VELOCITY PLOTTER DATA ***
C 
        WRITE(6,2000)
C
        IF(NAXY.EQ.0)GOTO 20
C
C     * TERM, AXIS AND CHARACTER *
C
        WRITE(6,2010)
        WRITE(6,2020)XOR,YOR
        WRITE(6,2030)XLEN,NUX
        WRITE(6,2020)XOR-110.0D0,YOR-60.0D0
        WRITE(6,2040)0.0D0
        WRITE(6,2020)XOR+XLEN-100.0D0,YOR-60.0D0
        WRITE(6,2040)XVAL
        WRITE(6,2020)XOR+XLEN-220.0D0,YOR-120.0D0
        WRITE(6,2050)
        WRITE(6,2020)XOR,YOR
        WRITE(6,2060)YLEN,NUY
        WRITE(6,2020)XOR-40.0D0,YOR-80.0D0
        WRITE(6,2070)0.0D0
        WRITE(6,2020)XOR-40.0D0,YOR+YLEN-80.0D0
        WRITE(6,2080)YVAL/1000.D0
        WRITE(6,2020)XOR-100.0D0,YOR+YLEN-200.0D0
        WRITE(6,2090)
C
 2000   FORMAT(2HJ1)
 2010   FORMAT(3H=<<)
 2020   FORMAT(1HM,F6.1,1H,F6.1,1H,)
 2030   FORMAT(3HX3,,F6.1,1H,I2,2H<<)
 2040   FORMAT(7HR0,S30,,1HP,F5.1,2H<<)
 2050   FORMAT(7HR0,S30,,12HP(RAD/SEC)<<)
 2060   FORMAT(3HX2,,F6.1,1H,I2,2H<<)
 2070   FORMAT(9HR900,S30,,1HP,F4.1,2H<<)
 2080   FORMAT(9HR900,S30,,1HP,F4.1,2H<<)
 2090   FORMAT(9HR900,S30,,11HP(KM/SEC)<<)
C
   20   DO 100 MODE=NSMD,NCFOM
C
C     * PHASE VELOCITY DATA PLOT *
C
          YV=VS(NL)
          IF(MODE.EQ.0)YV=CCMAX
C
          WRITE(6,2100)XOR+CFOM(MODE)*XFAC,YOR+YV*YFAC
          DO 110 IOM=1,NOM
            OMX=DOM*DFLOAT(IOM)
            IF(OMX.LT.CFOM(MODE))GOTO 110
            WRITE(6,2110)XOR+OMX*XFAC,YOR+CC(IOM,MODE)*YFAC
            IF(OMX.GE.XVAL)GO TO 30
  110     CONTINUE
C
 2100     FORMAT(1HM,F6.1,1H,F6.1,1H,)
 2110     FORMAT(1HD,F6.1,1H,F6.1,1H,)
C
C       * GROUP VELOCITY DATA PLOT *
C
   30     WRITE(6,2200)XOR+CFOM(MODE)*XFAC,YOR+YV*YFAC
          DO 120 IOM=1,NOM
            OMX=DOM*DFLOAT(IOM)
            IF(OMX.LT.CFOM(MODE))GOTO 120
            WRITE(6,2110)XOR+OMX*XFAC,YOR+UU(IOM,MODE)*YFAC
            IF(OMX.GE.XVAL)GO TO 100
  120     CONTINUE
C
 2200     FORMAT(1HM,F6.1,1H,F6.1,1H,)
 2210     FORMAT(1HD,F6.1,1H,F6.1,1H,)
C
  100   CONTINUE
C
C    ** OUTPUT OF AMPLITUDE W SPECTRA PLOTTER DATA **
C
        XORA=XOR+XVAL*XFAC+200.D0
        IF(NAXY.EQ.0)GOTO 40
C
C     * MAXIMUM AMPLITUDE SPECTRA *
C
        IOMMAX=0
        MDMAX=0
        ARKMAX=0.D0
        DO 210 MODE=NSMD,NCFOM
          DO 200 IOM=1,NOM
            ARKI=DABS(ARK(IOM,MODE))
            IF(ARKI.GT.ARKMAX)THEN
              IOMMAX=IOM
              MDMAX=MODE
              ARKMAX=ARKI
            END IF
  200     CONTINUE
  210   CONTINUE
        ARKMAX=ARK(IOMMAX,MDMAX)
        AMPFAC=YAMP/ARKMAX
C
C     * PEN, TERM, AXIS AND CHARACTER *
C
        WRITE(6,2500)
        WRITE(6,2510)
        WRITE(6,2520)XORA,YOR
        WRITE(6,2530)XLEN,NUX
        WRITE(6,2520)XORA-110.0D0,YOR-60.0D0
        WRITE(6,2540)0.0D0
        WRITE(6,2520)XORA+XLEN-100.0D0,YOR-60.0D0
        WRITE(6,2540)XVAL
        WRITE(6,2520)XORA+XLEN-220.0D0,YOR-120.0D0
        WRITE(6,2550)
        WRITE(6,2520)XORA,YOR
        WRITE(6,2560)YAMP,1
        WRITE(6,2520)XORA-40.0D0,YOR-80.0D0
        WRITE(6,2570)0.0D0
        WRITE(6,2520)XORA-40.0D0,YOR+YAMP-250.0D0
        WRITE(6,2580)ARKMAX
C
 2500   FORMAT(2HJ1)
 2510   FORMAT(3H=<<)
 2520   FORMAT(1HM,F6.1,1H,F6.1,1H,)
 2530   FORMAT(3HX3,,F6.1,1H,I2,2H<<)
 2540   FORMAT(7HR0,S30,,1HP,F5.1,2H<<)
 2550   FORMAT(7HR0,S30,,12HP(RAD/SEC)<<)
 2560   FORMAT(3HX2,,F6.1,1H,I2,2H<<)
 2570   FORMAT(9HR900,S30,,1HP,F4.1,2H<<)
 2580   FORMAT(9HR900,S30,,1HP,D11.3,2H<<)
C
C     * AMPLITUDE SPECTRA DATA PLOT ( W ) *
C
   40   DO 300 MODE=NSMD,NCFOM
C
          WRITE(6,2700)XORA+CFOM(MODE)*XFAC,YOR
          DO 310 IOM=1,NOM
            OMX=DOM*DFLOAT(IOM)
            IF(OMX.LT.CFOM(MODE))GOTO 310
            WRITE(6,2710)XORA+OMX*XFAC,YOR+ARK(IOM,MODE)*AMPFAC
            IF(OMX.GE.XVAL)GO TO 300
  310     CONTINUE
C
 2700     FORMAT(1HM,F6.1,1H,F6.1,1H,)
 2710     FORMAT(1HD,F6.1,1H,F6.1,1H,)
C
  300   CONTINUE
C
C    ** OUTPUT OF AMPLITUDE U SPECTRA PLOTTER DATA **
C
        IF(IYAMP.EQ.0)THEN
          WRITE(6,2800)
          GOTO 50
        END IF
 2800   FORMAT(2HJ2)
C
        XORA=XORA+XVAL*XFAC+200.D0
        IF(NAXY.EQ.0)GOTO 50
C
C     * PEN, TERM, AXIS AND CHARACTER *
C
      
        WRITE(6,3000)
        WRITE(6,3010)
        WRITE(6,3020)XORA,YOR
        WRITE(6,3030)XLEN,NUX
        WRITE(6,3020)XORA-110.0D0,YOR-60.0D0
        WRITE(6,3040)0.0D0
        WRITE(6,3020)XORA+XLEN-100.0D0,YOR-60.0D0
        WRITE(6,3040)XVAL
        WRITE(6,3020)XORA+XLEN-220.0D0,YOR-120.0D0
        WRITE(6,3050)
        WRITE(6,3020)XORA,YOR
        WRITE(6,3060)YAMP,1
        WRITE(6,3020)XORA-40.0D0,YOR-80.0D0
        WRITE(6,3070)0.0D0
        WRITE(6,3020)XORA-40.0D0,YOR+YAMP-250.0D0
        WRITE(6,3080)ARKMAX
C
 3000   FORMAT(2HJ1)
 3010   FORMAT(3H=<<)
 3020   FORMAT(1HM,F6.1,1H,F6.1,1H,)
 3030   FORMAT(3HX3,,F6.1,1H,I2,2H<<)
 3040   FORMAT(7HR0,S30,,1HP,F5.1,2H<<)
 3050   FORMAT(7HR0,S30,,12HP(RAD/SEC)<<)
 3060   FORMAT(3HX2,,F6.1,1H,I2,2H<<)
 3070   FORMAT(9HR900,S30,,1HP,F4.1,2H<<)
 3080   FORMAT(9HR900,S30,,1HP,D11.3,2H<<)
C
C     * AMPLITUDE SPECTRA DATA PLOT ( U ) *
C
   50   DO 500 MODE=NSMD,NCFOM
C
          WRITE(6,3210)XORA+CFOM(MODE)*XFAC,YOR
          DO 510 IOM=1,NOM
            OMX=DOM*DFLOAT(IOM)
            IF(OMX.LT.CFOM(MODE))GOTO 510
            ARKU=DABS(ARK(IOM,MODE)*AMPFAC*UW0(IOM,MODE))
            WRITE(6,3220)XORA+OMX*XFAC,YOR+ARKU
            IF(OMX.GE.XVAL)GO TO 500
  510     CONTINUE
C
 3210     FORMAT(1HM,F6.1,1H,F6.1,1H,)
 3220     FORMAT(1HD,F6.1,1H,F6.1,1H,)
C
  500   CONTINUE
C
C     * SKIP OF PLOTTER PEN *
C
        WRITE(6,3300)XORA+XVAL*XFAC+500.D0,YOR
 3300   FORMAT(1HM,F6.1,1H,F6.1,1H,)
C
C *** POTTER DATA FILE CLOSE ***
C
      CLOSE(6)
C
C *** AXIS DATA OUTPUT TO FILE ***
C
      IF(NAXY.NE.0)THEN
        OPEN(5,FILE='AX.DAT')
        WRITE(5,3500)XOR,YOR,XFAC,YFAC,AMPFAC,XVAL,IYAMP
        CLOSE(5)
 3500   FORMAT(6D15.7,I5)
      END IF
C
      WRITE(*,3700)
 3700 FORMAT(/' USE FOR PLOTTER : A> COPY PS.DAT PRN')
      WRITE(*,*)
C
      RETURN
      END
