The WebbPSF module

This module provides the primary interface, both for programmers and for interactive non-GUI use. It provides five classes corresponding to the JWST instruments, with consistent interfaces. See below for the detailed API; for now let’s dive into some example code.

Usage and Examples

Simple PSFs are easily obtained:

>>> from webbpsf import *
>>> nc = NIRCam()
>>> nc.filter =  'F200W'
>>> psf = nc.calcPSF(oversample=4)      # returns a pyfits.HDUlist containing PSF and header
>>> pylab.imshow(psf[0].data]           # display it on screen yourself, or
>>> display_PSF(psf)                    # use this convenient function to make a nice log plot with labeled axes
>>>
>>> psf = nc.calcPSF(filter='F470N', oversample=4)    # this is just a shortcut for setting the filter, then computing a PSF
>>>
>>> nc.calcPSF("myPSF.fits", filter='F480M' )         # you can also write the output directly to disk if you prefer.

For interactive use, you can have the PSF displayed as it is computed:

>>> nc.calcPSF(display=True)                          # will make nice plots with matplotlib.
Sample PSF image

More complicated instrumental configurations are available by setting the instrument’s attributes. For instance, one can create an instance of MIRI and configure it for coronagraphic observations, thus:

>>> miri = MIRI()
>>> miri.filter = 'F1065C'
>>> miri.image_mask = 'FQPM1065'
>>> miri.pupil_mask = 'MASKFQPM'
>>> miri.calcPSF('outfile.fits')
Sample PSF image

Input Source Spectra

WebbPSF attempts to calculate realistic weighted broadband PSFs taking into account both the source spectrum and the instrumental spectral response.

The default source spectrum is, if pysynphot is installed, a G2V star spectrum from Castelli & Kurucz 2004. Without pysynphot, the default is a flat spectrum in F_\nu such that the same number of photons are detected at each wavelength.

You may choose a different illuminating source spectrum by specifying a source parameter in the call to calcPSF(). The following are valid sources:

  1. A pysynphot.Spectrum object. This is the best option, providing maximum ease and accuracy, but requires the user to have pysynphot installed. In this case, the Spectrum object is combined with a pysynphot.ObsBandpass for the selected instrument and filter to derive the effective stimulus in detected photoelectrons versus wavelength. This is binned to the number of wavelengths set by the nlambda parameter.

  2. A dictionary with elements source["wavelengths"] and source["weights"] giving the wavelengths in meters and the relative weights for each. These should be numpy arrays or lists. In this case, the wavelengths and weights are used exactly as provided, without applying the instrumental filter profile.

    >>> src = {'wavelengths': [2.0e-6, 2.1e-6, 2.2e-6], 'weights': [0.3, 0.5, 0.2]}
    >>> nc.calcPSF(source=src, outfile='psf_for_src.fits')
    
  3. A tuple or list containing the numpy arrays (wavelength, weights) instead.

To calculate a monochromatic PSF, just use the monochromatic parameter. Wavelengths are always specified in meters. This is just a shorthand for a single-element source dict.

>>> miri.calcPSF(monochromatic=9.876e-6)

As a convenience, webbpsf includes a function to retrieve an appropriate pysynphot.Spectrum object for a given stellar spectral type from the PHOENIX or Castelli & Kurucz model libraries.

>>> src = webbpsf.specFromSpectralType('G0V', catalog='phoenix')
>>> psf = miri.calcPSF(source=src)

Array sizes, star positions, and centering

Output array sizes may be specified either in units of arcseconds or pixels. For instance,

>>> mynircam = NIRCam()
>>> result = mynircam.calcPSF(fov_arcsec=7, oversample=2, filter='F250M')
>>> result2= mynircam.calcPSF(fov_pixels=512, oversample=2, filter='F250M')

By default, the PSF will be centered at the exact center of the output array. This means that if the PSF is computed on an array with an odd number of pixels, the PSF will be centered exactly on the central pixel. If the PSF is computed on an array with even size, it will be centered on the “crosshairs” at the intersection of the central four pixels. If one of these is particularly desirable to you, set the parity option appropriately:

>>>  instrument.options['parity'] = 'even'
>>>  instrument.options['parity'] = 'odd'

