C **********************************************************************
C                       Hankel Function 1                              *
C **********************************************************************

C Subroutine FNHANK: This computes Hankel functions of the first kind
C  and of order zero and one.

C Input parameter
C complex Z    : The argument of the Hankel function (Im(Z)=0).

C Output parameter
C H(0:1)       : The Hankel functions of order zero and one respectively.

C External routines
C S17AEF: From the NAG library.
C S17AFF: From the NAG library.
C S17ACF: From the NAG library.
C S17ADF: From the NAG library.

      SUBROUTINE FNHANK(Z,H)
      COMPLEX*16 Z,J0,J1,Y0,Y1,H(0:1)
      REAL*8 X,REH0,REH1,IMH0,IMH1
      REAL*8 S17AEF,S17AFF,S17ACF,S17ADF
      INTEGER IFAIL
      IF (ABS(DIMAG(Z)).GT.1.0D-6) THEN
        WRITE(*,*) 'ERROR(FNHANK) - Im(Z) must be zero'
        STOP
      END IF
      IFAIL=0
      X=DREAL(Z)
      J0=S17AEF(X,IFAIL)
      J1=S17AFF(X,IFAIL)
      Y0=S17ACF(X,IFAIL)
      Y1=S17ADF(X,IFAIL)
      REH0=DREAL(J0)-DIMAG(Y0)
      IMH0=DIMAG(J0)+DREAL(Y0)
      REH1=DREAL(J1)-DIMAG(Y1)
      IMH1=DIMAG(J1)+DREAL(Y1)
      H(0)=CMPLX(REH0,IMH0)
      H(1)=CMPLX(REH1,IMH1)
      END