THE LEECH EXOPLANET IMAGING SURVEY: CHARACTERIZATION OF THE COLDEST DIRECTLY IMAGED EXOPLANET, GJ 504 b, AND EVIDENCE FOR SUPERSTELLAR METALLICITY^*

When brown dwarfs cool below ∼1200 K, their atmospheres transition from cloudy to clear, and methane becomes a dominant absorber in their spectral energy distributions (SEDs). Although gas giants were thought to be analogs to brown dwarfs, the first generation of directly imaged exoplanets had cloudy, methane-free atmospheres, even though their temperatures are well below the temperature where field brown dwarfs have had their "L $\to$ T" transition (Chauvin et al. 2004; Marois et al. 2008, 2010; Skemer et al. 2011). GJ 504 b, discovered by the SEEDS survey (Strategic Explorations of Exoplanets and Disks with Subaru; Tamura 2009), is the first example of an exoplanet that is cold enough (∼500 K) to be relatively cloud-free and have strong methane absorption features (Janson et al. 2013; Kuzuhara et al. 2013). Another T-dwarf exoplanet was recently discovered around 51 Eri (Macintosh et al. 2015).

At a separation of 25 (43.5 AU) from its G-star host, GJ 504 b is easily accessible to most high-contrast imaging systems (Kuzuhara et al. 2013). Its $H-{{K}}_{{s}}$ color (0.63 ± 0.15) is highly discrepant with similar-luminosity field brown dwarfs (approximately −0.2; Dupuy & Liu 2012). At the same time, it has strong methane absorption at 1.66 μm, clearly placing it in a different class than other directly imaged exoplanets (Janson et al. 2013). GJ 504 also has supersolar metallicity ([M/H] = 0.10–0.28, although most determinations are toward the lower end of this range; Edvardsson et al. 1993; Mishenina et al. 2004; Valenti & Fischer 2005; Takeda 2007; Gonzalez et al. 2010; Maldonado et al. 2012; Ramírez et al. 2013). It is therefore relatively likely to have a gas giant planet (Fischer & Valenti 2005).

The age of GJ 504 A, and thus the mass and planetary status of GJ 504 b, is uncertain. Kuzuhara et al. (2013) find consistency among multiple age indicators, such as X-ray activity, rotation rate, chromospheric activity, and H-R diagram location, indicating an age for GJ 504 A of 0.1–0.5 Gyr. However, a reanalysis of GJ 504 A's stellar properties by Fuhrmann & Chini (2015) suggests that the star lies above the main sequence on an H-R diagram with a corresponding age of  4.5 Gyr. Fuhrmann & Chini (2015) argue that the rapid rotation and other signs of youth arise because a massive planet has fallen into the star, carrying its orbital angular momentum with it. This leaves the presence of strong lithium absorption (Kuzuhara et al. 2013) unexplained (Soderblom 2010). Some of the lithium could have been replenished by the planet (Carlberg et al. 2012), but an usually massive planet would be required. Since there is no consensus on the age of the system, we consider both estimates in the following discussion; the younger age range implies a planet mass of ∼3–9 M_Jup, while the older one would suggest a mass of ∼30 M_Jup.

The LBTI Exozodi Exoplanet Common Hunt (LEECH) is a ∼100-night survey with the Large Binocular Telescope (LBT) to search for and characterize exoplanets in the mid-infrared (Skemer et al. 2014a; Maire et al. 2015). In this work, we confirm the detection of GJ 504 b and obtain photometry of the first T-dwarf exoplanet in three narrow L-band filters (3.71, 3.88, and 4.00 μm). For T dwarfs, L-band photometry can probe the broad methane fundamental absorption feature (centered at 3.3 μm) and put strong constraints on the luminosity of the planet, which peaks at ∼4 μm. The overall benefit of this additional photometry is to improve our ability to constrain GJ 504 b's bulk properties with atmospheric modeling. In Section 2, we present our observations and reductions, which comprise some of the deepest images taken from the ground at these wavelengths. In Section 3, we present our new photometry and adjust the literature photometry onto a common photometric system. In Section 4, we fit the photometry with a grid of models and discuss the physical nature of GJ 504 b. Finally, we present our conclusions in Section 5.

