EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 641 (September 2020) A&A, 641 (2020) A26 Full HTML
Free Access
Issue
A&A
Volume 641, September 2020
Article Number A26
Number of page(s) 16
Section Stellar structure and evolution
DOI https://doi.org/10.1051/0004-6361/202038446
Published online 01 September 2020

© ESO 2020

1. Introduction

Stars with initial masses Mi ≳ 8 M⊙ are expected to end their lives as core-collapse supernovae (SNe), forming compact stellar remnants such as neutron stars and black holes (BHs) (Heger et al. 2003). Through their strong stellar winds, intense radiation and powerful explosions, they dominate the transport of energy and momentum into the interstellar medium (cf. Mac Low et al. 2005). A small subset of massive stars appear as Wolf-Rayet (WR) stars. They are a spectroscopic class of stars dominated by broad emission lines that are a tell-tale sign of strong, radiation driven winds (see review by Crowther 2007). WR stars are subdivided into three main classes: nitrogen-rich (WN), carbon-rich (WC), and the extremely rare oxygen-rich (WO) sequences. They reflect the chemical composition of their atmospheres - N-rich (products of the CNO cycle) or C/O-rich (products of He-burning).

Due to the absence of a hydrogen envelope and their chemical composition, most WR stars are believed to be evolved, post-main sequence objects. These objects are called “classical WR stars” (cf. Massey 1981) and they are thought to evolve spectroscopically from O-type stars to WN, WC and, finally, WO stars, possibly through an intermediate luminous blue variable or red supergiant phase (Conti 1976; Crowther 2007). However, very massive stars can also appear as WR stars during their main-sequence phase (de Koter et al. 1997). In contrast to the classical WRs, young WR stars contain and display hydrogen in their atmospheres.

The role of multiplicity in the formation of WR stars is a long-standing question that is still widely debated (Paczyński 1967; Vanbeveren et al. 1998; Shenar et al. 2019, 2020). It is still not clear whether the dominant channel is through stripping via stellar winds or through a binary companion. This knowledge and, hence, the study of WR stars, is important to improve our understanding of observed SNe (Zapartas et al. 2017; Eldridge et al. 2018) and gravitational-wave (Eldridge et al. 2019) rates, the properties of the resulting compact objects (e.g. de Mink et al. 2014; Marchant et al. 2016) and the evolution of massive stars (e.g. Meynet & Maeder 2003; Hamann et al. 2006; Shenar et al. 2016; Sander et al. 2019; Woosley 2019).

Large-scale multiplicity surveys in the past decade have shown that the majority of massive stars live in binaries or higher multiple systems (Mason et al. 2009; Sana et al. 2012, 2013, 2014; Kobulnicky et al. 2014; Barbá et al. 2014; Dunstall et al. 2015; Maíz Apellániz et al. 2016, 2019). These systems are so close that the majority will interact (Fryer et al. 2007; Sana et al. 2012, 2013; Dunstall et al. 2015), fundamentally impacting the evolution of most of these systems (Pols 1994). However, most modern radial velocity (RV) surveys of massive stars in the Galaxy have focussed on main-sequence OB-stars, with the exception of WN stars from the OWN survey (Barbá et al. 2014). The state of affairs for extragalactic surveys is much better, with numerous surveys of the Magellanic Clouds (Bartzakos 1999; Foellmi et al. 2003a,b; Schnurr et al. 2008) and the like for M31 and M33 (Neugent & Massey 2014).

This is not without good reason, as surveys of WR stars face numerous challenges. First and foremost, the classical WR phase is quite short-lived in the evolution of massive stars (∼10% of their total lifetime). Combined with the fact that most WR stars that have been discovered in the past two decades are mainly infra-red bright, the sample size for an optical RV survey is severely limited. Second, spectral lines for WR stars are formed far out in the winds. They are therefore broad and affected by variability in the winds, rendering traditional tools to measure RVs unsuitable. Additionally, the complexity of each WR system along with the diverse physics that they present (e.g. wind-variability, clumping, dust production, etc.) has prompted a system-to-system modelling, making it difficult to apply a single method in a systematic fashion to a large sample. Finally, with input on the initial conditions of the systems, recent simulations have shown that post-interaction binary systems most likely exist with periods of the order of a year (Langer et al. 2019).

The Galactic Catalogue of WR stars1 (henceforth GCWR) is maintained by Paul Crowther (originally, Appendix 1 in Rosslowe & Crowther 2015), with 666 objects currently (v1.24). In the literature, the observed multiplicity fraction of WR stars is reported to be ∼0.4 (van der Hucht 2001). This value includes both spectroscopic and visual WR binaries. However, this fraction has not been corrected for observational biases and hence the intrinsic multiplicity of WR stars is yet to be constrained. An even worse situation stands for the distributions of mass-ratios and orbital periods, which are important resources that can be used to calibrate the physics in binary stellar evolution codes (e.g. efficiency of mass-transfer and angular momentum transport) as well as what is used in stellar population synthesis (e.g. BPASS, Eldridge et al. 2017). Accurate RV measurements and well-constrained uncertainties are of critical importance in order to obtain a reliable correction for observational biases.

Binary WR candidates are usually identified through Doppler shifts, the dilution of emission lines (the companion OB-star can sometimes be more luminous in the optical), the presence of absorption lines in the spectra, photometric variability, X-ray spectroscopy for hints of RV variability and wind-wind collision, and with astrometry. Modern techniques to detect periodicities in WR stars include measuring the equivalent-width and change in the full-width half-maxima of emission lines (e.g. WR1; Morel et al. 1999). Other state-of-the-art techniques to measure RVs include measuring the centroid of the emission (Fahed & Moffat 2012), fitting Gaussians to line profiles (e.g. WR127; de La Chevrotière et al. 2011), using (Fourier) cross-correlation (Lefèvre et al. 2005; David-Uraz et al. 2012; Shenar et al. 2019; Chené et al. 2019), or even modelling excess emission (e.g. Hill et al. 2000).

The variable nature of these strong, broad emission lines of WR stars have traditionally made it challenging to detect velocity shifts of the order of a few km s−1. These shifts are of the order of a fraction of a resolution element in low-resolution spectra and so they are extremely difficult to detect even with a high signal-to-noise (S/N) ratio. In terms of the parameter space for binary systems, these RV shifts correspond to systems that have long periods, and those with extreme mass-ratios.

Clumping and time-variability studies, however, have shown that emission lines are variable on the scale of days to weeks. They have further shown that this variability is different for different ions and is mostly limited to the centers of emission lines (Lépine et al. 2000; St-Louis et al. 2009; Chené & St-Louis 2011). Moreover, studies with high-resolution and high-S/N spectra have reported RV measurement accuracies of WR stars with formal uncertainties as small as ∼1 km s−1 (e.g. WR 113; David-Uraz et al. 2012). This is mainly due to the fact that high-resolution spectra contain more information, and this allows for a statistically more accurate solution. Furthermore, the high-ionisation “weak” lines can also be used to derive reliable RV measurements (e.g. R144, WR21a and R145; Sana et al. 2013; Tramper et al. 2016; Shenar et al. 2017, respectively). This is an advantage for systems with wind-wind collision, as the weaker lines are usually unaffected. Finally, high-resolution spectra make it possible to detect any narrow absorption lines that could be superimposed on top of the broad emission lines.

In this context, we initiated a new monitoring program aiming at improving the observational constraints on the WR multiplicity properties. In this first paper, we present and test out RV measurement method and apply it to a magnitude-limited sample of Northern, Galactic WC stars. The sample, observations and data reduction are presented in Sect. 2. Section 3 elaborates on the RV measurements and relevant masks used for cross-correlation. In Sect. 4, the correction for observational biases is discussed, along with the spectroscopic binary fraction and intrinsic period and mass-ratio distributions. Section 5 presents our conclusions.

2. Overview of the sample and data reduction

2.1. The sample

Our sample was selected from the GCWR and contains all WC stars with V ≤ 12 and declination δ ≥ −30° that were selected from the GCWR. For objects with missing V band magnitudes, narrow-band v magnitudes (Smith 1968; Massey 1984) were used with the condition v ≤ 13. This resulted in a final sample of 12 WC stars. The spectral type distribution among our sample is shown in Fig. 1.

