© ESO, 2011

1. Introduction

The characterization of the over 500 detected exoplanets has now begun to take place. Most of the studies are carried out at optical and infrared wavelengths, because this is where the planetary reflected light and thermal emission peak, respectively. The discovery of transiting planets (Charbonneau et al. 2000; Henry et al. 2000) allowed astronomers to constrain new physical parameters such as the radii and masses of the planets, which are not measurable by the radial velocity method alone. It is on these systems that in the last years the planetary atmosphere characterization has achieved the most exciting progress through the use of spectroscopy and broadband photometry with space telescopes. Examples are the identification of molecules such as water absorption (e.g. Tinetti et al. 2007) or methane (Swain et al. 2008), or the observation of the thermal emission variation with orbital phase (Knutson et al. 2007).

Although great improvements in characterizing the composition of transiting Hot-Jupiters have been achieved, they only represent about 20% of the known extrasolar planets¹. The characterization of non transiting planets would require the direct detection of their light, but the very low flux ratios between the planets and their host stars makes a direct detection a very challenging goal. Secondary eclipse observations from Spitzer show that planet-to-star flux ratios can be as high as 2.5    ×    10^-3 between 3.6 and 24 μm (e.g. Knutson et al. 2008). At 2.14 μm the expected flux should be less than these values. Many authors have attempted a direct detection of the Doppler-shifted signature in high-resolution spectroscopy from ground-based telescopes. In the optical Cameron et al. (1999) tried to observe the starlight reflected from the giant exoplanet Tau-Boötis b, they found an upper limit to the albedo and radius using a least-squares deconvolution method that is well described in the appendices of Collier Cameron et al. (2002), later the author repeated the analysis on υ Andromeda b (Collier Cameron et al. 2002). Recently Rodler et al. (2008, 2010) searched in the visible spectra of HD 75289Ab and Tau-Boötis b and found upper limits for their albedos using a model synthesis method. They constructed a model of the observation composed by a stellar template plus a shifted and scaled-down version of the stellar template to simulate the starlight reflected from the planet, these models were compared to the data by means of χ². In the near-infrared, several attempts have been made to detect Hot-Jupiters by trying to distinguish the planetary thermal emission from the starlight (Wiedemann et al. 2001; Lucas & Roche 2002; Barnes et al. 2007, 2008, 2010), they also found upper limits for the emitted flux of the planets. All these authors have used their own variation of a method based on the same principle of separating the planetary and stellar spectra given their relative Doppler shifts. Only recently, Snellen et al. (2010) claimed the detection of carbon monoxide from the transmission spectrum of HD 209458 b during a transit observation by using high-resolution spectra; nonetheless, his technique required a transiting system.

In this work we present an effort to constrain new physical parameters of the non-transiting Hot-Jupiter HD 217107 b. We attempt to trace its Doppler-shifted signature (estimated to be  ~10^-4 times dimmer than the star flux) with a correlation function between high-resolution data and models of its atmospheric spectrum. With positive detections this method would provide new information on its characteristics, such as its temperature, chemical composition, and the presence of chemical tracers associated with life. At the same time, the method enables the calibration of high-resolution spectroscopic models for a larger sample of planets that do not necessarily transit their parent star.

In Sect. 2 we review the planetary system HD 217107; in Sect. 3 we describe the observations, data reduction, and calibration procedures; in Sect. 4 we detail the theoretical atmospheric spectrum of the planet and the method used to extract and analyze the planetary signal and present the results of our data; in Sect. 5 we develop a strategy for the ideal data acquisition situation and simulate observations of other planetary systems; and in Sect. 6 we give the conclusions of our work.

2. The planetary system HD 217107

2.1. HD 217107 b discovery

HD 217107 is a main-sequence star that is similar to the Sun in mass, radius, and effective temperature; its spectral type, G8 IV, indicates that it is starting to evolve into the red-giant phase (Wittenmyer et al. 2007). The presence of HD 217107 b was first reported by Fischer et al. (1999) through radial velocity measurements of the star, the detection was then confirmed by Naef et al. (2001). Later, Fischer et al. (2001) identified a trend in the residuals of the fit, and Vogt et al. (2005) postulated the existence of a third companion in an external orbit with a period of 8.6 ± 2.7 yr. The presence of this third object promoted the study of this system in subsequent surveys (Butler et al. 2006; Wittenmyer et al. 2007; Wright et al. 2009), constraining more precisely the companions’ orbital parameters. Table 1 summarizes the parameters used in this work.

