© ESO, 2017

1. Introduction

Binary stars play a fundamental role in astrophysics since they allow direct measurements of masses, radii, and luminosities that put constraints on stellar physics, Galactic archaeology, high-energy physics, etc. Binary systems are found at all evolutionary stages, and after strong interaction, some may end up as double degenerate systems or merged compact objects.

Spectroscopic binaries (SBs) exist in different flavours. On the one hand, SB1 (SBs with one observable spectrum) can only be detected from the Doppler shift of the stellar spectral lines. On the other hand, SBn (n ≥ 2) are characterized by a composite spectrum made out of n stellar components, and are detected either from the composite nature of the spectrum or from the Doppler shift of the spectral lines. Spectroscopic binaries are certainly the binaries that cover the widest range of masses (from brown dwarfs to massive twins) and all ranges of periods (from hours to hundreds of years as observed so far, e.g. Pourbaix 2000). To date, more than 3500 SBs with orbital elements have been catalogued, of which about 1126 are SB2 (Pourbaix et al. 2004, and the latest online version of the SB9 catalogue). The Geneva-Copenhagen Survey catalogue (Nordström et al. 2004; Holmberg et al. 2009) contains approximately 4000 SB1, 2100 SB2, and 60 SB3 out of 16700 F and G dwarf stars in the solar neighbourhood, most without orbits. In the vast majority of cases, these binaries have not yet been confirmed but correspond to an overall binary fraction in the Milky Way of almost 40%. A census of binary fraction is also available from the Hipparcos catalogue (Frankowski et al. 2007), though the binary fraction per spectral type is probably affected by selection biases in the Hipparcos entry catalogue. New recent Galactic surveys like APOGEE (Majewski et al. 2017) or LAMOST (Luo et al. 2015) allow new investigations of binarity over a large sample of stars (see e.g. Gao et al. 2014; Troup et al. 2016; Fernandez et al. 2017). For instance, the RAVE survey led to the detection of 123 SB2 candidates out of 26 000 objects (Matijevič et al. 2010, 2011). We refer the reader to Duchêne & Kraus (2013) for a recent review of the physical properties of multiplicity among stars and more specifically to Raghavan et al. (2010) for a complete volume-limited sample of solar-type stars in the solar neighbourhood (distances closer than 25 pc).

The Gaia-ESO Survey (GES) is an ongoing ground-based high-resolution spectroscopic survey of 10⁵ stellar sources (Gilmore et al. 2012; Randich et al. 2013) covering the main stellar populations of the Galaxy (bulge, halo, thin and thick discs) as well as a large number of open clusters spanning large metallicity and age ranges. All evolutionary stages are encountered within the GES, from pre-main sequence objects to red giants. It aims to complement the spectroscopy of the Gaia ESA space mission (Wilkinson et al. 2005). The GES uses the FLAMES multi-fibre back end at the high-resolution UVES (R ~ 50 000) and moderate-resolution GIRAFFE (R ~ 20 000) spectrographs. The visual magnitude of the faintest targets reaches V ~ 20. The spectral coverage spans the optical wavelengths (from 4030 to 6950 Å) and the near-infrared around the Ca ii triplet and the Paschen lines (from 8490 to 8900 Å including the wavelength range of the Radial Velocity Spectrometer of the Gaia mission). The median signal-to-noise ratio (S/N) per pixel is similar for UVES and GIRAFFE single exposures (~30), whereas the most frequent values are around 20 and 5, respectively.

The motivation of the present work is to take advantage of a very large sample to automatically detect SBs with more than one visible component¹ that are not always detected by the GES single-star main analysis pipelines. Spectroscopic binaries may be a potential source of error when deriving atmospheric parameters and detailed abundances. This project presents a new method to automatically identify the number of velocity components in each cross-correlation function (CCF) using their successive derivatives and the analysis of about 51 000 stars available within the GES internal data release 4 (iDR4).

In Sect. 2, we describe the iDR4 stellar observations, their associated CCFs, and the selection criteria applied to them. The method used to detect the velocity components in a CCF, its parameters, and the formal uncertainty are presented in Sect. 3. In Sect. 4, the set of SBn (n ≥ 2) detected in iDR4 using this method is discussed, organized according to the stellar populations they belong to.

2. Data selection

2.1. Observations and CCF computation

Our analysis was performed on the iDR4 consisting of ~260 000 single exposures (corresponding to ~100 000 stacked spectra) of about 51 000 distinct stars observed with the FLAMES instrument feeding the optical spectrographs GIRAFFE (with setups HR3, HR5A, HR6, HR9B, HR10, HR14A, HR15N, HR15, HR21) and UVES (with setups U520 and U580) covering the optical and near-IR wavelength ranges given in Table 1.

The classical definition of a CCF function applied to the stellar spectra is (1)where f is a normalized spectrum, g a normalized template spectrum, and h is the lag expressed in km s^-1. The computation of the CCFs is performed by pipelines at CASU (Cambridge Astronomy Survey Unit²) for GIRAFFE spectra (Lewis et al., in prep.) and at INAF-Arcetri for UVES spectra (Sacco et al. 2014). For UVES CCFs, spectral templates from the library produced by de Laverny et al. (2012), and based on MARCS models (Gustafsson et al. 2008), are used. For GIRAFFE CCFs, spectral templates from the library produced by Munari et al. (2005), and based on Kurucz’s models (Kurucz 1993; Castelli & Kurucz 2003), are used. We stress that for a given spectrum, CCFs are calculated for all the templates and the CCF with the highest peak is selected. For UVES spectra, Hα and Hβ are masked in the observations. As illustrated in Fig. 1, CCFs are characterized by a maximum value (CCF peak), a minimum value (lowest point of the CCF tail) and a full amplitude (maximum – minimum). The constant velocity steps of GIRAFFE and UVES CCFs are 2.75 (mainly) and 0.50 km s^-1 (for a sampling of 401 and 4000 velocity points), respectively.

Examples of spectra and CCFs in the setups mentioned above are displayed in Figs. 2 and 3. These figures are built from the solar and Aldebaran spectra. The CCFs are represented over the same velocity range to allow an easy comparison between the various setups. When a lot of weak absorption lines are present (as in setups HR6 and HR10), the CCF peak is narrow and well-defined with a width smaller than for setups with strong features like Hδ (HR3), Hγ (HR5A), the Mg b triplet (HR9B), Hα (HR14A and HR15N), and the Ca II triplet (HR21). For HR15, the presence of telluric lines from 685 nm onwards reduces the maximum amplitude of the CCF to a value as low as 0.25, even with a S/N higher than 1000.

For the UVES setups, Aldebaran (α Tau, spectral type K5III) spectra and corresponding CCFs are presented in Fig. 3. Each setup is composed of two spectral chunks. In the present case, the lower chunk comes with S/N ≈ 70 and the upper one with S/N> 100. For the setup U520 low, the leftward CCF tail is negative, probably as a result of poor spectrum normalization due to the co-existence of many weak and strong lines. Since the wavelength range of the UVES setups is 2 or 3 times wider than those of GIRAFFE, the UVES setups are well suited to detecting SBn candidates.

The final GES spectrum of a given object is a stack of all individual exposures, wavelength calibrated, sky subtracted, and heliocentric radial velocity corrected. This could be a source of confusion in the case of composite spectra where the radial velocity of the different components changes between exposures. Moreover, a double-lined CCF coming from stacked spectra (and mimicking an SB2) can be the result of the SB1 combination taken at different epochs and stacked. To avoid this problem, we performed the binarity detection on the individual exposures (rather than on the stacked ones). This choice avoids spurious spectroscopic binary detection at the expense of using spectra with lower S/N, which will be shown not to be detrimental as long as S/N> 5 (see Sect. 3.4).

The number of individual observations per target is plotted in Fig. 4. Most of the stars observed with GIRAFFE have 2 or 4 observations because generally observed with HR10 and HR21 setups, whereas there are 4 or 8 observations in the case of UVES, due to the presence of two spectral chunks per setup. Moreover, the time span between consecutive observations is very often less than three days, as shown in Fig. 5. Benchmark stars (i.e. a sample of stars with well-determined parameters to be used as reference; see Heiter et al. 2015a) are the most observed objects, some having more than 100 observations.

2.2. Data selection in iDR4

Our sample has been drawn from the individual spectra database of the GES iDR4³, covering observations until June 2014, to which the following selection criteria were applied:

• S/N higher than 5;

• CCF maximum larger than 0.15;

• CCF minimum larger than −1;

• CCF full amplitude larger than 0.10;

• left CCF continuum − right CCF continuum smaller than 0.15.

These criteria were empirically determined thanks to a visual inspection of a representative sample of CCFs. We allow negative values for the CCF minimum to keep CCFs computed on unperfectly normalized spectra (without allowing spectra with a completely incorrect normalization). Criteria on the S/N and on the CCF maximum are presented in Fig. 6 for setups HR10 and HR21 which contain the most numerous observations. This figure clearly shows the impact of the S/N of a spectrum on its associated CCF: the higher the S/N, the higher the CCF maximum. For a given S/N, the interval spanned by the CCF maximum is mainly due to spectrum – template mismatch. For HR10, the over-density located at 30 < S/N < 200 and CCF max < 0.15 is mainly due to NGC 6705 members. In HR21, the clump located at 1000 < S/N < 2000 and 0.80 < CCF max < 0.85 is due to repeated observations of the solar spectrum.

