Discovery of a White Dwarf Companion to HD 159062

There are more than 200 known white dwarfs within 25 pc of the Sun (Sion et al. 2014), most of which are kinematically consistent with the thin disk population. Recent Gaia DR2 catalog queries report 139 white dwarfs within 20 pc (Hollands et al. 2018) and 153 within 25 pc (Jiménez-Esteban et al. 2018) based on strict cuts in Gaia color–magnitude space. The DA spectral type, which is distinguished by Balmer lines in the spectra for ${T}_{\mathrm{eff}}$ >  5000 K, is the most common, with DA white dwarfs approximately twice as abundant as all other spectral types (Giammichele et al. 2012; Sion et al. 2014; Kilic et al. 2018).

Approximately one-quarter of the known nearby white dwarfs reside in binary systems, most with a less-evolved companion but several in double-degenerate binaries (Holberg et al. 2016). This fraction is notably significantly lower than the field binary fraction. However, Toonen et al. (2017) performed population synthesis modeling including stellar multiplicity and evolution, and determined that the known WD/MS binary population is actually in reasonably good agreement with models. They determined that as high a fraction as 10%–30% of single white dwarfs originate as binary systems in which the components merge during post-main-sequence evolution. Low detection sensitivity to faint white dwarfs near bright, nearby main-sequence stars also likely plays a role in the observed rate of WD/MS binaries (Toonen et al. 2017).

Understanding the population of nearby white dwarfs, especially those in binary systems, is important for constraining the star formation history of the local neighborhood, as well as low-mass stellar evolution and white dwarf cooling models. Additionally, studying white dwarfs in binary systems can provide insights into the possible progenitors for SNe 1a. Several nearby benchmark white dwarf–main-sequence binary systems have been discovered with the combination of imaging and radial velocity (RV) data (e.g., Crepp et al. 2013, 2018), and a new detection of such a binary system is reported here.

HD 159062 is a bright nearby G9 dwarf star with low metallicity ([Fe/H] = −0.31 dex). It is likely fairly old; using its observed Ca ii H&K emission diagnostic ${R}_{\mathrm{HK}}^{{\prime} }$, and the empirical relation between ${R}_{\mathrm{HK}}^{{\prime} }$, rotation, and age from Mamajek & Hillenbrand (2008) results in an age estimate of ≈7 Gyr. HD 159062 is in the solar neighborhood at 21.7 pc (Gaia Collaboration et al. 2018), and as expected based on its age estimate, its kinematics are inconsistent with any nearby young clusters (Gagné et al. 2018). HD 159062 is consistent with an old main sequence G–K dwarf based on its $\mathrm{log}g$ and ${T}_{\mathrm{eff}}$. Its properties are detailed in Table 1, most of which derive from the literature and are indicated with superscripts on each parameter value.

Table 1.  Properties of HD 159062

HD 159062 Properties
R.A. (J2000)	17 30 16.4238
Decl. (J2000)	+47 24 07.922
J mag	5.804 ± 0.018a
${K}_{s}$ mag	5.392 ± 0.021a
W1 mag	5.397 ± 0.164b
π (mas)	46.123 ± 0.024c
d (pc)	21.68 ± 0.01
μα (mas yr−1)	169.747 ± 0.061c
μδ (mas yr−1)	77.164 ± 0.058c
${T}_{\mathrm{eff}}$ (K)	5283 ± 100d
$\mathrm{log}g$ (dex)	4.4 ± 0.1d
[Fe/H] (dex)	−0.31 ± 0.06d
Mass (M⊙)	0.76 ± 0.03d
Radius (R⊙)	0.76 ± 0.04d
$\mathrm{log}{R}_{\mathrm{HK}}^{{\prime} }$	−4.97
${S}_{\mathrm{HK}}$	0.170 ± 0.002
$V\sin i$ ($\mathrm{km}\,{{\rm{s}}}^{-1}$)	2.06 ± 0.5e

Note.

^a2MASS; Cutri et al. (2003). ^bWISE; Wright et al. (2010). ^cGaia DR2; Gaia Collaboration et al. (2018). ^dSpecMatch; Petigura et al. (2017). ^eSME; Piskunov & Valenti (2017).

Download table as:  ASCIITypeset image

Based on its [Fe/H] and [Mg/Fe] abundances as determined from high-resolution spectroscopy, Fuhrmann et al. (2017a) demonstrated that HD 159062 is most likely a Population II star. They argued that its age is most likely greater than the estimates from chromospheric activity or standard isochronal analysis would indicate, τ ≥ 12 Gyr. The authors list HD 159062 as a candidate blue straggler, indicating that rejuvenation from winds from a more evolved stellar companion may be complicating the standard age diagnostics.

Fuhrmann et al. (2017b) noted that HD 159062 has an anomalously high barium abundance of [Ba/Fe] = +0.4 ± 0.01 dex for a star of its type and metallicity. They argue that this may be further evidence of accretion of material onto HD 159062 via winds from an AGB companion that later become a white dwarf. Other abundance surveys of the solar neighborhood reported more typical values for the Ba abundance of HD 159062; Mishenina et al. (2008) reported [Ba/Fe] = +0.15 ± 0.1 and Reddy et al. (2006) reported [Ba/Fe] = +0.17 ± 0.12.

There exist many previously studied examples of binary star systems in which the more massive component has evolved off the main sequence through the AGB phase and into a white dwarf. These systems are typically distinguished by enhanced abundances of s-process heavy elements in the spectrum of the less-evolved star, which are typically thought to be a result of contamination from winds or accretion from the evolved companion during its AGB phase. Barium stars and CH stars are two subsets of this category of binary systems.

Barium stars are typically observed to be G and K giants, in binary systems with a white dwarf. These binaries have typical orbital periods of 500–10⁴ days, indicating that accretion from stellar winds (and not Roche-lobe overflow) is the dominant mechanism for mass transfer (Boffin & Jorissen 1988; Izzard et al. 2010; Van der Swaelmen et al. 2017).

