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 premainsequence 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 crosscorrelation. Together with a Monte Carlo bootstrap technique, it is possible to measure radial velocities with an accuracy of \(13\,\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 hypervelocity stars an option is given that the wavelength limits are automatically adjust for each iteration \(n\) based on the solution \(n1\) 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 hypervelocity 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 restframe
 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 restframe
 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 reinitiating 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 hypervelocity 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 sigmalevel 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 pep8
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.