thumbnail Fig. 1.

Spectral type distribution of the stars in our sample.

2.2. Observations

The WR RV monitoring campaign began in 2017 using the High-Efficiency and high-Resolution Mercator Echelle Spectrograph (HERMES, Raskin et al. 2011) mounted on the 1.2 m Flemish Mercator Telescope on La Palma. HERMES covers the optical wavelength range from 3800 Å to 9000 Å, with a spectral resolving power of R ∼ 85 000. As we show below, the combination of high-resolution with good S/N over a long period of time results in a powerful tool to search for binarity.

For most objects, at least four epochs were obtained over the course of 15 months. Table 1 gives the an oerview of the campaign (including archival HERMES data). The number of epochs for WR 137 is larger because we obtained an additional 16 epochs with a high cadence over 15 days, with three observations over the same night to investigate the effect of short-term atmospheric/wind variability on the measured RVs (Sect. 3.3).

Table 1.

Sample of WC stars in the RV monitoring campaign with the number of spectra, time coverage, and average S/N per resolution element at 5200 Å.

2.3. Data reduction and normalisation

The HERMES spectra were first reduced by the HERMES pipeline which includes the bias, flat-field corrections, cosmic removal and order merging of the spectra (Raskin et al. 2011). However, the resulting spectra are still affected by atmospheric extinction, telluric lines, interstellar reddening and the instrumental response. In order to achieve a normalisation procedure that is homogeneous over multiple epochs, we developed an algorithm that corrects for these effects. This minimises the effect of variable normalisation on the RV measurements.

We first correct the spectra for atmospheric extinction, applying a wavelength and airmass-dependent correction. Molecfit (Smette et al. 2015; Kausch et al. 2015) is then used to remove telluric lines with meteorological conditions as input from the Mercator meteo-station. Standard calibration stars that are used for the order merging are then used to derive the instrumental response function for the night. Once corrected for these three effects, the only remaining effect is due to that of interstellar (or even circumstellar) reddening.

To normalise these spectra, pseudo-continuum regions were chosen based on normalised model spectra for WC stars from the Potsdam Wolf-Rayet code (PoWR: Gräfener et al. 2002; Hamann & Gräfener 2003; Sander et al. 2012). Two continuum regions were adopted: in the red, around 8100 Å, and in the blue, around 5200 Å. By anchoring the respective continuum models to the red point and then reddening them based on Fitzpatrick (2004) we fit the slope and normalise the spectra. For this purpose, we used an RV value of 3.1, the average value for the Galaxy. The various steps in the data handling are illustrated in Fig. 2. We note that we only rely on the above-mentioned reddening law to obtain a homogeneous normalisation process and therefore refrain from providing the reddening values obtained though the continuum-fitting procedure.

thumbnail Fig. 2.

Main steps in our data reduction process. From top to bottom: (a) flatfield corrected spectrum produced by the HERMES pipeline. (b) Spectrum corrected for the instrumental response function and telluric contamination. (c) Continuum flux models from PoWR to the pseudo-continuum region around 5200 Å for different values of E(B − V). (d) Resulting normalised spectrum.

3. Radial velocity determination

In order to derive RVs for broad WR emission lines in a systematic manner, one needs a method that can be applied to most, if not all, of the objects in our sample and that yields meaningful errors. The uncertainty on each measurement will have contributions from different sources: statistical (based on the method and S/N), systematic (instrumental and normalisation), and intrinsic (variability in spectral lines). It is important to have reliable estimates of the uncertainties affecting the RV measurements as these will directly propagate into the observational bias correction. For binary systems, an upper limit on the RV measurement uncertainties can be estimated by measuring the dispersion in the residuals after fitting an orbit. However, the very long periods involved do not allow us to use this.

In this study, we focus solely on RV measurements of WR stars. In the case of SB2 binaries (i.e. double-lined spectroscopic binaries), it is, in principle, also possible to measure the RVs of the secondary star (typically an OB-type star). However, as the goal of our current study is binary detection rather than orbital analysis, there is little advantage in doing so. The WR component, being typically the lighter component, would exhibit larger RV amplitudes than its OB-type companion. Moreover, applying additional binary-criteria to SB2 binaries would require us to introduce non-trivial adjustments to our bias correction (see Sect. 4.3.2).

In all cases, the estimation of uncertainties depends on the method used to measure the RVs. In the case of Gaussian fitting, the errors on the RV measurements are usually derived using the co-variance matrix of the least-squares fit. Cross-correlation relies on applying a peak-finding algorithm to the cross-correlation function (CCF), which can be a Gaussian, Lorentzian, parabola or, a sinc function (e.g. Tody 1993). Another approach is to find the bisector of the CCF. The errors are derived based on the width, number of pixels, height of the CCF peak, and the antisymmetric noise (Tonry & Davis 1979). When measured for WR stars with low-resolution spectra, the errors are magnified by the low number of pixels around a broad CCF peak. In this work, we use a different statistical formulation for the CCF that represents it as a log-likelihood function. The latter is then maximised and fit with a parabola. This is explained in detail below.

3.1. Cross-correlation

The method that we used for the cross-correlation is based on Zucker (2003). The method, assumptions and limitations are briefly discussed here. For a more thorough derivation we refer to Zucker (2003). The observed spectrum is denoted by f(n) and is assumed to be Doppler shifted with respect to the “template” of zero-shift, denoted by g(n). Both the spectrum and the template are described as a function of pixel, n, such that a Doppler shift results in a uniform, linear shift in the spectrum (Tonry & Davis 1979).

In this work, the observed spectrum is assumed to be reproducible by scaling the template with a constant (a0) and shifting it by a certain number of pixels (p0). The noise in the observed spectrum is assumed to be Gaussian with a fixed standard deviation (σ0), so that:

f ( n ) = a 0 g ( n − p 0 ) + d n , with $$ \begin{aligned} f(n) = a_0 g(n-p_0) + d_n \mathrm{, with} \end{aligned} $$(1)

d n ∼ N ( 0 , σ 2 ) . $$ \begin{aligned} d_n \sim N(0,\sigma ^2) {.} \end{aligned} $$(2)

The spectra (or specific masks) are assumed to be “continuum-subtracted”, such that they have a zero mean:

∑ n f ( n ) = 0 , and $$ \begin{aligned} \sum _n f(n) = 0 \mathrm{, and} \end{aligned} $$(3)

∑ n g ( n ) = 0 . $$ \begin{aligned} \sum _n g(n) = 0 {.} \end{aligned} $$(4)

The last assumption that the number of pixels, N, is constant breaks down. For very small velocity shifts, the edge effects are negligible. However, for velocities of the order of 1000 km s−1, the assumption is not applicable and hence, we use a rolling function to carry out the shift in pixels. This does not have any impact on the computation of the CCF or the determination of the peak.

In order to express the CCF in terms of a log-likelihood function, for a given shift of p pixels we first define the cross-covariance function,

R ( p ) = 1 N ∑ n f ( n ) g ( n − p ) . $$ \begin{aligned} R(p) =\frac{1}{N}\sum _n f(n) g(n-p){.} \end{aligned} $$(5)

The standard deviation of the spectrum (sf) and template (sg) can be expressed as

s f 2 = 1 N ∑ n f 2 ( n ) , and $$ \begin{aligned} s_f^2 = \frac{1}{N} \sum _n f^2 (n) \mathrm{, and} \end{aligned} $$(6)

s g 2 = 1 N ∑ n g 2 ( n ) . $$ \begin{aligned} s_g^2 = \frac{1}{N} \sum _n g^2 (n) {.} \end{aligned} $$(7)

Following this, we define the CCF of f and g as

C ( p ) = R ( p ) s f s g · $$ \begin{aligned} C(p) = \frac{R(p)}{s_f s_g}{\cdot } \end{aligned} $$(8)

This CCF can be written in the form of a log-likelihood function. We denote a shift p ̂ $ \hat{p} $ that maximises the CCF. The errors can then be derived from the covariance matrix following maximum log-likelihood theory and can be expressed as

σ p 2 = − [ N C ″ ( p ̂ ) C ( p ̂ ) C 2 ( p ̂ ) 1 − C 2 ( p ̂ ) ] − 1 , $$ \begin{aligned} \sigma _p^2 = - \left[N \frac{C{\prime \prime }(\hat{p})}{C(\hat{p})} \frac{C^2(\hat{p})}{1 - C^2(\hat{p})} \right]^{-1} {,} \end{aligned} $$(9)

where C″(p) is the second derivative of the CCF. Equation (9) shows that the errors depend on the number of pixels N, the sharpness of the peak C ″ ( p ̂ ) / C ( p ̂ ) $ C{{\prime\prime}}(\hat{p}) / C(\hat{p}) $ and the line S/N ratio C 2 ( p ̂ ) / ( 1 − C 2 ( p ̂ ) ) $ C^2(\hat{p}) / (1 - C^2(\hat{p})) $.

The success of the cross-correlation depends on the accuracy and stability of the template. This is due to the fact that the template is implicitly assumed to represent the data accurately. In practice, the template is chosen to be the observations with the highest S/N. A parabola is fit to the top of the CCF in order to determine the peak and the related derivatives. The errors from the least-squares fitting of the parabola agree with those derived from Eq. (9), giving confidence in the validity of the method. Once the RVs and measurement uncertainties are derived, a weighted average of the shifted spectra is used to derive a co-added spectrum that has a higher S/N. This is then used as a template to re-derive the relative RVs and, thus, in an iterative fashion, an ideal mask and a final measurement are converged upon. We note that the RVs derived are relative to a particular observation in time, and that the velocities measured are, by no means, absolute values. In order to compare with those in the literature or to the rest frame, a uniform shift must be imposed.

3.2. Mask selection

After addressing the choice of the template, the question of which spectral lines or wavelength ranges (henceforth “masks”) to be used for the RV measurements arises. Most traditional methods of measuring RVs are applied to the strongest lines, as they contain the most signal (see Sect. 1). On the other hand, different lines originate in different ions and form in different regions in the wind. Therefore, by comparing the RVs measured using different sets of lines, we can both probe different parts of the wind as well as identify the most reliable mask to represent the motion of the star.

However, with high-resolution and good S/N, weak lines also contain significant amounts of information that can be used. As the WR spectra contains much broader (and stronger) spectral features, using the entire spectrum in the cross-correlation process can mitigate the problem of contamination from the secondary (David-Uraz et al. 2012). To investigate which masks are the most appropriate and if individual lines are better than the full spectrum, RVs are first measured using a variety of masks. These masks are created by selecting lines of different strengths and different ions, as well as with the full spectrum.

For WC stars, most lines originate from carbon (C II, C III and C IV), oxygen (O II, O III, O V, O VI) and helium (He I, He II). Naturally, the strength of the lines of various ions depends on the temperature (among other things) of the WR star, and are identified star-by-star using PoWR models with parameters from Sander et al. (2012). For each star in the sample, this results in “strong”, “moderate” and “weak” masks for different ions. An example of C IV masks for WR 137 is shown in Fig. 3.

thumbnail Fig. 3.

Spectrum of WR 137, rebinned for clarity, with different masks of C IV indicated (see legend). This process is repeated for all other ions, producing a plethora of masks.

The motivation for using lines of the same ion in the same mask is straightforward, while that for using only strong or weak lines is a bit more complex. Lines of the same ion are thought to have formed in the same region of the wind and should thus provide the same RV measurement. Moreover, it is useful to use only strong lines to measure RVs in order to compare the results of this work to what is available in the literature. The errors from the CCF (Eq. (9)), additionally, depend inversely on the S/N in the line, which is higher for strong lines. However, weak lines of highly ionised ions are formed closer to the surface and thus are thought to be less susceptible to wind variability, making them a logical choice for a mask. It is further essential to isolate and quantify the effects of wind variability on WR spectra. Once this is accounted for, we can decide which masks are the most appropriate to detect genuine RV variability. This is addressed in the following section.

3.3. Uncertainties from systematic sources and wind-variability: A proto-study of the binary WR 137

In this section, we investigate the effect of systematic sources and wind variability on RV measurements. To understand the extent to which wind variability can affect our measured RVs, WR 137 (WC7pd + O9 V) is an ideal candidate to examine. It is a well-studied long-period binary system with an orbital period of 13.05 years (Lefèvre et al. 2005; Richardson et al. 2016). We obtained 16 spectra over the course of 15 days during August 2019, significantly away from its periastron passage (which will happen in 2022). Therefore, we can neglect the effects of orbital motion and wind-wind collision (e.g. Lührs et al. 2002) due to binary interaction on the RV measurements for the purpose of this short-cadence study.

The RVs were measured as described in the previous subsection with different masks, with the additional condition that lines contaminated by absorption from the companion were excluded. The RV measurements using different masks with strong C IV lines, weak He II lines and the full spectrum can be seen in Fig. 4. Given RV measurements (vi ± σi), we calculated the weighted mean (μRV) and standard deviation (σRV) as:

w i = 1 σ i 2 , $$ \begin{aligned} { w}_i = \frac{1}{\sigma _i^2}{,} \end{aligned} $$(10)

μ RV = ∑ i w i v i ∑ i w i , $$ \begin{aligned} \mu _{\rm RV} = \frac{\sum _i { w}_i { v}_i}{\sum _i { w}_i} {,} \end{aligned} $$(11)

σ RV = ( ∑ i w i ( v i − μ RV ) 2 ∑ i w i ) 1 2 , $$ \begin{aligned} \sigma _{\rm RV} = \left( \frac{\sum _i { w}_i ({ v}_i - \mu _{\rm RV})^2}{\sum _i { w}_i} \right)^{\frac{1}{2}} \text{,} \end{aligned} $$(12)

thumbnail Fig. 4.

Phased RV plot for WR 137. Archival RV measurements are shown in blue and are taken from Lefèvre et al. (2005). The RV measurements from HERMES using the strong C IV lines are shown in red. Along with these measurements, the inset shows the short cadence measurements with the full spectrum (green) and weak He II lines (purple). All HERMES RV measurements have been shifted by −23.5 km s−1 in order to convert relative RV measurements to absolute values.

where σRV is a measure of the real dispersion in the data. Since we have already ruled out orbital and wind-wind collision effects, this dispersion can be assumed to represent the statistical error (σp) along with some contribution from systematic sources and wind variability (σw). These errors are not expected to be correlated and are therefore added in quadrature, such that they can be expressed as follows:

σ RV = σ p 2 + σ w 2 . $$ \begin{aligned} \sigma _{\rm RV} = \sqrt{\sigma _p^2 + \sigma _{\rm w}^2} {.} \end{aligned} $$(13)

Accounting for the wind variability provides us with an estimate on the effective uncertainty on the RV measurements (σRV). As σp is already derived, σw can be calculated for different ions and different line strengths, giving us a direct probe of how wind variability affects these lines. For now, however, it is sufficient to inspect σRV for different ions. Figure 5 depicts this comprehensively, plotting the RVs measured by each mask over the short cadence period, with σRV indicated at the bottom of the plot.

thumbnail Fig. 5.

RV measurements for WR 137 during the short cadence study. At the bottom of the plot are the values of (σw, σRV) (km s−1) for each mask. It is observed that the helium lines (and weak lines in general) have large dispersions. The masks least affected by wind variability are “full spec”, “C4 strong”, “C4 moderate”, and “C4 all” (determined by their values of σw).

Upon examination, two observations can be made based on Fig. 5. The first is that for masks with He II lines, a large scatter is seen independent of the strength of the lines. The second is that across other ions, masks with strong lines perform better than those with weak lines. It is not trivial to explain these observations, as there are multiple factors that are at play for each line (e.g. oscillator strength, number populations of the particular state, line-formation region, etc.).

In order to ascertain that the statistical errors do not bias our mask selection, we calculated values for σw using Eq. (13). The results are twofold – masks with ions of a higher ionisation state are the least affected by wind variability and strong lines are affected less than weak lines. Finally, by comparing the absolute values of σw in Fig. 5, the full spectrum mask is least affected by wind variability. For strong C IV lines, the average σp is larger than the corresponding value of σRV, implying that the wind variability is within the statistical error (and hence σw = 0). Given that the winds of all WC stars do not necessarily exhibit the same characteristics, this exercise has to be conducted for each object individually in order to select a specific ion to use. Upon inspection, diagnostic plots such as Fig. 5 for each object indicate that the full spectrum shows the least scatter. Considering the fact that the WR component dominates the spectrum for all of our objects, we decided to use the full spectrum to measure RVs.

