C ELPDBE    SOURCE    CHAT      05/01/12    23:36:17     5004           
      SUBROUTINE ELPDBE(CZ,N,NP,CH10,CH11,XM)
C     
C =====================================================================
C     
C     CALCUL DES FONCTIONS DE HANKEL PAR SERIE
C     
C     
C =====================================================================
C     
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-B,D-H,O-Z)
      IMPLICIT COMPLEX*16(C)
C    Include contenant quelques constantes dont XPI :
-INC CCREEL
C     
C     
      DIMENSION CZ(*)
      DIMENSION CH10(*)
      DIMENSION CH11(*)
C     
      CI = (0d0,1d0)
      CXPI = XPI
C     Constante d'Euler
      CGAM = .5772156649015328606d0
C     
C     
      XERR = 0D0
      DO  1000  I=1,N

         CJ0 = 0d0
         CJ1 = 0d0

C     

         IF ( ABS ( CZ(I)) .LT. XM) THEN

            CT = CZ(I) **2 / 4.d0
            CAL0 = 1d0
            CAL1 = CZ(I) /2.d0

            DO  1100  IK=1,NP
               K= IK - 1
               CJ0 = CJ0 + CAL0
               CAL0 = -1.d0*CT* CAL0/ ( ( K+1) **2)
               CJ1    = CJ1    + CAL1
               CAL1 = -1.d0*CT* CAL1/ ( ( K+1) *(K + 2))
 1100       CONTINUE
C     
            XERR = MAX (XERR,ABS(CAL0/CJ0))
            XERR = MAX (XERR,ABS(CAL1/CJ1))

            CT = CZ(I) **2 / 4.d0
            CY0 = LOG ( CZ(I) / 2.d0)
            CY0 = CJ0  * ( CY0  + CGAM)
            CS0   = 1d0
            CS0P  = 0d0
            CBET0 = CT

            CY1 = LOG ( CZ(I) / 2.d0)
            CY1 = CJ1 *  CY1  - 1.d0/CZ(I)
            CS1   = -2.d0*CGAM + 1.d0
            CS1P  = 0d0
            CBET1 =-1.d0*CZ(I) *CS1/4.d0


            DO  1110  IK=1,NP

               CY0 = CY0 + CBET0
               CS0P =  CS0 + 1.d0/( IK + 1)
               CBET0= -1.d0*CT* CBET0/ ( (IK+1) *(IK + 1)) * ( CS0P/CS0)
               CS0  =  CS0P

               K= IK - 1
               CY1 = CY1 + CBET1
               CS1P =  CS1 +( 1.D0/( K + 1)) + (1.D0/( K + 2 ))
               CBET1= -1.d0*CT* CBET1/ ( ( K+1) *( K+ 2)) * ( CS1P/CS1)
               CS1  =  CS1P

 1110       CONTINUE
            XERR = MAX (XERR,ABS(CBET0/CY0))
            XERR = MAX (XERR,ABS(CBET1/CY1))
            CY0 = CY0 *2.d0 / CXPI
            CY1 = CY1 *2.D0 / CXPI
            CH10(I) = CJ0 + CI*CY0
            CH11(I) = CJ1 + CI*CY1
         ELSE
            C8 = CZ(I) *8d0
            CAL0 = 1d0
            CAL1 = 1d0
            CP0= CAL0
            CP1= CAL1
            CQ0=0.D0
            CQ1=0.D0

            DO  1200  K=1,NP
               I1= 2*K -1
               I2= 2*K
               CBET0= -1.D0*CAL0*(2*I1 - 1)**2/ (I1 *C8)
               CBET1=     CAL1*(4 - (2*I1 - 1)**2)/ (I1 *C8)
               CAL0=     CBET0*(2*I2 - 1)**2/ (I2 *C8)
               CAL1= -1.D0*CBET1*(4 - (2*I2 - 1)**2)/ (I2 *C8)
               CP0    = CP0 +  CAL0
               CQ0    = CQ0 + CBET0
               CP1    = CP1 +  CAL1
               CQ1    = CQ1 + CBET1
 1200       CONTINUE
            XERR = MAX (XERR,ABS(CBET0/CQ0))
            XERR = MAX (XERR,ABS(CBET1/CQ1))
            XERR = MAX (XERR,ABS(CAL0/CP0))
            XERR = MAX (XERR,ABS(CAL1/CP1))

            CCOEF=(2 / ( CXPI * CZ(I))) ** .5D0
            CKI0= CZ(I) - CXPI /4.d0
            CH10(I)=CCOEF* (CP0 + CI * CQ0) * EXP ( CI * CKI0)
            CKI1= CZ(I) - (CXPI*3.d0 /4.d0)
            CH11(I)=CCOEF* (CP1 + CI * CQ1) * EXP ( CI * CKI1)
         ENDIF

 1000 CONTINUE

C     
      RETURN
      END