Setting one of these options will ensure that a field of view specified in arcseconds is properly rounded to either odd or even when converted from arcsec to pixels. Alternatively, you may also just set the desired number of pixels explicitly in the call to calcPSF():

>>>  instrument.calcPSF(fov_npixels = 512)

Offset sources

The PSF may also be shifted off-center by adjusting the offset of the stellar source. This is done in polar coordinates:

>>> instrument.options['source_offset_r'] = 0.3         # offset in arcseconds
>>> instrument.options['source_offset_theta'] = 45.     # degrees counterclockwise from instrumental +Y in the science frame

If these options are set, the offset is applied relative to the central coordinates as defined by the output array parity.

For coronagraphic modes, the coronagraph occulter is always assumed to be at the center of the output array. Therefore, these options let you offset the source away from the coronagraph.

Pixel scales, sampling, and oversampling

The derived instrument classes all know their own instrumental pixel scales. You can change the output pixel scale in a variety of ways, as follows. See the JWInstrument.calcPSF documentation for more details.

  1. Set the oversample parameter to calcPSF(). This will produce a PSF with a pixel grid this many times more finely sampled. oversample=1 is the native detector scale, oversample=2 means divide each pixel into 2x2 finer pixels, and so forth. You can automatically obtain both the oversampled PSF and a version rebinned down onto the detector pixel scale by setting rebin=True in the call to calcPSF:

    >>> hdulist = instrument.calcPSF(oversample=2, rebin=True)     # hdulist will contain a primary HDU with the
    >>>                                                            # oversampled data, plus an image extension
    >>>                                                            # with the PSF rebinned down to regular sampling.
    
  2. For coronagraphic calculations, it is possible to set different oversampling factors at different parts of the calculation. See the calc_oversample and detector_oversample parameters. This is of no use for regular imaging calculations (in which case oversample is a synonym for detector_oversample). Specifically, the calc_oversample keyword is used for Fourier transformation to and from the intermediate optical plane where the occulter (coronagraph spot) is located, while detector_oversample is used for propagation to the final detector. Note that the behavior of these keywords changes for coronagraphic modeling using the Semi-Analytic Coronagraphic propagation algorithm (not fully documented yet - contact Marshall Perrin if curious).

    >>> miri.calcPSF(calc_oversample=8, detector_oversample= 2)    # model the occulter with very fine pixels, then save the
    >>>                                                           # data on a coarser (but still oversampled) scale
    
  3. Or, if you need even more flexibility, just change the instrument.pixelscale attribute to be whatever arbitrary scale you require.

    >>> instrument.pixelscale = 0.0314159
    

Note that the calculations performed by WebbPSF are somewhat memory intensive, particularly for coronagraphic observations. All arrays used internally are double-precision complex floats (16 bytes per value), and many arrays of size (npixels*oversampling)^2 are needed (particularly if display options are turned on, since the Matplotlib graphics library makes its own copy of all arrays displayed). Your average laptop with a couple GB of RAM will do perfectly well for most computations so long as you’re not too ambitious with setting array size and oversampling. If you’re interested in very high fidelity simulations of large fields (e.g. 1024x1024 pixels oversampled 8x) then we recommend a large multicore desktop with >16 GB RAM.

Output format options for sampling

As just explained, WebbPSF can easily calculate PSFs on a finer grid than the detector’s native pixel scale. You can select whether the output data should include this oversampled image, a copy that has instead been rebinned down to match the detector scale, or optionally both. This is done using the options['output_mode'] parameter.

>>> nircam.options['output_mode'] = 'oversampled'
>>> psf = nircam.calcPSF()       # the 'psf' variable will be an oversampled PSF, formatted as a FITS HDUlist
>>> nircam.options['output_mode'] = 'detector sampled'
>>> psf2 = nircam.calcPSF()      # now 'psf2' will contain the resul as resampled onto the detector scale.
>>> nircam.options['output_mode'] = 'both'
>>> psf3 = nircam.calcPSF()      # 'psf3' will have the oversampled image as primary HDU, and
>>>                              # the detector-sampled image as the first image extension HDU.

The JWInstrument generic class