The low-metallicity counterparts of Ba stars are CH stars, which are typically Population II stars with enhanced heavy-metal abundances, also due to accretion from an evolved binary companion. These types of stars were named for the lines arising from the CH molecule observed in their spectra. They are also formed via accretion from an evolved companion, and the observed binary period distribution for this class of stars is similar to that of Ba stars, with periods nearly always below 10⁴ days (Jorissen et al. 2016).

Recently, Escorza et al. (2019) published RV observations and orbit fits for a compilation of 60 Ba and CH binary systems. Of these, 27 had well-measured orbits with periods ranging from hundreds to >8000 days. Other systems had incomplete orbital coverage, so orbital periods were not estimated, but it can be inferred that these systems may have significantly longer orbital periods.

The inconsistent measurements of the metallicity of HD 159062 raise questions about whether or not the system could be an example of a mild Ba or CH binary. Either way, we report here that HD 159062 does definitively reside in a binary with an evolved companion.

In this paper, we report the discovery of the white dwarf companion to HD 159062, using 14 years of precise RV data from Keck/HIRES, as well as multi-epoch, multi-band imaging observations from the ShaneAO system at Lick Observatory, PHARO at Palomar Observatory and NIRC2 at Keck Observatory. By combining the RV and imaging data, we can constrain the orbit and cooling age of HD 159062 A and its white dwarf companion, HD 159062 B. We describe the observations and data reduction in Section 2, and the astrometry and common proper motion of the companion in Section 3. We demonstrate that HD 159062 B is not consistent with a brown dwarf or main-sequence companion based on its dynamics and color in Section 4. We detail our joint Bayesian analysis of the combined data set in Section 5, and present the best orbital parameters. Finally, we discuss the implications for the system's evolution, especially in the context of typical Ba star systems, in Section 6.

2.1. Radial Velocity Observations

Using HIRES at the 10 m Keck Telescope (Vogt et al. 1994), we have obtained 45 RV observations of HD 159062 since 2003 as part of the California Planet Survey (CPS). Several of these measurements were taken as double exposures, and were later binned on a 1 day cadence.

For HD 159062, we collected spectra through the B5 or C2 decker, with a width of 086 and a spectral resolution of R ≈ 55,000. The exposures were timed to yield a per-pixel signal-to-noise ratio (S/N) of ≈160 at 550 nm. A template spectrum (without iodine) was taken through the B3 decker, with a width of 057 and spectral resolution of R ≈ 66,000, and with a per-pixel S/N of ≈220.

To reduce the spectra and extract RVs, we used the standard CPS pipeline (e.g., Howard et al. 2010; Howard & Fulton 2016). In brief, starlight passes through a cell containing molecular iodine, imprinting a dense forest of absorption lines on the stellar spectrum. These lines are used as a high-fidelity wavelength calibration and point-spread function (PSF) reference. The template spectrum is obtained without the iodine cell to use as a stellar spectral reference. This reference spectrum is deconvolved against the instrumental PSF, and the wavelength solution, instrumental PSF, and RV are forward-modeled for each observing epoch. Each spectrum is divided into ∼700 "chunks," and the RV is extracted individually for each. The RVs and uncertainties reflect the weighted average and scatter of the results for the ensemble of spectral chunks.

The RV time series of HD 159062 shows a strong RV acceleration of −13.29 ± 0.12 ${\rm{m}}\,{{\rm{s}}}^{-1}\,{\mathrm{yr}}^{-1}$. Careful inspection of the RV time series shows slight curvature about this linear trend, which provides greater leverage on the full orbit. Figure 1 displays the RV time series, as well as a linear fit to the data. The RV data set is provided in Table 2, along with the internal RV uncertainties.

Zoom In Zoom Out Reset image size

Table 2.  Keck/HIRES Radial Velocity Data

Epoch (BJD)	RV (${\rm{m}}\,{{\rm{s}}}^{-1}$)	σ (${\rm{m}}\,{{\rm{s}}}^{-1}$)	${S}_{\mathrm{HK}}$
2452832.91668	89.22	1.62	⋯
2453074.09155	83.23	1.68	⋯
2453196.81631	76.88	1.58	⋯
2453430.15779	65.68	0.93	0.166
2453547.94513	60.98	0.84	0.167
2453604.84622	57.70	1.03	0.168
2453807.13748	50.30	1.00	0.168
2453926.45063	43.33	0.57	0.168
2453927.88530	44.96	0.83	0.167
2454248.03230	31.71	1.44	0.170
2454249.95458	32.00	1.49	0.170
2454252.04380	32.81	1.02	0.172
2454255.93510	32.28	1.01	0.172
2454277.86337	31.99	1.33	0.171
2454278.90748	32.22	1.22	0.171
2454279.94547	31.42	1.29	0.164
2454285.90800	30.58	1.34	0.170
2454294.90713	30.60	1.30	0.173
2454634.43562	19.74	1.03	0.170
2454635.94255	19.00	1.43	0.170
2454636.93316	17.04	1.53	0.171
2454638.33247	18.76	0.99	0.171
2454640.38633	14.57	1.06	0.171
2454641.94950	19.55	1.47	0.172
2454644.08325	15.11	1.52	0.170
2454688.88580	17.38	1.47	0.175
2454689.89920	9.37	1.58	0.174
2454957.13845	3.77	1.67	0.178
2455401.78282	−16.12	1.40	0.172
2455722.86450	−22.42	1.56	0.176
2455782.89029	−26.95	1.55	0.170
2456164.72278	−40.66	1.54	0.172
2456451.97582	−50.77	1.03	0.169
2456475.78721	−49.97	1.43	0.173
2456709.08945	−54.30	1.43	0.172
2456883.74892	−65.35	1.49	0.171
2456910.83550	−66.37	1.62	0.171
2457061.17682	−66.74	1.43	0.169
2457211.95513	−76.17	1.40	0.170
2457671.70374	−91.12	1.72	0.167
2457831.16397	−94.97	1.59	0.167
2457853.14530	−91.51	1.92	0.170
2457970.94501	−103.25	1.86	0.169
2458263.03029	−108.99	1.22	0.168
2458349.72648	−110.69	1.58	0.167

