c # yennie.f         				Lu, 6/24/1991
c # subroutine using yennie's trick to sum point coulomb phase.
c # It is powerful to sum series of leg pol. 
	subroutine yennie (nth,k0,xgam,scoul,ar,ai,index,
     +ldmax,xcos,zcul)
	implicit real*8 (a-h,o-y)
	implicit complex*16 (z)
	real*8 k0
	dimension pl(0:99),ar(0:99),ai(0:99),br(0:99),bi(0:99)
	dimension scoul(100)
c # if index = 0, use input coef; 
c              1, for pure coul;
c 	     666, return analytic coulomb amp.
c
	hbarc = 197.3289d0
	pi = 3.141592654d0
cp0 = 1030.13340377d0
cgama = 0.53175388408238d-1
c
	lmax = ldmax - 1
	if (index .eq. 666) then
		sinhalf2 = .5 * (1. - xcos)
		fcoul = - xgam /2./k0/sinhalf2
		phc = 2. * scoul(1) - xgam * dlog(sinhalf2)
		fr = fcoul * dcos(phc)
		fi = fcoul * dsin(phc)
		zcul = cmplx(fr,fi)
		return
	else if ( index .eq. 1 ) then
		do j = 0, lmax
	  		if (nth .eq. 0) then
	   			del = scoul(j+1) - scoul(lmax+1)
	  		else if (nth .eq. 1) then
	   			del = scoul(j+1) 
	  		endif
           	ar(j) = dcos(del) * dsin(del) * (2.*j + 1.) 
           	ai(j) = dsin(del)**2 * (2.*j + 1.)
		enddo 
	endif
c
c ---- transform of coeff ---
	call trans (nth, ar, lmax, br)
	call trans (nth, ai, lmax, bi)
c
	leff = lmax - nth
	  call legpol (xcos,pl,lmax)
	  fr = 0.
	  fi = 0.
	     if( nth .eq. 0 ) then
		  do j = 0, leff
		    fr = fr + ar(j) * pl(j)
		    fi = fi + ai(j) * pl(j)
		  enddo
	     else
	         do j = 0, leff
	            fr = fr + br(j) * pl(j)
 	            fi = fi + bi(j) * pl(j)
	         enddo
		factor = 1./( 1. - xcos )**nth
	        xr =  factor * fr/k0
		xi =  factor * fi/k0
	 	zcul = CMPLX(xr, xi)
	     endif
	return
	end
c # Lu: 6/26/1991, make mtime reductions by yennie's trick
	subroutine trans (mtime, a, lmax, b)
	implicit real*8 (a-h, o-z)
	dimension a(0:99), b(0:99)
c
	m = 1
	if (mtime .eq. 0) return
	la = lmax
	lb = lmax - 1
c
101     call trans0 (a, la, b, lb)
	if (m .eq. mtime) return
	   do l = 0, lb
	     a(l) = b(l)
	   enddo
           m = m + 1
	   la = lb
	   lb = la - 1
	goto 101
	end
c # rearrangement of coeff of leg poly
	subroutine trans0 (a, la, b, lb)
	implicit real*8 (a-h, o-z)
	dimension a(0:99), b(0:99)
c
	b(0) = a(0) - a(1) * 1.d0/3.d0
	do j = 1, lb
	  b(j) = a(j) - a(j+1)*(j+1.)/(2.*j+3.) - a(j-1)*j/(2.*j-1.)
	enddo
	return
	end