digraph inheritance93c8e25c8f { rankdir=LR; size="8.0, 12.0"; "webbpsf.webbpsf.NIRISS" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "webbpsf.webbpsf.JWInstrument" -> "webbpsf.webbpsf.NIRISS" [arrowsize=0.5,style="setlinewidth(0.5)"]; "webbpsf.webbpsf.FGS" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "webbpsf.webbpsf.JWInstrument" -> "webbpsf.webbpsf.FGS" [arrowsize=0.5,style="setlinewidth(0.5)"]; "webbpsf.webbpsf.JWInstrument" [style="setlinewidth(0.5)",URL="#webbpsf.webbpsf.JWInstrument",fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25,shape=box,fontsize=10]; "poppy.instrument.Instrument" -> "webbpsf.webbpsf.JWInstrument" [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]; "webbpsf.webbpsf.MIRI" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "webbpsf.webbpsf.JWInstrument" -> "webbpsf.webbpsf.MIRI" [arrowsize=0.5,style="setlinewidth(0.5)"]; "webbpsf.webbpsf.NIRCam" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "webbpsf.webbpsf.JWInstrument" -> "webbpsf.webbpsf.NIRCam" [arrowsize=0.5,style="setlinewidth(0.5)"]; "webbpsf.webbpsf.NIRSpec" [shape=box,style="setlinewidth(0.5)",fontsize=10,fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,height=0.25]; "webbpsf.webbpsf.JWInstrument" -> "webbpsf.webbpsf.NIRSpec" [arrowsize=0.5,style="setlinewidth(0.5)"]; }

class webbpsf.webbpsf.JWInstrument(name='')

A generic JWST Instrument class.

Note: Do not use this class directly - instead use one of the specific_instrument subclasses!

This class provides a simple interface for modeling PSF formation through the JWST instruments, with configuration options and software interface loosely resembling the configuration of the instrument mechanisms.

This module currently only provides a modicum of error checking, and relies on the user being knowledgable enough to avoid trying to simulate some physically impossible or just plain silly configuration (such as trying to use a FQPM with the wrong filter).

The instrument constructors do not take any arguments. Instead, create an instrument object and then configure the filter or other attributes as desired. The most commonly accessed parameters are available as object attributes: filter, image_mask, pupil_mask, pupilopd. More advanced configuration can be done by editing the JWInstrument.options dictionary, either by passing options to __init__ or by directly editing the dict afterwards.

Methods

pupil = None

Filename or pyfits.HDUList for JWST pupil mask. Usually there is no need to change this.

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.

force_coron : bool

Set this to force full coronagraphic optical propagation when it might not otherwise take place (e.g. calculate the non-coronagraphic images via explicit propagation to all optical surfaces, FFTing to intermediate pupil and image planes whether or not they contain any actual optics, rather than taking the straight-to-MFT shortcut)

no_sam : bool

Set this to prevent the SemiAnalyticMethod coronagraph mode from being used when possible, and instead do the brute-force FFT calculations. This is usually not what you want to do, but is available for comparison tests. The SAM code will in general be much faster than the FFT method, particularly for high oversampling.

filter_list = None

List of available filters

pupilopd = None

Filename or pyfits.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 pyfits.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.

image_mask_list = None

List of available image_masks

pupil_mask_list = None

List of available pupil_masks

pixelscale = None

Detector pixel scale, in arcsec/pixel

filter

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

image_mask

Currently selected coronagraphic image plane mask, or None for direct imaging

pupil_mask

Currently selected Lyot pupil mask, or None for direct imaging

