mangadap.proc.emissionlinemodel module

A class hierarchy that fits the emission lines.


License

Copyright © 2019, SDSS-IV/MaNGA Pipeline Group


class mangadap.proc.emissionlinemodel.EmissionLineModel(method, binned_spectra, stellar_continuum=None, emission_line_moments=None, redshift=None, dispersion=None, minimum_error=None, method_list=None, artifact_path=None, ism_line_path=None, emission_line_db_path=None, output_path=None, output_file=None, hardcopy=True, overwrite=None, checksum=False, loggers=None, quiet=False)[source]

Bases: object

Class that holds the emission-line model fits.

Parameters:
  • method (EmissionLineModelDef) – Object that defines the main parameters used for the emission-line model fitting.

  • loggers (list, optional) – List of logging.Logger objects to log progress; ignored if quiet=True. Logging is done using log_output().

  • quiet (bool, optional) – Suppress all terminal and logging output.

_add_method_header(hdr, model_binid=None)[source]

Add fitting method information to the header.

_assign_image_arrays()[source]

Set image_arrays, which contains the list of extensions in hdu that are on-sky image data.

_assign_input_kinematics(emission_line_moments, redshift, dispersion, default_dispersion=100.0)[source]

Set the initial redshift and velocity dispersion for each spectrum.

In terms of precedence, directly provided redshifts and dispersions override those in the EmissionLineMoments model, which overrides those in the StellarContinuumModel object. The default_dispersion does not take precedence over any provided disperison.

If emission_line_moments, redshift, and stellar_continuum are all None, the redshift is set to 0.0. If dispersion is also None, the dispersion is set to default_dispersion.

To get the stellar kinematics, the function calls mangadap.proc.stellarcontinuummodel.StellarContinuumModel.matched_kinematics(). In this function, the provided redshift and dispersion must be a single value or None; therefore, the means of any vectors provided as redshift or disperison are passsed to this function instead of the full vector.

Parameters:
_assign_spectral_arrays()[source]
_construct_2d_hdu(good_snr, model_flux, model_base, model_mask, model_eml_par, model_fit_par=None, model_binid=None)[source]

Construct hdu that is held in memory for manipulation of the object. See construct_3d_hdu() if you want to convert the object into a DRP-like datacube.

_finalize_cube_mask(mask)[source]

Finalize the mask after the 2D mask has been reconstructed into a 3D cube.

This mostly handles the masks for regions outside the datacube field of view.

Note that the input mask is both edited in-place and returned.

Todo

  • This needs to be abstracted for non-DRP datacubes.

  • Describe MAPMASK usage

Parameters:

mask (numpy.ndarray) – 3D array with the current bitmask data.

Returns:

Edited bitmask data.

Return type:

numpy.ndarray

_get_line_fit_metrics(model_flux, model_base, model_mask, model_eml_par, model_binid=None, window=15, debug=False)[source]

Compute fit-quality metrics near each fitted line.

Primarily a wrapper of EmissionLineFit.line_metrics().

Parameters:
  • model_flux (numpy.ndarray) – The best-fitting emission-line model spectra.

  • model_base (numpy.ndarray) – The “baseline” of the emission-line model defined as the difference between the best-fitting stellar-continuum model and the best-fitting continuum from the emission-line modeling procedure.

  • model_mask (numpy.ndarray) – A boolean or integer (bitmask) array with the mask for the model data.

  • model_eml_par (EmissionLineFitDataTable) – A record array with the best-fitting results and parameters for each emission line. The format is expected to be given by EmissionLineFit._per_emission_line_dtype(). This object is edited in place and returned.

  • model_binid (numpy.ndarray, optional) – A 2D array with the bin ID assigned to each spaxel with a fitted emission-line model. The number of non-negative bin IDs must match the number of fitted models. If None, this is generated from the BINID column in the model_eml_par object.

  • window (int, optional) – The window size in pixels around each line to use for calculation of the figures of merit.

Returns:

Return the input model_eml_par after filling the LINE_PIXC, AMP, ANR, LINE_NSTAT, LINE_CHI2, LINE_RMS, and LINE_FRMS columns.

