Very Low-Mass Stellar and Substellar Companions to Solar-like Stars From MARVELS VI: A Giant Planet and a Brown Dwarf Candidate in a Close Binary System HD 87646
Bo Ma, Jian Ge, Alex Wolszczan, Matthew W. Muterspaugh, Brian Lee, Gregory W. Henry, Donald P. Schneider, Eduardo L. Martin, Andrzej Niedzielski, Jiwei Xie, Scott W. Fleming, Neil Thomas, Michael Williamson, Zhaohuan Zhu, Eric Agol, Dmitry Bizyaev, Luiz Nicolaci da Costa, Peng Jiang, A.F. Martinez Fiorenzano, Jonay I. Gonzalez Hernandez, Pengcheng Guo, Nolan Grieves, Rui Li, Jane Liu, Suvrath Mahadevan, Tsevi Mazeh, Duy Cuong Nguyen, Martin Paegert, Sirinrat Sithajan, Keivan Stassun, Sivarani Thirupathi, Julian C. van Eyken, Xiaoke Wan, Ji Wang, John P. Wisniewski, Bo Zhao, Shay Zucker
aa r X i v : . [ a s t r o - ph . E P ] A ug Very Low-Mass Stellar and Substellar Companions to Solar-likeStars From MARVELS VI: A Giant Planet and a Brown DwarfCandidate in a Close Binary System HD 87646
Bo Ma( 馬 波 ) , Jian Ge , Alex Wolszczan , Matthew W. Muterspaugh , , Brian Lee ,Gregory W. Henry , Donald P. Schneider , , Eduardo L. Mart´ın , Andrzej Niedzielski ,Jiwei Xie , , Scott W. Fleming , , Neil Thomas , Michael Williamson , , Zhaohuan Zhu Eric Agol , Dmitry Bizyaev , Luiz Nicolaci da Costa , , Peng Jiang , , A.F. MartinezFiorenzano , Jonay I. Gonz´alez Hern´andez , , Pengcheng Guo , Nolan Grieves Rui Li ,Jane Liu , Suvrath Mahadevan , , Tsevi Mazeh Duy Cuong Nguyen , Martin Paegert ,Sirinrat Sithajan , Keivan Stassun , Sivarani Thirupathi , Julian C. van Eyken , XiaokeWan , Ji Wang , John P. Wisniewski , Bo Zhao , Shay Zucker [email protected] Department of Astronomy, University of Florida, 211 Bryant Space Science Center, Gainesville, FL,32611-2055, USA Department of Astronomy and Astrophysics, The Pennsylvania State University, University Park, PA16802, USA Department of Mathematical Sciences, College of Life and Physical Sciences, Tennessee State University,Boswell Science Hall, Nashville, TN 37209, USA Center of Excellence in Information Systems Engineering and Management, Tennessee State University,3500 John A. Merritt Blvd., Box No. 9501, Nashville, TN 37209-1561, USA Department of Physics & Astronomy, Vanderbilt University, Nashville, TN 37235 USA Center for Exoplanets and Habitable Worlds, The Pennsylvania State University, University Park, PA16802, USA Centro de Astrobiolog´ıa (INTA-CSIC), Carretera de Ajalvir km 4, E-28550 Torrej´on de Ardoz, Madrid,Spain Toru´n Centre for Astronomy, Nicolaus Copernicus University in Toru´n, Grudziadzka 5, 87-100 Toru´n,Poland Department of Astronomy & Key Laboratory of Modern Astronomy and Astrophysics in Ministry ofEducation, Nanjing University, Nanjing 210093, China Apache Point Observatory and New Mexico State University, P.O. Box 59, Sunspot, NM, 88349-0059,USA Laborat´rio Interinstitucional de e-Astronomia (LIneA), Rio de Janeiro, RJ, 20921-400, Brazil Observat´orio Nacional, Rua General Jos´e Cristino, 77, 20921-400 S˜ao Crist´ov˜ao, Rio de Janeiro, RJ,Brazil Department of Physics, UC Santa Barbara, Santa Barbara, CA 93106-9530, USA Space Telescope Science Institute, 3700 San Martin Dr., Baltimore, MD, USA 21218 Computer Sciences Corporation, 3700 San Martin Dr., Baltimore, MD, USA 21218 Key Laboratory for Research in Galaxies and Cosmology, The University of Science and Technology ofChina, Hefei, Anhui 230026, China School of Physics and Astronomy, Raymond and Beverly Sackler Faculty of Exact Sciences, Tel AvivUniversity, Tel Aviv 69978, Israel Dunlap Institute for Astronomy and Astrophysics, University of Toronto, Toronto, ON, M5S 3H4,Canada Department of Astronomy, Yale University, New Haven, CT 06511, USA Department of Geosciences, Raymond and Beverly Sackler Faculty of Exact Sciences, Tel Aviv University,6997801 Tel Aviv, Israel HL Dodge Department of Physics & Astronomy, University of Oklahoma, 440 W Brooks St, Norman,
ABSTRACT
We report the detections of a giant planet (MARVELS-7b) and a brown dwarfcandidate (MARVELS-7c) around the primary star in the close binary system,HD 87646. It is the first close binary system with more than one substellarcircum-primary companion discovered to the best of our knowledge. The detec-tion of this giant planet was accomplished using the first multi-object Dopplerinstrument (KeckET) at the Sloan Digital Sky Survey (SDSS) telescope. Subse-quent radial velocity observations using ET at Kitt Peak National Observatory,HRS at HET, the “Classic” spectrograph at the Automatic Spectroscopic Tele-scope at Fairborn Observatory, and MARVELS from SDSS-III confirmed thisgiant planet discovery and revealed the existence of a long-period brown dwarfin this binary. HD 87646 is a close binary with a separation of ∼
22 AU betweenthe two stars, estimated using the Hipparcos catalogue and our newly acquiredAO image from PALAO on the 200-inch Hale Telescope at Palomar. The pri-mary star in the binary, HD 87646A, has T eff = 5770 ± g =4.1 ± − . ± .
08. The derived minimum masses of the two substellar com-panions of HD 87646A are 12.4 ± Jup and 57.0 ± . Jup . The periods are13.481 ± ± ± ±
1. INTRODUCTION
One of the most surprising astronomical developments of the last 25 years has beenthe discovery of an abundant population of extra-solar planets and brown dwarfs (BDs)(Wolszczan & Frail 1992; Mayor & Queloz 1995; Rebolo et al. 1995; Nakajima et al. 1995).
OK 73019, USA Instituto de Astrof´ısica de Canarias, E-38205 La Laguna, Tenerife, Spain Universidad de La Laguna, Dpto. Astrof´ısica, E-38206 La Laguna, Tenerife, Spain Fundaci´on Galileo Galilei-INAF, Rambla Jos´e Ana Fernandez P´erez, 738712 Bre˜na Baja, Tenerife, Spain Department of Astronomy, University of Washington, Box 351580, Seattle, WA 98195-1580, USA Department of Astrophysics, Princeton University, Princeton, NJ, 08544, USA
2. HD 87646
The target star, HD 87646, is a bright (V=8) G-type star with a fainter K-type stellarcompanion at a separation of 0.213 arcseconds and a position angle (PA) of 136 degreesaccording to the Hipparcos and Tycho Catalogues (Perryman 1997). Its Hipparcos parallaxof 13 . ± .
58 milliarcsecond (mas) places it at a distance of 73 . ± .
68 pc. Photom-etry and high resolution spectroscopic observation of HD 87646A have been obtained byFeltzing & Gustafsson (1998). They obtained a effective temperature of T eff = 5961 K fromphotometry, and spectroscopically derived log g = 4 .
41 and [Fe / H] = 0 .
3. This star is quitemetal rich, prompting Gonzalez et al. (2001) to explicitly recommend that it be observedwith precise radial velocity instruments due to the significantly higher probability of discov-ering hot Jupiter planets around metal rich stars (Fischer & Valenti 2005). While HD 87646was observed as part of the Geneva Copenhagen Survey (Nordstrom et al. 2008) we arepresently unaware of any high precision radial velocity observation, or the star being part ofany ongoing surveys.Such binary systems are challenging for precision radial velocity detection due to thecomplexity of analyzing spectra from two different stars. While detections of exoplanets inunresolved stellar binaries have been reported before (Konacki 2005), higher precision ob-servations and better cadence have not detected the same signal (Eggenberger et al. 2007).We speculate that this difficulty in detection is the reason this target was not observed insome ongoing surveys, like N2K (Robinson et al. 2006), that target high-metallicity stars.Multi-object surveys do not need to be as selective due to their inherent multiplicity advan-tage. Binaries may be excluded, but the existence of a few binaries among 60 stars observedsimultaneously is not a significant problem. In addition, once observations are well underwaythere is little advantage gained in removing the target since any replacement target wouldthen only be observed for a few epochs. We will study the impact of spectral contaminationfrom a faint companion star on the RV measurements for our target in Section 6.2. Ourstudy shows that the only substellar companions that can be detected in such close binariesare those massive enough to generate RV signals much larger than the noise induced by thespectral contamination.
3. The Multi-object KeckET pilot survey
The design of KeckET is based on a single object Exoplanet Tracker (ET) design forthe KPNO 2.1m telescope (Ge et al. 2003, 2006; Mahadevan et al. 2008). This instrumentadopts the dispersed fixed-delay interferometry (DFDI) approach for Doppler measurements(Erskine & Ge 2000; Ge 2002; Ge et al. 2002). Instead of the line centroid shifts in the 6 –high resolution cross-dispersed echelle spectrograph approach, the DFDI method measuresthe Doppler motion by monitoring the fringe shifts of stellar absorption lines created by aMichelson-type interferometer with a fixed-delay between the two interferometer arms. Themeasurement of this fixed-delay is described in Wang et al. (2012b,a).The KeckET instrument consists of 8 subsystems–a multi-object fiber feed, an iodinecell, a fixed-delay interferometer system, a slit, a collimator, a grating, a camera, and a 4k × µ m core diameters, whichare coupled to 180 µ m core diameter short fibers from the SDSS telescope, correspondingto 3 arcsec on the sky at f / ∼
900 ˚A, centered at 5400 ˚A. Details of theinstrument design can be found in Ge et al. (2006b), Wan et al. (2006), and Zhao & Ge(2006). KeckET has one spectrograph and one 4k ×
4k CCD camera that captures one ofthe two interferometer outputs, and has a 5.5% detection efficiency from the telescope tothe detector without the iodine cell under the typical APO seeing conditions ( ∼ ± − .The corresponding average photon-limit error is 5.5 ± − . The instrument’s precisionover longer time intervals has been measured with repeated observations of sky scatteredlight over a period of 45 days in Fall 2006, and 150 days in Winter/Spring 2007. The rmsdispersion of RV measurements of the sky over these periods, after subtracting the photonlimiting errors in quadrature, are 11.7 ± − and 11.3 ± − , respectively.The instrumental contributions to random measurement errors are mainly caused by in-homogeneous illumination of the slit, image aberration, and the interferometer comb aliasing(sampling on the detector). However, the dominant measurement RV error is produced bythe mathematical approximation used for extracting iodine and stellar Doppler signals in themixed stellar and iodine fringing spectra, which is on the order of 50 m s − (van Eyken, Ge 7 –& Mahadevan 2010), which is included in the RV errors showing in the data table. Althoughthis error has largely limited our capability of detecting relatively low mass planets, it doesnot affect the Doppler detection of massive giant planets, brown dwarfs and binaries.
4. Survey Data Processing and RV Results
The pipeline processing steps are described in detail in van Eyken et al. (2004), Ge etal. (2006), Mahadevan et al. (2008) and van Eyken et al. (2010). The data were processedusing standard IRAF procedures (Tody 1993), as well as software written in IDL. The imageswere corrected for biases, dark current, and scattered light and then trimmed, illuminationcorrected, slant corrected and low-pass filtered. The visibilities ( V ) and the phases ( θ )of the fringes were determined for each channel by fitting a sine wave to each column ofpixels in the slit direction. To determine differential velocity shifts the star+iodine data canbe considered as a summation of the complex visibilities ( V = V e iθ ) of the relevant star( V S e iθ S ) and iodine ( V I e iθ I ) templates ( Erskine 2003; van Eyken et al. 2004, 2010). Forsmall velocity shifts the complex visibility of the data, for each wavelength channel, can bewritten as V D e iθ D = V S e iθ S e iθ S − iθ S + V I e iθ I e iθ I − iθ I , (1)where V D , V S , and V I are the fringe visibilities for a given wavelength in the star+iodine data,star template, and iodine template respectively, and θ D , θ S , and θ I are the correspondingmeasured phases. In the presence of velocity shifts of the star and instrument drift, thecomplex visibilities of the star and iodine template best match the data with phase shifts of θ S − θ S and θ I − θ I , respectively. The iodine is a stable reference and the iodine phase shifttracks the instrument drift. The difference between star and iodine shifts is the real phaseshift of the star, ∆ φ , corrected for any instrumental drifts∆ φ = ( θ S − θ S ) − ( θ I − θ I ) . (2)This phase shift can be converted to a velocity shift, ∆ v , using a known phase-to-velocityscaling factor: ∆ v = cλ πd ∆ φ, (3)where c is the speed of light, λ is the wavelength and d is the optical delay in the Michelsoninterferometer. The KeckET data analysis pipeline identifies the shift in phase of the starand iodine templates that are the best match for the data, and uses these phase shifts tocalculate the velocity shift of the star relative to the stellar template. Since HD 87646 is aclose binary system, we need an additional complex visibility term in equation 1 to accountfor the contamination of the secondary star. Mathematically this is equivalent to adding a 8 –small noise term δθ S to the phase of the primary star θ S (van Eyken, Ge & Mahadevan 2010),which will translate to the measured RV according to equation 3. This noise term dependsmainly on the flux ratio of the two stars collected through the fiber and the radial velocityoffset between the two stars. Thus this noise term varies slowly between observations. Tosimplify the RV fitting process, we treat this noise as a constant value and will study itsimpact on the RV measurements of HD 87646 in section 6.2. In practice, one can try tomodel the visibility and phase of the secondary star if both star spectra and their flux ratiovariations with wavelength are known precisely. This method is not very practical in ourcurrent case because of the lack of these information.We have obtained a total of 16 observations of HD 87646 using KeckET from 2006December to 2007 June. The radial velocities obtained are listed in Table 1.
5. Follow-up Observations5.1. KPNO ET RV Observations
Subsequent observations were performed using the Exoplanet Tracker (ET) instrumentat KPNO (Ge et al. 2006b). Initial follow-up was performed in November of 2007, whichconfirmed the variability seen in the KeckET data. Additional data points were obtained atKPNO in 2008 January, 2008 February and 2008 May. The integration time was 35-40 minsin 2007 November and 20 mins in 2008 January, 2008 February and 2008 May.The data were reduced using software described in Mahadevan et al. (2008) and refer-ences therein. See van Eyken, Ge & Mahadevan (2010) for the theory behind the technique.A total of 40 data points were obtained from 2007 November to 2008 May and are listed inTable 2. The observations confirmed the linear trend shown in the KeckET data, which islater to be found due to another substellar companion.
Follow-up observations of HD 87646 were conducted with the fiber-fed High ResolutionSpectrograph (HRS, Tull 1998) of the Hobby Eberley telescope (HET, Ramsey et al. 1998).The observations were executed in queue scheduled mode (Shetrone et al. 2007), and used a 2arcsecond fiber and with the HRS slit set to yield a spectral resolution of R ∼ , − m s − −
592 nm)and 24 orders on the red one (602 −
784 nm). The spectral data used for RV measurementswere extracted from the 17 orders (505 −
592 nm) in which the I cell superimposed strongabsorption lines. The radial velocities obtained are listed in Table 3. HD 87646 was selected as an RV survey target by the MARVELS preselection criterion(Paegert et al. 2015). The star has been monitored at 23 epochs using the MARVELS instru-ment mounted on the SDSS 2.5m Telescope at APO between 2009 May and 2011 December(Ge et al. 2008, Ma et al. 2013). The MARVELS instrument is a fiber-fed dispersed fixed-delay interferometer instrument capable of observing 60 objects simultaneously and coversa wavelength range of 5000 − ∼ To investigate the nature of the linear RV trend found in previous RV data, we haveobtained additional observations of HD 87646 with a fiber-fed echelle spectrograph situatedat the 2 m Automatic Spectroscopic Telescope (AST) in Fairborn Observatory (Eaton &Williamson 2004, 2007). The robotic nature of the AST allowed for high cadence obser-vations, which removed orbital period degeneracies and helped solidify the longer-periodcompanion’s orbit. Through 2011 June the detector was a 2048 × − ×
4K 15 micron pixels, which requireda new readout electronics package, and a new dewar with a Cryotiger refrigeration system.The echelle spectrograms that were obtained with this new detector have 48 orders, coveringthe wavelength range 3800 − − m s − − m s − A total of nine high-resolution spectra ( R = 164 , eff and v sin i ), and to search for evidence of a secondset of lines in the system. The typical signal-to-noise ratio (S/N) for each spectrum is about150 per resolution element around 5500˚A. We also obtained spectroscopic observations of HD 87646 from the 2.1m telescope atKitt Peak National Observatory using the R = 30 ,
000 Direct Echelle Mode of the EXPERTspectrograph (Ge et al. 2010). A total of seven EXPERT spectra were acquired between2014 Feb to 2014 June. The exposure time for each observation ranged from 20-40 minutes,yielding an
S/N ∼
250 per resolution element around 5500˚A.
On 2008 May 29, high angular resolution lucky images of HD 87646 were obtained withthe FastCam instrument (Oscoz et al. 2008) on the 2.5-meter Nordic Optical Telescope atthe Roque de los Muchachos Observatory in La Palma (Spain). Five data cubes of 1000images each were obtained using the I-band filter. Individual exposure times were 30 ms foreach image. High spatial resolution was obtained by combining the best 1% of the images.Fig. 1 shows the processed image. The image scale was 30.95 ± − . In thisfigure the secondary star HD 87646B in this binary system is visible. The point spreadfunction (PSF) of the star has a full width half maximum (FWHM) value of 0 . − m s − − ) Velocity Error (ms − )2454903.7521 21690 2122454904.8249 20810 2892454905.6178 21580 1872454906.7753 21470 2582454908.7537 21290 1422454909.7757 21530 1792454910.7670 21600 3422454911.7909 22720 1712454913.7546 22740 2142454914.7619 22140 217 14 – On 2009 June 4, we acquired high resolution AO images of the binary star systemHD 87646 from PALAO on the 200-inch Hale Telescope at Palomar. Data sets were takenin the J and K bands. The AO system was running at 500 Hz. The seeing was roughly1.3” in K band. We also observed a calibrator star and subtracted the K-band’s PSF, whichimproved sensitivity by a factor of a few. The utility of PSF subtraction is limited, in thiscase, by the difference in stellar spectral types, since the filters are broad.Images were flat-fielded, background subtracted and cleaned; the final images are dis-played on a logarithmic scale in Fig. 2 with a scale of 25 mas pixel − . The binary systemis well-resolved in both bands. The angular separation is measured to be 401 ±
12 mas,nearly twice that quoted from the Hipparcos and Tycho Catalogues (Perryman 1997). Theposition angle is 69 . ◦ ± . ◦ . There is no evidence in the high resolution images for a tertiary(stellar) companion. The J and K-band brightness ratios of the two stars are 6 . ± . . ± .
10, respectively. We can not use these ratios and their error bars to put ameaningful constraint on the optical band flux ratio because the spectral energy distribution(SED) curve slope is basically flat in the J and K-bands, but very sharp in the optical band.
Photometry of HD 87646 was obtained in the Stromgren b and y bands between 2007December and 2015 June with the T12 0.8-m Automatic Photoelectric Telescope (APT) atFairborn Observatory in Arizona. Our primary goal with photometry is to detect if thecompanions transit the primary star. The data were processed using software described inHenry (1999). The measurements have a typical accuracy of ∼
6. Results6.1. Stellar Parameters
The SARG spectra taken at TNG without the iodine cell were used to derive the stel-lar parameters. HD 87646 is flagged as a binary in the Hipparcos and Tycho Catalogues(Perryman 1997), with a Hipparcos magnitude (a broad band V filter) difference betweenthe primary and secondary to be 2.66 ± − and the star PSF has a FWHM of 0 . ± − . The FWHMs are 70 mas and 120 mas forJ and K band images respectively. 16 –the normalized spectra show minor changes, only affecting the wings of strong lines (e.g theMgb lines, Fig. 3). The equivalent widths of most weak lines are essentially unchanged, sowe used Fe I and FeII lines with equivalent widths below 140m˚A and performed traditionalspectroscopic analysis.We use the latest MARCS model atmospheres (Gustafsson et al. 2008) for the analysis.Generation of synthetic spectra and the line analysis were performed using the turbospectrumcode (Alvarez & Plez 1998), which employs line broadening according to the prescription ofBarklem & O’Mara (1998). The line lists used are drawn from a variety of sources. Atomiclines are taken mainly from the VALD database (Kupka et al. 1999). The molecular speciesCH, CN, OH, CaH and TiO are provided by B. Plez (see Plez & Cohen 2005), while theNH, MgH and C molecules are from the Kurucz line lists. The solar abundances used hereare the same as Asplund (2005). We use FeI excitation equilibrium and derived an effectivetemperature T eff = 5770 ± α and H β wings also agree better forthis lower T eff value. We find log g =4.1 ± − is derived by forcing weak and strong FeI lines to yieldthe same abundances. We are not able to confirm the super solar metallicity of this object(Feltzing & Gustafsson 1998); we derived [Fe/H] = − . ± .
08. When we adopt the same T eff and microturbulence as Feltzing & Gustafsson (1998), we obtain the same metallicityvalue [Fe/H]=0.3, but with a large slope in the excitation potential versus FeI abundanceand reduced equivalent width versus FeI abundance. The derived stellar parameters aresummarized in Table 6.We attempted to place constraints on the secondary star by fitting the Balmer lineand Mgb line wings. Based on the Hipparcos data, which suggest a Hipparcos magnitudeTable 6: Parameters of the Star HD 87646AParameter Value T eff ±
80 Klog( g ) 4.1 ± / H] -0.17 ± V sin i − ξ t − Mass 1.12 ± M ⊙ Radius 1 . ± . R ⊙
17 –difference between the primary and secondary 2.66 ± T eff , log( g ), and [Fe/H].These relations were derived from a sample of eclipsing binaries with precisely measuredmasses and radii. We estimate the uncertainties in M ∗ and R ∗ by propagating the uncertain-ties in T eff , log( g ), and [Fe/H] using the covariance matrices of the Torres, Andersen & Gim´enez(2010) relations (kindly provided by G. Torres). Since the polynomial relations of Torres, Andersen & Gim´enez(2010) were derived empirically, the relations were subject to some intrinsic scatter, which weadd in quadrature to the uncertainties propagated from the stellar parameter measurements( σ log m = 0 .
027 and σ log r = 0 . M ∗ = 1 . ± . M ⊙ and R ∗ = 1 . ± . R ⊙ . HD 87646 is a binary system, and contamination of the primary star’s spectrum by thesecondary star leads to an increased RV jitter by interfering with the analysis pipeline. Inthis section we investigate the possible systematic RV errors caused by the blended binaryspectra using simulations. Since our RV observations were produced by two different kindsof spectrographs, we decided to perform two simulations, one for the DFDI instruments,including KeckET, KPNO ET, and MARVELS, and the other for traditional echelle spec-trographs, including HRS at HET and the AST fiber-fed echelle spectrograph at Fairborn.In both simulations, we first create a set of stellar spectra by combining a G-type star (forthe primary) and a K-type (for the secondary) star spectra with varying radial velocitiesfor both stars. Then we calculate the differential radial velocities for the G-star from thesimulated spectra. The differences between the output G-star RVs and input G-star RVsare the RV errors caused by secondary star spectra contamination. Both simulations yieldsimilar RV errors on the order of 200 m s − . We expect to see this level of systematic errorand will include it in the RV ‘jitter’ term when we perform the RV curve fitting in the nextsection. Traditionally the ‘jitter’ term used to denote any RV noise caused by stellar activity;our ‘jitter’ term also contains the RV noise caused by blended binary spectra. 18 – We have performed a Markov Chain Monte Carlo (MCMC) analyses of the combinedRV observations from KeckET, ET, HET, MARVELS and Fairborn instruments. In thisanalyses, we initially used a one planet RV model to fit our RV observations, and later foundthat there is another strong periodic RV signal in the RV residuals. We then adopted a twoobject (a planet and a brown dwarf) RV model to fit our RV data. The RV model detailsare presented in § ~θ = { P , K , e , ω , M , P , K , e , ω , M , C i , σ jitter } , (4)where P and P are orbital periods, K and K are the radial velocity semi-amplitudes, e and e are the orbital eccentricities, ω and ω are the arguments of periastron, M and M are the mean anomalies at chosen epoch ( τ ), C i is constant velocity offset between thedifferential RV data shown in Tables 1, 2, 3, 4, and 5 and the zero-point of the KeplerianRV model ( i = 1 for KeckET observations, i = 2 for KPNO ET observations, i = 3 for HETobservations, i = 4 for MARVELS observations, and i = 5 for Fairborn observations), and σ jitter is the “jitter” parameter. The jitter parameter describes any excess noise, includingboth astrophysical noise (e.g. stellar oscillation, stellar spots; Wright 2005), any instrumentnoise not accounted for in the quoted measurement uncertainties and systematic RV errorsfrom analyzing blended binary spectra discussed in the last section. We use standard priorsfor each parameter (see Gregory 2007). The prior is uniform in the logarithm of the orbitalperiod ( P and P ) from 1 to 5000 days. For K , K and σ jitter we use a modified Jefferysprior which takes the form of p ( x ) = ( x + x o ) − [log(1 + x max /x o ] − , where x o = 0 . − and x m ax = 2128 m s − (Gregory 2005). Priors for e and e are uniform between zero andunity. Priors for ω , ω , M and M are uniform between zero and 2 π . For C i , the priorsare uniform between min( v i )-5 km s − and max( v i )+5 km s − , where v i are the set of radialvelocities obtained from each of the four instruments ( i = 1 for KeckET observations, i = 2for KPNO ET observations, i = 3 for HET observations, i = 4 for MARVELS observation,and i = 5 for Fairborn observations). We verified that the chains did not approach thelimiting value of P , P , K , K and σ jitter . 19 –Following Ford (2006), we adopt a likelihood (i.e., conditional probability of making thespecified measurements given a particular set of model parameters) of p ( v | ~θ, M ) ∝ Y k exp[ − ( v k,θ − v k ) / σ k,obs + σ )] p σ k,obs + σ jitter2 , (5)where v k is observed radial velocity at time t k , v k,θ is the model velocity at time t k giventhe model parameters ~θ , and σ k,obs is the measurement uncertainty for the radial velocityobservation at time t k .We combine the Markov chains described above to estimate the joint posterior prob-ability distribution for the orbital model for HD 87646. In Table 7 we report the medianvalue and an uncertainty estimate for each model parameter based on the marginal posteriorprobability distributions. The uncertainties are calculated as the standard deviation aboutthe mean value from the combined posterior sample. Since the shape of the marginal poste-rior distribution is roughly similar to a multivariate normal distribution, the median valueplus or minus the reported uncertainty roughly corresponds to a 68.3% credible interval. Inthe same table, we also reported the rms of the RV fitting residuals for the five different RVinstruments used here, which are rms = 245 m s − for KeckET, rms = 248 m s − for KPNOET, rms = 261 m s − for HET, rms = 270 m s − for MARVELS, and rms = 312 m s − for Fairborn RV observations.This two-Keplerian orbital solution is shown in Figs 4 and 5 together with the KeckET,HET, KPNO ET, Fairborn and MARVELS RV data. The residuals in these two plots can notbe explained only by the errors in our RV data. A stellar jitter term σ jitter = 240 ±
12m s − is required in our fitting to explain these residuals. As discussed in the last section, most ofthe ‘jitter’ noise arises from our data pipeline handling the blended binary spectra insteadof a single star spectra. We also did an MCMC analysis using two ‘jitter’ noise terms,one for DFDI instruments and the other one for echelle spectrographs, and find the orbitalparameters for the giant planet candidate and the BD candidate are barely changed withinthe error bars. So to keep it simple, we choose to use one ‘jitter’ term for all our RVobservations from different instruments.HD 87646 is a binary system, so we have done another MCMC analysis by including alinear RV trend (v trend × (t − t )) to the two-objects RV model used above. This linear RVtrend is used to account for the perturbation of the primary star induced by the gravitationalforce of the secondary star. We note here that the offsets between different data sets willhinder the modeling of this linear trend as there is expected strong correlation between theoffsets and this linear trend. Our new RV fitting yields orbital parameters for the two sub-stellar companions in addition to a linear RV trend of v trend = − ±
18 m s − yr, whichare summarized in Table 8. Since all the main orbital parameters of the two sub-stellar 20 –Fig. 3.— High resolution TNG spectra of HD 87646 centered around H α , H β and Mgb lines(black lines). The three other lines correspond to synthetic spectra for a binary with a Gdwarf primary ( T eff = 5770 K, log g =4.1, [Fe/H] =-0.17, Microturbulence=1.8km s − and V sin i =7.5km s − ) and a K dwarf secondary ( T eff = 4000 K, log g =5.0). The green, red andblue lines correspond to a 0%, 10% and 50% flux contribution from the secondary to thewhole binary.Fig. 4.— Top:
Expanded section of the middle panel plot (the dotted rectangular region) toshow the short-period giant planet signal in the two-Keplerian RV model.
Middle:
Radialvelocity observations of HD 87646 with the two-Keplerian model.
Bottom:
RV residuals ofthe two-Keplerian orbit model. Each panel shows radial velocity observations from KeckET(yellow stars), HET (black triangles), KPNO ET (red squares), Fairborn (blue diamonds),and MARVELS (red cross). 21 –Fig. 5.— Phased RV curves for the two signals in the two-Keplerian RV model. In eachcase, the contribution of the other signal was subtracted. Each panel shows radial velocityobservations from KeckET (yellow stars), HET (black triangles), KPNO ET (red squares),Fairborn (blue diamonds), and MARVELS (red cross).Table 7: Orbital Parameters for HD 87646b and HD 87646cParameter HD 87646b HD 87646cMinimum Mass 12.4 ± M Jup ± M Jup a ± ± K ±
25 m s − ±
54 m s − P ± ± e ± ± ω (radians) 5.20 ± ± T prediction for transit (JD UTC ) 2454093.85 ± T periastron (JD UTC ) 2454088.3 ± . ± σ jitter ±
12 m s − C ± − C ± − C ± − C -0.306 ± − C ± − rms
245 m s − rms
248 m s − rms
261 m s − rms
270 m s − rms
312 m s −
22 –Table 8: Orbital Parameters for HD 87646b and HD 87646c with a Linear RV TrendParameter HD 87646b HD 87646cMinimum Mass 12.4 ± M Jup ± M Jup a ± ± K ±
24 m s − ±
56 m s − P ± ± e ± ± ω (radians) 5.18 ± ± T prediction for transit (JD UTC ) 2454093.86 ± T periastron (JD UTC ) 2454088.2 ± . ± trend − ±
18 m s − yr − σ jitter ±
13 m s − C ± − C ± − C ± − C -0.262 ± − C ± − rms
246 m s − rms
248 m s − rms
261 m s − rms
269 m s − rms
312 m s −
23 –companions are barely changed within their respect error bars and the strong correlationbetween this RV trend and telescopes RV offsets, we decide to keep using the numberspresent in Table 7 throughout this paper. This linear trend is not significant, which meansit is more likely that either the secondary star is close to its ascending or descending nodeduring 2008-2013, or this binary is on a relatively low-inclination (face on) orbit. It is notpossible for us to distinguish these two scenarios using our current data. Future high precisionastrometry observations, like GAIA, will help to solve this binary orbital problem. We alsowant to note here that this linear trend is not exact the real RV trend of the primary starinduced by the gravitational perturbation of the secondary star because of the secondarystar’s spectral contamination. It is close to ∼ ∼
10 in the optical band and mass ratio is ∼ Santos et al. (2002) found small radial velocity variations and line asymmetries for thestar HD 41004, which is a visual binary and is unresolved at the spectrograph. It was initiallythought to have a planetary companion around the primary star, but from the line bisectoranalysis they were able to infer a possible brown dwarf orbiting the secondary star insteadof a planet orbiting the primary star. Their conclusions were subsequently corroborated byZucker et al. (2003).Similar to HD 41004, HD 87646 is also a binary system with a small angular separa-tion (0 . ′′ ), which renders the spectrum a blended spectrum of the two stellar components.Following the same philosophy of Santos et al. (2002), we performed a bisector analysis forHD 87646 to determine from which star in the binary system the RV signal was produced.We have analyzed spectra taken at the Kitt peak 2m telescope using EXPERT (Ge et al.2010). Spectra were reduced using an IDL pipeline modified from an early version describedin Wang (2012). Frames were trimmed, bias subtracted, flat-field corrected, aperture-traced,and extracted. Cross-Correlation Functions (CCFs) are derived by cross-correlation with aspectral mask from the wavelength range 4900 − The top panel of Fig. 7 presents all 1077 photometric observations plotted against thelatest transit ephemeris of HD 87646a: T c = 2454093 .
85 d, P = 13 .
481 d. The differentialmagnitudes are measured against the mean of three comparison stars to improve precision.The standard deviation of these data from their mean is 0.0014 mag. A least-squares sine-curve fit to the phased data yields a full amplitude of 0.000089 ± ± . ∼ .
015 units of phase and adepth of 0 .
5% or ∼ .
005 mag (Kane & von Braun 2008). The ± σ uncertainty in the transitwindow timing is indicated by the two vertical dotted lines. There are 1005 observation thatlie outside the predicted transit window, which have a mean of 0.99998 ± ± . ± . ∼ .
005 mag are excluded by the photometry at the predicted transit time.
The mass function is related to the observed period, eccentricity and radial velocitysemi-amplitude as: ( m sin i ) ( M ∗ + m ) = P (1 − e ) K πG (6)where M ∗ is the mass of the primary and m the mass of the companions. Since the firstcompanion is known not to transit the star, we cannot break the degeneracy of mass andsin i with radial velocity observations alone. Using the derived stellar mass (1.12 M ⊙ ) forthe primary with the orbital parameters determined from the radial velocity (Table 7) wedetermine that the minimum mass of the inner companion for an edge-on orbit (sin i = 1) 25 –Fig. 6.— Measured radial-velocity vs. BIS from EXPERT spectroscopic data for HD 87646.The solid and dotted lines show simulation results when assuming a giant planet orbitingthe primary star and the secondary star, respectively.Fig. 7.— Top: the 1077 differential magnitudes of HD 87646 phased to the period of thegiant planet HD 87646Ab, taken by the 0.8m APT from 2008-2015. The horizontal dashedline corresponds to the mean brightness level of the 1077 observations. The vertical dashedline marks the expected time of mid-transit. Bottom: an expanded portion of the topplot, centered on the predicted central transit window. The solid curve shows the predictedcentral transit, with a depth of 0.005 mags and duration of 0.015 units of phase. The ± σ uncertainty in the transit window timing is indicated by the two vertical dotted lines. 26 –is 12.4 ± Jup . This mass is quite close to the deuterium burning limit, and the detectedcompanion is likely burning deuterium, although its minimum mass places it in the giantplanet regime. The second companion’s minimum mass when assuming an edge-on orbit is57.0 ± . Jup , which falls right into the brown dwarf regime.
7. Summary and Discussion7.1. Summary of the main results
Our SDSS MARVELS pilot survey and additional observations at the HET, KPNO2.1m telescope, and Fairborn observatory confirm the detection of two massive substellarcompanions in a close binary system HD 87646. The first companion, HD 87646Ab, has aminimum mass of 12.4 ± Jup , period of 13.481 ± ± ± . Jup , period of 674 ± ± > Jup ) BDs tend tohave higher eccentricities. This new BD is consistent with this trend.This is the eleventh detection of a substellar companion(s) in a binary system with sep-aration of only about 20 AU. The other ten systems are Gliese 86 (Queloz et al. 2000; Mu-grauer & Newhauser 2005; Lagrange et al. 2006), γ Cephei (Hatzes et al. 2003), HD 41004(Zucker et al. 2003, 2004), HD 188753 (Konacki 2005), HD 176051 (Muterspaugh et al.2010), HD 126614 (Howard et al. 2010), α Centauri (Dumusque et al. 2012), HD 196885(Correia et al.2008), OGLE-2013-BLG-0341 (Gould et al. 2014) and HD 59686 (Ortiz et al. 2016). How-ever, Eggenberger et al. (2007) did not confirm the planet in HD 188753, and Rajpaul, Aigrain, & Roberts(2016) suggest the planet signal discovered from α Centauri B is not from a real planet, butfrom the observation window function. To the best of our knowledge, HD 87646A is the firstmultiple planet/BD system detected in such close binaries.
In this section we will discuss the dynamical stability of the binary system HD 87646.First we have collected observational data of HD 87646 from the literature (Horch et al. 2008,Hartkopf & Mason 2009, Horch et al. 2010, Balega et al. 2013) and combined them with our 27 –AO data to constrain the binary orbit of HD 87646. Our best fitting binary orbital solution(solid line) and the observational data (black dots) are shown in Fig. 8. The best fittingparameters are P = 51 . e = 0 .
54 and a = 0 .
26 arcsec. At a distance of 73 . ± .
68 pc,this angular separation corresponds to a binary semi-major axis of 19 ± B ) and eccentricity (e B ) of the binarysystem up to a million years. We have assumed the giant planet, brown dwarf and binaryto be coplanar. The results are shown in Fig. 9. There is a stable zone in the binary a B -e B diagram. We have scaled the error bars of these astrometry data to force the best orbital fitto have a reduced chi-squared χ = 1 and then over-plotted the binary orbital parameterfit from astrometry data with max χ = 2 in Fig. 9. The big uncertainty of the binaryorbital fitting from astrometry data arises from the big error bar ( ∼ . χ = 1 can be problematic. The main conclusion from the simulation study andastrometry data fitting is, with a large binary semi-major axis (a B >
17 AU) and a relativelylow binary eccentricity (e B < (a B − × .
57 + 0 . HD 87646 is the first known system to have two massive substellar objects orbitinga star in a close binary system. Interestingly, the masses of these two substellar objectsare close to the minimum masses for burning deuterium ( ∼
13 M
Jup , Spiegel et al. 2011) 28 –Fig. 8.— Orbital data for binary HD 87646. The plus symbol marks the location of theprimary (HD 87646A), filled circles are measured position of HD 87646B from literatureand this paper, and line segments are drawn from the ephemeris prediction to the observedlocation of the secondary in each case.Fig. 9.— Dynamical simulation results for HD 87646. The contour lines show the time thesystem will remain stable according to our Mecury simulation. The unit of the color baris years. A stable zone is found in the eccentricity-semimajor axis diagram, which showsthat if the binary orbit has a large semi-major axis and low eccentricity, the system willremain stable. The white triangle symbol shows the best binary orbital fit from the currentastrometry data with χ = 1 after we rescale the error bars for astrometric data (seealso Fig. 8). The white ellipse shows shows the distribution of binary orbital parametersfrom fitting the astrometric data fitting with max χ = 2. There is an overlap regionbetween the distribution of binary orbital parameters and the dynamically stable region,which demonstrates the binary system has stable orbital solutions. 29 –and hydrogen ( ∼ Jup , Chabrier et al. 2000) which are generally assumed to be thegeneral mass boundaries between planet and brown dwarf and between brown dwarf andstar, respectively. All these peculiarities raise a question: how could such a system beformed? Here we briefly discuss this intriguing issue.The large masses of these two substellar objects suggest that they could be formedas stars with their binary hosts: a large molecular cloud collapsed and fragmented intofour pieces; the larger two successfully became stars and formed the HD 87646 binary,and the other smaller ones failed to form stars and became the substellar objects in thissystem (Chabrier et al. 2014). This scenario might be relevant for the binary stars butseems problematic for the two substellar objects on orbits within ∼ ∼
20 AU (seethe review by Thebault & Haghighipour 2014, and the references therein). A commonlyrecognized issue is that the binary perturbations generally inhibit the growth of planetesi-mals in the disk (Thebault 2011). Even if their growth could proceed in certain favorableconditions (Xie & Zhou 2008, 2009), it would be significantly slowed, requiring a time scaleof 10 yr or even longer (Xie et al. 2010). This result raises a problem for the formationof a gaseous giant planet, as it would not form a planetary core (via planetesimal growth)to accrete gas before the gas disk dissipation, which takes a timescale as short as 10 -10 yr for such close binaries (Cieza et al. 2009; Kraus et al. 2012). Following the above logic,Xie, Zhou & Ge (2010) found that Jupiter-like planets are unlikely to form around AlphaCentauri B. As for the case of HD 87646, the formation of the two massive substellar objectsvia the core accretion model would be more problematic because it requires a more massivedisk with mass larger than 68 M Jup . Such a massive disk is seldom observed in close binaries,which indicates that should such a massive disk exist, it would dissipate much faster than anormal lighter one.Conversely, the disk instability model could circumvent most of the above barriers.First, disk instability usually requires very high disk mass, which is in line with the massesof the two detected substellar objects. Second, planet formation via disk instability requiresa short timescale, which is also consistent with the short disk dissipation timescale observedin close binaries. In addition, the model of disk instability is recently advocated by Duchˆene(2010), who argued that planet formation might be dominated by disk instability in close 30 –binaries based on the fact that exoplanets within close binaries (separation <
100 AU) aresignificantly more massive than those within wide binaries or single stars. We could use theplanet and brown dwarf mass to estimate the minimum surface density of the primordialdisk, and test if such a disk is gravitationally unstable. We adopt the similarity solutionof the evolving viscous disk (Hartmann et al. 1998) where the surface density follows R − from the star to the disk edge. Since the binary separation is 19 AU, the tidal truncationradius for the circumstellar disk around the primary star is ∼ . At 6 AU, the temperature is around 90 K with 1solar luminosity. The sound speed is 0.56 km/s. Then the Toomre Q parameter with a 1.12solar mass star is 1.5. At 1 AU, the temperature is around 220 K. The surface density is15610 g/cm REFERENCES
Alam, S., Albareti, F. D., Allende Prieto, C., et al. 2015, ApJS, 219, 12Alvarez M. et al., 1998, A&A, 340, 149Antonini, F., Hamers, A. S., & Lithwick, Y. 2016, arXiv:1604.01781D’Angelo, G., Durisen, R. H., & Lissauer, J. J. 2010, Exoplanets, 319Asplund M., 2005, ARA&A, 43, 481Balega I. I., Balega Y. Y., Gasanova L. T., Dyachenko V. V., Maksimov A. F., MalogolovetsE. V., Rastegaev D. A., Shkhagosheva Z. U., 2013, Astrophysical Bulletin, 68, 53Barklem P. S., O’Mara B. J., 1998, MNRAS, 300, 863 32 –Boss A. P., 2006, ApJ, 643, 501Butler R. P., Marcy G. W., Fischer D. A., Vogt S. S., Tinney C. G., Jones H. R. A., PennyA. J., Apps K., 2004, in IAU Symposium, Vol. 202, Planetary Systems in the Universe,Penny A., ed., p. 3Chabrier G., Baraffe I., 2000, ARA&A, 38, 337Chabrier G., Johansen A., Janson M., Rafikov R., 2014, Protostars and Planets VI, 619Chambers, J. E. 1999, MNRAS, 304, 793Chatterjee, S., Ford, E. B., Matsumura, S., & Rasio, F. A. 2008, ApJ, 686, 580-602Cieza L. A. et al., 2009, ApJ, 696, L84Correia, A. C. M., Udry, S., Mayor, M., et al. 2008, A&A, 479, 271Duchˆene G., 2010, ApJ, 709, L114Dumusque X. et al., 2012, Nature, 491, 207Eaton J. A., Williamson M. H., 2004, in Society of Photo-Optical Instrumentation Engineers(SPIE) Conference Series, Vol. 5496, Advanced Software, Control, and Communica-tion Systems for Astronomy, Lewis H., Raffi G., eds., pp. 710–717Eaton J. A., Williamson M. H., 2007, PASP, 119, 886Eggenberger A., Udry S., Mazeh T., Segal Y., Mayor M., 2007, A&A, 466, 1179Eisenstein D. J. et al., 2011, AJ, 142, 72Erskine D. J., 2003, PASP, 115, 255Erskine D. J., Ge J., 2000, in Astronomical Society of the Pacific Conference Series, Vol.195, Imaging the Universe in Three Dimensions, van Breugel W., Bland-HawthornJ., eds., p. 501Fabrycky, D., & Tremaine, S. 2007, ApJ, 669, 1298Fekel F. C., Rajabi S., Muterspaugh M. W., Williamson M. H., 2013, AJ, 145, 111Feltzing S., Gustafsson B., 1998, A&AS, 129, 237Fischer D. A., Valenti J., 2005, ApJ, 622, 1102 33 –Fleming S. W. et al., 2010, ApJ, 718, 1186Ford E. B., 2006, ApJ, 642, 505Perryman, M. A. C., de Boer, K. S., Gilmore, G., et al. 2001, A&A, 369, 339Ge J., 2002, ApJ, 571, L165Ge J., Erskine D. J., Rushford M., 2002, PASP, 114, 1016Ge J., Mahadevan S., van Eyken J., Dewitt C., Shaklan S., 2003, in Astronomical Society ofthe Pacific Conference Series, Vol. 294, Scientific Frontiers in Research on ExtrasolarPlanets, Deming D., Seager S., eds., pp. 573–580Ge J. et al., 2006a, ApJ, 648, 683Ge J. et al., 2006b, in Society of Photo-Optical Instrumentation Engineers (SPIE) Confer-ence Series, Vol. 6269, Society of Photo-Optical Instrumentation Engineers (SPIE)Conference Series, p. 2Ge J. et al., 2010, in Society of Photo-Optical Instrumentation Engineers (SPIE) Confer-ence Series, Vol. 7735, Society of Photo-Optical Instrumentation Engineers (SPIE)Conference Series, p. 0Gonzalez A. H., Zaritsky D., Dalcanton J. J., Nelson A., 2001, ApJS, 137, 117Gould A. et al., 2014, Science, 345, 46Gregory P. C., 2005, ApJ, 631, 1198Gregory P. C., 2007, MNRAS, 381, 1607Gunn J. E. et al., 2006, AJ, 131, 2332Gustafsson B., Edvardsson B., Eriksson K., Jørgensen U. G., Nordlund ˚A., Plez B., 2008,A&A, 486, 951Han, E., Wang, S. X., Wright, J. T., et al. 2014, PASP, 126, 827Hartkopf W. I., Mason B. D., Worley C. E., 2001, AJ, 122, 3472Hartkopf, W. I., & Mason, B. D. 2009, AJ, 138, 813Hartmann, L., Calvet, N., Gullbring, E., & D’Alessio, P. 1998, ApJ, 495, 385 34 –Hatzes A. P., 2013, ApJ, 770, 133Helled, R., Bodenheimer, P., Podolak, M., et al. 2014, Protostars and Planets VI, 643Henry G. W., 1999, PASP, 111, 845Horch E. P., Falta D., Anderson L. M., DeSousa M. D., Miniter C. M., Ahmed T., vanAltena W. F., 2010, AJ, 139, 205Horch E. P., van Altena W. F., Cyr, Jr. W. M., Kinsman-Smith L., Srivastava A., Zhou J.,2008, AJ, 136, 312Howard, A. W., Johnson, J. A., Marcy, G. W., et al. 2010, ApJ, 721, 1467Howard A. W. et al., 2014, ApJ, 794, 51Howard A. W. et al., 2010, Science, 330, 653Kane S. R., von Braun K., 2008, ApJ, 689, 492Konacki M., 2005, Nature, 436, 230Kratter K. M., Murray-Clay R. A., 2011, ApJ, 740, 1Kraus A. L., Ireland M. J., Hillenbrand L. A., Martinache F., 2012, ApJ, 745, 19Kupka F., Piskunov N., Ryabchikova T. A., Stempels H. C., Weiss W. W., 1999, A&AS,138, 119Lagrange A.-M., Beust H., Udry S., Chauvin G., Mayor M., 2006, A&A, 459, 955Lovis C. et al., 2006, Nature, 441, 305Ma B., et al., 2013, AJ, 145, 20Ma B., Ge J., 2014, MNRAS, 439, 2781Mahadevan S., van Eyken J., Ge J., DeWitt C., Fleming S. W., Cohen R., Crepp J., VandenHeuvel A., 2008, ApJ, 678, 1505Marcy G. W., Butler R. P., 1992, PASP, 104, 270Mayer L., Wadsley J., Quinn T., Stadel J., 2005, MNRAS, 363, 641Mayor, M., & Queloz, D. 1995, Nature, 378, 355 35 –Marzari, F., Weidenschilling, S. J., Barbieri, M., & Granata, V. 2005, ApJ, 618, 502McArthur B. E. et al., 2004, ApJ, 614, L81Mugrauer M., Neuh¨auser R., 2005, MNRAS, 361, L15Muterspaugh M. W. et al., 2010, AJ, 140, 1657Nakajima, T., Oppenheimer, B. R., Kulkarni, S. R., et al. 1995, Nature, 378, 463Nelson A. F., 2000, ApJ, 537, L65Nordstrom B. et al., 2008, VizieR Online Data Catalog, 5117, 0Ortiz, M., Reffert, S., Trifonov, T., et al. 2016, arXiv:1608.00963Oscoz, A., Rebolo, R., L´opez, R., et al. 2008, Proc. SPIE, 7014, 701447Paegert, M., Stassun, K. G., De Lee, N., et al. 2015, AJ, 149, 186Perryman M. A. C., 1997, in ESA Special Publication, Vol. 402, Hipparcos - Venice ’97,Bonnet R. M., Høg E., Bernacca P. L., Emiliani L., Blaauw A., Turon C., KovalevskyJ., Lindegren L., Hassan H., Bouffard M., Strim B., Heger D., Perryman M. A. C.,Woltjer L., eds., pp. 1–4Plez B., Cohen J. G., 2005, A&A, 434, 1117Queloz D. et al., 2001, A&A, 379, 279Queloz D. et al., 2000, A&A, 354, 99Rajpaul V., Aigrain S., Roberts S., 2016, MNRAS, 456, L6Ramsey L. W. et al., 1998, in Society of Photo-Optical Instrumentation Engineers (SPIE)Conference Series, Vol. 3352, Advanced Technology Optical/IR Telescopes VI, SteppL. M., ed., pp. 34–42Rebolo, R., Zapatero Osorio, M. R., & Mart´ın, E. L. 1995, Nature, 377, 129Rivera E. J. et al., 2005, ApJ, 634, 625Robinson S. E., Strader J., Ammons S. M., Laughlin G., Fischer D., 2006, ApJ, 637, 1102Santos N. C. et al., 2004, A&A, 426, L19Santos N. C. et al., 2002, A&A, 392, 215 36 –Shetrone M. et al., 2007, PASP, 119, 556Spiegel D. S., Burrows A., Milsom J. A., 2011, ApJ, 727, 57Thebault P., 2011, Celestial Mechanics and Dynamical Astronomy, 111, 29Thebault P., Haghighipour N., 2014, ArXiv e-printsThomas, N., Ge, J., Grieves, N., Li, R., & Sithajan, S. 2016, PASP, 128, 045003Tody, D. 1993, Astronomical Data Analysis Software and Systems II, 52, 173Torres G., Andersen J., Gim´enez A., 2010, A&A Rev., 18, 67Tull R. G., 1998, in Society of Photo-Optical Instrumentation Engineers (SPIE) ConferenceSeries, Vol. 3355, Optical Astronomical Instrumentation, D’Odorico S., ed., pp. 387–398Udry S. et al., 2006, A&A, 447, 361van Eyken J. C., Ge J., Mahadevan S., 2010, ApJS, 189, 156van Eyken J. C., Ge J., Mahadevan S., DeWitt C., 2004, ApJ, 600, L79Wan X., Ge J., Guo P., Zhao B., Hariharan L., McDavitt D., 2006, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 6269, Society ofPhoto-Optical Instrumentation Engineers (SPIE) Conference Series, p. 2Wang J., Ge J., Wan X., De Lee N., Lee B., 2012a, PASP, 124, 1159Wang J., Ge J., Wan X., Lee B., De Lee N., 2012b, PASP, 124, 598Wolszczan, A., & Frail, D. A. 1992, Nature, 355, 145Wright J. T., 2005, PASP, 117, 657Wu, Y., & Murray, N. 2003, ApJ, 589, 605Xie J.-W., Zhou J.-L., 2008, ApJ, 686, 570Xie J.-W., Zhou J.-L., 2009, ApJ, 698, 2066Xie J.-W., Zhou J.-L., Ge J., 2010, ApJ, 708, 1566 37 –Zhao B., Ge J., 2006, in Society of Photo-Optical Instrumentation Engineers (SPIE) Con-ference Series, Vol. 6269, Society of Photo-Optical Instrumentation Engineers (SPIE)Conference Series, p. 2Zucker S., Mazeh T., Santos N. C., Udry S., Mayor M., 2003, A&A, 404, 775Zucker S., Mazeh T., Santos N. C., Udry S., Mayor M., 2004, A&A, 426, 695