mangadap.proc.spectralstack module
Stack some spectra!
Copyright © 2019, SDSS-IV/MaNGA Pipeline Group
- class mangadap.proc.spectralstack.SpectralStack[source]
Bases:
object
Class whose primary function is to stack a set of spectra while treating covariance between spectra.
See
covariance_mode_options()
for available methods of accounting for covariance.Class also approximates the resulting spectral resolution of the stacked data.
Todo
List attributes
- static _check_covariance_shape(covariance_mode, covar, nwave, nspec)[source]
Check that the input covariance object has the correct shape for the given mode.
- Parameters:
covariance_mode (
str
) – Covariance handling mode; seecovariance_mode_options()
.covar (None,
float
,mangadap.util.covariance.Covariance
) – The object to check the type against the covariance handling mode.nwave (
int
) – Number of wavelength channels.nspec (
int
) – Number of spectra.
- Returns:
Flag that the covariance data has the correct shape.
- Return type:
bool
- static _check_covariance_type(covariance_mode, covar, ivar)[source]
Check that the covariance variable has the correct type for the given mode.
- Parameters:
covariance_mode (
str
) – Covariance handling mode; seecovariance_mode_options()
.covar (None,
float
,mangadap.util.covariance.Covariance
) – The object to check the type against the covariance handling mode.ivar (None, object) – Inverse variance. Only check performed is whether or not this is None.
- Returns:
Flag that type is correct.
- Return type:
bool
- static _check_input_sres(sres, nspec)[source]
Check the shape and type of the input spectral resolution.
If the input is a masked array, interpolate over the mask.
- Parameters:
sres (numpy.ndarray, numpy.ma.MaskedArray) – Input spectral resolution data. If None, None is returned.
nspec (
int
) – Expected number of spectra.
- Returns:
The interpolated spectral resolution vectors. If the input
sres
is a single vector, the spectral resolution is repeated for each spectrum. Shape is \((N_{\rm spec}, N_{\rm wave})\). If the inputsres
is None, instead None is returned.- Return type:
- _covar_in_mean()[source]
Compute the covariance in the mean spectrum by propagating the division by the number of pixels through the covariance matrix.
- Returns:
Covariance in the mean spectrum. Returns None if
covar
is None.- Return type:
- static _get_input_mask(flux, ivar=None, mask=None, dtype=<class 'bool'>)[source]
Construct the baseline mask using the input arrays.
- Parameters:
flux (numpy.ndarray, numpy.ma.MaskedArray) – Input flux array.
ivar (numpy.ndarray) – Inverse variance array. Mask excludes any inverse variance values that are not greater than 0.
mask (numpy.ndarray) – Baseline boolean mask array.
dtype (data type) – Output data type for array.
- Returns:
Boolean mask array converted to
dtype
.- Return type:
- _get_stack_mean()[source]
Convert the summed stack to the mean stacked spectra using the internal data.
- Returns:
See the return statement for
stack()
.- Return type:
tuple
- _recalibrate_ivar_figure(ratio, ofile=None)[source]
Make a plot showing the computed inverse variance recalibration based on the covariance matrices.
- Parameters:
ratio (numpy.ndarray) – Ratio between the error with and without covariance.
ofile (
str
, optional) – File for the plot. If None, output to the screen.
- _set_rebin_transfer_matrix(binid, binwgt=None)[source]
Construct the transfer matrix that rebins the spectra.
The shape of the transfer matrix is \((N_{\rm bin} \times N_{\rm spec})\) with the expectation that the spectrum flux array has shape \((N_{\rm spec} \times N_{\rm wave})\).
The binned spectra are calculated by matrix multiplication, \(\mathbf{B} = \mathbf{T} \times \mathbf{F}\) such that the covariance matrix can be calculated as \(\mathbf{C} = \mathbf{T} \times \mathbf{\Sigma} \times \mathbf{T}^{\rm T}\), where \(\mathbf{\Sigma}\) is the covariance matrix in the flux array, \(\mathbf{F}\).
If weighting, the sum of the weights is normalized to the number of points included in the bin.
- Parameters:
binid (numpy.ndarray) – Index, one per spectrum in the flux array, for the binned spectrum. Indices of less than one are ignored.
binwgt (numpy.ndarray, optional) – List of weights for the spectra. If not provided, the weights are uniform.
- _stack_with_covariance(flux, covariance_mode, covar, ivar=None, sres=None)[source]
Stack the spectra and incorporate covariance.
- Parameters:
flux (numpy.ma.MaskedArray) – Flux array.
covariance_mode (
str
) – Covariance handling mode; seecovariance_mode_options()
.covar (None,
float
,mangadap.util.covariance.Covariance
) – The relevant covariance object that must match the needs of the covariance handling mode.ivar (
numpy.ma.MaskedArray
, optional) – The inverse variance array. Must not be None ifcovariance_mode
is ‘channels’ or ‘wavelengths’.sres (numpy.ma.MaskedArray, optional) – 1D or 2D spectral resolution as a function of wavelength for all or each input spectrum. Default is to ignore any spectral resolution data.
- _stack_without_covariance(flux, ivar=None, sres=None, linear_sres=False)[source]
Stack the spectra, ignoring any covariance that may or may not exist.
This method calculates the sum of the flux, flux^2, and determines the number of pixels in the sum. Sets
flux
,fluxsqr
,npix
,ivar
,sres
. The stored data is always based on the sum of the spectra in the stack.- Parameters:
flux (numpy.ma.MaskedArray) – Flux array.
ivar (
numpy.ma.MaskedArray
, optional) – The inverse variance array.sres (numpy.ma.MaskedArray, optional) – 1D or 2D spectral resolution as a function of wavelength for all or each input spectrum. Default is to ignore any spectral resolution data.
linear_sres (
bool
, optional) – Construct the combined resolution as the mean of the input instrumental FWHM (1/sres). If False, the quadratic mean is used (i.e., rms = sqrt(mean(square))); this mode is formally correct for the second moment of the combined LSF (however, note that this not the same as the accuracy of a Gaussian fit to the combined line profile; cf. Law et al. 2020).
- static build_covariance_data(cube, covariance_mode, covariance_par)[source]
Construct the covariance data needed by the specified covariance mode.
- Parameters:
cube (
mangadap.datacube.datacube.DataCube
) – Datacube object to bin.covariance_mode (
str
) – Mode to use. Must be an allowed mode; seecovariance_mode_options()
.covariance_par (
int
,float
,list
) – Parameters required by the specified mode.
- Returns:
None,
float
,mangadap.util.covariance.Covariance
: Covariance data relevant to the selected mode.- Raises:
ValueError – Raised if the covariance mode is not valid or if a covariance matrix cannot be computed for the provided datacube.
TypeError – Raised if the covariance parameters do not have the correct type.
- static covariance_mode_options(par_needed=False)[source]
Return the list of allowed covariance options.
The two parameters
covar_mode
andcovar_par
(seeSpectralStackPar
) set how covariance is accounted for in the stacking procedure. The valid options are:none
: The noise in the stacked spectrum is a nominal propagation of the error assuming no covariance. No parameters needed.calibrate
: Where \(N_{\rm bin}\) is the number of binned spaxels and \(\alpha\) as a provided parameter, the spectral noise is calibrated following:\[n_{\rm calib} = n_{\rm nominal} (1 + \alpha \log\ N_{\rm bin})\]channels
: The noise vector of each stacked spectrum is adjusted based on the mean ratio of the nominal and formally correct calculations of the noise measurements over a number of spectral channels. The channels are drawn from across the full spectral range. The number of channels to use is a defined parameter. The covariance matrix must be provided tostack()
.wavelengths
: Functionally equivalent tochannels
; however, the channels to use are set by a list of provided wavelengths. The covariance matrix must be provided tostack()
.approx_correlation
: Approximate the covariance matrix using a Gaussian description of the correlation between pixels. Seemangadap.datacube.datacube.DataCube.approximate_correlation_matrix()
. The value of \(\sigma\) provides for the Gaussian desciption of \(\rho_{ij}\) in the correlation matrix. The covariance matrix must be provided tostack()
.full
: The full covariance cube is calculated and the noise vectors are constructed using the formally correct calculation. No parameters needed. The covariance matrix must be provided tostack()
.
- Returns:
List of the allowed modes.
- Return type:
list
- static min_max_wave(wave, cz)[source]
Determine the minimum and maximum of all shifted wavelength ranges.
- Parameters:
wave (array-like) – Original wavelengths. Should be 1D.
cz (
float
, array-like) – The redshift of one or more spectra to be removed. Each element is applied to the wavelength vector to determine the maximum wavelength range required to full sample all deshifted spectra.
- Returns:
Two floats with the minimum and maximum wavelengths.
- Return type:
tuple
- Raises:
ValueError – Raised if either
wave
orcz
have the wrong dimensionality.
- static operation_options()[source]
Return the allowed stacking operations.
- Current operations are:
mean
: Construct the mean of the spectrasum
: Construct the spectrum sum.
- Returns:
List of available operations.
- Return type:
list
- static parse_covariance_parameters(mode, par)[source]
Parse the parameters needed for the treatment of the covariance when stacking spectra.
- Parameters:
mode (
str
) – Mode to use. Must be an allowed mode; seecovariance_mode_options()
.par (
str
) – String representation of the parameters for the specified mode.
- Returns:
Parameters parsed from the input string for the designated covariance mode.
- Return type:
float
,list
- Raises:
TypeError – Raised if the input parameter could not be converted to a float as needed by the specified mode.
ValueError – Raised if the mode is not recognized.
- static register(wave, cz, flux, ivar=None, mask=None, sres=None, log=False, base=10.0, keep_range=False, flim=0.5)[source]
Register a set of spectra to the same wavelength range given a set of measured velocities.
Todo
Allow for correction for deredshifting flux.
- Parameters:
wave (numpy.ndarray) – Single wavelength vector for all input spectra. Must be 1D with shape \((N_{\rm wave},)\).
cz (
float
, array-like) – The measured \(cz\) velocities of all spectra or each spectrum individually. The “registration” deshifts all the spectra such that \(cz = 0\).flux (numpy.ndarray, numpy.ma.MaskedArray) – Flux array to register. Must be 2D with shape \((N_{\rm spec},N_{\rm wave})\).
ivar (numpy.ndarray, numpy.ma.MaskedArray, optional) – Inverse variance in the spectrum fluxes. Shape must match
flux
.mask (numpy.ndarray, optional) – Boolean array with values in
flux
to mask. Default assumes nothing is masked. Ifflux
and/orivar
are masked array, this is included in a union with those masks. Shape must matchflux
.sres (numpy.ndarray, numpy.ma.MaskedArray, optional) – 1D or 2D spectral resolution as a function of wavelength for all or each input spectrum. Default is to ignore any spectral resolution data. If a masked array, the masked pixels are interpolated.
log (
bool
, optional) – Flag that the wavelength vector is sampled geometrically in wavelength.base (
float
, optional) – If the wavelength vector is geometrically sampled, this is the base of the logarithmic step.keep_range (
bool
, optional) – When registering the wavelengths of the shifted spectra, keep the identical spectral range as input.flim (
float
, optional) – Mask any pixels in the output flux arrays that are covered by less than this fraction by the input masked flux array.
- Returns:
Returns four arrays: (1) the new wavelength vector common to all spectra; (2) the new masked flux array; (3) the new masked inverse variance array, which will be None if no inverse variances are provided to the function; and (4) the new spectral resolution vectors, which will be None if no spectral resolution vectors are provided.
- Return type:
tuple
- Raises:
ValueError – Raised if the wavelength or velocity offset vectors are not one-dimensional, if the flux array is not two-dimensional, if the inverse variance or mask arrays do not have the same shape as the flux array, or if the number of wavelengths does not match the second axis of the flux array.
- stack(wave, flux, operation='mean', binid=None, binwgt=None, ivar=None, mask=None, sres=None, cz=None, log=False, base=10.0, covariance_mode=None, covar=None, keep_range=False)[source]
Stack a set of spectra.
The method also sets all the returned variables to the internal attributes; however, they are always kept as the sum, not the mean, of the spectra, regardless of the value of
operation
.- Parameters:
wave (numpy.ndarray) – Single wavelength vector for all input spectra.
flux (numpy.ndarray, numpy.ma.MaskedArray) – Spectrum flux values. Shape must be \((N_{\rm spec}, N_{\rm wave})\).
operation (
str
, optional) – Stacking operation to perform. Seeoperation_options()
.binid (numpy.ndarray, optional) – Indices of the bin in which to place each spectrum. Shape is \((N_{\rm spec},)\). If not provided, all the spectra will be combined.
binwgt (numpy.ndarray, optional) – Weights for each of the spectra. All weights must be positive. A binwgt of 0 is considered identical to masking the spectrum. Shape is \((N_{\rm spec},)\). If None, weights are uniform.
ivar (numpy.ndarray, numpy.ma.MaskedArray, optional) – Inverse variance in the spectrum fluxes. Shape must match
flux
.mask (numpy.ndarray, optional) – Binary mask values for the spectrum fluxes (False=unmasked; True=masked). Shape must match
flux
. Default assumes no pixel mask.sres (numpy.ndarray, numpy.ma.MaskedArray, optional) – 1D or 2D spectral resolution as a function of wavelength for all or each input spectrum. Shape must be \((N_{\rm wave},)\) or \((N_{\rm spec},N_{\rm wave})\). If provided as a masked array, masked pixels are replaced with the interpolated resolution at that location.
cz (numpy.ndarray, optional) – Vector with measured \(cz\) velocities used to deshift the spectra before stacking. See
register()
.log (
bool
, optional) – Flag that the wavelength vector is geometrically stepped in wavelength.base (
float
, optional) – If the wavelength vector is geometrically stepped, this is the base of the logarithmic step.covariance_mode (
str
, optional) – Keyword for method to use for dealing with covariance; seecovariance_mode_options()
.covar (
float
,mangadap.util.covariance.Covariance
, optional) – Covariance object to use, which must match the expectation from the covariance mode. Seecovariance_mode_options()
and_check_covariance_type()
.keep_range (
float
, optional) – When registering the wavelengths of the shifted spectra, keep the identical spectral range as input.
- Returns:
Returns seven objects: the wavelength vector, the stacked flux, the standard deviation in the stacked flux, the number of spectra in each stacked pixel, the stacked inverse variance, the stacked spectral resolution, and the stacked covariance.
- Return type:
tuple
- Raises:
ValueError – Raised if the wavelength vector is not one-dimensional, if the flux array is not two-dimensional, if the ancillary arrays (inverse variance, weights, mask) do not have the same shape as the flux array, if any of the weights are negative, or if the number of wavelengths does not match the second axis of the flux array. Also raised if the covariance mode is not recognized; see
covariance_mode_options()
.TypeError – Raised if the object type of
covar
does not match what is expected givencovariance_mode
.NotImplementedError – Raised if the requesting to both velocity register the spectra and correct the error vectors for spatial covariance.
- stack_datacube(cube, binid, par=None)[source]
Wrapper function for
stack()
that accepts a datacube object.Todo
List the required elements of
par
.
- Parameters:
cube (
mangadap.datacube.datacube.DataCube
) – Datacube to stack.binid (
numpy.ndarray
) – Indices of the bin in which to place each spectrum. Shape must be the flattened spatial shape of the datacube; i.e., \((N_{\rm spec},)\).par (
SpectralStackPar
, optional) – Set of parameters used to define how to stack the spectra. Seestack()
. Does not need to be provided on initialization, but a parameter set is required for all the stacking routines.
- Returns:
See the return statement for
stack()
.- Return type:
tuple
- class mangadap.proc.spectralstack.SpectralStackPar(operation='mean', register=False, cz=None, covar_mode='none', covar_par=None)[source]
Bases:
KeywordParSet
Class with parameters used to set how to stack a set of spectra. See
mangadap.par.parset.ParSet
for attributes.Todo
Allow for sigma rejection.
The defined parameters are:
Key
Type
Options
Default
Description
operation
str
mean
,sum
mean
Operation to perform for the stacked spectrum. See
SpectralStack.operation_options()
for the available operation options.register
bool
False
Flag to register the spectra by deshifting them based on their observed cz velocities. This is done before adding them based on a provided prior measurement of the velocities.
cz
ndarray, list
List of measured cz velocities used to register the spectra.
covar_mode
str
calibrate
,approx_correlation
,channels
,wavelengths
,none
,full
none
Describes how to incorporate covariance into the spectral stacking. See
SpectralStack.covariance_mode_options()
for the available options.covar_par
int, float, ndarray, list
The parameter(s) needed to perform a given method of handling the covariance. See
SpectralStack.covariance_mode_options()
for the available options.Warning
Velocity registration is currently not implemented.
- fromheader(hdr)[source]
Copy the information from the header
- Parameters:
hdr (astropy.io.fits.Header) – Header object to read from.
- toheader(hdr)[source]
Copy the information to a header.
- Parameters:
hdr (astropy.io.fits.Header) – Header object to write to.
- Returns:
Edited header object.
- Return type: