| Re: LeGendre Polynomials? [message #3742] |
Sun, 26 February 1995 11:51 |
thompson
Messages: 584
Registered: August 1991
|
Senior Member |
|
|
chris@mercury.sfsu.edu (Christopher McCarthy) writes:
> To My knowledge IDL does not do Legendre polynomials--an obvious
> oversite if true, given that it does Bessel Functions (except
> beselK!), etc.
> Has anyone written a Legendre Polynomial program?
> Thanks,
> Chris
Try this. Hope it meets your needs.
Bill Thompson
============================================================ ===================
FUNCTION LEGENDRE,X,L,M
;
IF N_PARAMS(0) LT 2 THEN BEGIN
PRINT,'*** LEGENDRE must be called with 2-3 parameters:'
PRINT,' X, L [, M ]'
RETURN,X
END ELSE IF N_PARAMS(0) EQ 2 THEN M = 0
;
IF M LT 0 THEN BEGIN
PRINT,'*** M must not be less than 0, routine LEGENDRE.'
RETURN,X
END ELSE IF M GT L THEN BEGIN
PRINT,'*** M must not be greater than L, routine LEGENDRE.'
RETURN,X
END ELSE BEGIN
S = SIZE(X)
IF S(0) EQ 0 THEN TEST = ABS(X) ELSE TEST = MAX(ABS(X))
IF TEST GT 1 THEN BEGIN
PRINT,'*** X must be in the range -1 to 1, routine LEGENDRE.'
RETURN,X
ENDIF
ENDELSE
;
PMM = 0.*X + 1.
IF M GT 0 THEN BEGIN
SOMX2 = SQRT( (1. - X) * (1. + X) )
FACT = 1.
FOR I = 1,M DO BEGIN
PMM = -PMM*FACT*SOMX2
FACT = FACT + 2.
ENDFOR
ENDIF
;
IF L EQ M THEN PLGNDR = PMM ELSE BEGIN
PMMP1 = X * (2*M + 1) * PMM
IF L EQ M+1 THEN PLGNDR = PMMP1 ELSE BEGIN
FOR LL = M+2,L DO BEGIN
PLL = (X*(2*LL-1)*PMMP1 - (LL+M-1)*PMM) / (LL-M)
PMM = PMMP1
PMMP1 = PLL
ENDFOR
PLGNDR = PLL
ENDELSE
ENDELSE
;
RETURN,PLGNDR
END
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
