I am working on building a very simple optical simulator for my workflow.
I am stuck at a point where I am trying to simulate the impact of diffraction on a lens system that has geometric aberrations.
Textbooks such as Goodman and others specify the wavefront aberrations at the exit pupil essentially add a phase-shift to the perfect spherical wave. This phase-shift term $e^{jW(x, y)}$ where $W(x, y)$ are the Seidel or Zernike polynomials that represents the optical path difference (OPD) between the actual and ideal wavefronts.
I have implemented a detailed ray tracing code that helps me generate spot diagrams. I am even able to calculate optical path length for every ray that traverses through the system. However, I am unclear on how I can (a) compute the actual wavefront OR (b) Determine the function $W(x, y)$ that represents the OPD at the exit pupil.
I came across this post on Zemax communities where they mention something about "subtracting a chief-ray centered reference sphere phase from the optical path lengths computed via ray tracing".
Any help on explaining this or pointing to resources that can help with this will be really appreciated.