FresnelWavefront¶
-
class
poppy.
FresnelWavefront
(beam_radius, units=Unit("m"), rayleigh_factor=2.0, oversample=2, **kwargs)[source]¶ Bases:
poppy.poppy_core.BaseWavefront
Wavefront for Fresnel diffraction calculation.
This class inherits from and extends the Fraunhofer-domain poppy.Wavefront class.
- beam_radius : astropy.Quantity of type length
- Radius of the illuminated beam at the initial optical plane. I.e. this would be the pupil aperture radius in an entrance pupil.
- units : astropy.units.Unit
- Astropy units of input parameters
- rayleigh_factor:
- Threshold for considering a wave spherical.
- oversample : float
- Padding factor to apply to the wavefront array, multiplying on top of the beam radius.
- Lawrence, G. N. (1992), Optical Modeling, in Applied Optics and Optical Engineering., vol. XI,
- edited by R. R. Shannon and J. C. Wyant., Academic Press, New York.
- https://en.wikipedia.org/wiki/Gaussian_beam
- IDEX Optics and Photonics(n.d.), Gaussian Beam Optics,
- [online] Available from: https://marketplace.idexop.com/store/SupportDocuments/All_About_Gaussian_Beam_OpticsWEB.pdf
- Krist, J. E. (2007), PROPER: an optical propagation library for IDL,
- vol. 6675, p. 66750P-66750P-9. [online] Available from: http://dx.doi.org/10.1117/12.731179
- Andersen, T., and A. Enmark (2011), Integrated Modeling of Telescopes, Springer Science & Business Media.
Attributes Summary
angular_coordinates
Should coordinates be expressed in arcseconds instead of meters at the current plane? divergence
Half-angle divergence of the gaussian beam fov
FOV in arcseconds, if applicable param_str
Formatted string of gaussian beam parameters. pixelscale
Pixelscale, in meters by default or in arcseconds if angular_coordinates is True waists
each [z_w_0,w_0] for each waist generated by an optic z_r
Rayleigh distance for the gaussian beam, based on current beam waist and wavelength. Methods Summary
apply_lens_power
(optic[, ignore_wavefront])Adds lens wavefront curvature to the wavefront corresponding to the lens’ focal length f_l, and updates the Gaussian beam parameters of the wavefront. coordinates
()Return Y, X coordinates for this wavefront, in the manner of numpy.indices() display
(*args, **kwargs)Display wavefront on screen from_wavefront
(wavefront)Convert a Fraunhofer type wavefront object to a Fresnel one planar_range
(z)Returns True if the input range z is within the Rayleigh range of the waist. propagate_direct
(z)Implements the direct propagation algorithm as described in Andersen & Enmark (2011). propagate_fresnel
(delta_z[, display_intermed])Top-level routine for Fresnel diffraction propagation propagate_to
(optic, distance)Propagates a wavefront object to the next optic in the list, after some separation distance (which might be zero). pupil_coordinates
(x, y, pixelscale)Utility function to generate coordinates arrays for a pupil plane wavefront r_c
([z])The gaussian beam radius of curvature as a function of distance z spot_radius
([z])radius of a propagating gaussian wavefront, at a distance z Attributes Documentation
-
angular_coordinates
= False¶ Should coordinates be expressed in arcseconds instead of meters at the current plane?
-
divergence
¶ Half-angle divergence of the gaussian beam
I.e. the angle between the optical axis and the beam radius (at a large distance from the waist) in radians.
-
fov
¶ FOV in arcseconds, if applicable
-
param_str
¶ Formatted string of gaussian beam parameters.
-
pixelscale
¶ Pixelscale, in meters by default or in arcseconds if angular_coordinates is True
-
waists
¶ each [z_w_0,w_0] for each waist generated by an optic
-
z_r
¶ Rayleigh distance for the gaussian beam, based on current beam waist and wavelength.
I.e. the distance along the propagation direction from the beam waist at which the area of the cross section has doubled. The depth of focus is conventionally twice this distance.
Methods Documentation
-
apply_lens_power
(optic, ignore_wavefront=False)[source]¶ Adds lens wavefront curvature to the wavefront corresponding to the lens’ focal length f_l, and updates the Gaussian beam parameters of the wavefront.
- optic : QuadraticLens
- An optic
- ignore_wavefront : boolean
- If True then only gaussian beam propagation parameters will be updated and the wavefront surface will not be calculated. Useful for quick calculations of gaussian laser beams
-
coordinates
()[source]¶ Return Y, X coordinates for this wavefront, in the manner of numpy.indices()
This function knows about the offset resulting from FFTs. Use it whenever computing anything measured in wavefront coordinates.
The behavior for Fresnel wavefronts is slightly different from Fraunhofer wavefronts, in that the optical axis is not the exact center of an array (the corner between pixels for an even number of pixels), but rather is a specific pixel (e.g. pixel 512,512 for a 1024x1024 array). This is for consistency with the array indexing convention used in FFTs since this class depends on FFTs rather than the more flexible matrix DFTs for its propagation.
For Fresnel wavefronts, this depends on the focal length to get the image scale right.
- Y, X : array_like
- Wavefront coordinates in either meters or arcseconds for pupil and image, respectively
-
display
(*args, **kwargs)[source]¶ Display wavefront on screen
- what : string
- What to display. Must be one of {intensity, phase, best}. ‘Best’ implies to display the phase if there is nonzero OPD, or else display the intensity for a perfect pupil.
- nrows : int
- Number of rows to display in current figure (used for showing steps in a calculation)
- row : int
- Which row to display this one in? If set to None, use the wavefront’s self.current_plane_index
- vmin, vmax : floats
- min and maximum values to display. When left unspecified, these default to [0, intens.max()] for linear (scale=’linear’) intensity plots, [1e-6*intens.max(), intens.max()] for logarithmic (scale=’log’) intensity plots, and [-0.25, 0.25] waves for phase plots.
- scale : string
- ‘log’ or ‘linear’, to define the desired display scale type for intensity. Default is log for image planes, linear otherwise.
- imagecrop : float, optional
- Crop the displayed image to a smaller region than the full array. For image planes in angular coordinates, this is given in units of arcseconds. The default is 5, so only the innermost 5x5 arcsecond region will be shown. This default may be changed in the POPPY config file. If the image size is < 5 arcsec then the entire image is displayed. For planes in linear physical coordinates such as pupils, this is given in units of meters, and the default is no cropping (i.e. the entire array will be displayed unless this keyword is set explicitly).
- showpadding : bool, optional
- For wavefronts that have been padded with zeros for oversampling, show the entire padded arrays, or just the good parts? Default is False, to show just the central region of interest.
- colorbar : bool
- Display colorbar
- crosshairs : bool
- Display a crosshairs indicator showing the axes centered on (0,0)
- ax : matplotlib Axes, optional
- axes to display into. If not set, will create new axes.
- use_angular_coordinates : bool, optional
- Should the axes be labeled in angular units of arcseconds? This is used by FresnelWavefront, where non-angular coordinates are possible everywhere. When using Fraunhofer propagation, this should be left as None so that the coordinates are inferred from the planetype attribute. (Default: None, infer coordinates from planetype)
- figure : matplotlib figure
- The current figure is modified.
-
classmethod
from_wavefront
(wavefront)[source]¶ Convert a Fraunhofer type wavefront object to a Fresnel one
Note, for now this function only works if the input wavefront is at a pupil plane, so the Fraunhofer wavefront has pixelscale in meters/pix rather than arcsec/pix. Conversion from image planes may be added later.
- wavefront : Wavefront
- The (Fraunhofer-type) wavefront to be converted
-
planar_range
(z)[source]¶ Returns True if the input range z is within the Rayleigh range of the waist.
- z : float
- distance from the beam waist
-
propagate_direct
(z)[source]¶ Implements the direct propagation algorithm as described in Andersen & Enmark (2011). Works best for far field propagation. Not part of the Gaussian beam propagation method.
- z : float or Astropy.Quantity length
- the distance from the current location to propagate the beam.
-
propagate_fresnel
(delta_z, display_intermed=False)[source]¶ Top-level routine for Fresnel diffraction propagation
- Each spherical wavefront is propagated to a waist and then to the next appropriate plane
- (spherical or planar).
- delta_z : float
- the distance from the current location to propagate the beam.
- display_intermed : boolean
- If True, display the complex start, intermediates waist and end surfaces.
-
propagate_to
(optic, distance)[source]¶ Propagates a wavefront object to the next optic in the list, after some separation distance (which might be zero). Modifies this wavefront object itself.
Transformations between most planes use Fresnel propagation. If the target plane is an image plane, the output wavefront will be set to provide its coordinates in arcseconds based on its focal length, but it retains its internal dimensions in meters for future Fresnel propagations. Transformations to a Detector plane are handled separately to allow adjusting the pixel scale to match the target scale. Transformations from any frame through a rotation plane simply rotate the wavefront accordingly.
- optic : OpticalElement
- The optic to propagate to. Used for determining the appropriate optical plane.
- distance : astropy.Quantity of dimension length
- separation distance of this optic relative to the prior optic in the system.
-
static
pupil_coordinates
(x, y, pixelscale)[source]¶ Utility function to generate coordinates arrays for a pupil plane wavefront
- x, y : array_like
- pixel indices
- pixelscale : float or 2-tuple of floats
- the pixel scale in meters/pixel, optionally different in X and Y
- Y, X : array_like
- Wavefront coordinates in either meters or arcseconds for pupil and image, respectively
-
r_c
(z=None)[source]¶ The gaussian beam radius of curvature as a function of distance z
- z : float, optional
- Distance along the optical axis. If not specified, the wavefront’s current z coordinate will be used, returning the beam radius of curvature at the current position.
Astropy.units.Quantity of dimension length
-
spot_radius
(z=None)[source]¶ radius of a propagating gaussian wavefront, at a distance z
- z : float, optional
- Distance along the optical axis. If not specified, the wavefront’s current z coordinate will be used, returning the beam radius at the current position.
Astropy.units.Quantity of dimension length