Note. RV uncertainties reported here are statistical uncertainties only, calculated from the weighted variance of the RV extracted in each spectral chunk (see Section 2.1). Jitter is not included.

Download table as:  ASCIITypeset image

2.2. High-resolution Imaging

2.2.1. Lick/ShaneAO

We observed HD 159062 using the ShaneAO system on the 3 m Shane Telescope at Lick Observatory (Gavel et al. 2014; Srinath et al. 2014). We obtained images in the ${K}_{s}$ band on UT 2016 April 20 and UT 2016 August 17. These initial images were taken as part of a uniform survey of stars in the solar neighborhood, in which the achievable contrast was improved by allowing the bright nearby primary stars to saturate the detector. With this technique, contrasts of 6–10 mag were frequently obtained at separations of 1''–3'' in under one minute of exposure time, at the expense of a saturation region including the core of the stellar PSF and occasionally the first Airy ring. Additionally, unsaturated follow-up data were required to characterize any newly discovered companions.

The saturated ${K}_{s}$ images both had total integration times of 19 s, composed of 13 frames of 1.46 s each, captured in a 4-point dither pattern with a single centered setup frame included and a dither throw of 2''. The resolution at the ${K}_{s}$ band of the Shane Telescope is 018. The inner working angle of the observations, outside of the saturation region, was ≈037, and the field of view extended to 10'' from the primary star in all directions. The plate scale was 32.6 ± 0.13 mas/pixel.

An additional image of HD 159062 was obtained on UT 2018 May 30 through the J and narrowband CH₄ − 1.2 μm filters. The image had a total integration time of 118.3 s, composed of 81 1.46 s integrations captured in a 4-point dither pattern with dither throws of 2'', 3'', and 4''. The field of view and plate scale were the same as those in the saturated ${K}_{s}$ images. The typical ShaneAO correction at J is significantly worse than the correction at ${K}_{s}$.

All images were flat-fielded and then sky-subtracted using 5 frames of dithered blank sky data of the same integration time as the science exposures. For the ${K}_{s}$ images, we used a KLIP reference differential imaging algorithm to suppress the PSF of the primary star, using all ${K}_{s}$ images of other stars obtained on the same night as the target exposure as a reference library (Wang et al. 2015). However, the companion was sufficiently separated from the primary star's PSF at 27, such that the local image statistics were not significantly affected by the PSF subtraction.

In the initial images, we detected a faint companion at ≈27 from the primary star. In the first image from 2016 April, the S/N of the detection was low, S/N = 3.5. In the second image from 2016 August, S/N = 5.6. The persistence of this low S/N companion over several months indicated its physical nature, but precise relative photometry could not be calculated from the K-band data due to the primary star's saturation. Therefore, we pursued follow-up imaging observations.

In the unsaturated J + CH₄ image, the companion detection had S/N = 8.4. We performed aperture photometry to obtain the J contrast, using an aperture radius of 1 FWHM = 3.9 pixels = 013 to extract relative flux values for the primary star and faint companion. A sky annulus spanning 3–5 FWHM was used.

Because of the brightness of the primary star, we extracted the flux of the companion using a high-pass filtered version of the image, to suppress the background flux at the location of the companion due to the bright primary.

We calibrated the companion flux measurement by injection/recovery. We injected 100 companions with similar contrast to the true companion into our image at 27 from the primary and at various position angles, yielding synthetic companion images, which were then high-pass-filtered with the same window size as was used for the true companion photometry. The flux ratios for the injected companions were then extracted and compared to the injected flux ratios. We determined the median multiplicative correction required, and this correction was applied to the aperture photometry for the true companion. The scatter in measured flux for the injected companions was added in quadrature to the uncertainty on our actual photometric data point.

This methodology yielded a measured contrast of ΔJ = 10.09 ± 0.38 mag. Here, we treat the contrast measured in the narrowband CH₄ − 1.2 μm filter as a proxy for the J-band contrast. The ShaneAO images are shown in Figure 2. The J-band image was not processed with the KLIP algorithm due to a lack of other J-band images from the night of the target observations to serve as a PSF reference library.

Zoom In Zoom Out Reset image size

2.2.2. Hale/PHARO

We observed HD 159062 with the PHARO adaptive optics system on the 200'' Hale Telescope at Palomar Observatory (Hayward et al. 2001) on UT 2017 June 7. The observations were carried out in the ${K}_{s}$ band with the narrow H₂ 2 − 1 filter in the grism wheel. The total integration time was 21 s, composed of 15 frames of 1.416 s each, captured in a 5-point dither pattern with a dither throw of 3''–5''. The pixel scale of the PHARO observations was 25.0 mas/pixel and the diffraction limit at the wavelength of observation was 011.

The bright primary star was unsaturated in these integrations. The frames were flat-fielded, and a sky background image was constructed from the median of the dithered images and subtracted. The companion was recovered at S/N = 6.5 in these Palomar data, and the unsaturated core of the primary PSF allowed us to calculate the relative brightness for this companion.

Aperture photometry was calculated on the unsaturated Palomar image using the same methodology as in the ShaneAO J-band image, with an aperture radius of 1 FWHM = 3.65 pixels = 009. This resulted in a contrast of Δ ${K}_{s}$ = 10.06 ± 0.22 for the companion, which is located approximately 30 FWHM away from the primary.

The imaging detection from PHARO is displayed in Figure 2, with the faint companion indicated by the white circle.

2.2.3. Keck/NIRC2