4. Results and discussion

4.1. Observed spectroscopic binary fraction

For all the objects in our sample, RVs were measured using multiple masks. The tests described in the previous section reveal that the measurements using the full spectrum are the most reliable and we adopt those in the following section. These are given in Appendix B. To identify objects with significant RV variation, we look at the maximum RV variation for each star, given by

Δ RV = max ( | v i − v j | ) , $$ \begin{aligned} \Delta \mathrm{RV} = \mathrm{max}(|{ v}_i - { v}_j|), \end{aligned} $$(14)

where vi and vj are the RV measurements for two epochs i − j for a given star. One then imposes a minimum amplitude threshold C for RV variations computed from any individual pair of epochs to be “significant”.

Δ RV ≥ C . $$ \begin{aligned} \Delta \mathrm{RV} \ge C. \end{aligned} $$(15)

Unlike in Sana et al. (2013), we do not enforce a significance criteria for each pair of RV measurements (Eq. (4) in Sana et al. 2013 compared to Eq. (15) above). This is due to the fact that the values for σp are typically very small (sub-km s−1) and underestimate the true error (see Sect. 3.3).

Given this criterion, the observed fraction of RV variable objects (fobs) can be calculated as a function of threshold C (Fig. 6). It is encouraging to see all the known binary systems group together on the right side of the plot. We can now take advantage of the fact that a subset of our sample are confirmed binaries to constrain a suitable value for C.

thumbnail Fig. 6.

Observed binary fraction (fobs) as a function of the threshold C. The adopted threshold is the vertical dotted-dashed line at 10 km s−1, resulting in fobs = 0.58. The objects are labelled as and when they are classified as binaries based on Δ RV. Names in bold are known binaries.

thumbnail Fig. 7.

Top: phased RV plot for WR 140. Archival RV measurements are shown in blue (Fahed et al. 2011) and green (Thomas et al. 2020). The RV measurements from HERMES (red) using C IV strong lines have been shifted by 19 km s−1 in order to convert relative RV measurements to absolute ones. Bottom: same for WR 113. Archival RV measurements are shown in blue (Hill et al. 2018) and green (Massey & Niemela 1981). The RV measurements from HERMES (red and yellow) have been shifted by −135 km s−1 in order to convert relative RV measurements to absolute values.

Two kinks can be seen, in the ranges of 5 to 10 km s−1 and of 12 to 19 km s−1. A threshold value of 10 km s−1 leads to an observed binary fraction of fobs = 0.58, and that of 15 km s−1 results in fobs = 0.33. To include the known and candidate binaries in our binary fraction, we choose the former threshold of 10 km s−12.

4.2. Comparison with the literature

Although the magnitude-limited sample of WC stars in this study is small, the even distribution across spectral types (Fig. 1) omits any selection bias. With this in mind, two important conclusions can be drawn from Fig. 6. The first is that our value of fobs is comparable to the value of 0.56 observed for main-sequence O-type stars by Sana et al. (2012).

The second is that there seems to be a lack of short-period (P < 100 d) binaries compared to the overabundance of short-period systems observed for main-sequence O-type stars (Sana et al. 2012). This is especially surprising given that the observational bias favours the detection of short-period WR binaries, which are expected to exhibit high RV amplitudes. Of the known binaries in our sample, orbits are derived only for WR 137, WR 140 and WR 113, with only WR 113 being part of a short-period binary system (P ∼ 29.7 d). This low fraction of detected short-period binary systems (0.08) is in stark contrast to what has been found for main-sequence O-type stars.

In order to investigate the observed lack of short-period WC binary systems, we scoured the literature in order to expand our small sample. The GCWR reports 237 WC stars in the Galaxy, of which 87 were included in van der Hucht (2001). For the sake of considering WC stars that were thoroughly studied for binarity, we focus on this subset of 87 stars.

Of this subset, 37 WCs are classified as detected or probable binaries, resulting in an observed Galactic WC binary fraction of 0.41. However, this fraction is severely biased towards short-period systems. This is because historical studies were mainly sensitive to short-period systems, given typical RV measurement uncertainties of the order of ∼20 km s−1. Despite the large bias towards detecting short-period systems, 15 report periods shorter than 100 d, seven of which have derived orbits. Therefore, the fraction of Galactic WC stars that are conclusively in short-period binary systems is only 0.09, which is consistent with the findings of this study.

4.3. Instrinsic multiplicity properties

4.3.1. Binary detection probability

To estimate the sensitivity domain of our survey, we perform Monte-Carlo (MC) simulations following the method described in Patrick et al. (2019). In short, for a grid of orbital period and secondary masses ranging from 1 to 105 d and 1 to 20 M⊙ respectively, we adopt a flat eccentricity distribution between 0.0 and 0.9, a random orientation of the orbital planes in the three-dimensional space, random times of periastron passage, and RV measurements uncertainties of 1.6 km s−1 that are taken from the single stars in Table 2. Simulated RV time series of 2500 systems are computed at each log P − M2 grid cells, with steps of Δ log P = 0.1 and ΔM2 = 1 M⊙, allowing us to compute the probability that the given orbital configuration would yield a RV signal passing our detection threshold, i.e. max(ΔRV) > 10 km s−1, given the observational sampling corresponding to each specific object. WC primary masses of 10 M⊙ are adopted throughout the simulations (cf. Massey 1981). An example of the resulting detection probability maps is given in Fig. 8. Plots for all other sources in our sample are included in Fig. B.1. Detection probabilities well over 90% are reached up to periods of 102.5 to 103.5 d (i.e., 1 to 10 years) depending on the object, except at very low-mass companion masses. Figure 9 shows the combined detection probability as a function of the orbital period assuming a flat mass-ratio distribution between 0.5 and 2.0 (e.g., as in Langer et al. 2019) and indicates excellent detection rates close to 100% up to periods of 100 d, above 80% up to 1000 d, then dropping to 20% at 10 000 d.

thumbnail Fig. 8.

Detection probability for WR 137 as a function of the orbital period and secondary mass. A representative primary mass 10 M⊙ has been adopted for the WR star.

Table 2.

Results of the RV measurements for the sample of WC stars.

4.3.2. Intrinsic multiplicity properties

The observed binary fraction is a minimum estimate on the number of true WC binaries in our sample. The number of undetected binary systems depends on our time sampling, our RV measurement accuracy, and the intrinsic distribution of orbital parameters. Given our limited sample size, it is beyond the scope of this paper (and probably unrealistic) to perform a self-consistent derivation of the distributions of orbital parameters, and of the binary fraction (as in, e.g., Sana et al. 2012, 2013). Nevertheless, it is still possible to obtain useful information on the multiplicity properties of the parent population. In this context, we adopt two different and complementary approaches, both using MC simulations and the same underlying assumptions as in Sect. 4.3.1.

As a first step, we attempt to constrain the true binary fraction among our sample of WC stars. A simple approach would be to correct the number of detected binaries using a detection probability factor that depends on the orbital period of the detected binaries. This probability can be readily read from Fig. 9. Unfortunately, only three of the objects above our binary threshold have known orbital periods (see Table 2). Lacking further information, we conservatively adopt the detection probability corresponding to period of 103.5 d (approximately 10 years) for the remaining objects (shorter periods can be excluded on the basis of the previous step, as they would have yield a RV-amplitude larger than 30 km s−1 in over 90% of the cases). Under these hypotheses, the obtained correction factor is found to be of the order of 3/2, suggesting an intrinsic binary fraction of about 0.9 for the objects in our sample in the considered parameter space (P up to 105 d and mass-ratios up to 3.0).

thumbnail Fig. 9.

Average detection probability curve of our campaign as a function of the orbital period and assuming a flat mass-ratio distribution between 0.5 and 2.0. Known orbital periods of WC stars in our sample are indicated by vertical lines.

