CCDData class¶
Getting started¶
Getting data in¶
The tools in ccdproc
accept only CCDData
objects, a
subclass of NDData
.
Creating a CCDData
object from any array-like data is easy:
>>> import numpy as np
>>> import ccdproc
>>> ccd = ccdproc.CCDData(np.arange(10), unit="adu")
Note that behind the scenes, NDData
creates references to
(not copies of) your data when possible, so modifying the data in ccd
will
modify the underlying data.
You are required to provide a unit for your data. The most frequently used
units for these objects are likely to be adu
, photon
and electron
, which
can be set either by providing the string name of the unit (as in the example
above) or from unit objects:
>>> from astropy import units as u
>>> ccd_photon = ccdproc.CCDData([1, 2, 3], unit=u.photon)
>>> ccd_electron = ccdproc.CCDData([1, 2, 3], unit="electron")
If you prefer not to use the unit functionality then use the special unit
u.dimensionless_unscaled
when you create your CCDData
images:
>>> ccd_unitless = ccdproc.CCDData(np.zeros((10, 10)),
... unit=u.dimensionless_unscaled)
A CCDData
object can also be initialized from a FITS file:
>>> ccd = ccdproc.CCDData.read('my_file.fits', unit="adu")
If there is a unit in the FITS file (in the BUNIT
keyword), that will be
used, but a unit explicitly provided in read
will override any unit in the
FITS file.
There is no restriction at all on what the unit can be – any unit in
astropy.units
or that you create yourself will work.
In addition, the user can specify the extension in a FITS file to use:
>>> ccd = ccdproc.CCDData.read('my_file.fits', hdu=1, unit="adu")
If hdu
is not specified, it will assume the data is in the primary
extension. If there is no data in the primary extension, the first extension
with data will be used.
Metadata¶
When initializing from a FITS file, the header
property is initialized using
the header of the FITS file. Metadata is optional, and can be provided by any
dictionary or dict-like object:
>>> ccd_simple = ccdproc.CCDData(np.arange(10), unit="adu")
>>> my_meta = {'observer': 'Edwin Hubble', 'exposure': 30.0}
>>> ccd_simple.header = my_meta # or use ccd_simple.meta = my_meta
Whether the metadata is case sensitive or not depends on how it is initialized. A FITS header, for example, is not case sensitive, but a python dictionary is.
Getting data out¶
A CCDData
object behaves like a numpy array (masked if the
CCDData
mask is set) in expressions, and the underlying
data (ignoring any mask) is accessed through data
attribute:
>>> ccd_masked = ccdproc.CCDData([1, 2, 3], unit="adu", mask=[0, 0, 1])
>>> 2 * np.ones(3) * ccd_masked # one return value will be masked
masked_array(data = [2.0 4.0 --],
mask = [False False True],
fill_value = 1e+20)
>>> 2 * np.ones(3) * ccd_masked.data # ignores the mask
array([ 2., 4., 6.])
You can force conversion to a numpy array with:
>>> np.asarray(ccd_masked)
array([1, 2, 3])
>>> np.ma.array(ccd_masked.data, mask=ccd_masked.mask)
masked_array(data = [1 2 --],
mask = [False False True],
fill_value = 999999)
A method for converting a CCDData
object to a FITS HDU list
is also available. It converts the metadata to a FITS header:
>>> hdulist = ccd_masked.to_hdu()
You can also write directly to a FITS file:
>>> ccd_masked.write('my_image.fits')
Masks and flags¶
Although not required when a CCDData
image is created you
can also specify a mask and/or flags.
A mask is a boolean array the same size as the data in which a value of
True
indicates that a particular pixel should be masked, i.e. not be
included in arithmetic operations or aggregation.
Flags are one or more additional arrays (of any type) whose shape matches the
shape of the data. For more details on setting flags see
astropy.nddata.NDData
.
WCS¶
The wcs
attribute of CCDData
object can be set two ways.
- If the
CCDData
object is created from a FITS file that has WCS keywords in the header, thewcs
attribute is set to aastropy.wcs.WCS
object using the information in the FITS header. - The WCS can also be provided when the
CCDData
object is constructed with thewcs
argument.
Either way, the wcs
attribute is kept up to date if the
CCDData
image is trimmed.
Uncertainty¶
Pixel-by-pixel uncertainty can be calculated for you:
>>> data = np.random.normal(size=(10, 10), loc=1.0, scale=0.1)
>>> ccd = ccdproc.CCDData(data, unit="electron")
>>> ccd_new = ccdproc.create_deviation(ccd, readnoise=5 * u.electron)
See Gain correct and create deviation image for more details.
You can also set the uncertainty directly, either by creating a
StdDevUncertainty
object first:
>>> from astropy.nddata.nduncertainty import StdDevUncertainty
>>> uncertainty = 0.1 * ccd.data # can be any array whose shape matches the data
>>> my_uncertainty = StdDevUncertainty(uncertainty)
>>> ccd.uncertainty = my_uncertainty
or by providing a ndarray
with the same shape as the data:
>>> ccd.uncertainty = 0.1 * ccd.data
INFO: array provided for uncertainty; assuming it is a StdDevUncertainty. [...]
In this case the uncertainty is assumed to be
StdDevUncertainty
. Using StdDevUncertainty
is required to enable error propagation in CCDData
If you want access to the underlying uncertainty use its .array
attribute:
>>> ccd.uncertainty.array
array(...)
Arithmetic with images¶
Methods are provided to perform arithmetic operations with a
CCDData
image and a number, an astropy
Quantity
(a number with units) or another
CCDData
image.
Using these methods propagates errors correctly (if the errors are
uncorrelated), take care of any necessary unit conversions, and apply masks
appropriately. Note that the metadata of the result is not set if the operation
is between two CCDData
objects.
>>> result = ccd.multiply(0.2 * u.adu)
>>> uncertainty_ratio = result.uncertainty.array[0, 0]/ccd.uncertainty.array[0, 0]
>>> round(uncertainty_ratio, 5)
0.2
>>> result.unit
Unit("adu electron")
Note
In most cases you should use the functions described in Reduction toolbox to perform common operations like scaling by gain or doing dark or sky subtraction. Those functions try to construct a sensible header for the result and provide a mechanism for logging the action of the function in the header.
The arithmetic operators *
, /
, +
and -
are not overridden.
Note
If two images have different WCS values, the wcs on the first
CCDData
object will be used for the resultant object.