A RADIAL VELOCITY STUDY OF COMPOSITE-SPECTRA HOT SUBDWARF STARS WITH THE HOBBY–EBERLY TELESCOPE^*

We estimated the spectral types of the cool companions by comparing the observed HET spectra to high-resolution templates in the ELODIE 3.1 stellar spectral library (Prugniel & Soubiran 2001). Spectral types F0–M9 were considered in luminosity classes IV and V. Using our own IDL code, we compared the best individual spectrum for each target to the template spectra over the wavelength range 5075–5200 Å, which includes important discriminators such as the Mg ib triplet (+MgH) and Fe lines. Before any comparisons were made, we degraded the resolution of the ELODIE templates (R = 42,000) by convolving them with Gaussians to match the resolution of the HET/MRS observations (R = 10,000). We constructed a library of possible "composite-spectrum" models, diluting each ELODIE template spectrum by a flat continuum (representing the contribution from the hot subdwarf) using a range of dilution factors (5%–95%). We then normalized the observed MRS spectra with high-order polynomials using the continuum function in IRAF and shifted each spectrum to its rest frame.

For each object–template pair we subtracted the model from the object and computed the sum of squared residuals (SSR); the adopted best-fitting dilution is the one that minimizes the SSR. We repeated this process for all ELODIE templates and assigned a final spectral type for each cool companion target corresponding to the global SSR minimum. To test the accuracy of our fitting technique, we also classified several twilight-sky spectra and a dozen RV standards with known spectral types (drawn from the Astronomical Almanac and Nidever et al. 2002). Our results matched the known values reasonably well, falling within two subtypes of the published classification 80% of the time. The largest disagreement encountered in our consistency check was four subtypes.

The second column in Table 1 lists our classifications for the target stars, which range from F4 V to G8 V. With 80% confidence, they are accurate to within two subtypes. They are also reasonably consistent with previously published values (e.g., Stark 2005). Although we report all companions as dwarfs, two of the targets, PG 0039+049 and PG 1317+123, are matched equally well with subgiant templates. Classifications could not be made for PG 1206+165 (companion lines too weak) or PG 1629+081 (resolved binary, companion not on fiber). In these cases, we adopt the spectral types reported by Stark (2005).

We used iraf routines to find RVs of the targets. Our primary strategy was to exploit the abundance of spectral features from the cool companion by cross-correlating each spectrum against a twilight-sky (i.e., solar) spectrum obtained on the same night. Before this was done, we created sampling region files for each spectrum that defined acceptable wavelength bins to use for the velocity fitting. In addition to ignoring bins compromised by cosmic rays and sky emission features, we also avoided spectral regions heavily influenced by absorption lines from the hot subdwarf component, namely, the Balmer lines (Hβ and Hγ), He i lines (4472, 4713, 5015, and 5876 Å), and the He ii 4686 Å line. Using only these data samples, we used the iraf continuum task to fit the continuum and normalize it to unity. We then cross-correlated the normalized spectra against the solar spectra, order-by-order, using fxcor and rvcorrect to determine heliocentric RVs.

Figure 1 shows an example target–sky pair and its cross-correlation function. Most targets matched the solar spectrum with a correlation nearing 0.8, although systems with companions deviating significantly from G2 V had lower values. Our fitting process resulted in ∼14 RV measurements per target observation, depending on the exact number of echelle orders used. RV errors were computed from the standard deviation of these values about the mean. We find typical errors of 600–1100 m s⁻¹.

Zoom In Zoom Out Reset image size