2.2. Radial velocity

The radial velocity of the planet, v_psini, around the center of mass of the system is given by the reflex motion of the star: ${\mathit{v}}_{\mathrm{p}}\mathrm{\left(}\mathit{t}\mathrm{\right)}\sin \mathit{i}\mathrm{=}\mathrm{-}{\mathit{v}}_{\mathrm{s}}\mathrm{\left(}\mathit{t}\mathrm{\right)}\sin \mathit{i}\frac{{\mathit{m}}_{\mathrm{s}}}{{\mathit{m}}_{\mathrm{p}}\sin \mathit{i}}\mathrm{\times }\sin \mathit{i}\mathit{.}$(1)It depends on the mass of the star, m_s; the minimum mass of the planet, m_psini; the projected radial velocity curve of the host star, v_s(t)sini (which in turn depends on the parameters T_p,  P,  e,  ω,  and  K_s); and the inclination of the orbit, i, and also on the velocity of the center of mass of the system, v_g, when measured from Earth. Thus, the radial velocity curve of the planet is a distinctive curve in time, parameterized by the values summarized in Table 1, where the only unknown parameter is the inclination of the orbit. Figure 1 shows the radial velocity curve of the star owing to the interaction with HD 217107 b, phased over one orbit, with the origin in phase (φ = 0) at the time of periastron. The radial velocity of the planet is proportional to this radial velocity curve (Eq. (1)).

Fig. 1 Radial velocity curve of HD 217107 vs. orbital phase. The crosses mark the observations of Wittenmyer et al. (2007), which we used to compute this orbital solution. The boxes over the curve indicate the coverage of our observations, the filled boxes represent the runs utilized in the analysis, while the open boxes represent the discarded runs (details in Sect. 4.2).

2.3. Flux estimate

By simulating the spectra of the planet and its host star as black bodies, we can estimate the order of magnitude of the planet-to-star flux ratio as a function of wavelength. The black body emission, F_λ(T), is determined by the surface temperature of the object. While for the star the temperature is well known from models (see Table 1), for the planet our best approximation is the equilibrium temperature ${\mathit{T}}_{\mathrm{eq}}\mathrm{=}{\left(\frac{\mathrm{1}\mathrm{-}\mathit{A}}{\mathrm{4}}\right)}^{\mathrm{1}\mathit{/}\mathrm{4}}{\left(\frac{{\mathit{R}}_{\mathrm{s}}}{\mathit{a}}\right)}^{\mathrm{1}\mathit{/}\mathrm{2}}{\mathit{T}}_{\mathrm{eff}}\mathit{.}$(2)For a reference value of the bond albedo of A = 0, we found an equilibrium temperature for HD 217107 b of T_eq = 1040 ± 19 K. Figure 2 shows the black body spectrum of the star and the planet assuming a radius between one and two Jupiter radii, which is the range of the radii for giant extrasolar planets measured to date.

Fig. 2 Black body emission of HD 217107 and HD 217107 b assuming a planet radius of 1.0 and 2.0 Jupiter radii. The vertical dashed line marks the waveband of our data (2.14 μm). The directly reflected light component has little contribution in the infrared and is thus omitted.

The planet-to-star flux ratio is given by $\mathrm{Flux}\mathrm{ratio}\mathrm{=}\frac{{\mathit{F}}_{\mathit{\lambda }}\mathrm{\left(}{\mathit{T}}_{\mathrm{planet}}\mathrm{\right)}}{{\mathit{F}}_{\mathit{\lambda }}\mathrm{\left(}{\mathit{T}}_{\mathrm{star}}\mathrm{\right)}}\mathrm{=}\frac{{\mathit{B}}_{\mathit{\lambda }}\mathrm{(}\mathit{T}\mathrm{=}\mathrm{1040}\mathrm{K}\mathrm{)}}{{\mathit{B}}_{\mathit{\lambda }}\mathrm{(}\mathit{T}\mathrm{=}\mathrm{5646}\mathrm{K}\mathrm{)}}{\left(\frac{{\mathit{R}}_{\mathrm{p}}}{{\mathit{R}}_{\mathrm{s}}}\right)}^{\mathrm{2}}\mathrm{·}$(3)At 2.14 μm, the flux ratio varies between 3    ×    10^-5 and 1.5    ×    10^-4 from one to two Jupiter radii of the planet’s radii, respectively. For shorter wavelengths the flux ratio decreases, because the star light dominates the emission spectrum. For longer wavelengths the net fluxes and thus the signal to noise ratio are lower.

