MultiBeam

class grizli.multifit.MultiBeam(beams, group_name=None, fcontam=0.0, psf=False, polyx=[0.3, 2.5], MW_EBV=0.0, min_mask=0.01, min_sens=0.08, sys_err=0.0, mask_resid=True, verbose=True, replace_direct=None, **kwargs)[source]

Bases: grizli.fitting.GroupFitter

Tools for dealing with multiple BeamCutout instances

Parameters:
beams : list

List of BeamCutout objects.

group_name : str, None

Rootname to use for saved products. If None, then default to ‘group’.

fcontam : float

Factor to use to downweight contaminated pixels. The pixel inverse variances are scaled by the following weight factor when evaluating chi-squared of a 2D fit,

weight = np.exp(-(fcontam*np.abs(contam)*np.sqrt(ivar)))

where contam is the contaminating flux and ivar is the initial pixel inverse variance.

psf : bool

Fit an ePSF model to the direct image to use as the morphological reference.

MW_EBV : float

Milky way foreground extinction.

min_mask : float

Minimum factor relative to the maximum pixel value of the flat f-lambda model where the 2D cutout data are considered good. Passed through to BeamCutout.

min_sens : float

See BeamCutout.

sys_err : float

Systematic error added in quadrature to the pixel variances:

var_total = var_initial + (beam.sci*sys_err)**2

Attributes:
TBD : type

Attributes Summary

N
Ngrism dictionary containing number of exposures by grism
PA Available PAs in each grism
grisms Available grisms
id

Methods Summary

apply_trace_shift([set_to_zero]) Set beam.yoffset back to zero
check_for_bad_PAs([poly_order, …])
compute_exptime() Compute number of exposures and total exposure time for each grism
compute_model([id, spectrum_1d, is_cgs, …]) Compute the dispersed 2D model for an assumed input spectrum
compute_model_psf([id, spectrum_1d, is_cgs]) Compute the dispersed 2D model for an assumed input spectrum and for ePSF morphologies
drizzle_fit_lines(fit, pline[, force_line, …]) TBD
drizzle_grisms_and_PAs([size, fcontam, …]) Make figure showing spectra at different orients/grisms
drizzle_segmentation([wcsobj, kernel, …]) Drizzle segmentation image from individual MultiBeam.beams.
drizzle_segmentation_id([id, wcsobj, …]) Drizzle segmentation image for a single ID
eval_poly_spec(coeffs_full) Evaluate polynomial spectrum
eval_trace_shift(shifts, self, indices, …) TBD
extend(new[, verbose]) Concatenate MultiBeam objects
fit_at_z([z, templates, fitter, …]) TBD
fit_redshift([prior, poly_order, fwhm, …]) TBD
fit_stars([poly_order, fitter, …]) TBD
fit_trace_shift([split_groups, max_shift, …]) TBD
flag_with_drizzled(hdul[, sigma, update, …]) Update MultiBeam masks based on the blotted drizzled combined image
init_poly_coeffs([flat, poly_order]) TBD
load_beam_fits(beam_list[, conf, verbose]) TBD
load_master_fits(beam_file[, verbose]) Load a “beams.fits” file.
make_simple_hdulist() Make a`~astropy.io.fits.HDUList` object with just a simple PrimaryHDU
oned_spectrum([tfit, get_contam, …]) Compute full 1D spectrum with optional best-fit model
oned_spectrum_to_hdu([sp, outputfile, units]) Generate 1D spectra fits HDUList
parse_fit_outputs(z, templates, coeffs_full, A) Parse output from fit_at_z.
redshift_fit_twod_figure(fit[, …]) Make figure of 2D spectrum
replace_direct_image_cutouts([ref_image, …]) Replace “REF” extensions in a beams.fits file
replace_segmentation_image_cutouts([ref_image]) Replace “REF” extensions in a beams.fits file
reshape_flat(flat_array) TBD
run_full_diagnostics([pzfit, pspec2, pline, …]) TBD
run_individual_fits([z, templates]) Run template fits on each exposure individually to evaluate variance in line and continuum fits.
show_redshift_fit(fit_data[, plot_flambda, …]) TBD
write_beam_fits([verbose]) TBD
write_master_fits([verbose, get_hdu, strip, …]) Store all beams in a single HDU TBD

Attributes Documentation

N
Ngrism

dictionary containing number of exposures by grism

PA

Available PAs in each grism

grisms

Available grisms

id

Methods Documentation

apply_trace_shift(set_to_zero=False)[source]

Set beam.yoffset back to zero