In some cases, it was also possible to estimate the sdB star's RV from the available He i lines and, for the hotter stars, He ii 4686 Å. We avoided using Hβ for this purpose since its profile spans almost an entire echelle order; also, its core is blended with the Hβ profile of the cool companion. For each of these targets, after selecting sampling regions around the sdB He lines, we cross-correlated each spectrum against a self-template, chosen to be the highest S/N spectrum obtained for the target. This process only provides relative RV measurements with respect to this template; converting them to absolute velocities requires determining the heliocentric velocity of the template. Fitting simple Gaussians or other symmetric functions to the sdB lines cannot be used for this purpose since the He i lines are asymmetric and would skew the measured centroids. Instead, we cross-correlated the observed profiles against a model sdB spectrum with T_eff = 30,000 K, log g = 5.6, and log N(He)/N(H) = −2. We then adjusted the relative measurements appropriately. Naturally, the RV precision achieved from this process is not ideal, but it was sufficient to detect and measure sdB accelerations in many of the systems. RVs determined from He i lines may show a systematic offset, however, as some of the lines are sensitive to pressure shifts in the high-gravity atmosphere (especially the 5876 Å D₃ line). Additionally, the gravitational redshift difference between the cool companion and the compact hot subdwarf can give the appearance of a ∼2 km s⁻¹ offset. These issues have been encountered before in studies of hot subdwarfs (e.g., Orosz et al. 1997) and also DB WDs (e.g., Koester 1987), where the offset can approach tens of km s⁻¹. We will discuss their effects on our results object-by-object, as needed.

To investigate the zero-point accuracy of the MRS velocities, we cross-correlated all of the observed RV standards against twilight-sky spectra (usually taken on the same night) and compared our results to the Nidever et al. (2002) values. A histogram of the RV differences is described reasonably well by a Gaussian distribution with a dispersion of σ ≈ 340 m s⁻¹; its centroid is within 100 m s⁻¹ of zero. We performed a similar test using all of the available twilight-sky spectra, cross-correlating them against the highest S/N twilight-sky spectrum. In this case, we find a scatter with σ = 250 m s⁻¹. Our RV standard star scatter is roughly consistent with the twilight-sky result, once the additional scatter of ≈150 m s⁻¹ from the Nidever et al. (2002) measurements is added in quadrature.

We fitted each RV curve with a model of the form:

$$\begin{equation} V_{\rm rad}= \gamma + K [\cos (\theta + \omega) + e \cos (\omega)], \end{equation} \tag{ 1 }$$

where ω is the periastron angle, K is the RV semiamplitude, γ is the systemic velocity, and e is the eccentricity. The true anomaly (θ) is defined in terms of the period (P), eccentricity, and time of periastron passage (T₀) through the relation

$$\begin{equation} \tan \left(\frac{\theta }{2}\right) = \sqrt{\frac{1+e}{1-e}} \tan \left(\frac{E}{2}\right) \end{equation} \tag{ 2 }$$

and Kepler's equation

$$\begin{equation} E - e \sin E = M = \frac{2 \pi }{P}(t-T_0), \end{equation} \tag{ 3 }$$

where E is the eccentric anomaly and M is the mean anomaly. We used the idl-based rvlin package⁴ (Wright & Howard 2009), which uses a Levenberg–Marquardt algorithm for the nonlinear parts of this optimization problem. Orbital parameter uncertainties were estimated using a bootstrapping technique driven by rvlin (Wang et al. 2012). We determined the best solution with all parameters left free, and we also fitted the data assuming a circular orbit (e = 0 fixed), with two fewer degrees of freedom.

We also computed the floating-mean periodogram (e.g., Cumming et al. 1999) for each RV data set and investigated χ² as a function of orbital period over the range P = 0.1–10,000 days. We consider the orbital period "solved" if (1) the best-fitting orbit has a χ² value at least 10 times smaller than the second best alias and (2) the uncertainty in P for the best fit is less than 5%. Once the first condition was met, we used a range of initial parameter guesses in RVLIN to be certain we converged on the correct solution, sometimes changing the inputs only slightly.

