GroupFitter

class grizli.fitting.GroupFitter[source]

Bases: object

Base class for StackFitter and MultiBeam spectrum fitting objects

Methods Summary

check_tfit_coeffs(tfit, templates[, …]) Compare emission line fluxes fit at each grism/PA to the combined value.
compute_D4000(z[, fit_background, fit_type, …]) Compute D4000 with step templates
compute_scale_array(pscale, wave) Return the scale array given the coefficients
fit_line_width([bl, nl, z0, max_nfev, tol, …]) Fit for emisson line width
get_SDSS_photometry([bands, templ, radius, …]) Try too get SDSS photometry from astroquery
get_flat_background(bg_params[, apply_mask]) Generate background array the same size as the flattened total science array.
get_flat_model(spectrum_1d[, id, …]) Generate model array based on the model 1D spectrum in spectrum_1d
initialize_masked_arrays([seg_ids]) Initialize flat masked arrays for faster likelihood calculation
objfun_scale(pscale, AxT, data, self, retval) Objective function for fitting for a scale term between photometry and spectra
oned_figure([bin, wave, show_beams, minor, …]) Make a figure showing the 1D spectra
optimal_extract([data, bin, wave, ivar, …]) Binned optimal extractions by grism with algorithm from Horne 1984
scale_to_photometry([tfit, tol, order, …]) Compute scale factor between spectra and photometry
set_photometry([flam, eflam, filters, …]) Set photometry attributes
template_at_z([z, templates, fwhm, …]) Get the best-fit template at a specified redshift
unset_photometry() Unset photometry-related attributes
xfit_at_z([z, templates, fitter, …]) Fit the 2D spectra with a set of templates at a specified redshift.
xfit_redshift([prior, templates, fwhm, …]) Two-step procedure for fitting redshifts
xfit_star([tstar, spline_correction, …]) Fit stellar templates
xmake_fit_plot(fit, tfit[, show_beams, bin, …]) Make a diagnostic plot of the redshift fit

Methods Documentation

check_tfit_coeffs(tfit, templates, refit_others=True, fit_background=True, fitter='nnls', bounded_kwargs={'method': 'bvls', 'tol': 1e-08, 'verbose': 0})[source]

Compare emission line fluxes fit at each grism/PA to the combined value. If refit_others=True, then compare the line fluxes to a fit from a new object generated excluding that grism/PA.

Returns:
max_line : str

Line species with the maximum deviation

max_line_diff : float

The maximum deviation for max_line (sigmas).

compare : dict

The full comparison dictionary

compute_D4000(z, fit_background=True, fit_type='D4000', fitter='lstsq')[source]

Compute D4000 with step templates

Parameters:
z : float

Redshift where to evaluate D4000

fit_background : bool

Include background in step template fit

fit_type : ‘D4000’, ‘Dn4000’
Definition to use:

D4000 = f_nu(3750-3950) / f_nu(4050-4250) Dn4000 = f_nu(3850-3950) / f_nu(4000-4100)

fitter : str

Least-squares method passed to template_at_z.

Returns:
w_d4000, t_d4000 : ndarray, dict

Step wavelengths and template dictionary

tfit : dict

Fit dictionary returned by template_at_z.

d4000, d4000_sigma : float

D4000 estimate and uncertainty from simple error propagation and step template fit covariance.

static compute_scale_array(pscale, wave)[source]

Return the scale array given the coefficients

Parameters:
pscale : array-like

Coefficients of the linear model normalized by factors of 10 per order, i.e, pscale = [10] is a constant unit scaling. Note that parameter order is reverse that expected by numpy.polyval.

wave : array-like

Wavelength grid in Angstroms. Scaling is normalized to (wave - 1e4)/1000.

Returns:
wscale : array-like

Scale factor

