Represent a wavefront as a sum of Zernike polynomials using a matrix inversion.
This function attempts to solve the a_i's in equation, phi(rho,theta) = SUM(a_i * Z_i(rho,theta), i=1 to M),
where the Z_i(rho,theta)'s are the Zernike polynomials from the zernfun.m file, phi is the wavefront to be represented as a sum of Zernike polynomials, the a_i's are the Zernike coefficients, and M is the number of Zernike polynomials to use.
Input: phi - Phase to be represented as a sum of Zernike polynomials that must be an nXn array (square)
(optional) M - Number of Zernike polynomials to use (Default = 12)
Output: a - Zernike coefficients (a_i's) as a vector Note: zernfun.m is required for use with this file. It is available here: http://www.mathworks.com/matlabcentral/fileexchange/7687
I suggest you make the output a structure array, with each coefficient corresponding to a field name the specifies the aberration.
@Eric Thanks for the suggestions! I'll work array size options in. The Zernike cell, I have kept because it can be useful to have this function return them to reconstruct the wavefront without having to calculate them a second time in another script. Mainly was just for me.
Also, thats a great suggestion for the non-square case, I hadn't really decided what to do about that case, but your suggestion sounds like a good idea.
This is a very good function. Some comments:
1. The cell array "Zernike" seems unnecessary. The for loop on lines 51-54 can be wrapped into the later loop for i=1:M. This saves some memory and a for loop.
2. You should mention somewhere in the documentation that the Zernike polynomials are only calculated on a 256x256 grid, regardless of the size of the input wavefront. I understand this will speed calculations when larger wavefronts are decomposed, but some users may see errors associated with this. It might be worthwhile to have this as an option (i.e., calculate Zernike polynomials over the entire array when accuracy is critical or use the 256x256 option when speed is critical).
3. You might modify the function to more elegantly handle non-square wavefront data. You could pad the smaller dimensions with zeros to make it square, for instance. Or you could solve the equation using a different technique that does not impose this constraint. You may need to force the user to supply (x,y) position data for that to make sense.