3. Observations and data reduction

3.1. Observations

We observed the planetary system HD 217107 in 11 nights between 2007 August 14 and November 28 using Phoenix (Hinkle et al. 2003), a high-resolution near-infrared spectrometer at the Gemini South Observatory.

The spectrograph has a 256    ×    1024 InSb Aladdin II array with a resolving power of μm per pixel, the slit covers 14 arcsec in length. Its gain is 9.2 e⁻/ADU and it has a readout noise of 40 electrons. An argon hollow cathode wavelength calibration source is supplied with the instrument. Over 950 frames of the system were obtained in service mode, using the standard ABBA nodding sequences to easily remove sky emission, they cover a portion in the infrared spectral range from 2.136 to 2.145 μm (see Table 2). We tuned the data acquisition after receiving the data from the first runs since the instrument was not fully characterized for use on the Gemini Telescope. For the first two nights, the exposure time was set to 25 s, whereas for the rest of the nights it was set to 80 s. We requested arc-lamp calibration exposures as well.

3.2. Reduction

We wrote our own interactive data language (IDL)² routines for the data reduction and analysis, processing each night and slit position as an independent data set to minimize systematics caused by different atmospheric conditions or instrumental set-up. We used the flat-field images to identify hot pixels, marking a pixel as bad if it had a value beyond 3.5 sigma from the median of the values of the nine subsequent pixels in its neighborhood. Bad pixels were masked in all further processing stages. Then, we divided the frames by a per-night master flat-field and subtracted their corresponding opposite A or B frame to remove bias and sky. Finally, we extracted the spectra from the frames with an IDL implementation³ of the optimal spectrum extraction algorithm described in Horne (1986), this algorithm identified cosmic ray hits, which were also masked from subsequent processing.

3.3. Wavelength calibration

First, we calibrated the wavelength dispersion using the ThAr lamps, identifying the line positions and strengths in a high-resolution ThAr line atlas (Hinkle et al. 2001). Because there was only one calibration lamp for each night, this solution represented only a rough wavelength calibration, because there are (sub pixel) offsets in wavelength in the data. To reach the high precision needed for this work, we fine-tuned the calibration with a high-resolution spectrum of the Sun⁴ to identify the telluric lines (identified as those present both in the solar spectrum and in an average spectrum of our data set).

We constructed an average spectrum to increase the S/N ratio by aligning and adding the spectra of each night. To determine the relative shifts, we selected within each set the first spectrum as reference, while the rest were shifted (using spline interpolations) to calculate the shift that minimized the root-mean-square of the correlation with the reference. The centers of fifteen common absorption lines were identified in wavelength values for the solar spectrum and in pixel position for our average spectrum. The wavelength solution is obtained by fitting a second-order polynomial (λ = c₀ + c₁p + c₂p²) to the solar wavelength vs. the pixel position. Typical fitting coefficients are c₀ = 2.145407, c₁ = −1.0305    ×    10^-5, and c₂ = −1.8496    ×    10^-10. The dispersion of the residuals is rms = 4.21    ×    10^-6μm. No pattern is seen in the residuals. Pixels at wavelengths dominated by the identified telluric absorption lines were discarded from subsequent processing owing to their highly variable nature. About 65% of the pixels remained for the next analysis steps.

4. Data analysis and results

4.1. Correlation

Because it is impossible to directly distinguish the planet’s signature from the stellar one in a single spectrum, following the idea of Deming et al. (2000) and Wiedemann et al. (2001), we searched for the planetary Doppler-shift signature through a correlation method between the (stellar-subtracted) residual data and a synthetic model of the planet’s spectrum. To remove the stellar flux, we aligned the spectra for each set (Doppler-shifting them and using a spline interpolation) in a reference system in which the star remains at rest, and constructed a stellar template from the average of the set. Then, the stellar templates and the spectra are normalized dividing by their respective medians. Finally, the wavelengths of the stellar templates are shifted according to the orbital phase of the star in each individual spectrum, and then the stellar template is subtracted from them. We avoided combining the different nights to obtain the stellar template, because it is highly probable that other systematics would be introduced.

