Spectrum, Image and Cube format

Attributes

Spectrum, Image and Cube objects all have the following items, where O denotes the name of the object:

Component Description
O.filename A FITS filename if the data were loaded from a FITS file
O.primary_header A FITS primary header instance
O.wcs World coordinate spatial information (WCS object)
O.wave World coordinate spectral information (WaveCoord object)
O.shape An array of the dimensions of the cube
O.data A numpy masked array of pixel values
O.data_header A FITS data header instance
O.unit The physical units of the data values
O.dtype The data-type of the data array (int, float)
O.var An optional numpy masked array of pixel variances
O.mask An array of the masked state of each pixel
O.ndim The number of dimensions in the data, variance and mask arrays

Masked arrays are arrays that can have missing or invalid entries. The numpy.ma module provides a nearly work-alike replacement for numpy that supports data arrays with masks. See the DataArray documentation for more details.

When an object is constructed from a MUSE FITS file, O.data will contain the DATA extension, O.var will contain the STAT extension and O.mask will contain the DQ extension if it exists. A DQ extension contains the pixel data quality. By default all bad pixels are masked. But it is possible for the user to create his mask by using the Python interface for euro3D convention.

Indexing

The format of each numpy array follows the indexing used by Python to handle 2D and 3D arrays. For an MPDAF cube C, the pixel in the bottom-lower-left corner is C[0,0,0] and pixel C[k,p,q] refers to the horizontal position q, the vertical position p, and the spectral position k, as follows:

_images/gridcube.jpg

In total, this cube C contains nq pixels in the horizontal direction, np pixels in the vertical direction and nk channels in the spectral direction. In the cube, each numpy masked array has 3 dimensions, array[k,p,q], where k is the spectral axis, and p and q are the spatial axes:

In [1]: from mpdaf.obj import Cube

# data array is read from the file (extension number 0)
In [2]: cube = Cube(filename='../data/sdetect/minicube.fits')

# The 3 dimensions of the cube:
In [3]: cube.shape
Out[3]: (3681, 40, 40)

In [4]: cube.data.shape
Out[4]: (3681, 40, 40)

In [5]: cube.var.shape
Out[5]: (3681, 40, 40)

In [6]: cube.mask.shape
Out[6]: (3681, 40, 40)

Cube[k,p,q] returns the value of Cube.data[k,p,q]:

In [7]: cube[3659, 8, 28]
Out[7]: 493.97263

Similarly Cube[k,p,q] = value sets the value of Cube.data[k,p,q], and Cube[k1:k2,p1:p2,q1:q2] = array sets the corresponding subset of Cube.data. Finally Cube[k1:k2,p1:p2,q1:q2] returns a sub-cube, as demonstrated in the following example:

In [8]: cube.info()
[INFO] 3681 x 40 x 40 Cube (../data/sdetect/minicube.fits)
[INFO] .data(3681 x 40 x 40) (1e-20 erg / (Angstrom cm2 s)), .var(3681 x 40 x 40)
[INFO] center:(10:27:56.3962,04:13:25.3588) size in arcsec:(8.000,8.000) step in arcsec:(0.200,0.200) rot:-0.0 deg
[INFO] wavelength: min:4749.89 max:9349.89 step:1.25 Angstrom

In [9]: cube[3000:4000,10:20,25:40].info()
[WARNING] The wavelength coordinates and data have different dimensions. Modifying the shape of the WaveCoord object
[INFO] 681 x 10 x 15 Cube (../data/sdetect/minicube.fits)
[INFO] .data(681 x 10 x 15) (1e-20 erg / (Angstrom cm2 s)), .var(681 x 10 x 15)
[INFO] center:(10:27:55.3962,04:13:25.1893) size in arcsec:(2.000,3.000) step in arcsec:(0.200,0.200) rot:-0.0 deg
[INFO] wavelength: min:8499.89 max:9349.89 step:1.25 Angstrom

Likewise, Cube[k,:,:] returns an Image, as demonstrated below:

In [10]: ima1 = cube[1000, :, :]

In [11]: plt.figure()
Out[11]: <matplotlib.figure.Figure at 0x7fb7475b8a10>

In [12]: ima1.plot(colorbar='v', title = '$\lambda$ = %.1f (%s)' %(cube.wave.coord(1000), cube.wave.unit))
Out[12]: <matplotlib.image.AxesImage at 0x7fb745e1f410>

In [13]: ima2 = cube[3000, :, :]

In [14]: plt.figure()
Out[14]: <matplotlib.figure.Figure at 0x7fb735eebfd0>

