Physical Optics Propagation in PYthon (POPPY)

Introduction:

The module poppy implements an object-oriented system for modeling physical optics propagation with diffraction, particularly for telescopic and coronagraphic imaging. Right now only image and pupil planes are supported; intermediate planes are a future goal. It also

This code makes use of the python standard module logging for output information. Top-level details of the calculation are output at level logging.INFO, while details of the propagation through each optical plane are printed at level logging.DEBUG. See the Python logging documentation for an explanation of how to redirect the poppy logger to the screen, a textfile, or any other log destination of your choice.

Poppy also includes a system for modeling a complete instrument (including optical propagation, synthetic photometry, and pointing jitter), and a variety of useful utility functions for analysing and plotting PSFs, documented below.

Key Concepts:

To model optical propagation, Poppy implements an object-oriented system for representing an optical train. There are a variety of OpticalElement classes representing both physical elements as apertures, mirrors, and apodizers, and also implicit operations on wavefronts, such as rotations or tilts. Each OpticalElement may be defined either via analytic functions (e.g. a simple circular aperture) or by numerical input FITS files (e.g. the complex JWST aperture with appropriate per-segment WFE). A series of such OpticalElements is chained together in order in an OpticalSystem class. That class is capable of generating Wavefronts (another class) suitable for propagation through the desired elemeents (with correct array size and sampling), and then performing the optical propagation onto the final image plane.

There is an even higher level class Instrument which adds support for selectable instrument mechanisms (such as filter wheels, pupil stops, etc). In particular it adds support for computing via synthetic photometry the appropriate weights for multiwavelength computations through a spectral bandpass filter, and for PSF blurring due to pointing jitter (neither of which effects are modeled by OpticalSystem). Given a specified instrument configuration, an appropriate OpticalSystem is generated, the appropriate wavelengths and weights are calculated based on the bandpass filter and target source spectrum, the PSF is calculated, and optionally is then convolved with a blurring kernel due to pointing jitter. All of the WebbPSF instruments are implemented by subclassing poppy.Instrument.

Poppy presently assumes that optical propagation can be modeled using Fraunhofer diffraction (far-field), such that the relationship between pupil and image plane optics is given by two-dimensional Fourier transforms. Fresnel propagation is not currently supported.

Two different algorithmic flavors of Fourier transforms are used in Poppy. The familiar FFT algorithm is used for transformations between pupil and image planes in the general case. This algorithm is relatively fast (O(N log(N))) but imposes strict constraints on the relative sizes and samplings of pupil and image plane arrays. Obtaining fine sampling in the image plane requires very large oversized pupil plane arrays and vice versa, and image plane pixel sampling becomes wavelength dependent. To avoid these constraints, for transforms onto the final Detector plane, instead a Matrix Fourier Transform (MFT) algorithm is used (See Soummer et al. 2007 Optics Express). This allows computation of the PSF directly on the desired detector pixel scale or an arbitrarily finely subsampled version therof. For equivalent array sizes N, the MFT is slower than the FFT(O(N^3)), but in practice the ability to freely choose a more appropriate N (and to avoid the need for post-FFT interpolation onto a common pixel scale) more than makes up for this and the MFT is faster.

List of Classes