Because the planet is approximately a thousand times less massive than its host star, the planetary Doppler wobble is greater by the same order of magnitude (see Eq. (1)), consequently the planetary signature will not be added coherently in the stellar template and thus appear blurred. The stellar template subtraction leaves a residual spectrum that consists of the signature of the planet, which is slightly attenuated in the averaging process and immersed in Poisson noise. The blurring of the planet signature (see Fig. 3) is determined by the planetary velocity span, which in turn depends on the time span of an observation and the orbital phase at the time of the observation. Observations near inferior or superior conjunction provide the greatest radial velocity spans, while observations close to the greater elongation of the planet’s orbit produce the smallest radial velocity spans, rendering the data useless. The rejected data sets in Table 2 were observed near greater elongation.

Fig. 3 Planetary spectrum blurring in the stellar template. Using our synthetic spectra of HD 217107 b we simulated the smearing of the planetary spectrum over one observing run. The dark-gray and light-gray lines denote the first and last spectra of HD 217107 b during an observing night. The relative shift owing to the orbital motion of HD 217107 b is 40 km s-1 in this simulation. The bottom black line shows averaged the planetary spectra.

For the high-resolution synthetic planetary spectra of HD 217107 b we used customized theoretical thermal emission models of its atmosphere (model described in Fortney et al. 2005, 2006, 2008) at three different distances from the star to account for the non-negligible eccentricity of the planet. The models are cloud-free, with solar metallicity, gravity g = 20 m s^-2, and the molecular abundances are those appropriate for chemical equilibrium. At these effective temperatures, the main absorbing molecules are H₂O, CH₄, CO, and CO₂. The chemistry is described in detail in Lodders & Fegley (2002) and Visscher et al. (2006). We empirically characterized the instrumental resolution through the analysis of the emission lines in the calibration lamps. We convolved the model spectra by the instrumental resolution, which we determined to be . Then, for an assumed value of sini, the synthetic spectra are Doppler-shifted to mimic the radial velocity of the planet at the time that the data frame was obtained. The correlation degree, C(i), is calculated according to the formula $\mathit{C}\mathrm{\left(}\mathit{i}\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}}{\mathit{N}}\sum _{\mathit{k}\mathrm{=}\mathrm{1}}^{\mathit{N}}\frac{{\sum }_{\mathit{j}\mathrm{=}\mathrm{1}}^{{\mathit{N}}_{\mathit{k}}}\mathrm{(}{\mathit{r}}_{\mathit{kj}}\mathrm{-}\mathit{r̅}{\mathit{k}}^{\mathrm{)}}\mathrm{\left(}{\mathit{\tau }}_{\mathit{kj}}\mathrm{\right(}\mathit{i}\mathrm{)}\mathrm{-}\mathit{\tau ̅}\mathit{k}\mathrm{(}\mathit{i}{\mathrm{\left)}}^{\mathrm{\right)}}}{\sqrt{\left\{{\sum }_{\mathit{j}\mathrm{=}\mathrm{1}}^{{\mathit{N}}_{\mathit{k}}}{\mathrm{(}{\mathit{r}}_{\mathit{kj}}\mathrm{-}\mathit{r̅}{\mathit{k}}^{\mathrm{)}}}^{\mathrm{2}}\right\}\left\{{\sum }_{\mathit{j}\mathrm{=}\mathrm{1}}^{{\mathit{N}}_{\mathit{k}}}{\mathrm{(}{\mathit{\tau }}_{\mathit{kj}}\mathrm{\left(}\mathit{i}\mathrm{\right)}\mathrm{-}\mathit{\tau ̅}\mathit{k}\mathrm{(}\mathit{i}{\mathrm{)}}^{\mathrm{)}}}^{\mathrm{2}}\right\}}}\mathrm{·}$(4)In this equation we used the notation f_kj for the value of the function at the pixel j of the spectrum k, and for the mean value of the function in the spectrum k. Here, “r” refers to the residual spectrum while “τ(i)” to the shifted planetary model spectrum, with N_k the number of pixels in spectrum k and N the total number of spectra. The denominator in the expression normalizes the correlation, and thus a value of 1.0 would indicate a perfect correlation.

We thus produce a curve of the correlation degree vs. the inclination of the orbit, evaluated in the range 0 < i < π/2. A positive value of this function indicates that the data spectrum resembles that of the model, while a negative one suggests anti-correlation. As consequence of the random nature of the Poisson noise, the value of the correlation between the residual spectra and the models should be close to zero, except when the adopted i matches that of the planetary system. Therefore, an appreciable peak in the correlation curve would represent a successful detection of the planetary signature and immediately indicates the value of i. By constraining the inclination with this method, the mass of the planet would be immediately determined via Eq. (1).

4.2. Data results

