: [Z, dZx, dZy] = zernikes_and_derivarives_cartesian_OSA (x, y, n)

: [Z, dZx, dZy] = zernikes_and_derivarives_cartesian_OSA (x, y, n, nan_zero)

Return the cartesian Zernike’s pollynomials and its partial derivatives up to radial degree n, i.e. until Z[n,n]

x is a matrix with the X coordinates of the points where the Zernike’s polynomials and its derivatives are computed. y is a matrix with the Y coordinates of the same points. n is an integer with the maximum radial degree desired. nan_zero is a string that determines the values of polynomial and derived values outside the radio unit circle.

Strictly, the polynoms are only defined for 0 <= X²+Y² <= 1. If variable nan_zero = ’nan’, the values of the polynomials for which it is verified that (X²+Y²)>1 are set = NaN. If variable nan_zero = ’zero’, the values of the polynomials for which it is verified that (X²+Y²)>1 are set = 0.

Z is a 3D matrix. Each page contains a 2D matrix with the values of a Zernike’s polynomial at the points defined by x and y.

dZx is a 3D matrix. Each page contains the values of the partial derivative in x.

dZy is a 3D matrix. Each page contains the values of the partial derivative in y.

It should be noted that in standard OSA/ANSI the simple-index j starts at 0, but in octave the indices of the vectors and matrices start at 1. So that page 1 of the 3D Z, dZx and dZy matrices corresponds to the single-index j = 0, and therefore to the double-index m = 0 and n = 0. Page 2 corresponds to j = 1, page 3 –> j = 2, etc.

Example

x = linspace(-1,1,101);
[X,Y] = meshgrid(x,x);
[Z,dZx,dZy] = zernikes_and_derivatives_cartesian_OSA (X,Y,7,'zero');
Z_00 = Z(:,:,1);
   # Z_00 is a 2D matrix with the values of Zernike's polynomial
   # with simple-index j = 0, and double-index m = 0 & n = 0. 
dZx_-24 = dZx(:,:,11);
   # Z_-44 is a 2D matrix with the values of the partial 
   # derivative in x of Zernike's polynomial with 
   # simple-index j = 10, and double-index m = -4 & n = 4.

Run the demo to see a more complete example.

Size of x must be equal size of y.

References:

1. Andersen T.B., "Efficient and robust recurrence relations for the Zernike circle polynomials and their derivatives in Cartesian coordinates". Optic Express 26(15), 18878-18896 (2018).
2. Thibos, L.N, Applegate, R.A., Schwiegerling, J.T. & Webb, R., Standards for reporting the optical aberrations of eyes. Journal of refractive surgery, 18(5), S652-S660 (2002).

See also: zernike_osa_ansi_to_nm, zernike_cartesian, zernike_name, zernike_polar, zernike_R_poly.

Demonstration 1

 Nmgr = 7;                        # Maximum radial degree number, default = 7
 Nmodos = ((Nmgr+3)*Nmgr/2)+1;    # Mode number, default = 36

 pr = 20;                         # Pupil radius, mm
 # Let us suppose that we obtain the heights of a biconical surface
 # for a series of points in the xy plane
 Rx = 0.002;			          # mm⁻¹
 px = 0.8;
 Ry = 0.001;			          # mm⁻¹
 py = 0.5;
 x = pr*(2*rand(200,1) - 1);	     # mm
 y = pr*(2*rand(200,1) - 1);	     # mm
 surface = (Rx*x.^2 + Ry*y.^2)/(1+sqrt(1-(px*(Rx*x).^2+py*(Ry*y).^2)));
 sag = max(max(surface));
 surface = sag - surface;         # surface height, mm

 # Coordinates normalization
 xnorm=x/pr;
 ynorm=y/pr;
 usefuldata = find(sqrt(x.^2+y.^2)
Produces the following figure


Figure 1

Package: optics