We followed up HD 159062 in $L^{\prime}$ band (3.8 μm) with Keck NIRC2 (Wizinowich et al. 2000, 2006; van Dam et al. 2006) on UT 2018 April 27 with natural guide star adaptive optics and the vortex coronagraph (Serabyn et al. 2017). The angular resolution was 008 and the plate scale was 9.942 ± 0.05 mas per pixel (Service et al. 2016). The field rotator was set to vertical angle mode, such that the telescope pupil tracks the elevation axis, to enable angular differential imaging (Marois et al. 2006). We took 85 frames with a discrete integration time of 0.5 s and 30 coadds resulting in a total integration time of 21.3 minutes over a 90 minute observation. During that time, we also took five images of the off-axis PSF with a discrete integration time of 0.008 s and 100 coadds as well as five images of the sky background with integration times matching that of the science and off-axis PSF frames. The sky off-axis PSF and background frames were taken every 15–20 minutes during the observing sequence. The alignment of the star and the center of the vortex focal plane mask was maintained by the QACITS tip-tilt control algorithm (Huby et al. 2017).

Bad pixels identified in dark frames and sky flats were replaced by the median of neighboring values. The sky flat was the median of 10 images of a blank patch of sky with the coronagraph focal plane mask removed (0.75 s discrete integrations, 10 coadds each). We subtracted the sky background frames from each frame individually using a scale factor to account for background variability. The frames were centered based on the position of the optical vortex core in the median of the science frames and each of the individual frames was co-registered using the peak of the cross-correlation with the median frame.

We applied principal component analysis (PCA, Soummer et al. 2012) to estimate and subtract the stellar contribution (on-axis PSF) from the images using the Vortex Image Processing (VIP) software package (Gomez Gonzalez et al. 2017). The final data product is a cube of median science frames with the starlight subtracted using 1–50 principal components (PCs). The significance of the imaged companion was maximized by subtracting the model of the on-axis PSF using 8 PCs. The PCs were computed within an annulus about the star centered on the companion with a width of 4 times the FWHM of the off-axis PSF (FWHM = 8 pixels).

To compute the photometry and astrometry of the companion, we subtracted a copy of the off-axis PSF at the location of the white dwarf in each science frame, varied its location and brightness, and repeated the PCA reduction with eight PCs until the values were minimized in the final reduced images in a 2× FWHM radius about the companion's position using a downhill simplex algorithm. To estimate the error, we re-injected the best-fit model of the companion into the pre-processed science frames at the same separation and brightness, but varied the parallactic angle by 10° and retrieved the photometry and astrometry using the same method 36 times, tracing a full circle about the star. We used the standard deviation of the measured flux and position of the injected companions as the uncertainty on each parameter. The NIRC2 narrow camera distortion correction was then applied by determining the location of the companion in each frame, interpolating the x-axis and y-axis distortion maps to the companion position, and then recalculating the astrometry of the primary star and companion. The scatter in the distortion map at the companion locations through the ADI sequence were added to the astrometric uncertainties in quadrature.

The companion was detected at an S/N of 17.5 with a flux ratio and angular separation between the companion and star of ${\rm{\Delta }}L^{\prime} =9.65\pm 0.06$ and 2673 ± 0003, respectively. The parallactic angle from the north toward the east was 30127 ± 005. The final reduced image is displayed in Figure 3.

Zoom In Zoom Out Reset image size

In addition to the new NIRC2 image obtained for this analysis, archival data on HD 159062 from UT 2012 June 26 and UT 2014 October 12 were discovered in the Keck Observatory archive. These data were collected as part of the TRENDS high-contrast imaging survey, and both epochs included a high-S/N detection of the faint companion.

Zoom In Zoom Out Reset image size

Table 3.  Summary of Imaging Observations

Epoch (JD)	Instrument	Filter	Separation ('')	P.A. (deg)	Contrast (Δ-mag)
2456102.89193	NIRC2	K'	2.594 ± 0.014	298.6 ± 0.1	⋯
2456942.70077	NIRC2	${K}_{s}$	2.637 ± 0.014	299.5 ± 0.1	⋯
2457498.99002	ShaneAO	${K}_{s}$	2.66 ± 0.03	301.5 ± 0.56	⋯
2457617.69501	ShaneAO	${K}_{s}$	2.68 ± 0.02	301.7 ± 0.47	⋯
2457911.93061	PHARO	${K}_{s}$	2.671 ± 0.004	301.3 ± 0.2	10.06 ± 0.22
2458236.09489	NIRC2	$L^{\prime}$	2.673 ± 0.003	301.27 ± 0.04	9.67 ± 0.08
2458268.86233	ShaneAO	J	2.67 ± 0.01	301.0 ± 0.24	10.09 ± 0.38

Download table as:  ASCIITypeset image

HD 159062 was observed through the TRENDS survey (Crepp et al. 2012) on 2012 June 24 UT and 2014 October 12 UT. Observations were taken with the narrow mode camera setting with plate scale 9.952 ± 0.002 mas pix⁻¹ (Yelda et al. 2010) using the full 1024'' × 1024'' field of view. Images were processed using standard flat-fielding, sky background subtraction, and focal plane distortion correction applied as prescribed in Yelda et al. (2010).

The UT 2012 June 24 observations consisted of 10 frames taken with integration times of 1 s per coadd and 10 coadds, collected in position angle mode. These observations were taken under clear conditions and at an airmass of 1.1 with the K' filter. The star was obscured by the 300 mas corona300 coronograph.

The UT 2014 October 12 observations consisted of 76 usable frames with integration times of 4 s per coadd and 6 coadds per frame. These data were taken under clear conditions and at an airmass of 1.5 with the ${K}_{s}$ filter and the instrument in vertical angle mode. In this epoch the 600 mas corona600 coronograph was used to obscure the bright primary star. Prior to measuring accurate separation and position angle, each image was de-rotated using the prescription of Yelda et al. (2010).

Astrometric information was extracted from the archival NIRC2 images, but due to duplication of filters and lack of unocculted images, photometry was not calculated from these observations. The reduced images are displayed in Figure 4.