For this analysis, we excluded the nights where the velocity span of the star was less than 10 m s^-1 (Col. 4 of Table 2) since they do not represent any significant improvement in the results, because the shift of the planet (~8.3 km s^-1) is not significantly higher than the instrumental resolution (~7.5 km s^-1). Figure 4 (Top panel) shows the correlation curve derived from our data as a function of sini. The degree of correlation found was close to zero at all inclinations, and we do not distinguish any identifiable positive peak that could indicate an atmosphere with absorption features resembling those of the models.

4.3. Planet-to-star flux ratio fitting

While the inclination determines the maximum of the correlation curve, the planet-to-star flux ratio () is the main physical parameter bounded to the magnitude of the correlation. In this section we determine the most probable values in the parameter space [ , sini], which gave rise to our result, and estimate the statistical significance of the value of the correlation reached. We searched for the best fitting values comparing our data results (Fig. 4 Top) with “synthetic” correlation curves. We generated the synthetic correlation curves by recreating our observations, adding a synthetic planetary spectrum, with known inclination and planet-to-star flux ratio, according to the following scheme: Step 1: we rearranged the order of the data set with random permutations within each night, but kept the original order of the dates of the observations. As a consequence, any real planet signature disappeared, but the noise level of the data was conserved. Step 2: using the atmospheric models of the planet, we injected a synthetic spectrum in the scrambled data set, Doppler-shifted and with a relative flux according to specific values of sini and , respectively. For simplicity, we adopted a constant along the orbit. Step 3: we processed these synthetic data through the same routines as in our original data (Sect. 4.1). We then iterated for a grid of values in the ranges: 0 ≤ i ≤ π/2 and 10^-5 ≤ F_p/F_s ≤ 10^-2, obtaining a set of synthetic correlations for sini and .

Once we obtained these models, we searched for the best-fit parameters through a χ² minimization between the data correlation curve and the synthetic correlation curves, generating a goodness-of-fit map (Fig. 4, bottom panel).

Fig. 4 Top: correlation result for our data set as a function of sini. The correlation remains flat along every inclination without any distinctive peak, the maximum value is reached at sini = 0.71. Bottom: goodness of fit, χ2-map, of the correlation models to the data. The horizontal and vertical axes refer to the fitting parameters sini and Fp/Fs respectively, at which the synthetic planetary spectrum was added in the correlation models creation, from which we calculated the minimum squares (${\mathit{\chi }}_{\mathit{i,}\mathrm{fr}}^{\mathrm{2}}$). We plot χ2 relative to the best fit ($\mathrm{\Delta }{\mathit{\chi }}_{\mathit{i,}\mathrm{fr}}^{\mathrm{2}}\mathrm{=}{\mathit{\chi }}_{\mathit{i,}\mathrm{fr}}^{\mathrm{2}}\mathrm{-}{\mathit{\chi }}_{\min }^{\mathrm{2}}$) using the function exp(−αΔχ2). The gray scale denotes the goodness of fit, from black for the best fit (at ${\mathit{\chi }}_{\min }^{\mathrm{2}}$), to white for the poorest fit. The plotting parameter, α, just enables a good contrast in the plots (the same value was used for all plots). Additionally, we determined with bootstrap procedures the solid lines (bottom to top) that mark the (1, 2, 3 and 4-σ) levels of false-alarm probability. The white cross marks the best fit at sini = 0.84 and Fp/Fs = 3.6    ×    10-3, situated below the 3-σ confidence level.

In addition, we used a bootstrap procedure to calculate false-alarm-probability limits for this map. Following Collier Cameron et al. (2002), we determined the frequency with which the correlation degree exceeds a given value as a result of noise in the absence of a planet signal. The routine consists of performing a random permutation of the data sets and the subsequent data analysis (steps 1 and 3 of previous paragraphs) which we repeated with a large number of times (~5000), recording the correlation curve after each trial. This set of correlation curves represents the correlation found in the absence of a planetary signal, and, because it is created from the data themselves, defines an empirical probability distribution that includes both the photon statistics and instrumental systematics.

Then, at each inclination, we stacked and sorted the values of the correlation in increasing order. We determined the 1, 2, 3, and 4-σ false-alarm confidence levels as the value of the correlation degree at the 65, 90, 99 and 99.9 percentiles of the trials. They represent the signal strengths at which spurious detections occur with 35, 10, 1, and 0.1 percent false-alarm probability respectively, at each value of the inclination. This allows us to assess the probability of obtaining a certain correlation degree in the absence of planetary emission.