These criteria allow us to avoid detecting spurious (noise-induced) CCF peaks. Over the 260 000 individual science spectra (corresponding to the 100 000 stacked spectra) within the iDR4, 9.3% have a S/N lower than 5, 1.0% have a null CCF (data processing issues), 7.8% have a CCF maximum lower than 0.15, 0.2% have a CCF minimum lower than −1.0, and 0.02% have a CCF full amplitude lower than 0.10. We ended up with about 205 000 CCFs (77.7%), corresponding to ~51 000 different stars.

3. Methods

3.1. Detection of extrema (DOE) code

The detection of extrema (doe) code has been designed to identify the (local and global) extrema in a given signal even when these extrema are strongly blended. By using successive derivatives of a function, it is possible to characterize it in a powerful way. Applied on spectral-line profiles for instance, the method makes it possible to identify all contributing blends (Sousa et al. 2007). Here we apply it to the CCFs. The method is inspired from signal-processing techniques (Foster 2013) which convolve the signal (here the CCF) with the derivatives of a Gaussian kernel to smooth and calculate the derivative of the CCF in a single operation. In other words, the first, second, and third derivatives of the Gaussian kernel are used to obtain the smoothed derivatives of the CCFs. Indeed, one of the interesting properties of the convolution of two generalized functions is defined as: (2)where f′ and g′ are the first derivatives of the generalized functions f and g. Convolving the CCF with the derivative of a Gaussian kernel is equivalent to computing the derivative of the CCF and to convolving (i.e. smoothing) it by a Gaussian kernel. We use the routine gaussian_filter1d of the sub-module ndimage of the scipy module (Jones et al. 2001) in Python. The routine first calculates the derivative of the Gaussian kernel before correlating it with the CCF function. The width of the Gaussian kernel controls the amount of smoothing.

A zero in the descending part of the first derivative obviously provides the position of the maximum of the CCF. However, in the case of a CCF composed of two or more peaks, the zeros of the first derivative will only provide the positions of well-separated peaks, i.e. peaks with a local minimum between them. Blended peaks might thus be missed. However, this difficulty can be circumvented by using the third derivative, whose zeros occurring in an ascending part provide the positions of all the peaks including the blended ones. Figure 7 shows that the use of the first derivative only does not allow a satisfactory detection of the CCF components. Indeed, although the CCF in the middle panel clearly exhibits two peaks, the first derivative has only one descending zero-crossing, thus resulting in the detection of only one component. However, the second derivative shows two local minima corresponding to the two CCF velocity components. The position of these two minima can be found by detecting the ascending zero-crossing of the third derivative. By using the third derivative, the different CCF components can thus be identified as regions where the CCF curvature is sufficiently negative (minima of the second derivative, or ascending zeros of the third derivative), separated by a region of larger curvature. To get the velocities of the various components, the CCF third derivative is simply interpolated to find its intersection with the x-axis. Some detection thresholds had to be set to automate the process in order to match the results obtained from a visual inspection of multiple-component CCFs.

The procedure is illustrated on simulated CCFs with one and two peaks (Figs. 1 and 7, respectively). We first test the operation of the doe code on single peaks at the lowest numerical resolution, i.e. peaks defined with only six velocity points (left panel of Fig. 1). The doe code applied on a more realistic (more noisy) simulated single-peaked CCF (as shown in the right panel of Fig. 1) also provides satisfactory results, with an accuracy on the radial velocity of the order of 0.20 km s^-1. We will show in Sect. 3.4 that the doe code has a small internal error of 0.25 km s^-1.

The first threshold (THRES0), expressed as a fraction of the full CCF amplitude, defines the considered velocity range: the doe code is applied only in the region where the CCF is higher than THRES0. The THRES0 threshold is represented by the horizontal red line in the top panels of Fig. 1 and subsequent figures. However, if several well-defined peaks are identified in the CCF, the THRES0 criterion is overridden, and all data points between the CCF peaks are included in the analysis of the derivatives, even though the CCF may be lower than THRES0.

A second threshold, THRES2, is set on the second CCF derivative. The THRES2 parameter is expressed as a fraction of the full amplitude of the CCF second derivative. This negative threshold is represented by the horizontal red line in the “2nd derivative” panel in Fig. 1 (and subsequent figures) such that only minima lower than this threshold are selected for the final peak detection (vertical black lines), whereas second-derivative minima higher than this threshold are not considered to be related to real components (vertical light grey lines in e.g. Fig. 9).

The width of the Gaussian kernel for the convolution of the CCF, SIGMA, is the third parameter. It is a smoothing parameter and aims to make the successive derivatives of the CCF less sensitive to the data noise.

The three parameters of doe (THRES0, THRES2, and SIGMA) have to be set by the user. Their value may have an impact on the number of detected peaks and the radial velocities associated with them. These three parameters need to be adjusted in order to give meaningful results (i.e. matching the efficiency of a by-eye detection) on all CCFs, but once fixed for each instrumental setup (see Table 1 and Sect. 3.3), they are kept constant to ensure homogeneous detection efficiency over the whole GES sample.

The parameter values result from a compromise between antagonistic requirements:

• the THRES0 parameter must not be too low to avoid anunrealistically large velocity range, or too high in order to be ableto detect real, albeit low, secondary peaks;

• the THRES2 parameter must be calibrated on extreme cases (two very close or very separated peaks). The choice of this parameter is important: it ensures that the second derivative (i.e. the curvature) of the CCF is negative enough, therefore corresponding to real components;

• the SIGMA parameter must not be too large, which would result in excessive smoothing and endanger the detection of close peaks, and not too small to reduce the impact of the numerical noise induced by the successive derivatives.

The empirical method used to set these parameters is described in Sect. 3.3.

3.2. Detection of peaks on simulated CCFs

We tested the efficiency of the doe code on simulated double-peak CCFs. Using the radiative transfer code turbospectrum (Plez 2012; de Laverny et al. 2012), the MARCS library of model atmospheres with spherical geometry (Gustafsson et al. 2008), and the GES atomic linelist (Heiter et al. 2015b), we computed the synthetic spectrum of a star with the following stellar parameters: T_eff = 5000 K, log g = 1.5, [Fe/H] = 0.0, and ξ_t = 1.5 km s^-1, between 5330 Å and 5610 Å for a resolution of R ~ 20 000, i.e. to reproduce an HR10 spectrum (see Sect. 2.1). Then, we shifted this spectrum so that the radial velocity of this simulated star is v_rad,0 = 72 km s^-1.

We also add Gaussian noise to reproduce spectra with S/N = 20. Then we combine the spectra shifted at different radial velocities to simulate a composite spectrum. Assuming a flux ratio between the two components of 2/3, we set the main peak at a fixed velocity of 72.0 km s^-1, whereas the position of the second peak is set at 36.0, 48.0, or 54.0 km s^-1. The cross-correlation function between the composite and the initial spectrum is calculated and then normalized by the maximum value of the mask auto-correlation (auto-correlation of the initial spectrum).

The three simulated CCFs and their derivatives are shown in Fig. 7, the value of SIGMA being 2.1 km s^-1. From the first derivative, only one crossing of the x-axis leads to the detection of one single peak. From the second derivative, we see clearly two minima in the left and middle panels, whereas we see only one minimum in the right panel. This leads to the conclusion that the detection limit between two components is 18 km s^-1. This detection limit depends on the typical width of absorption lines in the tested spectrum and also depends on the SIGMA parameter. However, reducing the SIGMA parameter too much could increase false peak detections for bumpy CCFs. The compromise adopted is described in Sect. 3.3.

3.3. Choice of the DOE parameters for the different setups

The three parameters of the doe code described in Sect. 3.1 have to be adjusted to optimize the detection of the CCF components. These parameters were adjusted by performing individual calibrations for the different setups (GIRAFFE HR10, HR15N, HR21, and UVES U520 and U580) using examples of single-, double-, and triple-peak CCFs with different separations between the components, and different component widths (i.e. different degrees of blending). For the remaining GIRAFFE setups, a standard value of the SIGMA parameter (3 km s^-1) was adopted. The adopted values are listed in Table 1. The parameter adjustment aims to obtain the same detection efficiency on the test CCFs as through visual inspection, especially in the extreme cases (blended CCFs). Figures 8 and 9 illustrate favourable and extreme cases. The value of THRES0 is higher for the GIRAFFE CCFs than for the UVES data because the correlation noise (i.e. the signal level in the CCF continuum) was observed to be higher in GIRAFFE CCFs.

Depending on the setup resolution and the number and strength of lines, the minimum separation for peak detection was empirically found to be in the range [20−60] km s^-1 for GIRAFFE setups (15 km s^-1 for UVES). As an example, in Sect. 3.2 and Fig. 7, we showed with simulated CCFs that the detection limit is reached for a minimum separation of 18 km s^-1 at R ~ 20 000 for slowly rotating stars. The spectrograph resolution and the CCF sampling are not the only relevant parameters here since the intrinsic line broadening (macroturbulence and stellar rotation) also has an impact on the CCF width.

The doe code includes a procedure to compare the number of valleys in the second derivative with the number of detected peaks. When these numbers are not identical, iteration on the detection occurs after increasing the SIGMA parameter. This procedure prevents false detections because in these situations the wide CCF often exhibits inflexion points which cause zeros in the third derivative (see left panels of Fig. 10). The number of valleys, defined as regions where the second derivative is continuously negative, is assessed first. For example, in the left “second derivative” panel of Fig. 10, one valley is detected. For low values of the SIGMA parameter, the number of detected velocity components is systematically higher than the number of valleys (left panels of Fig. 10). As long as the number of valleys is lower than the number of velocity components detected from the third derivative, the SIGMA parameter is increased by 2 km s^-1 until the number of detected velocity components equals the number of valleys. The iterative process is then stopped and the radial velocities of the detected velocity components are identified.

