____ _ _________________________________________ \/ \ _/\ / \/ \ / oigtFit \/
VoigtFit Interface¶
class DataSet¶
-
class
VoigtFit.
DataSet
(redshift, name='')[source]¶ Main class of the package
VoigtFit
. The DataSet handles all the major parts of the fit. Spectral data must be added using theadd_data
method. Hereafter the absorption lines to be fitted are added to the DataSet using theadd_line
oradd_many_lines
methods. Lastly, the components of each element is defined using theadd_component
method. When all lines and components have been defined, the DataSet must be prepared by calling theprepare_dataset
method and subsequently, the lines can be fitted using thefit
method.Attributes
- redshift : float
- Systemic redshift of the absorption system.
- name : str [default = ‘’]
- The name of the DataSet, this will be used for saving the dataset to a file structure.
- verbose : bool [default = True]
- If False, the printed information statements will be suppressed.
- data : list(data_chunks)
- A list of data chunks defined for the dataset. A data chunk is
a dictionary with keys ‘wl’, ‘flux’, ‘error’, ‘res’, ‘norm’.
See
DataSet.add_data
. - lines : dict
- A dictionary holding pairs of defined (line_tag :
dataset.Line
) - all_lines : list(str)
- A list of all the defined line tags for easy look-up.
- molecules : dict
- A dictionary holding a list of the defined molecular bands and Jmax
for each molecule:
{molecule* : [[band1, Jmax1], [band2, Jmax2], etc...]}
- regions : list(
regions.Region
) - A list of the fitting regions.
- cheb_order : int [default = -1]
- The maximum order of Chebyshev polynomials to use for the continuum fitting in each region. If negative, the Chebyshev polynomials will not be included in the fit.
- norm_method : str [default = ‘linear’]
- Default normalization method to use for interactive normalization if Chebyshev polynomial fitting should not be used.
- components : dict
- A dictionary of components for each ion defined:
(ion : [z, b, logN, options]). See
DataSet.add_component
. - velspan : float [default = 500]
- The default velocity range to use for the definition of fitting regions.
- ready2fit : bool [default = False]
- This attribute is checked before fitting the dataset. Only when
the attribute has been set to True can the dataset be fitted.
This will be toggled after a successful run of
DataSet.prepare_dataset
. - best_fit : lmfit.Parameters [default = None]
- Best-fit parameters from lmfit. This attribute will be None until the dataset has been fitted.
- pars : lmfit.Parameters [default = None]
- Placeholder for the fit parameters initiated before the fit.
The parameters will be defined during the call to
DataSet.prepare_dataset
based on the defined components.
Methods
activate_all
()Activate all lines defined in the DataSet. activate_fine_lines
(line_tag[, levels])Activate all lines associated to a given fine-structure complex. activate_line
(line_tag)Activate a given line defined by its line_tag activate_molecule
(molecule, band)Activate all lines for the given band of the given molecule. add_component
(ion, z, b, logN[, var_z, …])Add component for a given ion. add_component_velocity
(ion, v, b, logN[, …])Same as for add_component
but input is given as relative velocity instead of redshift.add_data
(wl, flux, res[, err, normalized, …])Add spectral data to the DataSet. add_fine_lines
(line_tag[, levels, …])Add fine-structure line complexes by providing only the main transition. add_line
(line_tag[, velspan, active])Add an absorption line to the DataSet. add_lines
(line_tags[, velspan])Alias for self.add_many_lines. add_many_lines
(tags[, velspan])Add many lines at once. add_molecule
(molecule, band[, Jmax, …])Add molecular lines for a given band, e.g., “AX(0-0)” of CO. add_spectrum
(filename, res[, normalized, …])Add spectral data from ASCII file (three or four columns accepted). all_active_lines
()Returns a list of all the active lines defined by their line_tag. clear_mask
(line_tag[, idx])Clear the mask for the Region
containing the given line_tag.copy_components
([to_ion, from_ion, logN, …])Copy velocity structure to ion from another ion. deactivate_all
()Deactivate all lines defined in the DataSet. deactivate_fine_lines
(line_tag[, levels, …])Deactivate all lines associated to a given fine-structure complex. deactivate_line
(line_tag)Deactivate a given line defined by its line_tag. deactivate_molecule
(molecule, band)Deactivate all lines for the given band of the given molecule. delete_component
(ion, index)Remove component of the given ion with the given index. find_ion
(ion)Return a list of all line tags for a given ion. find_line
(line_tag)Look up the fitting Region
for a given line tag.fit
([verbose])Fit the absorption lines using chi-square minimization. fix_structure
([ion])Fix the velocity structure, that is, the redshifts and the b-parameters. free_structure
([ion])Free the velocity structure, that is, the redshifts and the b-parameters. get_lines_for_ion
(ion)Return a list of Line
objects corresponding to the given ion.get_name
()Returns the name of the DataSet. get_resolution
(line_tag[, verbose])Return the spectral resolution for the fitting Region
where the line with the given line_tag is defined.has_ion
(ion[, active_only])Return True if the dataset has lines defined for the given ion. has_line
(line_tag[, active_only])Return True if the given line is defined. interactive_components
(line_tag[, velocity])Define components interactively for a given ion. load_components_from_file
(fname[, fit_pars])Load best-fit parameters from an output file fname. mask_line
(line_tag[, reset, mask, telluric, …])Define exclusion masks for the fitting region of a given line. mask_range
(line_tag, x1, x2[, idx])Define mask in a range from x1 to x2 in velocity space. normalize_line
(line_tag[, norm_method, velocity])Normalize or re-normalize a given line plot_fit
([rebin, fontsize, xmin, xmax, …])Plot all the absorption lines and the best-fit profiles. plot_line
(line_tag[, index, plot_fit, loc, …])Plot a single fitting Region
containing the line corresponding to the given line_tag.prepare_dataset
([norm, mask, verbose, …])Prepare the data for fitting. print_cont_parameters
()Print Chebyshev coefficients for the continuum fit. print_metallicity
(logNHI[, err])Print the total column densities for each element relative to HI in Solar units. print_results
([velocity, elements, systemic])Print the best fit parameters. print_total
()Print the total column densities of all components. remove_fine_lines
(line_tag[, levels])Remove lines associated to a given fine-structure complex. remove_line
(line_tag)Remove an absorption line from the DataSet. remove_molecule
(molecule, band)Remove all lines for the given band of the given molecule. reset_all_regions
()Reset the data in all Regions
defined in the DataSet to use the raw input data.reset_components
([ion])Reset component structure for a given ion. reset_region
(reg)Reset the data in a given regions.Region
to use the raw input data.save
([filename, verbose])Save the DataSet to file using the HDF5 format. save_cont_parameters_to_file
(filename)Save the best-fit continuum parameters to ASCII table output. save_fit_regions
([filename, individual, path])Save the fitting regions to ASCII table output. save_parameters
(filename)Save the best-fit parameters to ASCII table output. set_name
(name)Set the name of the DataSet. set_resolution
(res[, line_tag, verbose, nsub])Set the spectral resolution in km/s for the Region
containing the line with the given line_tag.set_systemic_redshift
(z_sys)Update the systemic redshift of the dataset show_components
([ion])Show the defined components for a given ion. show_lines
()Print all defined lines to terminal. sum_components
(ions[, components])Calculate the total column density for the given components of the given ion. velocity_plot
(**kwargs)Create a velocity plot, showing all the fitting regions defined, in order to compare different lines and to identify blends and contamination. get_NHI remove_all_lines -
activate_fine_lines
(line_tag, levels=None)[source]¶ Activate all lines associated to a given fine-structure complex.
Parameters: - line_tag : str
The line tag of the ground state transition to activate.
- levels : str, list(str) [default = None]
The levels of the fine-structure complexes to activate, with the string “a” referring to the first excited level, “b” is the second, etc… Several levels can be given at once as a list: [‘a’, ‘b’] or as a concatenated string: ‘abc’. By default, all levels are included.
-
activate_molecule
(molecule, band)[source]¶ Activate all lines for the given band of the given molecule.
- Ex:
activate_molecule('CO', 'AX(0-0)')
- Ex:
-
add_component
(ion, z, b, logN, var_z=True, var_b=True, var_N=True, tie_z=None, tie_b=None, tie_N=None)[source]¶ Add component for a given ion. Each component defined will be used for all transitions defined for a given ion.
Parameters: - ion : str
The ion for which to define a component: e.g., “FeII”, “HI”, “CIa”, etc.
- z : float
The redshift of the component.
- b : float
The effective broadening parameter for the component in km/s. This parameter is constrained to be in the interval [0 - 1000] km/s.
- logN : float
The 10-base logarithm of the column density of the component. The column density is expected in cm^-2.
- var_z : bool [default = True]
If False, the redshift of the component will be kept fixed.
- var_b : bool [default = True]
If False, the b-parameter of the component will be kept fixed.
- var_N : bool [default = True]
If False, the column density of the component will be kept fixed.
- tie_z, tie_b, tie_N : str [default = None]
Parameter constraints for the different variables.
Notes
The ties are defined relative to the parameter names. The naming is as follows: The redshift of the first component of FeII is called “z0_FeII”, the logN of the second component of SiII is called “logN1_SiII”. For more information about parameter ties, see the documentation for lmfit.
-
add_component_velocity
(ion, v, b, logN, var_z=True, var_b=True, var_N=True, tie_z=None, tie_b=None, tie_N=None)[source]¶ Same as for
add_component
but input is given as relative velocity instead of redshift.
-
add_data
(wl, flux, res, err=None, normalized=False, mask=None, nsub=1)[source]¶ Add spectral data to the DataSet. This will be used to define fitting regions.
Parameters: - wl : ndarray, shape (n)
Input vacuum wavelength array in Ångstrøm.
- flux : ndarray, shape (n)
Input flux array, should be same length as wl
- res : float or string
Spectral resolution either given in km/s (c/R), which is assumed to be constant over the whole spectrum, or as a string referring to a file containing the detailed line-spread function for the given spectrum. See details in the data section of the VoigtFit Parameter Language.
- err : ndarray, shape (n) [default = None]
Error array, should be same length as wl If None is given, an uncertainty of 1 is given to all pixels.
- normalized : bool [default = False]
If the input spectrum is normalized this should be given as True in order to skip normalization steps.
- mask : array, shape (n)
Boolean/int array defining the fitting regions. Only pixels with mask=True/1 will be included in the fit.
- nsub : integer
Kernel subsampling factor relative to the data. This is only used if the resolution is given as a LSF file.
-
add_fine_lines
(line_tag, levels=None, full_label=False, velspan=None)[source]¶ Add fine-structure line complexes by providing only the main transition. This function is mainly useful for the CI complexes, where the many lines are closely located and often blended.
Parameters: - line_tag : str
Line tag for the ground state transition, e.g., “CI_1656”
- levels : str, list(str) [default = None]
The levels of the fine-structure complexes to add, starting with “a” referring to the first excited level, “b” is the second, etc.. Several levels can be given at once: [‘a’, ‘b’]. Note that the ground state transition is always included. If levels is not given, all levels are included.
- full_label : bool [default = False]
If True, the label will be translated to the full quantum mechanical description of the state.
-
add_line
(line_tag, velspan=None, active=True)[source]¶ Add an absorption line to the DataSet.
Parameters: - line_tag : str
The line tag for the transition which should be defined. Ex: “FeII_2374”
- velspan : float [default = None]
The velocity span around the line center, which will be included in the fit. If None is given, use the default self.velspan defined (default = 500 km/s).
- active : bool [default = True]
Set the
Line
as active (i.e., included in the fit).
Notes
This will initiate a
Line
class with the atomic data for the transition, as well as creating a fittingRegion
containing the data cutout around the line center.
-
add_many_lines
(tags, velspan=None)[source]¶ Add many lines at once.
Parameters: - tags : list(str)
A list of line tags for the transitions that should be added.
- velspan : float [default = None]
The velocity span around the line center, which will be included in the fit. If None is given, use the default
velspan
defined (500 km/s).
-
add_molecule
(molecule, band, Jmax=0, velspan=None, full_label=False)[source]¶ Add molecular lines for a given band, e.g., “AX(0-0)” of CO.
Parameters: - molecule : str
The molecular identifier, e.g., ‘CO’, ‘H2’
- band : str
The vibrational band of the molecule, e.g., for CO: “AX(0-0)” These bands are defined in
molecules
.- Jmax : int [default = 0]
The maximal rotational level to include. All levels up to and including J will be included.
- velspan : float [default = None]
The velocity span around the line center, which will be included in the fit. If None is given, use the default
velspan
defined (500 km/s).- full_label : bool [default = False]
If True, the label will be translated to the full quantum mechanical description of the state.
-
add_spectrum
(filename, res, normalized=False, mask=None, continuum=None, nsub=1)[source]¶ Add spectral data from ASCII file (three or four columns accepted). This is a wrapper of the method add_data.
Parameters: - filename : string
Filename of the ASCII spectrum containing at least 3 columns: (wavelength, flux, error). A fourth column may be included containing a spectral mask.
- res : float or string
Spectral resolution either given in km/s (c/R), which is assumed to be constant over the whole spectrum, or as a string referring to a file containing the detailed line-spread function for the given spectrum. See details in the documentation.
- normalized : bool [default = False]
If the input spectrum is normalized this should be given as True in order to skip normalization steps.
- mask : array, shape (n)
Boolean/int array defining the fitting regions. Only pixels with mask=True/1 will be included in the fit.
- continuum : array, shape (n)
Continuum spectrum. The input spectrum will be normalized using this continuum spectrum.
- nsub : integer
Kernel subsampling factor relative to the data. This is only used if the resolution is given as a LSF file.
-
clear_mask
(line_tag, idx=None)[source]¶ Clear the mask for the
Region
containing the given line_tag. If more regions are defined for the same line (when fitting multiple spectra), the given region can be specified by passing an index idx.
-
copy_components
(to_ion='', from_ion='', logN=0, ref_comp=None, tie_z=True, tie_b=True)[source]¶ Copy velocity structure to ion from another ion.
Parameters: - to_ion : str
The new ion to define.
- from_ion : str
The base ion which will be used as reference.
- logN : float
If logN is given, the starting guess is defined from this value following the pattern of the components defined for anchor relative to the ref_comp (default: the first component).
- ref_comp : int
The reference component to which logN will be scaled.
- tie_z : bool [default = True]
If True, the redshifts for all components of the two ions will be tied together.
- tie_b : bool [default = True]
If True, the b-parameters for all components of the two ions will be tied together.
-
deactivate_all
()[source]¶ Deactivate all lines defined in the DataSet. This will not remove the lines.
-
deactivate_fine_lines
(line_tag, levels=None, verbose=True)[source]¶ Deactivate all lines associated to a given fine-structure complex.
Parameters: - line_tag : str
The line tag of the ground state transition to deactivate.
- levels : str, list(str) [default = None]
The levels of the fine-structure complexes to deactivate, with the string “a” referring to the first excited level, “b” is the second, etc… Several levels can be given at once as a list: [‘a’, ‘b’] or as a concatenated string: ‘abc’. By default, all levels are included.
-
deactivate_line
(line_tag)[source]¶ Deactivate a given line defined by its line_tag. This will exclude the line during the fit but will not remove the data.
-
deactivate_molecule
(molecule, band)[source]¶ Deactivate all lines for the given band of the given molecule. To see the available molecular bands defined, see the manual pdf or
VoigtFit.molecules
.
-
find_line
(line_tag)[source]¶ Look up the fitting
Region
for a given line tag.Parameters: - line_tag : str
The line tag of the line whose region will be returned.
Returns: - regions_of_line : list of
Region
A list of the fitting regions containing the given line. This can be more than one region in case of overlapping or multiple spectra.
-
fit
(verbose=True, **kwargs)[source]¶ Fit the absorption lines using chi-square minimization.
Parameters: - verbose : bool [default = True]
This will print the fit results to terminal.
- plot : bool [default = False]
This will make the best-fit solution show up in a new window.
- **kwargs
Keyword arguments are derived from the `scipy.optimize`_ minimization methods. The default method is leastsq_, used in lmfit. This can be changed using the method keyword. See documentation in lmfit and scipy.optimize_.
- rebin : int [default = 1]
Rebin data by a factor rebin before fitting. Passed as part of kwargs.
- sampling : int [default = 3]
Subsampling factor for the evaluation of the line profile. This is only used if the kernel is constant in velocity along the spectrum. Passed as part of kwargs.
Returns: - popt : lmfit.MinimizerResult_
The minimzer results from lmfit containing best-fit parameters and fit details, e.g., exit status and reduced chi squared. See documentation for lmfit.
- chi2 : float
The chi squared value of the best-fit. Note that this value is not the reduced chi squared. This value, and the number of degrees of freedom, are available under the popt object.
- .. _scipy.optimize: https://docs.scipy.org/doc/scipy/reference/tutorial/optimize.html
- .. _leastsq: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.leastsq.html
- .. _lmfit.MinimizerResult: https://lmfit.github.io/lmfit-py/fitting.html#lmfit.minimizer.MinimizerResult
-
fix_structure
(ion=None)[source]¶ Fix the velocity structure, that is, the redshifts and the b-parameters.
Parameters: - ion : str [default = None]
The ion for which the structure should be fixed. If None is given, the structure is fixed for all ions.
-
free_structure
(ion=None)[source]¶ Free the velocity structure, that is, the redshifts and the b-parameters.
Parameters: - ion : str [default = None]
The ion for which the structure should be released. If None is given, the structure is released for all ions.
-
get_lines_for_ion
(ion)[source]¶ Return a list of
Line
objects corresponding to the given ion.Parameters: - ion : str
The ion for the lines to get.
Returns: - lines_for_ion : list(
Line
) List of Lines defined for the given ion.
-
get_resolution
(line_tag, verbose=False)[source]¶ Return the spectral resolution for the fitting
Region
where the line with the given line_tag is defined.Parameters: - line_tag : str
The line-tag for the line to look up: e.g., “FeII_2374”
- verbose : bool [default = False]
If True, print the returned spectral resolution to std out.
Returns: - resolutions : list of float
A list of the spectral resolution of the fitting regions where the given line is defined.
-
has_ion
(ion, active_only=False)[source]¶ Return True if the dataset has lines defined for the given ion.
-
interactive_components
(line_tag, velocity=False)[source]¶ Define components interactively for a given ion. The components will be defined on the basis of the given line for that ion. If the line is defined in several spectra then the interactive window will show up for each. Running the interactive mode more times for different transitions of the same ion will append the components to the structure. If no components should be added, then simply click enter to terminate the process for the given transition.
Parameters: - line_tag : str
Line tag for the line belonging to the ion for which components should be defined.
- velocity : bool [default = False]
If a True, the region is displayed in velocity space relative to the systemic redshift instead of in wavelength space.
Notes
This will launch an interactive plot showing the fitting region of the given line. The user can then click on the positions of the components which. At the end, the redshifts and estimated column densities are printed to terminal. The b-parameter is assumed to be unresolved, i.e., taken from the resolution.
-
load_components_from_file
(fname, fit_pars=True)[source]¶ Load best-fit parameters from an output file fname. If fit_pars is True, then update the best_fit parameters.
-
mask_line
(line_tag, reset=True, mask=None, telluric=True, velocity=False)[source]¶ Define exclusion masks for the fitting region of a given line. Note that the masked regions are exclusion regions and will not be used for the fit. If components have been defined, these will be shown as vertical lines.
Parameters: - line_tag : str
Line tag for the
Line
whoseRegion
should be masked.- reset : bool [default = True]
If True, clear the mask before defining a new mask.
- mask : array_like, shape (n) [default = None]
If the mask is given, it must be a boolean array of the same length as the region flux, err, and wl arrays. Passing a mask this was supresses the interactive masking process.
- telluric : bool [default = True]
If True, a telluric absorption template and sky emission template is shown for reference.
- velocity : bool [default = False]
If a True, the regions are displayed in velocity space relative to the systemic redshift instead of in wavelength space when masking and defining continuum normalization interactively.
-
mask_range
(line_tag, x1, x2, idx=None)[source]¶ Define mask in a range from x1 to x2 in velocity space.
-
normalize_line
(line_tag, norm_method='spline', velocity=False)[source]¶ Normalize or re-normalize a given line
Parameters: - line_tag : str
Line tag of the line whose fitting region should be normalized.
- norm_method : str [default = ‘spline’]
Normalization method used for the interactive continuum fit. Should be on of: [“spline”, “linear”]
- velocity : bool [default = False]
If a True, the regions are displayed in velocity space relative to the systemic redshift instead of in wavelength space when masking and defining continuum normalization interactively.
-
plot_fit
(rebin=1, fontsize=12, xmin=None, xmax=None, max_rows=4, ymin=None, ymax=None, filename=None, subsample_profile=1, npad=50, loc='left', highlight_props=None, residuals=True, norm_resid=False, default_props={}, element_props={}, legend=True, label_all_ions=False, xunit='vel', line_props=None, hl_line_props=None)[source]¶ Plot all the absorption lines and the best-fit profiles. For details, see
VoigtFit.output.plot_all_lines()
.
-
plot_line
(line_tag, index=0, plot_fit=False, loc='left', rebin=1, nolabels=False, axis=None, fontsize=12, xmin=None, xmax=None, ymin=None, ymax=None, show=True, subsample_profile=1, npad=50, residuals=True, norm_resid=False, legend=True, default_props={}, element_props={}, highlight_props=None, label_all_ions=False, xunit='velocity', line_props=None, hl_line_props=None)[source]¶ Plot a single fitting
Region
containing the line corresponding to the given line_tag. For details, seeoutput.plot_single_line()
.
-
prepare_dataset
(norm=True, mask=False, verbose=True, active_only=False, force_clean=True, velocity=False, check_lines=True, f_lower=0.0, f_upper=100.0, l_lower=0.0, l_upper=10000.0)[source]¶ Prepare the data for fitting. This function sets up the parameter structure, and handles the normalization and masking of fitting regions.
Parameters: - norm : bool [default = True]
Opens an interactive window to let the user normalize each region using the defined
norm_method
.- mask : bool [default = True]
Opens an interactive window to let the user mask each fitting region.
- verbose : bool [default = True]
If this is set, the code will print small info statements during the run.
- force_clean : bool [default = False]
If this is True, components for inactive elements will be removed.
- velocity : bool [default = False]
If True, the regions are displayed in velocity space relative to the systemic redshift instead of in wavelength space when masking and defining continuum normalization interactively.
- check_lines : bool [default = True]
If True, all available lines covered by the data will be checked. The user will be propmted if lines are available for ions that have already been defined.
- f_lower : float [default = 0.]
Lower limit on oscillator strengths for transitions when verifying all transitions for defined ions.
- f_upper : float [default = 100.]
Upper limit on oscillator strengths for transitions when verifying all transitions for defined ions.
- l_lower : float [default = 0.]
Lower limit on rest-frame wavelength for transitions when verifying all transitions for defined ions.
- l_upper : float [default = 1.e4]
Upper limit on rest-frame wavelength for transitions when verifying all transitions for defined ions.
Returns: - bool
The function returns True when the dataset has passed all the steps. If one step fails, the function returns False. The
ready2fit
attribute of the dataset is also updated accordingly.
-
print_metallicity
(logNHI, err=0.1)[source]¶ Print the total column densities for each element relative to HI in Solar units.
-
print_results
(velocity=True, elements='all', systemic=None)[source]¶ Print the best fit parameters.
Parameters: - velocity : bool [default = True]
If True, show the relative velocities of each component instead of redshifts.
- elements : list(str) [default = ‘all’]
A list of elements for which to show parameters.
- systemic : float [default = None]
The systemic redshift used as reference for the relative velocities.
-
remove_fine_lines
(line_tag, levels=None)[source]¶ Remove lines associated to a given fine-structure complex.
Parameters: - line_tag : str
The line tag of the ground state transition to remove.
- levels : str, list(str) [default = None]
The levels of the fine-structure complexes to remove, with “a” referring to the first excited level, “b” is the second, etc.. Several levels can be given at once as a list: [‘a’, ‘b’] or as a concatenated string: ‘abc’. By default, all levels are included.
-
remove_line
(line_tag)[source]¶ Remove an absorption line from the DataSet. If this is the last line in a fitting region the given region will be eliminated, and if this is the last line of a given ion, then the components will be eliminated for that ion.
Parameters: - line_tag : str
Line tag of the transition that should be removed.
-
reset_all_regions
()[source]¶ Reset the data in all
Regions
defined in the DataSet to use the raw input data.
-
reset_components
(ion=None)[source]¶ Reset component structure for a given ion.
Parameters: - ion : str [default = None]
The ion for which to reset the components: e.g., ‘FeII’, ‘HI’, ‘CIa’, etc. If `None`is given, all components for all ions will be reset.
-
reset_region
(reg)[source]¶ Reset the data in a given
regions.Region
to use the raw input data.
-
save_cont_parameters_to_file
(filename)[source]¶ Save the best-fit continuum parameters to ASCII table output.
Parameters: - filename : str [default = None]
If None, the
name
attribute will be used.
-
save_fit_regions
(filename=None, individual=False, path='')[source]¶ Save the fitting regions to ASCII table output. The format is as follows: (wavelength , normalized flux , normalized error , best-fit profile , mask)
Parameters: - filename : str [default = None]
Filename for the fitting regions. If None, the
name
attribute will be used.- individual : bool [default = False]
Save the fitting regions to individual files. By default all regions are concatenated into one file.
- path : str [default = ‘’]
Specify a path to prepend to the filename in order to save output to a given directory or path. Can be given both as relative or absolute path. If the path doesn’t end in / it will be appended automatically. The final filename will be:
path/ + filename [+ _regN] + .reg
-
save_parameters
(filename)[source]¶ Save the best-fit parameters to ASCII table output.
Parameters: - filename : str
Filename for the fit parameter file.
-
set_name
(name)[source]¶ Set the name of the DataSet. This parameter is used when saving the dataset.
-
set_resolution
(res, line_tag=None, verbose=True, nsub=1)[source]¶ Set the spectral resolution in km/s for the
Region
containing the line with the given line_tag. If multiple spectra are fitted simultaneously, this method will set the same resolution for all spectra. If line_tag is not given, the resolution will be set for all regions, including the raw data chunks defined inVoigtFit.DataSet.data
!Note – If not all data chunks have the same resolution, this method should be used with caution. It is advised to check the spectral resolution beforehand and only update the appropriate regions using a for-loop.
-
show_components
(ion=None)[source]¶ Show the defined components for a given ion. By default, all ions are shown.
-
show_lines
()[source]¶ Print all defined lines to terminal. The output shows whether the line is active or not and the number of components for the given ion.
-
sum_components
(ions, components=None)[source]¶ Calculate the total column density for the given components of the given ion.
Parameters: - ions : str or list(str)
List of ions or a single ion for which to calculate the total column density.
- components : list(int)
List of integers corresponding to the indeces of the components to sum over.
Returns: - total_logN : dict()
Dictionary containing the log of total column density for each ion.
- total_logN_err : dict()
Dictionary containing the error on the log of total column density for each ion.
class Component¶
-
class
components.
Component
(z, b, logN, var_z=True, var_b=True, var_N=True, tie_z=None, tie_b=None, tie_N=None)[source]¶ Component object which contains the parameters for each velocity component of an absorption profile: redshift, z; broadning parameter, b; and column density in log(cm^-2). Options can control whether the components parameters are variable during the fit or whether they are tied to another parameter using the var_ and tie_ options.
Attributes: - b
- logN
- tie_N
- tie_b
- tie_z
- var_N
- var_b
- var_z
- z
Methods
get_option
(key)Return the value for the given option, key must be either tie_ or var_. get_pars
()Unpack the physical parameters [z, b, logN] set_option
(key, value)Set the value for a given option, must be either tie_ or var_.
class Line¶
-
class
lines.
Line
(tag, active=True)[source]¶ Line object containing atomic data for the given transition. Only the line_tag is passed, the rest of the information is looked up in the atomic database.
Attributes
- tag : str
- The line tag for the line, e.g., “FeII_2374”
- ion : str
- The ion for the line; The ion for “FeII_2374” is “FeII”.
- element : str
- Equal to
ion
for backwards compatibility. - l0 : float
- Rest-frame resonant wavelength of the transition. Unit: Angstrom.
- f : float
- The oscillator strength for the transition.
- gam : float
- The radiation damping constant or Einstein coefficient.
- mass : float
- The atomic mass in atomic mass units.
- active : bool [default = True]
- The state of the line in the dataset. Only active lines will be included in the fit.
Methods
get_properties
()Return the principal atomic constants for the transition: l0, f, and gam. set_active
()Set the line active; include the line in the fit. set_inactive
()Set the line inactive; exclude the line in the fit.
class Region¶
-
class
regions.
Region
(velspan, specID, line=None)[source]¶ A Region contains the fitting data, exclusion mask and line information. The class is instantiated with the velocity span, velspan, and a spectral ID pointing to the raw data chunk from DataSet.data, and can include a
dataset.Line
instance for the first line belonging to the region.Attributes
- velspan : float
- The velocity range to used for the fitting region.
- lines : list(
dataset.Line
) - A list of Lines defined in the region.
- label : str
- A LaTeX label describing the lines in the region for plotting purposes.
- res : float
- Spectral resolution of the region in km/s.
- wl : array_like, shape (N)
- Data array of wavelengths in Ångstrøm.
- flux : array_like, shape (N)
- Data array of fluxes (normalized if
normalized
is True). - err : array_like, shape (N)
- Array of uncertainties for each flux element.
- normalized : bool
- True if the data in the region are normlized.
- mask : array_like, shape (N)
- Exclusion mask for the region: 0/False = pixel is not included in the fit. 1/True = pixel is included in the fit.
- new_mask : bool
- Internal parameter for
VoigtFit.DataSet.prepare_dataset()
. If True, an interactive masking process will be initiated in the preparation stage. - cont_err : float
- An estimate of the uncertainty in the continuum fit.
- specID : str
- A spectral identifier to point back to the raw data chunk.
Methods
add_data_to_region
(data_chunk, cutout)Define the spectral data for the fitting region. add_line
(line)Add a new dataset.Line
to the fitting region.clear_mask
()Clear the already defined mask in the region. define_mask
([z, dataset, telluric, z_sys])Use an interactive window to define the mask for the region. generate_label
([active_only, ignore_finelines])Automatically generate a descriptive label for the region. has_active_lines
()Return True is at least one line in the region is active. has_line
(line_tag)Return True if a line with the given line_tag is defined in the region. is_normalized
()Return True if the region data is normalized. normalize
([plot, norm_method, z_sys])Normalize the region if the data are not already normalized. remove_line
(line_tag)Remove absorption line with the given line_tag from the region. set_label
(text)Set descriptive text label for the given region. unpack
()Return the data of the region (wl, flux, error, mask) set_mask -
add_data_to_region
(data_chunk, cutout)[source]¶ Define the spectral data for the fitting region.
Parameters: - data_chunk : dict()
A data_chunk as defined in the data structure of
DataSet.data
.- cutout : bool array
A boolean array defining the subset of the data_chunk which makes up the fitting region.
-
define_mask
(z=None, dataset=None, telluric=True, z_sys=None)[source]¶ Use an interactive window to define the mask for the region.
Parameters: - z : float [default = None]
If a redshift is given, the lines in the region are shown as vertical lines at the given redshift.
- dataset :
VoigtFit.DataSet
[default = None] A dataset with components defined for the lines in the region. If a dataset is passed, the components of the lines in the region are shown.
- telluric : bool [default = True]
Show telluric absorption and sky emission line templates during the masking.
- z_sys : float [default = None]
If a systemic redshift is given, the region is displayed in velocity space relative to the given systemic redshift instead of in wavelength space.
-
generate_label
(active_only=True, ignore_finelines=True)[source]¶ Automatically generate a descriptive label for the region.
-
normalize
(plot=True, norm_method='linear', z_sys=None)[source]¶ Normalize the region if the data are not already normalized. Choose from two methods:
- 1: define left and right continuum regions
- and fit a linear continuum.
- 2: define the continuum as a range of points
- and use spline interpolation to infer the continuum.
If z_sys is not None, show the region in velocity space using instead of wavelength space.
module voigt¶
The module contains functions to evaluate the optical depth, to convert this to observed transmission and to convolve the observed spectrum with the instrumental profile.
-
voigt.
Voigt
(wl, l0, f, N, b, gam, z=0)[source]¶ Calculate the optical depth Voigt profile.
Parameters: - wl : array_like, shape (N)
Wavelength grid in Angstroms at which to evaluate the optical depth.
- l0 : float
Rest frame transition wavelength in Angstroms.
- f : float
Oscillator strength.
- N : float
Column density in units of cm^-2.
- b : float
Velocity width of the Voigt profile in cm/s.
- gam : float
Radiation damping constant, or Einstein constant (A_ul)
- z : float
The redshift of the observed wavelength grid l.
Returns: - tau : array_like, shape (N)
Optical depth array evaluated at the input grid wavelengths l.
-
voigt.
convolve_numba
[source]¶ Define convolution function for a wavelength dependent kernel.
Parameters: - P : array_like, shape (N)
Intrinsic line profile.
- kernel : np.array, shape (N, M)
Each row of the kernel corresponds to the wavelength dependent line-spread function (LSF) evaluated at each pixel of the input profile P. Each LSF must be normalized!
Returns: - P_con : np.array, shape (N)
Resulting profile after performing convolution with kernel.
Notes
This function is decorated by the jit decorator from `numba`_ in order to speed up the calculation.
-
voigt.
evaluate_continuum
(x, pars, reg_num)[source]¶ Evaluate the continuum model using Chebyshev polynomials. All regions are fitted with the same order of polynomials.
Parameters: - x : array_like, shape (N)
Input wavelength grid in Ångstrøm.
- pars : dict(lmfit.Parameters)
An instance of lmfit.Parameters containing the Chebyshev coefficients for each region.
- reg_num : int
The region number, i.e., the index of the region in the list
VoigtFit.DataSet.regions
.
Returns: - cont_model : array_like, shape (N)
The continuum Chebyshev polynomial evaluated at the input wavelengths x.
-
voigt.
evaluate_profile
(x, pars, z_sys, lines, components, kernel, sampling=3, kernel_nsub=1)[source]¶ Evaluate the observed Voigt profile. The calculated optical depth, tau, is converted to observed transmission, f:
\[f = e^{-\tau}\]The observed transmission is subsequently convolved with the instrumental broadening profile assumed to be Gaussian with a full-width at half maximum of res. The resolving power is assumed to be constant in velocity space.
Parameters: - x : array_like, shape (N)
Wavelength array in Ångstrøm on which to evaluate the profile.
- pars : dict(lmfit.Parameters)
An instance of lmfit.Parameters containing the line parameters.
- lines : list(
Line
) List of lines to evaluate. Should be a list of
Line
objects.- components : dict
Dictionary containing component data for the defined ions. See
VoigtFit.DataSet.components
.- kernel : np.array, shape (N, M) or float
The convolution kernel for each wavelength pixel. If an array is given, each row of the array must specify the line-spread function (LSF) at the given wavelength pixel. The LSF must be normalized! If a float is given, the resolution FWHM is given in km/s (c/R). In this case the spectral resolution is assumed to be constant in velocity space.
- sampling : integer [default = 3]
The subsampling factor used for defining the input logarithmically space wavelength grid. The number of pixels in the evaluation will be sampling * N, where N is the number of input pixels. The final profile will be interpolated back onto the original wavelength grid defined by x.
- kernel_nsub : integer
Kernel subsampling factor relative to the data. This is only used if the resolution is given as a LSF file.
Returns: - profile_obs : array_like, shape (N)
Observed line profile convolved with the instrument profile.