The position of the companion in the saturated ShaneAO images, and the positions of both the primary and companion in the unsaturated ShaneAO, PHARO, and archival NIRC2 images, were extracted using the python implementation of DAOStarFinder in the photutils package (Bradley et al. 2017).

For the PHARO data, we adopted an uncertainty on the position of the primary star of 0.1 pixel, consistent with previous astrometric efforts on this instrument (D. Ciardi 2018, private communication). Since the companion was detected at low S/N, we use a heuristic estimate for uncertainty, σ ≈ FWHM/2.355/S/N = 0.24 pixel. The plate scale of the Palomar AO instrument has been measured to be 25.09 ± 0.04 mas/pixel by Metchev (2006), and the rotational position of the instrument has been shown to be stable to within 012. These uncertainties are added in quadrature to the astrometric data points.

For the unsaturated ShaneAO J-band image, we assume an astrometric precision of 0.2 pixel for both the primary and companion. The plate scale for the ShaneAO instrument is 32.6 ± 0.13 mas/pixel and the field rotation has recently been measured to be 187 ± 013 (G. Duchêne 2017, private communication). Once the separation and position angle were calculated from the ShaneAO data, 187 were added to the P.A. measurements, and the rotational and plate scale uncertainties were added in quadrature.

For the saturated ShaneAO images, the position of the saturated primary star was determined by rotational symmetry. In brief, we calculated the residuals between the reduced image and a rotated version of the image, for several different rotations about each prospective center pixel. We interpolated this residual map to find the optimal rotational centroid for the primary star. This method was adapted from Morzinski et al. (2015). We adopt 0.5 pixel uncertainty for the primary star due to its saturation, and 0.4 pixel for the secondary based on the heuristic diagnostic. The derivation of the astrometric position of the companion in the new NIRC2 image is described in Section 2.2.3. For the archival NIRC2 data, we adopted a 0.1 pixel uncertainty on the position of the companion, and due to the difficulty of centroiding the primary star behind the coronographic masks, we adopted a 0.5 pixel uncertainty on the position of the primary in both epochs. We added the uncertainty in the plate scale and true north angle (Yelda et al. 2010) in quadrature to our separation and P.A. measurements. Our complete astrometric data set, along with measured photometry from each image, is listed in Table 3.

Our astrometry of the imaged companion to HD 159062 spans six years, allowing us to check for common proper motion between the primary and imaged companion. HD 159062 has a proper motion of ΔR.A. = 169.747 ± 0.061 $\mathrm{mas}\,{\mathrm{yr}}^{-1}$ and Δdecl. = 77.164 ± 0.058 $\mathrm{mas}\,{\mathrm{yr}}^{-1}$ (Gaia Collaboration et al. 2018), large enough to allow us to detect the relative motion of a background source. Figure 5 shows the expected relative motion of a background source due to the proper motion of HD 159062, as well as the observed astrometry. We conclusively rule out the scenario that the imaged companion is a distant background stellar source.

Zoom In Zoom Out Reset image size

Since it is likely bound, HD 159062 B must be located at the same distance from Earth as its primary star. Combining this distance with its observed photometry rules out a main-sequence dwarf scenario for the companion; for reference, an M8V star at the bottom of the main sequence would have an expected contrast of only ${\rm{\Delta }}{K}_{s}\approx 6.5$ mag with the G9 primary, much less than the measured contrast of HD 159062 B. Instead, this companion must be an intrinsically faint object, either a brown dwarf or white dwarf.

A simple argument invoking acceleration rules out the brown dwarf scenario for HD 159062 B based on the mass constraint from the RV acceleration of HD 159062 A and the projected separation between the components.

Rearranging Newton's second law and the Newtonian law of gravity, for a companion B orbiting a primary star A we can equate the magnitude of the force vector and the magnitude of the acceleration, such that

$$\begin{eqnarray}&&{M}_{B}=\displaystyle \frac{{a}_{A}{r}^{2}}{G},\end{eqnarray} \tag{ 1 }$$

where M_B is the mass of the companion, a_A is the instantaneous acceleration of the primary, r is the physical separation between the two components, and G is the gravitational constant.

The lower limit for the instantaneous acceleration of the primary star comes from the slope of its radial velocities, dRV/dt = −13.29 ± 0.12 ${\rm{m}}\,{{\rm{s}}}^{-1}\,{\mathrm{yr}}^{-1}$. Because this describes only the acceleration in the radial direction, and does not account for acceleration along the plane of the sky, $\left|{dRV}/{dt}\right|$ is a lower limit for the full acceleration, or ${a}_{A}\geqslant \left|{dRV}/{dt}\right|$.

Likewise, the lower limit for physical separation is the projected separation of the two components, ${\rho }_{\mathrm{proj}}=57.78\pm 0.07\,\mathrm{au}$. This is also a projection of the full separation, so r ≥ ρ_proj.

By substitution into Equation (1), we find:

$$\begin{eqnarray}&&{M}_{B}\geqslant \left|\displaystyle \frac{{dRV}}{{dt}}\right|{\rho }_{\mathrm{proj}}^{2}{G}^{-1}\geqslant 0.24\,{M}_{\odot }.\end{eqnarray} \tag{ 2 }$$

This mass, derived from a simple acceleration argument, represents a dynamical lower limit for the companion, placing it securely above the maximum mass for a brown dwarf.

HD 159062 B's infrared colors confirm that it is not a brown dwarf. Since $L^{\prime}$ photometry was not readily available for HD 159062 A, we use the WISE W1 measurement of HD 159062 A as a proxy for $L^{\prime}$. This is reasonable because the W1 − $L^{\prime}$ color has been shown to be zero for stars earlier than M0 (De Rosa et al. 2016). From the apparent magnitudes of the primary star and the measured contrasts, we find that HD 159062 B has apparent magnitudes of J = 15.89 ± 0.38 mag, ${K}_{s}=15.45\pm 0.22$ mag, and $L^{\prime} =15.07\pm 0.18$ mag. Its colors are therefore $J-{K}_{s}=0.44\pm 0.44$ mag and ${K}_{s}-L^{\prime} =0.38\pm 0.28$ mag.