calcPSF(outfile=None, source=None, filter=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 :

filter : string, optional

Filter name. Setting this is just a shortcut for setting the object’s filter first, then calling calcPSF afterwards.

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.

Notes on Specific Instruments

NIRCam

class webbpsf.NIRCam

A class modeling the optics of NIRCam.

Relevant attributes include filter, image_mask, and pupil_mask.

The NIRCam class is smart enough to automatically select the appropriate pixel scale for the short or long wavelength channel based on whether you request a short or long wavelength filter.

Methods

See methods under JWInstrument

NIRSpec

class webbpsf.NIRSpec

A class modeling the optics of NIRSpec, in imaging mode.

This is not a substitute for a spectrograph model, but rather a way of simulating a PSF as it would appear with NIRSpec in imaging mode (e.g. for target acquisition). NIRSpec support is relatively simplistic compared to the other instruments at this point.

Relevant attributes include filter. In addition to the actual filters, you may select ‘IFU’ to indicate use of the NIRSpec IFU, in which case use the monochromatic attribute to set the simulated wavelength. Note: IFU to be implemented later

Methods

See methods under JWInstrument

MIRI

class webbpsf.MIRI

A class modeling the optics of MIRI, the Mid-InfraRed Instrument.

Relevant attributes include filter, image_mask, and pupil_mask.

In addition to the actual filters, you may select ‘MRS-IFU Ch1’ to indicate use of the MIRI IFU in Channel 1, and so forth. In this case, the monochromatic attribute controls the simulated wavelength. Note that the pixel scale varies with channel, which is why they are implemented separately. Note: IFU to be implemented later

Methods

See methods under JWInstrument

Sample PSF image for MIRI

An example MIRI PSF in F1000W.

Note that the MIRI imager field of view is rotated by 4.56 degrees relative to the JWST pupil; the coronagraph optics are correspondingly counterrotated to align them with the pupil. For direct imaging PSF calculations, this is most simply handled by rotating the pupil mask and OPD file prior to the Fourier propagation. For MIRI coronagraphy on the other hand, the rotation is performed as the last step prior to the detector.

Technical aside: Note that for computational reasons having to do with accurately simulating PSF centering on an FQPM, MIRI corongraphic simulations will include two ‘virtual optics’ called ‘FQPM FFT aligners’ that will show up in the display window for such calculations. These can be ignored by most end users of this software; interested readers should consult the POPPY documentation for more detail.

NIRISS

class webbpsf.NIRISS
A class modeling the optics of the Near-IR Imager and Slit Spectrograph
(formerly nTFI)

Relevant attributes include image_mask, and pupil_mask.

Methods

See methods under JWInstrument

FGS

class webbpsf.FGS

A class modeling the optics of the FGS.

Not a lot to see here, folks: There are no selectable options, just a great big detector-wide bandpass.

Methods

See methods under JWInstrument

TFI

Deprecated in favor of NIRISS. The TFI class is still included in this version of WebbPSF for compatibility with existing code, just to be on the safe side, but is likely to go away in the next release.

Utility Functions for Display and Plotting

webbpsf.Instrument(name)

This is just a convenience function, allowing one to access instrument objects based on a string. For instance,

>>> t = Instrument('NIRISS')
Parameters :

name : string

Name of the instrument class to return. Case insensitive.

Display Functions

webbpsf.display_PSF(HDUlist_or_filename=None, ext=0, vmin=1e-08, vmax=0.1, scale='log', cmap=<matplotlib.colors.LinearSegmentedColormap instance at 0x10388bcb0>, title=None, imagecrop=None, adjust_for_oversampling=False, normalize='None', crosshairs=False, markcentroid=False, colorbar=True, colorbar_orientation='vertical', pixelscale='PIXELSCL', ax=None, return_ax=False)

Display nicely a PSF from a given HDUlist or filename

This is extensively configurable. In addition to making an attractive display, for interactive usage this function provides a live display of the pixel value at a given (x,y) as you mouse around the image.

Parameters :

HDUlist_or_filename : pyfits.HDUlist or string

FITS file containing image to display.

ext : int

FITS extension. default = 0

vmin, vmax : float

min and max for image display scaling

scale : str

‘linear’ or ‘log’, default is log

cmap : matplotlib.cm.Colormap instance

Colormap to use. Default is matplotlib.cm.jet

ax : matplotlib.Axes instance

Axes to display into.

title : string, optional

imagecrop : float

size of region to display (default is whole image)

normalize : string

set to ‘peak’ to normalize peak intensity =1, or to ‘total’ to normalize total flux=1. Default is no normalization.

adjust_for_oversampling : bool

rescale to conserve surface brightness for oversampled PSFs? (making this True conserves surface brightness but not total flux) default is False, to conserve total flux.

markcentroid : bool

Draw a crosshairs at the image centroid location? Centroiding is computed with the JWST-standard moving box algorithm.

colorbar : bool

Draw a colorbar?

colorbar_orientation : str

either ‘horizontal’ or ‘vertical’; default is vertical.

pixelscale : str or float

if str, interpreted as the FITS keyword name for the pixel scale in arcsec/pixels. if float, used as the pixelscale directly.

webbpsf.display_PSF_difference(HDUlist_or_filename1=None, HDUlist_or_filename2=None, ext1=0, ext2=0, vmax=0.0001, title=None, imagecrop=None, adjust_for_oversampling=False, crosshairs=False, colorbar=True, colorbar_orientation='vertical', print_=False, ax=None, return_ax=False, vmin=None, normalize=False, normalize_to_second=False)

Display nicely the difference of two PSFs from given files

Parameters :

HDUlist_or_filename1,2 : pyfits.HDUlist or string

FITS files containing image to difference

ext1, ext2 : int

FITS extension. default = 0

vmax : float

for the scaling

title : string, optional

imagecrop : float

size of region to display (default is whole image)

normalize : bool

Display (difference image)/(mean image) instead of just the difference image.

adjust_for_oversampling : bool

rescale to conserve surface brightness for oversampled PSFs? (making this True conserves surface brightness but not total flux) default is False, to conserve total flux.

webbpsf.display_profiles(HDUlist_or_filename=None, ext=0, overplot=False)

Metrics of PSF Quality

webbpsf.radial_profile(HDUlist_or_filename=None, ext=0, EE=False, center=None, stddev=False, binsize=None, maxradius=None)

Compute a radial profile of the image.

This computes a discrete radial profile evaluated on the provided binsize. For a version interpolated onto a continuous curve, see measure_radial().

Code taken pretty much directly from pydatatut.pdf

Parameters :

HDUlist_or_filename : string

what it sounds like.

ext : int

Extension in FITS file

EE : bool

Also return encircled energy (EE) curve in addition to radial profile?

center : tuple of floats

Coordinates (x,y) of PSF center, in pixel units. Default is image center.

binsize : float

size of step for profile. Default is pixel size.

stddev : bool

Compute standard deviation in each radial bin, not average?

Returns :

results : tuple

Tuple containing (radius, profile) or (radius, profile, EE) depending on what is requested. The radius gives the center radius of each bin, while the EE is given inside the whole bin so you should use (radius+binsize/2) for the radius of the EE curve if you want to be as precise as possible.

webbpsf.measure_EE(HDUlist_or_filename=None, ext=0, center=None, binsize=None)

Returns a function object which when called returns the Encircled Energy inside a given radius.

Parameters :

HDUlist_or_filename : string

what it sounds like.

ext : int

Extension in FITS file

center : tuple of floats

Coordinates (x,y) of PSF center. Default is image center.

binsize: :

size of step for profile. Default is pixel size.

Returns :

encircled_energy: function :

A function which will return the encircled energy interpolated to any desired radius.

Examples

>>> EE = measure_EE("someimage.fits")
>>> print "The EE at 0.5 arcsec is ", EE(0.5)
webbpsf.measure_fwhm(HDUlist_or_filename=None, ext=0, center=None, level=0.5)

Measure FWHM* by interpolation of the radial profile (* or full width at some other fraction of max.)

Parameters :

HDUlist_or_filename, ext : string, int

Same as above

center : tuple

center to compute around. Default is image center.

level : float

Fraction of max to compute for; default is 0.5 for Half Max. You can also measure widths at other levels e.g. FW at 10% max by setting level=0.1

webbpsf.measure_sharpness(HDUlist_or_filename=None, ext=0)

Compute image sharpness, the sum of pixel squares.

See Makidon et al. JWST-STScI-001157 for a discussion of this image metric and its relationship to noise equivalent pixels.

Parameters :

HDUlist_or_filename, ext : string, int

Same as above

webbpsf.measure_centroid(HDUlist_or_filename=None, ext=0, slice=0, boxsize=20, print_=False, units='pixels', relativeto='origin', **kwargs)

Measure the center of an image via center-of-mass

Parameters :

HDUlist_or_filename, ext : string, int

Same as above

boxsize : int

Half box size for centroid

relativeto : string

either ‘origin’ for relative to pixel (0,0) or ‘center’ for relative to image center. Default is ‘origin’

units : string

either ‘pixels’ for position in pixels or ‘arcsec’ for arcseconds. Relative to the relativeto parameter point in either case.

Returns :

CoM : array_like

[Y, X] coordinates of center of mass.


Documentation last updated on June 27, 2012