3. Spectral line fitting and radial velocity measurements¶
The main purpose of MUSEreduce is the spectral line fitting in order to create template spectra from observations, which are used to measure the radial velocities of stars and gas, especially in the absence of working spectral template libraries, e.g., for pre-main-sequence stars in the optical.
The radial_velocities.RV_spectrum
class and the spectral fitting module line_fitter
are the heart of this package using pyspeckit as basic spectral fitting routine and ppxf for the spectral cross-correlation. Together with a Monte Carlo bootstrap technique, it is possible to measure radial velocities with an accuracy of \(1-3\,\rm{km}\,\rm{s}^{-1}\) in the absence of any template library. For a detailed description on how this technique works, we refer to Zeidler et al. 2019.
3.1. Spectral line fitter¶
This module is the work horse for fitting spectral lines from an input catalog to a given 1D spectrum. The line fitting itself is performed using the module pyspeckit.
pyspeckit uses an iterative method to fit the spectral lines and the continuum simultaneously. line_fitter
automatically adjusts the input parameters to pyspeckit in the manner to optimize the fitting (continuum order and wavelegth limits)
The line_fitter
determines whether the fit fullfils the set parameters and was succesfull. The line_fitter
also allows to fit blends by fixing a maximum ratio between the primary line (the one the user is interested in) and the secondary line (blend). This ensures that the blend does not become the dominant line.
In order to accomodate hyper-velocity stars an option is given that the wavelength limits are automatically adjust for each iteration \(n\) based on the solution \(n-1\) of the primary line taking into account \(\Delta \lambda / \lambda\).
If the spectra with high radial velocities, or if the object is red shifted an initial guess of the radial velocity needs to be provided. This may be done by the Kwargs rv_sys of the radial_velocities.RV_spectrum
. This is important if the shift is much larger than the wavelength limits, where the automatic adjustment fails. Additionally, the more accurate this rv_sys is provided the faster the fit converges.
-
line_fitter.
line_fitter
(self, linecat, line_idx, niter, input_resid_level, max_contorder, max_ladjust, adjust_preference, input_continuum_deviation, llimits, max_exclusion_level, blends, autoadjust, fwhm_block)[source]¶ This is module fits the spectral lines using pyspeckit. It automatically determines the goodness of fit and decides the best solution for the line and continuum fit. It handles blended lines in a way a maximum lines ratio can not be exceeded to not make the weaker blend the dominant lines. Additionally, the limits in wavelength range in iteration \(n\) can be adjusted automatically based on iteration \(n\) to accomodate for larger wavelength shifts.
In this way it is possible to fit spectral lines using a line catalog as only input.
- Args:
- linecat :
numpy.array()
- Array with the spectral lines and their wavelengths
- line_idx :
str
- Name of the primary line
- niter :
int
- Number of iterations
- input_resid_level :
float
- The maximum MAD for the fit residuals for a succesfull fit
- max_contorder :
int
- The maximum polynominal order of the continuum to have
- max_ladjust :
str
- The maximum number of wavelength range adjustments in steps of 5 Angstrom
- adjust_preference :
str
contorder: continuum order is adjusted first
wavelength: wavelength range is adjusted first
- input_continuum_deviation :
float
- by how much the continuum is allowed to deviate from a running median estimate. This is set to prevent lines mimicking a continuum
- llimits :
list
- the limits for the wavelength fit as set in
ppxf
- max_exclusion_level :
float
- The exclusion level for lines to be excluded from the next baseline
estimate as set in
pyspeckit
- blends :
ascii
-file orNone
- A file with primary lines that contain blends to provide a maximum amplitude ratio of the primary and the blend to prevent that the blend becomes the dominant line in the fit.
- autoadjust :
bool
True
: the wavelength limitsllimit
will be adjusted to the fit of the previous iteration. All other wavelength range are adjusted accordingly taking into account the proper velocity corrected shift \(\Delta \lambda/\lambda\). This is especially important to detect hyper-velocity stars.False
: no adjustment to the limits done- fwhm_block :
bool:obj:
True
: The minimum fwhm of the voigt profiles of the fitted lines is the instrument’s dispersionFalse
: The minimum fwhm of the voigt profiles of the fitted lines is zero
- linecat :
- Returns:
- line_idx :
str
- Name of the primary line
- temp_l :
float
- fitted wavelength of the primary line
- temp_a :
float
- fitted amplitude of the primary line
- temp_sl :
float
- fitted Lorentzian gamma of the primary line
- temp_sg :
float
- fitted Gaussian sigma of the primary line
- spec_select_idx :
numpy.array()
- indices of the used part of the spectrum
- template_f :
numpy.array()
- template spectrum shifted to the rest-frame
- continuum :
numpy.array()
- the fitted continuum
- lstart :
float
- first used wavelength bin (might differ from input if it was adjusted during fitting)
- lend :
float
- last used wavelength bin (might differ from input if it was adjusted during fitting)
- contorder :
int
- The order of the polynominal used for the continuum
- fit_f : array
- the fitted spectrum
- significance :
float
- the line strength over the continuum
- fit_failed :
bool
True
: if the fit failed for some reason. This line will be excluded from further analyses- fit_f_highres :
float
- the fitted spectrum and the oversampling resolution
- spec_select_idx_highres :
numpy.array()
- indices of the used part of the over sampled spectrum
- template_f_highres :
numpy.array()
- over sampled template spectrum shifted to the rest-frame
- continuum_highres :
numpy.array()
- the fitted over sampled continuum
- line_idx :
Warning
The automatic handling of absorption and emission lines implemented in vers. 0.1.2 has not been tested yet.
3.2. Radial Velocities¶
This is the main class for measuring the radial velocities. A basic example is descibed in Examples. A detailed description can be found in Zeidler et al. 2019.
Throughout the document the primary line is the spectral line of interest for which the line fitting will be executed. Secondary lines are spectral lines in blends.
A detailed demonstration on how to use the radial velocity fitter is provided in Examples including template files.
3.2.1. The radial velocity fitter¶
-
class
radial_velocities.
RV_spectrum
(spec_id, spec_f, spec_err, spec_lambda, loglevel='INFO', templatebins=100000, specbinsize=1.25, dispersion=2.4, linetype='absorption', rv_sys=0.0)[source]¶ This class contains the 1D spectrum, as well as all fitting parameters and output values including the radial velocity catalog.
- Args:
- spec_id :
str
- unique identifier of the spectrum
- spec_f :
numpy.array()
- flux array of the spectrum in erg/s/cm:math:^2/Angstrom
- spec_err :
numpy.array()
- flux uncertainty array of the spectrum
- spec_lambda :
numpy.array()
- wavelength_array of the spectrum in Angstrom
- spec_id :
- Kwargs:
- loglevel :
str
(optional, default:INFO
) DEBUG
: all functions run on a single core to obtain proper output in the correct order for proper debugging.- templatebins :
float
(optional, default: 100000) - Number of wavelength bins for the oversampled spectrum used to fit the spectral lines.
- specbinsize :
float
(optional, default: 1.25) - The spectral bin size. The default is set to fit the MUSE dataset
- dispersion :
float
(optional, default: 2.4) - The dispersion of the spectrograph. The default is set to the nominal MUSE instrument dispersion
- linetype:
str
(optional, default: absorption) absorption: All spectral lines are absorption lines
emission: All spectral lines are emission lines
both : spectral lines can be absorption or emission lines. This mode has not been tested yet !!!
- rv_sys :
float
(optional, default: 0) - systematic RV shift in km/s Should be provided for large velocity offsets or redshifted objects
- loglevel :
-
catalog
(initcat=None, save=False, load=None, printcat=False)[source]¶ This module handles the catalog that holds the fit results of the primary lines.
- Kwargs:
- initcat :
ascii
(optional, default:None
) ascii
file containing the primary lines format: name, lambda, start, end- save :
bool
(optional, default:False
) True
: The catalog is written to a file.- load :
str
(optional, default:None
) - Loads a catalog from catalog file.
- printcat :
bool
(optional, default:False
) True
: The catalog is printed to the terminal
- initcat :
-
clean
()[source]¶ This modules resets all of the output. This must be executed before repeating the spectral fit without re-initiating the class.
-
line_fitting
(input_cat, line_idxs, niter=5, n_CPU=-1, resid_level=None, max_contorder=2, max_ladjust=4, adjust_preference='contorder', input_continuum_deviation=0.05, llimits=[-2.0, 2.0], max_exclusion_level=0.3, blends=None, autoadjust=False, fwhm_block=False)[source]¶ Initializing the line fitting by calling
line_fitter
- Args:
- input_cat:
numpy.array()
- input spectral line catalog with wavelengths in Angstrom
- line_idxs:
numpy.array()
str
numpy.array()
: of the line names for which RV fits should be performed, must be identical to “name” in init_cat
- input_cat:
- Kwargs:
- niter :
int
(optional, default: 5) - The maximum number of iteration if convergence is not reached
- n_CPU :
float
(optional, default: -1) - Setting the number of CPUs used for the parallelization. If set to -1 all available system resources are used. Maximum number of CPUs is the number of spectral lines the fit is performed on.
- resid_level:
float
orNone
(optional, default: None) - The maximum MAD for the fit residuals for a succesfull fit
- max_contorder :
int
(optional, default: 2) - The maximum polynominal order of the continuum
- max_ladjust :
int
(optional, default: 4) - Ahe maximum number of wavelength range adjustments in steps of 5 Angstrom
- adjust_preference:
str
(optional, default: contorder) contorder
: continuum order is adjusted firstwavelength
: wavelength range is adjusted first- input_continuum_deviation
float
(optional, default: 0.05) - Fraction by how much the continuum is allowed to deviate from a running median estimate. This is set to prevent lines mimicking a continuum
- llimits:
list
(optional, default: [-2., 2.]) - the limits for the wavelength fit as set in ppxf
- max_exclusion_level
float
(optional, default: 0.3) - The exclusion level for lines to be excluded from the next baseline estimate as set in pyspeckit
- blends:
ascii
orNone
(optional, default: None) - A file with primary lines that contain blends to provide a maximum amplitude ratio of the primary and the blend to prevent that the blend becomes the dominant line in the fit.
- autoadjust
bool
(optional, default:False
) True
: the wavelength limitsllimit
will be adjusted to the fit of the previous iteration. All other wavelength range are adjusted accordingly taking into account the proper velocity corrected shift \(\Delta \lambda/\lambda\). This is especially important to detect hyper-velocity stars.- fwhm_block
bool
(optional, default:False
) True
: The minimum fwhm of the voigt profiles of the fitted lines ais the instrument’s dispersion. This prevents unphysical lines.False
: The minimum fwhm of the voigt profiles of the fitted lines is zero.
- niter :
-
plot
(oversampled=False)[source]¶ Plots the spectrum and the regions of the primary lines to a file including the spectral fit and the template.
-
rv_fit
(guesses, niter=10000, line_sigma=3, n_CPU=-1, line_significants=5, RV_guess_var=0.0)[source]¶ The module runs the radial velocity fit using ppxf and the Monte Carlo radial velocity determination.
- Args:
- guesses :
numpy.array()
- The initial guesses for the the radial velocity fit guesses in the form [RV,sepctral_dispersion]
- guesses :
- Kwargs:
- niter :
int
(optional, default: 10000) - number of iterations to bootstrap the spectrum
- line_sigma:
int
(optional, default: 3): - sigma for the RV clipping for the individual lines
- n_CPU :
float
(optional, default: -1) - Setting the number of CPUs used for the parallelization. If set to -1 all available system resources are used. Maximum number of CPUs is the number of spectral lines the fit is performed to.
- line_significants:
int
(optional, default: 5) - The sigma-level for the spectral line to be above the continuum in order to be considered valid
- RV_guess_var :
float
(optional, default: 0) - The maximum variation the RV guess will be varied using a uniform distribution.
- niter :
3.2.2. Monte Carlo radial velocity determination¶
This is the module that performs the Monte Carlo boot strapping to measure radial velocities. It is automatically called by the radial_velocities.RV_spectrum.rv_fit
attribute of radial_velocities.RV_spectrum
class and there is no need to repeat this step manually. Nevertheless, we show this part of the code since it may be useful for other applications to measure RVs.
-
ppxf_MC.
ppxf_MC
(log_template_f, log_spec_f, log_spec_err, velscale, guesses, nrand=100, degree=4, goodpixels=None, moments=4, vsyst=0, sigma=5, spec_id=None, RV_guess_var=0.0, n_CPU=-1)[source]¶ This module runs the Monte Carlo ppxf runs, which is needed for the RV measurements. Most of the input parameters are similar to the standard ppxf parameters (see Cappellari and Emsellem 2004) for a more detailed explanation).
- Args:
- log_template_f :
numpy.array()
- The logarithmically binned template spectrum
- log_spec_f :
numpy.array()
- The logarithmically binned source spectrum
- log_spec_err :
numpy.array()
- The logarithmically source spectrum uncertainties
- velscale:
float
- The velocity scale of the source spectrum
- guesses :
numpy.array()
- The initial guesses for the the radial velocity fit guesses in the form [RV,sepctral_dispersion]
- log_template_f :
- Kwargs:
- nrand :
int
(optional, default: 5) - The maximum number of iteration if convergence is not reached
- degree :
int
(optional, default: 4) - The degree of the additive polynomial to fit offsets in the continuum. A low order polynominal may be needed to get better results
- goodpixels :
numpy.array()
(optional, default: 5) - A
numpy.array()
: of pixels that are used by ppxf to fit the template to the source spectrum - moments :
int
(optional, default: 4) - The moments to be fit (v, sigma, h3, h4)
- vsyst :
float
(optional, default: 0.) - A systematic velocity. This may be needed if the system move at high velocities compared to the rest frame. If the guess of vsyst is good, the fit runs faster and more stable.
- sigma :
int
(optional, default: 5) - The sigma used to clip outliers in the histogram determination.
- spec_id :
str
(optional, default: None) - An ID number of the source spectrum. This becomes handy when fitting many individual sources because the output files will be named with the ID.
- RV_guess_var :
float
(optional, default: 0) - The maximum variation the RV guess will be varied using a uniform distribution.
- n_CPU :
int
(optional, default: -1) - The number of cores used for the Monte Carlo velocity fitting. If n_CPU=-1 than all available cores will be used.
- nrand :
3.3. History¶
New in version 0.1.0: working radial velocity fitting package
New in version 0.1.0: no RV calculation for lines that are below the significance level
New in version 0.1.0: instrument dispersion added as keyword
New in version 0.1.1: moved to pep-8
New in version 0.1.2: now handles absorption and emission lines. Not tested yet, though
New in version 0.1.3: During the boot strap procedure, the initial guesses of the velocities are varied between certain limits, to ensure more stability to the
fit in the case of an “unlucky choice of initial parameters. Added Kwarg to RV_spectrum.rv_fit
is RV_guess_var, which is by default 0.
It describes the min/max variation of the initial RV guess for each fit.
New in version 1.0: The release version as originally published in Zeidler et al. 2019.
New in version 1.1: Introducing the Kwargs rv_sys, which should be used for large initial RV offsets from rest frame or for redshifted spectra. rv_sys should provide an estimate of that shift.