To test these intriguing results, we set up a last experiment where we test different period distributions and intrinsic binary fractions with the same assumptions on q, i, and e as in Sect. 4.3.1. We reject combinations of orbital configurations and binary fractions that cannot reproduce a simple observational fact: six objects are observed with values of ΔRV between 10 and 30 km s−1 within 1-σ (Fig. 6)3. We constructed two experiments based on:

  • a uniform log P distribution with variable minimum and maximum boundaries, and

  • a normal log P distribution with variable central period and 1σ-dispersion.

For these experiments, we removed WR 113 from the sample as its significantly different orbital period would require more complicated functional forms than those adopted here. In that case, we also request that no system with an RV variation amplitude above 30 km s−1 be found.

In both cases, parent distributions with intrinsic binary fractions larger than 90% performed significantly better than those with lower binary fractions. This is in line with our above estimate based on the detection probability.

4.3.3. Intrinsic binary fraction of the parent WC population

Overall, the results of these MC experiments suggest that WC stars in our sample have a large binary fraction and are dominated by long-period systems. The present sample was selected based on declination and magnitude limits and should not contain any specific selection bias. Therefore, it should remain representative of the entire WC population. Under these considerations, a test-of-hypothesis at the 10%-significance level allows us to reject parent binary fractions smaller than 0.724. This suggests that, with a 90% confidence limit, the true binary fraction of the carbon-sequence WR population in the Milky Way lays in the range 0.72 to 1.00. In contrast to the results in Sect. 4.3.2, the latter value does account for the sample size.

5. Evolutionary implications and conclusions

With a small, magnitude-limited sample of northern Galactic WC stars, we have demonstrated that using cross-correlation with high-resolution spectra allows us to derive accurate RV measurements of WR stars, and to even quantify the wind-variability. With a sample of 12 stars spread evenly over different spectral-types, we find an observed binary fraction of 0.58. Using MC simulations, we attempt to gain information about the parent WC population and conclude that the fraction of Galactic WC stars that are part of a binary system is at least 0.72, significantly larger than the literature value of ∼0.4.

Our observations clearly show that there is a deficiency in the number of short-period WC binaries in the Galaxy. Furthermore, our MC simulations indicate that Galactic WC stars primarily exist in long-period binaries (P > 100 d). Comparing these results to the observed population of main-sequence O-type binaries can help us put constraints on the viable evolutionary channels for WC stars and on the forward evolution of massive binaries. We briefly discuss these evolutionary channels below.

O-type stars are generally expected to evolve through the modified Conti scenario (O→LBV→WN→WC→WO, Crowther 2007), although the dominant mechanism leading to the formation of WR stars is still debated. Keeping this in mind, we briefly look at the outcome of possible evolutionary channels for short-period O-type binaries. These binaries can interact through Roche Lobe overflow (RLOF), common-envelope evolution, or even through mergers. We explore possible mechanisms that can explain the formation of long-period WC binaries when starting from short-period O-type binaries.

A first viable scenario that produces long-period binaries is the RLOF channel. Depending on the initial orbital configuration, the primary may undergo Case-A (on the main sequence) or early Case-B (post main sequence) mass transfer. The companion will, thus, strip the primary of its hydrogen envelope, leading to the formation of a WN star. Except in the case of a low-mass companion, this results in a mass-ratio inversion, irrespective of whether the mass transfer is conservative or non-conservative. This then leads to the widening of the orbit. Further loss of mass through winds during the entire WN and the early WC phase adds to the widening of the orbit.

The simulation study by Langer et al. (2019) explores the RLOF channel, albeit with the aim of constraining the period distribution of OB+BH binary systems (top panel of their Fig. 6). Given the fact that WC stars are thought to be the direct progenitors of BHs and assuming that the orbit does not undergo drastic change during the short-lived WC phase, this distribution should be valid for WC binary systems. The distribution has a small peak at short periods (log P ∼ 0.7) and a larger peak at long periods (log P ∼ 2.2). It is important to note that the simulations in Langer et al. (2019) are performed with MESA binary evolution tracks at lower metallicity, for the Large Magellanic Cloud. At Galactic metallicity, the corresponding evolutionary tracks will have higher mass-loss rates and hence, shift the distribution towards longer periods. Coupled with the fact that the theoretical distribution of Langer et al. (2019) can contain BHs formed from WN stars, it is plausible that our observations sample the tail of the period distribution.

Another evolutionary scenario focuses on common-envelope evolution and the possibility of a merger. Common-envelope evolution leads to the shrinking of the orbit, resulting in a very short-period binary system that may merge. In the case of a merger, the binary system ceases to exist and may result in an apparently single WC star. This channel may partly be responsible for the formation of the apparently-single WC stars in our sample. However, it is also possible that these stars formed as single stars.

A final channel to consider that produces long-period binaries invokes the existence of a hierarchical triple system (e.g. Sana et al. 2014; Maíz Apellániz et al. 2016). If the outer component can evolve into a WR due to wind-stripping then it can be observed as a component in a wide binary system. Alternatively, the merger of the inner binary followed by strong mass loss can also lead to the formation of a WR star, leading to multiple possibilities. These scenarios invoke the formation of a WN star by wind-stripping, in contrast to companion-stripping. Lastly, the expansion of the inner system might lead to dynamical interaction with the tertiary companion, potentially resulting in the expulsion of one of the stars in the system.

In conclusion, these results exhibit that there is a dire need to increase the sample by expanding RV surveys to the southern sky and to fainter magnitudes. This will enable us to derive the binary fraction and the intrinsic distributions of the orbital parameters of massive stars at the end stages of their lives in a self-consistent manner. This combined with what is known for main sequence O-type stars, is the key to advancing our knowledge of massive binary evolution.

2

David-Uraz.

3

Similar results were obtained when trying to reproduce the number of stars with RV amplitude between 15 and 30 km s−1.

4

Specifically, binomial chances to draw 11 or more binaries in a sample of 12 (corresponding the retrieved intrinsic binary fraction in the sample of > 0.9) drops below 10% for parent binary fraction smaller than 0.72.

Acknowledgments

KD, TS and HS acknowledge support from the European Research Council (ERC) innovation programme under the European Union’s DLV-772225-MULTIPLES Horizon 2020 research and innovation programme. PM acknowledges support from the FWO junior postdoctoral fellowship No. 12ZY520N.