We observed GJ 504 on UT 2013 April 21 and UT 2014 March 11–13 with the Large Binocular Telescope Interferometer (LBTI; Hinz et al. 2012) and its 1–5 μm imager, L/M Infrared Camera (LMIRcam; Skrutskie et al. 2010; Leisenring et al. 2012). The LBT has twin deformable secondary adaptive optics (AO) systems (Esposito et al. 2011; Bailey et al. 2014), which provide excellent sensitivity in the thermal infrared (≳2 μm) compared to traditional AO systems (Lloyd-Hart 2000). The diffraction-limited beams from the AO systems are fed into LBTI, which can overlap or separate the two images on LMIRcam. For contrast-limited observations, such as the LEECH planet search (Skemer et al. 2014a), we typically separate the beams to allow independent and redundant observations of the inner, speckle-noise-limited regime. For sensitivity-limited observations, we overlap the beams (incoherently) so that the faint astronomical source can be extracted from a single bright sky background aperture, rather than two. GJ 504 b, at a separation of ∼25, falls into the latter category, where sensitivity is a greater priority than contrast. For the UT 2013 April 21 observations, the LBT's right-side adaptive optics system was unavailable, so we acquired data using just the left side of the telescope. For the UT 2014 March 11–13 observations, both sides of the telescope were operable, and we overlapped the two images of GJ 504 b to increase our sensitivity.

We observed GJ 504 in narrowband filters: L_NB6 (3.61–3.80 μm), L_NB7 (3.76–3.99 μm), and L_NB8 (3.97–4.03 μm). Basic properties for these filters are tabulated in Table 1, and transmission profiles are shown in Figure 1. The weather was photometric on the nights we obtained L_NB6 (UT 2014 March 12) and L_NB7 (UT 2013 April 21) data. The first night of L_NB8 observations (UT 2014 March 11) was nonphotometric, which prompted us to repeat this filter on UT 2014 March 13. UT 2014 March 13 had patchy clouds away from the telescope, which cleared early in the observations. Integration times were chosen to be long enough that the off-star data were sky background noise limited, rather than read noise limited. This choice meant that the star was saturated in the frames used to detect GJ 504 b. Additional frames with shorter integration times were used to measure the brightness of the host star. A summary of our observations is presented in Table 2.

Zoom In Zoom Out Reset image size

Table 2.  Observations

Date	Filter	Aperture	Frame Time (Sat/Unsat)	Int Time	Conditions
			(s)	(minutes)
UT 2013 Apr 21	LNB7	8.4 m	0.524/0.058	55	photometric, 09 seeing
UT 2014 Mar 11	LNB8	2 × 8.4 m	0.990/0.087	30	patchy clouds, 14 seeing
UT 2014 Mar 12	LNB6	2 × 8.4 m	0.291/0.029	101	photometric, 09 seeing
UT 2014 Mar 13	LNB8	2 × 8.4 m	0.873/0.087	44	patchy clouds and then clear, 10 seeing

Download table as:  ASCIITypeset image

We reduced the data using a custom LMIRcam pipeline developed at MPIA (Bonnefoy et al. 2014). The pipeline (1) replaces bad detector pixels with the median of their adjacent neighbors, based on a table of outlier pixels cataloged from off-sky calibration frames; (2) removes the detector bias and background sky/telescope emission by subtracting the median of images from chronologically neighboring nod subsequences; (3) determines the subpixel centroid of the star point-spread function (PSF) in each image by fitting a Gaussian, masking the inner saturated pixels, and shifts and crops the image to a common, aligned frame; (4) flags images with peak star fluxes below a specified threshold and images with abnormal background levels to exclude data contaminated by clouds and poor seeing; (5) removes residual detector bias from columns and rows, based on overscan pixels; and (6) forms the cube of reduced, photometric-quality images, for inspection and PSF subtraction. Angular differential imaging (ADI) and principal component analysis (PCA) PSF subtraction is carried out with a separate program described later in this section.

Figure 2 shows the result of aligning, de-rotating, and co-adding all of the photometric quality images acquired of GJ 504 b in the L_NB7 (3.95 μm) filter—a composite with an effective integration time of 3121 s. Owing to the relatively wide angular separation of GJ 504 b, the planet is visible without subtracting the star, in the northwest corner of the co-added image. The image flux scale is given in units of the star's peak intensity, determined from the short-exposure (unsaturated PSF) image sequence.