Determining the eccentricity of long-period sdB binaries is an astrophysical imperative. As already mentioned in Section 1, if the MS companions are responsible for the formation of the sdB (and the sdB progenitor filled its Roche lobe), one would expect their orbits to have e = 0 due to the circularizing nature of the RLOF and CE processes. If, however, the hot subdwarf formed from the merger of a binary originally part of a triple-star system, the cool companion currently observed would have had no part in the formation of the hot subdwarf, in which case the orbit could be noncircular. All of the long-period systems studied to date have sufficient room to accommodate an inner binary, even after allowing for moderate "adiabatic expansion" of the orbit (Eggleton et al. 1989; Debes & Sigurdsson 2002) resulting from mass loss in the sdB formation process.

In order to determine whether eccentric-orbit fits are preferred for our targets, we performed an F-test similar to Lucy & Sweeney (1971) using the statistic defined by

$$\begin{equation} \mathcal {F} = \frac{\big(\chi ^2_{\rm circ} - \chi ^2_{\rm ecc}\big)}{({\rm DOF}_{\rm ecc}-{\rm DOF}_{\rm circ})}\frac{{\rm DOF}_{\rm ecc}}{\chi ^2_{\rm ecc}} \end{equation} \tag{ 4 }$$

where DOF_ecc and DOF_circ are the degrees of freedom in the eccentric- and circular-orbit models, respectively, and χ²_ecc and χ²_circ the chi-squared values of their fits to the data. Under the assumption that $\mathcal {F}$ follows an F-distribution, this statistic provides a quantitative means for deciding whether we can justify fitting the data with two additional degrees of freedom (e, ω) compared to our null hypothesis, a circular-orbit solution. As $\mathcal {F}$ increases, so does the likelihood that one can discredit the null hypothesis. We reject a circular-orbit solution if the probability of obtaining the observed $\mathcal {F}$ (the p-value) falls below our false-rejection probability, which we conservatively set to 10⁻⁴ due to the considerable phase gaps in the RV measurements. If the p-value is above this limit, we do not rule out a circular-orbit solution using the current data.

For those systems for which we can measure the RVs of both components, we can calculate the mass ratio q of the binary. It is straightforward to show that the instantaneous sdB and MS companion RVs at time t_i are related by

$$\begin{equation} V_{\rm MS}(t_i) = \gamma (1+q) - qV_{\rm sdB}(t_i), \end{equation} \tag{ 5 }$$

where q ≡ K_MS/K_sdB = M_sdB/M_MS and γ_sdB = γ_MS has been assumed. Thus, q can be measured from the slope of a linear fit to a plot of RV_MS versus RV_sdB. The resulting q is unaffected by any systematic offsets between the component velocities, such as that caused by different gravitational redshifts. It is also independent of the inclination angle, which remains unknown for all systems studied (none show eclipses). We used the IDLfitexy⁵ function, an orthogonal distance regression routine (procedure from Press et al. 1992), to determine q from the best straight-line fit to the data, taking into account the errors from both RV measurements.

Previous classifications for the cool companion of PG 1104+243 have ranged from F8 V (Stark 2005) to K3.5 V (Allard et al. 1994). Comparisons of the HET/MRS spectra to the ELODIE templates show a best-fitting spectral type near G2 V. The dilution of the companion spectrum is D ≡ L_comp/L_total ≈ 0.5, where the luminosity ratio refers to the 5075–5200 Å spectral region. As noted by Ø12, the hot subdwarf in PG 1104+243 might be described best as an "sdOB" star since its spectrum shows the He ii 4686 Å line in addition to strong H Balmer and He i features. Previous RV measurements have been published by Orosz et al. (1997), who found no significant variations over a three-day timespan. An orbital solution with P = 752 ± 14 days was recently published by Ø12. Combining our data with those of both Ø12 and Orosz et al. (1997), we find an orbital period of P = 758.8 ± 0.13 days, in excellent agreement with the Ø12 result. This value falls frustratingly close to two years, as clearly seen in the RV curve shown in Figure 2. By chance, current observing seasons align with quadrature phases; otherwise, all points would have been obtained near a conjunction, a less-than-ideal situation for measuring the mass ratio from spectroscopy. At the same time, information near conjunction is needed to improve constraints on the geometry of the orbit. Several more orbital cycles must elapse before complete phase coverage can be obtained. At this time, we have no reason to prefer an eccentric orbit over a circular one.