digraph inheritancec987e26b1b { rankdir=LR; size="8.0, 12.0"; "poppy.poppy.OpticalElement" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.IdealFieldStop" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.IdealFieldStop" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.instrument.Instrument" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.CircularAperture" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.CircularAperture" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.FITSOpticalElement" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.OpticalElement" -> "poppy.poppy.FITSOpticalElement" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.IdealCircularOcculter" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.IdealCircularOcculter" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.HexagonAperture" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.HexagonAperture" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.CompoundAnalyticOptic" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.CompoundAnalyticOptic" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.Rotation" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.OpticalElement" -> "poppy.poppy.Rotation" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.IdealBarOcculter" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.IdealBarOcculter" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.Detector" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.OpticalElement" -> "poppy.poppy.Detector" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.AnalyticOpticalElement" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.OpticalElement" -> "poppy.poppy.AnalyticOpticalElement" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.BandLimitedCoron" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.BandLimitedCoron" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.SquareAperture" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.SquareAperture" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.FQPM_FFT_aligner" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.FQPM_FFT_aligner" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.IdealFQPM" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.AnalyticOpticalElement" -> "poppy.poppy.IdealFQPM" [arrowsize=0.5,style="setlinewidth(0.5)"]; "poppy.poppy.Wavefront" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "poppy.poppy.OpticalSystem" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; }

Wavefront

class poppy.Wavefront(wavelength=2e-06, npix=1024, dtype=<type 'numpy.complex128'>, diam=8.0, oversample=2, pixelscale=None)

A class representing a monochromatic wavefront that can be transformed between pupil and image planes (but not to intermediate planes, yet).

In a pupil plane, a wavefront object wf has
  • wf.diam, a diameter in meters
  • wf.pixelscale, a scale in meters/pixel
In an image plane, it has
  • wf.fov, a field of view in arcseconds
  • wf.pixelscale, a scale in arcsec/pixel

Use the wf.propagateTo() method to transform a wavefront between conjugate planes. This will update those properties as appropriate.

By default, Wavefronts are created in a pupil plane. Set pixelscale=# to make an image plane instead.

Parameters :

wavelength : float

Wavelength of light in meters

npix : int

Size parameter for wavefront array to create, per side.

diam : float, optional

For _PUPIL wavefronts, sets physical size corresponding to npix. Units are meters. At most one of diam or pixelscale should be set when creating a wavefront.

pixelscale : float, optional

For _IMAGE PLANE wavefronts, use this pixel scale.

oversample : int, optional

how much to oversample by in FFTs. Default is 2. Note that final propagations to Detectors use a different algorithmg and, optionally, a separate oversampling factor.

dtype : numpy.dtype, optional

default is double complex.

Methods

wavelength = None

Wavelength in meters

diam = None

Diameter in meters. Applies to a pupil plane only.

fov = None

Field of view in arcsec. Applies to an image plane only.

pixelscale = None

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

history = None

List of strings giving a descriptive history of actions performed on the wavefront. Saved to FITS headers.

location = None

A descriptive string for where a wavefront is instantaneously located (e.g. ‘before occulter’). Used mostly for titling displayed plots.

copy()

Return a copy of the wavefront as a different object.

normalize()

Set this wavefront’s total intensity to 1

asFITS(what='intensity', includepadding=False)

Return a wavefront as a pyFITS HDUList object

Parameters :

what : string

what kind of data to write. Must be one of ‘parts’, ‘intensity’, ‘complex’. The default is to write a file containing intensity.

writeto(filename, clobber=True, **kwargs)

Write a wavefront to a FITS file.

Parameters :

filename : string

filename to use

what : string

what to write. Must be one of ‘parts’, ‘intensity’, ‘complex’

clobber : bool, optional

overwhat existing? default is True

Returns :

outfile: file on disk :

The output is written to disk.

display(what='intensity', nrows=1, row=1, showpadding=False, imagecrop=None, colorbar=False, crosshairs=True, ax=None, title=None)

Display wavefront on screen

Parameters :

what : string

What to display. Must be one of {intensity, phase, best}. ‘Best’ implies ‘display the phase if there is 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?

imagecrop : float, optional

For image planes, set the maximum # of arcseconds to display. Default is 5, so only the innermost 5x5 arcsecond region will be shown.

showpadding : bool, optional

Show the entire padded arrays, or just the good parts? Default is False

colorbar : bool

Display colorbar

ax : matplotlib Axes

axes to display into

Returns :

figure : matplotlib figure

The current figure is modified.

propagateTo(optic)

Propagates a wavefront object to the next optic in the list. Modifies this wavefront object itself.

tilt(Xangle=0.0, Yangle=0.0)

Tilt a wavefront in X and Y.

Recall from Fourier optics (although this is straightforwardly rederivable by drawing triangles...) that for a wavefront tilted by some angle theta in radians, for a point r meters from the center of the pupil has

extra_pathlength = sin(theta) * r extra_waves = extra_pathlength/ wavelength = r * sin(theta) / wavelength

So we calculate the U and V arrays (corresponding to r for the pupil, in meters from the center) and then multiply by the appropriate trig factors for the angle.

The sign convention is chosen such that positive Yangle tilts move the star upwards in theg array at the focal plane. (This is sort of an inverse of what physically happens in the propagation to or through focus, but we’re ignoring that here and trying to just work in sky coords)

Parameters :

Xangle, Yangle : float

tilt angles, specified in arcseconds

rotate(angle=0.0)

Rotate a wavefront by some amount

Parameters :

angle : float

Angle to rotate, in degrees counterclockwise.

coordinates()

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 measures in wavefront coordinates.

Returns :

Y, X : array_like

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

Optical System

class poppy.OpticalSystem(name='unnamed system', verbose=True, oversample=2)

A class representing a series of optical elements, either Pupil, Image, or Detector planes, through which light can be propagated.

The difference between Image and Detector planes is that Detectors have fixed pixels in terms of arcsec/pixel regardless of wavelength (computed via MFT) while Image planes have variable pixels scaled in terms of lambda/D. Pupil planes are some fixed size in meters, of course.

Parameters :

name : string

descriptive name of optical system

oversample : int

Either how many times above Nyquist we should be (for pupil or image planes), or how many times a fixed detector pixel will be sampled. E.g. oversample=2 means image plane sampling lambda/4*D (twice Nyquist) and detector plane sampling 2x2 computed pixels per real detector pixel. Default is 2.

verbose : bool

whether to print stuff while computing

Methods

addPupil(optic=None, function=None, **kwargs)

Add a pupil plane optic from file(s) giving transmission or OPD

  1. from file(s) giving transmission and/or OPD

    [set arguments transmission=filename and/or opd=filename]

  2. from an analytic function

    [set function=’Circle’, ‘Square’ and set additional kwargs to define shape etc.

  3. from an already-created OpticalElement object

    [set optic=that object]

Parameters :

optic : poppy.OpticalElement, optional

An already-created OpticalElement you would like to add

function: string, optional :

Name of some analytic function to add. Optional kwargs can be used to set the parameters of that function. Allowable function names are Circle, Square

opd, transmission : string, optional

Filenames of FITS files describing the desired optic.

Note: Now you can use the optic arguement for either an OpticalElement or a string function name, :

and it will do the right thing depending on type. Both existing arguments are left for back compatibility for now. :

Any provided parameters are passed to :ref:`OpticalElement`. :

addImage(optic=None, function=None, **kwargs)

Add an image plane optic, either

  1. from file(s) giving transmission or OPD

    [set arguments transmission=filename and/or opd=filename]

  2. from an analytic function

    [set function=’circle, fieldstop, bandlimitedcoron, or FQPM’ and set additional kwargs to define shape etc.

  3. from an already-created OpticalElement object

    [set optic=that object]

Parameters :

optic : poppy.OpticalElement

An already-created OpticalElement you would like to add

function: string :

Name of some analytic function to add. Optional kwargs can be used to set the parameters of that function. Allowable function names are CircularOcculter, fieldstop, BandLimitedCoron, FQPM

opd, transmission : string

Filenames of FITS files describing the desired optic.

Note: Now you can use the optic arguement for either an OpticalElement or a string function name, :

and it will do the right thing depending on type. Both existing arguments are left for back compatibility for now. :

addDetector(pixelscale, oversample=None, **kwargs)

Add a Detector object to an optical system. By default, use the same oversampling as the rest of the optical system, but the user can override to a different value if desired by setting oversample.

Other arguments are passed to the init method for Detector().

Parameters :

pixelscale : float

Pixel scale in arcsec/pixel

oversample : int, optional

Oversampling factor for this detector, relative to hardware pixel size. Optionally distinct from the default oversampling parameter of the OpticalSystem.

inputWavefront(wavelength=2e-06)

Create a Wavefront object suitable for sending through a given optical system, based on the size of the first optical plane, assumed to be a pupil.

If the first optical element is an Analytic pupil (i.e. has no pixel scale) then an array of 1024x1024 will be created (not including oversampling).

Uses self.source_offset to assign an off-axis tilt, if requested.

Parameters :

wavelength : float

Wavelength in meters

Returns :

wavefront : poppy.Wavefront instance

A wavefront appropriate for passing through this optical system.

propagate_mono(wavelength=2e-06, normalize='first', save_intermediates=False, display_intermediates=False, intermediate_fn='wave_step_%03d.fits', poly_weight=None)

Propagate a wavefront through some number of optics. Returns a pyfits.HDUList object.

Parameters :

wavelength : float

Wavelength in meters

normalize : string, {‘first’, ‘last’}

how to normalize the wavefront? * ‘first’ = set total flux = 1 after the first optic, presumably a pupil * ‘last’ = set total flux = 1 after the entire optical system. * ‘first=2’ = set total flux = 2 after the first optic (used for debugging only)

display_intermediates : bool

Should intermediate steps in the calculation be displayed on screen? Default False

save_intermediates : bool

Should intermediate steps in the calculation be saved to disk? Default False. If this is True, then setting poly_weight controls whether intermediate optical planes are actually saved to disk by this routine (for the monochromatic case) or are passed back up via memory and handled in calcPSF (for the polychromatic case).

poly_weight : float

is this being called as part of a polychromatic calculation? if not, set this to None. if so, set this to the weight for that wavelength. (This is used only for properly normalizing the multiwavelength FITs file written to disk if save_intermediates is set.)

calcPSFmultiproc(source, nprocesses=4, save_intermediates=False, **kwargs)

Calculate a multi-wavelength PSF over some weighted sum of wavelengths, across multiple processors

This version uses Python’s multiprocessing package to span tasks across available processor cores.

Any additional kwargs will be passed on to propagate_mono()

Parameters :

source : dict

a dict containing ‘wavelengths’ and ‘weights’ list. TBD - replace w/ pysynphot observation object

nprocesses : int

Number of processes. Don’t make this too large, or you will overfill your memory and things will grind painfully to a half as everything swaps to disk and back.

save_intermediates : bool

whether to output intermediate optical planes to disk. Default is False

Returns :

outfits : :

a pyfits.HDUList

Notes

Don’t set the number of processes too high, or your machine will exhaust its memory and grind to a halt. Figure about 600 MB per process, due to the large arrays used in the FFTing, and multiple copies etc. Yes, this is surprisingly high overhead - but the raw wavefront array is typically 2048^2 * 16 (since it’s type dcomplex) = 67 MB, so this is just ~ 8-10 copies of the array floating arround. TODO: be more clever and efficient with all this somehow? TODO: write an auto tool to optimize the number of processes automatically?

calcPSF(wavelength=1e-06, weight=None, save_intermediates=False, save_intermediates_what='all', display=False, return_intermediates=False, source=None, **kwargs)

Calculate a PSF, either - multi-wavelength PSF over some weighted sum of wavelengths (if you provide a source parameter) - monochromatic (if you provide just a wavelen= parameter)

Any additional kwargs will be passed on to propagate_mono()

Parameters :

source : dict

a dict containing ‘wavelengths’ and ‘weights’ list. TBD - replace w/ pysynphot observation object

wavelen : float, optional

wavelength in meters for monochromatic calculation.

save_intermediates : bool

whether to output intermediate optical planes to disk. Default is False

save_intermediate_what : string

What to save - phase, intensity, amplitude, complex, parts, all. default is all.

return_intermediates: bool :

return intermediate wavefronts as well as PSF?

display : bool

whether to display when finished or not.

Returns :

outfits : :

a pyfits.HDUList

display(**kwargs)

Display all elements in an optical system on screen.

Any extra arguments are passed to the optic.display() methods of each element.

Optical Elements

class poppy.OpticalElement(name='unnamed optic', verbose=True, planetype=None, oversample=1, opdunits='meters')

Base class for all optical elements, whether from FITS files or analytic functions.

If instantiated on its own, this just produces a null optical element (empty space, i.e. an identity function on transmitted wavefronts.) Use one of the many subclasses to create a nontrivial optic.

The OpticalElement class follows the behavoior of the Wavefront class, using units of meters/pixel in pupil space and arcsec/pixel in image space.

Parameters :

name : string

descriptive name for optic

verbose : bool

whether to print stuff while computing

planetype : int

either poppy._IMAGE or poppy._PUPIL

oversample : int

how much to oversample beyond Nyquist.

Methods

name = None

string. Descriptive Name of this optic

getPhasor(wave)

Compute a complex phasor from an OPD, given a wavelength.

The returned value should be the complex phasor array as appropriate for multiplying by the wavefront amplitude.

Parameters :

wave : float or obj

either a scalar wavelength or a Wavefront object

display(nrows=1, row=1, what='intensity', phase=False, wavelength=None, crosshairs=True, ax=None, colorbar=True, colorbar_orientation=None, title=None)

Display plots showing an optic’s transmission and OPD

Parameters :

what : str

What do display: ‘intensity’, ‘phase’, or ‘both’

ax : matplotlib.Axes

Axes to display into


class poppy.FITSOpticalElement(name='unnamed optic', transmission=None, opd=None, opdunits='meters', shift=None, rotation=None, pixelscale=None, planetype=None, **kwargs)

Defines an arbitrary optic, based on amplitude transmission and/or OPD FITS files.

This optic could be a pupil or field stop, an aberrated mirror, a phase mask, etc. The FITSOpticalElement class follows the behavior of the Wavefront class, using units of meters/pixel in pupil space and arcsec/pixel in image space.

Parameters :

name : string

descriptive name for optic

transmission, opd : string or pyfits HDUList

Either FITS filenames or actual pyfits.HDUList objects for the transmission (from 0-1) and opd (in meters)

opdunits : string

units for the OPD file. Default is ‘meters’. can be ‘meter’, ‘meters’, ‘micron(s)’, ‘nanometer(s)’, or their SI abbreviations

verbose : bool

whether to print stuff while computing

planetype : int

either _IMAGE or _PUPIL

oversample : int

how much to oversample beyond Nyquist.

shift : tuple of floats, optional

2-tuple containing X and Y fractional shifts for the pupil.

rotation : float

Rotation for that optic, in degrees

pixelscale : optical str or float

By default, poppy will attempt to determine the appropriate pixel scale by examining the FITS header, checking keywords “PUPLSCAL” and ‘PIXSCALE’ for pupil and image planes respectively. If you would like to override and use a different keyword, provide that as a string here. Alternatively, you can just set a floating point value directly too (in meters/pixel or arcsec/pixel, respectively, for pupil or image planes).

*NOTE:* All mask files must be *squares*. :

Methods

class poppy.AnalyticOpticalElement(**kwargs)

Bases: poppy.poppy.OpticalElement

Defines an abstract analytic optical element, i.e. one definable by some formula rather than by an input OPD or pupil file.

This class is useless on its own; instead use its various subclasses that implement appropriate getPhasor functions. It exists mostly to provide some behaviors & initialization common to all analytic optical elements.

Parameters :

name, verbose, oversample, planetype : various

Same as for OpticalElement

transmission, opd : string

These are not allowed for Analytic optical elements, and this class will raise an error if you try to set one.

Methods

class poppy.CircularAperture(name=None, radius=1.0, pad_factor=1.5, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal circular pupil aperture

Parameters :

name : string

Descriptive name

radius : float

Radius of the pupil, in meters. Default is 1.0

pad_factor : float, optional

Amount to oversize the wavefront array relative to this pupil. This is in practice not very useful, but it provides a straightforward way of verifying during code testing that the amount of padding (or size of the circle) does not make any numerical difference in the final result.

Methods

class poppy.HexagonAperture(name=None, flattoflat=None, side=None, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal hexagonal pupil aperture

Specify either the side length (= corner radius) or the flat-to-flat distance.

Parameters :

name : string

Descriptive name

flattoflat : float, optional

Flat-to-flat distance of the pupil, in meters. Default is 1.0

side : float, optional

side length of hexagon, in meters.

Methods

class poppy.SquareAperture(name=None, size=1.0, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal square pupil aperture

Parameters :

name : string

Descriptive name

size: float :

half width of the square, in meters. Default is 1.0

Methods

class poppy.IdealFieldStop(name='unnamed field stop', size=20.0, angle=0, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal square field stop

Parameters :

name : string

Descriptive name

size : float

Size of the field stop, in arcseconds. Default 20.

angle : float

Position angle of the field stop sides relative to the detector +Y direction, in degrees.

Methods

class poppy.IdealCircularOcculter(name='unnamed occulter', radius=1.0, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal circular occulter (opaque circle)

Parameters :

name : string

Descriptive name

radius : float

Radius of the occulting spot, in arcseconds. Default is 1.0

Methods

class poppy.IdealBarOcculter(name='bar occulter', width=1.0, angle=0, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal bar occulter (like in MIRI’s Lyot coronagraph)

Parameters :

name : string

Descriptive name

width : float

width of the bar stop, in arcseconds. Default is 1.0

angle : float

position angle of the bar, rotated relative to the normal +y direction.

Methods

class poppy.BandLimitedCoron(name='unnamed BLC', kind='circular', sigma=1, wavelength=None, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal band limited coronagraph occulting mask.

Parameters :

name : string

Descriptive name

kind : string

Either ‘circular’ or ‘linear’. The linear ones are custom shaped to NIRCAM’s design with flat bits on either side of the linear tapered bit. Also includes options ‘nircamcircular’ and ‘nircamwedge’ specialized for the JWST NIRCam occulters, including the off-axis ND acq spots and the changing width of the wedge occulter.

sigma : float

The numerical size parameter, as specified in Krist et al. 2009 SPIE

wavelength : float

Wavelength this BLC is optimized for, only for the linear ones.

Methods

class poppy.IdealFQPM(name='unnamed FQPM ', wavelength=1.065e-05, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Defines an ideal 4-quadrant phase mask coronagraph, with its retardance set perfectly to 0.5 waves at one specific wavelength and varying linearly on either side of that.

Parameters :

name : string

Descriptive name

wavelength : float

Wavelength in meters for which the FQPM was designed

Methods

class poppy.FQPM_FFT_aligner(name='FQPM FFT aligner', direction='forward', **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Helper class for modeling FQPMs accurately

Adds (or removes) a slight wavelength- and pixel-scale-dependent tilt to a pupil wavefront, to ensure the correct alignment of the image plane FFT’ed PSF with the desired quad pixel alignment for the FQPM.

This is purely a computational convenience tool to work around the pixel coordinate restrictions imposed by the FFT algorithm, not a representation of any physical optic.

Parameters :

direction : string

‘forward’ or ‘backward’

Methods

class poppy.CompoundAnalyticOptic(name='unnamed', verbose=True, opticslist=None, **kwargs)

Bases: poppy.poppy.AnalyticOpticalElement

Define a compound analytic optical element made up of the combination of two or more individual optical elements.

This is just a convenience routine for semantic organization of optics. It can be useful to keep the list of optical planes cleaner, but you can certainly just add a whole bunch of planes all in a row without using this class to group them.

All optics should be of the same plane type (pupil or image); propagation between different optics contained inside one compound is not supported.

Parameters :

opticslist : list

A list of AnalyticOpticalElements to be merged together.

Methods

class poppy.Rotation(angle=0.0, units='degrees', **kwargs)

Bases: poppy.poppy.OpticalElement

Performs a rotation of the axes in the optical train.

This is not an actual optic itself, of course, but can be used to model a rotated optic by appling a Rotation before and/or after light is incident on that optic.

This is basically a placeholder to indicate the need for a rotation at a given part of the optical train. The actual rotation computation is performed in the Wavefront object’s propagation routines.

Parameters :

angle : float

Rotation angle, counterclockwise. By default in degrees.

units : ‘degrees’ or ‘radians’

Units for the rotation angle.

Methods


class poppy.Detector(pixelscale, fov_pixels=None, fov_arcsec=None, oversample=1, name='Detector', offset=None)

Bases: poppy.poppy.OpticalElement

A Detector is a specialized type of OpticalElement that forces a wavefront onto a specific fixed pixelization.

Parameters :

name : string

Descriptive name

pixelscale : float

Pixel scale in arcsec/pixel

fov_pixels, fov_arcsec : float

The field of view may be specified either in arcseconds or by a number of pixels. Either is acceptable and the pixel scale is used to convert as needed. You may specify a non-square FOV by providing two elements in an iterable.

oversample : int

Oversampling factor beyond the detector pixel scale

offset : tuple (X,Y)

Offset for the detector center relative to a hypothetical off-axis PSF. Specifying this lets you pick a different sub-region for the detector to compute, if for some reason you are computing a small subarray around an off-axis source. (Has not been tested!)

Methods

Instrument

class poppy.Instrument(name='', *args, **kwargs)

A generic astronomical instrument, composed of (1) an optical system implemented using POPPY, and (2) some defined spectral bandpass(es), implemented using pysynphot. This provides the capability to model both the optical and spectral responses of a given system. PSFs may be calculated for given source spectral energy distributions and output as FITS files, with substantial flexibility.

It also provides capabilities for modeling some PSF effects not due to wavefront aberrations, for instance blurring caused by pointing jitter.

This is a base class for Instrument functionality - you will likely not want to use this directly, but rather should subclass it for your particular instrument of interest. Some of the complexity of this class is due to splitting up functionality into many separate routines to allow users to subclass just the relevant portions for a given task.

You will at a minimum want to override the following class methods:

_getOpticalSystem _getFilterList _getDefaultNLambda _getDefaultFOV _getFITSHeader
For more complicated systems you may also want to override:
_validateConfig _getSynphotBandpass _applyJitter

Methods

pupil = None

Aperture for this optical system. May be a FITS filename, pyfits.HDUList, or poppy.OpticalElement

pupilopd = None

Pupil OPD for this optical system. May be a FITS filename, or pyfits.HDUList. If the file contains a datacube, you may set this to a tuple (filename, slice) to select a given slice, or else the first slice will be used.

options = None

A dictionary capable of storing other arbitrary options, for extensibility. The following are all optional, and may or may not be meaningful depending on which instrument is selected.

Parameters :

source_offset_r : float

Radial offset of the target from the center, in arcseconds

source_offset_theta : float

Position angle for that offset

pupil_shift_x, pupil_shift_y : float

Relative shift of a coronagraphic pupil in X and Y, expressed as a decimal between 0.0-1.0 Note that shifting an array too much will wrap around to the other side unphysically, but for reasonable values of shift this is a non-issue.

rebin : bool

For output files, write an additional FITS extension including a version of the output array rebinned down to the actual detector pixel scale?

jitter : string

Type of jitter model to apply. Currently not implemented

parity : string “even” or “odd”

You may wish to ensure that the output PSF grid has either an odd or even number of pixels. Setting this option will force that to be the case by increasing npix by one if necessary.

pixelscale = None

Detector pixel scale, in arcseconds/pixel

filter

Currently selected filter name (e.g. “F200W”)

calcPSF(outfile=None, source=None, nlambda=None, monochromatic=None, fov_arcsec=None, fov_pixels=None, oversample=None, detector_oversample=None, fft_oversample=None, rebin=True, clobber=True, display=False, save_intermediates=False, return_intermediates=False)

Compute a PSF. The result can either be written to disk (set outfile=”filename”) or else will be returned as a pyfits HDUlist object.

Output sampling may be specified in one of two ways:

  1. Set oversample=<number>. This will use that oversampling factor beyond detector pixels for output images, and beyond Nyquist sampling for any FFTs to prior optical planes.
  2. set detector_oversample=<number> and fft_oversample=<other_number>. This syntax lets you specify distinct oversampling factors for intermediate and final planes.

By default, both oversampling factors are set equal to 2.

Parameters :

source : pysynphot.SourceSpectrum or dict

specification of source input spectrum. Default is a 5700 K sunlike star.

nlambda : int

How many wavelengths to model for broadband? The default depends on how wide the filter is: (5,3,1) for types (W,M,N) respectively

monochromatic : float, optional

Setting this to a wavelength value (in meters) will compute a monochromatic PSF at that wavelength, overriding filter and nlambda settings.

fov_arcsec : float

field of view in arcsec. Default=5

fov_pixels : int

field of view in pixels. This is an alternative to fov_arcsec.

outfile : string

Filename to write. If None, then result is returned as an HDUList

oversample, detector_oversample, fft_oversample : int

How much to oversample. Default=4. By default the same factor is used for final output pixels and intermediate optical planes, but you may optionally use different factors if so desired.

rebin : bool, optional

If set, the output file will contain a FITS image extension containing the PSF rebinned onto the actual detector pixel scale. Thus, setting oversample=<N> and rebin=True is the proper way to obtain high-fidelity PSFs computed on the detector scale. Default is True.

clobber : bool

overwrite output FITS file if it already exists?

display : bool

Whether to display the PSF when done or not.

save_intermediates, return_intermediates : bool

Options for saving to disk or returning to the calling function the intermediate optical planes during the propagation. This is useful if you want to e.g. examine the intensity in the Lyot plane for a coronagraphic propagation.

Returns :

outfits : pyfits.HDUList

The output PSF is returned as a pyfits.HDUlist object. If outfile is set to a valid filename, the output is also written to that file.

Notes

More advanced PSF computation options (pupil shifts, source positions, jitter, ...) may be set by configuring the .options dictionary attribute of this class.

display()

Display the currently configured optical system on screen



Documentation last updated on June 25, 2012