Figure 4 (Bottom panel) shows the probability map for HD 217107 b. The best fit occurs at sini = 0.838 and F_p/F_s = 3.6    ×    10^-3, although the relative improvement in χ² against the surrounding parameters is shallow. This value disagrees with the maximum value of the correlation curve (Fig. 4 Top), and furthermore, the bootstrap results indicate that this value is below the 3-σ confidence limit of the signal not being a false positive. Also, this is much higher than the predicted value from Sect. 2.3 (between 3    ×    10^-5 and 1.5    ×    10^-4). The disagreement of the results of the top and bottom panel in Fig. 4 suggests that systematics remain after the data reduction, while the strength of the correlation value, two orders of magnitude above the expected flux ratio, indicates that this result is not realistic. The 3-σ confidence limit only allows us to establish an upper limit in the flux ratio at 4–5     ×    10^-3 for inclinations greater than sini = 0.6. In conclusion, we cannot state the detection of HD 217107 b.

5. Future prospects

5.1. Observational strategy

Although our current data do not enable us to claim the detection of HD 217107 b, we identified a strategy to maximize the chances of a successful detection. This involves selecting suitable candidate systems and precisely choosing the phasing and span of the observations. To exemplify the advantages, we simulated realistic observations of other planetary systems.

We limited our sample to the currently known extrasolar planets without transits⁵, observables from Gemini South, for the first semester of the year. Even though we constrained ourselves to the instrument we characterized and to a fixed time span, the purpose of these simulations is to provide one successful detection with our method. It is plausible that considering the full extent of possibilities, stronger signals can be acquired. The improvement in a detection, limited by purely photon noise, can be quantified by the planet-to-star flux ratio and the stellar flux, according to the expression (fluxes in number of photons) ${\mathit{F}}_{\mathrm{p}}\mathit{/}\mathrm{Noise}\mathrm{=}\frac{\mathrm{\left(}{\mathit{F}}_{\mathrm{p}}\mathit{/}{\mathit{F}}_{\mathrm{s}}\mathrm{\right)}{\mathit{F}}_{\mathrm{s}}}{\sqrt{\mathrm{\left|}{\mathit{F}}_{\mathrm{s}}\mathrm{\right|}\mathrm{+}\mathrm{\left|}{\mathit{F}}_{\mathrm{p}}\mathrm{\right|}}}\mathrm{\approx }\mathrm{\left(}{\mathit{F}}_{\mathrm{p}}\mathit{/}{\mathit{F}}_{\mathrm{s}}\mathrm{\right)}{\sqrt{\mathit{F}}}_{\mathrm{s}}\mathit{.}$(5)Then, for example, the spectrometer CRIRES with four times the wavelength coverage of that of Phoenix, has twice the sensitivity of Phoenix. We decided to simulate the Phoenix instrument, since it is well characterized by our group, while other instruments should present their own systematics, which are hard to quantify.

We simulated the planets as if the strength of the high-resolution absorption features were the same as that of our models for HD 217107 b, but with the corresponding (estimated as in Sect. 2.3). A caveat for this assumption is that the strength of the lines is not very clear in planets that exhibit thermal inversions. Burrows et al. (2008) and Fortney et al. (2008) suggest that the emission features should be weaker.

Fig. 5 Correlation curves and χ2-maps of synthetic data of HD 179949. A synthetic planetary signal was injected in the spectra with parameters: sini = 0.77 and Fp/Fs = 0.003. Left: using our observing strategy. Right: random distributed observing dates. The routine successfully recovers the signal at sini = 0.78 and Fp/Fs = 2.8    ×    10-3 in both cases, although, when using our strategy, the correlation degree is stronger, and the parameters are better determined compared with the right panel.

The target selection criteria are based first on the radial velocity span of the planet, where we set a lower limit cutoff of 7.5 km s^-1 (equivalent to the FWHM of the instrument spectral resolution) for a three-hour observing run if the orbit was at i = 30°. Second, we look for higher apparent brightnesses of the stars for better signal-to-noise ratios. Table 3 lists two of the better suited selected planetary systems (HD 217107 listed for comparison). A brighter K-band magnitude of the star improves the signal-to-noise ratio, while a smaller semi-major axis involves a higher radial velocity span, which enables a greater Doppler shift of the planet spectra during the runs and at the same time favors higher planet-to-star flux ratios.