Low-mass stars and brown dwarfs spanning spectral types from M0 through T6 have ${K}_{s}$ − $L^{\prime}$ colors ranging from 0.5 to 2.0 mag, marginally consistent with the measured ${K}_{s}$ − $L^{\prime}$ color of HD 159062 B. However, these objects also have J − ${K}_{s}$ colors ranging from 0.8 to 1.8 mag (Leggett et al. 2001), significantly redder than the observed J − ${K}_{s}$ color of HD 159062 B.

We conclude that HD 159062 B is a white dwarf. Both the mass lower limit and the infrared color of the faint companion are consistent with this scenario. HD 159062 B is too massive to be a brown dwarf, and too faint to be a main-sequence dwarf. Additionally, the predictions from Fuhrmann et al. (2017b) regarding the enhanced Ba abundance of HD 159062 A, outlined in Section 1, strengthen the white dwarf interpretation.

We performed a joint MCMC analysis of the 45 RV data points and a subset of the astrometric data. Since the imaging observations were taken on three different instruments, uncertainties and variations in the field rotation and plate scale of the instruments make combining the astrometric data difficult. We therefore chose to limit our analysis to the three NIRC2 data points, which we expected to be the most consistent in plate scale and rotation, and which conveniently covered the largest time baseline of the orbit.

The likelihood function we used is

$$\begin{eqnarray}\begin{array}{rcl}\mathrm{log}{ \mathcal L } & = & -\displaystyle \sum _{i}\left[\displaystyle \frac{{\left({v}_{i}-{v}_{m}({t}_{i})\right)}^{2}}{2({\sigma }_{i}^{2}+{\sigma }_{\mathrm{jit}}^{2})}+\mathrm{log}\sqrt{2\pi ({\sigma }_{i}^{2}+{\sigma }_{\mathrm{jit}}^{2})}\right]\\ & & -\displaystyle \sum _{j}\left[\displaystyle \frac{{\left({x}_{j}-{x}_{m}({t}_{j})\right)}^{2}}{2{\sigma }_{j}^{2}}+\mathrm{log}\sqrt{2\pi {\sigma }_{j}^{2}}\right].\end{array}\end{eqnarray} \tag{ 3 }$$

Here, v_i (t_i) and x_j (t_j) correspond to the RV and astrometric measurements (epochs), where a zero-point offset γ has been added to the RV data as a free parameter in the model. Due to a HIRES CCD upgrade in 2004, we use one γ parameter for the first three RV data points, which were taken prior to the upgrade (γ_k), and another for the post-upgrade data (γ_j). v_m(t_i) and x_m(t_j) are the model velocities and positions at the RV and astrometric epochs, respectively. σ_i is the internal RV precision, and σ_jit is the additional uncertainty added in quadrature due to stellar jitter and instrumental uncertainties not included in the internal precision estimated from the RV chunks. σ_j is the astrometric uncertainty.

We implemented some standard priors on the sample parameters for this fit; they are detailed in Table 4. We sampled from a uniform mass distribution between our dynamical lower limit, 0.24 ${M}_{\odot }$ from Section 4, and the Chandrasekhar upper limit for white dwarf mass, 1.44 ${M}_{\odot }$. We chose to sample in $\mathrm{log}P$ as a variant of a Jeffreys prior for the orbital period because it is a scale parameter. We chose our lower bound on the orbital period to be several times lower than the minimum orbital period calculated by assuming the semimajor axis is equivalent to the measured projected separation. The projected physical separation we measured was approximately 60 au, which would correspond to a period of ≈400 yr assuming a total system mass of 1.3 ${M}_{\odot }$. The upper bound on period was conservatively chosen to be 10⁴ yr.

Table 4.  MCMC Priors on Orbital Parameters

Parameter	Prior
mA	Gaussian, 0.80 ± 0.05M⊙
mB	Uniform, [0.24, 1.44] M⊙
$\mathrm{log}\,P$	Uniform, [102, 104] yr
$\sqrt{e}\cos \omega$	Uniform, [0.0, 1.0]
$\sqrt{e}\sin \omega$	Uniform, [0.0, 1.0]
Ω	Uniform, [0, 2π]
M0(t0)	Uniform, [0, 2π]
$\cos i$	Uniform, [ − 1.0, 1.0]
π	Gaussian, 46.123 ± 0.024 mas
γj − γk	Gaussian, $0\pm 2\,{\rm{m}}\,{{\rm{s}}}^{-1}$

Download table as:  ASCIITypeset image

We sampled in $\sqrt{e}\cos \omega$ and $\sqrt{e}\sin \omega$ rather than in eccentricity and argument of periastron to improve convergence time while maintaining a uniform prior on the eccentricity. Our priors on longitude of the ascending node Ω and mean anomaly M₀ at the epoch of the first observation were uniform between 0 and 2π, as we have found that this improves convergence time over allowing these parameters to be unconstrained and later taking the modulus of 2π (R. De Rosa 2018, private communication). We sampled uniformly in $\cos i$ for uniform distribution of inclination over a sphere.

Previous studies have made use of RV observations of standard stars to place the HIRES pre-upgrade and post-upgrade data on a uniform scale. We therefore place a Gaussian prior on the difference between the pre-upgrade velocity zero-point γ_k and the post-upgrade zero-point γ_j such that ${\gamma }_{j}-{\gamma }_{k}=0\pm 2\,{\rm{m}}\,{{\rm{s}}}^{-1}$.

