A new method based on Markov chains for deriving SB2 orbits directly from their spectra
We present a new method to derive orbital elements of double-lined spectroscopic binaries (SB2). The aim is to have accurate orbital parameters of a selection of SB2 in order to prepare the exploitation of astrometric Gaia observations. Combined with our results, they should allow one to measure the mass of each star with a precision of better than 1 %.
The new method presented here consists of using the spectra at all epochs simultaneously to derive the orbital elements without templates. It is based on a Markov chain including a new method for disentangling the spectra.
A new method based on Markov chains for deriving SB2 orbits directly from their spectra
Keywords: binaries: spectroscopic
Double-lined spectroscopic binaries (SB2) form the basis of the least model–dependent methods used to derive stellar masses. Their orbital elements are employed to derive the products , where is the mass of a component and is the inclination of the orbital plane. Therefore, when the inclination may be obtained from another technique, such as eclipse observations or astrometric measurements, the accuracy of the masses depends on that of the SB2 orbital elements.
The derivation of the elements of a spectroscopic orbit is usually based on the derivation of the radial velocities (RV) of the components from the blended spectra of the binary star. This is generally done by fitting one (Hill 1993) or two spectral templates to each spectrum (Zucker & Mazeh 1994). The approach is justifiable when the templates are similar to the actual spectra of the components, but otherwise the results are degraded by template mis-match.
Other methods do not use templates, but are based on the disentangling of the spectra (Simon & Sturm 1994; González & Levato 2006). The method presented hereafter falls into this category. It is Bayesian, since it is assumed that the radial velocities obey Keplerian motions, and the orbital elements are derived at the same time that the spectra are disentangled. It is based on Markov chains and on an new disentangling technique. The method is described in Sect. 2 hereafter. In Sect. 3, it is applied to HD 89745.
2 A BiT MOrE with MCMC
A BiT MOrE with MCMC is an Algorithm for Binaries To Match Orbital Elements with Markov Chain Monte Carlo. It combines a MCMC routine and a disentangling method.
2.1 The MCMC method
Our method is based on a Markov chain Monte Carlo (MCMC) survey of model parameter space (see, e.g., MacKay 2003), which samples the free parameters and produces a posterior probability density function (PDF) for each dimension. To get that, the MCMC routine investigates randomly the space of parameters. Thus, it is well adapted for our case with 6 orbital elements (the systemic velocity is lost, see next section).
The MCMC algorithm requires an initial guess for the elements of the spectroscopic orbit: the amplitudes of the RV of the two components, the orbital period, the eccentricity, the epoch and the argument of the periastron. Subsequently, the MCMC routine solves the orbital equations to derive the RV for each time of observation.
2.2 The disentangling method
The goal of this step is to build a model from the orbital elements. For that, we developed a new spectral disentangling algorithm acting in Fourier space in order to extract each spectral component. It is different from iterative procedures like González & Levato (2006) or Hadrava (1995) because it does not use references line or least-square fitting.
Actually, the spectra () are the sum of the two components ( and ) which are shifted ( and ) according to their orbital RV at the time of observation (equation 2.1).
The Doppler shift is used to disentangle the two components in Fourier space where the gap between each epoch is just a phase shift. The Fourier transform is calculated using Fast Fourier Transform techniques (e.g., Press et al. 1992). So, it is possible to determine, right from the initialisation step, and . Then, the two disentangled spectra are summed to obtain the model spectra. Note that information on the values of continua are lost in passing to Fourier space. Moreover, the systemic velocity is not derived because only RV shifts are involved.
Finally, the model is compared to the actual spectrum to test the likelihood of the solution. It is injected again at the beginning of the algorithm to find the next parameters. This procedure is reiterated until the Markov chain converges.
To validate our method, synthetic spectra have been built assuming a set of orbital elements. In order to be sure that the method is reliable, the computation was started assuming initial elements quite different from those used in the simulation. The algorithm was run to attempt to recover the input parameter values and their uncertainties.
Our method was applied to a set of 11 spectra of HD 89745, taken over 1785 days, corresponding to about 77 % of the period. These spectra arise out of a large observation program at the Haute-Provence observatory with the T193/Sophie (Halbwachs & Arenou 2009), and also from another program carried out with the same instrument (Halbwachs et al. 2012). This star was chosen for the number of spectra, although it appears somewhat atypical. In the past, it received 27 Coravel observations, and an orbit was derived after adding 7 Sophie spectra (Halbwachs et al. 2012). For this calculation, all the RV were derived from the cross-correlation function (CCF) of the spectra with a single mask or a simple template, the 1-mask technique. However, when the TODMOR algorithm was applied to a subset of 9 Sophie spectra, the two-template CCF yielded RVs quite different from those obtained with a single mask. The consequence was that the RV semi-amplitudes, and , were quite different from the previous values: around 11 and 8 km/s instead of 7 and 11 km/s, respectively. The mass ratio of the binary is then changed from 0.62 to 1.30. If this discrepancy is due to a bias in the RV estimations obtained through a single template, one would expect to find it again with a disentangling technique.
A small wavelength range, just 3 orders, is taken into account between 5490.8 and 5787.9 Angströms to avoid telluric lines (for a the present assessment). All spectra have a good signal to noise ratio (higher than 50). The disentangling of the two components and their associated model are shown for one observation in Fig. 1. It was obtained after only 20000 iterations in A BiT MOrE with MCMC. From the PDF (see Fig. 2), the uncertainties at one sigma around the best likelihood value can be obtained (listed in table 1).
|Orbital elements||Halbwachs et al. (2012)||TODMOR||A BiT MOrE with MCMC|
The values obtained with A BiT MOrE with MCMC are roughly in agreement with the previous values. We have large uncertainties due to the small number of spectra, the incomplete period coverage and the fact that only three orders are taken into account. However, we confirm that the mass ratio could be inversed, as found with the TODMOR reduction. It is worth noticing that this does not come from a permutation between the primary component and the secondary one, since the orbital elements and are in agreement with those found with the 1-mask technique. Nevertheless, these results are only preliminary and more observations are needed to draw reliable conclusions.
The simulations have shown that A BiT MOrE with MCMC is a powerful new tool to perform orbital element calculations. This method should provide reliable orbital elements and recover the spectra of the individual stars at the same time.
We have also seen that the orbital elements obtained from a single template are only roughly consistent with those obtained when two templates are used or when the spectra are disentangled. For HD 89745, using two templates and disentangling both leads to orbital elements rather different from those obtained from a single template, and we plan to derive a new orbit when more spectra become available.
However, more detailed investigations are required to compare the reliability of the three techniques considered here. This will be done in the future in order to obtain the most accurate masses when Gaia measurements are available.
Thanks to AS Gaia for functioning support. The observations were successful thanks to the help of the staff of the OHP. It is a pleasure to thank Eric Gosset for valuable discussions.
- González & Levato (2006) González, J. F. & Levato, H. 2006, A&A, 448, 283
- Hadrava (1995) Hadrava, P. 1995, A&AS, 114, 393
- Halbwachs & Arenou (2009) Halbwachs, J.-L. & Arenou, F. 2009, in SF2A-2009: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics, ed. M. Heydari-Malayeri, C. Reyl’E, & R. Samadi, 53
- Halbwachs et al. (2012) Halbwachs, J.-L., Mayor, M., & Udry, S. 2012, MNRAS, 422, 14
- Hill (1993) Hill, G. 1993, in Astronomical Society of the Pacific Conference Series, Vol. 38, New Frontiers in Binary Star Research, ed. K.-C. Leung & I.-S. Nha, 127
- MacKay (2003) MacKay, D. 2003, in Information Theory, Inference, and Learning Algorithms (Cambridge University Press)
- Press et al. (1992) Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 1992, Numerical recipes in FORTRAN. The art of scientific computing
- Simon & Sturm (1994) Simon, K. P. & Sturm, E. 1994, A&A, 281, 286
- Zucker & Mazeh (1994) Zucker, S. & Mazeh, T. 1994, ApJ, 420, 806