Figure 10 shows an example of this procedure applied on the K1 pre-main sequence object 2MASS J06411542+0946396 (CNAME⁴ 06411542+0946396) member of the cluster NGC 2264 (Fűrész et al. 2006). The doe run starts with the standard SIGMA value of 5 km s^-1. Initially, the doe code detects three valleys in the second derivative and six velocity components from the third derivative, which are clearly spurious detections. After three iterations, one valley and three velocity components are identified (left panel of Fig. 10). After 11 iterations, SIGMA increases from 5 to 27 km s^-1 and the process results in one velocity component located at 18.31 km s^-1 (right panel of Fig. 10, compared with the velocity of 19.86 km s^-1 found by Fűrész et al. 2006). The case of CCF multiplicity that can be due to physical processes different from binarity (pulsating stars, nebular lines in spectra, etc.) is discussed in Sect. 4.7.

3.4. Estimation of the formal uncertainty of the method

In this section, we assess the choice of the SIGMA parameter and its effect on the derived radial velocities and their uncertainty. The uncertainty on the derived radial velocity for single-peak CCF depends mainly on the S/N of the spectrum used to compute the CCF, the normalization of this spectrum, and the mismatch between the spectrum and the mask (spectral type, atomic and molecular profiles, rotational velocity, etc.).

We performed Monte Carlo simulations to compute single-peak CCFs from spectra of different S/N but using the same atmospheric parameters defined in Sect. 3.2. We sliced this synthetic spectrum and degraded its resolution in order to match the following settings: GIRAFFE HR10 and HR21, UVES U520 and U580 (up and low). For each S/N level, we computed 251 realizations of our simulated GIRAFFE and UVES spectra by adding Gaussian noise and computed the corresponding CCFs using a mask made of a noise-free spectrum with a null radial velocity. We finally ran doe with different values of SIGMA (from 1 to 15 by step of 1 km s^-1). Figures 11 and 12 show the difference Δv_rad = v_rad,doe−v_rad,0, where v_rad,0 = 72.0 km s^-1, as a function of the doe parameter SIGMA (right panel) and the 251 CCFs (left panel) along with the noise-free CCF (labelled “+ ∞”). We show the results for the lowest S/N (i.e. the most unfavourable cases) for the setups GIRAFFE HR10 and HR21 and UVES U580 (low and up). The mean and standard deviation of Δv_rad are also superimposed with dark dots and error bars in the right panels.

The noise-free CCF (blue curve) in the left panels of Figs. 11 and 12 show striking differences from one setup to the other. This is directly related to the spectral information contained by the spectrum used in the CCF computation. For our simulated star, the HR10 and U580 (low) spectra are more crowded than the HR21 and the U580 (up) spectra. This results in a higher level of the CCF continuum. In addition, in HR21 the large wings of the CCFs are due to the strong Ca ii IR triplet that completely dominates this spectral range (see Fig. 2). Figures 11 and 12 also show that the spectral noise tends to shift downward the CCF in comparison to the noise-free CCF because the noisy spectra are less similar to the mask than the noise-free ones. In U580 (Fig. 12), we see that the distance between the noisy CCFs and the reference value is not similar in the upper and lower left panels, despite the same S/N. A greater distance is seen in U580 low than in U580 up because there are more weak lines in the low setup for our simulated star, and therefore they quickly vanish in the noise when the S/N drops.

The right panels of Figs. 11 and 12 show the effect of SIGMA on the derived radial velocity (uncertainty and/or bias). Our simulations clearly demonstrate that SIGMA has to be chosen in a specific range to ensure reliable results. While our simulated UVES CCFs show that the doe performance is very stable for any value of SIGMA, our simulated GIRAFFE CCFs show that only a limited range of SIGMA values can ensure reliable velocity measurements. Figure 11 suggests keeping SIGMA between ~2 and ~8 km s^-1 for HR10 and ~2 and ~7 km s^-1 for HR21, in agreement with our empirical calibration on a subsample of real GES CCFs (see Table 1). The behaviour of doe, while varying SIGMA, is different for GIRAFFE and UVES CCFs (Figs. 11 and 12). This is not due to the S/N, but rather to the sampling of the velocity grid onto which the GIRAFFE and UVES CCFs are computed, i.e. SIGMA is related to the velocity step of the CCFs. Indeed, in Sect. 2.1, we recall that the sampling frequency of the CCF is lower for GIRAFFE CCFs than for UVES CCFs: as SIGMA increases, a pronounced asymmetry on the second derivative appears for GIRAFFE CCFs, resulting in the high scatter displayed by Fig. 11.

Our simulations allow us to quantify the effect of the S/N of the spectra on the method. For U520 and U580, the standard deviation on the radial velocity at the recommended SIGMA goes from 0.05 km s^-1 at S/N = 5 to lower than 0.01 km s^-1 at S/N = 50. For GIRAFFE HR10, it goes from 0.20 km s^-1 at S/N = 5 to 0.02 km s^-1 at S/N = 50. For GIRAFFE HR21, the situation is the worst of all the setups with a standard deviation going from 0.25 km s^-1 at S/N = 10 to 0.06 km s^-1 at S/N = 50. The obvious conclusion is that the UVES setups tend to give more precise results for a given S/N compared to the GIRAFFE setups. This is understandable since a single UVES spectrum has a higher resolution and a larger wavelength coverage than any GIRAFFE spectrum. For our simulated star, the precision on the radial velocity derived by doe is up to five times higher for UVES setups than for GIRAFFE HR10 (this is even worse when compared with HR21).

This first approach of simulated CCFs shows that the method is quite robust with respect to the noise level in the GES spectra. Obviously, the presence of multiple components in the CCF may shift the detected radial velocities especially when the peaks blend with one another. In this case, the inaccuracy on the radial velocity can reach several km s^-1 (increasing as the blending degree increases). No quantitative calculations have been performed so far, but the middle panel of Fig. 7 shows a good example: the main peak is detected at 0.95 km s^-1 of its expected position and the second peak at 2.3 km s^-1, with a simulated distance of 24 km s^-1 between the two peaks. We conclude that the (conservative) random uncertainty on the radial velocity derived by doe is of the order of ± 0.25 km s^-1, while the systematic uncertainty is lower than 0.05 km s^-1 for single-peak CCFs and may reach a few km s^-1 for multi-peak CCFs. Other effects, like template mismatch or imperfect normalization, may have an effect on the uncertainty on the derived radial velocity. We also refer the reader to Jackson et al. (2015) where a discussion on the radial velocity uncertainties can be found, along with their empirical calibrations as a function of S/N, vsini, and the effective temperature of the source for GIRAFFE HR10, HR15N, and HR21 setups. As shown by Sacco et al. (2014) and Jackson et al. (2015), the errors on the GES radial velocities for most of the stars are dominated by the zero-point systematic errors of the wavelength calibration, which are not discussed here.

3.5. Detection efficiency as a function of S/N

Using Monte Carlo simulations, we assessed the impact of the S/N of GIRAFFE HR10 and HR21 spectra on the detection efficiency of the double-peaked CCF of an SB2. For that purpose we simulated synthetic SB2 spectra (a pair of twin stars) varying the S/N (from 1 to 100) and varying the difference in radial velocity of the two components Δv_rad (from 5 to 100 km s^-1). For each pair (Δv_rad, S/N), we computed as above 251 realizations of the spectra and their corresponding CCFs. We then applied doe with the parameters adapted to each setup (see Table 1).

The maps in Fig. 13 show the detection efficiency in HR10 and HR21. The green dots (respectively, the red triangles) indicate (Δv_rad, S/N) conditions when doe is able to detect the two expected peaks in more than 95% of cases (respectively, conditions when doe failed to detect the two expected peaks in more than 95% of cases). Blue plusses represent intermediate cases making detection efficiency dependent on the noise: due to the noise, spurious peaks may appear (i.e. detection failed) or the two peaks may have different heights (despite being twin stars) and become discernible to doe for small Δv_rad (i.e. detection succeeded; e.g. for HR21, at S/N = 10 and Δv_rad = 25 km s^-1).

These simulations demonstrate that even spectra with very low S/N carry sufficient information to reveal the binary nature of the targets. Specifically, in the HR10 setup, double peaks are detected in 95% of the cases when S/N ≥ 2 and Δv_rad ≥ 25 km s^-1, while in the HR21 setup, they are detected at the same rate when S/N ≥ 5 and Δv_rad ≥ 45 km s^-1. Thus, the S/N threshold that we adopted (i.e. analysis of CCFs for all spectra with S/N ≥ 5) protect us from mixed cases, which tend to happen for the lowest levels of S/N. This also shows that the HR10 setup is more able to detect SB2 than HR21 because HR21 is located around the IR Ca ii triplet whose lines have strong wings that decrease the detection efficiency. In Sect. 4.2, the histogram of the radial velocity separation of the effectively detected SB2 candidates is presented (Fig. 18). Observationally, HR10 spectra (respectively, HR21) allow us to detect SB2 with Δv_rad as low as ~25 km s^-1 (respectively, ~60 km s^-1): for both setups we are dealing with cases falling in the green dotted area of the maps. Thus, we expect in all cases an SB2 detection efficiency better than 95%.

