VETTING GALACTIC LEAVITT LAW CALIBRATORS USING RADIAL VELOCITIES: ON THE VARIABILITY, BINARITY, AND POSSIBLE PARALLAX ERROR OF 19 LONG-PERIOD CEPHEIDS

The Cepheid sample investigated here and listed in Table 1 was selected according to several criteria. As the primary goal is to determine parallax accurately using HST/WFC3 spatial scans (Riess et al. 2014; Casertano et al. 2016), the most crucial selection criteria were those centered on the astrometric measurement itself.

An optimal target for high-accuracy spatial scan parallax measurements has

• 1.  
  pulsation period longer than approximately 10 days as this reflects the periods of the predominant group of Cepheids found in other galaxies (due to higher luminosity) (S. L. Hoffmann et al. 2016) and avoids putative nonlinearities of the PLR intervening at 10 days (for differing recent opinions on the matter, see Inno et al. 2013; Bhardwaj et al. 2016; García-Varela et al. 2016);
• 2.  
  mean V-band magnitude fainter than 7.5 to avoid saturation (H-band > 5 mag);
• 3.  
  at least five, ideally more than ten, reference stars within six magnitudes of the Cepheid for which trails would be recorded simultaneously;
• 4.  
  an expected distance of d ≲ 3 kpc so that parallax can be determined to better than 10% for each individual Cepheid for a nominal parallax uncertainty of 20–40 μarcsec;
• 5.  
  no known companion star with ${P}_{\mathrm{orb}}$ on the order of the sparsely sampled HST/WFC3 spatial scan observations (typically five epochs);
• 6.  
  extinction (A_H ≲ 0.5 mag) to avoid excessive uncertainty in the inferred absolute magnitudes.

Although binaries were not strictly excluded from the sample (see, e.g., YZ Carinae below or visual binaries), we caution that the present sample is subject to selection effects concerning binarity and should not be considered random in this regard. Therefore, we stress that this sample alone should not be used to infer the properties of binary fractions unless this occurs in conjunction with further observational data capable of eliminating or reducing such selection effects (Evans et al. 2013).

Table 1.  Basic Information on the Sample of Cepheids Discussed Here

Cepheid	HD	R.A.(J2000)	Decl.(J2000)	$\langle {m}_{V}\rangle$	NCor	NHam	NHer	${\rm{\Delta }}{t}_{\mathrm{obs}}$	References
		(h:m:s)	(d:m:s)	(mag)				(years)
SY Aur	277622	05:12:39.20	42:49:54	9.1	0	78	31	2.6	⋯
SS CMa	HIP 36088	07:26:07.20	−25:15:26	9.9	24	42	14	3.0	(1), (2), (3)
VY Car	93203	10:44:32.70	−57:33:55	7.5	83	0	0	5.0	(1), (2), (4), (5), (6)
XY Car	308149	11:02:16.10	−64:15:46	9.3	70	0	0	1.9	⋯
XZ Car	305996	11:04:13.50	−60:58:48	8.6	118	0	0	4.2	Section 4.3.1
YZ Car	90912	10:28:16.80	−59:21:01	8.7	28	0	0	2.4	(7), (8), Section 4.2
AQ Car	89991	10:21:23.00	−61:04:27	8.9	59	0	0	1.9	Section 4.3.1
HW Car	92490	10:39:20.30	−61:09:09	9.2	68	0	0	5.0	⋯
DD Cas	HIP 118122	23:57:35.00	62:43:06	9.9	0	78	23	2.3	(9), (10)
KN Cen	HIP 66383	13:36:36.90	−64:33:30	9.9	70	0	0	2.1	(9), (11), (12), (13), (14), (15)
SZ Cyg	196018	20:32:54.30	46:36:05	9.4	0	66	41	2.3	(16), (17), Section 4.3.1
CD Cyg	227463	20:04:26.60	34:06:44	9.0	0	66	45	2.3	Section 4.3.1
VX Per	236948	02:07:48.50	58:26:37	9.4	0	78	40	4.0	⋯
X Pup	60266	07:32:47.00	−20:54:35	8.6	48	21	15	1.4	(18)
AQ Pup	65589	07:58:22.10	−29:07:48	8.7	51	41	6	3.0	(10), (19), (20)
WZ Sgr	167660	18:16:59.70	−19:04:33	8.1	48	6	29	4.0	(5), (15), (21), (22), (23)
RY Sco	162102	17:50:52.30	−33:42:20	8.0	52	0	0	2.1	(1), (22)
Z Sct	172902	18:42:57.30	−05:49:15	9.6	41	48	41	3.1	⋯
S Vul	338867	19:48:23.80	27:17:11	9.1	12	12	37	1.9	⋯