To simulate the observations, we recreated the same instrumental settings of our data, but carefully selected the observing schedule. For each one of the nights in the period and restricted to air masses under 1.5, we selected the three-hour range that gives the maximum velocity span. We recorded then, the radial velocity spans for each night, and chose those with the biggest spans. We used the solar spectrum to simulate the stellar spectrum, while for the planetary component we used the atmospheric models of HD 217107 b added with a given planet-to-star flux ratio and inclination. Each component is Doppler-shifted according to the orbital parameters. Finally we added Poisson noise to the spectra, according to the signal-to-noise corresponding to the magnitude of the target. The synthetic data were processed in exactly the same way as our original data.

5.2. Simulations

In our first test, we present two simulations of an observing campaign on a target with the physical parameters of HD 179949 to show the improvements of our observing strategy in contrast with a regular observation. The given parameters are sini = 0.77 and F_p/F_s = 3    ×    10^-3. Figure 5 left shows the simulation following our observing strategy, while Fig. 5 right shows the simulation selecting random observing dates. In both cases the correlation curves (top panels) mark the inclination of the synthetic orbit with a increment in the correlation degree near sini = 0.77, while the χ²-maps (bottom panels) effectively indicate the best fit at sini = 0.78 and F_p/F_s = 2.8    ×    10^-3.

We identify the main differences between these two simulations: first, given the larger radial velocity spans when implementing our strategy, the planetary spectrum is more blurred in the stellar template and consequently less reduced when the template is subtracted, the planetary spectrum signal is thus stronger in the residual spectrum, which increases the correlation degree. As consequence of these greater correlation degrees, all confidence levels are generally lower, because it is less probable to reach this correlation degree by chance in the no-planet case, and lastly, the χ²-map peak is much better determined. The improvement is reflected more in the distinction of the best fit against other values of the parameter space than in the distinction against the no-planet case.

In another simulation, we recreated the planetary system Tau Boo as close as possible to its real physical characteristics (Fig. 6). Tau Boo b had the parameters sini = 0.82 and F_p/F_s = 4    ×    10^-4, our routines returned the best fit: sini = 0.79 and F_p/F_s = 3.6    ×    10^-4, slightly underestimating the values. Nevertheless, the χ²-map shows an improvement in the region near the injected inclination and flux ratio. The bootstrap results set the 3–σ confidence limit near F_p/F_s 1.5     ×    10^-4 (for sini > 0.5), indicating a detection with 99% confidence.

6. Discussion and conclusions

Because the instrument was not well characterized at the time and our service-mode observational strategy had to be adapted after the first few observing windows, the data for HD 217107 were not as sensitive as expected. The correlation curve was featureless for all inclinations and with values close to zero, with a maximum at sini = 0.71. By fitting the sine of the inclination and the planet-to-star flux ratio through least-squares, we found the best-fit parameters of sini = 0.84 and F_p/F_s = 3.6    ×    10^-3 at a level below our 3-σ confident limit. As a consequence of the faint features in the results, the disagreement between the peak in the correlation (Fig. 4 Top panel) and the most probable value of sini (Fig. 4 Bottom panel), and the higher than predicted F_p/F_s, we could not claim a detection of HD 217107 b with our current data. Given the results of the bootstrap procedure, we reject the flux ratio of HD 217107 b to its host star to be over 5     ×    10^-3 (3-σ confidence). We attribute these results to the absence of an ideal strategy in the data acquisition at the time of the observations and a needed further treatment of the instrument systematics.

Fig. 6 Similar to Fig. 4, correlation curve (top panel) and χ2-map (bottom panel) of a simulation of the planetary system Tau Boo, with an injected companion at sini = 0.82 and Fp/Fs = 4    ×    10-4.

We could not detect HD 217107 b, but defined the outlines of future campaigns by carefully defining a candidate selection criterion and an observational strategy. We conclude that the best-suited candidates for this technique are those in very close orbits, which allow the planets to have high orbital velocities and higher planet-to-star flux ratios. We propose an observing strategy where for the period of observations we specifically select the nights with maximum radial velocity spans. To explore the capabilities of our routines, we simulated other planetary systems as observed by the Phoenix spectrograph, with the same number of hours and an appropriate schedule of observations. The system HD 179949 was recreated, contrasting the use of our observing strategy with a regular observation schedule. We recovered the planetary signature in both cases, but showing an improvement in the correlation degree, precision in the χ²-map, and lower σ limits when using our observing schedule. Finally, we performed a realistic simulation of the planetary system Tau Boo, and successfully detected its signature.