Zoom In Zoom Out Reset image size

The presence of strong He lines from the sdB (5876 and 4686 Å) permits us to measure accelerations for both components. As expected for stars in a binary, the sdB and MS velocity curves are 180° out of phase. The sdB RV curve shows a systemic velocity consistent with that measured from the companion (γ = −15.5 ± 0.03 km s⁻¹). From the relationship between RV_sdB and RV_MS, as seen in Figure 3, we compute a mass ratio of q = 0.68 ± 0.08. Table 5 summarizes the system parameters derived for PG 1104+243. If we assume a companion mass of 1.0 M_☉ and adopt a 10% error to include uncertainties in the spectral classification, we find M_sdB = 0.68 ± 0.11 M_☉, considerably larger than the "canonical" value of ∼0.48 M_☉, which corresponds to helium core ignition at the tip of the RGB. If we instead assume exactly the canonical mass for the hot subdwarf, we would compute a mass of 0.7 ± 0.1 M_☉ for the companion, which implies its spectral type might be significantly later than our result from the composite-spectrum analysis, or that the companion is evolved.

Zoom In Zoom Out Reset image size

Spectral classifications for the hot component in PG 1317+123 have varied over years. It is probably an sdO star as it displays a strong He ii 4686 Å line and a relatively weak He i 4472 Å line. Ulla & Thejll (1998) classify the cool companion as a G8 V, while Ø12 report it to be a "G" star. We find a spectral type near G1 V and dilution D ≈ 0.3 from our comparisons with ELODIE models. Ø12 reported an orbital period of P = 912 ± 52 days. After phasing the HET data with theirs, we find a significantly longer value, P = 1179 ± 12 days. Figure 2 shows the resulting RV curve, while Table 5 summarizes all of the orbital parameters derived. Our results for the systemic velocity and MS RV amplitudes also disagree. As pointed out by Ø12, the Mercator/HERMES measurements show much larger scatter than one would expect given the precision of their data. They briefly suggested γ Doradus-type pulsations or star spots on the cool companion as potential culprits. However, G-type classifications for the cool companion are inconsistent with those of γ Dor stars, which are typically early F or late A dwarfs. The other scenario, a bright spot on the companion's surface, seems more compatible with the observations. Assuming a solar-type dwarf and a star spot with a filling factor of a few percent (Saar & Donahue 1997), the amplitude of the scatter (1–2 km s⁻¹) would imply a rotational period around one week. As with PG 1104+243, we currently have no reason to prefer an eccentric-orbit solution. Even with additional data, constraining the eccentricity will be especially difficult until the rapid velocity variations are better characterized.

He i 5876 Å and He ii 4686 Å features provided RV estimates for the hot subdwarf that are again 180° out of phase with the cool companion. Figure 3 shows a strong linear relationship between RV_sdO and RV_MS, whose slope gives q = 0.42 ± 0.05. A circular orbit fitted to the sdO RVs with the period and phase fixed to the values in Table 5 gives a consistent systemic velocity and a semiamplitude in agreement with the expectation from combining K_MS and q. Assuming a companion mass from with our spectral classification, M_MS = 1.03 ± 0.11 M_☉, gives M_sdO = 0.43 ± 0.07 M_☉.

