JWInstrument

class webbpsf.JWInstrument(*args, **kwargs)[source]

Bases: SpaceTelescopeInstrument

Superclass for all JWST instruments

Attributes Summary

aperturename

SIAF aperture name for detector pixel to sky coords transformations

pupilopd

Filename or fits.HDUList for JWST pupil OPD.

telescope

Methods Summary

calc_psf([outfile, source, nlambda, ...])

Compute a PSF.

get_opd_file_full_path([opdfilename])

Return full path to the named OPD file.

get_optical_system([fft_oversample, ...])

Return an OpticalSystem instance corresponding to the instrument as currently configured.

get_wfe([kind, wavelength, plot])

Extract and return one component plane of the optical model for this instrument

interpolate_was_opd(array, newdim)

Interpolates an input 2D array to any given size.

load_wss_opd(filename[, output_path, ...])

Load an OPD produced by the JWST WSS into this instrument instance, specified by filename

load_wss_opd_by_date([date, choice, ...])

Load an OPD produced by the JWST WSS into this instrument instance, specified by filename.

set_position_from_aperture_name(aperture_name)

Set the simulated center point of the array based on a named SIAF aperture.

visualize_wfe_budget([slew_delta_time, ...])

Display a visual WFE budget showing the various terms that sum into the overall WFE for a given instrument

Attributes Documentation

aperturename

SIAF aperture name for detector pixel to sky coords transformations

pupilopd = None

Filename or fits.HDUList for JWST pupil OPD.

This can be either a full absolute filename, or a relative name in which case it is assumed to be within the instrument’s data/OPDs/ directory, or an actual fits.HDUList object corresponding to such a file. 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.

telescope = 'JWST'

Methods Documentation

calc_psf(outfile=None, source=None, nlambda=None, monochromatic=None, fov_arcsec=None, fov_pixels=None, oversample=None, detector_oversample=None, fft_oversample=None, overwrite=True, display=False, save_intermediates=False, return_intermediates=False, normalize='first', add_distortion=True, crop_psf=True)[source]

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

Output sampling may be specified in one of two ways:

  1. Set oversample=. 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= and fft_oversample=. This syntax lets you specify distinct oversampling factors for intermediate and final planes.

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

Parameters:
sourcesynphot.spectrum.SourceSpectrum or dict

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

nlambdaint

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

monochromaticfloat, optional

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

fov_arcsecfloat

field of view in arcsec. Default=5

fov_pixelsint

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

outfilestring

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

oversample, detector_oversample, fft_oversampleint

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.

overwritebool

overwrite output FITS file if it already exists?

displaybool

Whether to display the PSF when done or not.

save_intermediates, return_intermediatesbool

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.

normalizestring

Desired normalization for output PSFs. See doc string for OpticalSystem.calc_psf. Default is to normalize the entrance pupil to have integrated total intensity = 1.

add_distortionbool

If True, will add 2 new extensions to the PSF HDUlist object. The 2nd extension will be a distorted version of the over-sampled PSF and the 3rd extension will be a distorted version of the detector-sampled PSF.

crop_psfbool

If True, when the PSF is rotated to match the detector’s rotation in the focal plane, the PSF will be cropped so the shape of the distorted PSF will match it’s undistorted counterpart. This will only be used for NIRCam, NIRISS, and FGS PSFs.

Returns:
outfitsfits.HDUList

The output PSF is returned as a fits.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.

get_opd_file_full_path(opdfilename=None)[source]

Return full path to the named OPD file.

The OPD may be:
  • a local or absolute path,

  • or relative implicitly within an SI directory, e.g. $WEBBPSF_PATH/NIRCam/OPD

  • or relative implicitly within $WEBBPSF_PATH

This function handles filling in the implicit path in the latter cases.

get_optical_system(fft_oversample=2, detector_oversample=None, fov_arcsec=2, fov_pixels=None, options=None)[source]

Return an OpticalSystem instance corresponding to the instrument as currently configured.

When creating such an OpticalSystem, you must specify the parameters needed to define the desired sampling, specifically the oversampling and field of view.

Parameters:
fft_oversampleint

Oversampling factor for intermediate plane calculations. Default is 2

detector_oversample: int, optional

By default the detector oversampling is equal to the intermediate calculation oversampling. If you wish to use a different value for the detector, set this parameter. Note that if you just want images at detector pixel resolution you will achieve higher fidelity by still using some oversampling (i.e. not setting oversample_detector=1) and instead rebinning down the oversampled data.

fov_pixelsfloat

Field of view in pixels. Overrides fov_arcsec if both set.

fov_arcsecfloat