check_for_bad_PAs(poly_order=1, chi2_threshold=1.5, fit_background=True, reinit=True)[source]
compute_exptime()[source]

Compute number of exposures and total exposure time for each grism

compute_model(id=None, spectrum_1d=None, is_cgs=False, apply_sensitivity=True, scale=None, reset=True)[source]

Compute the dispersed 2D model for an assumed input spectrum

This is a wrapper around the grizli.model.GrismDisperser.compute_model method, where the parameters are described.

Nothing returned, but the model and modelf attributes are updated on the GrismDisperser subcomponents of the beams list.

compute_model_psf(id=None, spectrum_1d=None, is_cgs=False)[source]

Compute the dispersed 2D model for an assumed input spectrum and for ePSF morphologies

This is a wrapper around the grizli.model.GrismDisperser.compute_model_psf method, where the parameters are described.

Nothing returned, but the model and modelf attributes are updated on the GrismDisperser subcomponents of the beams list.

drizzle_fit_lines(fit, pline, force_line=['Ha+NII', 'Ha', 'OIII', 'Hb', 'OII'], save_fits=True, mask_lines=True, mask_sn_limit=3, mask_4959=True, verbose=True, include_segmentation=True, get_ir_psfs=True, min_line_sn=4)[source]

TBD

drizzle_grisms_and_PAs(size=10, fcontam=0, flambda=False, scale=1, pixfrac=0.5, kernel='square', usewcs=False, tfit=None, diff=True, grism_list=['G800L', 'G102', 'G141', 'F090W', 'F115W', 'F150W', 'F200W', 'F356W', 'F410M', 'F444W'], mask_segmentation=True, reset_model=True, make_figure=True, fig_args={'average_only': False, 'cmap': 'viridis_r', 'mask_segmentation': True, 'scale_size': 1}, **kwargs)[source]

Make figure showing spectra at different orients/grisms

TBD

drizzle_segmentation(wcsobj=None, kernel='square', pixfrac=1, verbose=False)[source]

Drizzle segmentation image from individual MultiBeam.beams.

Parameters:
wcsobj: `~astropy.wcs.WCS` or `~astropy.io.fits.Header`

Output WCS.

kernel: e.g., ‘square’, ‘point’, ‘gaussian’

Drizzle kernel, see drizzle.

pixfrac: float

Drizzle ‘pixfrac’, see drizzle.

verbose: bool

Print status messages.

Returns:
drizzled_segm: ndarray, type int64.

Drizzled segmentation image, with image dimensions and WCS defined in wcsobj.

drizzle_segmentation_id(id=None, wcsobj=None, kernel='square', pixfrac=1, verbose=True)[source]

Drizzle segmentation image for a single ID

eval_poly_spec(coeffs_full)[source]

Evaluate polynomial spectrum

static eval_trace_shift(shifts, self, indices, poly_order, lm, verbose, fit_with_psf)[source]

TBD

extend(new, verbose=True)[source]

Concatenate MultiBeam objects

Parameters:
new : MultiBeam

Beam object containing new beams to add.

verbose : bool

Print summary of the change.

fit_at_z(z=0.0, templates={}, fitter='nnls', fit_background=True, poly_order=0)[source]

TBD

fit_redshift(prior=None, poly_order=1, fwhm=1200, make_figure=True, zr=None, dz=None, verbose=True, fit_background=True, fitter='nnls', delta_chi2_threshold=0.004, zoom=True, line_complexes=True, templates={}, figsize=[8, 5], fsps_templates=False)[source]

TBD

fit_stars(poly_order=1, fitter='nnls', fit_background=True, verbose=True, make_figure=True, zoom=None, delta_chi2_threshold=0.004, zr=0, dz=0, fwhm=0, prior=None, templates={}, figsize=[8, 5], fsps_templates=False)[source]

TBD

fit_trace_shift(split_groups=True, max_shift=5, tol=0.01, verbose=True, lm=False, fit_with_psf=False, reset=False)[source]

TBD

flag_with_drizzled(hdul, sigma=4, update=True, interp='nearest', verbose=True)[source]

Update MultiBeam masks based on the blotted drizzled combined image

[in progress … xxx]

Parameters:
hdul : HDUList

FITS HDU list output from drizzle_grisms_and_PAs or read from a stack.fits file.

sigma : float

Residual threshold to flag.

update : bool

Update the mask.

interp : str

Interpolation method for ablot.

Returns:
Updates the individual fit_mask attributes of the individual beams
if update==True.
init_poly_coeffs(flat=None, poly_order=1)[source]

TBD