4. iDR4 results and discussion

The doe code is included in a specifically designed workflow to handle all the GES single-exposure spectra for all setups. The automated workflow includes three steps: first, the CCFs are selected using the set of criteria described in Sect. 2.2; second, the doe code is applied to the CCFs to identify the number of peaks and a confidence flag is assigned; third, the CCFs in a given setup are combined per star and a last criterion is applied: for a given star, if more than 75% of the CCFs in at least one setup show two peaks (respectively, three and four), then the star is classified as SB2 candidate (respectively, SB3 and SB4). This rather restrictive criterion (see Sect. 4.7) is adopted to prevent false positive SB detections (due to spectra normalization, cosmics, nebular lines, etc.).

After this automatic procedure, a visual inspection is performed to ensure that no false positive detection remains and that the confidence flag is relevant. We investigate the CCFs and the spectra of all the SBn candidates one by one. When a clear false detection is encountered, the SB candidate is removed from the list. When an SB is flagged by the automatic process as probable (A) or possible (B), but the visual inspection of the CCF series (all setups considered) casts doubts on this classification, the corresponding spectra for that object are inspected. The choice of the final flag for an object can be downgraded if other CCFs provide discrepant results. This procedure ensures that processes other than binarity moderately contaminate SB candidates flagged C, marginally contaminate SB candidates flagged B, and exceptionally contaminate those flagged A. Despite these difficulties, adopting clear classification criteria ensures the best possible consistency throughout the survey.

The SBn candidates reported in the present paper are much fainter on average than those listed in the Ninth Catalogue of Spectroscopic Binary Orbits (SB9; Pourbaix et al. 2004) (Fig. 14). The average visual magnitude of SB2 within the SB9 catalogue is around V ~ 8. For the GES SB2 candidates, the average is V ~ 15. The Gaia-ESO program targets both Milky Way field stars and stars in open and globular clusters. We refer the reader to Stonkutė et al. (2016) for the selection function of Milky Way field stars (excluding the bulge stars), to Bragaglia (2012) and Bragaglia et al. (in prep.) for the selection criteria in open clusters, and to Pancino et al. (2017) for the criteria in globular clusters and calibration open clusters. We note that the targets observed in regions like the bulge, Cha I (Sacco et al. 2017), and γ² Vel (Prisinzano et al. 2016) associations, as well as the ρ Oph (Rigliaco et al. 2016) molecular cloud, are selected on the basis of coordinates and photometry (VISTA and 2MASS), thus providing a rough membership criterion.

The list of the SB2 and SB3 candidates in the Milky Way field is given in Tables A.1 and A.2. The list of SB2 in the bulge, the Cha I, γ² Vel, and ρ Oph associations and the CoRoT field is given in Table B.1. Finally, the list of SBn in stellar clusters is given in Table C.1. The results (classifications and confidence flags) are included in the GES public releases (see footnote 3) using the nomenclature given in the GES outlier dictionary developed by the GES Working Group 14 (WG 14)⁶.

4.1. Binary classification

The binary classification⁷ was developed for the GES within WG 14. The following scheme is adopted: the peculiarity flag is built from the juxtaposition of a peculiarity index and a confidence flag letter. The peculiarity index is defined as 20n0, with n ≥ 2, where n is the number of distinct velocity components in the CCF. With this peculiarity index, an SB2 is classified as 2020, an SB3 2030, etc. Even though a star is flagged 2020 (i.e. SB2), a third component may be present but not visible during the observation or may be undetectable at the resolution and S/N of the considered exposure.

Moreover, the WG14 dictionary recommends the use of confidence flags (A: probable, B: possible, and C: tentative). Clearly, the closer the CCF peaks are, the less certain the detection is. The criteria to allocate these flags were defined as follows:

• A: the local minimum between peaks is deeper than 50% of the fullamplitude of the highest peak;

• B: the local minimum between peaks is higher than 50% of the full amplitude of the highest peak;

• C: no local minimum is detected between peaks, but the CCF slope changes.

With these definitions, the SB2 whose CCF is plotted in the left, middle, and right panels of Figs. 8 and 9 would be flagged as A, B, and C, respectively.

For triple-peak CCFs, the same type of criteria are applied to the second local minimum. If this second local minimum is lower than 70% of the full amplitude of the highest peak, then the confidence flag is set to A, else B. Examples of these two cases are shown on simulated CCFs in Fig. 15. The CCF in the middle panel is classified as 2030B because the leftmost local minimum is higher than 0.5 times the largest amplitude, but also as 2020A because the middle and leftmost peaks, taken as a whole, are well separated from the rightmost peak.

4.2. iDR4 SB2 candidates

Table 2 presents the breakdown of the detected SBn candidates in terms of confidence flags, whereas Table 3 provides the detailed results of the analysis per field in terms of automated detection (“doe”) and after visual checking (“confirmed”). A total of 1092 sources were identified as SB2 candidates by the automated procedure described in the previous section, 342 of which were confirmed after visual inspection, giving a success rate of about 30% similar to that of Matijevič et al. (2010) for the RAVE survey. Typical rejected cases include distorted CCFs caused by negatives fluxes or pulsating stars. Some confidence flags were also changed during the visual inspection phase (see Sect. 4.7). The largest number of stars has been observed with the GIRAFFE setup HR21 because it corresponds to the Gaia wavelength range of the radial velocity spectrometer. However, the rate of SBn detection in this setup is very low because it is dominated by the presence of the Ca ii triplet, which is a very strong feature in late-type stars, thus resulting in a broad CCF that can mask possible multiple peaks (Fig. 2, bottom panel). Moreover, emission in the line cores of this triplet induces fake double-peak CCFs because in the templates the lines are always in absorption. Consequently, it is very difficult to identify double peaks due to binarity based on HR21 CCFs (see Sect. 4.7 for more details). This explains why we only have two firm detections among the 31 970 stars observed with this setup alone. Hence, this setup is not well-suited to detecting stellar multiplicity at least in our situation (see Matijevič et al. 2010: although they were able to discover 123 SB2 out of 26 000 RAVE targets, they also had to deal with very broadened CCFs and could not detect binaries with Δv_rad ≤ 50 km s^-1).

The setup with the second largest number of observed objects is HR10. This setup covers the range [535−560] nm with many small absorption lines that result in a narrow CCF, suitable for the detection of stellar multiplicity (see Fig. 8). The largest number of probable SB2 candidates is indeed detected with this setup.

To illustrate how some setups are more adapted than others to detect SBn, we show spectra and CCFs in these setups for single stars (the Sun and Arcturus in Figs. 2 and 3) and for an SB2 candidate (NGC 6705 1936 observed in most of the GIRAFFE setups where the composite nature of the spectrum is clearly visible in Fig. 16).

Contrary to field stars, which are observed in HR10 and HR21 only, cluster stars were observed with many different setups. The number of SB2 candidates in the field is 185 out of 27786 stars (0.67%) whereas in the clusters it amounts to 127 out of 16468 (0.77%, see Table 3).

There are about 30 SB2 candidates detected with a double-peaked CCF in both GIRAFFE HR10 and HR21. For instance, the field star 02394731-0057248 (magnitude V = 13.8) is identified as an SB2 candidate with HR10 and HR21 (see Fig. 17). This new candidate has no entry in the Simbad database.

The histograms of the radial velocity separation of SB2 candidates for GIRAFFE HR10 and HR21 and for UVES U580 are shown in Fig. 18 (U520 is not represented, due to insufficient statistics). The smallest measured radial velocity separations are 23.3, 60.9, and 15.2 km s^-1 for HR10, HR21, and U580, respectively. This is well in line with the detection capabilities of the doe code as mentioned in Sect. 3.3 (~30 km s^-1 for GIRAFFE and ~15 km s^-1 for UVES setups). In U580, the high bin value around 72 km s^-1 is mainly due to the repeated observations of a specific object, the SB4 candidate 08414659-5303449 in IC 2391 (see Sect. 4.5).

Concerning the SB2 candidates in open clusters, not only did we check the cleanliness of the SB2 CCF profile, but we also compared the velocities of the two peaks with the cluster velocity. Assuming that most of the SB2 systems discovered by GES generally have components of about equal masses, then an SB2 that is member of the cluster should have a cluster velocity about midway between the two component velocities. This simple test allows us to assess the likelihood that the SB2 system is a cluster member. This method is applied for the SB2, SB3, and SB4 candidates analysed and full details are given in the present section and in Sects. 4.4 and 4.5. The results are shown in Table C.1. The column labelled “Member” in Table C.1 evaluates the likelihood of cluster membership based on the component velocities: if the cluster velocity falls in the range encompassed by the component velocities, we assume that the centre of mass of the system moves at the cluster velocity, which means that membership is likely. In that case, we put “y” in the column. On the contrary, if the CCF exhibits two well-defined peaks not encompassing the cluster velocity, the star is labelled as an SB2 non-member of the cluster (“n” in the column). Another possibility is that one component has a velocity close to that of the cluster and the second velocity is offset. In that case, the SB2 nature is questionable and the star is more probably a pulsating star (responsible for the secondary peak or bump) belonging to the cluster (“y” in column “Member”). The list of individual radial velocities based on iDR5 data will be given in a forthcoming paper. More extended remarks for each cluster are provided in Appendix C.