Note. Basic information on the sample of Cepheids discussed. Hipparcos identifiers are given where no HD number was available. Coordinates and mean magnitudes are based on the information from the GCVS (http://www.sai.msu.su/gcvs/gcvs/), average magnitudes are approximate. The number of observations obtained with Coralie, Hamilton, and Hermes are listed together with the total temporal baseline ${\rm{\Delta }}{t}_{\mathrm{obs}}$ of our new observations. Typically, we obtained three observations per pointing with the Hamilton spectrograph. The total number of new observations made available here is 1630.

References. References for work previously discussing the binarity or cluster membership of these objects are as follows (see also the binary Cepheids database by Szabados 2003): (1) Evans & Udalski (1994), (2) Szabados (1996), (3) Casertano et al. (2016), (4) Turner (1977), (5) Anderson et al. (2013), (6) Perryman & ESA (1997), (7) Coulson (1983), (8) Petterson et al. (2004), (9) Madore (1977), (10) Madore & Fernie (1980), (11) Walraven et al. (1964), (12) Lloyd Evans (1968), (13) Stobie (1970), (14) Pel (1978), (15) Szabados (1989), (16) Kurochkin (1966), (17) Szabados (1991), (18) Szabados et al. (2012), (19) Fernie et al. (1966), (20) Vinko (1991), (21) Bersier (2002), (22) Proust et al. (1981), (23) Turner et al. (1993). Sections discussing individual Cepheids in more detail are indicated in the column References.

A machine-readable version of the table is available.

Download table as:  DataTypeset image

We have secured time-series observations from three different high-resolution echelle spectrographs: Coralie (R ∼ 60,000) at the Swiss 1.2 m Euler telescope located at La Silla Observatory, Chile; Hermes (R ∼ 85,000) at the Flemish 1.2 m Mercator telescope⁸ located at the Roque de los Muchachos Observatory on La Palma, Canary Islands, Spain; Hamilton (R ∼ 60,000) at the 3 m Shane telescope located at Lick Observatory, California, USA.

Coralie and Hermes spectra were reduced using the dedicated pipelines available on site. Hamilton spectra were reduced using standard IRAF routines. All spectra were bias-corrected and flat-fielded, and cosmic ray hits were removed. ThAr (Coralie, Hermes) and TiAr (Pakhomov & Zhao 2013, Hamilton) lamps were used for wavelength calibration.

All RVs presented here were determined using the cross-correlation technique (Baranne et al. 1996; Pepe et al. 2002) using a numerical mask representative of a solar spectral type (G2 mask). An example of the data available is presented in Table 2.

Coralie RVs were corrected for temporal variations in the wavelength calibration using ThAr reference spectra that are interlaced with the science orders on the detector. Hermes spectra were corrected for such variations using frequent re-calibration of the wavelength solution and a model for estimating RV zero-point changes associated with changes in air pressure (as done in Anderson et al. 2015). Hamilton spectra are the most affected by temporal variations in the instrumental zero-point, which dominate the uncertainty for the associated RVs presented here. We use stable RV standard stars to track the RV variation due to intra-night changes of the wavelength solution and determine appropriate corrections for science exposures by interpolating the time sequence of offsets determined.

The precision of Coralie and Hermes measurements is typically on the order of 10–30 m s⁻¹, depending on the signal-to-noise ratio achieved. At this level of precision, the instrumental zero-points of Coralie and Hermes are compatible without adjustments. Hamilton RVs are significantly less precise due to the unstable zero-point; here we adopt 200 m s⁻¹ as a typical uncertainty for Hamilton RVs. This value includes the uncertainty associated with tracking the nightly zero-point variations using standard stars as well as (smaller) RV zero-point differences among instruments.

The time of observation for all newly observed spectra are given as solar system barycentric Julian dates minus 2,400,000 and all associated RV measurements are relative to the solar system barycenter.

Table 2.  Example RV Data Obtained for This Program

Cepheid	BJD–2.4M	${\phi }_{\mathrm{puls}}$	vr	σ(vr)	Instrument
	(days)		(km s−1)	(km s−1)
SY Aur	56402.68783	0.0580	−10.248	0.2	Hamilton
SY Aur	56519.01162	0.5236	6.516	0.2	Hamilton
SY Aur	56519.01798	0.5242	6.500	0.2	Hamilton
SY Aur	56581.00590	0.6341	10.463	0.2	Hamilton
SY Aur	56581.01220	0.6347	10.463	0.2	Hamilton
SY Aur	56581.01850	0.6354	10.496	0.2	Hamilton
SY Aur	56581.91266	0.7235	7.004	0.2	Hamilton
SY Aur	56581.91897	0.7241	6.930	0.2	Hamilton
SY Aur	56581.92527	0.7247	6.837	0.2	Hamilton
SY Aur	56609.97389	0.4894	4.237	0.2	Hamilton
...
S Vul	57500.899168	0.7406	15.663	0.032	Coralie
S Vul	57504.910489	0.7991	14.473	0.073	Coralie
S Vul	57507.895237	0.8426	10.209	0.035	Coralie
S Vul	57508.895744	0.8572	7.892	0.029	Coralie
S Vul	57511.898062	0.9010	−0.430	0.032	Coralie
S Vul	57526.868106	0.1192	−10.978	0.016	Coralie
S Vul	57528.902885	0.1489	−9.989	0.022	Coralie
S Vul	57529.903860	0.1635	−9.466	0.025	Coralie
S Vul	57534.837524	0.2105	−6.426	0.014	Coralie
S Vul	57536.886856	0.2399	−5.069	0.024	Coralie

Note. The full data set comprising 1630 observations of the 19 Cepheids obtained with the three spectrographs is published online via the journal and the CDS. Pulsation phase is defined such that ${\phi }_{\mathrm{puls}}$ = 0 at minimal RV and is computed using ephemerides listed in Table 3. Dates and RVs are relative to the solar system barycenter.

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

The observed RV curve of a (binary) Cepheid is a superposition of the systemic RV relative to the solar system barycenter, v_γ, the pulsational variability, ${v}_{r,\mathrm{puls}}$, and the orbital motion of the Cepheid relative to the center of gravity of the binary system, ${v}_{r,\mathrm{orb}}$. Thus,

$$\begin{eqnarray}&&{v}_{r}(t)={v}_{\gamma }+{v}_{r,\mathrm{puls}}+{v}_{r,\mathrm{orb}}.\end{eqnarray} \tag{ 1 }$$

We model the pulsational variability as a Fourier series with an appropriate (fixed) number of harmonics, N_FS, which is adopted during a preliminary inspection of the available RV data, see Table 3. The Fourier model of the pulsation is computed as:

$$\begin{eqnarray}&&{v}_{r,\mathrm{puls}}(t)=\displaystyle \sum _{i=1}^{N}{a}_{i}\sin 2\pi {\phi }_{\mathrm{puls}}+{b}_{i}\cos 2\pi {\phi }_{\mathrm{puls}},\end{eqnarray} \tag{ 2 }$$

with pulsation phase ${\phi }_{\mathrm{puls}}=(t-E)/{P}_{\mathrm{puls}}$, where t is time in Julian days, E is the reference epoch, and ${P}_{\mathrm{puls}}$ is determined by minimizing the internal scatter of our RV data using as starting point a reference value from the General Catalog of Variable Stars (Samus et al. 2009). We here employ a definition of the epoch E so that ϕ = 0.0 coincides with minimal RV near the mean (solar system) barycentric JD of the data considered. This choice of phase zero-point is arbitrary and not of particular importance to this work. ϕ = 0.0 is expected to be close to a time of maximum light although the values of E presented here are not necessarily comparable to times of maximum light determined from light curves. Similarly, the values of ${P}_{\mathrm{puls}}$ that minimize scatter of the present RV data are close, albeit not necessarily identical, to ${P}_{\mathrm{puls}}$ listed in the literature, depending on the number of harmonics used for the fit and how well the pulsations repeat over time. One exception to this procedure is the case of YZ Carinae, where the strong orbital RV signal complicates the determination of ${P}_{\mathrm{puls}}$ based on RV data alone. We therefore used V-band photometric data from the All Sky Automated Survey (ASAS; Pojmanski 2002) to derive a new best-fit ${P}_{\mathrm{puls}}$ and adopted this value for the RV modeling.

Table 3.  Results from RV Curve Modeling for All Program Cepheids

Cepheid	${P}_{\mathrm{puls}}$	E	NFS	vγ	Ap2p	A1	R21	ϕ21	R31	ϕ31	rms	${\sigma }_{{v}_{\gamma }}$
	(days)	BJD–2.4M		(km s−1)	(km s−1)	(km s−1)					(km s−1)	(km s−1)
SY Aur	10.145458	56990.536188	9	−3.172	22.403	10.524	0.39	4.33	0.09	6.99	0.216	0.051
SS CMa	12.352828	57062.269144	7	77.085	38.345	16.801	0.24	4.80	0.17	1.62	0.274	0.110
VY Car	18.882696	56502.102998	9	1.637	57.358	23.973	0.29	3.01	0.07	5.27	0.438	0.281
XY Car	12.436275	57145.616572	11	−5.698	48.171	20.987	0.03	3.61	0.13	2.78	0.176	0.161
XZ Car	16.652208	56554.510768	13	5.671	56.299	23.788	0.29	3.01	0.07	5.22	0.274	0.150
YZ Cara	18.1676	51928.9358	8	0.844	29.692	14.078	0.08	3.14	0.04	2.19	0.037	0.063
AQ Car	9.769452	57105.87848	7	1.577	32.932	13.976	0.30	5.08	0.17	1.75	0.167	0.121
HW Car	9.199135	56727.566552	7	−13.035	19.284	9.018	0.22	5.05	0.08	1.74	0.128	0.075
DD Cas	9.812156	56871.673758	7	−69.453	33.913	14.486	0.24	5.24	0.17	2.01	0.163	0.039
KN Cen	34.018969	57135.280429	9	−42.217	50.014	22.343	0.33	3.00	0.21	5.93	0.821	0.514
SZ Cyg	15.11133	56921.702556	9	−9.976	51.060	21.949	0.24	2.94	0.06	4.73	0.255	0.051
CD Cyg	17.076041	56946.237797	11	−8.709	58.892	24.724	0.27	3.01	0.05	5.00	0.352	0.083
VX Per	10.882827	56819.055602	7	−35.037	30.790	13.471	0.47	4.45	0.16	7.41	0.224	0.051
X Pup	25.959165	57262.157494	11	70.970	57.069	25.403	0.36	3.00	0.17	5.77	0.430	0.162
AQ Pup	30.182036	57123.479592	9	60.798	59.522	26.096	0.32	2.98	0.14	5.67	0.777	0.375
WZ Sgr	21.850992	57052.286612	11	−17.088	54.919	23.601	0.32	3.05	0.10	5.55	0.370	0.185
RY Sco	20.322084	57172.826932	7	−18.653	34.768	16.375	0.16	2.93	0.01	4.79	0.330	0.294
Z Sct	12.901867	56956.219712	7	29.924	52.156	22.614	0.06	4.09	0.14	2.46	0.590	0.120
S Vul	69.653841	57241.56046	5	1.137	29.766	12.962	0.35	3.12	0.15	6.24	0.070	0.033

Note. Pulsation periods and epochs of minimal RV determined for the number of harmonics indicated (N_FS). Mean velocity v_γ, peak-to-peak amplitude A_p2p, first harmonic amplitude A₁, Fourier amplitude and phase ratios R₂₁, ϕ₂₁, R₃₁, ϕ₃₁, fit rms, and standard mean error on v_γ.

^aYZ Carinae is a spectroscopic binary; see Table 5 and Section 4.2.

A machine-readable version of the table is available.

Download table as:  DataTypeset image

In general, RV orbital motion is modeled as a Keplerian with semi-amplitude K, eccentricity e, argument of periastron ω, and the true anomaly θ (see e.g., Hilditch 2001):

$$\begin{eqnarray}&&{v}_{r,\mathrm{orb}}=K[\cos (\omega +\theta )+e\cos \omega ].\end{eqnarray} \tag{ 3 }$$

However, the majority of Cepheids considered here do not exhibit evidence of orbital motion, and indeed one of the primary aims of this work is to set upper limits on undetected companions. We thus assume zero eccentricity unless required (and explicitly stated). This simplifies Equation (3) to an ordinary sinusoid with amplitude and phase, which we here model as:

$$\begin{eqnarray}&&{v}_{r,\mathrm{orb},e=0}={a}_{\mathrm{orb}}\sin 2\pi {\phi }_{\mathrm{orb}}+{b}_{\mathrm{orb}}\cos 2\pi {\phi }_{\mathrm{orb}},\end{eqnarray} \tag{ 4 }$$

with orbital semi-amplitude $K=\sqrt{{a}_{\mathrm{orb}}^{2}+{b}_{\mathrm{orb}}^{2}}$. Orbital period and semi-amplitude then yield the projected semimajor axis $a\sin i$ of the Cepheid's orbit around the common center of gravity:

$$\begin{eqnarray}&&a\sin i\,[\mathrm{au}]=9.192\times {10}^{-5}\cdot {\rm{K}}\,[\mathrm{km}\,{{\rm{s}}}^{-1}]\cdot {{\rm{P}}}_{\mathrm{orb}}\,[\mathrm{days}].\end{eqnarray} \tag{ 5 }$$

Highly eccentric or very long orbital period (${P}_{\mathrm{orb}}$ ≫ 5 years) systems may remain undetected by our RV measurements, depending on the geometry and which part of the orbit would be sampled by the observations. To this end, we also inspect the long-term stability of v_γ using a combination of our new data with published RVs from the literature, see Section 4.3. However, companions on such very long period orbits (${P}_{\mathrm{orb}}$ > 10 years) are not likely to affect the HST/WFC3 parallax measurements.

In the following subsections, we discuss the pulsational RV modeling of the program Cepheids. The search for spectroscopic binarity and estimation of parallax error due to companions is presented in Section 4 below.

Figure 1 presents the data for 18 of the program stars (YZ Car is described separately in Section 4.2) together with the fitted pulsation model. Figure 2 shows the corresponding residuals as a function of observation date. The figures illustrate that most Cepheids have very well-sampled RV curves, although a few cases could benefit from better phase-sampling. This includes in particular S Vul, which is observationally challenging due to its extremely long and unstable pulsation period (Makarenko 1978; Mahmoud & Szabados 1980). We here find a best-fit period of ${P}_{\mathrm{puls}}$ = 69.7 days, which is approximately 1.4% longer than previously reported values, compared to the range of periods in the literature (67.3–68.7 days). Nevertheless, our observations do not sample the complete pulsation curve (in particular the minimum RV), introducing a substantial systematic period uncertainty of ∼1 days.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The spectroscopic binarity of Cepheids is usually determined by investigating the long-term stability of the pulsation-averaged velocity v_γ. Specifically, all known Cepheid binaries have ${P}_{\mathrm{orb}}$ > 1 years and K > 1 km s⁻¹ (see Szabados 2003), possible smaller amplitude companions being masked by RV zero-point offsets among different instruments on the order of a few hundred m s⁻¹ (see Evans et al. 2015b) or noise intrinsic to Cepheid RV variability. The consistently flat residuals shown in Figure 2 thus provide no indication for spectroscopic binarity, leaving only the possibility of very low-amplitude (K ≲ 1 km s⁻¹) or very long timescale (${P}_{\mathrm{orb}}$ ≫ 5 years) orbital motion. These possibilities are investigated in detail in Section 4.

Table 3 lists the results of the pulsational modeling for all program Cepheids. Specifically, it includes best-fit pulsation periods and epochs, number of harmonics used for the Fourier series N_FS, systemic velocity v_γ (which will be of use for Gaia, see de Bruijne & Eilers 2012), peak-to-peak amplitude A_p2p, amplitude of the first harmonic A₁, Fourier ratios R₂₁, R₃₁, ϕ₂₁, and ϕ₃₁ (Simon & Lee 1981), fit rms and uncertainty on v_γ. Amplitude and phase of the ith harmonic are defined as ${A}_{i}=\sqrt{{a}_{i}^{2}+{b}_{i}^{2}}$ and $\tan {\phi }_{i}={b}_{i}/{a}_{i}$, and are computed using the coefficients obtained from Equation (2). Amplitude ratios among harmonics are defined as R_i1 = A_i/A₁, and phase ratios as ${\phi }_{i1}={\phi }_{i}-i\cdot {\phi }_{1}$. Figure 3 illustrates these results and their dependence on logarithmic ${P}_{\mathrm{puls}}$, ignoring S Vul for which the available RV data were insufficient to reliably determine these parameters.

Zoom In Zoom Out Reset image size

We find a dependence of Fourier amplitude and phase ratios on ${P}_{\mathrm{puls}}$ in broad agreement with previous observational (Kovacs et al. 1990) and simulation-based results (Aikawa & Antonello 2000). In particular, we recover the general morphology of increasing RV amplitudes that flatten off around 17 days, as well as the associated steep decline in ϕ₂₁. The period distribution of our program Cepheids nicely complements the sample presented by Kovacs et al. (1990), doubling the number of Cepheids in the ${P}_{\mathrm{puls}}$ range upward of 10 days. These parameters will be useful for hydrodynamical modeling of Cepheid variability, although such applications are outside the scope of this work. Here, we use these parameters to show that most Cepheids exhibit the RV variability behavior expected for their ${P}_{\mathrm{puls}}$.

We note that the RV amplitudes of RY Sco and YZ Car are outliers from the overall trend indicated by the other stars in the sample. Additionally, visual inspection of the data reveals that the (pulsational part of the) RV curve is more sinusoidal than that of other Cepheids; see Figure 1 and Section 4.2 for YZ Car. This simple RV curve shape is quantified as low-amplitude ratios between the first three harmonics; see the parameters R₂₁ and R₃₁ in Figure 3. However, it is known from photometric studies that light curve amplitudes can vary considerably at fixed ${P}_{\mathrm{puls}}$ according to the pulsation-average temperature of the Cepheid, i.e., its position in the instability strip (e.g., Antonello & Morelli 1996; Sandage et al. 2009; Kanbur et al. 2010). Inspection of the ASAS lightcurves (Pojmanski 2002) of RY Sco and YZ Car shows that both exhibit very similar, saw-tooth-shaped V-band variability with only a very minor bump-like feature near minimum light and similar peak-to-peak amplitude of ∼0.8 mag. Since other Cepheids in this ${P}_{\mathrm{puls}}$ range exhibit stronger bump features in their light curves and larger RV amplitudes, this suggests a connection between the low RV amplitudes of RY Sco and YZ Car and the weak bumps. We consider a detailed investigation of the connection between RV and light curve shapes to be out of the scope of this work. In the near future, our parallax measurements will help to clarify whether these differences in light and RV curve shapes among Cepheids are related to differences in luminosity.

Pulsation period changes due to secular evolution—i.e., linear variations on the order of 10–100 s yr⁻¹ for solar metallicity Cepheids in the period range 10–50 days (e.g., Anderson et al. 2016b)—are generally not an issue over the less than five year temporal baseline of our observations. However, Cepheids are also known to exhibit nonlinear pulsation period changes over shorter timescales, and this effect is particularly noticeable for ${P}_{\mathrm{puls}}$ ≳ 20 days. (e.g., Szabados 1989, 1991; Berdnikov et al. 2000, 2009; Poleski 2008; Anderson 2014). A well-known example is RS Pup (${P}_{\mathrm{puls}}$ ∼ 42 days), whose nonlinear ${P}_{\mathrm{puls}}$ variations can result in phase offsets of up to 20% over the course of 20 years (Berdnikov et al. 2009). For short-period overtone Cepheids, period fluctuations on similar timescales have been found using high-cadence photometry from Kepler and MOST (Derekas et al. 2012; Evans et al. 2015a).

The peak-to-peak RV amplitudes presented here are on the order of 20–60 km s⁻¹, see Table 3. Using instruments featuring extreme long-term instrumental stability and high RV precision on the order of a few m s⁻¹, it is now possible to detect pulsation irregularities on the order of 0.01% with confidence. This has led to the discovery of RV curve modulation (Anderson 2014), which is particularly erratic in long-period (${P}_{\mathrm{puls}}$ > 10 days) Cepheids, where cycle-to-cycle variations are found. For example, RS Pup's RV amplitude varies by approximately 1 km s⁻¹ from one pulsation cycle to the next, and up to 3 km s⁻¹ over the course of one year. Both effects are also present, albeit weaker, in the 35 day Cepheid ℓ Car.

Very dense time-sampling is required to model nonlinear period and amplitude fluctuations on a cycle-to-cycle basis (Anderson et al. 2016a; Anderson 2016). Cepheids with ${P}_{\mathrm{puls}}$ ≳ 20 days are most affected by these difficulties, since nonlinear period fluctuations are strongest for these stars and since achieving good phase-sampling is particularly challenging due to practical constraints such as telescope access, weather, and the Moon.

Modeling Cepheid RV variability using the adopted stable model (Equation (2)) thus fails to account for all (astrophysical) signals present in the data, leading to excess residuals and generally very high values of χ². Since these high χ² values are the result of model inadequacy, it would be incorrect to scale RV uncertainties, which have furthermore been shown to represent an adequate estimation of RV precision in the sense of the ability to reproduce a central value from multiple measurements (Anderson 2013).

The presence of additional signals dominates the reduction in χ² when a further model component is introduced in the fit, such as orbital motion (Equation (3) or (4)). We therefore caution that a detection of spectroscopic binarity should not be claimed based purely on a reduction in χ². Rather, additional (visual) inspection of the data is required and must be weighed against other indicators.

Binarity of variable stars can affect parallax measurements primarily in two ways: (1) via real positional modulation due to orbital motion that is not accounted for by the astrometric modeling and (2) via apparent positional modulation in phase with the variability (e.g., of the Cepheid). Here we are interested primarily in the former effect and how to constrain the related error using RVs. Objects affected by photocenter variations due to variability were denoted as variability induced movers (VIMs; Wielen 1996) in Hipparcos and will be considered in future work, since RVs cannot constrain this effect. It is worth noting that parallax error due to orbital motion can be strong even for low-mass companions that would not lead to VIM-type error; see the example of δ Cep (Anderson et al. 2015) whose companion is at least 5–6 H-band magnitudes fainter than the primary (Gallenne et al. 2016).

Parallax error due to orbital motion (henceforth: parallax error) is expected to be strongest for ${P}_{\mathrm{orb}}$ near one year due to the aliasing between orbital and parallactic motion of the HST spatial scan observations. In this section we estimate the possible hidden impact of binarity on the HST parallax measurements by determining upper limits on the potential astrometric impact of undetected binaries among the program Cepheids. Photometric effects are briefly discussed in Section 5.2. We first determine orbital configurations of spectroscopic binaries not ruled out by the high-precision RV data for a range of input orbital periods covering the range of possible binary orbits (here: ${P}_{\mathrm{orb}}$ > 222 days) up to the temporal baseline of the measurements.⁹ We then constrain the parallax error $\hat{\varpi }$ that could result from modeling the HST astrometry of the Cepheid as a single star if in reality it were a binary.

The temporal baseline of the RV data presented here is ideal for this purpose. The data have been recorded contemporaneously with (within approximately one year of) the HST spatial scan measurements and cover the range of orbital periods where the greatest impact of orbital motion on the parallax measurement is to be expected (longer ${P}_{\mathrm{orb}}$ would primarily affect proper motion). However, while RV data are highly sensitive, they can only measure the line-of-sight component of orbital motion. This limitation is explicitly described in the following.

We model the RV data as a sum of mean velocity, pulsation, and circular orbital motion (Equation (1) with Equation (4)) for a set of input (fixed) orbital periods ${P}_{\mathrm{orb}}$ using the pulsation ephemerides listed in Table 3. For each ${P}_{\mathrm{orb}}$, we obtain a best-fit solution for v_γ and the Fourier coefficients, as well as a semi-amplitude K and orbital phase ${\phi }_{\mathrm{orb}}$. The projected semimajor axis of each best-fit solution is calculated using Equation (5). We use $a\sin i$ in au, since the angular size of the orbit scales with the parallax of the Cepheid.

Figure 4 shows the results obtained using this procedure for the example of VX Persei. The left-hand panel shows the χ² map, which would tend to favor ${P}_{\mathrm{orb}}$ < 1 year with very small $a\sin i$, specified here in [au], i.e., equivalent to the fraction of the star's parallax.

Zoom In Zoom Out Reset image size

To estimate parallax error from these upper limits on spectroscopic binarity, we first compute positional offsets δx(t) due to orbital motion at times t of HST spatial scan observations using the best-fit values of $a\sin i$ and ${\phi }_{\mathrm{orb}}$ determined for each ${P}_{\mathrm{orb}}$. Although RV measurements are blind to inclination, positional offsets due to orbital motion intrinsically depend both on inclination i and the orientation of the line of nodes with respect to the astrometric resolution direction, θ. We therefore compute positional offsets for a two-dimensional grid of inclination and orientation angles for each orbital period ${P}_{\mathrm{orb}}$ used in the RV modeling. Each positional offset is computed as:

$$\begin{eqnarray}&&\delta x(t)=a\cos \phi \cos \theta +a\sin \phi \sin \theta \cos i,\end{eqnarray} \tag{ 6 }$$

where $a=a\sin i/\sin i$ is the semimajor axis of the Cepheid around the center of gravity of the hypothetical binary, and $\phi =2\pi \tfrac{t-E}{{P}_{\mathrm{orb}}}+{\phi }_{\mathrm{orb},0}$ with E the epoch of the RV modeling (see Table 3) and $\tan {\phi }_{\mathrm{orb},0}={B}_{\mathrm{orb}}/{A}_{\mathrm{orb}}$; see Equation (4).

Finally, we estimate the projected parallax error $\hat{\varpi }\sin i$ using the set of positional offsets computed for each orbital period using a least-squares procedure that takes into account the exact times t and parallax factors π_f  of each HST observation. Hence, we compute $\hat{\varpi }({P}_{\mathrm{orb}},i,\theta )$ using the δx(t) for all HST observations. As i and θ are unconstrained by the RV data, and since we are interested in conservative upper limits, we adopt the maximal unsigned $\hat{\varpi }$ for each ${P}_{\mathrm{orb}}$ and multiply this value by the sine of the inclination for which it was computed, i.e.,

$$\begin{eqnarray}&&\hat{\varpi }\sin i({P}_{\mathrm{orb}})=\max (| \hat{\varpi }({P}_{\mathrm{orb}},i,\theta )| )\cdot \sin {i}_{\max (| \hat{\varpi }| )}\end{eqnarray} \tag{ 7 }$$

Note that we here use the absolute value $| \hat{\varpi }({P}_{\mathrm{orb}},i,\theta )|$ to estimate the unsigned projected parallax error, since RV data do not constrain i and θ. Depending on the configuration, $\hat{\varpi }\sin i({P}_{\mathrm{orb}})$ could be positive or negative, resulting in an over- or underestimated parallax.

The quantity $\hat{\varpi }\sin i$ represents an upper limit in the sense that it reflects the maximal unsigned parallax error for a given orbital period. However, it preserves the notion that the basis for this upper limit, the modeling of RV data, cannot constrain inclination. We further note that $1/\sin i\lt 3$ for $i\gt 19.5\,\deg$ (94% of possible inclinations) and <10 for i > 5.7 deg (99.5%).

Table 4 lists the results thus obtained for all program Cepheids excluding YZ Car, whose orbit is updated in Section 4.2 below. For each Cepheid, we provide information for (1) the solution offering the weakest constraint on possible parallax error, and( 2) for the solution with minimal χ².

Table 4.  Upper Limits on $\hat{\varpi }\sin i$, the Maximum Unsigned Projected Parallax Error Estimated from Upper Limits on Spectroscopic Binarity

			Greatest Impact on Parallax					Overall Minimum χ2
Cepheid	Δt	DOF	${P}_{\mathrm{orb}}$	K	$a\sin i$	$\hat{\varpi }\sin i$	χ2	${P}_{\mathrm{orb}}$	K	$a\sin i$	$\hat{\varpi }\sin i$	χ2
	(years)		(years)	(km s−1)	(% ϖ)	(% ϖ)		(years)	(km s−1)	(% ϖ)	(% ϖ)
SY Aur	2.60	87	2.60	0.186 ± 0.001	1.63	0.0	125	2.07	0.179 ± 0.001	1.24	0.0	107
SS CMa	2.96	62	0.76	0.191 ± 0.022	0.49	0.2	867	0.77	0.193 ± 0.021	0.50	0.2	866
VY Car	5.01	61	0.95	0.242 ± 0.028	0.78	0.9	20618	0.64	0.103 ± 0.008	0.22	0.1	20207
XY Car	1.91	44	1.07	0.421 ± 0.025	1.52	1.7	3238	1.09	0.43 ± 0.029	1.58	1.6	3219
XZ Car	4.23	88	0.99	0.185 ± 0.005	0.61	0.8	10928	2.48	0.197 ± 0.001	1.64	0.6	8556
AQ Car	1.91	41	1.00	0.476 ± 0.012	1.60	1.9	709	1.01	0.472 ± 0.018	1.59	1.9	708
HW Car	5.01	50	0.96	0.149 ± 0.003	0.48	0.6	789	0.61	0.108 ± 0.001	0.22	0.1	544
DD Cas	2.28	82	1.05	0.099 ± 0.001	0.35	0.4	81	1.11	0.100 ± 0.001	0.37	0.4	80
KN Cen	2.06	47	0.97	0.858 ± 0.568	2.81	3.1	19831	0.99	0.763 ± 0.724	2.54	2.8	19791
SZ Cyg	2.28	85	0.96	0.076 ± 0.001	0.24	0.3	234	1.75	0.147 ± 0.001	0.87	0.1	204
CD Cyg	2.28	85	1.06	0.122 ± 0.002	0.43	0.6	552	1.96	0.142 ± 0.002	0.93	0.2	529
VX Per	3.99	100	3.99	0.073 ± 0.002	0.98	0.2	290	0.62	0.105 ± 0.001	0.22	0.0	262
X Pup	1.38	58	1.13	1.168 ± 0.098	4.43	4.2	2093	1.14	1.128 ± 0.093	4.31	3.9	2077
AQ Pup	2.96	76	1.01	1.567 ± 0.047	5.32	6.7	12079	1.66	1.038 ± 0.004	5.77	0.4	2299
WZ Sgr	4.02	57	0.93	0.222 ± 0.075	0.69	0.7	5126	1.90	0.406 ± 0.023	2.59	0.6	4913
RY Sco	2.06	33	0.99	1.901 ± 1.161	6.34	7.0	13759	0.67	0.410 ± 0.007	0.93	0.2	13424
Z Sct	3.10	111	0.94	0.444 ± 0.019	1.40	1.5	3156	0.93	0.450 ± 0.018	1.41	1.4	3148
S Vul	1.87	34	0.97	0.075 ± 0.005	0.24	0.3	32	0.61	0.192 ± 0.001	0.39	0.1	27

Note. For each star, we list the temporal baseline of our new observations, the degrees of freedom after the fit, as well as best-fit results for (a) the orbital solution leading to the maximum unsigned projected parallax error (typically near one year orbital period), and (b) the orbital solution with overall minimal χ² (see Section 3.3 for a related discussion). Uncertainties on K are based on the fit covariance matrix. $a\sin i$ is calculated using K and the corresponding fixed ${P}_{\mathrm{orb}}$ (see Equation (5)). $\hat{\varpi }\sin i$ is estimated as described in Section 4.1 and can have positive or negative sign. $\hat{\varpi }\sin i$ and $a\sin i$ are given in au, which is equivalent to percent of the parallax. For all solutions shown here, eccentricity e = 0. The orbital elements of YZ Carinae are given separately in Table 5.

A machine-readable version of the table is available.

Download table as:  DataTypeset image

As mentioned in Section 3.2 above, we find no indication for binarity for these 18 Cepheids for a range of orbital periods on the order of the observational baseline (two to five years, depending on the star). More importantly, the minimum-χ² solutions provide an estimated mean upper limit on parallax error of 0.8% for all 18 stars (<1% for 13, <4% for all 18), despite some imperfections in the sampling of some Cepheids. We caution, however, that the systematic uncertainty of ${P}_{\mathrm{puls}}$ due to incomplete phase coverage may affect the result for S Vul, see Section 3.2, although additional observations are required to determine whether this is the case.

The stars with the weakest constraints on $\max (| \hat{\varpi }\sin i| )$ are RY Sco (7%), AQ Pup (6.7%), X Pup (4.2%), and KN Cen (3.1%), all of which exhibit signs of cycle-to-cycle fluctuations of pulsation period and/or amplitude. The 14 remaining Cepheids have $\max (\hat{\varpi }\sin i)\lt 2 \%$, even for these solutions with maximal impact. As expected, most of the best-fit orbital solutions that would lead to the greatest parallax error are near the one year alias between orbital and parallactic motion. Both exceptions for which this is not the case, SY Aur and VX Per, yield the largest parallax error at ${P}_{\mathrm{orb}}$ corresponding to the baseline of the measurements. These results therefore strongly suggest that an astrometric modeling assuming a single star configuration is appropriate for all 18 Cepheids; see Section 4.2 for the exception of YZ Car.

Qualitatively, the flat pulsation-only residuals shown in Figure 2 already indicated that no large parallax error due to orbital motion was to be expected for these stars. The above results for $a\sin i$ and parallax error mirror and quantify this point. While a main limitation of this RV-based work is its insensitivity to inclination, this quantification of possible undetected configurations serves to increase confidence in the accuracy of the parallax measurements themselves and will be useful for future vetting of candidate high-accuracy calibrators of the Galactic Leavitt law.

YZ Carinae (${P}_{\mathrm{puls}}$ = 18.1676) is the only spectroscopic binary Cepheid among the program stars whose ${P}_{\mathrm{orb}}$ is shorter than our observational baseline. Its spectroscopic binary nature was discovered and originally reported by Coulson (1983) together with a preliminary orbital estimate of ${P}_{\mathrm{orb}}$ ∼ 850 days and low eccentricity. Petterson et al. (2004) obtained additional, higher-precision RV data and determined a significantly shorter orbital period of 657 days with similar eccentricity.

We here update and improve YZ Car's orbital solution by fitting a combination of a Fourier series and a Keplerian orbit to the new, highly precise, Coralie data presented here together with RVs published by Bersier (2002) and the post-1996 measurements by Petterson et al. (2004, Table A4). Figure 5 illustrates the quality of this solution. To verify our result, we also determined the orbit including older measurements by Pont et al. (1994) and Coulson (1983) in the fit, finding excellent agreement. However, these older data do not improve the quality of the solution due to larger measurement uncertainties and/or the possibility of pulsation period changes and we therefore prefer the solution based exclusively on RVs with uncertainties better than 300 m s⁻¹. The value of ${P}_{\mathrm{puls}}$ adopted for this modeling was determined using ASAS V-band photometry, since orbital motion significantly affects the measured RV on timescales of a month. Details of YZ Car's orbital solution and previous determinations are provided in Table 5.

Zoom In Zoom Out Reset image size

Table 5.  Orbital Elements for YZ Carinae

Parameter	C83a	P04b	This Work
${P}_{\mathrm{orb}}$ (days)	850 ± 11	657.3 ± 0.3	830.22 ± 0.34
e	0.13 ± 0.07	0.14 ± 0.03	0.041 ± 0.010
K (km s−1)	9.4 ± 0.5	10.0 ± 0.4	10.26 ± 0.82
vγ (km s−1)	1.0 ± 0.4	0.0 ± 0.2	0.844 ± 0.063
T0–2.4M	43575 ± 11	42250 ± 9	53422 ± 29
ω (deg)	239 ± 6	116 ± 5	195 ± 12
$a\sin i$ [106 km]	⋯	89	117.1 ± 9.4
fmass (${M}_{\odot }$)	0.071	0.066	0.093 ± 0.041
imin (deg)			22†
rms (km s−1)	⋯	⋯	0.39
$\hat{\varpi }\sin i\ ( \% \ \varpi )$	⋯	⋯	26
$\hat{\varpi }\sin i/\sin {i}_{\min }\ ( \% \ \varpi )$			76

Notes. Some uncertainties derived here are larger than in the literature due to the combined pulsation plus orbital fit. The minimum inclination (marked by ^†) is determined from the mass function assuming m₁ = 7 M_⊙. For m₁ = 9 M_⊙, min i = 20 deg.

^aCoulson (1983). ^bPetterson et al. (2004).

Download table as:  ASCIITypeset image

The value of ${P}_{\mathrm{orb}}$ ∼ 830 days determined here is nearly in agreement with the rough estimate provided by Coulson (1983, 850 days), and strongly disagrees with the solution presented by Petterson et al. (2004, 657.3 days), which is striking due to the small uncertainties quoted in the latter publication. Based on a visual inspection of the various available data sets, we conclude that the pre-1996 data employed by Petterson et al. (2004, Table A3) (pre-1996 Mount John University Observatory (MJUO) RVs) are not comparable with the other available RV data. This mismatch of data could be explained by the fact that the pre-1996 MJUO RVs were measured by different collaborators who may have employed nonstandard definitions of RV, such as bisector velocities, or measured velocities of Hα rather than metallic lines (Wallerstein et al. 1992). The fact that the updated result presented here is consistent with all other data spanning nearly 40 years, including Petterson et al.'s post-1996 data, strongly supports this conjecture. In addition, we point out the order of magnitude smaller residual scatter in Figure 5 compared to the scatter of residuals shown in Petterson et al. (2004, Figure 5).

Using this updated orbital solution together with the actual dates of the HST spatial scan observations, we determine $\hat{\varpi }\sin i=0.26\,\mathrm{au}$. Assuming an average inclination of 60 degrees, this would lead to a parallax error of up to 30% (e.g., ±100 μarcsec at 3 kpc distance), which should be clearly noticeable in the astrometric measurements. We therefore hope to obtain two additional epochs of spatial scan observations of YZ Carinae in order to improve the parallax measurement by accounting for orbital motion in the astrometric model.

Whereas the above sections focus on relatively short timescale orbital motion (${P}_{\mathrm{orb}}$ ≲ 5 years), we also investigate longer-term variations of v_γ, which are usually interpreted as evidence for spectroscopic binarity, by comparing our data with older measurements from the literature.

To this end, we divide the available data—literature RVs and our new RVs—into tranches that provide adequate phase sampling, balancing better phase coverage against temporal baseline per tranche. We then determine best-fit ${P}_{\mathrm{puls}}$, epoch of minimum RV, and v_γ for the data corresponding to each tranche using an RV template fitting approach. To achieve an accurate result it is crucial for the data belonging to a given tranche to sample both the rising and falling branch of the RV curve. The shaded regions in Figure 2 indicate how data tranches were selected for our new data; literature data were done analogously by inspection.

The RV templates used to fit each data tranche are created as Fourier series with harmonic coefficients resulting from the pulsational RV curve modeling described in Section 3.2. Each template fit determines two quantities for a fixed pulsation period ${P}_{\mathrm{puls}}$(t): v_γ(t), and a phase offset δϕ(t) relative to the mean observation date required to determine the time of minimum RV (E(t)) corresponding to this tranche and ${P}_{\mathrm{puls}}$(t). To account for period changes, we determine the globally best-fitting (minimum χ²) solution for a grid of input ${P}_{\mathrm{puls}}$ that lies within 0.1 days of the value listed in Table 3. We then repeat this procedure to within 0.01 days around the previous best-fit period to achieve a finer result. The final result of each fitted tranche is visually inspected to ensure a satisfactory result.

The main limitations of using time-variable v_γ as an indicator of spectroscopic binarity are (1) RV zero-point offsets among spectrographs and authors (up to several hundred m s⁻¹, see Evans et al. 2015b); (2) nonlinear period fluctuations preventing adequate phase-folding of a given tranche's data (which can be on the order of 1 km s⁻¹; see AQ Pup in Section 4.3.2); (3) apparent changes in in v_γ induced by cycle-to-cycle changes of RV variability (up to a few hundred m s⁻¹; see Anderson et al. 2016a). Determining the impact of (1) would require precision standard star RV time-series that are generally not available in the literature. Points (2) and (3) are particularly relevant for long-period (${P}_{\mathrm{puls}}$ ≳ 20 days) Cepheids as explained in Section 3.3. To avoid spurious detections, we therefore consider the overall behavior of v_γ(t) over all tranches and adopt a threshold of 1 km s⁻¹ as the minimum offset before concluding on spectroscopic binarity.

We thus investigated possible long-term variations of v_γ for 15 of our 18 Cepheids, excluding YZ Car (Section 4.2), HW Car and S Vul (both: lack of literature data). In the following, we report the discovery of three new candidate spectroscopic binaries (Section 4.3.1), followed by a critical investigation of Cepheids previously reported to be spectroscopic binaries (Section 4.3.2). Cepheids reported here or in the literature to be binaries are shown in Figure 6, whereas Figure 7 presents Cepheids for which no significant variations in v_γ were found and that have not previously been reported to be binaries.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

4.3.1. New Candidate Spectroscopic Binaries

Based on our RV template fitting approach, we report the discovery of three new spectroscopic binary candidates: XZ Car, AQ Car, and CD Cyg; Table 6 lists these results.

Table 6.  Time-variable v_γ of New Spectroscopic Binaries Based on RV Template Fitting

Cepheid	${P}_{\mathrm{puls}}$	Epoch	Δt	NRV	vγ	stdmer	σ(vγ)	References
	(days)	BJD–2.4M	(years)		(km s−1)	(km s−1)	(km s−1)
XZ Car	16.6525	44464.9312	2.21	38	1.07	0.22	0.16	CCG
XZ Car	16.6565	45164.0712	0.98	13	2.66	0.41	0.18	CCG; PBM
XZ Car	16.6573	55955.0183	0.05	15	5.49	0.36	0.003	h
XZ Car	16.6507	56021.6499	0.04	9	5.58	0.23	0.003	h
XZ Car	16.6495	56054.9312	0.04	11	5.55	0.26	0.003	h
XZ Car	16.6549	56088.2126	0.04	5	5.51	0.27	0.001	h
XZ Car	16.6523	56404.6109	0.33	12	5.72	0.18	0.001	h
XZ Car	16.6473	56671.0475	0.04	12	5.78	0.36	0.002	h
XZ Car	16.6473	56704.3392	0.04	14	5.53	0.39	0.004	h
XZ Car	16.6569	56787.6148	0.04	14	5.63	0.26	0.001	h
XZ Car	16.6549	57054.0691	0.04	19	5.69	0.21	<0.001	h
XZ Car	16.6541	57486.9476	0.03	7	5.54	0.46	0.011	h
AQ Car	9.7745	43975.7970	0.02	8	2.03	0.32	0.27	CCG
AQ Car	9.7717	44337.1307	0.58	18	2.61	0.49	0.70	CCG
AQ Car	9.7689	44698.7272	0.51	25	2.02	0.32	0.32	CCG
AQ Car	9.7679	45148.1286	1.32	19	1.43	0.17	0.018	CCG; PBM
AQ Car	9.7677	50599.3879	0.83	15	0.17	0.21	0.015	B02
AQ Car	9.7679	56793.2563	0.03	13	1.64	0.20	0.001	h
AQ Car	9.7705	57057.0354	0.04	21	1.54	0.14	<0.001	h
AQ Car	9.7673	57154.7358	0.03	13	1.69	0.15	0.001	h
AQ Car	9.7698	57457.5769	0.11	12	1.45	0.26	0.002	h
CD Cyg	17.075	25787.5216	5.19	12	−10.39	1.07	4.51	J37
CD Cyg	17.081	31301.5552	0.07	22	−11.90	0.63	2.54	S45
CD Cyg	17.0759	44686.6171	2.66	30	−11.29	0.22	0.020	I99; B88
CD Cyg	17.0758	46359.9408	3.54	10	−12.05	0.38	0.044	I99
CD Cyg	17.079	48852.8005	0.27	19	−11.98	0.37	0.051	G92
CD Cyg	17.071	49228.4290	0.15	6	−13.35	0.66	0.093	G92
CD Cyg	17.071	49569.8434	0.22	9	−13.29	0.36	0.055	G92
CD Cyg	17.081	49945.3098	0.12	12	−14.17	0.58	0.52	G92
CD Cyg	17.0792	56570.5614	0.36	30	−8.56	0.09	0.001	h
CD Cyg	17.0748	56877.9332	0.34	39	−8.74	0.13	0.001	h
CD Cyg	17.0742	57151.1361	0.19	17	−8.80	0.19	0.003	h
CD Cyg	17.0752	57304.8482	0.25	25	−8.57	0.13	0.001	h

Note. ${P}_{\mathrm{puls}}$, E, and v_γ are based on RV template fitting. Δt indicates the timespan of the measurement, N_RV the number of measurements fitted, "stdmer" the standard mean error based on residual scatter, σ(v_γ) the uncertainty from the fit covariance matrix.

References. CCG: Coulson et al. (1985), PBM: Pont et al. (1994), B02: Bersier (2002), J37: Joy (1937), I99: Imbert (1999), B88: Barnes et al. (1988), S45: Struve (1945), G92: Gorynya et al. (1992), h: this work.

A machine-readable version of the table is available.

Download table as:  DataTypeset image

Column "stdmer" quotes the standard mean error based on the residual scatter of the fit and should be compared to the fit uncertainty derived from the covariance matrix. It is interesting to note that "stdmer" tends to be smaller than σ(v_γ) for older, imprecise data, whereas the opposite is the case for new high-precision data. Specifically, the improvement of "stdmer" stagnates compared to the improvement in σ(v_γ) when using more precise (newer) data to determine v_γ. This is a consequence of the intrinsic astrophysical noise of Cepheid pulsations that manifest as fluctuations in period and RV curve shape (Anderson 2014; Anderson et al. 2016a). Inspection of the v_γ values derived for XZ Car shows that this astrophysical noise can lead to variations larger than a factor of several σ(v_γ). For XZ Car specifically, the (unweighted) mean v_γ inferred by template fitting of exclusively new measurements is v_γ = 5.603 ± 0.030 km s⁻¹, which is very close to the value of v_γ determined in a combined Fourier fit (5.671 km s⁻¹, see Table 3).

Adopting the above-stated threshold of 1 km s⁻¹, we find that XZ Car, AQ Car, and CD Cyg exhibit significantly time-dependent v_γ on timescales longer than a few years. This marks the first discovery of XZ Car's binarity, and our new data are decisive in demonstrating the likely binary nature of CD Cyg, which was previously considered not to be a spectroscopic binary (Evans et al. 2015b). While AQ Car's comparatively small v_γ variation (∼1.5 km s⁻¹) between that of Bersier (2002) and our RVs renders this evidence tentative, we note that our proposition that this is evidence of spectroscopic binarity is corroborated by the larger difference to older RVs (>2 km s⁻¹ Coulson et al. 1985; Pont et al. 1994) as well as the high quality of the Bersier (2002) data that provide a well-constrained fit.

We note that none of these three Cepheids are expected to have incurred significant parallax error due to orbital motion, Table 4 listing ≤1.6% for each of their respective $\max (\hat{\varpi }\sin i)$ solutions. However, proper motions estimated using long temporal baselines such as the Tycho–Gaia astrometric solution (Michalik et al. 2015; Lindegren 2016) may be affected by long-timescale orbital motion.

4.3.2. Revisiting Previously Reported Spectroscopic Binary Candidates

In addition to the new binaries presented above, Figure 6 also shows v_γ(t) of seven Cepheids that have previously been considered to be spectroscopic binaries.

We have recently discussed the binarity of SS Canis Majoris in light of our HST astrometric and recent high-precision RV measurements (Casertano et al. 2016). Using our RV template fitting technique, we here additionally investigate long-term variations of v_γ not discussed in our preceding paper. We find that the oldest data by Joy (1937) are not sufficiently accurate to determine v_γ with precision, although the central values of our fit results do reproduce the difference of ∼15 km s⁻¹ compared to RV data of Coulson & Caldwell (1985) as previously reported (Szabados 1996). Using Coulson & Caldwell's data, we determine a significant offset in v_γ of ∼3.8 km s⁻¹ compared to our new data, see Table 7, which would support a long-timescale spectroscopic binary nature of SS CMa.

Table 7.  Time-dependent v_γ Based on RV Template Fitting for Reported Binary Cepheids

Cepheid	${P}_{\mathrm{puls}}$	Epoch	Δt	NRV	vγ	stdmer	σ(vγ)	References
	(days)	BJD–2.4M	(years)		(km s−1)	(km s−1)	(km s−1)
SS CMa	12.3568	27554.7032	3.15	5	58.36	2.51	20.2	J37
SS CMa	12.3536	44373.6734	2.04	47	73.39	0.45	0.61	CC85
SS CMa	12.3528	56716.3689	1.58	27	77.16	0.10	0.002	h
SS CMa	12.3522	57062.2692	0.12	26	77.08	0.17	0.001	h
SS CMa	12.3478	57346.3907	0.06	19	77.07	0.21	0.002	h
SS CMa	12.3478	57482.2462	0.02	8	77.19	0.49	0.014	h
VY Car	18.8865	33999.9783	1.04	15	−2.31	0.59	1.28	S55
VY Car	18.8847	40567.9096	1.08	6	−2.77	0.64	0.95	L80
VY Car	18.8853	44483.7134	3.08	60	2.51	0.30	0.27	CC85
VY Car	18.8841	50551.7748	0.82	16	1.47	0.45	0.064	B02
VY Car	18.8831	55803.4453	1.21	41	1.64	0.41	0.002	h
VY Car	18.8825	57162.9905	1.19	42	1.61	0.26	0.001	h
KN Cen	34.0232	44545.6646	3.08	34	−39.05	0.41	0.50	CC85
KN Cen	34.0192	56965.1895	1.05	47	−42.15	0.48	0.008	h
KN Cen	34.0206	57509.4973	0.18	23	−42.32	0.98	0.036	h
SZ Cyg	15.1127	25826.3093	7.45	8	−16.97	2.26	19.3	J37
SZ Cyg	15.1063	31310.5534	0.07	17	−12.05	0.74	4.06	S45
SZ Cyg	15.1079	44622.2780	2.48	28	−10.91	0.26	0.038	I99; B88
SZ Cyg	15.1112	45589.2646	3.02	18	−11.38	0.17	0.012	I99
SZ Cyg	15.1113	47145.6313	3.29	14	−11.41	0.22	0.016	I99
SZ Cyg	15.1157	56574.1255	0.36	30	−10.10	0.11	0.002	h
SZ Cyg	15.1113	57057.7047	1.36	77	−9.97	0.05	<0.001	h
X Pup	25.9582	25896.3080	3.98	10	64.78	1.54	14.6	J37
X Pup	25.9628	44874.5896	2.33	32	67.19	0.49	0.82	B88; C01
X Pup	25.9608	50924.8055	2.10	33	71.40	0.33	0.034	B02; P05
X Pup	25.9605	54924.3983	1.23	42	72.15	0.37	0.017	S11
X Pup	25.9602	57235.7830	1.38	84	70.97	0.25	0.002	h
AQ Pup	30.1768	33078.6817	0.45	13	60.75	0.37	0.52	S55
AQ Pup	30.1768	44647.4321	0.49	30	56.95	0.43	0.63	CC85; B88
AQ Pup	30.1834	45159.8310	1.21	14	59.88	0.49	0.94	CC85; B88
AQ Pup	30.1806	54829.5366	1.74	38	61.25	0.70	0.068	S11
AQ Pup	30.1820	57093.2976	2.96	98	60.80	0.34	0.006	h
WZ Sgr	21.8521	25849.9613	5.09	9	−9.91	0.92	3.23	J37
WZ Sgr	21.8549	44553.0330	2.83	24	−16.97	0.37	0.46	CC85; B88
WZ Sgr	21.8461	49578.7275	0.17	26	−17.74	0.30	0.12	G92
WZ Sgr	21.8561	49950.1832	0.27	19	−18.14	0.39	0.05	G92
WZ Sgr	21.8463	50277.9982	0.22	38	−18.03	0.25	0.02	B02; G92
WZ Sgr	21.8535	50649.5795	0.49	32	−17.86	0.25	0.02	B02; G92
WZ Sgr	21.8521	56833.7761	0.23	33	−17.11	0.22	0.001	h
WZ Sgr	21.8507	57161.5552	0.18	29	−17.12	0.36	0.006	h
WZ Sgr	21.8516	57489.2777	0.18	13	−17.21	0.30	0.005	h

Note. See also Table 6. J37 data do not constrain the fit well for SS CMa, SZ Cyg, and X Pup.

References. S45: Struve (1945), CC85: Coulson & Caldwell (1985), L80: Lloyd Evans (1980), B88: Barnes et al. (1988), I99: Imbert (1999), S11: Storm et al. (2011), G92: Gorynya et al. (1992), B02: Bersier (2002), J37: Joy (1937), S55: Stibbs (1955), P05: Petterson et al. (2005), C01: Caldwell et al. (2001).

A machine-readable version of the table is available.

Download table as:  DataTypeset image

KN Centauri is a special case among binary Cepheids in that its hot main sequence companion has been detected and characterized using optical photometry (Madore 1977; Madore & Fernie 1980) as well as UV (Böhm-Vitense & Proffitt 1985; Evans 1994) and optical (Lloyd Evans 1980) spectra. However, the orbital signature on the RV curve had thus far escaped detection. We here report a detection of this signature based on the ∼3 km s⁻¹ offset in v_γ, see Table 7. The potential ability to detect the orbital motion of both components separately renders KN Cen an important target for a future model-independent mass measurement.

VY Carinae is a cluster Cepheid (Turner 1977; Anderson et al. 2013) whose previously reported spectroscopic binarity was based on evidence for low-amplitude variations of v_γ Szabados (1996, 2K < 5 km s⁻¹). The variation of v_γ for VY Car shown in Figure 6 is peculiar: the two oldest epochs indicate a rather constant velocity, and also the four newer epochs. At present, it is difficult to judge whether this pattern is caused by orbital motion or is due to systematics such as data quality or sampling.

SZ Cygni had previously been reported to exhibit time-variable v_γ by Struve (1945) and Szabados (1991). However, we find that the data upon which this evidence was based do not constrain v_γ at the crucial epochs. Notably, the data by Joy (1937) and Struve (1945) are imprecise and do not sample pulsation phase very well. We do, however, find that data by Imbert (1999) indicate a small difference of about 1.3 km s⁻¹ relative to v_γ determined using our new measurements.

X Puppis was first reported to exhibit time-variable v_γ by Szabados et al. (2012). Similarly to SZ Cyg, we find that the oldest RV data are insufficiently accurate and sampled to constrain the fit well. However, data published by Barnes et al. (1988) and Caldwell et al. (2001) yield a significantly lower v_γ than newer measurements, which furthermore exhibit a suspicious trend in v_γ(t). Additional high-precision observations of X Pup taken over the next few years will clarify whether this variation is consistent with orbital motion.

AQ Puppis exhibits very significant nonlinear changes of ${P}_{\mathrm{puls}}$ (Vinko 1991) in addition to exceptionally fast (300 s yr⁻¹) secular changes of ${P}_{\mathrm{puls}}$ (Berdnikov & Ignatova 2000) that approach values predicted for Cepheids on a first crossing of the instability strip (Turner et al. 2012; Anderson et al. 2016b). Fernie et al. (1966) reported a chance alignment with an OB association (see also Turner et al. 2012), which the HST spatial scans will further illuminate. Madore & Fernie (1980) presented evidence of a photometric companion based on both a phase-shift and color-loop methodology. AQ Pup's spectroscopic binary nature was reported by (Vinko 1991) based on a systematic offset in mean velocity between the data from Joy (1937) and those by Stibbs (1955), Barnes et al. (1988), and Coulson & Caldwell (1985). Ignoring the imprecise and extremely sparse (five measurements over four years) data by Joy (1937), we find that nearly all available RV data are consistent with a constant v_γ. However, RVs measured near epoch JD 44650 appear to be offset by ∼3 km s⁻¹ (see Table 7 and Figure 6) from measurements taken just one year later by the same authors (Coulson & Caldwell 1985; Barnes et al. 1988). Following visual inspection of the RV data and given that our new RV data do not indicate fast and significant variations in v_γ, we do not consider this apparent offset in v_γ to be a solid indication of spectroscopic binarity. Rather, it appears more likely that nonlinear changes of ${P}_{\mathrm{puls}}$ (see e.g., Turner et al. 2012, Figure 10) lead to problems in phase-folding the data, which is required to determine v_γ.

WZ Sagittarii is a member of the open cluster Turner 2 (Turner et al. 1993; Anderson et al. 2013), whose other members may aid in the determination of its accurate parallax. A spectroscopic binary nature of WZ Sgr has been argued both for (Szabados 1989) and against (Evans et al. 2015b). As Figure 6 shows, nearly all RV data are nicely consistent with a stable v_γ, the oldest data by Joy (1937, eight measurements) being the exception. We thus conclude that WZ Sgr is unlikely to be a spectroscopic binary.

The Cepheid binary fraction has been a topic of intense research for several decades (e.g., Lloyd Evans 1968; Madore 1977; Böhm-Vitense and Proffitt 1985; Szabados 1991; Evans & Udalski 1994; Evans et al. 2015b, 2016a, 2016b). Yet, inspection of the available literature data of previously reported candidate binaries (see Section 4.3) and the discovery of previously unknown binary systems among our relatively bright Cepheids suggests that far from everything is known for even relatively well-studied cases.

As mentioned in Section 2.1, the present sample of Cepheids is not random with respect to binarity and may therefore not be representative of the binary fraction of all Cepheids. Nevertheless, we summarize our investigation of binarity for the program Cepheids. Convincing evidence for spectroscopic binarity for five of the 19 Cepheids (SS CMa, YZ Car, XZ Car, KN Cen, and CD Cyg) exists. Four additional Cepheids (VY Car, AQ Car, SZ Cyg, X Pup) exhibit tentative evidence of variations in v_γ consistent with binarity, although imprecise literature data and sometimes heavily fluctuating pulsational variability render these results inconclusive. The literature further indicates DD Cas to have an unresolved photometric companion (this applies also to KN Cen), bringing the total number of binaries in this sample to between 6/19 (32%) and 10/19 (53%), depending on the inclusion of questionable candidates. This is broadly consistent with other recent estimates, e.g., by Neilson et al. (2015, 35%) and Evans et al. (2015, 29 ± 8% for ${P}_{\mathrm{orb}}$ < 20 years). These numbers do not include previously reported cases of visual binaries with extreme separations (greater than a couple of arcseconds) or Cepheids belonging to open clusters.

Furthermore, we stress that not all binary Cepheids are fundamentally unsuitable as high-accuracy Leavitt law calibrators, provided their photometry is not biased (see Section 5.2) and that parallax can be measured accurately (see Section 4.1). Additional photometric and interferometric observations will be useful to investigate these points.

The literature contains frequent references to binarity as being a significant source of photometric bias for the estimation of absolute magnitudes. For distance scale applications, this leads to two main questions: (1) what is the (pulsation-phase) average contrast between a Cepheid and a typical companion star? (2) How large of an effect on the distance scale could result from systematic differences in binary statistics among selected samples of Galactic and extragalactic Cepheids?

We therefore estimate the photometric contrast between Cepheids and typical, spatially unresolved, hot main sequence companions via isochrones computed using the Geneva stellar evolution group's (Ekström et al. 2012; Georgy et al. 2013) online tool.¹¹ We estimate the luminosity of the hot companion assuming a fixed mass ratio of q = m₂/m₁ = 0.7 as a conservative typical value based on the extensive work by N. Evans and collaborators (e.g., Evans 1995; Evans et al. 2013). For KN Cen and DD Cas, we use information on detected companions (Madore 1977) to estimate worst-case scenarios. Specifically, we assume a brighter B2 companion for KN Cen, despite IUE spectra indicating a B6 dwarf (Evans 1994). Masses are referred to here as lower case m to distinguish them from magnitudes M. The mass of the Cepheid is determined via a (pulsation) period-mass relation based on Geneva models (Anderson et al. 2016b) and thus sets the mass scale for both stars. The age of the isochrone is adopted based on period–age relations by Anderson et al. (2016b) using period change information measured or compiled by Turner et al. (2006) to inform the crossing number, where possible. The isochrone is computed for the inferred Cepheid's age, solar metallicity, and typical zero-age main sequence (ZAMS) rotation rate (ω = 0.5). The contrast in different pass-bands and filter combinations is estimated by forcing the Cepheid to be observed during blue loop evolution and looking up the properties of the companion with mass close to m₂ as per the adopted mass ratio. All hypothetical companions thus investigated are hot main sequence stars.

Using this approach, we estimate approximate magnitude differences in bolometric magnitudes, V-band, I-band, H-band, as well as reddening-free Wesenheit indices (Madore 1982) ${W}_{{VI}}=I-1.55(V-I)$ (Soszyński et al. 2008) and ${W}_{H,{VI}}\ =H-0.4(V-I)$ (Riess et al. 2011), which are known to be particularly useful for measuring Cepheid distances using PLRs thanks to reduced intrinsic PLR dispersion and reduced sensitivity to reddening.¹²

Table 8 lists the results obtained, including the adopted crossing number as well as inferred age and primary (Cepheid) mass m₁. The remaining columns list the quantities of the companion.

Table 8.  Estimating Typical Photometric Bias due to Main Sequence Companions

					Assumes m2/m1 = 0.7
Cepheid	Xing	age	$\mathrm{log}\,\mathrm{age}$	m1	m2	ΔMV	ΔMI	ΔMH	ΔWVI	${\rm{\Delta }}{W}_{{\rm{H}},\mathrm{VI}}$
SY Aur	3	81	7.91	6.0	4.2	−3.89	−4.58	−5.24	−5.65	−5.52
SS CMa	0	68	7.83	6.3	4.4	−4.04	−4.76	−5.44	−5.87	−5.73
VY Car	2	52	7.71	7.3	5.1	−4.31	−5.04	−5.75	−6.17	−6.05
XY Car	0	67	7.83	6.3	4.4	−4.04	−4.76	−5.44	−5.87	−5.73
XZ Car	0	57	7.75	6.9	4.9	−4.14	−4.87	−5.57	−6.01	−5.87
YZ Car	3	55	7.74	6.9	4.9	−4.14	−4.87	−5.57	−6.01	−5.87
AQ Car	0	78	7.89	6.0	4.2	−3.89	−4.58	−5.24	−5.65	−5.52
HW Car	0	80	7.91	6.0	4.2	−3.89	−4.58	−5.24	−5.65	−5.52
SZ Cyg	3	62	7.80	6.6	4.6	−3.99	−4.71	−5.40	−5.83	−5.69
CD Cyg	3	58	7.76	6.9	4.9	−4.14	−4.87	−5.57	−6.01	−5.87
VX Per	2	69	7.84	6.3	4.4	−4.04	−4.76	−5.44	−5.87	−5.73
X Pup	3	44	7.64	7.8	5.4	−4.42	−5.15	−5.87	−6.29	−6.17
AQ Pup	3	40	7.60	8.2	5.7	−4.48	−5.24	−5.98	−6.41	−6.29
WZ Sgr	3	49	7.69	7.3	5.1	−4.31	−5.04	−5.75	−6.17	−6.05
RY Sco	3	52	7.71	7.3	5.1	−4.31	−5.04	−5.75	−6.17	−6.05
Z Sct	3	69	7.84	6.3	4.4	−4.04	−4.76	−5.44	−5.87	−5.73
S Vul	3	24	7.37	10.2	7.2	−4.55	−5.76	−7.10	−7.63	−7.60
DD Cas (B4V)	3	82	7.92	6.0	4.95	−3.21	−3.93	−4.62	−5.05	−4.92
KN Cen (B2V)	3	37	7.57	8.6	6.9	−3.75	−4.74	−5.77	−6.28	−6.18

Note. Column X marks the instability crossing based on positive (assumed to be third crossings) and negative (second crossings) observed period changes (Turner et al. 2006). Ages are estimated using period-age relations for the appropriate crossing assuming average initial rotation rates (Anderson et al. 2016b). Cepheid masses m₁ are estimated using isochrones of solar metallicity and average ZAMS rotation (Ekström et al. 2012; Georgy et al. 2013) for a fixed adopted mean color $V-I=0.5$ during the blue loop phase. m₂ is the mass of a hypothetical companion, where m₂/m₁ = 0.7 for most cases (see Section 5.2) for the purpose of estimating the contrast in different filters and filter combinations. DD Cas and KN Cen are marked together with companion spectral-type estimates (Madore 1977). The Cepheids in this program are between 25 and 80 Myr old. Cepheids are much brighter than their main sequence companions, as expected, and this contrast increases when using longer-wavelength data and Wesenheit indices, as well as with pulsation period.

A machine-readable version of the table is available.

Download table as:  DataTypeset image

Figure 8 illustrates the comparison. It clearly shows that the contrast between Cepheid and companion increases with increasing wavelength. DD Cas and KN Cen stand out as spikes against the general trend due to higher mass ratios. Wesenheit indices amplify the contrast since Cepheids tend to be much redder than their hot main sequence companions. Figure 8 further suggests a ${P}_{\mathrm{puls}}$-dependence of this contrast that can be understood via the larger luminosity difference between stars on main sequence and blue loop evolutionary phases for younger (higher ${P}_{\mathrm{puls}}$) stars compared to older (lower ${P}_{\mathrm{puls}}$) stars.

Zoom In Zoom Out Reset image size

While a more detailed investigation of this effect would require population synthesis and examination of Cepheid colors and binary properties, we here note that the typical contrast using Wesenheit formulations or H-band photometry is on the order of six magnitudes at the typical period ($\mathrm{log}{P}_{\mathrm{puls}}\sim 1.3$) of extragalactic Cepheid samples (S. L. Hoffmann et al., submitted). The associated increase in brightness of −0.004 mag is much lower than the width of the instability strip (Riess et al. 2016, 0.08 mag in the H-band) and would lead to a distance error of merely 0.2% for an individual star, much less for an entire population. Binaries with lower contrast—such as KN Cen and DD Cas—are expected to be the exception and even for these cases, no strong photometric bias is expected in the H-band. Furthermore, selection criteria applied in the search for extragalactic Cepheids such as amplitude ratios are expected to remove very strong outliers (e.g., Hoffmann & Macri 2015). While photometric bias is comparatively stronger for shorter-period Cepheids and in optical pass-bands, we conclude that the light contributed by typical companion stars is on average negligible for distance scale applications where large numbers of Cepheids are observed using near-IR photometry and/or Wesenheit magnitudes.

An important additional consideration for the accuracy of parallax measurements is the phase-dependence of the contrast between a Cepheid and its spatially unresolved companion that shifts the measured photocenter in phase with the Cepheid's variability. This effect has been referred to as VIM (Wielen 1996) in the context of Hipparcos (Perryman & ESA 1997; van Leeuwen 2007). Whereas most Hipparcos VIM solutions were later found to be color-induced (Pourbaix et al. 2003), the correction for VIM-type effects is an integral part of Gaia data processing (CU4) due to Gaia's much increased astrometric precision. The same effect can also impact our HST spatial scanning measurements, depending on (i) the orientation of the scan direction with respect to the orientation of the binary; (ii) the (unresolved) angular separation of the two components; (iii) the average contrast between Cepheid and companion in the passband used to measure parallax; (iv) the pulsation phases sampled by the scan observations. We will study this effect in detail and constrain its impact on our Cepheid parallaxes in future work.

Cepheids with variability periods longer than ∼20 days have been shown in the literature to exhibit nonlinear variations of ${P}_{\mathrm{puls}}$ and other cycle-to-cycle modulations (e.g., Berdnikov et al. 2000; Anderson 2014). Among the present sample, the most affected stars are KN Cen, X Pup, AQ Pup, and S Vul, with AQ Pup having been discussed as a candidate for a first crossing Cepheid (Anderson et al. 2016b).

Nonlinear changes of ${P}_{\mathrm{puls}}$ complicate the inference of mean magnitudes from photometric measurements taken at random times using a known pulsation ephemeris. In the worst case, nonlinear period changes may lead to a complete loss of knowledge of the pulsation phase, leading to observations observed at random phase. Near-IR photometry can partially mitigate the scatter of the PLR determined by random-phase observations, since pulsation amplitudes decrease with increasing wavelength. For a given galaxy, the distance error contribution by this term is <0.12 mag, slightly larger than the intrinsic dispersion of the Leavitt law in the H-band (Riess et al. 2016).

At present, it is not clear what fraction of Cepheids exhibits such effects and how these phenomena are related to the ability to use affected Cepheids as precise standard candles. A characterization of nonlinear period changes in extragalactic Cepheids has thus far only been possible in the Magellanic Clouds (Poleski 2008; Soszyński et al. 2015, Süveges & Anderson, submitted). Gaia parallaxes and the ability to study Galactic Cepheids in great detail will enable a better understanding of pulsation irregularities and inform the vetting process of high-accuracy Leavitt law calibrators accordingly.