____ _ _________________________________________ \/ \ _/\ / \/ \ / 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
- redshiftfloat
Systemic redshift of the absorption system.
- namestr [default = ‘’]
The name of the DataSet, this will be used for saving the dataset to a file structure.
- verbosebool [default = True]
If False, the printed information statements will be suppressed.
- datalist(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
.- linesdict
A dictionary holding pairs of defined (line_tag :
dataset.Line
)- all_lineslist(str)
A list of all the defined line tags for easy look-up.
- moleculesdict
A dictionary holding a list of the defined molecular bands and Jmax for each molecule:
{molecule* : [[band1, Jmax1], [band2, Jmax2], etc...]}
- regionslist(
regions.Region
) A list of the fitting regions.
- cheb_orderint [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_methodstr [default = ‘linear’]
Default normalization method to use for interactive normalization if Chebyshev polynomial fitting should not be used.
- componentsdict
A dictionary of components for each ion defined: (ion : [z, b, logN, options]). See
DataSet.add_component
.- velspanfloat, Tuple(float, float) [default = 400]
The default velocity range to use for the definition of fitting regions.
- ready2fitbool [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_fitlmfit.Parameters [default = None]
Best-fit parameters from lmfit. This attribute will be None until the dataset has been fitted.
- parslmfit.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.- static_variableslmfit.Parameters
Parameter dictionary that holds fit variables other than those related to components of absorption lines and continuum parameters.
- Attributes
- velspan
Methods
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[, airORvac, …])Add spectral data from ASCII file (three or four columns accepted).
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 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.
equivalent_width_limit
(line_tag[, ref, …])Determine the equivalent width limit and corresponding limit on log(N).
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.
lines_of_ion
(ion)Return a list of all line tags 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 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 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
([active_only])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.
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.
add_variable
check_velspan
get_NHI
remove_all_lines
reset_static_variables
set_velspan
- activate_fine_lines(line_tag, levels=None)[source]¶
Activate all lines associated to a given fine-structure complex.
- Parameters
- line_tagstr
The line tag of the ground state transition to activate.
- levelsstr, 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)')
- 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
- ionstr
The ion for which to define a component: e.g., “FeII”, “HI”, “CIa”, etc.
- zfloat
The redshift of the component.
- bfloat
The effective broadening parameter for the component in km/s. This parameter is constrained to be in the interval [0 - 1000] km/s.
- logNfloat
The 10-base logarithm of the column density of the component. The column density is expected in cm^-2.
- var_zbool [default = True]
If False, the redshift of the component will be kept fixed.
- var_bbool [default = True]
If False, the b-parameter of the component will be kept fixed.
- var_Nbool [default = True]
If False, the column density of the component will be kept fixed.
- tie_z, tie_b, tie_Nstr [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, filename='')[source]¶
Add spectral data to the DataSet. This will be used to define fitting regions.
- Parameters
- wlndarray, shape (n)
Input vacuum wavelength array in Ångstrøm.
- fluxndarray, shape (n)
Input flux array, should be same length as wl
- resfloat 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.
- errndarray, 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.
- normalizedbool [default = False]
If the input spectrum is normalized this should be given as True in order to skip normalization steps.
- maskarray, shape (n)
Boolean/int array defining the fitting regions. Only pixels with mask=True/1 will be included in the fit.
- nsubinteger
Kernel subsampling factor relative to the data. This is only used if the resolution is given as a LSF file.
- filenamestring
The filename from which the data originated. Optional but highly recommended. Alternatively, use the method
DataSet.add_spectrum
.
- 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_tagstr
Line tag for the ground state transition, e.g., “CI_1656”
- levelsstr, 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_labelbool [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_tagstr
The line tag for the transition which should be defined. Ex: “FeII_2374”
- velspantuple(float, float) [default = None]
The velocity span around the line center, which will be included in the fit. Can either be given as a single number for a symmetric region, or as a range or (lower, upper). If None is given, use the default self.velspan (default = ±400 km/s).
- activebool [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
- tagslist(str)
A list of line tags for the transitions that should be added.
- velspanfloat [default = None]
The velocity span around the line center, which will be included in the fit. If None is given, use the default
velspan
(±400 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
- moleculestr
The molecular identifier, e.g., ‘CO’, ‘H2’
- bandstr
The vibrational band of the molecule, e.g., for CO: “AX(0-0)” These bands are defined in
molecules
.- Jmaxint [default = 0]
The maximal rotational level to include. All levels up to and including J will be included.
- velspanfloat [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_labelbool [default = False]
If True, the label will be translated to the full quantum mechanical description of the state.
- add_spectrum(filename, res, airORvac='vac', normalized=False, mask_array=None, mask=None, use_mask=True, continuum=None, nsub=1, verbose=False, ext=None)[source]¶
Add spectral data from ASCII file (three or four columns accepted). This is a wrapper of the method add_data.
- Parameters
- filenamestring
Filename of the ASCII spectrum containing at least 3 columns: (wavelength, flux, error). A fourth column may be included containing a spectral mask.
- resfloat 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.
- airORvacstring {‘vac’ or ‘air’}
Defines whether the input wavelengths are in ‘air’ or ‘vacuum’ units. If given as ‘air’ the wavelengths will be converted to ‘vacuum’.
- normalizedbool [default = False]
If the input spectrum is normalized this should be given as True in order to skip normalization steps.
- mask_arrayarray, shape (n)
Boolean/int array defining the fitting regions. Only pixels with mask=True/1 will be included in the fit.
- maskDeprecated, use mask_array
- use_maskbool [default = True]
If False, do not pass the mask array to the dataset. This may be useful if the mask defined in the spectrum is not appropriate for the fitting purposes.
- continuumarray, shape (n)
Continuum spectrum. The input spectrum will be normalized using this continuum spectrum.
- nsubinteger [default = 1]
Kernel subsampling factor relative to the data. This is only used if the resolution is given as a LSF file.
- verbosebool [default = False]
Print status messages.
- extint or string
Index or name of the extension in the HDU List. Only used if the input is a FITS 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_ionstr
The new ion to define.
- from_ionstr
The base ion which will be used as reference.
- logNfloat
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_compint
The reference component to which logN will be scaled.
- tie_zbool [default = True]
If True, the redshifts for all components of the two ions will be tied together.
- tie_bbool [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_tagstr
The line tag of the ground state transition to deactivate.
- levelsstr, 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
.
- equivalent_width_limit(line_tag, ref=None, nofit=False, sigma=3.0, verbose=True, threshold=1.5)[source]¶
Determine the equivalent width limit and corresponding limit on log(N).
- Parameters
- line_tagstr
The line ID of the absorption line for which the limit should be estimated. (Must be a line in the line-list of VoigtFit. Ex.: FeII_1611, TiII_1910, etc.)
- refstr or None
The reference line used to determine the integration range in velocity space, which contains 99.7% of the optical depth. By default, the strongest line of the same ionization state (e.g., I, II, IV, etc.) in the dataset is used. If the reference line has been fitted, the best-fit optical depth profile will be used, otherwise the range is determined by the observed apparent optical depth.
- nofitbool [default = False]
If True, always use the observed apparent optical depth of the reference line (ref, see above).
- sigmafloat [default = 3]
The significance level of the inferred upper limit, by default 3 sigma (99.7%) is used.
- verbosebool [default = True]
Print informative messages from the function?
- thresholdfloat [default = 15]
The continuum noise threshold used to infer the edge of the apparent optical depth of the observed profile. The noise is estimated as the median noise per pixel
- find_line(line_tag)[source]¶
Look up the fitting
Region
for a given line tag.- Parameters
- line_tagstr
The line tag of the line whose region will be returned.
- Returns
- regions_of_linelist 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.
- regions_of_linelist of
- fit(verbose=True, **kwargs)[source]¶
Fit the absorption lines using chi-square minimization.
- Parameters
- verbosebool [default = True]
This will print the fit results to terminal.
- plotbool [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.
- rebinint [default = 1]
Rebin data by a factor rebin before fitting. Passed as part of kwargs.
- samplingint [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
- poptlmfit.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.
- chi2float
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.
- fix_structure(ion=None)[source]¶
Fix the velocity structure, that is, the redshifts and the b-parameters.
- Parameters
- ionstr [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
- ionstr [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
- ionstr
The ion for the lines to get.
- Returns
- lines_for_ionlist(
Line
) List of Lines defined for the given ion.
- lines_for_ionlist(
- 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_tagstr
The line-tag for the line to look up: e.g., “FeII_2374”
- verbosebool [default = False]
If True, print the returned spectral resolution to std out.
- Returns
- resolutionslist 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_tagstr
Line tag for the line belonging to the ion for which components should be defined.
- velocitybool [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_tagstr
Line tag for the
Line
whoseRegion
should be masked.- resetbool [default = True]
If True, clear the mask before defining a new mask.
- maskarray_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.
- telluricbool [default = True]
If True, a telluric absorption template and sky emission template is shown for reference.
- velocitybool [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=True)[source]¶
Normalize or re-normalize a given line
- Parameters
- line_tagstr
Line tag of the line whose fitting region should be normalized.
- norm_methodstr [default = ‘spline’]
Normalization method used for the interactive continuum fit. Should be on of: [“spline”, “linear”]
- velocitybool [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=10, 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
- normbool [default = True]
Opens an interactive window to let the user normalize each region using the defined
norm_method
.- maskbool [default = True]
Opens an interactive window to let the user mask each fitting region.
- verbosebool [default = True]
If this is set, the code will print small info statements during the run.
- active_onlybool [default = False]
If True, only define masks for active lines.
- force_cleanbool [default = False]
If this is True, components for inactive elements will be removed.
- velocitybool [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_linesbool [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_lowerfloat [default = 0.]
Lower limit on oscillator strengths for transitions when verifying all transitions for defined ions.
- f_upperfloat [default = 100.]
Upper limit on oscillator strengths for transitions when verifying all transitions for defined ions.
- l_lowerfloat [default = 0.]
Lower limit on rest-frame wavelength for transitions when verifying all transitions for defined ions.
- l_upperfloat [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
- velocitybool [default = True]
If True, show the relative velocities of each component instead of redshifts.
- elementslist(str) [default = ‘all’]
A list of elements for which to show parameters.
- systemicfloat [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_tagstr
The line tag of the ground state transition to remove.
- levelsstr, 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_tagstr
Line tag of the transition that should be removed.
- reset_all_regions(active_only=True)[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
- ionstr [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.
- save_cont_parameters_to_file(filename)[source]¶
Save the best-fit continuum parameters to ASCII table output.
- Parameters
- filenamestr [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
- filenamestr [default = None]
Filename for the fitting regions. If None, the
name
attribute will be used.- individualbool [default = False]
Save the fitting regions to individual files. By default all regions are concatenated into one file.
- pathstr [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
- filenamestr
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
- ionsstr or list(str)
List of ions or a single ion for which to calculate the total column density.
- componentslist(int)
List of integers corresponding to the indeces of the components to sum over.
- Returns
- total_logNdict()
Dictionary containing the log of total column density for each ion.
- total_logN_errdict()
Dictionary containing the error on the log of total column density for each ion.
class Component¶
- class VoigtFit.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 VoigtFit.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
- tagstr
The line tag for the line, e.g., ‘FeII_2374’
- ionstr
The ion for the line; The ion for ‘FeII_2374’ is ‘FeII’.
- elementstr
Equal to ion for backwards compatibility.
- l0float
Rest-frame resonant wavelength of the transition. Unit: Angstrom.
- ffloat
The oscillator strength for the transition.
- gamfloat
The radiation damping constant or Einstein coefficient.
- massfloat
The atomic mass in atomic mass units.
- activebool [default = True]
The state of the line in the dataset. Only active lines will be included in the fit.
Methods
Return the principal atomic constants for the transition: l0, f, and gam.
Set the line active; include the line in the fit.
Set the line inactive; exclude the line in the fit.
class Region¶
- class VoigtFit.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
- velspanTuple(float, float)
The velocity ranges used for the fitting region. Given as a range by a tuple of (lower, upper).
- lineslist(
dataset.Line
) A list of Lines defined in the region.
- labelstr
A LaTeX label describing the lines in the region for plotting purposes.
- resfloat
Spectral resolution of the region in km/s.
- wlarray_like, shape (N)
Data array of wavelengths in Ångstrøm.
- fluxarray_like, shape (N)
Data array of fluxes (normalized if
normalized
is True).- errarray_like, shape (N)
Array of uncertainties for each flux element.
- normalizedbool
True if the data in the region are normlized.
- maskarray_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_maskbool
Internal parameter for
VoigtFit.DataSet.prepare_dataset()
. If True, an interactive masking process will be initiated in the preparation stage.- cont_errfloat
An estimate of the uncertainty in the continuum fit.
- specIDstr
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 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.
get_velocity
(z_sys[, line])Return the relative velocities of the region, with respect to systemtic redshift of the given line.
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.
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)
evaluate_region
set_mask
- add_data_to_region(data_chunk, cutout)[source]¶
Define the spectral data for the fitting region.
- Parameters
- data_chunkdict()
A data_chunk as defined in the data structure of
DataSet.data
.- cutoutbool 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
- zfloat [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.
- telluricbool [default = True]
Show telluric absorption and sky emission line templates during the masking.
- z_sysfloat [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.
- get_velocity(z_sys, line=None)[source]¶
Return the relative velocities of the region, with respect to systemtic redshift of the given line.
- 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.
- VoigtFit.funcs.voigt.H(a, x)[source]¶
Voigt Profile Approximation from T. Tepper-Garcia (2006, 2007).
- VoigtFit.funcs.voigt.Voigt(wl, l0, f, N, b, gam, z=0)[source]¶
Calculate the optical depth Voigt profile.
- Parameters
- wlarray_like, shape (N)
Wavelength grid in Angstroms at which to evaluate the optical depth.
- l0float
Rest frame transition wavelength in Angstroms.
- ffloat
Oscillator strength.
- Nfloat
Column density in units of cm^-2.
- bfloat
Velocity width of the Voigt profile in cm/s.
- gamfloat
Radiation damping constant, or Einstein constant (A_ul)
- zfloat
The redshift of the observed wavelength grid l.
- Returns
- tauarray_like, shape (N)
Optical depth array evaluated at the input grid wavelengths l.
- VoigtFit.funcs.voigt.convolve_numba(P, kernel)[source]¶
Define convolution function for a wavelength dependent kernel.
- Parameters
- Parray_like, shape (N)
Intrinsic line profile.
- kernelnp.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_connp.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.
- VoigtFit.funcs.voigt.convolve_profile(profile, width)[source]¶
Convolve with Gaussian kernel using constant width in pixels!
- Parameters
- profilearray, shape (N)
Input profile to convolve with a Gaussian kernel.
- widthfloat
Kernel width (or sigma) of the Gaussian profile in units of pixels. Note – this should not be the FWHM, as is usually used to denote the resolution element: width = FWHM / 2.35482
- Returns
- profile_obsarray, shape(N)
The convolved version of profile.
- VoigtFit.funcs.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
- xarray_like, shape (N)
Input wavelength grid in Ångstrøm.
- parsdict(lmfit.Parameters)
An instance of lmfit.Parameters containing the Chebyshev coefficients for each region.
- reg_numint
The region number, i.e., the index of the region in the list
VoigtFit.DataSet.regions
.
- Returns
- cont_modelarray_like, shape (N)
The continuum Chebyshev polynomial evaluated at the input wavelengths x.
- VoigtFit.funcs.voigt.evaluate_optical_depth(profile_wl, pars, lines, z_sys=None)[source]¶
Evaluate optical depth based on the component parameters in pars.
- Parameters
- profile_wlarray_like, shape (N)
Wavelength array in Ångstrøm on which to evaluate the profile.
- parsdict(lmfit.Parameters)
An instance of lmfit.Parameters containing the line parameters.
- lineslist(
Line
) List of lines to evaluate. Should be a list of
Line
objects.- z_sysfloat or None
The systemic redshift, used to determine an effective wavelength range within which to evaluate the profile. This is handy when fitting very large regions, such as HI and metal lines together.
- Returns
- tauarray_like, shape (N)
The resulting optical depth of all lines in the wavelength region.
- VoigtFit.funcs.voigt.evaluate_profile(x, pars, lines, kernel, z_sys=None, 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
- xarray_like, shape (N)
Wavelength array in Ångstrøm on which to evaluate the profile.
- parsdict(lmfit.Parameters)
An instance of lmfit.Parameters containing the line parameters.
- lineslist(
Line
) List of lines to evaluate. Should be a list of
Line
objects.- kernelnp.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.
- z_sysfloat or None
The systemic redshift, used to determine an effective wavelength range within which to evaluate the profile. This is handy when fitting very large regions, such as HI and metal lines together.
- samplinginteger [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_nsubinteger
Kernel subsampling factor relative to the data. This is only used if the resolution is given as a LSF file.
- Returns
- profile_obsarray_like, shape (N)
Observed line profile convolved with the instrument profile.
- VoigtFit.funcs.voigt.fwhm_to_pixels(wl, res_wl)[source]¶
Convert spectral resolution in wavelength to pixels R = lambda / dlambda, where dlambda is the FWHM
- Parameters
- wlarray, shape (N)
Input array of wavelength to evaluate
- res_wlfloat
The input spectral resolution element (FWHM) in Ångström.
- Returns
- wfloat
Kernel width converted to pixels (to be used with
convolve_profile
).
- VoigtFit.funcs.voigt.resvel_to_pixels(wl, res)[source]¶
Convert spectral resolution in km/s to pixels
- Parameters
- wlarray, shape (N)
Input array of wavelength to evaluate
- resfloat
The input spectral velocity resolution in km/s.
- Returns
- wfloat
Kernel width converted to velocity pixels (km/s). (to be used with convolution of logarithmically-spaced data/models).
module output¶
- VoigtFit.io.output.load_lsf(lsf_fname, wl, nsub=1)[source]¶
Load a Line-Spread Function table following format from HST: First line gives wavelength in Angstrom and the column below each given wavelength defines the kernel in pixel space:
wl1 wl2 wl3 ... wlN lsf11 lsf21 lsf31 ... lsfN1 lsf12 lsf22 lsf32 ... lsfN2 : : lsf1M lsf2M lsf3M ... lsfNM
- Parameters
- lsf_fnamestring
The filename containing the LSF data.
- wlarray like, shape (N)
The wavelength grid onto which the LSF will be evaluated
- nsubint [default = 1]
Kernel subsampling factor relative to the data. This is only used if the resolution is given as a LSF file.
- Returns
- kernelnp.array, shape(N, M)
A grid of interpolated LSF evaluated at each given input wavelength of the array wl of shape N, where M is the number of pixels in the LSF.
Notes
The output kernel is transposed with respect to the input format for ease of computation in the convolution since indexing is faster along rows than columns.
- VoigtFit.io.output.merge_two_dicts(default, x)[source]¶
Merge the keys of dictionary x into dictionary default.
- VoigtFit.io.output.plot_H2(dataset, n_rows=None, xmin=None, xmax=None, ymin=- 0.1, ymax=2.5, short_labels=False, rebin=1, smooth=0)[source]¶
Generate plot for H2 absorption lines.
- Parameters
- dataset
VoigtFit.DataSet
An instance of the class
VoigtFit.DataSet
containing the H2 lines to plot.- n_rowsint [default = None]
Number of rows to show in figure. If None, the number will be determined automatically.
- xminfloat
The lower x-limit in Å. If nothing is given, the extent of the fit region is used.
- xmaxfloat
The upper x-limit in Å. If nothing is given, the extent of the fit region is used.
- yminfloat [default = -0.1]
The lower y-limit in normalized flux units.
- ymaxfloat [default = 2.5]
The upper y-limit in normalized flux units.
- rebinint [defualt = 1]
Rebinning factor for the spectrum, default is no binning.
- smoothfloat [default = 0]
Width of Gaussian kernel for smoothing.
- dataset
- VoigtFit.io.output.plot_all_lines(dataset, plot_fit=True, rebin=1, fontsize=12, xmin=None, xmax=None, ymin=None, ymax=None, max_rows=4, filename=None, subsample_profile=1, npad=50, residuals=True, norm_resid=False, legend=True, loc='left', show=True, default_props={}, element_props={}, highlight_props=None, label_all_ions=False, xunit='vel', line_props=None, hl_line_props=None)[source]¶
Plot all active absorption lines. This function is a wrapper of
plot_single_line()
. For a complete description of input parameters, see the documentation forplot_single_line()
.- Parameters
- dataset
VoigtFit.DataSet
Instance of the
VoigtFit.DataSet
class containing the line regions to plot.- max_rowsint [default = 4]
The maximum number of rows of figures. Each row consists of two figure panels.
- filenamestr
If a filename is given, the figures are saved to a pdf file.
- dataset
- VoigtFit.io.output.plot_excitation(dataset, molecule)[source]¶
Plot the excitation diagram for a given molecule
- VoigtFit.io.output.plot_residual(dataset, line_tag, index=0, rebin=1, xmin=None, xmax=None, axis=None)[source]¶
Plot residuals for the best-fit to a given absorption line.
- Parameters
- datasetclass DataSet
An instance of DataSet class containing the line region to plot.
- line_tagstr
The line tag of the line to show, e.g., ‘FeII_2374’
- indexint [default = 0]
The line index. When fitting the same line in multiple spectra this indexed points to the index of the given region to be plotted.
- rebin: int
Integer factor for rebinning the spectral data.
- xminfloat
The lower x-limit in relative velocity (km/s). If nothing is given, the extent of the region is used.
- xmaxfloat
The upper x-limit in relative velocity (km/s). If nothing is given, the extent of the region is used.
- axismatplotlib.axes.Axes
The plotting axis of matplotlib. If None is given, a new figure and axis will be created.
- VoigtFit.io.output.plot_single_line(dataset, line_tag, index=0, plot_fit=False, loc='left', rebin=1, nolabels=False, axis=None, fontsize=12, subsample_profile=1, xmin=None, xmax=None, ymin=None, ymax=None, show=True, npad=50, legend=True, residuals=False, norm_resid=False, default_props={}, element_props={}, highlight_props=None, label_all_ions=False, xunit='velocity', line_props=None, hl_line_props=None, sort_f=True)[source]¶
Plot a single absorption line.
- Parameters
- dataset
VoigtFit.DataSet
Instance of the
VoigtFit.DataSet
class containing the lines- line_tagstr
The line tag of the line to show, e.g., ‘FeII_2374’
- indexint [default = 0]
The line index. When fitting the same line in multiple spectra this indexed points to the index of the given region to be plotted.
- plot_fitbool [default = False]
If True, the best-fit profile will be shown
- locstr [default = ‘left’]
Places the line tag (right or left).
- rebinint [default = 1]
Rebinning factor for the spectrum
- nolabelsbool [default = False]
If True, show the axis x- and y-labels.
- axismatplotlib.axes.Axes
The plotting axis of matplotlib. If None is given, a new figure and axis will be created.
- fontsizeint [default = 12]
The fontsize of text labels.
- xminfloat
The lower x-limit in relative velocity (km/s). If nothing is given, the extent of the region is used.
- xmaxfloat
The upper x-limit in relative velocity (km/s). If nothing is given, the extent of the region is used.
- yminfloat [default = None]
The lower y-limit in normalized flux units. Default is determined from the data.
- ymaxfloat [default = None]
The upper y-limit in normalized flux units. Default is determined from the data.
- showbool [default = True]
Show the figure.
- subsample_profileint [default = 1]
Subsampling factor to calculate the profile on a finer grid than the data sampling. By default the profile is evaluated on the same grid as the data.
- npadint [default = 50]
Padding added to the synthetic profile before convolution. This removes end artefacts from the FFT routine.
- residualsbool [default = False]
Add a panel above the absorption line to show the residuals of the fit.
- norm_residbool [default = False]
Show normalized residuals.
- legendbool [default = True]
Show line labels as axis legend.
- default_propsdict
Dictionary of transition tick marker properties. The dictionary is passed to matplotlib.axes.Axes.axvline to control color, linewidth, linestyle, etc.. Two additional keywords can be defined: The keyword text is a string that will be printed above or below each tick mark for each element. The keyword loc controls the placement of the tick mark text for the transistions, and must be one either ‘above’ or ‘below’.
- element_propsdict
Dictionary of properties for individual elements. Each element defines a dictionary with individual properties following the format for default_props.
- Ex:
element_props={'SiII': {'color': 'red', 'lw': 1.5},
'FeII': {'ls': '--', 'alpha': 0.2}}
This will set the color and linewidth of the tick marks of SiII transitions and the linestyle and alpha-parameter of FeII transitions.
- Ex:
- highlight_propsdict/list [default = None]
A dictionary of ions (e.g., “FeII”, “CIa”, etc.) used to calculate a separate profile for this subset of ions. Each ion as a keyword must specify a dictionary which can change individual properties for the given ion. Similar to element_props. If an empty dictionary is given, the default parameters will be used. Alternatively, a list of ions can be given to use default properties for all ions.
Ex:
highlight_props={'SiII':{}, 'FeII':{'color': 'blue'}}
This will highlight SiII transitions with default highlight properties, and FeII transistions with a user specified color.
Ex:
highlight_props=['SiII', 'FeII']
This will highlight SiII and FeII transitions using default highlight properties.
- label_all_ionsbool [default = False]
Show labels for all ions defined. The labels will appear above the component tick marks.
- xunitstring [default = ‘velocity’]
The unit of the x-axis, must be either ‘velocity’ or ‘wavelength’. Shortenings are acceptable too, e.g., ‘vel’/’v’ or ‘wave’/’wl’.
- line_propsdict [default = None]
A dictionary of keywords to change the default line properties of the best-fit profile; e.g., ‘color’, ‘lw’, ‘linestyle’. All keywords will be passed to the plot function of matplotlib.
- hl_line_propsdict [default = None]
A dictionary of keywords to change the default line properties of the best-fit profile for highlighted ions. All keywords will be passed to the plot function of matplotlib.
- sort_fbool [default = True]
If True, calculate velocities with respect to the line with the largest oscillator strength. Otherwise, use the given line_tag as reference.
- dataset
- VoigtFit.io.output.print_T_model_pars(dataset, thermal_model, filename=None)[source]¶
Print the turbulence and temperature parameters for physical model.
- VoigtFit.io.output.print_cont_parameters(dataset)[source]¶
Print the Chebyshev coefficients of the continuum normalization.
- VoigtFit.io.output.print_metallicity(dataset, params, logNHI, err=0.1)[source]¶
Print the metallicity derived from different species. This will add the column densities for all components of a given ion.
- Parameters
- dataset
VoigtFit.DataSet
An instance of the
VoigtFit.DataSet
class containing the definition of data and absorption lines.- paramslmfit.Parameters
Output parameter dictionary, e.g.,
VoigtFit.DataSet.best_fit
. See lmfit for details.- logNHIfloat
Column density of neutral hydrogen.
- errfloat [default = 0.1]
Uncertainty (1-sigma) on logNHI.
- .. _lmfit: https://lmfit.github.io/lmfit-py/index.html
- dataset
- VoigtFit.io.output.print_results(dataset, params, elements='all', velocity=True, systemic=0)[source]¶
Print the parameters of the best-fit.
- Parameters
- dataset
VoigtFit.DataSet
An instance of the
VoigtFit.DataSet
class containing the line region to plot.- paramslmfit.Parameters
Output parameter dictionary, e.g.,
dataset.Dataset.best_fit
. See lmfit for details.- elementslist(str) [default = ‘all’]
A list of ions for which to print the best-fit parameters. By default all ions are shown.
- velocitybool [default = True]
Show the components in relative velocity or redshift.
- systemicfloat [default = 0]
Systemic redshift used to calculate the relative velocities. By default the systemic redshift defined in the dataset is used.
- .. _lmfit: https://lmfit.github.io/lmfit-py/index.html
- dataset
- VoigtFit.io.output.print_total(dataset, verbose=True)[source]¶
Print the total column densities of all species. This will sum all the components of each ion. The uncertainty on the total column density is calculated using random resampling within the errors of each component.
- VoigtFit.io.output.rebin_bool_array(x, n)[source]¶
Rebin input boolean array x by an integer factor of n.
- VoigtFit.io.output.rebin_spectrum(wl, spec, err, n, method='mean')[source]¶
Rebin input spectrum by a factor of n. Method is either mean or median, default is mean.
- Parameters
- wlarray_like, shape (N)
Input wavelength array.
- specarray_like, shape (N)
Input flux array.
- errarray_like, shape (N)
Input error array.
- nint
Integer rebinning factor.
- methodstr [default = ‘mean’]
Rebin method, either ‘mean’ or ‘median’.
- Returns
- wl_rarray_like, shape (M)
Rebinned wavelength array, the new shape will be N/n.
- spec_rarray_like, shape (M)
Rebinned flux array, the new shape will be N/n.
- err_rarray_like, shape (M)
Rebinned error array, the new shape will be N/n.
- VoigtFit.io.output.save_EW(EW_list, filename, path='')[source]¶
Save the equivalent width data for each line in the input list.
- Parameters
- EW_listlist[namedtuple (EquivalentWidth)]
A list of namedtuple (EquivalentWidth), see implementation in VoigtFit.dataset.py
- filenamestr
Filename to be used. If the filename already exists, it will be overwritten.
- pathstr [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 directory does not exist, it will be created. The final filename will be:
path/ + filename
- VoigtFit.io.output.save_cont_parameters_to_file(dataset, filename, path='')[source]¶
Save Chebyshev coefficients to file.
- VoigtFit.io.output.save_fit_regions(dataset, filename, individual=False, path='')[source]¶
Save fit regions to ASCII file.
- Parameters
- filenamestr
Filename to be used. If the filename already exists, it will be overwritten. A .reg file extension will automatically be append if not present already.
- individualbool [default = False]
If True, save each fitting region to a separate file. The individual filenames will be the basename given as filename with _regN appended, where N is an integer referring to the region number.
- pathstr [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 directory does not exist, it will be created. The final filename will be:
path/ + filename [+ _regN] + .reg
- VoigtFit.io.output.save_individual_components(dataset, filename, path='')[source]¶
Save best-fit profile to ASCII file for each individual component of each ion. By default the profiles are calcuated on a logarithmically binned wavelength array over the entire range of the input data. If a sampling factor other than one is given, the wavelength grid is subsampled by this factor.
- Parameters
- datasetdataset.DataSet
Dataset with the fitted data and parameters
- filenamestr
Filename to be used. If the filename already exists, it will be overwritten.
- pathstr [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 directory does not exist, it will be created. The final filename will be:
path/ + filename
- VoigtFit.io.output.save_parameters_to_file(dataset, filename, path='')[source]¶
Save best-fit parameters to file.
- VoigtFit.io.output.show_H2_bands(ax, z, bands, Jmax, color='blue', short_labels=False)[source]¶
Add molecular H2 band identifications to a given matplotlib axis.
- Parameters
- axmatplotlib.axes.Axes
The axis instance to decorate with H2 line identifications.
- zfloat
The redshift of the H2 system.
- bandslist(str)
A list of molecular bands of H2. Ex: [‘BX(0-0)’, ‘BX(1-0)’]
- Jmaxint or list(int)
The highest rotational J-level to include in the identification. A list of Jmax for each band can be passed as well. Lines are currently only defined up to J=7.
- VoigtFit.io.output.show_components(dataset, ion=None)[source]¶
Show the defined components for a given ion. By default, all ions are shown.
- VoigtFit.io.output.sum_components(dataset, ion, components)[source]¶
Calculate the total abundance for given components of the given ion.
- Parameters
- dataset
VoigtFit.DataSet
An instance of the
VoigtFit.DataSet
class containing the definition of data and absorption lines.- ionstr
Ion for which to calculate the summed abundance.
- componentslist(int)
List of indeces of the components to sum over.
- dataset
- Returns
- total_logNfloat
The 10-base log of total column density.
- total_logN_errfloat
The error on the 10-base log of total column density.