4.3. Orbital elements of two confirmed SB2 in clusters

With the data collected so far, we were able to confirm the binary nature of two SB2 candidates in clusters by deriving reliable orbital solutions for the systems 06404608+0949173 (NGC 2264 92) and 19013257-0027338 (Berkeley 81, hereafter Be 81).

The first system 06404608+0949173 (magnitude V ~ 12) is a bona fide SB2 for which 24 spectra are available (20 GIRAFFE HR15N and 4 UVES U580) and an orbit can be computed, as shown in Fig. 19. Observations where only one velocity component is detected are not used to calculate the orbital solution because these velocities are not accurate (Fig. 19) since the two velocity components are blended. The orbital elements are listed in Table 4. The short period of 2.9637 ± 0.0002 d implies that neither of the components can be a giant, which is consistent with the classification of the system as K0 IV (Walker 1956). The centre-of-mass velocity of the system (14.6 km s^-1) is close to the cluster velocity (17.7 km s^-1), as it should be. The mass ratio is M_B/M_A = 1.10. Classified as FK Com in the GCVS (=V642 Mon), this source is chromospherically active with X-ray emission (ROSAT and XMM). This system thus adds to the two SB2 systems with available orbits (VSB 111 and VSB 126) already known in NGC 2264 (Karnath et al. 2013).

The second system 19013257-0027338 (magnitude V ~ 17) is a confirmed SB2 (2020 A) for which 18 spectra are available (8 GIRAFFE HR15N and 10 GIRAFFE HR9B). This source is not listed in the Simbad database. The orbital elements are given in Table 4 and the orbit is displayed in Fig. 20. Strangely enough, a good SB2 solution for this system could only be obtained by adding an extra parameter to the orbital elements, namely an offset between the systemic velocities derived from component A and from component B (see the ΔV_B term in Eq. (2) of Pourbaix & Boffin 2016). In most cases this offset is null, but there could be situations where it is not, like in the presence of gravitational redshifts or convective blueshifts that are different for components A and B (Pourbaix & Boffin 2016). Alternatively, if the spectrum of one of the components forms in an expanding wind (as in a Wolf-Rayet star), it would also lead to such an offset. However, what is puzzling in the considered case is the large value of the offset (24.8 ± 1.2 km s^-1) for which we could not find any convincing explanation. Indeed, no Wolf-Rayet stars are known in the Be 81 cluster according to the Simbad database. This very diffuse cluster of intermediate age lies towards the Galactic centre (Hayes & Friel 2014; Donati et al. 2014).

4.4. SB3 candidates

Tables 2 and 3 show that, in total, 11 SB3 candidates (7 probable: flag A, 1 possible: flag B, and 3 tentative: flag C) were detected. Five of these SB3 are found in the field (Fig. 21 and Table A.2) and six in clusters (Fig. 22 and within Table C.1). A total of 266 targets were initially labelled as SB3 candidates by the doe code, while only 11 were kept after visual inspection, giving a success rate of about 4% (compared to 30% for SB2 detection). The SB3 candidates are essentially detected in UVES setups and in GIRAFFE setups HR9B and HR10. The SB3 candidates in the stellar clusters were examined on a case-by-case basis, and the results are reported below.

NGC 2451.

The CCF of 07470917-3859003 exhibits three clear peaks (the CCF is classified as 2030A), at 25.0, 96.1, and 136.6 km s^-1. The first velocity is compatible with membership in NGC 2451A. The DSS⁸ image reveals the presence of a slightly fainter star about 12′′ south (a greater distance than the 1.2′′ size of the fibre, so no contamination is possible). Given the fact that the two fainter peaks are not located symmetrically with respect to the cluster velocity, it is doubtful that the system could be a physical triple system in the case of membership to NGC 2451.

NGC 2516.

NGC 2516 45 (system 07575737-6044162) is a star classified as A2 V (Hartoog 1976) with V = 9.9. The iDR4 recommended parameters (T_eff= 8500 K, log g = 4.1, and solar metallicity) suggest that it could be a δ Scu star. Its CCF is most likely associated with a fast rotator with a superimposed sharper central peak. The SB3 nature of this candidate is therefore doubtful and a follow-up of this source should be performed before drawing any firm conclusion.

NGC 6705.

In total, the doe routine finds 52 SB3 candidates in NGC 6705, one of the largest number of SB3 among all the targeted clusters (Table 3). After a first-pass analysis we discarded all of them but one, NGC 6705 1147 (system 18510286-0615250). The velocities corresponding to the three peaks observed in the CCF are listed in Table 5. They exhibit clear temporal variations. The cluster velocity is 29.5 km s^-1 (Cantat-Gaudin et al. 2014). This velocity is close to that of the middle peak in the CCF (C, i.e. the faintest). That central peak does not vary as much as the most extreme peaks, and moreover, the shape of peak C is not as sharp as are peaks A and B. Considering that the cluster NGC 6705 is a dense one, we believe that this third peak is from background contamination. We therefore conclude that the detection of NGC 6705 1147 as SB3 is spurious and should be downgraded to SB2. The SB2 analysis is presented in Table 5 where we computed the mass ratio, adopting 34 km s^-1 (Table 3) as the centre-of-mass (cluster) velocity. The observed velocity variations are consistent at all times with a mass ratio of the order of 1.32.

NGC 6005.

The CCF of 15553867-5724434 (classified as 2030B) shows three peaks, at −81.6, −14.4, and 32.7 km s^-1, compared with −25.2 km s^-1 for the cluster velocity (Carlberg 2014). The spectra are at the minimum required S/N. These data are compatible with 15553867-5724434 being an SB3 member of NGC 6005.

NGC 6802.

The CCF of 19302315+2013406 (classified as 2030C) shows three distinct peaks, at −22.4, 22.0, and 65.5 km s^-1, compared with 12.4 km s^-1 for the cluster velocity (Hayes & Friel 2014). These data are compatible with 19302315+2013406 being an SB3 member of NGC 6802.

Trumpler 20.

The CCF of 12391904-6035311 (classified as 2030C) shows three distinct peaks, at −85.78, −44.4, and 14.8 km s^-1, compared with −40.8 km s^-1 for the cluster velocity (Kharchenko et al. 2005). These data are compatible with 12391904-6035311 being an SB3 member of Trumpler 20. An extended analysis of the GES data for this cluster may be found in Donati et al. (2014).

4.5. The unique SB4 candidate HD 74438

We detected one SB4 candidate: the A2V star HD 74438 (CNAME 08414659-5303449, with V = 7.58) belonging to the open cluster IC 2391 (Platais et al. 2007).

The star was observed 45 times within 2.5 h on February 18, 2014, with the U520 and U580 setups. Its peculiarity was already noticed by Platais et al. (2007) since it lies 0.9 mag above the main sequence in a colour-magnitude diagram, and therefore was already thought to be a triple system (since the maximum deviation for a binary system with two components of equal brightness would amount to 2.5 × log 2 = 0.75 mag). It is nevertheless considered a bona fide member of the cluster by Platais et al. (2007). Therefore, the centre-of-mass velocity for the system can be considered identical to the cluster velocity, namely 14.8 ± 1 km s^-1 (Platais et al. 2007). A typical example of the CCF of HD 74438 is presented in Fig. 23 where its four distinct CCF peaks are clearly apparent. The velocities of the peaks at different times over the night of February 18, 2014, are collected in Table 6. In this table, we first notice that the velocities of components A and B (which correspond to the highest peaks) vary slowly and oppositely to each other. Their amplitude of variations is similar. If we compute the velocity variations with respect to the cluster velocity (which should correspond to the centre-of-mass velocity of the AB pair, neglecting the gravitational influence of components C and D – columns Δv_r(A) and Δv_r(B) in Table 6), we note that these variations obey a simple property: their ratio is almost constant. In a simple binary system, this property is expected since the ratio Δv_r(B)/Δv_r(A) equals the mass ratio M_A/M_B. Here we find M_A/M_B ~ 1.01. Thus, the brightest components in the system, which correspond to the most prominent peaks A and B, are almost twins since their masses differ by only 1%. We observe that the pair CD obeys the same property: Δv_r(C)/Δv_r(D) is almost constant, even though the amplitude of variations is larger than that of the AB pair. Again, assuming no perturbations from the AB pair, we get M_D/M_C = Δv_r(C)/Δv_r(D) ~ 1.12. It seems, therefore, that the observed variations do make sense and give credit for a physical nature of the ABCD system as a double pair AB/CD. We could nevertheless expect some perturbations of one pair on the other, at least in the form of a trend of the centre-of-mass velocities of each pair, if pair CD orbits around pair AB. The available observations do not span a time interval long enough to check that possibility.