Our analysis of PG 1338+061's spectrum reveals a G2–G7 V companion (best fit: G4 V with D ≈ 0.2), much earlier than previous classifications (K0 V, Stark 2005; K4.5 V, Allard et al. 1994). Two He i lines (5876 and 4472 Å) were used to measure the sdB RVs, which show the sdB and companion accelerating in anti-phase with one another. The last panel in Figure 2 presents the RV curve for PG 1338+061, combining the HET and Mercator observations. Unlike PG 1104+243 and PG 1317+123, we find significant evidence for an orbit with non-zero eccentricity. Figure 4 shows our best-fitting circular and eccentric fits to the cool companion's RV curve for comparison. Computation of the F-test statistic gives $\mathcal {F}$ = 12.6 with a corresponding p-value of 10⁻⁷, well below our criterion for claiming non-zero eccentricity. We therefore reject the circular-orbit hypothesis and report e = 0.15 ± 0.02 for this system. The best-fitting period in this case is P = 937 ± 10 days, in slight disagreement with Ø12, who assumed a circular-orbit solution and found P = 915 ± 16 days. In Section 6.2, we discuss the implications of non-zero eccentricity on formation scenarios of long-period sdB systems.

Zoom In Zoom Out Reset image size

The relationship between RV_sdB and RV_MS (Figure 3) is fitted well with q = 0.58 ± 0.07. Combining this value with K_MS predicts a velocity semiamplitude consistent with the best fit to the subdwarf RV curve (with period and phase fixed). The difference in the systemic velocities is γ_sdB − γ_MS = −2.1 ± 1.0 km s⁻¹. Assuming a mass of M_MS = 0.97 ± 0.1 M_☉ for the companion, we calculate M_sdB = 0.56 ± 0.09 M_☉ for the hot subdwarf, consistent with the canonical value.

One clear trend appears to emerge from the empirical distribution of orbital periods: closer binaries tend to have late-type MS or WD companions while the long-period systems host F/G/K dwarfs. This apparent dichotomy has been noted before (e.g., Green et al. 2008). Although the number of long-period binaries known is still relatively small, their addition allows us to begin testing the predictions of BPS models in both the short- and long-period regimes for the first time.

The Han et al. (2002, 2003) models reproduce the short-period end of the distribution quite well but predict that the sdB+MS (G/K) binaries should have periods shorter than ∼100 days, the outcome of CE evolution using the "α-formalism," in which orbital energy (and/or thermal energy stored in the envelope) is tapped to expel the envelope. Current observations suggest that the energy-based α-formalism for CE evolution in its simplest forms does not account for the long-period sdB binaries. We note that Han et al. (2002, 2003) deliberately excluded the GK-effect binaries when choosing their favored BPS models, focusing on the short-period M dwarf and WD systems (the only systems with measured periods available at the time). Thus, it is not surprising their model predictions are inconsistent with the current period distribution of sdB+F/G/KV binaries. Clausen et al. (2012) encountered a similar problem when modeling the sdB population using the α-formalism. This failing is not limited to the sdB binaries and has been discussed in other contexts by Webbink (2008), who points out the need in some cases for (1) additional energy sources to augment or substitute for some of the energy extracted from the orbit in expelling the envelope, or (2) a reduction in the binding energy of the envelope (possibly through more extensive mass loss prior to the CE phase), in order that orbital shrinkage not be so drastic.

An alternative scenario proposed by Nelemans et al. (2000, 2001) appears, at first glance, to be more successful. In this model, the long-period systems are the product of CE evolution, with a prescribed fractional angular momentum loss that is proportional to the fractional mass lost when the envelope is ejected. The current distribution of long-period systems overlaps surprisingly well with the densest region of the predicted population when the angular momentum parameter γ = 1.5 (Figure 2 of Nelemans 2010). Current observations hint at a dearth of binaries in the range P ≈ 30–110 days, and in this context we note that a period histogram drawn "by eye" from the sdB+MS population (with companion masses ∼1 M_☉) predicted by Nelemans (2010) shows a local minimum in the number of binaries with periods from ≈30 to 90 days. It would be well not to read too much into this apparent agreement, since the model is advertised as tentative, is not well matched to the PG sample in terms of its flux limit, and the transformation of sdB bolometric luminosities to V magnitudes is not described. In addition, the empirical diagram should be more thoroughly populated and the observational biases accounted for, before the significance of the gap is assessed.

