mangadap.proc.sasuke module¶
Implements an emission-line fitting class that largely wraps pPXF.
Copyright © 2019, SDSS-IV/MaNGA Pipeline Group
- class mangadap.proc.sasuke.Sasuke(bitmask, loggers=None, quiet=False)[source]¶
Bases:
mangadap.proc.spectralfitting.EmissionLineFit
Use ninja skills and pPXF to fit emission lines.
https://en.wikipedia.org/wiki/Sasuke_Uchiha
Effectively, nothing happens during the instantiation of this object. Instead, a typical usage of the class when fitting a set of emission lines would be:
# Read the emission-line database emldb = EmissionLineDB.from_key('ELPMILES') # Instantiate the emission-line fitter el_fitter = Sasuke(EmissionLineModelBitMask()) # Fit the spectra model_wave, model_flux, model_eml_flux, model_mask, model_fit_par, \ model_eml_par = el_fitter.fit(...)
See
fit()
for the arguments to the main fitting function.The class inherits attributes from
mangadap.proc.spectralfitting.EmissionLineFit
(fit_type, bitmask, par, fit_method). Other attributes are defined upon instantiation and set to None. This isn’t necessary (attributes can be defined elsewhere in the class methods), but it provides a collation of the class attributes for reference.- Parameters
bitmask (
BitMask
) –BitMask object use to flag fit results. This must be provided and should typically be an instantiation of
EmissionLineModelBitMask
; however, it can be any object withBitMask
as its base class. The flags set within the main fitting function (Sasuke.fit()
) are:DIDNOTUSE, INSUFFICIENT_DATA, FIT_FAILED, NEAR_BOUND, NO_FIT,
mangadap.proc.ppxffit.PPXFFit.rej_flag
(PPXF_REJECT), MIN_SIGMA, BAD_SIGMA, and MAXITER.
Also the DAP-specific calling function (
Sasuke.fit_SpatiallyBinnedSpectra()
) will also assign bits NON_POSITIVE_CONTINUUM during the equivalent width measurements (seemangadap.spectralfitting.EmissionLineFit.measure_equivalent_width()
).loggers (
list
, optional) – List of logging.Logger objects to log progress; ignored if quiet=True. Logging is done usingmangadap.util.log.log_output()
. Default is no logging. This can be reset in some methods.quiet (
bool
, optional) – Suppress all terminal and logging output. Default is False.
- loggers¶
See initialization argument.
- Type
list
- quiet¶
Suppress all terminal and logging output.
- Type
bool
- obj_wave¶
Wavelength vector for object spectra. Shape is \((N_{\rm pix},)\).
- Type
numpy.ndarray
- obj_flux¶
Object spectra to fit. Shape is \((N_{\rm spec},N_{\rm pix})\).
- Type
- obj_ferr¶
\(1\sigma\) errors in object spectra. Shape is \((N_{\rm spec},N_{\rm pix})\).
- Type
- obj_sres¶
Spectral resolution array for object spectra. Shape is \((N_{\rm spec},N_{\rm pix})\).
- Type
numpy.ndarray
- nobj¶
Number of object spectra (i.e., \(N_{\rm spec}\)).
- Type
int
- npix_obj¶
Number of pixels in each object spectrum (i.e., \(N_{\rm pix}\)).
- Type
int
- input_obj_mask¶
A copy of the input mask array of the object spectra (boolean array). Shape is \((N_{\rm spec},N_{\rm pix})\).
- Type
numpy.ndarray
- obj_to_fit¶
Flag to fit each object spectrum. Instantiating by fully masked spectra in
input_obj_mask
. Shape is \((N_{\rm spec},)\).- Type
numpy.ndarray
- input_cz¶
Input redshifts (in km/s) for each spectrum. Shape is \((N_{\rm spec},)\).
- Type
numpy.ndarray
- velscale¶
Velocity scale (km/s) of the object spectra.
- Type
float
- tpl_wave¶
Wavelength vector for template spectra. Shape is \((N_{\rm pix,t},)\).
- Type
numpy.ndarray
- tpl_flux¶
Template spectra to use in fit. Shape is \((N_{\rm tpl},N_{\rm pix,t})\).
- Type
numpy.ndarray
- tpl_to_use¶
Set of flags used to select the templates used to fit each spectrum. Shape is \((N_{\rm spec},N_{\rm tpl})\).
- Type
numpy.ndarray
- nstpl¶
Number of stellar templates.
- Type
int
- ntpl¶
Total number of templates (gas + stars).
- Type
int
- npix_tpl¶
Number of pixels in template spectra (i.e., \(N_{\rm pix,t}\)).
- Type
int
- tpl_npad¶
Nearest length for FFT, \(N_{\rm pad}\).
- Type
int
- tpl_rfft¶
The complex array with the real FFT of the template spectra. Shape is \((N_{\rm tpl}, N_{\rm pad}/2 + 1)\).
- Type
numpy.ndarray
- matched_resolution¶
The spectral resolution of the templates is matched to that of the galaxy data. WARNING: This functionality needs to be checked in relation to the gas templates!
- Type
bool
- velscale_ratio¶
The integer ratio between the velocity scale of the pixel in the galaxy data to that of the template data.
- Type
int
- emldb¶
Emission-line database that is parsed to construct the emission-line templates (see
EmissionLineTemplates
).
- neml¶
Number of emission lines in the database.
- Type
int
- fit_eml¶
Boolean array setting which emission lines are fit. Shape is \((N_{\rm eml,)\).
- Type
numpy.ndarray
- eml_tpli¶
Integer array with the template that includes each emission line. Shape is \((N_{\rm eml,)\).
- Type
numpy.ndarray
- eml_compi¶
Integer array with the kinematic component that includes each emission line. Shape is \((N_{\rm eml,)\).
- Type
numpy.ndarray
- ncomp¶
Total number of kinematic components to fit.
- Type
int
- tpl_comp¶
The integer kinematic component associated with each template. Shape is \((N_{\rm tpl},)\).
- Type
numpy.ndarray
- tpl_vgrp¶
The integer velocity group associated with each template. Shape is \((N_{\rm tpl},)\).
- Type
numpy.ndarray
- tpl_sgrp¶
The integer sigma group associated with each template. Shape is \((N_{\rm tpl},)\).
- Type
numpy.ndarray
- comp_moments¶
Number of moments for each component. Moments with negative numbers have fixed kinematics. Shape is \((N_{\rm comp},)\).
- Type
numpy.ndarray
- comp_start_kin¶
Array of lists where each list provdes the starting parameters for the kinematics of each component. Shape is \((N_{\rm comp},)\).
- Type
numpy.ndarray
- npar_kin¶
The total number of kinematic parameters, which is just the sum of the absolute value of
comp_moments
, \(N_{\rm kin}\).- Type
int
- nfree_kin¶
The total number of free kinematic parameters.
- Type
int
- velocity_limits¶
The upper and lower velocity limits imposed by pPXF. See
mangadap.proc.ppxffit.PPXFFit.losvd_limits()
.- Type
numpy.ndarray
- sigma_limits¶
The upper and lower velocity dispersion limits imposed by pPXF. See
mangadap.proc.ppxffit.PPXFFit.losvd_limits()
.- Type
numpy.ndarray
- gh_limits¶
The upper and lower limits on all higher order Gauss-Hermite moments imposed by pPXF. See
mangadap.proc.ppxffit.PPXFFit.losvd_limits()
.- Type
numpy.ndarray
- bias¶
pPXF bias parameter. (Currently irrelevant because gas is currently always fit with moments=2.)
- Type
float
- degree¶
pPXF degree parameter setting the degree of the additive polynomials to use, \(o_{\rm add}\).
- Type
int
- mdegree¶
pPXF mdegree parameter setting the degree of the multiplicative polynomials to use, \(o_{\rm mult}\).
- Type
int
- reddening¶
pPXF reddening parameter setting the initial \(E(B-V)\) to fit, based on a Calzetti law.
- Type
float
- reject_boxcar¶
Size of the boxcar to use when rejecting fit outliers.
- Type
int
- spectrum_start¶
Array with the starting index of the pixel in the object spectra to fit (inclusive). Shape is \((N_{\rm spec},)\).
- Type
numpy.ndarray
- spectrum_end¶
Array with the ending index of the pixel in the object spectra to fit (exclusive). Shape is \((N_{\rm spec},)\).
- Type
numpy.ndarray
- dof¶
Degrees of freedom in the fit.
- Type
int
Todo
velocity_limits
,sigma_limits
, andgh_limits
are not passed to ppxf during the fit. They are expected to match what’s in the pPXF code. Should change the code so that they are passed using pPXF’s BOUND keyword.
- _is_near_bounds(kin, kin_inp, vel_indx, sig_indx, lbound, ubound, tol_frac=0.01, fill_value=- 999.0)[source]¶
Check if the kinematics were fit and near the imposed limits.
Any value that is close to
fill_value
is flagged as not having been fit, presumably because too much/all of the data near the emission line is masked.The definition of “near” is that the velocity and higher moments cannot be closer than the provided fraction of the total width to the boundary. For the velocity dispersion, the fraction is done in log space.
- Parameters
kin (numpy.ndarray) – Best-fitting kinematic parameters. Shape should be \((N_{\rm spec},N_{\rm kin})\).
kin_inp (numpy.ndarray) – The initial guesses for the best-fitting kinematics. This is needed because the velocity limits are set relative to the input guess. Shape should be \((N_{\rm spec},N_{\rm kin})\).
vel_index (numpy.ndarray) – Boolean array setting if the parameter is a velocity. This is needed because the velocity limits are set relative to the input guess.
sig_index (numpy.ndarray) – Boolean array setting if the parameter is a velocity dispersion. This is needed because the proximity to the bounds is logarithmic for the velocity dispersion.
lbound (numpy.ndarray) – Lower bound on each parameter.
ubound (numpy.ndarray) – Upper bound on each parameter.
tol_frac (float) – (Optional) The fractional tolerance for classifying the returned parameter and near the boundary.
- Returns
Three boolean arrays that flag if (1) the parameter is near either boundary, (2) the parameter is near the lower boundary, and (3) the parameter is near the fill_value. The second array is important in case the fitted parameter is a velocity dispersion and that it’s near the lower boundary because it has hit the pixel sampling limit.
- Return type
numpy.ndarray
- _reset_to_fill_value(model_eml_par, fill_value=- 999.0)[source]¶
Reset all the emission lines masked as having insufficient data to fit to the provided fill value.
Note
The error values are currently reset to 0. by
mangadap.proc.emissionlinemodel.EmissionLineModel.
- _save_results(etpl, start, end, flux, ferr, spec_to_fit, model_flux, model_eml_flux, model_wgts, model_wgts_err, model_addcoef, model_multcoef, model_reddening, model_kin_inp, model_kin, model_kin_err, model_mask, model_fit_par, model_eml_par, fill_value=- 999.0)[source]¶
Save and assess the results of the ppxf fits.
The results are saved both as the direct output from pPXF and after parsing the data into the results for each emission line. The function also assesses the data to set any necessary flags. Much of this is the same as
mangadap.proc.ppxffit.PPXFFit._save_results
.- Parameters
etpl (
mangadap.proc.emissionlinetemplates.EmissionLineTemplates
) – The object used to construct and hold the emission-line templates.start (numpy.ndarray) – Starting pixel of the fit for each spectrum.
end (numpy.ndarray) – End pixel (exclusive) of the fit for each spectrum.
flux (numpy.ndarray) – The fitted spectra.
ferr (numpy.ndarray) – Errors in the fitted spectra.
spec_to_fit (numpy.ndarray) – Boolean flags for spectra with attempted fits.
model_flux (numpy.ndarray) – The best-fit model spectra.
model_eml_flux (numpy.ndarray) – The best-fit emission-line-only spectra.
model_wgts (numpy.ndarray) – The weights applied to each template.
model_wgts_err (numpy.ndarray) – The weight errors.
model_addcoeff (numpy.ndarray) – The coefficients of any fitted additive polynomials.
model_multcoeff (numpy.ndarray) – The coefficients of any fitted multiplicative polynomials.
model_reddening (numpy.ndarray) – The best-fitting \(E(B-V)\) if reddening was included in the fits.
model_kin_inp (numpy.ndarray) – The starting-guess input kinematics.
model_kin (numpy.ndarray) – The best-fitting kinematics.
model_kin_err (numpy.ndarray) – The errors in the kinematics.
model_mask (numpy.ndarray) – The array of bitmask values associated with spectral fitting flags. Shape is \((N_{\rm spec}, N_{\rm pix})\).
model_fit_par (
SasukeFitDataTable
) – Data table with the specific fit parameters. This is parsed into the main datatable by_save_results()
.model_eml_par (
EmissionLineDataTable
) – Data table with the parameters measured for each emission line.fill_value (scalar-like, optional) – The value used to fill masked measurements.
- Returns
Five numpy arrays are returned:
the best-fitting model spectra,
the best-fitting emission-line only spectra,
the bitmask values,
the per spectrum ppxf result, and
the per spectrum emission-line parameters.
The first 3 returned objects are of type numpy.ndarray and have shape \((N_{\rm spec}, N_{\rm pix})\); the last two are numpy.recarray instances with shape \((N_{\rm spec},)\).
- Return type
tuple
- _validate_dispersions(model_eml_par, rng=[0, 400])[source]¶
Check that the corrected velocity dispersion are in the provided range.
(FB) More flagging of kinematic paramters. DOES THIS APPLY TO EMISSION LINES? (KBW): Yes. The corrected dispsersion must be larger than 0 and less than 400 km/s. It’s easy to turn this off or change the limits.
(KBW) Currently not called…
- Parameters
model_eml_par (
EmissionLineDataTable
) – Data table with the parameters measured for each emission line.rng (
list
, optional) – Two-element list with the minimum and maximum allowed corrected velocity dispersion. Measurements outside this range are flagged asBAD_SIGMA
.
- Returns
Returns the input record array with any additional flags.
- Return type
- static construct_continuum_models(emission_lines, stpl_wave, stpl_flux, obj_wave, obj_flux_shape, model_fit_par, select=None, redshift_only=False, deredshift=False, dispersion_corrections=None, dvtol=1e-10)[source]¶
Construct the continuum models using the provided set of model parameters.
This is a wrapper for
PPXFModel
, and very similar tomangadap.proc.ppxffit.PPXFFit.contruct_models()
.The input velocities are expected to be cz, not “ppxf” (pixelized) velocities.
If redshift_only is true, the provided dispersion is set to 1e-9 km/s, which is numerically identical to 0 (i.e., just shifting the spectrum) in the tested applications. However, beware that this is a HARDCODED number.
- static deconstruct_bins_options()[source]¶
Methods allowed for deconstructing bins are:
ignore
: Do not deconstruct bins (default)nearest
: Use the on-sky proximity of spaxel and bin centers to do associate spaxels to bins (e.g., as done in DR15 hybrid binning).binid
: Use the known association between spaxels and bins based on which spaxels were used to construct the binned spectra.- Returns
List of allowed options.
- Return type
list
- static etpl_line_sigma_options()[source]¶
Keyword options for the mode used to set the standard deviation of the emission lines the the emission-line templates. Possible modes are:
default
: In order of precedence:Use the spectral resolution of the observed spectra, or
use the spectral resolution of the stellar templates, or
adopt a FWHM of 2 pixels and calculate the resolution assuming a Gaussian line profile.
zero
: The width of all lines is set to 0. The functionfrom
mangadap.util.lineprofiles
used to construct the line must be able to produce the line if the standard deviation is 0!offset
: Apply a quadrature offset of thewavelength-dependent trend of the instrumental dispersion resulting from the
default
mode such that the minimum instrumental dispersion is 0. Note that the minimum dispersion can be set to something other than 0 using either the etpl_line_sigma_min parameter inSasukePar
or the keyword argument inSasuke.fit()
. If the minimum dispersion is 0, the function frommangadap.util.lineprofiles
used to construct the line must be able to produce a line with a standard deviation of 0!
- Returns
List of allowed options.
- Return type
list
- fit(emission_lines, obj_wave, obj_flux, obj_ferr=None, obj_mask=None, obj_sres=None, guess_redshift=None, guess_dispersion=None, reject_boxcar=None, stpl_wave=None, stpl_flux=None, stpl_sres=None, stpl_to_use=None, stellar_kinematics=None, etpl_sinst_mode='default', etpl_sinst_min=0, remapid=None, remap_flux=None, remap_ferr=None, remap_mask=None, remap_sres=None, remap_skyx=None, remap_skyy=None, obj_skyx=None, obj_skyy=None, velscale_ratio=None, matched_resolution=True, waverange=None, bias=None, degree=- 1, mdegree=0, reddening=None, max_velocity_range=400.0, alias_window=None, dvtol=1e-10, loggers=None, quiet=False, plot=False, sigma_rej=3.0, ensemble=True)[source]¶
Fit a set of emission lines using pPXF to all provided spectra.
The input flux arrays are expected to be ordered with spectra along rows. I.e., to plot the first spectrum in the array, one would do:
from matplotlib import pyplot pyplot.plot(obj_wave, obj_flux[0,:]) pyplot.show()
The function will fit the spectra with or without any provided stellar templates.
The body of this function mostly deals with checking the input and setting up the template and object data for use with pPXF.
The function is meant to be general, but has been largely tested with MaNGA spectra.
If the spectral resolution is not matched between the templates and the object spectra, the provided stellar_kinematics are expected to include any resolution difference; i.e., it is not the astrophysical velocity dispersion.
Todo
guess_redshift and guess_dispersion are not actually optional. The function will fail without them!
Allow for moments != 2.
Allow for fixed components to be set from emission-line database
Allow for bounds to be set from emission-line database
Allow to ignore emission-line database tying and just tie all kinematics
Should probably put the emission-line construction in a separate function.
- Parameters
emission_lines (
mangadap.par.emissionlinedb.EmissionLineDB
) – Emission-line database that is parsed to construct the emission-line templates to fit (seemangadap.proc.emissionelinetemplates.EmissionLineTemplates
).obj_wave (numpy.ndarray) – Wavelength vector for object spectra. Shape is (\(N_{\rm pix}\),).
obj_flux (numpy.ndarray) – Object spectra to fit. Can be provided as a numpy.ma.MaskedArray. Shape is (\(N_{\rm spec},N_{\rm pix}\)).
obj_ferr (numpy.ndarray, optional) – \(1\sigma\) errors in object spectra. Can be provided as a numpy.ma.MaskedArray. Shape is (\(N_{\rm spec},N_{\rm pix}\)). If None, the quadrature sum of the fit residuals is used as the fit metric instead of chi-square (all errors are set to unity).
obj_mask (numpy.ndarray,
mangadap.util.pixelmask.SpectralPixelMask
, optional) – Boolean array ormangadap.util.pixelmask.SpectralPixelMask
instance used to censor regions of the spectra to ignore during fitting.obj_sres (numpy.ndarray, optional) – The spectral resolution of the object data. Can be a single vector for all spectra or one vector per object spectrum.
guess_redshift (array-like, optional) – The starting redshift guess. Can provide a single value for all spectra or one value per spectrum.
guess_dispersion (array-like, optional) – The starting velocity dispersion guess. Can provide a single value for all spectra or one value per spectrum.
reject_boxcar (
int
, optional) – Size of the boxcar to use when rejecting fit outliers. If None, no rejection iteration is performed.stpl_wave (numpy.ndarray, optional) – Wavelength vector for stellar template spectra. Shape is (\(N_{{\rm pix},\ast}\),).
stpl_flux (numpy.ndarray, optional) – Stellar template spectra to use in fit. Shape is (\(N_{{\rm tpl},\ast},N_{{\rm pix},\ast}\)).
stpl_sres (numpy.ndarray, optional) – Spectral resolution of the stellar template spectra. All templates are assumed to have the same spectral resolution. Shape is (\(N_{{\rm pix},\ast}\),). This is also used to set the spectral resolution of the emission-line templates. If not provided, the emission-line templates are constructed with an LSF that has a disperison of 1 pixel.
stpl_to_use (numpy.ndarray, optional) – Set of flags used to select the stellar templates used to fit each spectrum. Shape is (\(N_{\rm spec},N_{{\rm tpl},\ast}\)).
stellar_kinematics (numpy.ndarray, optional) – The kinematics to use for the stellar component. If the spectral resolution of the templates is different from the galaxy data, the provided stellar velocity dispersions must include the (assumed wavelength-independent) quadrature difference in the template and galaxy instrumental resolutions. The stellar kinematics are fixed for all calls made to ppxf. Shape is (\(N_{\rm spec},N_{{\rm kin},\ast}\)). The shape of this array determines the number of moments assigned to the stellar component.
etpl_sinst_mode (
str
, optional) – Mode used to set the instrumental dispersion of the emission-line templates; seeetpl_line_sigma_options()
for the options. Default mode is default, imagine that.etpl_sinst_min (
float
, optional) – Minimum allowed instrumental dispersion value. If the mode of constructing the instrumental dispersion of the templates results in values that are below this value, the whole function is offset in quadrature such that no value goes below the minimum. Default is 0 km/s.remapid (numpy.ndarray, optional) – The index of the input spectrum associated with each remapped spectrum. This is ignored if remapping is not done. If remapping spectra are provided, and this is None, coordinates must be provided and this vector is constructed by associating each remapped spectrum to the input spectra just based by proximity. If provided, this also associates each remapped spectrum to the spectrum to use for the stellar kinematics. Shape is (\(N_{\rm remap}\),).
remap_flux (numpy.ndarray, optional) – Ultimately the spectra to be modeled, if provided. The nominal usage is, e.g., if the object spectra (obj_flux) are binned spaxels, this can be flux arrays of all the spaxels that were used to construct the binned data. The function would then use a first fit to the binned spectra and a matching between the on-sky positions (see remap_skyx, remap_skyy, obj_skyx, obj_skyy) to ultimately fit the individual spaxel data. See
mangadap.contrib.xjmc.emline_fitter_with_ppxf()
. These spectra must be sampled identically to the object spectra (velscale) and have the same wavelength vector (obj_wave). Can be provided as a numpy.ma.MaskedArray. Shape is (\(N_{\rm remap},N_{\rm pix}\)).remap_ferr (numpy.ndarray, optional) – The 1-\(\sigma\) error in remap_flux. Can be provided as a numpy.ma.MaskedArray. Shape is (\(N_{\rm remap},N_{\rm pix}\)).
remap_mask (numpy.ndarray,
mangadap.util.pixelmask.SpectralPixelMask
, optional) – Boolean array ormangadap.util.pixelmask.SpectralPixelMask
instance used to censor regions of the remapping spectra to ignore during fitting.remap_sres (numpy.ndarray, optional) – The spectral resolution of the remapping spectra. Can be a single vector for all spectra or one vector per object spectrum.
remap_skyx (numpy.ndarray, optional) – On-sky x position of the remapping spectra. Shape is (\(N_{\rm remap}\),).
remap_skyy (numpy.ndarray, optional) – On-sky y position of the remapping spectra. Shape is (\(N_{\rm remap}\),).
obj_skyx (numpy.ndarray, optional) – On-sky x position of the object spectra. Shape is (\(N_{\rm spec}\),).
obj_skyy (numpy.ndarray, optional) – On-sky y position of the object spectra. Shape is (\(N_{\rm spec}\),).
velscale_ratio (
int
, optional) – The integer ratio between the velocity scale of the pixel in the galaxy data to that of the template data. If None, set to unity.matched_resolution (
bool
, optional) – The spectral resolution of the templates is matched to that of the galaxy data. WARNING: This functionality is never used!waverange (array-like, optional) – Lower and upper wavelength limits to include in the fit. This can be a two-element vector to apply the same limits to all spectra, or a N-spec x 2 array with wavelength ranges for each spectrum to be fit. Default is to use as much of the spectrum as possible.
bias (
float
, optional) – pPXF bias parameter. (Currently irrelevant because gas is currently always fit with moments=2.)degree (
int
, optional) – pPXF degree parameter setting the degree of the additive polynomials to use, \(o_{\rm add}\). Default is no polynomial included.mdegree (
int
, optional) – pPXF mdegree parameter setting the degree of the multiplicative polynomials to use, \(o_{\rm mult}\). Default is no polynomial included.reddening (
float
, optional) – pPXF reddening parameter used to fit \(E(B-V)\) assuming a Calzetti law. Cannot be fit simultaneously with multiplicative polynomial. Must be larger than 0 to start. Default is not fit.max_velocity_range (
float
, optional) – Maximum range (+/-) expected for the fitted velocities in km/s. Default is 400 km/s.alias_window (
float
, optional) – The window to mask to avoid aliasing near the edges of the spectral range in km/s. Default is six times max_velocity_range.dvtol (
float
, optional) – The velocity scale of the template spectra and object spectrum must be smaller than this tolerance. Default is 1e-10.plot (
bool
, optional) – Show the automatically generated pPXF fit plots for each iterations of each spectrum.loggers (
list
, optional) – List of logging.Logger objects to log progress; ignored if quiet=True. Logging is done usingmangadap.util.log.log_output()
.quiet (
bool
, optional) – Suppress all terminal and logging output.
- Returns
(1) the wavelength vector for the model spectra (should be the same as obj_wave),
the best-fitting model spectra,
the best-fitting emission-line only spectra,
the bitmask values,
the per spectrum ppxf result, and
the per spectrum emission-line parameters.
The first object is a numpy.ndarray instance with shape \((N_{\rm pix},)\), the next 3 objects are numpy.ndarray instances with shape \((N_{\rm spec}, N_{\rm pix})\), and the last two are numpy.recarray instances with shape \((N_{\rm spec},)\).
- Return type
The function returns 6 arrays
- Raises
ValueError – Raised if the length of the spectra, errors, or mask does not match the length of the wavelength array; raised if the wavelength, redshift, or dispersion arrays are not 1D vectors; and raised if the number of redshifts or dispersions is not a single value or the same as the number of input spectra.
- fit_SpatiallyBinnedSpectra(binned_spectra, par=None, loggers=None, quiet=False, debug=False)[source]¶
This DAP-specific function interprets the DAP-specific classes and constructs call(s) to the general
fit()
function to fit the spectra.The format of the calling sequence is dictated by the DAP interface: Any emission-line fitter (e.g.
mangadap.proc.elric.Elric
) must have the same function interface in order to be called within themangadap.proc.emissionlinemodel.EmissionLineModel
object. Because this is a DAP-specific function, it should not declare anything to self.The returned emission-line “baseline” array is set to be the difference between the best-fitting stellar continuum (passed to this function as par[‘stellar_continuum’]) and the best-fitting stellar continuum from the combined emission-line and stellar spectra produced by
fit()
. In this way, the best-fitting model for each spectrum is:best_fit_model = par['stellar_continuum']['FLUX'].data \ + model_eml_flux + model_eml_base
where model_eml_flux and model_eml_base are the 2nd and 3rd returned objects, respectively. These are written to the output DAP model LOGCUBE file in the extensions EMLINE and EMLINE_BASE, respectively.
- Parameters
binned_spectra (
mangadap.proc.spatiallybinnedspectra.SpatiallyBinnedSpectra
) – Spectra to fit.par (
SasukePar
) – Parameters provided from the DAP to the general Sasuke fitting algorithm (fit()
). Althought technically optional given that it is a keyword parameter, theSasukePar
parameter must be provided for proper execution of the function.loggers (
list
, optional) – List of logging.Logger objects to log progress; ignored ifquiet=True
. Logging is done usingmangadap.util.log.log_output()
. Default is no logging.quiet (
bool
, optional) – Suppress all terminal and logging output.
- Returns
The function returns:
wavelength vector of the fitted models, which should match the binned spectra wavelength vector (binned_spectra[‘WAVE’].data),
the best-fitting emission-line model model spectra,
the best-fitting emission-line baseline (see the description above),
the bitmask values,
the per spectrum ppxf result, and
the per spectrum emission-line parameters.
The first object is a numpy.ndarray instance with shape \((N_{\rm pix},)\) , the next 3 objects are numpy.ndarray instances with shape \((N_{\rm spec}, N_{\rm pix})\), and the last two are numpy.recarray instances with shape \((N_{\rm spec},)\).
- Return type
tuple
- static fit_figures_of_merit(model_fit_par)[source]¶
Parse the model fit parameters into only the set of fit-quality metrics per fitted spectrum.
This is primarily used to set the output data for the MAPS file.
- Parameters
model_fit_par (numpy.recarray) – The array with the fit parameters produced by
Sasuke.fit()
. Must have the datatype defined bySasuke._per_firt_dtype()
.- Returns
a list of names to assign each row (length is NFOM) and a 2D array with the figures of merit for each spectrum (shape is NSPEC x NFOM).
- Return type
Two objects are returned
- get_stellar_templates(par, cube, z=0.0, loggers=None, quiet=False)[source]¶
Return the stellar template library.
If fitting a different set of templates, the spectral resolution of the new templates must be matched to the galaxy data and the velocity dispersion used by the fit must be astrophysical (corrected for any resolution difference).
- Returns
Returns the template library, a flag if the resolution of the templates has been matched to the galaxy data, and the velocity sampling compared to the galaxy data.
- Return type
TemplateLibrary, bool, int
- class mangadap.proc.sasuke.SasukeFitDataTable(ntpl, nadd, nmult, nkin, mask_dtype, shape=None)[source]¶
Bases:
mangadap.util.datatable.DataTable
Holds results specific to how Sasuke fits the spectrum.
In
fit()
,BINID
andBINID_INDEX
are the same. They are not infit_SpatiallyBinnedSpectra()
.- Parameters
ntpl (
int
) – Number of spectral templates.nadd (
int
) – Number of coefficients in any additive polynomial included in the fit. Can be 0.nmult (
int
) – Number of coefficients in any multiplicative polynomial included in the fit. Can be 0.nkin (
int
) – Number of kinematic parameters (e.g., 2 for V and sigma)mask_dtype (
type
) – The data type used for the maskbits (e.g., numpy.int16). Typically this would be set bymangadap.util.bitmask.BitMask.minimum_dtype()
.shape (
int
,tuple
, optional) – The shape of the initial array. If None, the data array will not be instantiated; useinit()
to initialize the data array after instantiation.
- class mangadap.proc.sasuke.SasukePar(stellar_continuum=None, emission_lines=None, continuum_templates=None, etpl_line_sigma_mode=None, etpl_line_sigma_min=None, velscale_ratio=None, guess_redshift=None, guess_dispersion=None, minimum_snr=None, deconstruct_bins=None, pixelmask=None, reject_boxcar=None, bias=None, moments=None, degree=None, mdegree=None, reddening=None)[source]¶
Bases:
mangadap.par.parset.KeywordParSet
Hold the parameers necessary to run the Sasuke emission-line fitter.
This is the object that gets passed to
Sasuke.fit_SpatiallyBinnedSpectra()
. In the DAP, it is instantiated bymangadap.proc.emissionlinemodel.available_emission_line_modeling_methods()
and some of its components are filled bymangdap.proc.emissionlinemodel.EmissionLineModel._fill_method_par()
.The continuum templates can either be the string keyword used to construct the template library, or the template library itself. If the former, the
Sasuke.fit_SpatiallyBinnedSpectra()
method will construct the template library for later callback.When instantiated, the
mangadap.par.parset.KeywordParSet
objects test that the input objects match the provided dtypes. See documentation formangadap.par.parset.ParSet
for the list of attributes and exceptions raised.The defined parameters are:
Key
Type
Options
Default
Description
stellar_continuum
StellarContinuumModel
The result of the previous fit to the stellar continuum.
emission_lines
EmissionLineDB
Emission-line database with the details of the lines to be fit.
continuum_templates
str, TemplateLibrary
The new continuum template library to use during the emission-line fit.
etpl_line_sigma_mode
str
default
,zero
,offset
default
Mode used to set the instrumental dispersion of the emission-line templates. Mode options are explated by
Sasuke.etpl_line_sigma_options()
. Default isdefault
.etpl_line_sigma_min
int, float
0.0
Impose a minimum emission-line sigma by offsetting the nominal trend, in quadrature, to have this minimum value. Default is 0.
velscale_ratio
int
1
Integer ratio between the width of the spectral pixels in the galaxy and template spectra. If
velscale_ratio = 2
, the template pixels are half as wide as the galaxy pixels.guess_redshift
ndarray, list, int, float
Single or per-spectrum redshift to use as the initial velocity guess.
guess_dispersion
ndarray, list, int, float
Single or per-spectrum velocity dispersion to use as the initial guess.
minimum_snr
int, float
0.0
Minimum S/N of spectrum to fit (ignored when deconstructing bins).
deconstruct_bins
str
ignore
,nearest
,binid
ignore
Method to use for deconstructing binned spectra into individual spaxels for emission-line fitting. See
Sasuke.deconstruct_bin_options()
.pixelmask
SpectralPixelMask
Mask to apply to all spectra being fit.
reject_boxcar
int
Size of the boxcar to use when rejecting fit outliers.
bias
int, float
pPXF bias parameter. (Irrelevant because gas is currently always fit with moments=2.)
moments
int
2
pPXF moments parameter. (Irrelevant because gas is currently always fit with moments=2.)
degree
int
-1
pPXF degree parameter setting the degree of the additive polynomials to use.
mdegree
int
0
pPXF mdegree parameter setting the degree of the multiplicative polynomials to use.
reddening
int, float
pPXF reddening parameter setting the initial \(E(B-V)\) to fit, based on a Calzetti law.