Field of view, in arcseconds. Default is 2

Returns:
osyspoppy.OpticalSystem

an optical system instance representing the desired configuration.

get_wfe(kind='si', wavelength=2e-06, plot=False)[source]

Extract and return one component plane of the optical model for this instrument

This is a utility function for convenience, making it easier to access and plot various OPD maps. It doesn’t do anything unique which can’t be done otherwise, and in particular this isn’t used at all as part of the optical propagation calculations.

Note, all WFE terms are returned in OTE entrance pupil orientation (i.e. as if you were in front of the OTE and looking at it), regardless of pupil flips and orientations in the optical propagation.

Parameters:
kindstring

A type of WFE. Must be one of “SI”, “OTE”, “OTE_field_dep”, or other values TBD. Case insensitive.

plotbool

Make a quick plot of this WFE. Not very flexible or scriptable but useful for some interactive checks

interpolate_was_opd(array, newdim)[source]

Interpolates an input 2D array to any given size.

Parameters:
array: float

input array to interpolate

newdim: int

new size of the 2D square array (newdim x newdim)

Returns:
newopd: new array interpolated to (newdim x newdim)
load_wss_opd(filename, output_path=None, backout_si_wfe=True, verbose=True, plot=False, save_ote_wfe=False)[source]

Load an OPD produced by the JWST WSS into this instrument instance, specified by filename

This includes:
  • If necessary, downloading that OPD from MAST. Downloaded files are cached in $WEBBPSF_PATH/MAST_JWST_WSS_OPDs

  • calling import_wss_opd to load the OPD from the FITS file and

    perform some necessary format conversions

  • Subtract off the instrument WFE for the field point used in wavefront sensing, to get an

    OTE-only wavefront. WebbPSF will separately add back in the SI WFE for the appropriate

    field point, as usual.

  • Subtract off the modeled field dependence term in the OTE WFE for the sensing field point, to get

    an estimate of the OTE wavefront nominally at the master chief ray location (between the NIRCams). WebbPSF will automatically add back on top of this the OTE field dependent WFE for the appropriate field point. as usual.

Parameters:
filenamestr

Name of OPD file to load

output_pathstr

Downloaded OPD are saved in this location. This option is convinient for STScI users using /grp/jwst/ote/webbpsf-data/. Default is $WEBBPSF_PATH/MAST_JWST_WSS_OPDs

backout_si_wfebool

Subtract model for science instrument WFE at the sensing field point? Generally this should be true which is the default.

plotbool

Generate informative plots showing WFE, including the backout steps. Only works if backout_si_wfe is True.

save_ote_wfebool

Save OTE-only WFE model? This is not needed for calculations in WebbPSF, but can be used to export OTE WFE models for use with other software. The file will be saved in the WEBBPSF_DATA_PATH directory and a message will be printed on screen with the filename. Note that the exported OPD file will give the OTE estimated total WFE at the selected Instrument’s field point, not the OTE global at master chief ray, since it is the OTE WFE at the selected field point which is most of use for some other tool.

load_wss_opd_by_date(date=None, choice='closest', verbose=True, plot=False, **kwargs)[source]

Load an OPD produced by the JWST WSS into this instrument instance, specified by filename.

This does a MAST query by date to identify the relevant OPD file, then calls load_wss_opd.

Parameters:
date: string

Date time in UTC as ISO-format string, a la 2021-12-25T07:20:00 Note, if date is left unspecified as None, the most recent available measurement will be retrieved.

choicestring

Method to choose which OPD file to use, e.g. ‘before’, ‘after’

Further keyword parameters may be passed via **kwargs to load_wss_opd
set_position_from_aperture_name(aperture_name)[source]

Set the simulated center point of the array based on a named SIAF aperture. This will adjust the detector and detector position attributes.

visualize_wfe_budget(slew_delta_time=<Quantity 14. d>, slew_case='EOL', ptt_only=False, verbose=True)[source]

Display a visual WFE budget showing the various terms that sum into the overall WFE for a given instrument

Compares a WebbPSF instrument instance with the JWST optical budget for that instrument

Parameters:
instwebbpsf.JWInstrument

A JWST instrument instance

slew_delta_timeastropy.Quantity time

Time duration for thermal slew model

slew_casebasestring

‘BOL’ or ‘EOL’ for beginning of life or end of life thermal slew model. EOL is about 3x higher amplitude

ptt_onlybool

When decomposing wavefront into controllable modes, use a PTT-only basis? The default is to use all controllable pose modes. (This is mostly a leftover debug option at this point, not likely useful in general)

verbosebool

Be more verbose