function legendre,l,x

if(l gt 150) message,"l is too big -- factorial function will overflow!"

common legendre_cache,l_ok,coeffs

if not isvalid(l_ok) then l_ok=1
if not isvalid(coeffs) then coeffs = [1d0,0d0,0d0]

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; Check if the coeff table is large enough.  If not, then double its
;; size.
while nlm(l_ok) le l+1 do begin
    l_ok = [l_ok,l_ok * 0]
    c1 = coeffs * 0d0
    coeffs = [[coeffs,c1],[c1,c1]]
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; If necessary, calculate the coefficients.
if ( not l_ok(l) ) then begin
    
    ;;;;;;;;;;;;;;;;;;;;
    ; Calculate (x^2 - 1)^l coeffs
    c0  = coeffs(*,l)
    c0(0)=-1
    c0(2)=1
    for i=2,l do begin
        c1            =                  c0          * (-1) ; const
        c1(2:l*2+1+2) =  c1(2:l*2+1+2) + c0(0:l*2+1) * (1)  ; x^2
        c0 = c1
    end

    ;;;;;;;;;;;;;;;;;;;;
    ; Differentiate l times
    ii = indgen(l*2+1)+1
    for i=1,l do begin
        c0(0:l*2) = c0(1:l*2+1)*ii
    end

    ;;;;;;;;;;;;;;;;;;;;
    ; Multiply times the constant...
    c0 = c0 / 2^l / factorial(l)

    ;;;;;;;;;;;;;;;;;;;;
    ; Stuff the values into the coefficient matrix
    coeff(*,l) = c0
    l_ok(l) = 1
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; The coefficients are valid.  Do the multiplication...

out = 0
for i=0,l do out = coeffs(l-i,l) + x*out
return,out

end