The priors on primary mass m_A and parallax π were Gaussian. The parallax constraint is from the Gaia DR2 astrometric results (Gaia Collaboration et al. 2018). For the primary mass, we based our prior on both a new spectroscopic estimate of the mass of HD 159062 A and on literature values. We use the stellar spectroscopic analysis tool SpecMatch (Petigura et al. 2017), which performs fits to theoretical stellar spectra to obtain stellar parameters ${T}_{\mathrm{eff}}$, $\mathrm{log}g$, and [Fe/H]. These parameters are then converted to stellar mass and radius using the package isoclassify (Huber et al. 2017). This results in a mass estimate for HD 159062 A of 0.76 ± 0.03 ${M}_{\odot }$.

We also perform a literature search for measurements of the mass of HD 159062 A. Ramírez et al. (2012, 2013) reported a mass of m_A = 0.8 ^+ 0.02_−0.01 ${M}_{\odot }$; Brewer et al. (2016) reported masses of m_A = 1.01 ± 0.14 ${M}_{\odot }$ or m_A = 0.8 ± 0.02 ${M}_{\odot }$, depending on which mass derivation method is used; Casagrande et al. (2011) reported ${m}_{A}={0.84}_{-0.02}^{+0.04}\,{M}_{\odot };$ and Fuhrmann et al. (2017b) reported m_A = 0.84 ${M}_{\odot }$. All literature mass measurements agree within 3σ with our spectroscopic measurement, but all indicate slightly more massive solutions. Therefore, we adopt a prior for HD 159062 A of m_A = 0.80 ± 0.05 ${M}_{\odot }$.

We used the parallel-tempering MCMC sampler from Emcee (Foreman-Mackey et al. 2013). We used 10 temperatures with the default temperature scale, 300 walkers, and iterated for 10⁵ steps. We examined the walker position plots, autocorrelation functions, and the evolution of the median and 1σ values of each parameter to assess when the burn-in phase was complete. We present the constraints on the mass, eccentricity, and other orbital parameters in column (2) of Table 5.

Table 5.  MCMC Results

Parameter	Median & 68% CI
Model	RV/Ast	RV/Ast/Phot
mA (${M}_{\odot }$)	0.80 ± 0.05	0.80 ± 0.05
mB (${M}_{\odot }$)	${0.65}_{-0.04}^{+0.14}$	${0.65}_{-0.04}^{+0.12}$
$\mathrm{log}\,P$ (yr)	${2.38}_{-0.15}^{+0.19}$	${2.40}_{-0.16}^{+0.18}$
$\sqrt{e}\cos \omega$	${0.61}_{-0.42}^{+0.16}$	${0.58}_{-0.46}^{+0.18}$
$\sqrt{e}\sin \omega$	$-{0.27}_{-0.10}^{+0.07}$	$-{0.27}_{-0.07}^{+0.07}$
Ω (°)	${138}_{-5}^{+11}$	${137}_{-4}^{+9}$
M0(ti = 0) (°)	${144}_{-28}^{+59}$	${147}_{-27}^{+65}$
$\cos i$	${0.60}_{-0.13}^{+0.23}$	${0.58}_{-0.12}^{+0.22}$
π (mas)	46.12 ± 0.02	46.12 ± 0.02
τ (Gyr)	⋯	${8.2}_{-0.5}^{+0.3}$
Derived Parameters
P(yr)	${238}_{-68}^{+128}$	${250}_{-76}^{+130}$
e	${0.44}_{-0.31}^{+0.30}$	${0.40}_{-0.28}^{+0.31}$
ω (°)	$-{26}_{-29}^{+7}$	$-{26}_{-38}^{+7}$
i (°)	${53}_{-19}^{+9}$	${54}_{-18}^{+8}$
${T}_{\mathrm{eff},{\rm{B}}}$ (K)	⋯	${4580}_{-160}^{+440}$
Instrumental Parameters
γk (${\rm{m}}\,{{\rm{s}}}^{-1}$)	${934}_{-605}^{+492}$	${1018}_{-572}^{+431}$
γj (${\rm{m}}\,{{\rm{s}}}^{-1}$)	${933}_{-605}^{+492}$	${1018}_{-573}^{+431}$
σjit,k (${\rm{m}}\,{{\rm{s}}}^{-1}$)	${2.3}_{-1.5}^{+2.9}$	${2.3}_{-1.5}^{+2.9}$
σjit,j (${\rm{m}}\,{{\rm{s}}}^{-1}$)	${1.35}_{-0.30}^{+0.33}$	${1.36}_{-0.30}^{+0.32}$

Download table as:  ASCIITypeset image

5.1. White Dwarf Cooling Models

5.2. Results

We find that all of the constraints from the combined analysis are consistent to within 1σ with those from the RV and astrometric analysis. The only additional constraints resulting from inclusion of the photometric prior are on the cooling age and temperature of the white dwarf. We find that the photometry is consistent with an old white dwarf ($\tau ={8.2}_{-0.5}^{+0.3}$ Gyr), which has cooled significantly.

To demonstrate the full behavior of the posteriors, including correlations between the various orbital elements, mass, and age, a corner plot of the derived white dwarf parameters is displayed in Figure 6. It is clear from the covariance plots that many of the posteriors on orbital parameters are highly correlated, which is not surprising given the limited coverage of the full orbit available in the RV and imaging data.

Zoom In Zoom Out Reset image size

The white dwarf companion HD 159062 B has a mass constraint of ${m}_{B}={0.65}_{-0.04}^{+0.12}$ ${M}_{\odot }$ and a long orbital period of $P={250}_{-76}^{+130}$ yr, much longer than the 14 yr baseline of our RV observations. The period remains poorly constrained, with a 1σ credible interval ranging from 174 to 380 yr. The eccentricity posterior appears to be multi-modal, with a peak at low eccentricity near e = 0.1 and another at high eccentricity near e = 0.7.

A plot showing a random sampling of 50 orbits drawn from the posteriors of the joint MCMC run is provided in Figure 7. From this plot, it is clear that the orbital period is not fully constrained with the data available.

Zoom In Zoom Out Reset image size

