Spectrum

class mpdaf.obj.Spectrum(filename=None, hdulist=None, data=None, mask=False, var=None, ext=None, unit=Unit(dimensionless), copy=True, dtype=None, primary_header=None, data_header=None, convert_float64=True, **kwargs)[source]

Bases: mpdaf.obj.arithmetic.ArithmeticMixin, mpdaf.obj.DataArray

Spectrum objects contain 1D arrays of numbers, optionally accompanied by corresponding variances. These numbers represent sample fluxes along a regularly spaced grid of wavelengths.

The spectral pixel values and their variances, if any, are available as arrays[q that can be accessed via properties of the Spectrum object called .data and .var, respectively. These arrays are usually masked arrays, which share a boolean masking array that can be accessed via a property called .mask. In principle, these arrays can also be normal numpy arrays without masks, in which case the .mask property holds the value, numpy.ma.nomask. However non-masked arrays are only supported by a subset of mpdaf functions at this time, so masked arrays should be used where possible.

When a new Spectrum object is created, the data, variance and mask arrays can either be specified as arguments, or the name of a FITS file can be provided to load them from.

Parameters
filenamestring

An optional FITS file name from which to load the spectrum. None by default. This argument is ignored if the data argument is not None.

extint or (int,int) or string or (string,string)

The optional number/name of the data extension or the numbers/names of the data and variance extensions.

wavempdaf.obj.WaveCoord

The wavelength coordinates of the spectrum.

unitstr or astropy.units.Unit

The physical units of the data values. Defaults to astropy.units.dimensionless_unscaled.

datafloat array

An optional 1 dimensional array containing the values of each pixel of the spectrum, stored in ascending order of wavelength (None by default). Where given, this array should be 1 dimensional.

varfloat array

An optional 1 dimensional array containing the estimated variances of each pixel of the spectrum, stored in ascending order of wavelength (None by default).

Attributes
filenamestring

The name of the originating FITS file, if any. Otherwise None.

unitastropy.units.Unit

The physical units of the data values.

primary_headerastropy.io.fits.Header

The FITS primary header instance, if a FITS file was provided.

data_headerastropy.io.fits.Header

The FITS header of the DATA extension.

wavempdaf.obj.WaveCoord

The wavelength coordinates of the spectrum.

Attributes Summary

data

Return data as a numpy.ma.MaskedArray.

dtype

The type of the data.

mask

The shared masking array of the data and variance arrays.

ndim

The number of dimensions in the data and variance arrays : int

shape

The lengths of each of the data axes.

var

Return variance as a numpy.ma.MaskedArray.

Methods Summary

LSF_convolve(self, lsf, size, \*\*kwargs)

Convolve spectrum with LSF.

abmag_band(self, lbda, dlbda)

Compute AB magnitude corresponding to the wavelength band.

abmag_filter(self, lbda, eff)

Compute AB magnitude using array filter.

abmag_filter_name(self, name)

Compute AB magnitude using the filter name.

abs(self[, out])

Return a new object with the absolute value of the data.

add_asym_gaussian(self, lpeak, flux, …[, …])

Add an asymetric gaussian on spectrum in place.

add_gaussian(self, lpeak, flux, fwhm[, …])

Add a gaussian on spectrum in place.

clone(self[, data_init, var_init])

Return a shallow copy with the same header and coordinates.

convolve(self, other[, inplace])

Convolve a Spectrum with a 1D array or another Spectrum, using the discrete convolution equation.

copy(self)

Return a copy of the object.

correlate(self, other[, inplace])

Cross-correlate the spectrum with a other spectrum or an array.

crop(self)

Reduce the size of the array to the smallest sub-array that keeps all unmasked pixels.

fftconvolve(self, other[, inplace])

Convolve a Spectrum with a 1D array or another Spectrum, using the Fourier convolution theorem.

fftconvolve_gauss(self, fwhm[, nsig, unit, …])

Convolve the spectrum with a Gaussian using fft.

filter(self[, kernel])

Perform filtering on the spectrum.

fit_lines(self, redshift, \*\*kwargs)

Use pyplatefit to fit the spectrum lines.

fwhm(self, l0[, cont, spline, unit])

Return the fwhm of a peak.

gauss_asymfit(self, lmin, lmax[, lpeak, …])

Truncate the spectrum and fit it with an asymetric gaussian function.

gauss_dfit(self, lmin, lmax, wratio[, …])

Truncate the spectrum and fit it as a sum of two gaussian functions.

gauss_fit(self, lmin, lmax[, lpeak, flux, …])

Perform a Gaussian fit.

get_data_hdu(self[, name, savemask, …])

Return an ImageHDU corresponding to the DATA extension.

get_dq_hdu(self[, name, header])

Return an ImageHDU corresponding to the DQ (mask) extension.

get_end(self[, unit])

Return the wavelength of the last pixel of the spectrum.

get_range(self[, unit])

Return the wavelength range (Lambda_min, Lambda_max) of the spectrum.

get_start(self[, unit])

Return the wavelength value of the first pixel of the spectrum.

get_stat_hdu(self[, name, header, …])

Return an ImageHDU corresponding to the STAT extension.

get_step(self[, unit])

Return the wavelength step size.

get_wcs_header(self)

Return a FITS header containing coordinate descriptions.

info(self)

Print information.

integrate(self[, lmin, lmax, unit])

Integrate the flux over a specified wavelength range.

interp_mask(self[, spline])

Interpolate masked pixels.

mask_region(self[, lmin, lmax, inside, unit])

Mask spectrum pixels inside or outside a wavelength range, [lmin,lmax].

mask_selection(self, ksel)

Mask selected pixels.

mask_variance(self, threshold)

Mask pixels with a variance above a threshold value.

mean(self[, lmin, lmax, weight, unit])

Compute the mean flux over a specified wavelength range.

median_filter(self[, kernel_size, spline, …])

Perform a median filter on the spectrum.

new_from_obj(obj[, data, var, copy, unit])

Create a new object from another one, copying its attributes.

plot(self[, max, title, noise, snr, lmin, …])

Plot the spectrum.

poly_fit(self, deg[, weight, maxiter, nsig, …])

Perform polynomial fit on normalized spectrum and returns polynomial coefficients.

poly_spec(self, deg[, weight, maxiter, …])

Return a spectrum containing a polynomial fit.

poly_val(self, z)

Update in place the spectrum data from polynomial coefficients.

rebin(self, factor[, margin, inplace])

Combine neighboring pixels to reduce the size of a spectrum by an integer factor.

resample(self, step[, start, shape, unit, …])

Resample a spectrum to have a different wavelength interval.

set_wcs(self[, wcs, wave])

Set the world coordinates (spatial and/or spectral where pertinent).

sqrt(self[, out])

Return a new object with positive data square-rooted, and negative data masked.

subspec(self, lmin[, lmax, unit])

Return the flux at a given wavelength, or the sub-spectrum of a specified wavelength range.

sum(self[, lmin, lmax, weight, unit])

Obtain the sum of the fluxes within a specified wavelength range.

to_ds9(self[, ds9id, newframe, zscale, cmap])

Send the data to ds9 (this will create a copy in memory)

to_spectrum1d(self[, unit_wave])

Return a specutils.Spectrum1D object.

truncate(self[, lmin, lmax, unit])

Truncate the wavelength range of a spectrum in-place.

unmask(self)

Unmask the data (just invalid data (nan,inf) are masked).

wavelet_filter(self[, levels, sigmaCutoff, …])

Perform a wavelet filtering on the spectrum in 1 dimension.

write(self, filename[, savemask, checksum, …])

Save the data to a FITS file.

Attributes Documentation

data

Return data as a numpy.ma.MaskedArray.

The DataArray constructor postpones reading data from FITS files until they are first used. Read the data array here if not already read.

Changes can be made to individual elements of the data property. When simple numeric values or Numpy array elements are assigned to elements of the data property, the values of these elements are updated and become unmasked.

When masked Numpy values or masked-array elements are assigned to elements of the data property, then these change both the values of the data property and the shared mask of the data and var properties.

Completely new arrays can also be assigned to the data property, provided that they have the same shape as before.

dtype

The type of the data.

mask

The shared masking array of the data and variance arrays.

This is a bool array which has the same shape as the data and variance arrays. Setting an element of this array to True, flags the corresponding pixels of the data and variance arrays, so that they don’t take part in subsequent calculations. Reverting this element to False, unmasks the pixel again.

This array can be modified either directly by assignments to elements of this property or by modifying the masks of the .data and .var arrays. An entirely new mask array can also be assigned to this property, provided that it has the same shape as the data array.

ndim

The number of dimensions in the data and variance arrays : int

shape

The lengths of each of the data axes.

var

Return variance as a numpy.ma.MaskedArray.

If variances have been provided for each data pixel, then this property can be used to record those variances. Normally this is a masked array which shares the mask of the data property. However if no variances have been provided, then this property is None.

Variances are typically provided along with the data values in the originating FITS file. Alternatively a variance array can be assigned to this property after the data have been read.

Note that any function that modifies the contents of the data array may need to update the array of variances accordingly. For example, after scaling pixel values by a constant factor c, the variances should be scaled by c**2.

When masked-array values are assigned to elements of the var property, the mask of the new values is assigned to the shared mask of the data and variance properties.

Completely new arrays can also be assigned to the var property. When a masked array is assigned to the var property, its mask is combined with the existing shared mask, rather than replacing it.

Methods Documentation

LSF_convolve(self, lsf, size, **kwargs)[source]

Convolve spectrum with LSF.

Parameters
lsfpython function

mpdaf.MUSE.LSF object or function f describing the LSF.

The first three parameters of the function f must be lbda (wavelength value in A), step (in A) and size (odd integer).

f returns an np.array with shape=2*(size/2)+1 and centered in lbda

Example: from mpdaf.MUSE import LSF

sizeodd int

size of LSF in pixels.

kwargskwargs

it can be used to set function arguments.

Returns
outSpectrum
abmag_band(self, lbda, dlbda)[source]

Compute AB magnitude corresponding to the wavelength band.

Parameters
lbdafloat

Mean wavelength in Angstrom.

dlbdafloat

Width of the wavelength band in Angstrom.

Returns
outfloat, float

Magnitude value and its error

abmag_filter(self, lbda, eff)[source]

Compute AB magnitude using array filter.

Parameters
lbdafloat array

Wavelength values in Angstrom.

efffloat array

Efficiency values.

Returns
outfloat, float

Magnitude value and its error

abmag_filter_name(self, name)[source]

Compute AB magnitude using the filter name.

Parameters
namestring

‘U’, ‘B’, ‘V’, ‘Rc’, ‘Ic’, ‘z’, ‘R-Johnson’, ‘F606W’, ‘F775W’, ‘F814W’, ‘F850LP’

Returns
outfloat, float

Magnitude value and its error

abs(self, out=None)

Return a new object with the absolute value of the data.

Parameters
outmpdaf.obj.DataArray, optional

Array of the same shape as input, into which the output is placed. By default, a new array is created.

add_asym_gaussian(self, lpeak, flux, fwhm_right, fwhm_left, cont=0, peak=False, unit=Unit("Angstrom"))[source]

Add an asymetric gaussian on spectrum in place.

Parameters
lpeakfloat

Gaussian center.

fluxfloat

Integrated gaussian flux or gaussian peak value if peak is True.

fwhm_rightfloat

Gaussian fwhm on the right (red) side

fwhm_leftfloat

Gaussian fwhm on the right (red) side

contfloat

Continuum value.

peakbool

If true, flux contains the gaussian peak value.

unitastropy.units.Unit

Type of the wavelength coordinates. If None, inputs are in pixels.

add_gaussian(self, lpeak, flux, fwhm, cont=0, peak=False, unit=Unit("Angstrom"))[source]

Add a gaussian on spectrum in place.

Parameters
lpeakfloat

Gaussian center.

fluxfloat

Integrated gaussian flux or gaussian peak value if peak is True.

fwhmfloat

Gaussian fwhm.

contfloat

Continuum value.

peakbool

If true, flux contains the gaussian peak value

unitastropy.units.Unit

Type of the wavelength coordinates. If None, inputs are in pixels.

clone(self, data_init=None, var_init=None)

Return a shallow copy with the same header and coordinates.

Optionally fill the cloned array using values returned by provided functions.

Parameters
data_initcallable(), optional

An optional function to use to create the data array (it takes the shape as parameter). For example np.zeros or np.empty can be used. It defaults to None, which results in the data attribute being None.

var_initcallable(), optional

An optional function to use to create the variance array, with the same specifics as data_init. This default to None, which results in the var attribute being assigned None.

convolve(self, other, inplace=False)[source]

Convolve a Spectrum with a 1D array or another Spectrum, using the discrete convolution equation.

This function, which uses the discrete convolution equation, is usually slower than Image.fftconvolve(). However it can be faster when other.data.size is small, and it always uses much less memory, so it is sometimes the only practical choice.

Masked values in self.data and self.var are replaced with zeros before the convolution is performed, but are masked again after the convolution.

If self.var exists, the variances are propagated using the equation:

result.var = self.var (*) other**2

where (*) indicates convolution. This equation can be derived by applying the usual rules of error-propagation to the discrete convolution equation.

The speed of this function scales as O(Nd x No) where Nd=self.data.size and No=other.data.size.

Uses scipy.signal.convolve.

Parameters
otherSpectrum or numpy.ndarray

The 1D array with which to convolve the spectrum in self.data. This can be an array of the same size as self, or it can be a smaller array, such as a small gaussian profile to use to smooth the spectrum.

When other.data contains a symmetric filtering function, such as a gaussian profile, the center of the function should be placed at the center of pixel:

(other.shape - 1) // 2

If other is an MPDAF Spectrum object, note that only its data array is used. Masked values in this array are treated as zero. Any variances found in other.var are ignored.

inplacebool

If False (the default), return the results in a new Spectrum. If True, record the result in self and return that.

Returns
outSpectrum
copy(self)

Return a copy of the object.

correlate(self, other, inplace=False)[source]

Cross-correlate the spectrum with a other spectrum or an array.

Uses scipy.signal.correlate. self and other must have the same size.

Parameters
other1d-array or Spectrum

Second spectrum or 1d-array.

inplacebool

If True, replace the input spectrum with the correlation.

Returns
outSpectrum
crop(self)

Reduce the size of the array to the smallest sub-array that keeps all unmasked pixels.

This removes any margins around the array that only contain masked pixels. If all pixels are masked in the input cube, the data and variance arrays are deleted.

Returns
itemlist of slice

The slices that were used to extract the sub-array.

fftconvolve(self, other, inplace=False)[source]

Convolve a Spectrum with a 1D array or another Spectrum, using the Fourier convolution theorem.

This function, which performs the convolution by multiplying the Fourier transforms of the two arrays, is usually much faster than Spectrum.convolve(), except when other.data.size is small. However it uses much more memory, so Spectrum.convolve() is sometimes a better choice.

Masked values in self.data and self.var are replaced with zeros before the convolution is performed, but they are masked again after the convolution.

If self.var exists, the variances are propagated using the equation:

result.var = self.var (*) other**2

where (*) indicates convolution. This equation can be derived by applying the usual rules of error-propagation to the discrete convolution equation.

The speed of this function scales as O(Nd x log(Nd)) where Nd=self.data.size. This function temporarily allocates a pair of arrays that have the sum of the shapes of self.shape and other.shape, rounded up to a power of two along each axis. This can involve a lot of memory being allocated. For this reason, when other.shape is small, Spectrum.convolve() may be more efficient than Spectrum.fftconvolve().

Uses scipy.signal.fftconvolve.

Parameters
otherSpectrum or numpy.ndarray

The 1D array with which to convolve the spectrum in self.data. This can be an array of the same size as self.data, or it can be a smaller array, such as a small gaussian to use to smooth the spectrum.

When other contains a symmetric filtering function, such as a gaussian profile, the center of the function should be placed at the center of pixel:

(other.shape - 1) // 2

If other is an MPDAF Spectrum object, note that only its data array is used. Masked values in this array are treated as zero. Any variances found in other.var are ignored.

inplacebool

If False (the default), return the results in a new Spectrum. If True, record the result in self and return that.

Returns
outSpectrum
fftconvolve_gauss(self, fwhm, nsig=5, unit=Unit("Angstrom"), inplace=False)[source]

Convolve the spectrum with a Gaussian using fft.

Parameters
fwhmfloat

Gaussian fwhm.

nsigint

Number of standard deviations.

unitastropy.units.Unit

type of the wavelength coordinates

inplacebool

If True, convolve the original spectrum in-place, and return that.

Returns
outSpectrum
filter(self, kernel=<class 'astropy.convolution.kernels.Box1DKernel'>, **parameters)[source]

Perform filtering on the spectrum.

Uses astropy.convolution kernels and convolution.

Parameters
kernelastropy.convolution.Kernel1D

astropy kernel to use (see https://docs.astropy.org/en/stable/convolution/kernels.html#available-kernels)

  • Box1DKernel (default, parameters: width)

  • Gaussan1DKernel (parameters: stddev)

parameterskeywords

parameters to pass to the kernel

Returns
outSpectrum
fit_lines(self, redshift, **kwargs)[source]

Use pyplatefit to fit the spectrum lines.

This method uses the pyplatefit.fit_spec function to perform a line fitting. Refer to its documentation.

Parameters
redshiftfloat

Expected redshift of the spectrum.

**kwargskeyword arguments

Additional arguments passed to the fit_spec function.

Returns
See pyplatefit documentation.
fwhm(self, l0, cont=0, spline=False, unit=Unit("Angstrom"))[source]

Return the fwhm of a peak.

Parameters
l0float

Wavelength value corresponding to the peak position.

unitastropy.units.Unit

Type of the wavelength coordinates. If None, inputs are in pixels.

contint

The continuum [default 0].

splinebool

Linear/spline interpolation to interpolate masked values.

Returns
outfloat
gauss_asymfit(self, lmin, lmax, lpeak=None, flux=None, fwhm=None, cont=None, peak=False, spline=False, weight=True, plot=False, plot_factor=10, ax=None, unit=Unit("Angstrom"))[source]

Truncate the spectrum and fit it with an asymetric gaussian function.

Returns the two gaussian functions (right and left) as mpdaf.obj.Gauss1D objects.

From Johan Richard and Vera Patricio, modified by Jeremy Blaizot.

Parameters
lminfloat or (float,float)

Minimum wavelength value or wavelength range used to initialize the gaussian left value.

lmaxfloat or (float,float)

Maximum wavelength or wavelength range used to initialize the gaussian right value.

lpeakfloat

Input gaussian center. if None it is estimated with the wavelength corresponding to the maximum value in [max(lmin), min(lmax)].

fluxfloat

Integrated gaussian flux or gaussian peak value if peak is True.

fwhmfloat

Input gaussian fwhm, if None it is estimated.

peakbool

If true, flux contains the gaussian peak value .

contfloat

Continuum value, if None it is estimated by the line through points (max(lmin),mean(data[lmin])) and (min(lmax),mean(data[lmax])).

splinebool

Linear/spline interpolation to interpolate masked values.

weightbool

If weight is True, the weight is computed as the inverse of variance.

unitastropy.units.Unit

type of the wavelength coordinates. If None, inputs are in pixels.

plotbool

If True, the resulted fit is plotted.

plot_factordouble

oversampling factor for the overplotted fit

Returns
outmpdaf.obj.Gauss1D, mpdaf.obj.Gauss1D

Left and right Gaussian functions.

gauss_dfit(self, lmin, lmax, wratio, lpeak_1=None, flux_1=None, fratio=1.0, fwhm=None, cont=None, peak=False, spline=False, weight=True, plot=False, plot_factor=10, unit=Unit("Angstrom"))[source]

Truncate the spectrum and fit it as a sum of two gaussian functions.

Returns the two gaussian functions as mpdaf.obj.Gauss1D objects.

From Johan Richard and Vera Patricio.

Parameters
lminfloat or (float,float)

Minimum wavelength value or wavelength range used to initialize the gaussian left value.

lmaxfloat or (float,float)

Maximum wavelength or wavelength range used to initialize the gaussian right value.

wratiofloat

Ratio between the two gaussian centers

lpeak_1float

Input gaussian center of the first gaussian. if None it is estimated with the wavelength corresponding to the maximum value in [max(lmin), min(lmax)]

flux_1float

Integrated gaussian flux or gaussian peak value if peak is True.

fratiofloat

Ratio between the two integrated gaussian fluxes.

fwhmfloat

Input gaussian fwhm, if None it is estimated.

peakbool

If true, flux contains the gaussian peak value .

contfloat

Continuum value, if None it is estimated by the line through points (max(lmin),mean(data[lmin])) and (min(lmax),mean(data[lmax])).

splinebool

Linear/spline interpolation to interpolate masked values.

weightbool

If weight is True, the weight is computed as the inverse of variance.

plotbool

If True, the resulted fit is plotted.

plot_factordouble

oversampling factor for the overplotted fit

unitastropy.units.Unit

Type of the wavelength coordinates. If None, inputs are in pixels.

Returns
outmpdaf.obj.Gauss1D, mpdaf.obj.Gauss1D
gauss_fit(self, lmin, lmax, lpeak=None, flux=None, fwhm=None, cont=None, peak=False, spline=False, weight=True, plot=False, plot_factor=10, unit=Unit("Angstrom"), fix_lpeak=False)[source]

Perform a Gaussian fit.

Uses scipy.optimize.leastsq to minimize the sum of squares.

Parameters
lminfloat or (float,float)

Minimum wavelength value or wavelength range used to initialize the gaussian left value (in angstrom)

lmaxfloat or (float,float)

Maximum wavelength or wavelength range used to initialize the gaussian right value (in angstrom)

lpeakfloat

Input gaussian center (in angstrom), if None it is estimated with the wavelength corresponding to the maximum value in [max(lmin), min(lmax)]

unitastropy.units.Unit

Type of the wavelength coordinates. If None, inputs are in pixels.

fluxfloat

Integrated gaussian flux or gaussian peak value if peak is True.

fwhmfloat

Input gaussian fwhm (in angstrom), if None it is estimated.

peakbool

If true, flux contains the gaussian peak value .

contfloat

Continuum value, if None it is estimated by the line through points (max(lmin),mean(data[lmin])) and (min(lmax),mean(data[lmax])).

splinebool

Linear/spline interpolation to interpolate masked values.

weightbool

If weight is True, the weight is computed as the inverse of variance.

plotbool

If True, the Gaussian is plotted.

plot_factordouble

oversampling factor for the overplotted fit

Returns
outmpdaf.obj.Gauss1D
get_data_hdu(self, name='DATA', savemask='dq', convert_float32=True)

Return an ImageHDU corresponding to the DATA extension.

Parameters
namestr

Extension name, DATA by default.

savemaskstr

If ‘dq’, the mask array is saved in a DQ extension. If ‘nan’, masked data are replaced by nan in a DATA extension. If ‘none’, masked array is not saved.

convert_float32bool

By default float64 arrays are converted to float32, in order to produce smaller files.

Returns
astropy.io.fits.ImageHDU
get_dq_hdu(self, name='DQ', header=None)

Return an ImageHDU corresponding to the DQ (mask) extension.

get_end(self, unit=None)[source]

Return the wavelength of the last pixel of the spectrum.

Parameters
unitastropy.units.Unit

The units of the returned wavelength.

Returns
outfloat

The wavelength of the final pixel of the spectrum.

get_range(self, unit=None)[source]

Return the wavelength range (Lambda_min, Lambda_max) of the spectrum.

Parameters
unitastropy.units.Unit

The units of the returned wavelengths.

Returns
outfloat array

The minimum and maximum wavelengths.

get_start(self, unit=None)[source]

Return the wavelength value of the first pixel of the spectrum.

Parameters
unitastropy.units.Unit

The units of the returned wavelength.

Returns
outfloat

The wavelength of the first pixel of the spectrum.

get_stat_hdu(self, name='STAT', header=None, convert_float32=True)

Return an ImageHDU corresponding to the STAT extension.

Parameters
namestr

Extension name, STAT by default.

headerastropy.io.fits.Header

Fits Header to put in the extension, typically to reuse the same as in the DATA extension. Otherwise it is created with the wcs.

convert_float32bool

By default float64 arrays are converted to float32, in order to produce smaller files.

Returns
astropy.io.fits.ImageHDU
get_step(self, unit=None)[source]

Return the wavelength step size.

Parameters
unitastropy.units.Unit

The units of the returned step-size.

Returns
outfloat

The width of a spectrum pixel.

get_wcs_header(self)

Return a FITS header containing coordinate descriptions.

info(self)

Print information.

integrate(self, lmin=None, lmax=None, unit=Unit("Angstrom"))[source]

Integrate the flux over a specified wavelength range.

The units of the integrated flux depend on the flux units of the spectrum and the wavelength units, as follows:

If the flux units of the spectrum, self.unit, are something like Q per angstrom, Q per nm, or Q per um, then the integrated flux will have the units of Q. For example, if the fluxes have units of 1e-20 erg/cm2/Angstrom/s, then the units of the integration will be 1e-20 erg/cm2/s.

Alternatively, if unit is not None, then the unit of the returned number will be the product of the units in self.unit and unit. For example, if the flux units are counts/s, and unit=u.angstrom, then the integrated flux will have units counts*Angstrom/s.

Finally, if unit is None, then the units of the returned number will be the product of self.unit and the units of the wavelength axis of the spectrum (ie. self.wave.unit).

The result of the integration is returned as an astropy Quantity, which holds the integrated value and its physical units. The units of the returned number can be determined from the .unit attribute of the return value. Alternatively the returned value can be converted to another unit, using the to() method of astropy quantities.

Parameters
lminfloat

The minimum wavelength of the range to be integrated, or None (the default), to select the minimum wavelength of the first pixel of the spectrum. If this is below the minimum wavelength of the spectrum, the integration behaves as though the flux in the first pixel extended down to that wavelength.

If the unit argument is None, lmin is a pixel index, and the wavelength of the center of this pixel is used as the lower wavelength of the integration.

lmaxfloat

The maximum wavelength of the range to be integrated, or None (the default), to select the maximum wavelength of the last pixel of the spectrum. If this is above the maximum wavelength of the spectrum, the integration behaves as though the flux in the last pixel extended up to that wavelength.

If the unit argument is None, lmax is a pixel index, and the wavelength of the center of this pixel is used as the upper wavelength of the integration.

unitastropy.units.Unit

The wavelength units of lmin and lmax, or None to indicate that lmin and lmax are pixel indexes.

Returns
outastropy.units.Quantity, astropy.units.Quantity

The result of the integration and its error, expressed as a floating point number with accompanying units. The integrated value and its physical units can be extracted using the .value and .unit attributes of the returned quantity. The value can also be converted to different units, using the .to() method of the returned objected.

interp_mask(self, spline=False)[source]

Interpolate masked pixels.

Parameters
splinebool

False: linear interpolation (use scipy.interpolate.interp1d), True: spline interpolation (use scipy.interpolate.splrep and scipy.interpolate.splev).

mask_region(self, lmin=None, lmax=None, inside=True, unit=Unit("Angstrom"))[source]

Mask spectrum pixels inside or outside a wavelength range, [lmin,lmax].

Parameters
lminfloat

The minimum wavelength of the range, or None to choose the wavelength of the first pixel in the spectrum.

lmaxfloat

The maximum wavelength of the range, or None to choose the wavelength of the last pixel in the spectrum.

unitastropy.units.Unit

The wavelength units of lmin and lmax. If None, lmin and lmax are assumed to be pixel indexes.

insidebool

If True, pixels inside the range [lmin,lmax] are masked. If False, pixels outside the range [lmin,lmax] are masked.

mask_selection(self, ksel)

Mask selected pixels.

Parameters
kseloutput of np.where

Elements depending on a condition

mask_variance(self, threshold)

Mask pixels with a variance above a threshold value.

Parameters
thresholdfloat

Threshold value.

mean(self, lmin=None, lmax=None, weight=True, unit=Unit("Angstrom"))[source]

Compute the mean flux over a specified wavelength range.

Parameters
lminfloat

The minimum wavelength of the range, or None to choose the wavelength of the first pixel in the spectrum.

lmaxfloat

The maximum wavelength of the range, or None to choose the wavelength of the last pixel in the spectrum.

unitastropy.units.Unit

The wavelength units of lmin and lmax. If None, lmin and lmax are assumed to be pixel indexes.

weightbool

If weight is True, compute the weighted mean, inversely weighting each pixel by its variance.

Returns
out(float, float)

The mean flux and its error.

median_filter(self, kernel_size=1.0, spline=False, unit=Unit("Angstrom"), inplace=False)[source]

Perform a median filter on the spectrum.

Uses scipy.signal.medfilt.

Parameters
kernel_sizefloat

Size of the median filter window.

unitastropy.units.Unit

unit ot the kernel size

inplacebool

If False, return a filtered copy of the spectrum (the default). If True, filter the original spectrum in-place, and return that.

Returns
outSpectrum
classmethod new_from_obj(obj, data=None, var=None, copy=False, unit=None)

Create a new object from another one, copying its attributes.

Parameters
objmpdaf.obj.DataArray

The object to use as the template for the new object. This should be an object based on DataArray, such as an Image, Cube or Spectrum.

dataarray_like, optional

An optional data array, or None to indicate that obj.data should be used. The default is None.

vararray_like, optional

An optional variance array, or None to indicate that obj.var should be used, or False to indicate that the new object should not have any variances. The default is None.

copybool

Copy the data and variance arrays if True (default False).

plot(self, max=None, title=None, noise=False, snr=False, lmin=None, lmax=None, ax=None, stretch='linear', unit=Unit("Angstrom"), noise_kwargs=None, **kwargs)[source]

Plot the spectrum.

By default, the matplotlib drawstyle option is set to ‘steps-mid’. The reason for this is that in MPDAF integer pixel indexes correspond to the centers of their pixels, and the floating point pixel indexes of a pixel extend from half a pixel below the integer central position to half a pixel above it.

Parameters
maxfloat

If max is not None (the default), it should be a floating point value. The plotted data will be renormalized such that the peak in the plot has this value.

titlestring

The title to give the figure (None by default).

noisebool

If noise is True, colored extensions above and below the plotted points indicate the square-root of the variances of each pixel (if any).

snrbool

If snr is True, data/sqrt(var) is plotted.

lminfloat

The minimum wavelength to be plotted, or None (the default) to start the plot from the minimum wavelength in the spectrum.

lmaxfloat

The maximum wavelength to be plotted, or None (the default) to start the plot from the maximum wavelength in the spectrum.

axmatplotlib.Axes

The Axes instance in which the spectrum is drawn, or None (the default), to request that an Axes object be created internally.

unitastropy.units.Unit

The wavelength units of the lmin and lmax arguments, or None to indicate that lmin and lmax are floating point pixel indexes.

noise_kwargsdict

Properties for the noise plot (if noise=True). Default to color='0.75', facecolor='0.75', alpha=0.5.

kwargsdict

kwargs can be used to set properties of the plot such as: line label (for auto legends), linewidth, anitialising, marker face color, etc.

poly_fit(self, deg, weight=True, maxiter=0, nsig=(-3.0, 3.0), verbose=False)[source]

Perform polynomial fit on normalized spectrum and returns polynomial coefficients.

Parameters
degint

Polynomial degree.

weightbool

If weight is True, the weight is computed as the inverse of variance.

maxiterint

Maximum allowed iterations (0)

nsig(float,float)

The low and high rejection factor in std units (-3.0,3.0)

Returns
outndarray, shape.

Polynomial coefficients ordered from low to high.

poly_spec(self, deg, weight=True, maxiter=0, nsig=(-3.0, 3.0), verbose=False)[source]

Return a spectrum containing a polynomial fit.

Parameters
degint

Polynomial degree.

weightbool

If weight is True, the weight is computed as the inverse of variance.

maxiterint

Maximum allowed iterations (0)

nsig(float,float)

The low and high rejection factor in std units (-3.0,3.0)

Returns
outSpectrum
poly_val(self, z)[source]

Update in place the spectrum data from polynomial coefficients.

Uses numpy.poly1d.

Parameters
zarray

The polynomial coefficients, in increasing powers:

data = z0 + z1(lbda-min(lbda))/(max(lbda)-min(lbda)) + … + zn ((lbda-min(lbda))/(max(lbda)-min(lbda)))**n

rebin(self, factor, margin='center', inplace=False)[source]

Combine neighboring pixels to reduce the size of a spectrum by an integer factor.

Each output pixel is the mean of n pixels, where n is the specified reduction factor.

Parameters
factorint

The integer reduction factor by which the spectrum should be shrunk.

marginstring in ‘center’|’right’|’left’|’origin’

When the dimension of the input spectrum is not an integer multiple of the reduction factor, the spectrum is truncated to remove just enough pixels that its length is a multiple of the reduction factor. This sub-spectrum is then rebinned in place of the original spectrum. The margin parameter determines which pixels of the input spectrum are truncated, and which remain.

The options are:
‘origin’ or ‘center’:

The start of the output spectrum is coincident with the start of the input spectrum.

‘center’:

The center of the output spectrum is aligned with the center of the input spectrum, within one pixel.

‘right’:

The end of the output spectrum is coincident with the end of the input spectrum.

inplacebool

If False, return a rebinned copy of the spectrum (the default). If True, rebin the original spectrum in-place, and return that.

Returns
outSpectrum
resample(self, step, start=None, shape=None, unit=Unit("Angstrom"), inplace=False, atten=40.0, cutoff=0.25)[source]

Resample a spectrum to have a different wavelength interval.

Parameters
stepfloat

The new pixel size along the wavelength axis of the spectrum.

startfloat

The wavelength at the center of the first pixel of the resampled spectrum. If None (the default) the center of the first pixel has the same wavelength before and after resampling.

unitastropy.units.Unit

The wavelength units of the step and start arguments. The default is u.angstrom.

shapeint

The dimension of the array of the new spectrum (ie. the number of spectral pixels). If this is not specified, the shape is selected to encompass the wavelength range from the chosen start wavelength to the ending wavelength of the input spectrum.

inplacebool

If False, return a resampled copy of the spectrum (the default). If True, resample the original spectrum in-place, and return that.

attenfloat

The minimum attenuation (dB), of the antialiasing decimation filter at the Nyquist folding frequency of the new pixel size. Larger attenuations suppress aliasing better at the expense of worsened resolution. The default attenuation is 40.0 dB. To disable antialiasing, specify atten=0.0.

cutofffloat

Mask each output pixel of which at least this fraction of the pixel was interpolated from masked input pixels.

Returns
outSpectrum
set_wcs(self, wcs=None, wave=None)

Set the world coordinates (spatial and/or spectral where pertinent).

Parameters
wcsmpdaf.obj.WCS, optional

Spatial world coordinates. This argument is ignored when self is a Spectrum.

wavempdaf.obj.WaveCoord, optional

Spectral wavelength coordinates. This argument is ignored when self is an Image.

sqrt(self, out=None)

Return a new object with positive data square-rooted, and negative data masked.

Parameters
outmpdaf.obj.DataArray, optional

Array of the same shape as input, into which the output is placed. By default, a new array is created.

subspec(self, lmin, lmax=None, unit=Unit("Angstrom"))[source]

Return the flux at a given wavelength, or the sub-spectrum of a specified wavelength range.

A single flux value is returned if the lmax argument is None (the default), or if the wavelengths assigned to the lmin and lmax arguments are both within the same pixel. The value that is returned is the value of the pixel whose wavelength is closest to the wavelength specified by the lmin argument.

Note that is a wavelength range is asked for, a view on the original spectrum is returned and both will be modified at the same time. If you need to modify only the sub-spectrum, you’ll need to copy() it before.

Parameters
lminfloat

The minimum wavelength of a wavelength range, or the wavelength of a single pixel if lmax is None.

lmaxfloat or None

The maximum wavelength of the wavelength range.

unitastropy.units.Unit

The wavelength units of the lmin and lmax arguments. The default is angstroms. If unit is None, then lmin and lmax are interpreted as array indexes within the spectrum.

Returns
outfloat or Spectrum
sum(self, lmin=None, lmax=None, weight=True, unit=Unit("Angstrom"))[source]

Obtain the sum of the fluxes within a specified wavelength range.

Parameters
lminfloat

The minimum wavelength of the range, or None to choose the wavelength of the first pixel in the spectrum.

lmaxfloat

The maximum wavelength of the range, or None to choose the wavelength of the last pixel in the spectrum.

unitastropy.units.Unit

The wavelength units of lmin and lmax. If None, lmin and lmax are assumed to be pixel indexes.

weightbool

If weight is True, compute the weighted sum, inversely weighting each pixel by its variance.

Returns
outfloat, float

The total flux and its error.

to_ds9(self, ds9id=None, newframe=False, zscale=False, cmap='grey')

Send the data to ds9 (this will create a copy in memory)

Parameters
ds9idstr, optional

The DS9 session ID. If ‘None’, a new one will be created. To find your ds9 session ID, open the ds9 menu option File:XPA:Information and look for the XPA_METHOD string, e.g. XPA_METHOD:  86ab2314:60063. You would then calll this function as cube.to_ds9('86ab2314:60063')

newframebool

Send the cube to a new frame or to the current frame?

to_spectrum1d(self, unit_wave=Unit("Angstrom"))[source]

Return a specutils.Spectrum1D object.

truncate(self, lmin=None, lmax=None, unit=Unit("Angstrom"))[source]

Truncate the wavelength range of a spectrum in-place.

Parameters
lminfloat

The minimum wavelength of a wavelength range, or the wavelength of a single pixel if lmax is None.

lmaxfloat or None

The maximum wavelength of the wavelength range.

unitastropy.units.Unit

The wavelength units of the lmin and lmax arguments. The default is angstroms. If unit is None, then lmin and lmax are interpreted as array indexes within the spectrum.

unmask(self)

Unmask the data (just invalid data (nan,inf) are masked).

wavelet_filter(self, levels=9, sigmaCutoff=5.0, epsilon=0.05, inplace=False)[source]

Perform a wavelet filtering on the spectrum in 1 dimension.

Code contributed by Markus Rexroth (EPFL, 2016), and used in https://arxiv.org/abs/1703.09239 (with funding from ERC Advanced Grant LIDA, Swiss National Science Foundation, ERC Starting Grant 336736-CALENDS).

Parameters
levelsint

Highest scale level.

sigmaCutofffloat

Cleaning threshold. By default 5 for a 5 sigma cleaning in wavelet space.

epsilonfloat in ]0,1[

Residual criterion used to perform the cleaning

inplacebool

If False, return a filtered copy of the spectrum (the default). If True, filter the original spectrum in-place, and return that.

Returns
outSpectrum
write(self, filename, savemask='dq', checksum=False, convert_float32=True)

Save the data to a FITS file.

Overwrite the file if it exists.

Parameters
filenamestr

The FITS filename.

savemaskstr

If ‘dq’, the mask array is saved in a DQ extension If ‘nan’, masked data are replaced by nan in the DATA extension If ‘none’, masked array is not saved

checksumbool

If True, adds both DATASUM and CHECKSUM cards to the headers of all HDU’s written to the file.

convert_float32bool

By default float64 arrays are converted to float32, in order to produce smaller files.