Zoom In Zoom Out Reset image size

To reduce starlight contamination in the planet signal, we subtract an estimate of the star's PSF from each image using standard high-contrast imaging techniques. The data were taken in ADI mode (Marois et al. 2006), where the instrument does not rotate (LBTI is on a nonrotating mount), so that its diffractive pattern and static aberrations stay fixed, while the sky image rotates with parallactic angle. We then reduced the data with a PCA high-contrast algorithm (Amara & Quanz 2012; Soummer et al. 2012; Fergus et al. 2014), using a custom implementation that closely follows Soummer et al. (2012). At each separation, we optimize the number of principle components that are used in the subtraction (N_PC) by inserting 12 artificial planets at different azimuth angles and measuring their signal-to-noise ratios as a function of N_PC. We also optimize the minimum parallactic rotation gap (parameterized as N_FHWM) between the image being fit and the library of images used to do the fit (Lafrenière et al. 2007). Optimal values of N_PC vary from ∼15 to 30, and optimal values of N_FHWM vary from ∼0.5 to 1.0. Final reduced images are shown in Figure 3. We detect GJ 504 b at signal-to-noise ratios of 5.9, 5.3, and 5.9 in L_NB6, L_NB7, and L_NB8, respectively (as described below).

Zoom In Zoom Out Reset image size

Relative photometry between GJ 504 A and b is measured using the forward modeling approach described in Soummer et al. (2012), which uses an image of the star at the position of the planet to calibrate the planet self-subtraction that is common to high-contrast image processing techniques (Lafrenière et al. 2007). We optimize the source model over an annular sector centered on the planet and spanning 30° in position angle, conservatively fixing N_FHWM = 1.0, as we find that values less than this increase self-subtraction. We measure the error on our relative photometry by inserting artificial planets at the same separation as GJ 504 b, but at different position angles, and repeating the forward modeling. The standard deviation of our measurements of the artificial planets is adopted as our formal relative photometry uncertainty and our detection signal-to-noise ratio. Relative photometry between GJ 504 A and b is reported in Table 3. Note that the two nights of L_NB8 observations were combined in this analysis, after confirming that they produced similar photometry, within 1σ errors.

To convert contrast measurements (relative photometry) to apparent magnitudes, we adopt and calculate photometry for GJ 504 A. For J, H, and K_s, we adopt apparent photometry from Kidger & Martín-Luis (2003), converted to Two Micron All Sky Survey photometry (J = 4.13 ± 0.01, H = 3.88 ± 0.01, K_s = 3.81 ± 0.01) using the methodology and Vega spectrum of Rieke et al. (2008). For all other filters where GJ 504 b has been observed, we calculate GJ 504 A's apparent photometry by fitting a model stellar atmosphere (Castelli & Kurucz 2004) with parameters that best match GJ 504's temperature, gravity, and metallicity measurements (Valenti & Fischer 2005). In our fit, we use the previously quoted JHK_s photometry, as well as WISE (Cutri et al. 2013) W3 (3.831 ± 0.015) and W4 (3.757 ± 0.022) photometry (the WISE W1 and W2 filters were not used because they were flagged as saturated). Our fit to these five data points with a model atmosphere produced a reduced χ² of 0.70. We estimate GJ 504 A's apparent magnitude to be 3.87 in CH₄s, 3.85 in CH₄l, and 3.82 in L', L_NB6, L_NB7, and L_NB8.²⁰ This model-driven approach leads to a 1σ–2σ inconsistency between our CH₄ photometry and our H-band photometry (the two CH₄ bands span H band and should not be brighter than the H photometry), which we would like to correct in order to avoid propagating erroneous color information into GJ 504 b's photometry (incorrect colors could affect our derived model atmosphere parameters, whereas an incorrect overall luminosity will only affect radius). Therefore, we adjust our estimated CH₄ photometry to be fainter by 0.02 mag to be consistent with the broad H-band value. Similarly, we adjust the L photometry to be brighter by 0.02 mag to be consistent with the K_s photometry. Since inconsistencies in our photometric estimates appear to be at the ∼0.02 mag level, we adopt 0.02 mag uncertainties for all GJ 504 A photometry. These inconsistencies are independent of our exact choice of stellar parameters and are more likely the result of propagated photometric errors or stellar variability. Our adopted GJ 504 A photometry and the resulting GJ 504 b photometry are summarized in Table 3. For our atmosphere modeling, we convert to absolute magnitudes assuming a distance of 17.56 ± 0.08 pc (van Leeuwen 2007).