Return type:

numpy.recarray

_get_missing_models(unique_bins=None, debug=False)[source]
_initialize_primary_header(hdr=None)[source]

Construct the primary header for the reference file.

Parameters:

hdr (astropy.io.fits.Header, optional) – Input base header for added keywords. If None, uses the datacube header from binned_spectra (if there is one) and then cleans the header using mangadap.util.fitsutil.DAPFitsUtil.clean_dap_primary_header().

Returns:

Initialized header object.

Return type:

astropy.io.fits.Header

channel_names()[source]
construct_3d_hdu()[source]

Reformat the model spectra into a cube matching the shape of the DRP fits file.

construct_continuum_models(replacement_templates=None, redshift_only=False, deredshift=False, dispersion_corrections=None)[source]
copy_to_array(ext='FLUX', waverange=None, include_missing=False)[source]

Wrapper for mangadap.util.fitsutil.DAPFitsUtil.copy_to_array() specific for EmissionLineModel.

Return a copy of the selected data array. The array size is always \(N_{\rm models} \times N_{\rm wavelength}\); i.e., the data is always flattened to two dimensions and the unique spectra are selected.

Parameters:
  • ext (str, optional) – Name of the extension from which to draw the data. Must be allowed for the current mode; see data_arrays. Default is 'FLUX'.

  • waverange (array-like, optional) – Two-element array with the first and last wavelength to include in the computation. Default is to use the full wavelength range.

  • include_missing (bool, optional) – Create an array with a size that accommodates the missing models.

Returns:

A 2D array with a copy of the data from the selected extension.

Return type:

numpy.ndarray

copy_to_masked_array(ext='FLUX', flag=None, waverange=None, include_missing=False)[source]

Wrapper for mangadap.util.fitsutil.DAPFitsUtil.copy_to_masked_array() specific for EmissionLineModel.

Return a copy of the selected data array as a masked array. This is functionally identical to copy_to_array(), except the output format is a numpy.ma.MaskedArray. The pixels that are considered to be masked can be specified using flag.

Parameters:
  • ext (str, optional) – Name of the extension from which to draw the data. Must be allowed for the current mode; see data_arrays. Default is ‘FLUX’.

  • flag (str, list, optional) – (List of) Flag names that are considered when deciding if a pixel should be masked. The names must be a valid bit name as defined by bitmask. If not provided, ANY non-zero mask bit is omitted.

  • waverange (array-like, optional) – Two-element array with the first and last wavelength to include in the computation. Default is to use the full wavelength range.

  • include_missing (bool, optional) – Create an array with a size that accommodates the missing models.

Returns:

A 2D array with a copy of the data from the selected extension.

Return type:

numpy.ndarray

static default_paths(cube, method_key, rdxqa_method, binning_method, stelcont_method=None, output_path=None, output_file=None)[source]

Set the default directory and file name for the output file.

Parameters:
  • cube (mangadap.datacube.datacube.DataCube) – Datacube to analyze.

  • method_key (str) – Keyword designating the method used for the reduction assessments.

  • rdxqa_method (str) – The method key for the basic assessments of the datacube.

  • binning_method (str) – The method key for the spatial binning.

  • stelcont_method (str, optional) – The method key for the stellar-continuum fitting method. If None, not included in output file name.

  • output_path (str, Path, optional) – The path for the output file. If None, the current working directory is used.

  • output_file (str, optional) – The name of the output “reference” file. The full path of the output file will be directory_path/output_file. If None, the default is to combine cube.output_root and the method keys. The order of the keys is the order of operations (rdxqa, binning, stellar continuum, emission-line model).

Returns:

Returns a Path with the output directory and a str with the output file name.

Return type:

tuple

file_name()[source]

Return the name of the output file.

file_path()[source]

Return the full path to the output file.

fill_continuum_to_match(binid, replacement_templates=None, redshift_only=False, deredshift=False, corrected_dispersion=False, missing=None)[source]

Get the emission-line continuum model, if possible, that matches the input bin ID matrix. The output is a 2D matrix ordered by the bin ID; any skipped index numbers in the maximum of the union of the unique numbers in the binid and missing input are masked.