>>> pscale = [10]
>>> N = len(pscale)
>>> rescale = 10**(np.arange(N)+1)
>>> wscale = np.polyval((pscale/rescale)[::-1], (wave-1.e4)/1000.)
fit_line_width(bl=2.5, nl=1.1, z0=1.9367, max_nfev=100, tol=0.001, verbose=False)[source]

Fit for emisson line width

Returns:

width/(1000 km/s), z, nfev, (nfev==max_nfev)
get_SDSS_photometry(bands='ugriz', templ=None, radius=2, SDSS_CATALOG='V/147/sdss12', get_panstarrs=False)[source]

Try too get SDSS photometry from astroquery

(developmental)

get_flat_background(bg_params, apply_mask=True)[source]

Generate background array the same size as the flattened total science array.

Parameters:
bg_params : array with shape (self.N) or (self.N, M)

Background parameters for each beam, where the M axis is polynomial cofficients in the order expected by Polynomial2D. If the array is 1D, then provide a simple pedestal background.

Returns:
bg_model : Array with dimensions (self.fit_mask.sum(),)

Flattened, masked background array.

get_flat_model(spectrum_1d, id=None, apply_mask=True, is_cgs=True)[source]

Generate model array based on the model 1D spectrum in spectrum_1d

Parameters:
spectrum_1d : tuple, -1

Tuple of 1D arrays (wavelength, flux). If -1, then use the in_place model attributes of each beam.

id : int

Value that identifies pixels in the segmentation thumbnail with the desired object to model

apply_mask : bool

Return the model pixels applying the BeamCutout fit_mask attribute

is_cgs : bool

spectrum_1d flux array has CGS f-lambda flux density units.

Returns:
model : Array with dimensions (self.fit_mask.sum(),)

Flattened, masked model array.

initialize_masked_arrays(seg_ids=None)[source]

Initialize flat masked arrays for faster likelihood calculation

static objfun_scale(pscale, AxT, data, self, retval)[source]

Objective function for fitting for a scale term between photometry and spectra

oned_figure(bin=1, wave=None, show_beams=True, minor=0.1, tfit=None, show_rest=False, axc=None, figsize=[6, 4], fill=False, units='flam', min_sens_show=0.1, ylim_percentile=2, scale_on_stacked=False, show_individual_templates=False, apply_beam_mask=True, loglam_1d=True, trace_limits=None, show_contam=False, beam_models=None)[source]

Make a figure showing the 1D spectra

Parameters:
bin : float

Binning factor relative to nominal resolution (per pix) of each grism

wave : None, array

Fixed wavelength array for the sampled spectra

show_beams : bool

Show all individual beams

minor : float

Minor axis tick interval (microns)

tfit : dict

Fit information from template_at_z. If provided, then will include the best-fit models in the figure

show_rest : bool

Show rest-frame wavelengths

acx : AxesSubplot

If provided, then draw into existing axis without making a new figure

figsize : (float, float)

Figure size (inches)

fill : bool

plot filled spectra

show_individual_templates : bool

Show each individual template with its scaling along with the best-fit combination

units : str

Y-axis units

  • ‘flam’ = Scaled f-lambda cgs
  • ‘nJy’ = nanoJansky
  • ‘mJy’ = milliJansky
  • ‘eps’ = native detector units of electrons per second
  • ‘meps’ = “milli”-electrons per second
  • ‘spline[N]’ = Divide out a spline continuum
  • ‘resid’ = Residuals w.r.t. model in tfit
loglam_1d : bool

Plot as log wavelength

trace_limits : (float, float)

If provided, extract spectra relative to the (tilted) spectral trace

show_contam : bool

Include curves for contamination model

min_sens_show : float
ylim_percentile : float
Returns:
fig : Figure

Figure object

optimal_extract(data=None, bin=1, wave=None, ivar=None, trace_limits=None, loglam=True, **kwargs)[source]

Binned optimal extractions by grism with algorithm from Horne 1984

The spatial profile for each beam is the 2D model spectrum generated using its attached direct image thumbnail. The Horne (1984) algorithm is essentially a least-squares fit of the spatial model to the observed 2D spectrum, weighted by the uncertainties.