With the first discoveries of directly imaged planets, it was obvious that something about planetary atmospheres made them different from the atmospheres of similar-temperature field brown dwarfs (Chauvin et al. 2004; Marois et al. 2008). In near-infrared color–magnitude diagrams (see Figure 4), the HR 8799 planets fall on what appears to be an extension of the L-dwarf sequence. The implication is that these young planets have retained their dusty, methane-poor atmospheres at lower luminosities than old field brown dwarfs, which transition to methane-rich, cloud-free T dwarfs below M_H ∼ 14 mag. Surveys have now found young, dusty brown dwarfs that are more analogous to the HR 8799 planets (e.g., Faherty et al. 2013; Liu et al. 2013; Gauza et al. 2015), and there are theoretical justifications for why, in addition to effective temperature, an object's mass/surface gravity affects the clouds and chemistry of its photosphere (Marley et al. 2012; Zahnle & Marley 2014).

Zoom In Zoom Out Reset image size

Also plotted in Figure 4 is GJ 504 b, which is much less luminous than previously discovered exoplanets. In the $J-H$ versus J and ${{K}}_{{s}}-L^{\prime}$ versus K_s diagrams, GJ 504 b falls right on the late T-dwarf sequence. However, in the $H-{{K}}_{{s}}$ versus H diagram, GJ 504 b is much redder than the field T dwarfs. There is one other object near GJ 504 b in the color–magnitude diagrams: GJ 758 B, a ∼30–40 M_Jup companion to a ∼6 Gyr Sun-like star (Mamajek & Hillenbrand 2008; Thalmann et al. 2009; Janson et al. 2011).

Clearly, there is something about GJ 504 b and GJ 758 B that make them different from other objects with similar luminosities. For the HR 8799 planets, youth and low surface gravity are responsible for their unusual appearances (Marois et al. 2008). While GJ 504 b's age is uncertain and GJ 758 B is clearly old, their position in Figure 4 unambiguously demonstrates that at least one of their physical properties is unusual.

For GJ 504 b, Kuzuhara et al. (2013) and Janson et al. (2013) suggest that gravity or metallicity could be driving the planet's unusual near-infrared colors. In the rest of this section, we attempt to model GJ 504 b's atmosphere and directly constrain these properties.

4.1. Atmosphere Models

We attempt to fit the complete SED of GJ 504 b using the photometry described in Section 3. We use models similar to those described in Morley et al. (2012, 2014), which include opacities for T/Y-dwarf condensates. The methane line lists have been updated using Yurchenko et al. (2014), and the alkali line lists have been updated to use the results from Allard et al. (2005). Chemical equilibrium grids based on previous thermochemical models (Lodders 1999, 2002; Lodders & Fegley 2002, 2006; Visscher et al. 2006, 2010; Visscher 2012; Moses et al. 2013) have been revised and extended to include higher metallicities. These updates will be described in detail in a set of upcoming papers that focus on the new model grid. In addition to temperature, our model grid parameterizes surface gravity, metallicity, and cloud thickness.²¹ We allow radius to be a free parameter so that atmospheric properties rather than luminosity drive the fit. Our model grid contains 480 models, comprising temperatures of 450, 475, 500, 525, 550, 575, 600, and 625 K, surface gravities of 30, 100, 300, and 1000 m s⁻², metallicities of [M/H] = 0, 0.5, and 1.0, and cloud thicknesses of f_sed = 1, 2, 3, 5, and cloud-free. The parameters that we choose to vary are among the most fundamental to the bulk appearance of planetary atmospheres. However, we cannot rule out that additional parameters, such as nonequilibrium NH₃ chemistry (Zahnle & Marley 2014), might play an important role as well. We also note that systematic differences between model families can account for substantially disparate parameter estimates (Patience et al 2012). With these caveats in mind, our best-fit model is T = 550 K, [M/H] = 0.5, R = 0.94 R_Jup, g = 100 m s⁻², f_sed = 3, and log(L) = −6.13 L_⊙ (see Figure 5). The reduced χ² (counting only radius scaling as a free parameter) is 0.98 with 7 degrees of freedom. None of the other models provide a fit with a relative probability higher than 1%, based on the Bayesian probability: $P({\mathrm{model}}_{1}/{\mathrm{model}}_{2})={e}^{({\chi }_{{\mathrm{model}}_{2}}^{2}-{\chi }_{{\mathrm{model}}_{1}}^{2})/2}$ for Gaussian photometric errors and a uniform model prior. While this analysis demonstrates that there is at least one self-consistent model that adequately fits all of the data, our current grid is too sparse to sample the error distribution of each parameter. Without sampling the error distribution, our best-fit model may not be at the peak of the global probability distribution. Thus, in Section 4.2, we interpolate between the models to form a denser grid, which we use to adopt estimates for each parameter.

Zoom In Zoom Out Reset image size

To lend some intuition to the effects of varying individual parameters, Figure 6 shows the best-fit model in four panels, which individually vary temperature, surface gravity, metallicity, and cloud properties. As with our best fit, all models are scaled by radius. In this scheme (which is partially driven by the size of the error bars and the radius fit), temperature is primarily constrained by the L'- and narrow L-band photometry. Gravity is constrained by the photometry in the J and L bands. Metallicity is constrained by J, K_s, and L. Cloud properties are also constrained by J, K_s, and L. However, for metallicity, K_s and L move in opposite directions, while for cloud properties, they move in the same direction. No parameter is fully degenerate with a combination of other parameters, and all four parameters (plus radius) were necessary to obtain an adequate best fit. The L-band photometry, in particular, was critical for resolving degeneracies between temperature and the other parameters.

Zoom In Zoom Out Reset image size

4.2. Interpolated Atmosphere Models

While our model grid is able to produce a plausible fit to the available GJ 504 b data, it is too sparsely spaced to sample the error distribution of the model parameters. We cannot easily produce a much larger grid of models, so instead, we interpolate between the models using quad-linear interpolation (linear with temperature, metallicity, and cloud properties, logarithmic with surface gravity). We assign probabilities to each model using the Bayesian posterior, $P\propto {e}^{-{\chi }^{2}/2}$, with Bayesian priors set by grid spacing. We adopt uniform priors for temperature, f_sed, and [M/H] over their full model-allowed range. For gravity, we adopt uniform priors above a minimum surface gravity that is set by the radius of the planet and the minimum planet mass (3 M_jup) derived by Kuzuhara et al. (2013). For radius we adopt a uniform prior between 0.9 and 1.3 R_Jup, the plausible radius range for GJ 504 b, varying mass, core mass, initial entropy, and metallicity (Fortney et al. 2007, 2008). Marginal probabilities are shown in Figure 7. Gaussian fits to the marginalized probability distributions give the following marginalized parameter distributions for GJ 504 b: T_eff = 544 ± 10 K, [M/H] = 0.60 ± 0.12, R = 0.96 ± 0.07 R_⊙, log(L) = −6.13 ± 0.03 L_⊙. For cloudiness and gravity, whose distributions do not resemble Gaussians, we adopt f_sed(cloudiness) = 2–5 and g < 600 m s⁻². Parallax uncertainty has a negligible effect on the radius and luminosity errors. Two-dimensional probability contours, for each pair of model parameters, are shown in Figure 8. Correlations are evident between parameters.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Interpolation carries an intrinsic risk that the spectra change in nonlinear ways, which can vary by bandpass. For our particular model grid and photometry, the best example of this behavior is seen in the metallicity plot of Figure 6. The J-band photometry changes more quickly between [M/H] = 0.5 and 1.0 than between [M/H] = 0.0 and 0.5. At the same time, many of the other model photometry points move linearly with [M/H]. To test how this might affect our final marginalized probability distributions, we repeated our analysis interpolating exponentially (10^[M/H]) in the J band and linearly in the other bands. The net result is a shift in the probability distribution of 3 K for T_eff, 0.02 for [M/H], 0.02 R_Jup for radius, and negligible changes for the other parameters. In all cases, this shift is much smaller than our derived error bars.