Assuming that the ratio of the CCF amplitudes roughly scales with the luminosity ratio⁹, and adopting a ratio of 3 between the peak amplitudes of A and D (see Fig. 23), we get a magnitude difference between components A and D equal to Δm = 2.5log 3 = 1.2 mag. Consequently, the observed visual magnitude m_V = 7.58 is mainly due to the pair AB. With the parallax of the system π = 5.716 ± 0.298 mas provided by Gaia DR1 (Gaia Collaboration 2016), the distance of this system is 175 ± 9 pc. The absolute magnitude of AB pair is then M_V(AB) = 1.36. Assuming similar masses, we have M_V(A) = M_V(B) = 2.12. This corresponds to a spectral type A7 and to masses of M_A = M_B = 1.8M_⊙ if on the main sequence (luminosity class V). The absolute magnitudes of components C and D are consequently M_V(C) = M_V(D) = 3.31 corresponding to a spectral type F1, which correspond to a mass of about 1.5 M_⊙. Inserting these values in the defining relation for the orbital velocity semi-amplitude (expressed in km s^-1) (3)it is possible to derive an upper limit to the orbital period. Indeed, for the AB pair, we adopt e = 0, q = 1, M_A = 1.8M_⊙, and K_A> 36 km s^-1 (Table 6), and obtain an upper limit on the orbital period of the AB pair, P (d) < 93sin³i. The same method applied to the CD pair (with M_D = 1.5M_⊙, e = 0, q = 1.1, and K_D> 60 km s^-1) yields P (d) < 20sin³i, in agreement with the fast variation observed in Table 6 for the C and D velocities.

An even more constraining limitation on the orbital period can be derived from the fast variations exhibited by the D component over the 2.2 h time span covered by the observations (Table 6). We first assume that 74.55 km s^-1 corresponds to the maximum orbital velocity, from which we derive a semi-amplitude K_D = 59.75 km s^-1, corresponding to ωt₁ = π/ 2. It is then possible to find ω = 2π/P, and hence P, by assuming a sinusoidal velocity variation (in a circular orbit), reaching a velocity of 67.9 km s^-1 at time t₂<t₁, such that . From this, we derive ωt₁−ωt₂ = π/ 2−1.093 = 0.477 = 2π/P × (t₁−t₂) = 2π/P × 2.2 h, or P = 29 h as the tentative period of the CD pair.

To conclude, we note that the above arguments also allow us to estimate the deviation of HD 74438 in the colour−magnitude diagram for a system consisting of components with fluxes F_A = F_B, and F_D = F_C = 1/3F_A. The magnitude excess amounts to 2.5log (2 + 2/3) = 1.1 mag, not far from the 0.9 mag reported by Platais et al. (2007). The velocities of the components would definitely be worth monitoring over a few hundred days.

4.6. Multiplicity flagging by other GES working groups

It is worth mentioning that different nodes within the GES WGs have identified/detected spectroscopic systems for restricted subsamples of iDR4 data. Because we wanted to rely on a homogeneous detection process, we did not include the SBn detected by other WG in the present analysis. This detailed comparison will be performed for the next data release.

Working group 12, focusing on pre-main sequence stars in clusters, detected 176 SB2 (A: 168, B:2, C: 6), 1 SB3 and 2 SB4. The intersection with our list amounts to 66. In particular, the two SB4 detected by WG12 are classified as SB2 in our final list; we re-checked that only two peaks are visible on the CCFs computed from single exposures. Working group 12 developed a specific method of removing nebular contamination by masking the nebular lines in HR15N spectra for the clusters NGC 2264, NGC 6530, and Tr14. Indeed, these nebular lines can produce a double-peaked CCF that can be misclassified as an SB2 candidate; see Klutsch et al. (in prep.) for more details.

Working group 13, dedicated to OBA-star spectrum analyses, identified about 30 SB2 in clusters (NGC 2547, NGC 3293, NGC 6705 and Tr 14). They detected one SB3 candidate (system 10344470-5805229 in NGC 3293) that we have rejected. Indeed, the three peaks were detected by our method only in two CCFs and only in the HR5A setup, whereas ten CCFs of the same object displayed only one or two peaks in various other setups (HR3, HR6, HR9B, and HR14A). This SB3 detection was therefore not reliable enough considering our rejection criteria (discussed at the beginning of Sect. 4). However, we did not even select this object as an SB2 because the velocity difference between the two peaks is too large (>290 km s^-1), indicating possible spurious peak(s).

In summary, the GES working groups, which are very focused, will inevitably reach higher detection rates for specific types of objects, but their methods do not apply to the whole GES survey. The method presented here, on the contrary, provides homogeneous information for the whole survey, using all (GIRAFFE and UVES) individual spectra.

4.7. Multiple peak CCFs unrelated to binarity

Double- and triple-component CCFs may sometimes be mimicked by physical processes unrelated to binarity. To clearly establish the binary nature of field stars, multiple observations covering a complete orbital cycle are mandatory in order to derive the orbital elements that fit best the radial velocities. In the case of stars belonging to associations and clusters containing hot and cold gas, the situation is worse: emission lines, which are not masked prior to the CCF computations, may produce troughs in the CCFs that can be interpreted as multiple peaks. Moreover, hot and pulsating stars like δ Scu stars, or young hot stars with discs, may also produce bumps in the CCFs. It is beyond the scope of the present paper to study the specific signatures of such processes on the CCFs, which also depend on the considered setup. However, we provide below some examples of multiple peak CCFs probably unrelated to binarity. Furthermore, in order to remove some spectral signatures degrading the CCFs (emission lines, very strong lines, etc.), we plan to recompute consistently all GIRAFFE and UVES CCFs in a forthcoming paper.

For instance, NGC 6705 1507 (system 18505296-0617402) is classified as A0 (Cantat-Gaudin et al. 2014) and shows three peaks in its CCF (originally classified as 2030C; Fig. 24) for the setting HR6, at −25.1, 33.8, and 86.9 km s^-1. The central, highest peak is close to the cluster velocity, and the other two are almost symmetrically located from the central peak, at ± 50 km s^-1. The very edge of the CCF has a steep slope which is reminiscent of a fast rotator. Indeed, the full base width of the CCF is about 180 km s^-1, a value typical of the rotation velocities of A stars. Moreover, a spectrum in the HR9B setting, taken on the same night, confirms the above analysis, which leads us conclude that the triple-peak CCF of star 1507 in the cluster NGC 6705 is most probably caused by pseudo-absorptions superimposed on a rapid-rotator profile. A similar situation is encountered for the two other SB3 candidates 18510403-0616023 and 18511155-0606094. These three objects have been discarded from the final list.

An example of a star automatically classified as an SB2 with flag C and very likely to be instead a δ Scu star, i.e. a hot rapid rotator with pulsation and no emission in Hα, is 18503348-0619555 (NGC 6705 1916, V = 13.7). This star has recommended parameters of T_eff= 7821 K and log g = 3.96, compatible with a δ Scu-type star. The CCFs in different setups at different epochs are shown in Fig. 25. The first CCF has two components (SB2), one broader than the other. The asymmetry of the second CCF could potentially lead the doe code to identify three components (SB3). The last CCF is less ambiguous though it can be seen as an SB2 with close radial velocities. This SB2 candidate has been removed from the final list of SBn candidates.

In Trumpler 14, spectra are strongly contaminated by nebular lines around Hα. This may result from a reduction issue (inadequate sky subtraction in a nebular background). The nebular lines in emission, located at the cluster velocity, superimpose on the absorption lines of the star, also at the cluster velocity. Because the nebular lines in emission are narrower than the stellar lines in absorption, it results in a CCF with two clear peaks; sometimes the minimum between the two peaks goes even lower than the CCF continuum. Such false SB2 candidates could be unmasked (see Fig. 26 for the two examples 10433966-5935573 and 10444601-5935228) because in that cluster we found too many stars with radial velocities around −40 and 20 km s^-1, i.e. symmetrical with respect to the cluster velocity (~−10 km s^-1). They can be explained by an emission at the cluster velocity obliterating the Hα line resulting in a central absorption splitting the CCF (an emission line corresponds to absorption in the CCF).

4.8. Distribution in the (log g, T_eff) plane

The GES consortium provides recommended atmospheric parameters (T_eff, log g, and [Fe/H]) for 63% of the stars from the iDR4. They result from a delicate merging of atmospheric parameters obtained by different WGs using different methods, but all with the same model atmospheres and linelists. Among them, we identified a hundred of our confirmed SBn candidates, representing 30% of our detected SBn).

They are shown in the log g – T_eff plane (see Fig. 27). This figure reveals a sudden drop in the number of stars surveyed above 7000 K. This threshold corresponds to the transition between A and F stars, the latter being surveyed in a systematic way by the GES, the former being included only if they belong to specific clusters. For SBn, the atmospheric parameters provided by the GES pipeline are uncertain (or even wrong, as we show below) because (i) composite spectra cannot be fitted with single synthetic spectra, and (ii) spectra fitted by the automated pipelines are not individual exposures but rather stacked ones. Despite these shortcomings, the log g – T_eff diagram nevertheless allows us to identify systems of interest.

The two SB2 and SB3 systems on the warm end of the log g – T_eff diagram (with T_eff> 8000 K) are worth discussing. Their CNAMEs are 18280622+0642252 (NGC 6633 110, BD+06 3793, A3V), classified as 2020A, and 07575737-6044162 (NGC 2516 45, CD-60 1959, A2V), classified as 2030C.

System 18280622+0642252 shows two peaks of equal heights at −70 km s^-1 and 38 km s^-1. The first peak is particularly broad and is probably associated with a rapidly rotating star. Since the cluster velocity (−25.4 km s^-1) lies between the two peaks, and the double-peak CCF is very well-defined, we confirm the SB2 flag from the doe routine.

System 07575737-6044162 exhibits a broad CCF most likely associated with a fast rotator. It has a sharp central peak. It may be a physical double system, but certainly not a triple one (see Fig. 22).