It is interesting to consider the question of eccentricity of the observed wide orbits. PG 1338+061 (this paper) and possibly PG 1018−047 (Deca et al. 2012) have significant eccentricities. We also note that eccentric orbits have not been excluded for PG 1104+243 and PG 1317+123 based on the present data. If wide sdB binaries typically prove to have eccentric orbits, this would suggest one of three possibilities. Perhaps (1) these systems are or were hierarchical triple systems in which the sdB star was made in the inner binary, possibly by merger (see, e.g., Clausen & Wade 2011, who describe a "hydrogen merger channel" that leads to sdB's with canonical masses). We note that in PG 1338+061, the minimum periastron distance of 2 AU in the current orbit gives ample room to have once accommodated an inner system. If, however, the systems have always been binaries, then either (2) circularization was not achieved because the sdB progenitor never filled its Roche lobe or (3) eccentricity was pumped into the system via one of several possible mechanisms.

In the first case, the current long-period eccentric orbits are not directly comparable to the predictions of BPS models since the observed cool companions did not participate directly in the formation of the hot subdwarfs. They nevertheless provide some dynamical constraints on "binary" formation scenarios, since they limit the size of the inner orbit in the hierarchical triple. Studies of the sdB rotational velocities in these systems might shed light on whether a merger included two He-core WDs (fast rotation expected: Gourgouliatos & Jeffery 2006; although if there is a giant stage between merger and the sdB phase the rotation may be slowed; see Zhang & Jeffery 2012) or an He WD and an MS star (slow rotation, owing to angular momentum shedding near the tip of the RGB; Clausen & Wade 2011).

In the case where these systems have always been binaries, we suppose that a strong, tidally enhanced wind along the lines suggested by Tout & Eggleton (1988) can occasionally result in the complete loss of the envelope, just as the core of the sdB progenitor reaches the critical mass (or near it; see D'Cruz et al. 1996) to ignite He burning, without ever filling its Roche lobe. Some expansion of the orbit would occur owing to the wind mass loss, and the orbital eccentricity could in large part be preserved. If this timing coincidence does not occur, the extra mass loss can at least reduce the envelope binding energy and, after a CE episode, leave a (circular) orbit at longer period than would otherwise be attained, or even reduce the mass ratio to the point where RLOF mass transfer proceeds stably and a CE phase is avoided altogether, as in Tout & Eggleton (1988).

Finally, we consider briefly the case where these systems are in fact binaries in which the sdB progenitor once filled its Roche lobe. Several mechanisms exist that could work to increase the orbital eccentricity, even as mass is being lost or transferred, thereby preventing the binaries from circularizing. Most notably, these include pumping from a circumbinary disk (e.g., Dermine et al. 2012, and references therein) and asymmetric/phase-dependent mass loss (Bonačić Marinović et al. 2008; Soker 2000), both of which have been invoked to explain the large eccentricities observed in some post-asymptotic giant branch (AGB) systems having undergone strong tidal interactions. Long-period sdB binaries are post-RGB systems, and as such, their sdB progenitors had significantly different luminosities, radii, envelope structures, and mass-loss rates than AGB stars. Eccentricity-pumping mechanisms operating on RGB and AGB systems will neither act on the same timescales, nor will they have the same ability to compete with the effects of tides. We regard these scenarios, which have primarily been linked to post-AGB binaries, as uncertain in the case of hot subdwarf binaries; further discussion of them is beyond the scope of this paper. For completeness, we also recognize the action of a distant third body (the Kozai mechanism; Kozai 1962), which we regard as the "ultimate epicycle," as another possible avenue for increasing the eccentricity.