4.3. The Physical Properties of GJ 504 b

We obtained images of what is currently the coldest directly imaged exoplanet, GJ 504 b, in three narrow L-band filters, mapping out the 3.3 μm methane fundamental absorption feature, and putting further constraints on the planet's basic physical properties. With a best-fit temperature of 550 K, GJ 504 b is analogous to field brown dwarfs with a T spectral type. Indeed, the SED of GJ 504 b shows many similarities to a T-type brown dwarf: strong methane absorption features at 1.66 and 3.3 μm, blue $J-H$ colors that imply a relative lack of clouds compared to red L dwarfs, and an SED that is reasonably well fit by a water-vapor-dominated spectrum. On the other hand, GJ 504 b's SED is different from any currently known field brown dwarf, particularly with regard to its unusually red $H-{{K}}_{{s}}$ colors.

We constructed a model grid with radius, temperature, metallicity, surface gravity, and cloud types as free parameters to try to explain GJ 504 b's unusual SED. We find T_eff = 544 ± 10 K, f_sed (cloudiness) = 2–5, g < 600 m s⁻², [M/H] = 0.6 ± 0.12, R = 0.96 ± 0.07 R_Jup, and log(L) = −6.13 ± 0.03 L_⊙. If GJ 504 is young (0.1–0.5 Gyr), Baraffe et al. (2003) models predict that the companion is 3–8 M_jup. If GJ 504 is old (3–6.5 Gyr), Baraffe et al. (2003) models predict that the companion is 19–30 M_jup. Our estimate of the planet's surface gravity favors the low-mass interpretation, but not conclusively.

Of particular note, our planet atmosphere model requires a superstellar metallicity to explain GJ 504 b's complete SED, and particularly its $H-{{K}}_{{s}}$ colors shown in Figure 4. Various planet formation models predict that planet metallicities should differ from their host stars. Conversely, pairs of objects that form by binary fragmentation should have similar metallicities. Within the confines of our atmospheric models, this suggests that GJ 504 b formed like a planet, not like a binary companion. This result is independent of GJ 504 b's age and mass.

Additional photometry and spectroscopy at higher signal-to-noise ratio will further improve our ability to understand planets like GJ 504 b. In particular, observations over a broad wavelength range, including the mid-infrared (≳3 μm), can help break degeneracies in model parameters. Spectroscopy in the mid-infrared will be particularly valuable (Skemer et al. 2015). The James Webb Space Telescope (JWST) will discover and characterize a variety of new worlds (e.g., Beichman et al. 2010). The LEECH exoplanet imaging survey is searching for and characterizing cool exoplanets in the mid-infrared right now, with the goal of improving our theoretical understanding of their atmospheres so that we can take full advantage of JWST's limited life span.

The authors thank the anonymous referee for a very helpful report. We thank the HiCIAO team and Masayuki Kuzuhara for providing filter curves for the HiCIAO methane filters. A.S. is supported by the National Aeronautics and Space Administration through Hubble Fellowship grant HSTHF2-51349 awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS 5-26555. E.B. is supported by the Swiss National Science Foundation (SNSF). S.E., S.D., and A.L.M. acknowledge support from the "Progetti Premiali" funding scheme of the Italian Ministry of Education, University, and Research. C.V. is supported by NSF AST-1312305. LEECH is funded by the NASA Origins of Solar Systems Program, grant NNX13AJ17G. This material is based in part on work supported by the National Aeronautics and Space Administration under Agreement No. NNX15AD94G, Earths in Other Solar Systems, issued through the Science Mission Directorate interdivisional initiative Nexus for Exoplanet System Science. The Large Binocular Telescope Interferometer is funded by the National Aeronautics and Space Administration as part of its Exoplanet Exploration program. LMIRcam is funded by the National Science Foundation through grant NSF AST-0705296. This publication makes use of the Spex Prism Library Analysis Toolkit. This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, and NEOWISE, which is a project of the Jet Propulsion Laboratory/California Institute of Technology. WISE and NEOWISE are funded by the National Aeronautics and Space Administration.