In [15]: ima2.plot(colorbar='v', title = '$\lambda$ = %.1f (%s)' %(cube.wave.coord(3000), cube.wave.unit))
Out[15]: <matplotlib.image.AxesImage at 0x7fb735841810>

In [16]: plt.figure()
Out[16]: <matplotlib.figure.Figure at 0x7fb735e73890>

In [17]: ima2[5:25, 15:35].plot(colorbar='v',title = 'Zoom $\lambda$ = %.1f (%s)' %(cube.wave.coord(3000), cube.wave.unit))
Out[17]: <matplotlib.image.AxesImage at 0x7fb73610c3d0>
_images/ObjFormatIma1.png _images/ObjFormatIma2.png _images/ObjFormatZommIma2.png

In the Image objects extracted from the cube, Image[p1:p2,q1:q2] returns a sub-image, Image[p,q] returns the value of pixel (p,q), Image[p,q] = value sets value in Image.data[p,q], and Image[p1:p2,q1:q2] = array sets the corresponding part of Image.data.

Finally, Cube[:,p,q] returns a Spectrum:

In [18]: spe = cube[:, 8, 28]

In [19]: import astropy.units as u

In [20]: from mpdaf.obj import deg2sexa

In [21]: coord_sky = cube.wcs.pix2sky([8, 28], unit=u.deg)

In [22]: dec, ra = deg2sexa(coord_sky)[0]

In [23]: plt.figure()
Out[23]: <matplotlib.figure.Figure at 0x7fb744378910>

In [24]: spe.plot(title = 'Spectrum ra=%s dec=%s' %(ra, dec))

In [25]: plt.figure()
Out[25]: <matplotlib.figure.Figure at 0x7fb735f298d0>

In [26]: spe[1640:2440].plot(title = 'Zoom Spectrum ra=%s dec=%s' %(ra, dec))
_images/ObjFormatSpe.png _images/ObjFormatZoomSpe.png

Getters and setters

Cube.get_step, Image.get_step and Spectrum.get_step return the world-coordinate separations between pixels along each axis of a cube, image, or spectrum, respectively:

In [27]: cube.get_step(unit_wave=u.nm, unit_wcs=u.deg)
Out[27]: array([  1.25000000e-01,   5.55555556e-05,   5.55555556e-05])

In [28]: ima1.get_step(unit=u.deg)
Out[28]: array([  5.55555556e-05,   5.55555556e-05])

In [29]: spe.get_step(unit=u.angstrom)
Out[29]: 1.25

Cube.get_range returns the range of wavelengths, declinations and right ascensions in a cube. Similarly, Image.get_range returns the range of declinations and right ascensions in an image, and Spectrum.get_range returns the range of wavelengths in a spectrum, as demonstrated below:

In [30]: cube.get_range(unit_wave=u.nm, unit_wcs=u.deg)
Out[30]: 
array([ 474.9890625 ,   10.46458228,   63.35455983,  934.9890625 ,
         10.46674895,   63.35676316])

In [31]: ima1.get_range(unit=u.deg)
Out[31]: array([ 10.46458228,  63.35455983,  10.46674895,  63.35676316])

In [32]: spe.get_range(unit=u.angstrom)
Out[32]: array([ 4749.890625,  9349.890625])

The get_start and get_end methods of cube, image and spectrum objects, return the world-coordinate values of the first and last pixels of each axis:

In [33]: print cube.get_start(unit_wave=u.nm, unit_wcs=u.deg), cube.get_end(unit_wave=u.nm, unit_wcs=u.deg)
[ 474.9890625    10.46458228   63.35676315] [ 934.9890625    10.46674895   63.35455983]

In [34]: print ima1.get_start(unit=u.deg), ima2.get_end(unit=u.deg)
[ 10.46458228  63.35676315] [ 10.46674895  63.35455983]

In [35]: print spe.get_start(unit=u.angstrom), spe.get_end(unit=u.angstrom)
4749.890625 9349.890625

Note that when the declination axis is rotated away from the vertical axis of the image, the coordinates returned by get_start and get_end are not the minimum and maximum coordinate values within the image, so they differ from the values returned by get_range.

Cube.get_rot and Image.get_rot return the rotation angle of the declination axis to the vertical axis of the images within these objects:

In [36]: cube.get_rot(unit=u.deg)
Out[36]: -0.0

In [37]: ima1.get_rot(unit=u.rad)
Out[37]: -0.0

Updated flux and variance values can be assigned directly to the O.data and O.var attributes, respectively. Similarly, elements of the data can be masked or unmasked by assigning True or False values to the corresponding elements of the O.mask attribute. Changes to the spatial world coordinates must be performed using the set_wcs method: