zern_radial.pro

IDL语言编写的用于天文自适应光学仿真的软件CAOS V6.0的第一部分。

; $Id: zern_radial.pro,v 1.1.1.1 2002/03/12 11:53:46 riccardi Exp $           
;           
; A. Riccardi, Dipartimento di Astronomia di Firenze (Italy).         
; e-mail address: [e-mail]           
; Please, send me a message if you modify this code.    
    
    
function zern_radial, n, m, r
;+
; NAME:
;       ZERN_RADIAL
;
; PURPOSE:
;       ZERN_RADIAL returns the value of radial part of Zernike polynomial
;       of radial degree n and azimuthal frequency m in abs(r).
;       (We use Noll formulation of Zernike polynomials)
;
; CATEGORY:
;       Special polynomial
;
; CALLING SEQUENCE:
;
;       Result = ZERNKE_RADIAL(N, M, R)
;
; INPUTS:
;       N:      radial degree, integer, N >= 0
;       M:      azimuthal frequency, integer 0 <= M <= N and M-N even
;       R:      piont to evaluate
;
; OUTPUTS:
;       this function returns the value of radial part of Zernike polynomial
;       of radial degree n and azimuthal frequency m in r.
;       On error return 0.
;
; EXAMPLE:
;       Evaluate radial part of Zernike coma in rho. Enter:
;
;           radial = ZERN_RADIAL(3, 1, rho)
;
; MODIFICATION HISTORY:
;       Written by:     A. Riccardi; March, 1995.
;-

	if n lt 0 then begin
		print, 'zern_radial -- must have n>=0'
		return, 0.
	endif

	if m lt 0 or m gt n or not is_even(n-m) then begin
		print, 'zern_radial -- must have 0<=m<=n and n-m even'
		return, 0.
	endif
 
	npm2 = (n+m)/2
	nmm2 = (n-m)/2

	; calculate the coefficient of highest order of the radial polynomial

	; force the same data type of r for a
	
	a=0.*r(0)+1.
	temp=npm2+1  
	for i=1,nmm2 do begin  
		a = a/i*temp  
		temp = temp+1  
	endfor
	result = a*(1. + r*0.)

	; evaluate the polynomial
	for s=0,nmm2-1 do begin
		; calculate the coefficient of lower order
		a = -a/(n-s)*(npm2-s)/(s+1)*(nmm2-s)
		; evaluate the polynomial
		result = result*r*r+a
	endfor

	if (m eq 0) then $
		return, result $
	else $
		return, result*r^m
end