Along with the optimal extraction, this method also implements an option to extract an effective “aperture” within a specified region above and below the spectral trace.

While the traces may not be directly aligned with the x axis of the 2D spectra, both the optimal and trace extractions extract along y pixels at a given x.

Parameters:
data : ndarray, None

Data array with same dimensions as self.scif_mask (flattened & masked) 2D spectra of all beams. If None, then use self.scif_mask.

bin : bool

Binning factor relative to the grism-dependent resolution values, specified in GRISM_LIMITS.

wave : ndarray, None

Wavelength bin edges. If None, then compute from parameters in GRISM_LIMITS.

ivar : ndarray, None

Inverse variance array with same dimensions as self.scif_mask (flattened & masked) 2D spectra of all beams. If None, then use self.weighted_sigma2_mask.

trace_limits : [float, float] or None

If specified, perform a simple sum in cross-dispersion axis between trace_limits relative to the central pixel of the trace rather than the optimally-weighted sum. Similarly, the output variances are the sum of the input variances in the trace interval.

Note that the trace interval is evaluated with < >, as opposed to <= >=, as the center of the trace is a float rather than an integer pixel index.

loglam : bool

If True and wave not specified (see above), then output wavelength grid is log-spaced.

Returns:
tab : dict

Dictionary of Table spectra for each available grism.

scale_to_photometry(tfit=None, tol=0.0001, order=0, init=None, fit_background=True, Rspline=50, use_fit=True, **kwargs)[source]

Compute scale factor between spectra and photometry

Parameters:
tfit : dict

Template fit info at a specific redshift from template_at_z. If not specified, then makes and fits spline templates

tol : float

Fit tolerance passed to the minimizer

order : int

Order of the polynomial scaling to fit

init : None

Initial parameters

fit_background : bool

Include additive background

Rspline : float

Spectral resolution R of spline templates

use_spline : bool

Use spline templates

Returns:
res : object

Result from scipy.optimize.least_squares. The coefficients of the linear scaling are in res.x.

set_photometry(flam=[], eflam=[], filters=[], ext_corr=1, lc=None, force=False, tempfilt=None, min_err=0.02, TEF=None, pz=None, source='unknown', **kwargs)[source]

Set photometry attributes

Parameters:
flam, eflam : array-like

Flux densities and uncertainties in f-lambda cgs units

filters : list

List of FilterDefinition objects

ext_corr : float or array-like

MW extinction correction

lc : array-like

Precomputed filter central wavelengths. Will automatically be computed from filters if not specified

force : bool

Don’t try to set if already specified (Nphot > 0)

tempfilt : eazy.photoz.TemplateGrid

Precomputed grid of templates integrated through the filters bandpasses

min_err : float

minimum or systematic error to add in quadrature to eflam

TEF : eazy.templates.TemplateError

Template error function

pz : None, (array, array)

Precomputed (z, pz) pdf from, e.g., eazy

source : str

String to indicate the provenance of the photometry

Returns:
photom_flam : array_like

Flux densities from flam

photom_eflam : array-like

Uncertainties including min_err

photom_filters : list

filters

Nphot : int

Number of photometry bandpasses

The returned parameters above are not returned but are rather
attributes that are set. This function also updates the
sivarf, weightf, fit_mask attributes to include the
spectra + photometry
template_at_z(z=0, templates=None, fwhm=1400, get_uncertainties=2, draws=0, **kwargs)[source]

Get the best-fit template at a specified redshift

Parameters:
z : float

Redshift

templates : dict

Dictionary of SpectrumTemplate objects

fwhm : float

FWHM of line templates if templates generated in-place

get_uncertainties : int

Get coefficient uncertainties from covariance matrix

draws : int

Number of random draws from covariance matrix

kwargs : dict

Any additional keywords are passed to xfit_at_z

Returns:
tfit : dict

