ZernikePolynomials.jl
This package provides functionality to work with Zernike polynomials.
It features the following:
- functions to convert between common serial indices of Zernike polynomials
- functions to evaluate the polynomials on a grid
- functions to estimate Zernike coefficients in a least squares sense from a user-provided input.
Conversion to and between (Noll & OSA/ANSI) sequential indices
This package provides conversion utility functions between three common ways of specifying a Zernike polynomial.
The first is a tuple (m,n) of integers to specify the polynomial Z_n^m(ρ,θ), where n ≥ m, and (ρ,θ) are polar coordinates (radius and angle).
The second (sequential) index is the OSA/ANSI standard index (0,1,2,3,...).
The third is Noll's sequential index (1,2,3,...).
An invalid combination (m,n) returns NaN in conversion.
>> using ZernikePolynomials
>> [mn2OSA(n,m) for n in 0:4, m in -5:5]
>> [mn2Noll(n,m) for n in 0:4, m in -5:5]
>> [OSA2mn(j) for j in 0:5]
>> [Noll2mn(j) for j in 0:5]
>> [Noll2OSA(OSA2Noll(j)) for j = 0:5]
Evaluating Zernike polynomials
The Zernike polynomials (ρ,θ) -> Z_n^m(ρ,θ) or (x,y) -> Z_n^m(x,y) can be obtained as a function as follows:
>> Zernike(1,1)
>> Zernike(1,1,coord=:polar)
>> Z = Zernike(5,index=:Noll,coord=:cartesian) # 5th polynomial by Noll's numbering
>> Z = Zernike(1,1,coord=:cartesian)
>> Z(0.5,0.2) # evaluate the function at cartesian coordinate (0.5,0.2)
Polar coordinates are used by default. The functions return 0 for ρ > 1.
Zernike polynomials and affine combination thereof can easily be evaluated on a grid of points using the function evaluateZernike()
.
For example, to evaluate 0.5Z_2 + 0.3Z_3 (by OSA numbering) on a 256x256 grid, use the following
>> evaluateZernike(256, [2, 3], [0.5, 0.3], index=:OSA)
To evaluate it on a (square) grid with predefined coordinates, use for example
>> x = LinRange(-2,2,256)
>> ϕ = evaluateZernike(x, [2, 3], [0.5, 0.3], index=:OSA)
>> using Plots
>> heatmap(ϕ)
Here the range x
is a range that gives the x and y coordinates.
Zernike polynomial normalization
The Zernike polynomials are normalized according to Thibos et al. - "Standards for Reporting the Optical Aberrations of Eyes", i.e.
N_n^m = sqrt(2 (n+1) / (1 + δ(m,0)))
where δ(m,0) = 1 for m = 0 and 0 otherwise.
These normalization constants can be obtained by normalization(m,n)
:
>> [normalization(OSA2mn(i)...) for i in 0:5]
Note that the definition on the Wikipedia page on Zernike polynomials is different. Here the polynomials are normalized between [-1,1]. When using or reporting (estimated) Zernike coefficients, it is important to be aware of which normalization has been used.
Estimating Zernike coefficients in a least squares sense
A common use for Zernike polynomials is to approximate a given 2D input. Since Zernike polynomials are often used in optics to approximate a 2D phase, this is the term used in the function documentation.
The function Zernikecoefficients()
estimates (in a least-squares sense) the optimal coefficients for a sum of Zernike polynomials. A vector of (sequential) indices should be provided to specify which coefficients should be estimated
>> ϕ = evaluateZernike(256, [2, 3], [0.5, 0.3], index=:OSA)
>> Zernikecoefficients(ϕ, [2, 3]) # ≈ [0.5, 0.3]
It is also possible to specify on which coordinates the phase is defined:
>> x = LinRange(-2,2,256)
>> ϕ = evaluateZernike(x, [2, 3], [0.5, 0.3], index=:Noll)
>> Zernikecoefficients(x, ϕ, [1, 3], index=:Noll) # ≈ [0.0, 0.3]
Note that in this last example the coefficient for the Zernike polynomial with Noll index 1 is estimated, but this is not present in ϕ. Therefore the estimated coefficient is 0.
Plotting example
using Plots
using ZernikePolynomials
x = LinRange(-1,1,256)
ϕ = evaluateZernike(x, [2, 3], [0.5, 0.3], index=:Noll)
f(x,y) = (x^2 + y^2) <= 1 ? 1. : NaN # NaNs are not plotted
heatmap( ϕ .* [f(X,Y) for X in x, Y in x])