References

  1. Barbá, R., Gamen, R., Arias, J. I., et al. 2014, Rev. Mex. Astron. Astrofis. Conf. Ser., 44, 148 [Google Scholar]
  2. Bartzakos, P. 1999, Ph.D. Thesis, Université de Montreal (Canada) [Google Scholar]
  3. Chené, A. N., & St-Louis, N. 2011, ApJ, 736, 140 [NASA ADS] [CrossRef] [Google Scholar]
  4. Chené, A. N., St-Louis, N., Moffat, A. F. J., et al. 2019, MNRAS, 484, 5834 [NASA ADS] [CrossRef] [Google Scholar]
  5. Conti, P. 1976, Proc. 20th Colloq. Int. Astrophys., Univ. Liege, 193 [Google Scholar]
  6. Crowther, P. A. 2007, ARA&A, 45, 177 [NASA ADS] [CrossRef] [Google Scholar]
  7. David-Uraz, A., Moffat, A. F. J., Chené, A.-N., et al. 2012, MNRAS, 426, 1720 [NASA ADS] [CrossRef] [Google Scholar]
  8. de Koter, A., Heap, S. R., & Hubeny, I. 1997, ApJ, 477, 792 [Google Scholar]
  9. de La Chevrotière, A., Moffat, A. F. J., & Chené, A. N. 2011, MNRAS, 411, 635 [NASA ADS] [CrossRef] [Google Scholar]
  10. de Mink, S. E., Sana, H., Langer, N., Izzard, R. G., & Schneider, F. R. N. 2014, ApJ, 782, 7 [NASA ADS] [CrossRef] [Google Scholar]
  11. Desforges, S., St-Louis, N., Chené, A. N., de la Chevrotière, A., & Hénault-Brunet, V. 2017, MNRAS, 465, 1227 [NASA ADS] [CrossRef] [Google Scholar]
  12. Dougherty, S. M., Williams, P. M., van der Hucht, K. A., Bode, M. F., & Davis, R. J. 1996, MNRAS, 280, 963 [NASA ADS] [CrossRef] [Google Scholar]
  13. Dougherty, S. M., Williams, P. M., & Pollacco, D. L. 2000, MNRAS, 316, 143 [NASA ADS] [CrossRef] [Google Scholar]
  14. Dunstall, P. R., Dufton, P. L., Sana, H., et al. 2015, A&A, 580, A93 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Eldridge, J. J., Stanway, E. R., Xiao, L., et al. 2017, PASA, 34, e058 [NASA ADS] [CrossRef] [Google Scholar]
  16. Eldridge, J. J., Xiao, L., Stanway, E. R., Rodrigues, N., & Guo, N. Y. 2018, PASA, 35, 49 [NASA ADS] [CrossRef] [Google Scholar]
  17. Eldridge, J. J., Stanway, E. R., & Tang, P. N. 2019, MNRAS, 482, 870 [NASA ADS] [CrossRef] [Google Scholar]
  18. Fahed, R., & Moffat, A. F. J. 2012, MNRAS, 424, 1601 [NASA ADS] [CrossRef] [Google Scholar]
  19. Fahed, R., Moffat, A. F. J., & Bonanos, A. Z. 2009, MNRAS, 392, 376 [CrossRef] [Google Scholar]
  20. Fahed, R., Moffat, A. F. J., Zorec, J., et al. 2011, MNRAS, 418, 2 [NASA ADS] [CrossRef] [Google Scholar]
  21. Fitzpatrick, E. L. 2004, in Interstellar Extinction in the Milky Way Galaxy, ASP Conf. Ser., 309, 33 [NASA ADS] [Google Scholar]
  22. Foellmi, C., Moffat, A. F. J., & Guerrero, M. A. 2003a, MNRAS, 338, 360 [NASA ADS] [CrossRef] [Google Scholar]
  23. Foellmi, C., Moffat, A. F. J., & Guerrero, M. A. 2003b, MNRAS, 338, 1025 [NASA ADS] [CrossRef] [Google Scholar]
  24. Fryer, C. L., Mazzali, P. A., Prochaska, J., et al. 2007, PASP, 119, 1211 [NASA ADS] [CrossRef] [Google Scholar]
  25. Gräfener, G., & Hamann, W. R. 2005, A&A, 432, 633 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Gräfener, G., Koesterke, L., & Hamann, W. R. 2002, A&A, 387, 244 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Güdel, M., & Nazé, Y. 2009, A&ARv, 17, 309 [NASA ADS] [CrossRef] [Google Scholar]
  28. Hamann, W. R., & Gräfener, G. 2003, A&A, 410, 993 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Hamann, W. R., Gräfener, G., & Liermann, A. 2006, in The Galactic WN stars: Line-blanketed Analyses versus Evolutionary Models, ASP Conf. Ser., 353 [Google Scholar]
  30. Heger, A., Fryer, C. L., Woosley, S. E., Langer, N., & Hartmann, D. H. 2003, ApJ, 591, 288 [NASA ADS] [CrossRef] [Google Scholar]
  31. Hill, G. M., Moffat, A. F. J., St-Louis, N., & Bartzakos, P. 2000, MNRAS, 318, 402 [NASA ADS] [CrossRef] [Google Scholar]
  32. Hill, G. M., Moffat, A. F. J., & St-Louis, N. 2018, MNRAS, 474, 2987 [NASA ADS] [CrossRef] [Google Scholar]
  33. Kausch, W., Noll, S., Smette, A., et al. 2015, A&A, 576, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Kobulnicky, H. A., Kiminki, D. C., Lundquist, M. J., et al. 2014, ApJS, 213, 34 [NASA ADS] [CrossRef] [Google Scholar]
  35. Langer, N., Schürmann, C., Stoll, K., et al. 2019, A&A, 638, A39 [Google Scholar]
  36. Lefèvre, L., Marchenko, S. V., Lépine, S., et al. 2005, MNRAS, 360, 141 [NASA ADS] [CrossRef] [Google Scholar]
  37. Lépine, S., Moffat, A. F. J., St-Louis, N., et al. 2000, AJ, 120, 3201 [NASA ADS] [CrossRef] [Google Scholar]
  38. Lührs, S., & Moffat, A. F. J. 2002, in Interacting Winds from Massive Stars, eds. A. F. J. Moffat, & N. St-Louis, ASP Conf. Ser., 260, 555 [Google Scholar]
  39. Mac Low, M.-M., Balsara, D. S., Kim, J., & de Avillez, M. A. 2005, ApJ, 626, 864 [NASA ADS] [CrossRef] [Google Scholar]
  40. Maíz Apellániz, J., Sota, A., Arias, J. I., et al. 2016, ApJS, 224, 4 [NASA ADS] [CrossRef] [Google Scholar]
  41. Maíz Apellániz, J., Trigueros Páez, E., Negueruela, I., et al. 2019, A&A, 626, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Marchant, P., Langer, N., Podsiadlowski, P., Tauris, T. M., & Moriya, T. J. 2016, A&A, 588, A50 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Marchenko, S. V., Moffat, A. F. J., & Grosdidier, Y. 1999, ApJ, 522, 433 [NASA ADS] [CrossRef] [Google Scholar]
  44. Mason, B. D., Hartkopf, W. I., Gies, D. R., Henry, T. J., & Helsel, J. W. 2009, AJ, 137, 3358 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  45. Massey, P. 1981, ApJ, 246, 153 [NASA ADS] [CrossRef] [Google Scholar]
  46. Massey, P. 1984, ApJ, 281, 789 [NASA ADS] [CrossRef] [Google Scholar]
  47. Massey, P., & Niemela, V. S. 1981, ApJ, 245, 195 [NASA ADS] [CrossRef] [Google Scholar]
  48. Meynet, G., & Maeder, A. 2003, A&A, 404, 975 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Moffat, A. F. J., & Shara, M. M. 1986, AJ, 92, 952 [NASA ADS] [CrossRef] [Google Scholar]
  50. Monnier, J. D., Zhao, M., Pedretti, E., et al. 2011, ApJ, 742, L1 [NASA ADS] [CrossRef] [Google Scholar]
  51. Morel, T., Georgiev, L. N., Grosdidier, Y., et al. 1999, A&A, 349, 457 [Google Scholar]
  52. Neugent, K. F., & Massey, P. 2014, ApJ, 789, 10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Niemela, V. S., Morrell, N. I., Barba, R. H., & Bosch, G. L. 1996, in Revista Mexicana de Astronomia y Astrofisica Conference Series, eds. V. Niemela, N. Morrell, P. Pismis, & S. Torres-Peimbert, Rev. Mex. Astron. Astrofis. Conf. Ser., 5, 100 [NASA ADS] [Google Scholar]
  54. Niemela, V. S., Shara, M. M., Wallace, D. J., Zurek, D. R., & Moffat, A. F. J. 1998, AJ, 115, 2047 [NASA ADS] [CrossRef] [Google Scholar]
  55. Niemela, V. S., Gamen, R., Morrell, N. I., & Giménez Benítez, S. 1999, in Wolf-Rayet Phenomena in Massive Stars and Starburst Galaxies, eds. K. A. van der Hucht, G. Koenigsberger, & P. R. J. Eenens, IAU Symp., 193, 26 [Google Scholar]
  56. Paczyński, B. 1967, Acta Astron., 17, 355 [NASA ADS] [Google Scholar]
  57. Patrick, L. R., Lennon, D. J., Britavskiy, N., et al. 2019, A&A, 624, A129 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. Pols, O. R. 1994, A&A, 290, 119 [Google Scholar]
  59. Raskin, G., van Winckel, H., Hensberge, H., et al. 2011, A&A, 526, A69 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Richardson, N. D., Shenar, T., Roy-Loubier, O., et al. 2016, MNRAS, 461, 4115 [CrossRef] [Google Scholar]
  61. Rosslowe, C. K., & Crowther, P. A. 2015, MNRAS, 447, 2322 [NASA ADS] [CrossRef] [Google Scholar]
  62. Rustamov, D. N., & Cherepashchuk, A. M. 1989, Sov. Ast., 33, 36 [NASA ADS] [Google Scholar]
  63. Sana, H., de Mink, S. E., de Koter, A., et al. 2012, Science, 337, 444 [Google Scholar]
  64. Sana, H., de Koter, A., de Mink, S. E., et al. 2013, A&A, 550, A107 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Sana, H., Le Bouquin, J. B., Lacour, S., et al. 2014, ApJ, 215, 15 [NASA ADS] [Google Scholar]
  66. Sander, A., Hamann, W. R., & Todt, H. 2012, A&A, 540, A144 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Sander, A. A. C., Hamann, W. R., Todt, H., et al. 2019, A&A, 621, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. Schnurr, O., Moffat, A. F. J., St-Louis, N., Morrell, N. I., & Guerrero, M. A. 2008, MNRAS, 389, 806 [NASA ADS] [CrossRef] [Google Scholar]
  69. Shenar, T., Hainich, R., Todt, H., et al. 2016, A&A, 591, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Shenar, T., Richardson, N. D., Sablowski, D. P., et al. 2017, A&A, 598, A85 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  71. Shenar, T., Sablowski, D. P., Hainich, R., et al. 2019, A&A, 627, A151 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Shenar, T., Gilkis, A., Vink, J. S., Sana, H., & Sander, A. A. C. 2020, A&A, 634, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. Skinner, S., Güdel, M., Schmutz, W., & Zhekov, S. 2006, Astrophys. Space Sci., 304, 97 [NASA ADS] [CrossRef] [Google Scholar]
  74. Smette, A., Sana, H., Noll, S., et al. 2015, A&A, 576, A77 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Smith, L. F. 1968, MNRAS, 140, 409 [NASA ADS] [Google Scholar]
  76. St-Louis, N., Chené, A. N., Schnurr, O., & Nicol, M. H. 2009, ApJ, 698, 1951 [NASA ADS] [CrossRef] [Google Scholar]
  77. Thomas, J. D., Richardson, N. D., Sana, H. S., et al. 2020, MNRAS, submitted [Google Scholar]
  78. Tody, D. 1993, in IRAF in the Nineties, ASP Conf. Ser., 52, 173 [NASA ADS] [Google Scholar]
  79. Tonry, J., & Davis, M. 1979, AJ, 84, 1511 [NASA ADS] [CrossRef] [Google Scholar]
  80. Tramper, F., Sana, H., Fitzsimons, N. E., et al. 2016, MNRAS, 455, 1275 [NASA ADS] [CrossRef] [Google Scholar]
  81. Vanbeveren, D., De Donder, E., Van Bever, J., Van Rensbergen, W., & De Loore, C. 1998, New A, 3, 443 [NASA ADS] [CrossRef] [Google Scholar]
  82. van der Hucht, K. A. 2001, New Astron. Rev., 45, 135 [NASA ADS] [CrossRef] [Google Scholar]
  83. Varricatt, W. P., & Ashok, N. M. 2006, MNRAS, 365, 127 [NASA ADS] [CrossRef] [Google Scholar]
  84. Williams, P. M. 2014, MNRAS, 445, 1253 [NASA ADS] [CrossRef] [Google Scholar]
  85. Williams, P. M., Kidger, M. R., van der Hucht, K. A., et al. 2001, MNRAS, 324, 156 [NASA ADS] [CrossRef] [Google Scholar]
  86. Williams, P. M., van der Hucht, K. A., & Rauw, G. 2005, ArXiv e-prints [arXiv: astro-ph/0508157] [Google Scholar]
  87. Woosley, S. E. 2019, ApJ, 878, 49 [NASA ADS] [CrossRef] [Google Scholar]
  88. Zapartas, E., de Mink, S. E., Izzard, R. G., et al. 2017, A&A, 601, A29 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  89. Zhekov, S. A. 2017, MNRAS, 472, 4374 [NASA ADS] [CrossRef] [Google Scholar]
  90. Zucker, S. 2003, MNRAS, 342, 1291 [Google Scholar]