In conclusion high-resolution instruments like Phoenix are capable of detecting extrasolar planet Doppler-shifted signals with flux ratios as low as 10⁴ with this method if we perform a careful treatment of the systematics (approaching the photon noise limit), if we count with appropriate theoretical models, and if we follow an optimized scheme in the data acquisition. Furthermore, using other instruments like CRIRES or NIRSPEC could increase the confidence of the results. Since our simulations exclude systematics effects specific to the instrument, an adequate treatment to remove them would be necessary. Our optimal observing strategy tends to select observations at superior conjunctions of the planet’s orbit, capturing the highest amount of light possible from the planet and at the same time covering the highest radial velocity span for a determined time extent. Refinements of this technique will involve the optimization of the distribution of time designated to the length of an observing run vs. the number of nights of observation, while adding phase-dependent functions of the planet’s brightness to account for the changing observed portion of the day/night side of the planet, and for different amounts of irradiation in eccentric orbits, will increase the accurateness of the fitted parameters.

Acknowledgments

Patricio Cubillos and Patricio Rojo are supported by the FONDAP Center for Astrophysics 15010003, the center of excellence in Astrophysics and Associated Technologies (PFB06) and the FONDECYT project 11080271.

References

All Tables

All Figures

	Fig. 1 Radial velocity curve of HD 217107 vs. orbital phase. The crosses mark the observations of Wittenmyer et al. (2007), which we used to compute this orbital solution. The boxes over the curve indicate the coverage of our observations, the filled boxes represent the runs utilized in the analysis, while the open boxes represent the discarded runs (details in Sect. 4.2).
In the text

	Fig. 2 Black body emission of HD 217107 and HD 217107 b assuming a planet radius of 1.0 and 2.0 Jupiter radii. The vertical dashed line marks the waveband of our data (2.14 μm). The directly reflected light component has little contribution in the infrared and is thus omitted.
In the text

Fig. 3 Planetary spectrum blurring in the stellar template. Using our synthetic spectra of HD 217107 b we simulated the smearing of the planetary spectrum over one observing run. The dark-gray and light-gray lines denote the first and last spectra of HD 217107 b during an observing night. The relative shift owing to the orbital motion of HD 217107 b is 40 km s-1 in this simulation. The bottom black line shows averaged the planetary spectra.

In the text

Fig. 4 Top: correlation result for our data set as a function of sini. The correlation remains flat along every inclination without any distinctive peak, the maximum value is reached at sini = 0.71. Bottom: goodness of fit, χ2-map, of the correlation models to the data. The horizontal and vertical axes refer to the fitting parameters sini and Fp/Fs respectively, at which the synthetic planetary spectrum was added in the correlation models creation, from which we calculated the minimum squares (${\mathit{\chi }}_{\mathit{i,}\mathrm{fr}}^{\mathrm{2}}$). We plot χ2 relative to the best fit ($\mathrm{\Delta }{\mathit{\chi }}_{\mathit{i,}\mathrm{fr}}^{\mathrm{2}}\mathrm{=}{\mathit{\chi }}_{\mathit{i,}\mathrm{fr}}^{\mathrm{2}}\mathrm{-}{\mathit{\chi }}_{\min }^{\mathrm{2}}$) using the function exp(−αΔχ2). The gray scale denotes the goodness of fit, from black for the best fit (at ${\mathit{\chi }}_{\min }^{\mathrm{2}}$), to white for the poorest fit. The plotting parameter, α, just enables a good contrast in the plots (the same value was used for all plots). Additionally, we determined with bootstrap procedures the solid lines (bottom to top) that mark the (1, 2, 3 and 4-σ) levels of false-alarm probability. The white cross marks the best fit at sini = 0.84 and Fp/Fs = 3.6    ×    10-3, situated below the 3-σ confidence level.

In the text

Fig. 5 Correlation curves and χ2-maps of synthetic data of HD 179949. A synthetic planetary signal was injected in the spectra with parameters: sini = 0.77 and Fp/Fs = 0.003. Left: using our observing strategy. Right: random distributed observing dates. The routine successfully recovers the signal at sini = 0.78 and Fp/Fs = 2.8    ×    10-3 in both cases, although, when using our strategy, the correlation degree is stronger, and the parameters are better determined compared with the right panel.

In the text

	Fig. 6 Similar to Fig. 4, correlation curve (top panel) and χ2-map (bottom panel) of a simulation of the planetary system Tau Boo, with an injected companion at sini = 0.82 and Fp/Fs = 4    ×    10-4.
In the text