Dictionary of fit results, used in various other places like oned_figure, etc.

Keyword Description
cfit Dict of template normalizations and uncertainties
cont1d SpectrumTemplate of best-fit continuum
line1d SpectrumTemplate of best-fit continuum + emission line
coeffs Array of fit coefficients
chi2 (float) chi-squared of the fit
z (float) The input redshift
templates Copy of the input templates dictionary
line1d_err If draws > 0, this will be template draws with the same dimension as line1d
unset_photometry()[source]

Unset photometry-related attributes

xfit_at_z(z=0, templates=[], fitter='nnls', fit_background=True, get_uncertainties=False, get_design_matrix=False, pscale=None, COEFF_SCALE=1e-19, get_components=False, huber_delta=4, get_residuals=False, include_photometry=True, use_cached_templates=False, bounded_kwargs={'method': 'bvls', 'tol': 1e-08, 'verbose': 0}, apply_sensitivity=True)[source]

Fit the 2D spectra with a set of templates at a specified redshift.

Parameters:
z : float

Redshift.

templates : list

List of templates to fit.

fitter : str

Minimization algorithm to compute template coefficients.

Available options are:

For the last option, the line flux limits are set by the limits in the global grizli.fitting.LINE_BOUNDS list and bounded_kwargs are passed to lsq_linear.

fit_background : bool

Fit additive pedestal background offset.

get_uncertainties : bool, int

Compute coefficient uncertainties from the covariance matrix. If specified as an int > 1, then the covariance matrix is computed only for templates with non-zero coefficients

get_design_matrix : bool

Return design matrix and data, rather than nominal outputs.

huber_delta : float

Use the Huber loss function (scipy.special.huber) rather than direct chi-squared. If huber_delta < 0, then fall back to chi-squared.

Returns:
chi2 : float

Chi-squared of the fit

coeffs, coeffs_err : numpy.ndarray

Template coefficients and uncertainties.

covariance : numpy.ndarray

Full covariance

xfit_redshift(prior=None, templates={}, fwhm=1200, line_complexes=True, fsps_templates=False, zr=[0.65, 1.6], dz=[0.005, 0.0004], zoom=True, verbose=True, fit_background=True, fitter='nnls', bounded_kwargs={'method': 'bvls', 'tol': 1e-08, 'verbose': 0}, delta_chi2_threshold=0.004, poly_order=3, make_figure=True, figsize=[8, 5], use_cached_templates=True, get_uncertainties=True, Rspline=30, huber_delta=4, get_student_logpdf=False)[source]

Two-step procedure for fitting redshifts

  1. polynomial, spline template fits
  2. redshift grids
Parameters:
prior : None, (array, array)

Redshift prior (z, pz). Will be interpolated to the redshift fit grid

templates : dict

Dictionary the SpectrumTemplate objects to use for the fits

fwhm, line_complexes, fsps_templates : float, bool, bool

Parameters passed to load_templates if templates is empty.

make_figure, fig_size : bool, (float, float)

Make the diagnostic figure with dimensions fig_size

zr : (float, float)

Redshift limits of the logarithmic (1+z) redshift grid

dz : (float, float)

Step size of the grid. The second value will be used to “zoom in” on the peaks found in the coarse grid step from the first value.

zoom : bool

Do the second pass with the dz[1] step size

verbose : bool

Some verbosity control

fit_background : bool

Include contribution of additive background

fitter, bounded_kwargs : str, dict

Least-squares optimization method. See xfit_at_z

delta_chi2_threshold : float

Not used

poly_order : int

Order of polynomials for the “uninformative” polynomial fit. The parameters of the polynomial and full template fits are computed to evaluate the extent to which the galaxy / stellar templates improve the fit

Rspline : float

Spectral resolution, R, of spline templates for another “uninformative” fit.

use_cached_templates : bool

Try to used cached versions of dispersed template models for templates that don’t depend on redshift (polynomials, splines)