Appendix A: Comments on specific objects

WR 4. According to the seventh catalogue of WR stars (van der Hucht 2001), WR 4 is classified as a WC5+?, an SB1 with no observable line dilution or absorption effects. Photometric variability with a period of 4.3 days (Moffat & Shara 1986) was determined. Rustamov & Cherepashchuk (1989) found both a photometric and spectroscopic period of 2.4096 days with a peak-to-peak amplitude of ∼60 km s−1, but we can safely reject the presence of such large amplitude RV variation in our data set (Table 2) questioning the fact that the 2.4 d period traces an orbital motion. We therefore do not include WR 4 in the list of short-period binaries in our literature study (Sect. 4.2).

In a study of co-rotating interaction regions (CIRs), St-Louis et al. (2009) classified WR 4 as “Not Variable” in the weak C IV 5016 Å line. Chené & St-Louis (2011) further elaborated that the maximum line-profile variability was 1.2% of the level of flux in the line, hence making it ineligible to be a CIR candidate. In this work, we obtained 22 epochs over a time coverage of 890 days, out of which about 12 were observed within two nights, attempting a probe on shorter periods. Although we do not report a short period, the adopted choice of our threshold allows us to classify it as a binary (Table 2). Further monitoring of this system is required to determine the orbital period of this system.

WR 5. van der Hucht (2001) reports WR 5 as a single star with spectral type of WC6. No emission was found in the X-rays (Skinner et al. 2006; Güdel & Nazé 2009). In this work, WR 5 is found to have a Δ RV of 4.9 km s−1 over 718 days with a σRV of 1.7 km s−1 for RVs measured using the full spectrum. It does not pass our binary detection criterion and is, therefore, classified as presumably single.

WR 111. WR 111 is a presumably single WC5 star (Niemela et al. 1999; van der Hucht 2001). It has been studied with spectral analyses and hydrodynamical studies (Gräfener et al. 2002; Gräfener & Hamann 2005; Sander et al. 2012), and was used to show for the first time that WC winds can be explained through radiation-driven mass loss. In this work, WR 111 shows a peak-to-peak variability of 4.5 km s−1 over 742 days with a σRV of 1.5 km s−1. It does not pass our binary detection criterion and is therefore classified as a presumably single star.

WR 113. WR 113 is a WC8d+O8-9 binary system first identified by Massey & Niemela (1981) with a period of about 29.707 days in a circular orbit. A debate about the eccentricity arose when Niemela et al. (1996) derived similar orbital parameters with a significant eccentricity (e = 0.19). The matter was resolved by David-Uraz et al. (2012), who showed that the two were indistinguishable solutions and adopted a circular orbit and was later confirmed by Hill et al. (2018) with a period of 29.700 ± 0.001 d. In this work, WR 113 has the largest Δ RV and σRV of 282.6 and 93.9 km s−1 respectively, passing our binary detection criterion. The RV measurements phased with the orbital period are shown in Fig. 7.

WR 117. According to the GCWR, WR 117 is classified as a WC9d, which is a dust producing star. Sander et al. (2012) find it to be a bit too hotter than most WC9 counterparts, even suggesting it might be in the WC8 parameter regime. Williams et al. (2005) did not find any evidence for an OB-type companion, and hypothesize that WR 117 might be an evolved WC9 star. A large scatter in photometric magnitudes was observed by Williams (2014) with possible 0.5 mag eclipses due to dust. In this work, WR 117 shows a Δ RV of 4.9 km s−1 and a σRV of 1.9 km s−1 over a coverage of 488 days and is hence classified as a presumably single star.

WR 119. WR 119 is classified as a single WC9d star according to the GCWR. Fahed et al. (2009) studied it using V and I band photometry, finding a larger than normal variability in these bands. Desforges et al. (2017) performed a temporal variance spectrum analysis over timescales of days to weeks and found significant variability in the C III 5696 Å line. In this work, we find values of Δ RV and σRV of 11.0 and 3.2 km s−1 respectively, resulting in WR 119 being classified as a binary.

WR 135. WR 135 is an apparently single WC8 star (GCWR). Previous RV studies have failed to determine any variability (Niemela et al. 1999), nor even photometric variability (Moffat & Shara 1986). In this work, we find a Δ RV of 4.7 km s−1 and σRV of 1.6 km s−1 over a baseline of 2660 days, classifying it as a presumably single star.

