mangadap.proc.ppxffit module
Implements a wrapper class for pPXF.
Todo
Allow new iteration mode that iteratively fits until the velocity is not up against one of the +/- 2000 km/s limits? Could be useful for poor redshift guesses.
Copyright © 2019, SDSS-IV/MaNGA Pipeline Group
- class mangadap.proc.ppxffit.PPXFFit(bitmask, par=None)[source]
Bases:
StellarKinematicsFit
Use pPXF to measure the stellar kinematics. Although it can also fit the composition and emission lines, for now just force it to be a
StellarKinematicsFit
objec.- bitmask
Bitmask to use for masking spectral pixels and parameter values. For
fit_SpatiallyBinnedSpectra()
, must have bit for ‘LOW_SNR’. Forfit()
must have bits for ‘TPL_PIXELS’, ‘TRUNCATED’, ‘PPXF_REJECT’, ‘LARGE_CHI2’, ‘LARGE_RESID’, ‘INSUFFICIENT_DATA’, ‘FIT_FAILED’, ‘NEAR_BOUND’.- Type:
- _fit_all_spectra(templates, templates_rfft, tpl_to_use, plot=False, plot_file_root=None)[source]
Fit all spectra provided.
Get an initial fit
Reject
Mask and smooth templates and objects
Fit ratio
Mask and smooth
Fit ratio
Fit unmasked with fixed kinematics to get models (with emission lines?)
- _fit_dispersion_correction(templates, templates_rfft, result, baseline_dispersion=None)[source]
- Calculate the dispersion correction:
Construct the optimized, redshifted template without the convolution with the best-fitting LOSVD.
Convolve it with the resolution difference in the data.
Use pPXF to fit the matched version to the unmatched version; this should be the dispersion correction.
Todo
Decide how to deal with regions below 2-pixel resolution in synthetic spectrum being fit. How does this alter the correction?
- _fit_global_spectrum(obj_to_include=None, plot=False)[source]
- Fit the global spectrum. This:
Sets the base-level good pixel mask for the fits to the individual bins
Gives the template weights for the global template
Todo
Only include spectra above a given S/N in global spectrum?
Allow for a number of iterations as input.
- _is_near_bounds(result, guess_velocity, tol_frac=0.01)[source]
Check if the fitted kinematics are near the imposed limits.
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.
- _nominal_dispersion_correction(obj_sres, gpm, cz)[source]
Calculate the dispersion corrections as the quadrature difference between the spectral resolution of the template and object spectra. Returns a masked array!
- _run_fit_iteration(obj_flux, obj_ferr, start, end, base_velocity, tpl_flux, tpl_rfft, guess_kin, fix_kinematics=False, obj_to_fit=None, tpl_to_use=None, degree=None, mdegree=None, dof=None, weight_errors=False, plot=False)[source]
Fit all the object spectra in obj_flux.
Todo
Calculate DOF in this function?
Explicitly set the bounds to use instead of using the pPXF defaults?
- Parameters:
obj_flux (numpy.ma.MaskedArray) – Size is \(N_{\rm spec}\times N_{\rm chan}\), object spectra.
obj_ferr (numpy.ma.MaskedArray) – Size is \(N_{\rm spec}\times N_{\rm chan}\), object errors
start (array) – Size is \(N_{\rm spec}\), starting pixel for each spectrum
end (array) – Size is \(N_{\rm spec}\), ending pixel (+1) for each spectrum
base_velocity (array) – Size is \(N_{\rm spec}\), base velocity offset between each object spectrum and the template spectra.
tpl_flux (array) – Size is \(N_{\rm tpl}\times N_{\rm tpl chan}\), template spectra
tpl_rfft (array) – Size is \(N_{\rm tpl}\times N_{\rm tpl pad}\), real FFT of the template spectra
guess_kin (array) – Initial guess for kinematics. Size is \(N_{\rm spec}\times N_{\rm moments}\).
fix_kinematics (bool) – (Optional) Flag to fix the kinematics to the input values during the fit.
obj_to_fit (array) – (Optional) Size is \(N_{\rm spec}\), boolean flag to fit object spectrum
tpl_to_use (array) – (Optional) Size is \(N_{\rm spec}\times N_{\rm tpl}\), boolean flag to use a template for the fit to each object spectrum
plot (bool) – (Optional) Produce the default ppxf fit plot.
degree (int) – (Optional) Additive polynomial order. Default is to use the internal attribute
degree
.mdegree (int) – (Optional) Multiplicative polynomial order. Default is to use the internal attribute
mdegree
.dof (int) – (Optional) Number of degrees of freedom in the fit. Default is to use the internal attribute
dof
.
- Returns:
Array with \(N_{\rm spec}\) instances of
PPXFFitResult
.- Return type:
numpy.ndarray
- _validate_kinematics(model_mask, model_par)[source]
Validate the returned kinematics.
- Checks:
corrected velocity dispersion must be in the range 50-400 km/s
- static check_objects(obj_wave, obj_flux, obj_ferr=None, obj_sres=None)[source]
Check that the input object data is valid for use with pPXFFit.
- Parameters:
obj_wave (numpy.ndarray) – 1D vector of object wavelengths in angstroms. Does NOT need to be same as the template wavelengths.
obj_flux (numpy.ndarray) – \(N_{\rm spec}\times N_{\rm pix}\) array of object spectra to fit. Can be a numpy.ma.MaskedArray.
obj_ferr (numpy.ndarray) – (Optional) \(N_{\rm spec}\times N_{\rm pix}\) array with the errors in the object spectra. Can be a numpy.ma.MaskedArray.
obj_sres (numpy.ndarray) – (Optional) 1D or 2D array with the spectral resolution (\(R = \lambda/\Delta\lambda\)) at each wavelength for (each of) the object spectra. Default is the resolution is not provided and assumed to be same as the template resolution.
- Returns:
Returns four arrays: the object wavelength vector, the object flux array, the object flux error and the object spectral resolution. The latter two objects can be None if the relevant object is None on input. All arrays have the same shape as the input
obj_flux
. Ifobj_sres
is input as a 1D vector, the spectral resolution is repeated for all input flux vectors. All input arrays are copies of the input.- Return type:
tuple
- static check_pixel_scale(tpl_wave, obj_wave, velscale_ratio=None, dvtol=1e-10)[source]
Confirm that the pixel scale of the template and object spectra are identical within a certain tolerance, accounting for an input pixel-scale ratio. Returns the velocity scale of the object spectra and the velocity scale ratio wrt the template spectra.
- static check_templates(tpl_wave, tpl_flux, tpl_sres=None, velscale_ratio=None)[source]
Check that the input template data is valid for use with pPXFFit.
- static compile_model_flux(obj_flux, ppxf_result, rescale=False)[source]
Return the model flux but pulling the models out of the ppxf results. The model can be rescaled to the data based on the fitted pixels if rescale is True.
The output array is masked in the spectral regions below and above the fitted wavelength range; any intervening pixels are not masked, even if they’re not included in the fit.
- static construct_models(tpl_wave, tpl_flux, obj_wave, obj_flux_shape, model_par, select=None, redshift_only=False, deredshift=False, corrected_dispersion=False, dvtol=1e-10)[source]
Construct models using the provided set of model parameters. This is a wrapper for
PPXFModel
.Only the shape of the object data is needed, not the data itself.
Allows for a replacement template library that must have the same shape as
tpl_flux
.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.
To convolve the model to the corrected dispersion, instead of the uncorrected dispersion, set corrected_dispersion=True. Correction always uses SIGMACORR_EMP data.
- static convert_velocity(v, verr)[source]
Convert kinematics from pPXF from pixel shifts to redshifts. pPXF determines the velocity offset by making the approximation that every pixel (logarithmically binned in wavelength) is a constant change in velocity. Over large velocity shifts, this approximation can become poor. An e-mail summary from Michele Cappellari:
The velocity scale per pixel is input as
\[\delta v = c \delta\ln\lambda = c (\ln\lambda_1 - \ln\lambda_0)\]The velocites output by pPXF are:
\[V = \delta v N_{\rm shift} = c \ln(\lambda_{N_{\rm shift}}/\lambda_0)\]which implies that the relation between PPXF output velocity and redshift is
\[1 + z = exp(V/c),\]which reduces z~vel/c in the low-redshift limit. This function converts the \(V\) values provided by pPXF to \(cz\) velocities.
Note
IMPORTANT: After this conversion, you must revert the velocities back to the “pixel-based velocities” (using
_revert_velocity()
) before using the velocities to reconstruct the pPXF fitted model.
- fit(tpl_wave, tpl_flux, obj_wave, obj_flux, obj_ferr, guess_redshift, guess_dispersion, iteration_mode='global_template', reject_boxcar=100, ensemble=True, velscale_ratio=None, mask=None, usetpl=None, matched_resolution=True, tpl_sres=None, obj_sres=None, waverange=None, bias=None, degree=4, mdegree=0, moments=2, loggers=None, quiet=False, max_velocity_range=400.0, alias_window=None, dvtol=1e-10, plot=False, plot_file_root=None)[source]
Wrapper for pPXF with some additional convenience functions. Limited implementation at the moment.
- Parameters:
tpl_wave (numpy.ndarray) – 1D vector of template wavelengths at rest in angstroms.
tpl_flux (numpy.ndarray) – N-templates x N-wavelengths array of template spectra to fit.
obj_wave (numpy.ndarray) – 1D vector of object wavelengths in angstroms. Does NOT need to be same as the template wavelengths.
obj_flux (numpy.ndarray) – N-spec x N-wavelengths array of object spectra to fit. Can be a numpy.ma.MaskedArray.
obj_ferr (numpy.ndarray) – N-spec x N-wavelengths array with the errors in the object spectra.
guess_redshift (
float
, numpy.ndarray) – Single or spectrum-specific redshift used to set the initial guess kinematics.guess_dispersion (
float
, numpy.ndarray) – Single or spectrum-specific velocity dispersion used to set the initial guess kinematics.iteration_mode (
str
, optional) – Iteration sequence to perform. Seeiteration_modes()
.reject_boxcar (
int
, optional) – Size of the boxcar to use during the rejection iteration. Default is 100. If None, rejection uses the entire residual spectrum.ensemble (
bool
, optional) – Treat the list of input spectra as an ensemble. Currently, this only affects how the spectra are masked. Default is to treat them as an ensemble. When not treated as an ensemble, each spectrum is masked individually according to the input wavelength range and velocity offsets. It does not make sense to set the iteration_mode to something that will include a fit to the global spectrum if you’re not treating the list of object spectra as an ensemble.velscale_ratio (
int
, optional) – Ratio of velocity scale per pixel in the object spectra to that for the template spectra. Default is that they should be identical.mask (numpy.ndarray,
mangadap.util.pixelmask.SpectralPixelMask
, optional) – A baseline pixel mask to use during the fitting. Other pixels may be masked via the convenience functions, but these pixels will always be masked.matched_resolution (
bool
, optional) – Flag that the object and template spectra have identical spectral resolution. Default is True.tpl_sres (numpy.ndarray, optional) – One-dimensional vector with the spectral resolution (\(R = \lambda/\Delta\lambda\)) at each wavelength of the template spectra. Default is the resolution is not provided and assumed to be same as the object resolution.
obj_sres (numpy.ndarray, optional) – One- or Two-dimensional array with the spectral resolution (\(R = \lambda/\Delta\lambda\)) at each wavelength for (each of) the object spectra. Default is the resolution is not provided and assumed to be same as the template resolution.
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) – From the pPXF documentation: This parameter biases the (h3, h4, …) measurements towards zero (Gaussian LOSVD) unless their inclusion significantly decreases the error in the fit. Set this to BIAS=0.0 not to bias the fit: the solution (including [V, sigma]) will be noisier in that case. The default BIAS should provide acceptable results in most cases, but it would be safe to test it with Monte Carlo simulations. This keyword precisely corresponds to the parameter lambda in the Cappellari & Emsellem (2004) paper. Note that the penalty depends on the relative change of the fit residuals, so it is insensitive to proper scaling of the NOISE vector. A nonzero BIAS can be safely used even without a reliable NOISE spectrum, or with equal weighting for all pixels.degree (
int
, optional) – From the pPXF documentation: degree of the additive Legendre polynomial used to correct the template continuum shape during the fit (default: 4). Set DEGREE = -1 not to include any additive polynomial.mdegree (
int
, optional) –From the pPXF documentation: degree of the multiplicative Legendre polynomial (with mean of 1) used to correct the continuum shape during the fit (default: 0). The zero degree multiplicative polynomial is always included in the fit as it corresponds to the weights assigned to the templates. Note that the computation time is longer with multiplicative polynomials than with the same number of additive polynomials.
Note
IMPORTANT: Multiplicative polynomials cannot be used when the REDDENING keyword is set.
moments (
int
, optional) –From the pPXF documentation: Order of the Gauss-Hermite moments to fit. Set this keyword to 4 to fit [h3, h4] and to 6 to fit [h3, h4, h5, h6]. Note that in all cases the G-H moments are fitted (non-linearly) together with [V, sigma].
If MOMENTS=2 or MOMENTS is not set then only [V, sigma] are fitted and the other parameters are returned as zero.
If MOMENTS is negative then the kinematics of the given COMPONENT are kept fixed to the input values.
EXAMPLE: We want to keep fixed component 0, which has an LOSVD described by [V, sigma, h3, h4] and is modelled with 100 spectral templates; At the same time we fit [V, sigma] for COMPONENT=1, which is described by 5 templates (this situation may arise when fitting stellar templates with pre-determined stellar kinematics, while fitting the gas emission). We should give in input to ppxf() the following parameters: component = [0]*100 + [1]*5 # –> [0, 0, …, 0, 1, 1, 1, 1, 1] moments = [-4, 2] start = [[V, sigma, h3, h4], [V, sigma]]
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.quiet (
bool
, optional) – Suppress all terminal and logging output. Default is False.max_velocity_range (
float
, optional) – Maximum range (+/-) expected for the fitted velocities in 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.plot (
bool
, optional) – Show the automatically generated pPXF fit plots.
- Returns:
Returns 4 numpy.ndarray objects:
The wavelengths of the best fitting model spectra. Nominally the same as the wavelengths of the input object spectra (obj_wave).
The fluxes of the best-fitting model spectra.
A mask for the best-fitting models spectra, following from the internal bitmask.
A record array with the fitted model parameters; see
spectralfitting.StellarKinematicsFit._per_stellar_kinematics_dtype
.
- Return type:
tuple
- Raises:
ValueError – Raised if the input arrays are not of the correct shape or if the pixel scale of the template and object spectra is greater than the specified tolerance.
- fit_SpatiallyBinnedSpectra(binned_spectra, gpm=None, par=None, loggers=None, quiet=False, debug=False)[source]
This is a basic interface that is geared for the DAP that interacts with the rest of the, more general, parts of the class.
This should not declare anything to self!
- static fitting_mask(tpl_wave, obj_wave, velscale, velscale_ratio=None, waverange=None, velocity_offset=None, max_velocity_range=400.0, alias_window=None, loggers=None, quiet=False)[source]
Return a list of pixels in the object spectrum to be fit using pPXF.
Be clear between velocity (ppxf) vs. redshift (cz) !
The limits applied to the fitted pixels are:
Apply the provided wavelength range limit (waverange).
pPXF will only allow a fit when the number of template pixels is the same as or exceeds the number of pixels in the object spectrum. The first step toward limiting template spectra that are too long is to truncate the blue and red edges that likely won’t be used given the provided velocity offsets (velocity_offset) and the expected velocity range (max_velocity_range).
Remove leading and trailing pixels that will cause alias problems during the convolution with the LOSVD (alias_window).
- Parameters:
obj_wave (array) – Wavelength vector of the object spectrum to be fit.
tpl_wave (array) – Wavelength vector of the template library to fit to the object spectrum.
velscale (float) – Velocity scale of the pixel.
velscale_ratio (int) – (Optional) Ratio of the object velscale to the template velscale. Default is 1 (i.e. the two have the same pixel scale).
waverange (array) – (Optional) Lower and upper wavelength limits to include in the analysis. The array can either define a single wavelength range – shape is (2,) – or a set of wavelength ranges – shape is (n,2). Default is to apply no wavelength range limitation.
velocity_offset (array) – (Optional) Vector with the velocity offset (expected or actual) between the template and the object spectrum in km/s. Used to estimate which wavelengths can be removed from the template. This can be a single offset or a set of offsets. If both waverange and velocity_offset are 2D arrays, the number of wavelength ranges and offsets must be the same. Default is that there is no velocity offset.
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.
loggers (list) – (Optional) List of logging.Logger objects to log progress; ignored if quiet=True. Logging is done using
mangadap.util.log.log_output()
. Default is no logging.quiet (bool) – (Optional) Suppress all terminal and logging output. Default is False.
- Returns:
Four boolean vectors are returned:
flags for pixels to include in the fit
flags for pixels that were excluded because they were outside the designated wavelength range
flags for pixels that were excluded to ensure the proper length of the template spectrum wrt the object spectrum.
flags for pixels that were truncated to avoid convolution aliasing
- Return type:
numpy.ndarray
- Raises:
ValueError – Raised if (1) no pixels are valid for the fit, (2) the template and object do not have overlapping spectral regions given the expected velocity offset between the two, or (3) the turncation to deal with aliasing removes all remaining pixels.
- static initialize_model_mask(obj_wave, obj_flux, mask=None, bitmask=None, velocity_offset=None)[source]
- static initialize_pixels_to_fit(tpl_wave, obj_wave, obj_flux, obj_ferr, velscale, velscale_ratio=None, waverange=None, mask=None, bitmask=None, velocity_offset=None, max_velocity_range=400.0, alias_window=None, ensemble=True, loggers=None, quiet=False)[source]
obj_flux and obj_ferr must be masked arrays.
Todo
This function can alter obj_flux and obj_ferr !! Return them instead?
- static iteration_modes()[source]
Possible iteration methods:
none
: Fit all bins with all templates with a single call to pPXF.no_global_wrej
: Do not fit the global spectrum first, but include a rejection iteration. All templates are fit in each step.global_template
: Fit the global spectrum with all templates and include a single rejection iteration. The pixel mask for this fit is the base mask for all fits to the individual bins. A single rejection iteration is done for each bin. Only the global template is used when fitting each bin.nonzero_templates
: Fit the global spectrum with all templates and include a single rejection iteration. The pixel mask for this fit is the base mask for all fits to the individual bins. A single rejection iteration is done for each bin. Only the templates with non-zero weights are used when fitting each bin.all_templates
: Fit the global spectrum with all templates and include a single rejection iteration. The pixel mask for this fit is the base mask for all fits to the individual bins. A single rejection iteration is done for each bin. All templates are used when fitting each bin.
- Returns:
List of allowed options.
- Return type:
list
- static losvd_limits(velscale)[source]
Return the limits on the LOSVD parameters used by pPXF.
Velocity limits are \(\pm 2000\) km/s
Velocity-disperison limits are from 1/10 pixels to 1000 km/s
Limits of the higher orders moments are from -0.3 to 0.3
- static obj_tpl_pixelmatch(velscale, tpl_wave, velscale_ratio=None, dvtol=1e-10)[source]
Confirm that the pixel scale of the template and object data are the same within some tolerance, accounting for and input ratio of the two.
- static ppxf_tpl_obj_voff(tpl_wave, obj_wave, velscale, velscale_ratio=None)[source]
Determine the pseudo offset in velocity between the template and object spectra, just due to the difference in the starting wavelengths.
This calculation is independent of the base of the logarithm used the sampling of the spectra.
Assumes wavelengths are logarithmically binned.
- Parameters:
tpl_wave (numpy.ndarray) – Wavelength vector for the template library to fit to the object spectrum.
obj_wave (numpy.ndarray) – Wavelength vector for the object spectrum to be fit.
velscale (float) – Velocity step per pixel in km/s for the object spectrum.
velscale_ratio (int) – (Optional) The integer ratio between the velocity scale of the pixel in the galaxy data to that of the template data. This is used only when constructing the template library. Default is None, which is the same as assuming that the velocity scales are identical.
- Returns:
Velocity offset in km/s between the initial wavelengths of the template and object spectra.
- Return type:
float
Todo
Implement a check that calculates the velocity ratio directly?
- rej_flag = 'PPXF_REJECT'
- static reject_model_outliers(obj_flux, ppxf_result, rescale=False, local_sigma=False, boxcar=None, nsigma=3.0, niter=9, loggers=None, quiet=False)[source]
- static revert_velocity(v, verr)[source]
Revert the velocity back to the “pixelized” velocity returned by pPXF.
Order matters here. The error computation is NOT a true error propagation; it’s just the inverse of the “convert” operation
- rng_flag = 'OUTSIDE_RANGE'
- snr_flag = 'LOW_SNR'
- tpl_flag = 'TPL_PIXELS'
- trunc_flag = 'TRUNCATED'
- class mangadap.proc.ppxffit.PPXFFitPar(guess_redshift=0.0, guess_dispersion=100.0, iteration_mode='nonzero_templates', template_library='MILESHC', matched_resolution=False, pixelmask=None, reject_boxcar=100, velscale_ratio=1, bias=None, degree=8, mdegree=0, moments=2)[source]
Bases:
KeywordParSet
Define a parameter set used by the pPXF fitting method.
The defined parameters are:
Key
Type
Options
Default
Description
guess_redshift
ndarray, list, int, float
0.0
Initial guess for the redshift (\(cz\)) of each binned spectrum.
guess_dispersion
ndarray, list, int, float
100.0
Initial guess for the velocity dispersion for each binned spectrum.
iteration_mode
str
none
,no_global_wrej
,global_template
,nonzero_templates
,all_templates
nonzero_templates
Iteration mode to use; see
PPXFFit.iteration_modes()
.template_library
str, TemplateLibraryDef
MILESHC
Keyword identifier or definition object used to build the spectral template library used during the fit.
matched_resolution
bool
False
Flag that the spectral resolution of the templates are (should be) matched to the galaxy data.
pixelmask
PixelMask
Pixel mask to include during the fitting; see
SpectralPixelMask()
.reject_boxcar
int
100
Number of pixels in the boxcar used to determine the local sigma for rejecting outliers.
velscale_ratio
int
1
The integer ratio between the velocity scale of the pixel in the galaxy data to that of the template data.
bias
int, float
ppxf
bias
parameter used to penalize low S/N spectra toward a Gaussian LOSVD.degree
int
8
ppxf
degree
parameter used to set the order of the additive polynomial to include in the fit.mdegree
int
0
ppxf
mdegree
parameter used to set the order of the multiplicative polynomial to include in the fit.moments
int
2, 4, 6
2
ppxf
moments
parameter used to set the number of moments of the LOSVD to fit. The DAP has not been well tested for fits that include any more than \(V\) and \(\sigma\).- fill(binned_spectra, guess_vel=None, guess_sig=None, **kwargs)[source]
Use the provided object to “fill” the parameters by constructing the template library spectra used during the fit and to set the initial guess kinematics for each spectrum.
kwargs are passed directly to the template library construction
- class mangadap.proc.ppxffit.PPXFFitResult(degree, mdegree, start, end, tpl_to_use, ppxf_fit, ntpl, weight_errors=False, component_fits=False)[source]
Bases:
object
A basic utility to save the critical parts of the pPXF model.
- class mangadap.proc.ppxffit.PPXFModel(templates, galaxy, velscale, velscale_ratio=None, vsyst=None, sigma_diff=0, moments=2, degree=4, mdegree=0, component=0, gas_component=None, lam=None, reddening=None, gas_reddening=None, reddening_func=None, sky=None, templates_rfft=None, trig=False, quiet=False)[source]
Bases:
object
Class that reconstructs a pPXF model given a set of templates and the model parameters.
This pulls functions from M. Cappellari’s ppxf class, version 6.7.0.
Input common to ppxf() input should have the same format as pPXF input.
The only reason the galaxy spectrum (or spectra) are provided is to set if the data are to fit as reflection-symmetric and how many spectral pixels there are.
- construct(kin, tplwgts, addpoly=None, multpoly=None, reddening=None, gas_reddening=None)[source]
Construct a pPXF model spectrum provided the input parameters. Mostly a copy of ppxf._linear_fit().
- Parameters:
kin (list, numpy.ndarray) – Must have the kinematics of each component. Must be a list or array of vectors with a shape (NCOMP,). The length of each vector in the list/array must be the same as MOMENTS.
tplwgts (numpy.ndarray) – Weights to apply to each template. Shape must be (NTEMP,).
addpoly (numpy.ndarray) – (Optional) Coefficients of the additive polynomials. Shape must be (NSPEC*(DEGREE+1),). Exception is raised if coefficients are expected (self.degree > -1), but no coefficiencts are provided.
multpoly (numpy.ndarray) – (Optional) Coefficients of the multiplicative polynomials. Shape must be (NSPEC*MDEGREE,). Exception is raised if coefficients are expected (self.mdegree > 0), but no coefficiencts are provided.
reddening (float) – (Optional) E(B-V) value for the continuum fit. Exception is raised if value is expected (self.reddening is not None), but none is provided.
gas_reddening (float) – (Optional) E(B-V) value for the gas components. Exception is raised if value is expected (self.reddening is not None), but none is provided.
- Return type:
numpy.ndarray