Cube object¶
MPDAF Cube objects contain a series of images at neighboring wavelengths. These
are stored as a 3D array of pixel values. The wavelengths of the images increase
along the first dimension of this array, while the Y and X axes of the images
lie along the second and third dimensions. In addition to the data array, Cubes
contain a WaveCoord
object that describes the wavelength
scale of the spectral axis of the cube, and a WCS
object that
describes the spatial coordinates of the spatial axes. Optionally a 3D array of
variances can be provided to give the statistical uncertainties of the
fluxes. Finally a mask array is provided for indicating bad pixels.
The fluxes and their variances are stored in numpy masked arrays, so virtually all numpy and scipy functions can be applied to them.
Preliminary imports:
In [1]: import numpy as np
In [2]: import matplotlib.pyplot as plt
In [3]: import astropy.units as u
In [4]: from mpdaf.obj import Cube, WCS, WaveCoord
Cube creation¶
There are two common ways to obtain a Cube
object:
- A cube can be created from a user-provided 3D array of pixel values, or from both a 3D array of pixel values, and a corresponding 3D array of variances. These can be simple numpy arrays, or they can be numpy masked arrays in which bad pixels have been masked. For example:
In [5]: wcs1 = WCS(crval=0, cdelt=0.2)
In [6]: wave1 = WaveCoord(cdelt=1.25, crval=4000.0, cunit=u.angstrom)
# numpy data array
In [7]: MyData = np.ones((400, 30, 30))
# cube filled with MyData
In [8]: cube = Cube(data=MyData, wcs=wcs1, wave=wave1) # cube 400X30x30 filled with data
In [9]: cube.info()
[INFO] 400 x 30 x 30 Cube (no name)
[INFO] .data(400 x 30 x 30) (no unit), no noise
[INFO] spatial coord (pix): min:(0.0,0.0) max:(5.8,5.8) step:(0.2,0.2) rot:-0.0 deg
[INFO] wavelength: min:4000.00 max:4498.75 step:1.25 Angstrom
- Alternatively, a cube can be read from a FITS file. In this case the flux and variance values are read from specific extensions:
# data and variance arrays read from the file (extension DATA and STAT)
In [10]: obj1 = Cube('obj/CUBE.fits')
In [11]: obj1.info()
[INFO] 1595 x 10 x 20 Cube (obj/CUBE.fits)
[INFO] .data(1595 x 10 x 20) (1e-40 erg / (Angstrom cm2 s)), .var(1595 x 10 x 20)
[INFO] center:(-30:00:00.4494,01:20:00.4376) size:(2.000",4.000") step:(0.200",0.200") rot:-0.0 deg frame:FK5
[INFO] wavelength: min:7300.00 max:9292.50 step:1.25 Angstrom
By default, if a FITS file has more than one extension, then it is expected to have a ‘DATA’ extension that contains the pixel data, and possibly a ‘STAT’ extension that contains the corresponding variances. If the file doesn’t contain extensions of these names, the “ext=” keyword can be used to indicate the appropriate extension or extensions.
The spatial and spectral world-coordinates of a Cube are recorded within
WCS
and WaveCoord
objects, respectively. When a cube
is read from a FITS file, these are automatically generated, based on FITS
header keywords. Alternatively, when an cube is extracted from another cube,
they are derived from the world coordinate information of the original cube.
In the above example, information about a cube is printed using the info
method. In this case the information is as follows:
- The cube has dimensions, 1595 x 10 x 20, which means that there are 1595 images of different wavelengths, each having 10 x 20 spatial pixels.
- In addition to the pixel array (.data(1595 x 10 x 20)), there is also a variance array of the same dimensions (.var(1595 x 10 x 20)).
- The flux units of the pixels are 10-20 erg/s/cm2/Angstrom.
- The center of the field of view is at DEC: -30° 0’ 0.45” and RA: 1°20‘0.437”.
- The size of the field of view is 2 arcsec x 4 arcsec. The pixel dimensions are 0.2 arcsec x 0.2 arcsec.
- The rotation angle of the field is 0°, which means that North is along the positive Y axis.
- The wavelength range is 7300-9292.5 Angstrom with a step of 1.25 Angstrom
The indexing of the cube arrays follows the Python conventions for indexing a 3D array. The pixel in the bottom-lower-left corner is referenced as [0,0,0] and the pixel [k,p,q] refers to the horizontal position q, the vertical position p, and the spectral position k, as follows:
(see Spectrum, Image and Cube format for more information).
The following example computes a reconstructed white-light image and displays it. The white-light image is obtained by summing each spatial pixel of the cube along the wavelength axis. This converts the 3D cube into a 2D image. The cube in this examples contains an observation of a single galaxy.
In [12]: ima1 = obj1.sum(axis=0)
In [13]: plt.figure(figsize=(8, 4))
Out[13]: <Figure size 800x400 with 0 Axes>
In [14]: ima1.plot(scale='arcsinh', colorbar='v')