WR 137. WR 137 is a well-studied, long-period binary system with a late-type WC star and a late-type O star (WC7pd + O9 V; p: peculiar, d: dust producing). Marchenko et al. (1999) reported dust production in 1997 September – 1998 May and Williams et al. (2001) derived a period of 13.05 ± 0.15 years by tracking the periodic dust formation. This period was verified by the orbital analysis using RVs by Lefèvre et al. (2005), who found the maximum emission from dust to occur just after the periastron passage. The orbit was recently resolved using interferometry with CHARA (Richardson et al. 2016), reiterating the importance of long term RV monitoring campaigns in identifying WR binaries. The broad emission lines of the WR star are superimposed by a narrow emission and an absorption component, making the companion star a prime candidate to be an Oe star. The RVs measured in this campaign are in agreement with those obtained from the literature (Fig. 4).

WR 140. WR 140 is classified as a (WC7pd + O5.5fc) system with a well established spectroscopic (Fahed et al. 2011) and visual (interferometric; Monnier et al. 2011) orbit (P = 7.938 years, e = 0.8996). These parameters were derived by combining archival data along with measurements during the 2009 periastron passage. An improved derivation of the masses and orbital parameters using spectroscopic and interferometric data from the 2016 December periastron passage makes it one of the most accurately constrained masses for a WR binary system to date (Thomas et al. 2020). In this work, WR 140 has a Δ RV of 21.6 km s−1 and a σRV of 6.5 km s−1, which passes our binary detection criterion. The RV measurements along with archival data are shown in Fig. 7.

WR 143. WR 143 is classified as a binary system (WC4 + Be; Varricatt & Ashok 2006), who detected the presence of hydrogen emission lines and He I lines in the infrared spectrum. In the optical spectra, a hint of this is seen in Hα, with features similar to those exhibited by WR 137. We detect significant RV variability over our 76 d baseline and thus classify WR 143 as a binary (Δ RV = 11.2 km s−1, σRV = 4.2 km s−1).

WR 146. WR 146 is a visual binary with a spectral classification of WC5 + O8. The components show non-thermal radio emission (Dougherty et al. 1996) and were resolved with HST images (Niemela et al. 1998). The computed mass loss rate for WR 146 is much higher than prototypical WC5 stars, prompting a debate on the multiplicity of the O-star companion (Dougherty et al. 2000). An X-ray study by Zhekov (2017) shows a discrepancy in the mass-loss rate by a factor of 8 to 10. Given the angular separation, Niemela et al. (1998) suggest binary periods of hundreds of years. However, our 847 day baseline with values for Δ RV and σRV of 18.9 and 5.9 km s−1 suggest periods of the order of 3000 d (Fig. B.1) for a massive companion. Further monitoring is required to constrain the orbit.

WR 154. WR 154 is classified as a single WC6 star (GCWR). In this work, we do not detect significant RV variation. However, given the short baseline (24 days), we cannot rule out the possibility of it being a long period binary as the detection probability of our campaign drops significantly for periods above 100 d (Fig. B.1).

Appendix B: Relative RV measurements

Relative RVs for the objects in our sample, measured with the full spectrum as a mask. The reference epoch is marked with a “*”. The Barycentric Julian Date (BJD) is given as the middle of the exposure. The average S/N is given in Table 1.

thumbnail Fig. B.1.

Binary detection probability of the observational campaign for the various stars in our sample, adopting a minimum RV variability threshold of 10 km s−1.

thumbnail Fig. B.1.

continued.

Table B.1.

Journal of the HERMES observations for WR 4.

Table B.2.

Journal of the HERMES observations for WR 5.

Table B.3.

Journal of the HERMES observations for WR 111.

Table B.4.

Journal of the HERMES observations for WR 113.

Table B.5.

Journal of the HERMES observations for WR 117.

Table B.6.

Journal of the HERMES observations for WR 119.

Table B.7.

Journal of the HERMES observations for WR 135.

Table B.8.

Journal of the HERMES observations for WR 137.

Table B.9.

Journal of the HERMES observations for WR 140.

Table B.10.

Journal of the HERMES observations for WR 143.

Table B.11.

Journal of the HERMES observations for WR 146.

Table B.12.

Journal of the HERMES observations for WR 154.

All Tables

Table 1.

Sample of WC stars in the RV monitoring campaign with the number of spectra, time coverage, and average S/N per resolution element at 5200 Å.

In the text
Table 2.

Results of the RV measurements for the sample of WC stars.

In the text
Table B.1.

Journal of the HERMES observations for WR 4.

In the text
Table B.2.

Journal of the HERMES observations for WR 5.

In the text
Table B.3.

Journal of the HERMES observations for WR 111.

In the text
Table B.4.

Journal of the HERMES observations for WR 113.

In the text
Table B.5.

Journal of the HERMES observations for WR 117.

In the text
Table B.6.

Journal of the HERMES observations for WR 119.

In the text
Table B.7.

Journal of the HERMES observations for WR 135.

In the text
Table B.8.

Journal of the HERMES observations for WR 137.

In the text
Table B.9.

Journal of the HERMES observations for WR 140.

In the text
Table B.10.

Journal of the HERMES observations for WR 143.

In the text
Table B.11.

Journal of the HERMES observations for WR 146.

In the text
Table B.12.

Journal of the HERMES observations for WR 154.

In the text

All Figures

thumbnail Fig. 1.

Spectral type distribution of the stars in our sample.
In the text
thumbnail Fig. 2.

Main steps in our data reduction process. From top to bottom: (a) flatfield corrected spectrum produced by the HERMES pipeline. (b) Spectrum corrected for the instrumental response function and telluric contamination. (c) Continuum flux models from PoWR to the pseudo-continuum region around 5200 Å for different values of E(B − V). (d) Resulting normalised spectrum.
In the text
thumbnail Fig. 3.

Spectrum of WR 137, rebinned for clarity, with different masks of C IV indicated (see legend). This process is repeated for all other ions, producing a plethora of masks.
In the text
thumbnail Fig. 4.

Phased RV plot for WR 137. Archival RV measurements are shown in blue and are taken from Lefèvre et al. (2005). The RV measurements from HERMES using the strong C IV lines are shown in red. Along with these measurements, the inset shows the short cadence measurements with the full spectrum (green) and weak He II lines (purple). All HERMES RV measurements have been shifted by −23.5 km s−1 in order to convert relative RV measurements to absolute values.
In the text
thumbnail Fig. 5.

RV measurements for WR 137 during the short cadence study. At the bottom of the plot are the values of (σw, σRV) (km s−1) for each mask. It is observed that the helium lines (and weak lines in general) have large dispersions. The masks least affected by wind variability are “full spec”, “C4 strong”, “C4 moderate”, and “C4 all” (determined by their values of σw).
In the text
thumbnail Fig. 6.

Observed binary fraction (fobs) as a function of the threshold C. The adopted threshold is the vertical dotted-dashed line at 10 km s−1, resulting in fobs = 0.58. The objects are labelled as and when they are classified as binaries based on Δ RV. Names in bold are known binaries.
In the text
thumbnail Fig. 7.

Top: phased RV plot for WR 140. Archival RV measurements are shown in blue (Fahed et al. 2011) and green (Thomas et al. 2020). The RV measurements from HERMES (red) using C IV strong lines have been shifted by 19 km s−1 in order to convert relative RV measurements to absolute ones. Bottom: same for WR 113. Archival RV measurements are shown in blue (Hill et al. 2018) and green (Massey & Niemela 1981). The RV measurements from HERMES (red and yellow) have been shifted by −135 km s−1 in order to convert relative RV measurements to absolute values.
In the text
thumbnail Fig. 8.

Detection probability for WR 137 as a function of the orbital period and secondary mass. A representative primary mass 10 M⊙ has been adopted for the WR star.
In the text
thumbnail Fig. 9.

Average detection probability curve of our campaign as a function of the orbital period and assuming a flat mass-ratio distribution between 0.5 and 2.0. Known orbital periods of WC stars in our sample are indicated by vertical lines.
In the text
thumbnail Fig. B.1.

Binary detection probability of the observational campaign for the various stars in our sample, adopting a minimum RV variability threshold of 10 km s−1.
In the text
thumbnail Fig. B.1.

continued.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.