        COMPLEX E1,ZIP,TEXP,DEN,Y      
        DI=1.33
        F=1.0
        RHO=7.28E-04
        RHOW=9.32E-05
        DKP=0.83
        C0=60000.0
        G=1.15E+07  
        CD=0.0
        AR=4.0
        E=3.0E+07
        N1=1
        N2=20
        N3=1
        IOMG1=1
        IOMG2=2
        IOMG3=3
        IX1=10.0
        IX2=400.0
        IX3=10. 
        IT1=1
        IT2=2
        IT3=3
        PI=3.1416
        XL=200.0
        X0=100.0
        W=2.0
        K2=20
        A=DI*E
        B=RHO*AR
        C=(RHO*DI)+(E*DI*RHO/(DKP*G))
C        C=0
        D=(RHO**2)*(DI/(DKP*G))
C        D=0
C        DO 902 IOMG=IOMG1,IOMG2,IOMG3
C   902     CONTINUE
        DO 4000 IOMG=IOMG1,IOMG2,IOMG3
        OMEG=IOMG*942.0
        DO 3000 IX=IX1,IX2,IX3 
        X=IX*0.5
        DO 2000 IT=IT1,IT2,IT3
        T=0.
        Y=0
        DO 1000 N=N1,N2,N3
C        GO TO 20
10      IF(N.GT.1) GO TO 11
        ALPHA=0.7340
        BETNL=1.8751
        GO TO 12
11      ALPHA=1.0
        BETNL=(2.0*N-1.0)*1.5708
        GO TO 13 
12      BETN=BETNL/XL
        ZR1=COSH(BETN*X)-COS(BETN*X)
        ZR=ZR1-ALPHA*(SINH(BETN*X)-SIN(BETN*X))
        ZR01=COSH(BETN*X0)-COS(BETN*X0)
        ZR0=ZR01-ALPHA*(SINH(BETN*X0)-SIN(BETN*X0))
        GO TO 14
13      BETN=BETNL/XL
        ZR=EXP(-BETN*X)-COS(BETN*X)+SIN(BETN*X)
        ZR0=EXP(-BETN*X0)-COS(BETN*X0)+SIN(BETN*X0)
        GO TO 14        
14      XI=2.0*ALPHA/BETN
        XIS=XL
        XII=ALPHA*(BETN)*(2.0-ALPHA*BETNL)
        XIV=(BETN**4)*XIS
        GO TO 100
20      BETNL=N*PI
        BETN=BETNL/XL
        ZR=SIN(BETN*X)
C        PRINT *,BETNL
        ZR0=SIN(BETN*X0)
        XI=1.0
        XIS=XL/2.0
        XII=-(BETN**2)*XL/2.0
        XIV=(BETN**4)*XL/2.0
        GO TO 100
30      BETNL=(2.0*N+1)*1.5708
        ALPHA=1.0
        BETN=BETNL/XL
        ZR=EXP(-BETN*X)+COS(BETN*X)-SIN(BETN*X)
        ZR0=EXP(-BETN*X0)+COS(BETN*X0)-SIN(BETN*X0)
        XIS=XL
        XII=BETN*(2.0-BETNL)
        XIV=(BETN**4)*XL
100        CHI=-0.48*(OMEG/C0)*SQRT(XL*W)
C        CHI=0.0
        THETA=0.13*((OMEG/C0)**2)*XL*W
C        THETA=0.00
        ZIP=CMPLX(THETA,CHI)
        E1=CMPLX(CD+W*RHOW*C0*THETA,W*RHOW*C0*CHI)
C        GO TO 100
        TEXP=CMPLX(COS(OMEG*T),SIN(OMEG*T))
        AP=A*XIV
        BP=-B*(OMEG**2)*XIS
        CP=C*(OMEG**2)*XII
C        CP=0
        DP=D*(OMEG**4)*XIS
C        DP=0
        REP=(CD+W*RHOW*C0*THETA)*OMEG*XIS
C        REP=0
        AEP=(W*RHOW*C0*CHI)*OMEG*XIS
C        AEP=0 
        DEN=CMPLX(AP+BP+CP+DP-AEP,REP)
C        PRINT *,N,AP,BP,CP,DP,AEP,REP
        Y=Y+(F*ZR0*ZR*TEXP)/DEN
1000    CONTINUE
C        DO 923 IOMG=IOMG1,IOMG2,IOMG3
        Y1=ABS(Y)
C    923     CONTINUE        
C        PRINT *,OMEG,X,T,Y1
2000    CONTINUE
        PRINT *,OMEG,X,T,Y 
3000    CONTINUE
C        PRINT *,OMEG,X,T,Y1
4000    CONTINUE
C        PRINT *,OMEG,X,T,Y1
C         DO 962 IOMG=IOMG1,IOMG2,IOMG3
C         OMEG=IOMG*.1
C         Y1(IOMG)=ABS(Y(IOMG))
C         PRINT *,IOMG,OMEG,Y1(IOMG)
C     962      CONTINUE         
C         DO 965 IOMG=IOMG1+1,IOMG2,IOMG3
C         OMEG=IOMG*.1 
C         IF(Y1(IOMG-1)<Y1(IOMG)>Y1(IOMG+1)) THEN IB=IOMG
C         PRINT *,IB
C         GO TO 965
C    216      PRINT *,OMEG,IOMG,Y1(IOMG-1),Y1(IOMG),Y1(IOMG+1)
C    965      CONTINUE
        END 
       
            