The two eccentricity modes can be seen in the possible RV curves, with some solutions showing sharp peaks indicative of a high eccentricity, while others are more sinusoidal in shape. However, both classes of orbit have the same astrometric signature, due to an anti-correlation between eccentricity and inclination. As might be expected, the low-eccentricity solutions have more edge-on orientations, while the higher-eccentricity solutions are closer to face-on, such that both types of orbit still fit the astrometric motion well. Additionally, eccentricity is highly anti-correlated with orbital period, with shorter orbital periods corresponding to higher eccentricities. A correlation between mass and period is also observed, thus the range of possible orbits follows a continuum from lower-mass, shorter-period, circular, edge-on cases to higher-mass, longer-period, eccentric, face-on solutions.

The maximum probability orbit model is plotted against the RV data in Figure 8, with residuals shown. A periodogram of these residuals is also displayed. Together, these plots demonstrate that the RV trend and slight curvature are fully accounted for in the single-companion Keplerian orbital fit. The residuals do not show correlated structure, and no significant shorter-period peaks are seen in the periodogram.

Zoom In Zoom Out Reset image size

Using the posteriors on mass and cooling age and the Montreal white dwarf cooling models, we calculate the posterior on the temperature of the white dwarf as well. This can be compared to the ${T}_{\mathrm{eff}}$ implied by the white dwarf spectrum if spectroscopy is obtained for this target. We obtain a derived constraint on HD 159062 B's temperature of ${T}_{\mathrm{eff}}={4580}_{-160}^{+440}$ K, slightly cooler than the main-sequence primary star, HD 159062 A.

Finally, for future observations, we predict the expected location of the white dwarf over the next two decades. We convert the posteriors on orbital parameters into posteriors on the angular offsets between HD 159062 A and B at several epochs spanning 2020 through 2040, and plot these 2D posterior distributions, along with the astrometric data, in Figure 9. The plot demonstrates that orbital motion on the scale of a few hundred mas will be detectable in the next decade. As expected, the uncertainty on the position of the white dwarf increases over time due to uncertainties in the orbital parameters. More astrometric data from Keck/NIRC2 will help to place new constraints on the orbital parameters after several years.

Zoom In Zoom Out Reset image size

5.3. Total System Age

The indirect suggestion of a white dwarf companion to HD 159062 A has been detailed in Fuhrmann et al. (2017a, 2017b). Specifically, the overabundance of barium measured by these studies indicates contamination from an AGB companion sometime during the evolution of HD 159062 A. We do note, however, that other abundance surveys of the local neighborhood by Reddy et al. (2006) and Mishenina et al. (2013) found lower and more consistent [Ba/Fe] abundance values of +0.17 and +0.15 dex, respectively.

In comparison with typical Ba or CH stars, HD 159062 A differs in several ways. First, its $\mathrm{log}g$ measured from high-resolution spectroscopy indicates that HD 159062 A is a main-sequence star. Typical Ba stars are observed to be G and K giants (Pols et al. 2003; Izzard et al. 2010; Mahanta et al. 2016), although dwarf Ba systems do exist. CH stars are more often main-sequence companions, but the metallicity of HD 159062 is too high to qualify as a typical CH star, which have [Fe/H] < −0.5 (Escorza et al. 2019). Ba binaries typically have orbital periods ranging from a few to a few tens of years. The observed separation of HD 159062 B makes it significantly wider than known Ba or CH star binaries, and its orbital period is P = 250 ^+ 130₋₇₆ yr, several times longer than typical Ba or CH binary orbital periods (Jorissen et al. 1998, 2016; Van der Swaelmen et al. 2017; Escorza et al. 2019).

The wide separation of HD 159062 B raises some doubts about its capability of contaminating HD 159062 A sufficiently via a stellar wind to explain the observed Ba abundance. Izzard et al. (2010) propose that the eccentricity-orbital period relation for Ba star binaries could be better explained by a process that "kicks" the white dwarf as it forms, possibly due to asymmetric mass loss or magnetic fields. In their model, a tail of high eccentricity, very-long-period Ba binaries is formed, possibly consistent with the orbit of HD 159062 B.

HD 159062 is a newly detected binary system with one old main-sequence component, and one initially more massive component that has since evolved into a white dwarf. We have shown that HD 159062 B cannot be a main-sequence or brown dwarf companion, and is also inconsistent with a background source. Based on its mass and photometry, it is consistent with a white dwarf companion. The system's old age means that the white dwarf has cooled significantly, with a temperature of ${{T}_{\mathrm{eff}}}_{,B}={4580}_{-160}^{+440}$ K and cooling age of ${t}_{\mathrm{cool}}={8.2}_{-0.5}^{+0.3}$ Gyr. Continued monitoring of the radial velocities of HD 159062 A via high-resolution spectroscopy, and especially additional measurements of the astrometric position of HD 159062 B through high-contrast imaging, will help to further constrain the long-period binary orbit and improve our understanding of the evolutionary history of this system.

We would like to thank the anonymous referee for many helpful and detailed suggestions that have helped to improve this manuscript. We thank Pierre Bergeron for providing updated white dwarf cooling models including $L^{\prime}$-band synthetic photometry. L.A.H. would like to thank James Graham for helpful discussion on this project, and Erik Petigura for carrying out HIRES observations and performing a spectroscopic analysis of the properties of HD 159062 A. M.R.K acknowledges support from the NSF Graduate Research Fellowship, grant No. DGE 1339067. Some of the data presented herein were obtained at the W. M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California and the National Aeronautics and Space Administration. The Observatory was made possible by the generous financial support of the W. M. Keck Foundation. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Maunakea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain.

Facilities: Keck:I (HIRES) - , Keck:II (NIRC2) - , Lick:Shane (ShARCS) - , Palomar:Hale (PHARO). -

Software: Emcee (Foreman-Mackey et al. 2013), VIP (Gomez Gonzalez et al. 2017), isoclassify (Huber et al. 2017), pyKLIP (Wang et al. 2015), photutils (Bradley et al. 2017).