Use replacement_templates only if the number of templates is identical to the number used during the fitting procedure. Use redshift_only to produce the best-fitting model without and velocity dispersion. Use corrected_dispersion to produce the model using the corrected velocity dispersion.

fill_to_match(binid, include_base=False, missing=None)[source]

Get the emission-line model that matches the input bin ID matrix. The output is a 2D matrix ordered by the bin ID; any skipped index numbers in the maximum of the union of the unique numbers in the binid and missing input are masked.

Use include_base to include the baseline in the output emission-line models.

fit(binned_spectra, stellar_continuum=None, emission_line_moments=None, redshift=None, dispersion=None, minimum_error=None, output_path=None, output_file=None, hardcopy=True, overwrite=None, loggers=None, quiet=False)[source]

Fit the emission lines.

fit_figures_of_merit()[source]

Use the fitting class to set the figures-of-merit of the full fit.

This is primarily used to set the output data for the MAPS file.

Returns:

Returns (1) a list of names to assign each column (length is NFOM), (2) the units of the data in this column (length is NFOM), and (3) a 2D array with the figures of merit for each spectrum (shape is NSPEC x NFOM). If the fitting class used does not have a fit_figures_of_merit function, the name is set to ‘null’, the unit is an empty string, and the data is a (NMOD,1) array of zeros.

Return type:

tuple

matched_kinematics(binid, line_name, redshift=None, dispersion=100.0, constant=False, cz=False, corrected=False, min_dispersion=None, nearest=False, missing=None)[source]

Return the redshift and velocity dispersion for all the spectra based on the emission-line model fit to a single line.

For spectra the were not fit, use the median (unmasked) kinematic measurement if no default value is provided as an argument.

Parameters:
  • binid (numpy.ndarray) – 2D array with the bin ID associate with each spaxel in the datacube.

  • line (str) – Line to use for the kinematics. Must match one of the lines created by channel_names(). Function will raise an exception if the line_name is None or if the channel names cannot be constructed.

  • redshift (float) – (Optional) The default redshift to use for spectra without an emission-line fit. Default is to use the median of the unmasked measurements.

  • dispersion (float) – (Optional) The default velocity dispersion to use for spectra without an emission-line fit. Default is 100 km/s.

  • constant (bool) – (Optional) Force the function to return a constant redshift and dispersion for each spectrum, regardless of any fitted kinematics.

  • cz (bool) – (Optional) Return the redshift as cz in km/s, as opposed to unitless z.

  • corrected (bool) – (Optional) By default, the returned velocity dispersions are the measured values, which will include the instrumental resolution. Setting corrected to True will return the corrected dispersions.

  • min_dispersion (float) – (Optional) Impose a minimum dispersion.

  • nearest (bool) – (Optional) Instead of the median of the results for the spectra that were not fit, use the value from the nearest bin.

  • missing (list) – (Optional) Any bin ID numbers missing from the input binid image needed for constructing the output matched data.

Returns:

Returns arrays with a redshift (or cz) and dispersion for each spectrum.

Return type:

numpy.ndarray

Raises:

TypeError – Raised if the input redshift or dispersion values are not single numbers.

read(ifile=None, strict=True, checksum=False, debug=False)[source]

Read an existing file with an existing set of emission-line models.

write(match_DRP=False, overwrite=False)[source]

Write the hdu object to the file.

class mangadap.proc.emissionlinemodel.EmissionLineModelBitMask[source]

Bases: DAPBitMask

Derived class that specifies the mask bits for the emission-line modeling. The maskbits defined are:

Key

Bit

Description

DIDNOTUSE

0

Pixel was ignored because it was flagged as DONOTUSE or FORESTAR by the DRP, or as LOW_SPECCOV, LOW_SNR, or NONE_IN_STACK in the binning step.

FORESTAR

1

Pixel was ignored because it was flagged as FORESTAR by the DRP.

LOW_SNR

2

Pixel was ignored because the S/N estimate of the spectrum was below the set threshold; see header keyword ELFMINSN.

ARTIFACT

3