load_beam_fits(beam_list, conf=None, verbose=True)[source]

TBD

load_master_fits(beam_file, verbose=True)[source]

Load a “beams.fits” file.

make_simple_hdulist()[source]

Make a`~astropy.io.fits.HDUList` object with just a simple PrimaryHDU

oned_spectrum(tfit=None, get_contam=True, get_background=False, masked_model=None, **kwargs)[source]

Compute full 1D spectrum with optional best-fit model

Parameters:
bin : float / int

Bin factor relative to the size of the native spectral bins of a given grism.

tfit : dict

Output of template_at_z.

Returns:
sp : dict

Dictionary of the extracted 1D spectra. Keys are the grism names and the values are Table objects.

oned_spectrum_to_hdu(sp=None, outputfile=None, units=None, **kwargs)[source]

Generate 1D spectra fits HDUList

Parameters:
sp : optional, dict

Output of oned_spectrum. If None, then run that function with **kwargs.

outputfile : None, str

If a string supplied, then write the HDUList to a file.

Returns:
hdul : HDUList

FITS version of the 1D spectrum tables.

parse_fit_outputs(z, templates, coeffs_full, A)[source]

Parse output from fit_at_z.

Parameters:
z : float

Redshift at which to evaluate the fits.

templates : list of SpectrumTemplate objects

Generated with, e.g., load_templates.

coeffs_full : ndarray

Template fit coefficients

A : ndarray

Matrix generated for fits and used for computing model 2D spectra:

>>> model_flat = np.dot(coeffs_full, A)
>>> # mb = MultiBeam(...)
>>> all_models = mb.reshape_flat(model_flat)
>>> m0 = all_models[0] # model for mb.beams[0]
Returns:
line_flux : dict

Line fluxes and uncertainties, in cgs units (erg/s/cm2)

covar : ndarray

Covariance matrix for the fit coefficients

cont1d, line1d, model1d : SpectrumTemplate

Best-fit continuum, line, and full (continuum + line) templates

model_continuum : ndarray

Flat array of the best fit 2D continuum

redshift_fit_twod_figure(fit, spatial_scale=1, dlam=46.0, NY=10, figsize=[8, 3.5], **kwargs)[source]

Make figure of 2D spectrum

TBD

replace_direct_image_cutouts(ref_image='gdn-100mas-f160w_drz_sci.fits', ext=0, interp='poly5', cutout=200, background_func=<function median>, thumb_labels=None)[source]

Replace “REF” extensions in a beams.fits file

Parameters:
ref_image : str or HDUList

Filename or preloaded FITS file.

interp : str

Interpolation function to use for do_blot.

cutout : int

Make a slice of the ref_image with size [-cutout,+cutout] around the center position of the desired object before passing to blot.

background_func : function, None

If not None, compute local background with value from background_func(ref_image[cutout]).

Returns:
beams_image : HDUList

Image object with the “REF” extensions filled with the new blotted image cutouts.

replace_segmentation_image_cutouts(ref_image='gdn-100mas-f160w_seg.fits')[source]

Replace “REF” extensions in a beams.fits file

Parameters:
ref_image : str, HDUList, ImageHDU

Filename or preloaded FITS file.

Returns:
beams_image : HDUList

Image object with the “REF” extensions filled with the new blotted image cutouts.

reshape_flat(flat_array)[source]

TBD

run_full_diagnostics(pzfit={}, pspec2={}, pline={}, force_line=['Ha+NII', 'Ha', 'OIII', 'Hb', 'OII'], GroupFLT=None, prior=None, zoom=True, verbose=True)[source]

TBD

size=20, pixscale=0.1, pixfrac=0.2, kernel=’square’

run_individual_fits(z=0, templates={})[source]

Run template fits on each exposure individually to evaluate variance in line and continuum fits.

Parameters:
z : float

Redshift at which to evaluate the fit

templates : list of SpectrumTemplate objects

Generated with, e.g., load_templates.

Returns:
line_flux, line_err : dict

Dictionaries with the measured line fluxes and uncertainties for each exposure fit.

coeffs_list : ndarray [Nbeam x Ntemplate]

Raw fit coefficients

chi2_list, DoF_list : ndarray [Nbeam]

Chi-squared and effective degrees of freedom for each separate fit

show_redshift_fit(fit_data, plot_flambda=True, figsize=[8, 5])[source]

TBD

write_beam_fits(verbose=True)[source]

TBD

write_master_fits(verbose=True, get_hdu=False, strip=True, include_model=False, get_trace_table=True)[source]

Store all beams in a single HDU TBD