The three giant SB2 candidates (19262489+0137506, 22180319-5834560, and 11265745-4100160) appearing in the log g – T_eff diagram (with log g < 2) are surprising since they should have a mass ratio very close to 1. Their CCFs are displayed in Fig. 28. To our knowledge, there are only a few SB2 systems known so far involving two giant stars: (i) HD 172481 (more precisely an F2Ia post-AGB star and an M giant; Reyniers & Van Winckel (2001), Jorissen et al. (2009); (ii) HD 187669 (a double-line eclipsing binary; Hełminiak et al. 2015a); (iii) TYC 6861-523-1/ASAS J182510-2435.5 (Ratajczak et al. 2013); (iv) KIC 09246715 (a double-lined spectroscopic and eclipsing binary; Hełminiak et al. 2015b).

System 19262489+0137506 (a CoRoT target with T_eff = 4300 K, log g = 1.0), classified as 2020A, has two peaks well separated by 117 km s^-1, of almost equal intensities, implying a rather short period for a pair of giants (Fig. 28, middle). Adopting K = 117/2 km s^-1, q = 1, sini = 1, e = 0, and M₁ = 1M_⊙, Eq. (3) predicts a period of the order of 7.5 d for the associated binary. This is rather short considering the giant nature of the two components. For instance, the minimum orbital period in the large sample of binaries with a K giant component in open clusters (Mermilliod et al. 2007) is just above 25 d. The situation is even worse for the sample of field M giants from Jorissen et al. (2009) where the shortest orbital period is above 200 d. This trend reflects the increase in the stellar radius along the giant branch. Independently, the spectral type of the system was estimated to be M2III from broad-band photometry (Exo-Dat, Deleuil et al. 2009). In any case, this system is worth a follow-up investigation, especially looking for signs of mass-transfer activity (e.g. possible Hα emission in its spectrum, but the two spectra available in HR15 are too noisy to see any such sign of activity).

System 22180319-5834560, classified as 2020C (and T_eff = 4100 K, log g = 1.8), exhibits a very broad CCF coming from the strong Ca ii triplet in the HR21 setup, with two bumps responsible for the SB2 classification (Fig. 28, bottom). Observations in HR10 one day later does not show any sign of binarity. Inspection of the HR21 spectra reveals that the bumps observed in the CCF may be due to emission in the Ca ii triplet line cores, making the SB nature doubtful.

System 11265745-4100160 (V = 13), classified as 2020B (with T_eff = 4400 K, log g = 1.9, top CCF in Fig. 28), exhibits two close velocity components in HR10 (separated by about 32 km s^-1) but not visible in HR21. The validity of the atmospheric parameters may have been disturbed by the SB2 nature of the star.

4.9. Comparison with other catalogues

To estimate the proportion of new SBn candidates, we cross-checked our 352 distinct candidates with published online catalogues of stars. The intersection with the Simbad database (Wenger et al. 2000) provides 96 matches. Among them one is classified as a double or multiple star (WDS J08513+1150, CNAME 08511901+1150056 belonging to M67), anf four as spectroscopic binary stars: 2MASS J06435849-0100515 (CNAME 06435847-0100516) in the Corot field, CD-52 2472 (CNAME 08385566-5257516) in the cluster IC2391, 2MASS J08512291+1148493 (CNAME 08512291+1148493) in M67, and NGC 2682 165 (CNAME 08512940+1154139) also in M67. Two are classified as eclipsing binary stars: 2MASS J18505933-0622051 (CNAME 18505933-0622051) in NGC 6705 and CoRoT 101129018 (CNAME 19263739+0152562). All these previously known binaries have been attributed a “A” confidence flag by our workflow.

We cross-matched our detections with various other catalogues, using the X-Match and the Vizier Search online tools from the CDS¹⁰ by uploading the J2000 coordinates built from the CNAME of our SB candidates. For each catalogue, we set the matching area within a radius of 3 arcsec.

The comparison with the Ninth Catalogue of Spectroscopic Binary Orbits (Pourbaix et al. 2004, SB9) leads to three systems in common, namely 08511868+1147026, 08511901+1150056, and 08512291+1148493, which are members of the M67 cluster (NGC 2682) with a visual magnitude of about 12.5.

The comparison with the Washington visual Double Star catalogue (WDS, Mason et al. 2016) leads to an intersection of two systems, namely WDS J08513+1150 in M 67 and WDS J11088-7519 in Cha I (CNAME 08511901+1150056 and 11085326-7519374, respectively).

Cross-matches with the Geneva-Copenhagen Survey of the solar neighbourhood III (Holmberg et al. 2009), with the bibliographic catalogue of stellar radial velocities (Malaroda et al. 2006), with the RAVE catalogue of SB2 candidates (Matijevič et al. 2010), and with the Multiple Star Catalogue (MSC) (Tokovinin 1997) resulted in empty intersections. We note that the limiting magnitudes of all these catalogues are much brighter than that of the GES (V ~ 19); therefore, we expected a small intersection.

Four of our SBn candidates are known in the WEBDA cluster database, with available orbital parameters (see Table 7). We also found that two SB2 known in WEBDA are observed in iDR4, but were discarded by the workflow. M67 111 has been observed (08511799+1145541), but the second peak is too low to be automatically detected. The same issue occurs with NGC 2264 134 (06405934+0955201), known to be an SB2 in WEBDA and known to be a pre-main sequence star. It has been observed eight times and seems to be an SB3 candidate because four CCFs have one peak, two CCFs have two peaks, and two CCFs have three peaks (see Fig. 29).

For the sake of completeness, we also checked whether the doe algorithm retrieved the known SBn candidates from the Geneva-Copenhagen Survey and from SB9. It turns out that only one SB2 (08511799+1145541 in M67) present in SB9 was not found by doe: ten observations in U580 were performed, but the second peak is only visible and detected in two of them. Because only stars with more than 75% of multiple peaks detection in a given setup were flagged as SB2 candidates, 08511799+1145541 was rejected. This shows that the 75% criterion, chosen to be conservative, might be too restrictive in some cases, although it prevents many false positive detections.

Previously known SB2 systems flagged as such by the GES are listed in Table 7. We note that the analysis of the GES data provides a substantial number of new SB2 and SB3 candidates because the SB detection was performed on a huge data sample (~51 000 stars) characterized by a faint limiting magnitude with respect to previous surveys. The new SB2, SB3, and SB4 candidates clearly deserve more observations in order to derive their orbital elements.

5. Conclusion

We present a method for identifying multiple-lined spectroscopic binaries (SBn, n ≥ 2) based on the successive derivatives of the CCFs. A list of SBn among the GES iDR4, in the Galactic field and in the stellar clusters, is presented. In addition, orbital solutions for binary systems belonging to the open clusters NGC 2264 and Be 81 have been calculated.

The detection method has been tested on all the setups of the GIRAFFE and UVES spectrographs available within the GES. It turns out that UVES U580 and GIRAFFE HR10 are the most appropriate setups for detecting multiplicity with velocity differences as low as 15 km s^-1 and 23 km s^-1, respectively. Simulations show that the doe algorithm reliably derives radial velocities (with a formal error of the order of 0.20 km s^-1 at a typical S/N of 10 for GIRAFFE and lower than 0.01 km s^-1 at S/N = 50 for UVES setups; for multi-component CCFs the formal error will be slightly higher and, in addition, the systematic error may reach a few km s^-1 at the detection limit).

The detection method leads to a number of false positive detections in stellar clusters. Using physical properties of the clusters and combining information from the spectra and CCFs of different setups, we discussed and discarded a fraction of candidates. A confusing SB2-like signature could be imprinted to the CCF by pulsations in δ Scuti variable stars, by Hα emission in circumstellar discs or interstellar absorption by cold clouds along the line of sight. In such cases, spurious peaks or bumps appear in the CCF.

We discovered 340 SB2, 11 SB3, and one SB4 out of 51 000 stars with more than 205 000 single exposures. The most confident binary candidates (“A” flag) most often show very clear composite spectra. Incidentally, we warn against the use of the GES recommended atmospheric parameters for these SBn candidates. Indeed, one-third of SBn candidates have GES recommended parameters, but the presence of multiple components in spectral lines can potentially lead to incorrect atmospheric parameters.

The frequency of SBn (n ≥ 2) found by our method in the GES iDR4 sample is 0.7%. Most of the SBn candidates are new because they belong to a sample of stars much fainter than was covered by previous catalogues. If we extrapolate this percentage of 0.7% SBn binaries to the final GES pool of 10⁵ stars, we expect to reach about 1000 new SBn systems in the upcoming data releases because the number of observed stars will increase by a factor of two and because we plan to further fine-tune our detection criteria. Indeed the aim of the present analysis was to detect binaries and minimize the number of false positive detections (i.e. stars incorrectly classified as SBn). The method presented in this paper can be readily applied to the ESA Gaia mission spectra.

---

Acknowledgments

T.M., M.V.d.S., and S.v.E. are supported by a grant from the Fondation ULB. This work has been partly funded by an Action de recherche concertée (ARC) from the Direction générale de l’Enseignement non obligatoire et de la Recherche scientifique – Direction de la recherche scientifique – Communauté française de Belgique. T.M. is supported by the FNRS-F.R.S. as temporary post-doctoral researcher under grant No. 2.4513.11. This work was supported by the Fonds de la Recherche Scientifique FNRS under grant No. T.0198.13. C.A. acknowledges to the Spanish grant AYA2015-63588-P within the European Founds for Regional Development (FEDER). M.T.C. acknowledges the financial support from the Spanish Ministerio de Economía y Competitividad, through grant AYA2013-40611-P. R.S. acknowledges support from the Polish Ministry of Science and Higher Education. This work was partly supported by the European Union FP7 programme through ERC grant number 320360 and by the Leverhulme Trust through grant RPG-2012-541. We acknowledge the support from INAF and Ministero dell’Istruzione, dell’Università e della Ricerca (MIUR) in the form of the grant “Premiale VLT 2012”. The results presented here benefit from discussions held during the Gaia-ESO workshops and conferences supported by the ESF (European Science Foundation) through the GREAT Research Network Programme. This research has made use of the Washington Double Star Catalogue maintained at the US Naval Observatory. This research has made use of the WEBDA database, operated at the Department of Theoretical Physics and Astrophysics of the Masaryk University. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This research has made use of Python, in particular the Python module pyfits.py which is a product of the Space Telescope Science Institute, which is operated by AURA for NASA. This research has made used of the Digitized Sky Surveys which were produced at the Space Telescope Science Institute under US Government grant NAG W-2166. The images of these surveys are based on photographic data obtained using the Oschin Schmidt Telescope on Palomar Mountain and the UK Schmidt Telescope. The plates were processed into the present compressed digital form with the permission of these institutions. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. The authors thank the referee for his comments which helped to improve the manuscript.

References

Appendix A: SB2 and SB3 candidates in the field

Appendix B: SB2 candidates in selected fields

Appendix C: SB2, SB3, and SB4 candidates in stellar clusters

IC 2391.

This open cluster includes the unique SB4 candidate 08414659-5303449 in the current iDR4 GES data. At some epochs, the two weakest components are hardly visible, which is why we classified this source with both 2020A and 2040A flags (see Fig. 23). This SB4 candidate is analysed in detail in Sect. 4.5.

IC 2602.

All three systems are consistent with cluster membership. System 10403116-6416249 seems to be a pair of rapidly rotating stars.

IC 4665.

System 17452506+0540233 has a broad CCF with a secondary bump in its tail, but the velocities are not centred on the cluster velocity.

M67.

All four SB2s are confirmed through visual inspection, having composite spectra and having membership confirmed.

NGC 2243.

Only one clear SB2 candidate with a composite spectrum is retained. Two other candidates (06292559-3116070 and 06294409-3116276) are not retained since the major peak of the CCF is at the cluster velocity, with a secondary bump offset by −60 km s^-1 and −100 km s^-1, respectively.

NGC 2264.

The doe procedure has flagged a lot of stars as SB2 in this (and in all other) young clusters. Many of these stars have broad CCFs with a secondary bump, as illustrated in Fig. C.1. As the centre of this very broad CCF is close to the cluster velocity, these stars are thought to be both rapidly rotating and pulsating (δ Scu variables), and this combination is responsible for the peculiar and specific CCFs observed in young clusters, whose turn-off is located higher up on the main sequence to allow the presence of δ Scu stars. Not all of them are A stars though, and therefore we suggest the alternative hypothesis that this peculiar CCF profile is related to the disc still surrounding these young stars. In that case, the CCFs offer an interesting diagnostic to study/detect these discs (see Rebull et al. 2002). 06405650+0911389 (HD 261905) has its main peak at the cluster velocity, and a clearly defined, well-separated second peak at a velocity of 71.9 km s^-1. Although the field is not especially crowded, the DSS image¹¹ reveals that the stellar image might be not perfectly round and seems contaminated by a nearby source. 06421531+0942581 has a secondary peak close to the cluster velocity, but the main peak is totally offset (99 km s^-1). That peak might be due to a somewhat brighter star (NGC 2264 SBL 560) located about 4 arcsec west of the target (probably not a member, given its high offset velocity).

NGC 2451.

The situation for this cluster is special since there are in fact two different clusters, located at different distances, superimposed at the same location on the sky (Dias et al. 2002). These authors report a velocity of +22.7 km s^-1 for the nearest NGC 2451A cluster and 14.0 km s^-1 for the farthest NGC 2451B cluster. 07401559-3735416, a genuine SB2 system, cannot be a member of NGC 2451. On JD 2 456 634 the CCF exhibits peaks at 21 and 62 km s^-1, while at JD 2 456 638 and JD 2 456 677, the peaks are located at around −2 and 85 km s^-1, implying a centre-of-mass velocity of the order of 40 km s^-1, significantly offset with respect to the velocities of NGC 2451A and NGC 2451B. 07422055-3833429 bears similarities with the cases discussed in relation with NGC 2264, namely a very broad CCF (base width of about 400 km s^-1), a main peak at 105 km s^-1, well offset with respect to the cluster velocity, and another bump at 20 km s^-1, close to the cluster velocity. The spectrum seems to show Hα emission. This star has been discarded from the final list.

NGC 2516.

07593671-6021483 is probably a genuine SB2, with the peaks (22 and 50 km s^-1) centred on the cluster velocity (23.6 km s^-1). 07594121-6109251 has a broad CCF (base width 100 km s^-1), with two bumps (−23 and −5 km s^-1) not centred on the cluster velocity and may be contaminated by nebular absorption lines.

NGC 2547.

08081564-4908244 is a genuine SB2, but is probably not a member of NGC 2547, since the component velocities (52 and 122 km s^-1) do not encompass the cluster velocity (15.7 km s^-1).

NGC 3293.

10361099-5814310, classified as 2020C, is probably a δ Scu star (the recommended parameters are T_eff = 8985 K, log g = 4.01) and shows emission in Hα. Rather than SB2 systems, 10353288-5813498 and 10353397-5813178 are rapidly rotating (and probably pulsating) stars (because their CCFs are distorted). They are pre-main sequence star candidates (Delgado et al. 2007).

NGC 3532.

The source 11085927-5849560 is identified for the first time as an SB2 candidate.

NGC 6005.

Three SB2 candidates have only been observed with the HR9B setup (around the Mg I b triplet) where it is difficult to assess if the spectra are composite or not. 15555518-5725349 shows a broad CCF due to Hα in HR15N with a main peak at −69 km s^-1, and a bump at −27 km s^-1 close to the cluster velocity.

NGC 6530.

Numerous spurious detections of SB candidates are due to the presence of nebular lines in emission in HR15. Also, there are thin and deep absorption lines around 6678, 6715, and 6730 Å. Nebular lines are present around 6717 and 6730 Å in HR15 and have led to some reduction issues since negative fluxes are observed at these wavelengths in some stars (18044420-2415380, 18045889-2415261, 18043887-2427164). Their associated CCFs are thus not reliable and have been discarded from the final list. There is strong and deep Hα absorption in several stars. Surprisingly, many discarded SB2 components have velocities close to −60 km s^-1. This raises the question of the presence of another possible velocity component for that cluster.

NGC 6633.

18280622+0642252 (NGC 6633 110, V = 10.1, A3) is an interesting case of a fast rotator which could be a δ Scu-type star according to the iDR4 recommended parameters (T_eff = 9600 K, log g = 4.80, and solar metallicity). Only the upper grating of U580 is available showing very thin and deep absorption lines superimposed on the less deep and rotationally broadened Na I D doublet probably caused by nebular line contamination.

NGC 6705.

The composite spectra and the associated CCFs of one of the five SB2 A candidates are presented in Fig. 16. This is an illustration of a very favourable case: 18503230-0617112 (NGC 6705 1936) has been observed in eight setups and shows a two-component CCF in all of them. 18511434-0617090 has four observations with the HR15N (Hα) setup. In all cases, the main peak is around 33 km s^-1, thus close the cluster velocity, whereas the CCF exhibits a secondary bump around −50 km s^-1. The contrast of that bump is variable, however suggesting that its origin may be related to stellar variability (but B−V ~ 1, suggesting that the star is a red giant, and Hα variability is not expected; Cantat-Gaudin et al. 2014). On the contrary, if the system is an SB2, its kinematics is not compatible with membership in the cluster.

NGC 6752.

19105940-5957059 is the star A13 in Moni Bidin et al. (2006) which has not been detected as binary. The CCFs show clearly the presence of two peaks (flagged 2020 B). Twenty-three observations covering more than 1500 days are available, but we unsuccessfully tried to fit an orbit. Indeed the radial velocities of the components remain constant within few km s^-1. Moreover, the star is located in a very dense region of this globular cluster and we conclude that this is an “optical” SB2.

NGC 6802.

Two SB2 and one SB3 candidates has been found in this cluster.

Trumpler 14.

A large number of false SB2 detections were identified, due to the presence of very strong nebular lines and reduction issues in HR15N where spectra have Hα with negative flux and core emission (see Sect. 4.7 and figures within for discussion). Nebular emission in Trumpler 14 and more generally in the Carina nebula has been investigated in detail by Damiani et al. (2016).

Trumpler 20.

12384378-6037077 has one component located at the cluster velocity. The unique CCF of 12393764-6038190 could either be indicative of a rapidly rotating star with some asymmetries in the line profile or of a cluster member (−36 km s^-1) blended with a non-member (−77 km s^-1). The same remark holds true for 12393362-6041446. The secondary peak of 12391767-6036083 is probably from a non-member. The main peak in the CCF of 12391992-6029552 at −3.4 km s^-1 is probably from a non-member.

Trumpler 23.

The radial velocity of this cluster has just been assessed at −61.3 ± 1.9 km s^-1 within the GES consortium (Overbeek et al. 2017). Therefore, the SB2 candidate 16004521-5332044 can be considered a member of this cluster.

All Tables

All Figures