FresnelWavefront

class poppy.FresnelWavefront(beam_radius, units=Unit('m'), rayleigh_factor=2.0, oversample=2, **kwargs)[source]

Bases: BaseWavefront

Wavefront for Fresnel diffraction calculation.

This class inherits from and extends the Fraunhofer-domain poppy.Wavefront class.

Parameters:
beam_radiusastropy.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.

unitsastropy.units.Unit

Astropy units of input parameters

rayleigh_factor:

Threshold for considering a wave spherical.

oversamplefloat

Padding factor to apply to the wavefront array, multiplying on top of the beam radius.

References

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_image_plane_fftmft(optic)

Apply a focal plane mask using fft and mft methods to highly sample at the focal plane.

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

Pixel scale, in arcsec/pixel or meters/pixel depending on plane type

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_image_plane_fftmft(optic)[source]

Apply a focal plane mask using fft and mft methods to highly sample at the focal plane.

Parameters:
opticFixedSamplingImagePlaneElement
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.

Parameters:
opticQuadraticLens

An optic

ignore_wavefrontboolean

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.

Returns:
Y, Xarray_like

Wavefront coordinates in either meters or arcseconds for pupil and image, respectively

display(*args, **kwargs)[source]

Display wavefront on screen

Parameters:
whatstring

What to display. Must be one of {intensity, phase, wfe, best, ‘both’}. ‘intensity’ shows the wavefront intensity, ‘wfe’ shows the wavefront error in meters or microns, ‘phase’ is similar to ‘wfe’ but shows wavefront phase in radians at the given wavelength. ‘Best’ implies to display the phase if there is nonzero OPD, or else display the intensity for a perfect pupil. ‘both’ will show two panels, for the wavefront intensity and wavefront error.

nrowsint

Number of rows to display in current figure (used for showing steps in a calculation)

rowint

Which row to display this one in? If set to None, use the wavefront’s self.current_plane_index

vmin, vmaxfloats

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.

vmax_wfefloat

Max value to display for phase or wfe, if distinct from vmax for intensity for the case of what=’both’. In other words you can use this to separately set the min/max scales for the two plots. The minimum scale for wfe will be set to the negative of this to ensure a balanced display scale symmetric around zero. This parameter is ignored if the ‘what’ parameter is not equal to ‘both’.

scalestring

‘log’ or ‘linear’, to define the desired display scale type for intensity. Default is log for image planes, linear otherwise.

imagecropfloat, 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 (unless you specify an angular_coordinate_unit parameter to choose another unit). 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).

showpaddingbool, 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.

colorbarbool

Display colorbar

crosshairsbool

Display a crosshairs indicator showing the axes centered on (0,0)

axmatplotlib Axes, optional

axes to display into. If not set, will create new axes.

use_angular_coordinatesbool, optional

Should the axes be labeled in angular units, e.g. 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)

angular_coordinate_unitastropy unit

Unit to use for angular coordinates display; default is arcsecond.

Returns:
figurematplotlib 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.

Parameters:
wavefrontWavefront

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.

Parameters:
zfloat

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.

This algorithm was implemented early on as a cross-check for the main propagate_fresnel pathway which uses the angular spectrum method. In most cases you should use that; calling this method is not recommended in general for most use cases, unless you have a particular reason to do so.

Parameters:
zfloat 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).

Parameters:
delta_zfloat

the distance from the current location to propagate the beam.

display_intermedboolean

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.

Parameters:
opticOpticalElement

The optic to propagate to. Used for determining the appropriate optical plane.

distanceastropy.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

Parameters:
x, yarray_like

pixel indices

pixelscalefloat or 2-tuple of floats

the pixel scale in meters/pixel, optionally different in X and Y

Returns:
Y, Xarray_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

Parameters:
zfloat, 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.

Returns:
Astropy.units.Quantity of dimension length
spot_radius(z=None)[source]

radius of a propagating gaussian wavefront, at a distance z

Parameters:
zfloat, 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.

Returns:
Astropy.units.Quantity of dimension length