get_uncertainties : bool

Get template fit coefficient uncertainties from the fit covariance matrix

huber_delta : float

Parameter for Huber loss function (see xfit_at_z)

get_student_logpdf : bool

Get logpdf for likelihood assuming Student-t distribution rather than standard normal assumption

Returns:
fit : Table

Table with fit information on the redshift grid and metadata on some fit characteristics.

Table metadata
Meta Description
N Number of spectrum extensions / beams
polyord Order of the polynomial fit
chi2poly \(\chi^2\) of the polynomial fit
chi2spl \(\chi^2\) of the spline fit
Rspline Spectral resolution of the spline templates
kspl Effective number of parameters of spline fit
huberdel huber_delta
splf[i] Flux of spline fit at fixed wavelengths
sple[i] Unc. of spline fit at fixed wavelengths
NTEMP Number of templates
DoF Degrees of freedom of the fit (total number of unmasked pixels in all 2D beams)
ktempl N parameters of the template fit
chimin Minimum \(\chi^2\) of the template fit
chimax Maximum \(\chi^2\) of the template fit
fitter Least squares method
as_epsf Fit was done as EffectivePSF
bic_poly Bayesian Information Criterion (BIC) of the polynomial fit. BIC = log(DoF)*k + min(chi2) + C
bic_spl BIC of the spline fit
bic_temp BIC of the template (redshift) fit
st_df Student-t df of spline fit
st_loc Student-t loc of spline fit
st_scl Student-t scale of spline fit
Z02 Integrated cdf(<z) = 0.025
Z16 Integrated cdf(<z) = 0.16
Z50 Integrated cdf(<z) = 0.50
Z84 Integrated cdf(<z) = 0.84
Z97 Integrated cdf(<z) = 0.975
ZWIDTH1 Z84 - Z16
ZWIDTH2 Z97 - Z02
z_map Redshift at Max(PDF)
zrmin Start of the redshift grid zr
zrmax End of the redshift grid zr
z_risk Redshift at minimum risk
min_risk Minimum risk
gam_loss gamma parameter of risk
Column Description
zgrid Redshift grid NZ
chi2 \(\chi^2(z)\) NZ
student_logpdf log PDF of student-t likelihood NZ
coeffs Template coefficients (NZ,NTEMP)
covar Template covariance (NZ,NTEMP,NTEMP)
pdf Full likelihood including optional prior
risk “Risk” parameter from Tanaka et al. (arXiv/1704.05988)
xfit_star(tstar=None, spline_correction=True, fitter='nnls', fit_background=True, spline_args={'Rspline': 5}, oned_args={})[source]

Fit stellar templates

Parameters:
tstar : dict

Dictionary of stellar SpectrumTemplate objects

spline_correction : bool

Include spline scaling correction for template mismatch

fitter : str

Least-squares method passed to template_at_z.

fit_background : bool

Fit for additive background component

spline_args : dict

Parameters passed to split_spline_template for generating the spline correction arrays

oned_args : dict

Keywords passed to oned_figure

Returns:
fig : Figure

Figure object

line : str

Line of text describing the best fit

tfit : dict

Fit information from template_at_z

xmake_fit_plot(fit, tfit, show_beams=True, bin=1, minor=0.1, scale_on_stacked_1d=True, loglam_1d=True, zspec=None)[source]

Make a diagnostic plot of the redshift fit

Parameters:
fit : Table

Redshift fit results from xfit_redshift

tfit : dict

Template fit at best redshift from template_at_z

show_beams : bool

Show 1D spectra of all individual “beams”

bin : float

Binning factor relative to nominal wavelength resolution (1 pix) of each grism

minor : float

Minor axis ticks, microns

scale_on_stacked_1d : bool

Set y limits based on stacked spectrum

loglam_1d : bool

Show log wavelengths

zspec : float, None

Spectroscopic redshift that will be indicated on the figure

Returns:
fig : Figure

Figure object