Pixel was ignored during the emission-line model fit because it was designated as containing an artifact.

OUTSIDE_RANGE

4

Pixel was ignored during the emission-line model fit because it was outside of any emission-line fitting window.

NOCONTINUUM

5

Pixel did not have any model of the stellar-continuum subtracted.

INSUFFICIENT_DATA

6

There were insufficient data over the relevant wavelength range to fit the line profile(s).

FIT_FAILED

7

Emission-line fit failed according to status returned by scipy.optimize.least_squares.

NEAR_BOUND

8

Emission-line parameter(s) are within at or near the imposed boundaries in parameter space to within the error of the parameter.

UNDEFINED_COVAR

9

Emission-line fit resulted in a covariance matrix with negative values along the diagonal, meaning that the formal error is undefined. (Need to test, but this should mean that the fit is bad.)

EXCLUDED_FROM_MODEL

10

Emission-line fit excluded from the emission-line flux model because the input emission-line database requested it be one of the lines in the fitting window be excluded or one of the following bits were set: INSUFFICIENT_DATA, FIT_FAILED, NEAR_BOUND, UNDEFINED_COVAR.

UNDEFINED_SIGMA

11

Emission-line velocity dispersion is undefined because it was found to be smaller than the instrumental dispersion.

NON_POSITIVE_CONTINUUM

12

Equivalent width measurements were not computed because the continuum in either the blue or red sidebands was not positive.

NO_FIT

13

Emission-line or spectrum was not fit.

TPL_PIXELS

14

These pixels were removed to ensure that the number of template spectral pixels was >= the number of fitted object pixels during the stellar-continuum fit.

TRUNCATED

15

This region was truncated to avoid convolution errors at the edges of the spectral range during the stellar-continuum fit.

PPXF_REJECT

16

Pixel was rejected during the pPXF fit via the clean parameter.

MIN_SIGMA

17

The fitted velocity dispersion is at the minimum allowed by pPXF (1/100th of a pixel)

BAD_SIGMA

18

Corrected velocity dispersion is below 0 km/s or above 400 km/s

MAXITER

19

The fit optimizer reached the maximum number of iterations during the fit, which may imply failure to converge.

INVALID_ERROR

20

Pixel ignored because the flux error was invalid.

EML_REGION

21

Pixel was ignored during the stellar-continuum fit because it contains an emission-line.

cfg_root = 'emission_line_model_bits'
class mangadap.proc.emissionlinemodel.EmissionLineModelDef(key='EFITMPL11SSPDB', minimum_snr=0.0, mom_vel_name='Ha-6564', mom_disp_name=None, fitpar=None, fitclass=None, fitfunc=None, overwrite=False)[source]

Bases: KeywordParSet

A class that holds the parameters necessary to perform the emission-line profile fits.

The defined parameters are:

Key

Type

Options

Default

Description

key

str

EFITMPL11SSPDB

Keyword used to distinguish between different emission-line moment databases.

minimum_snr

int, float

0.0

Minimum S/N of spectrum to fit

mom_vel_name

str

Ha-6564

Name of the emission-line moments band used to set the initial velocity guess for each spaxel.

mom_disp_name

str

Name of the emission-line moments band used to set the initial velocity dispersion guess for each spaxel.

fitpar

ParSet, dict

Fitting function parameters

fitclass

Undefined

Class used to perform the fit.

fitfunc

Undefined

Function or method that performs the fit; must be callable.

overwrite

bool

False

If the output file already exists, redo all the calculations and overwrite it.

_validate()[source]

Validate the modeling method.

classmethod from_dict(d)[source]

Instantiate the object from a nested dictionary.

The dictionary is expected to have keys for all the class keywords; see the class instantiation method. It should also include the 'fit_method' keyword used to select the type of fit to perform, and a nested dictionary selected by the keyword 'fit' with the parameters required by the fitting method. If these are not provided, the default fitter, Sasuke and parameters will be used.

Warning

Currently, the only accepted value for 'fit_method' is 'sasuke', and the 'fit' dictionary must provide the parameters defined by SasukePar.

Parameters:

d (dict) – Dictionary used to instantiate the class.