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.
This package provides conversion utility functions between three common ways of specifying a Zernike polynomial.
Numbering conventions are implemented as subtypes of the abstract type ZernikeIndex.
The first is to use NM(n, m) to specify the polynomial Zₙᵐ(ρ,θ), where n and m are integers with n ≥ abs(m) (further, n - abs(m) must be even), and (ρ,θ) are polar coordinates (radius and angle).
The second (sequential) index is the OSA/ANSI standard index (0,1,2,3,...), specified as OSA(j).
The third is Noll's sequential index (1,2,3,...), specified as Noll(j).
Invalid inputs to the ZernikeIndex constructors result in an ArgumentError, which ensures that you can only construct valid indexes. You can convert between types by calling the constructor or using convert:
julia> using ZernikePolynomials
julia> nm = NM(3, -3)
NM(3, -3)
julia> OSA(nm)
OSA(6)
julia> Noll(nm)
Noll(9)
julia> NM(OSA(6))
NM(3, -3)
julia> [NM(n, m) for n in 0:4, m in -5:5]
ERROR: ArgumentError: Invalid Zernike index pair (n,m)=(0,-5).
...
julia> Noll(0)
ERROR: ArgumentError: Invalid Noll index 0.
...
julia> supertype(NM)
ZernikeIndexThe Zernike polynomials (ρ,θ) -> Zₙᵐ(ρ,θ) or (x,y) -> Zₙᵐ(x,y) can be obtained as a function as follows:
>> zernike(NM(1,1))
>> zernike(NM(1,1),coord=:polar)
>> Z = zernike(Noll(5),coord=:cartesian) # 5th polynomial by Noll's numbering
>> Z = zernike(NM(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, OSA.([2, 3]), [0.5, 0.3])To evaluate it on a (square) grid with predefined coordinates, use for example
>> x = LinRange(-2,2,256)
>> ϕ = evaluatezernike(x, OSA.([2, 3]), [0.5, 0.3])
>> using Plots
>> heatmap(ϕ)Here the range x is a range that gives the x and y coordinates.
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 normalization(m,n):
>> [normalization(NM(OSA(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.
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, OSA.([2, 3]), [0.5, 0.3])
>> zernikecoefficients(ϕ, OSA.([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, Noll.([2, 3]), [0.5, 0.3])
>> zernikecoefficients(x, ϕ, Noll.([1, 3])) # ≈ [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.
using Plots
using ZernikePolynomials
x = LinRange(-1,1,256)
ϕ = evaluatezernike(x, Noll.([2, 3]), [0.5, 0.3])
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])