      SUBROUTINE CDFUN(M,EGVL,FNORM,H1,ZR,NDEP,CDF,MX,
     1                 MD,IMAT,BEGVL,BMATL,ALPHA,FAC,RAT,BETA)
C
C     THIS SUBROUTINE COMPUTES THE VALUES OF THE COMPLEX*16 DEPTH FUNCTION
C     CDF AT EACH OF THE NDEP OUTPUT DEPTHS ZR IN THE REGION WITH
C     EIGENVALUES GIVEN BY EGVL
C
      IMPLICIT REAL*8 (A-H,O-Z)
      DIMENSION ZR(MD)
      COMPLEX*16 BEGVL(MX),BMATL(MX,MX)
      COMPLEX*16 EGVL(MX),FNORM(MX),CDF(MX,MD)
      COMPLEX*16 ALPHA(MX), FAC(MX), RAT(MX), BETA(MX)
      COMPLEX*16 CDSRT,ONE,CI
      COMPLEX*16 SUM,ZERO
      COMPLEX*16 CSDS, BDF1
C
      COMMON /BLKEVN/ HB,CW,CB,FKW,FKB,ROHW,ROHB,ATEN
C
      ONE=DCMPLX(1.0,0.0)
      ZERO=DCMPLX(0.0,0.0)
      CI=DCMPLX(0.0,1.0)
      FKW2=FKW**2
      FKB2=FKB**2
C
      IF(IMAT .LT. 1)   THEN
C
       DO 100   I=1,M
       BETA(I)=CDSRT(FKB2*(EGVL(I)+CI*ATEN))
       ALPHA(I)=CDSRT(FKW2+FKB2*(EGVL(I)-ONE))
       FAC(I)=FKB2/(ALPHA(I)*CDSRT(FNORM(I)))
       RAT(I)=FAC(I)*CDSIN(ALPHA(I)*H1)
  100  CONTINUE
C
       DO 300 I=1,M
C
       DO 200 J=1,NDEP
       IF(ZR(J) .LT. H1)   THEN
        CDF(I,J)=FAC(I)*CDSIN(ALPHA(I)*ZR(J))
       ELSE
        CDF(I,J)=RAT(I)*CSDS(ZR(J),H1,BETA(I))
       END IF
  200  CONTINUE
C
  300  CONTINUE
C
      ELSE
C
       DO 420   I=1,M
       ALPHA(I)=CDSRT(FKW2+FKB2*(BEGVL(I)-ONE))
       FAC(I)=FKB2/(ALPHA(I)*CDSRT(FNORM(I)))
       RAT(I)=FAC(I)*CDSIN(ALPHA(I)*H1)
       BETA(I)=CDSRT(FKB2*BEGVL(I))
  420  CONTINUE
C
       DO 900 J=1,NDEP
       IF(ZR(J) .LE. H1)   THEN
C
        DO 600 I=1,M
        SUM=ZERO
        DO 500 II=1,M
        BDF1=FAC(II)*CDSIN(ALPHA(II)*ZR(J))
        SUM=SUM+BMATL(II,I)*BDF1
  500   CONTINUE
        CDF(I,J)=SUM
  600   CONTINUE
C
       ELSE
C
        DO 800 I=1,M
        SUM=ZERO
        DO 700 II=1,M
        BDF1=RAT(II)*CSDS(ZR(J),H1,BETA(II))
        SUM=SUM+BMATL(II,I)*BDF1
  700   CONTINUE
        CDF(I,J)=SUM
  800   CONTINUE
C
       END IF
C
  900  CONTINUE
C
      END IF
C
      RETURN
      END
