Background & Theory
To determine PSF we must know a) the electric field, b) the aperture, and c) the observation point. First, let us begin with the ‘point’ in the point spread function. As an optical scientist, we care to investigate how a point is transferred through our optical system. Thus, we begin by selecting a point of interest in our field. This point will radiate rays outward, with the output rays being normal to a sphere surrounding the point (i.e. a spherical wavefront). Not all of these rays enter the optical system of course; all we care about is selecting for the rays that enter the entrance pupil of our system, which is determined by tracing the rays from the point source to the entrance pupil. From the entrance pupil, the rays must then be propagated through the system, using geometrical ray propagation. Unlike a traditional design process however, we only trace the rays to the conjugate of the entrance pupil, known as the exit pupil. This makes sense, from an ideal, infinitely small (i.e. only happens in textbooks and software) point source, way determine where rays enter the system in the entrance pupil. We then propagate the rays through the system and re-evaluate them as the exit of the system, the exit pupil.
This is where the system gets slightly more complex (and more interesting!). Thus far, we treated the light from our point source as geometrical rays. This is fine for first order approximations, however, we know light is more complex than just rays. This is where we begin to consider the wave nature of light from our point source. With a wave description of light, the light propagates as a wave and also can interfere with other waves of light (or even itself). At the exit pupil, we consider the propagated rays of light as sampling points along our electric field wavefront. We can treat the electric field at the exit pupil, which we will refer to as Ui, as a collection of point sources (Hyugen wavelets) [3]. These Hyugen wavelets are treated as infinitely small point sources which radiate spherical electromagnetic waves outward from the point source. The electric field component of the wavelets retains the same frequency as the source wavefront, Ui. As the wavelets propagate, they interfere with one another. This is known as diffraction; for a robust discussion on diffraction we recommend Dr. Tom Milster’s lectures, or Dr. Mahajan’s discussion on Wave Diffraction Optics /[1,2/].
For comparison, a ‘perfect’ system without aberrations with form an Airy disk, which is a 2D Bessel function. Thus, we see that while if we had retained our more naive picture of purely geometrical rays, we may be led to think a perfect system would transfer an idealized point source into an idealized point in image space. In reality, due to the diffraction that occurs, the best transfer we can hope for is an Airy disk. The ideal Airy disk diameter, and thus the best resolution that can possibly be obtained for an optical system via geometrical optics is defined as 2.44 \lambda F/#. Of course, most of the time we aren’t so lucky as to design a perfect system. Instead, when we look at our PSF, we may notice the differences it has from an ideal Airy disk, and what this implies for the performance. Key questions to ask yourself may include the following: Is the PSF symmetric? What is the width of the PSF, and what implications does this have on the system resolution? Does the PSF change across the field, and if so how?
Algorithm
For the interested scientist and user, it is essential to know what exactly the software is doing to determine the PSF. The following section lays out the algorithms used:
Step 1: User selects a field coordinate for the PSF process
This takes in the field coordinate and determines the x and y position at the object plane of the point source.
Step 2: Ray transfer
The rays are geometrically propagated through the optical system using the following algorithms
First, find intercept between ray and object:
For a spherical surface we use Line Sphere Intersection:
Where \vec{l} is the ray vector, o is the origin of the ray, c is the sphere center, r is the radius of curvature of the sphere, and d is the distance from the ray origin, along it’s pointing vector, to the intercept. Thus, we solve for d by using the quadratic formula.
For a plane surface we use Line Plane Intersection:
Where p_0 is a point on the plane, o is the origin of the ray, \vec{n} is the normal vector of the plane, \vec{l} is the ray vector, and d is the distance from the ray origin, along it’s pointing vector, to the intercept. We solve for d directly in this case.
Step 3: Ray Deviation
After the intercept at the next surface is determined, the ray is deviated by solving the vector form of Snell’s Law. The following algorithm is utilized:
Given the incident ray vector \vec{l}_{inc}, and the normal of the incident point, \vec{n}, we can define the cosine of the angle of incidence \theta _{inc} as
\cos \theta _{inc} must be positive, if it is not then \vec{n} sign is changed to positive in eq. 3.
Next, we determine the cosine of the angle of transmission, \cos \theta _{trans}, as:
if \cos \theta _{trans} is less than zero, total internal reflection occurred. The Fresnel coefficients are next calculated as:
Where R_s is the reflectance for s-polarized light, while R_p is the reflectance for p-polarized light, and R is the effective reflectivity.
We approximate the interaction at this stage, and state that if R is greater than 0.5, the ray reflects entirely, with the power scaled by the effective reflectivity, otherwise, it transmits, with the power scaled by 1 - R. The direction for the reflected ray, \vec{l} _{refl}, is given as:
While the direction for the refracted ray, \vec{l} _{trans}, is given as:
Step 4: Exit Pupil Ray Phase
The rays from Step 1 are traced through the system to find the image rays in the image plane, following the algorithms in Steps 2-4 (the ray intersection is determined, then the ray deviation is calculated, through surfaces).
The rays are then directly traced from the image to the Exit Pupil (XP) plane, using equation (2). The optical path difference (OPD) from the real rays in the exit pupil as compared to the ideal reference sphere is determined. The ideal reference sphere is defined as a sphere whose origin as located at the paraxial image plane, and whose radius of curvature is defined as the distance from the exit pupil plane to the paraxial image plane. The OPD algorithm is :
Where n is the refractive index the ray is in at the exit pupil (image space refractive index), OPL_{real} is the optical path length of the ray from the image plane to the exit pupil plane, and OPL_{ideal} is the optical path length from the paraxial image point in the image plane to the exit pupil plane.
We then determine the OPD phase across the exit pupil for all rays as:
Where \phi_{ray} is the phase difference for each ray, whose sampling occurs at every ray intercept point at the exit pupil, between the ideal wavefront for a paraxial image and our actual wavefront that forms our real image of the traced point source object.
Finally, a check is performed to determine if the points in the exit pupil where the phase was determined actually pass through the exit pupil or if they are obscured, either by missing the pupil or hitting a physical mask. If they do not pass, the phase difference is defined as zero. Otherwise, the real component of the phase is scaled to correctly define irradiance. The scaling is defined as:
Where C_{irr} is the irradiance scaling factor, c is the speed of light in a vacuum, n is the refractive index in image space, and \epsilon _0 is the permittivity of the vacuum. Then, the real component of the each rays optical path difference phase is defined as:
Where \vec{S} is the Poynting vector of the ray. The collection of phase differences define a field, termed A(x,y), which complex monochromatic field that has illuminated the exit pupil and diffracts to the image plane. We are making the assumption currently that the diffraction can be approximated as Fraunhofer diffraction, although, more options will be added in due time.
Step 5: Field Propagation
Following the determination of the phase difference field in Step 4, the field is propagated to the image plane via Fraunhofer propagation. The field in the image plane is thus determined as:
Where U(x',y',z') is the propagated field in the image plane, defined at points x', y', and z' = z_{img}. The integration is performed over the limits of the exit pupil, termed XP, A(x,y) is the field at the exit pupil from step 4, `i` is the imaginary unit, defined as \sqrt{-1}, \lambda is the wavelength of the light, and l and m are the direction cosines of the point (x',y') with respect to the origin of the exit pupil.
Finally, the irradiance of our transferred point source object point in the image plane is given as:
If you have any questions regarding these algorithms, or find any errors in our algorithms, please let us know by contacting us at mail@eleoptics.com
References:
[1] Milster, Tom. “OPTI 505R - Interference and Diffraction” Lecture, The University of Arizona, Tucson, AZ, 2015.
[2] https://en.wikipedia.org/wiki/Fraunhofer_diffraction
[3] C. Huygens, “Treatise on light,” rendered into English by S. P. Thompson, Dover Publications, Inc., New York, 1962.
[4] Mahajan, Virendra. “Chapters 2 - Image Systems with Circular Pupils” in Optical Imaging and Aberrations, Part II. Wave Diffraction Optics, Second Edition : SPIE Press, 2011
[5] Smith, Warren J. Modern Optical Engineering : The Design of Optical Systems. 3rd ed. New York: McGraw Hill, 2000. Print. Optical and Electro-optical Engineering Ser.
[6] Smith, Warren J. Modern Lens Design. 2nd ed. New York: McGraw-Hill, 2005. Print. McGraw-Hill Professional Engineering. Modern Lens Design.
[7] Zimmer, Hans-Georg. Geometrical Optics. Berlin, Heidelberg: Springer Berlin Heidelberg, 1970. Print. Applied Physics and Engineering, An International Ser., 9.