# Subaru/CHARIS High-Resolution Mode Spectroscopy of the Brown Dwarf Companion
HD 33632 Ab

Aidan B. Gibbs Department of Physics & Astronomy, University of California,
Los Angeles, CA 90095, USA  Briley L. Lewis Department of Physics & Astronomy,
University of California, Los Angeles, CA 90095, USA  Michael P. Fitzgerald
Department of Physics & Astronomy, University of California, Los Angeles, CA
90095, USA  Timothy D. Brandt Department of Physics, University of California,
Santa Barbara, Santa Barbara, CA 93106, USA  Minghan Chen Department of
Physics, University of California, Santa Barbara, Santa Barbara, CA 93106, USA
Yiting Li Department of Astronomy, University of Michigan, Ann Arbor, MI
48109, USA  Rachel Bowens-Rubin Department of Astronomy, University of
California, Santa Cruz, Santa Cruz, CA 95064, USA  Rebecca Jensen-Clem
Department of Astronomy, University of California, Santa Cruz, Santa Cruz, CA
95064, USA  Benjamin A. Mazin Department of Physics, University of California,
Santa Barbara, Santa Barbara, CA 93106, USA

(Received May 20, 2024; Revised July 15, 2024; Accepted July 16, 2024)

###### Abstract

Brown dwarfs (BD) are model degenerate in age and mass. High-contrast imaging
and spectroscopy of BD companions to host stars where the mass and age can be
independently constrained by dynamics and stellar age indicators respectively
provide valuable tests of BD evolution models. In this paper, we present a new
epoch of Subaru/CHARIS H\- and K-band observations of one such previously
discovered system, HD 33632 Ab. We reanalyze the mass and orbit using our new
epoch of extracted relative astrometry, and fit extracted spectra to the
newest generation of equilibrium, disequilibrium, and cloudy spectral and
evolution models for BDs. No spectral model perfectly agrees with evolutionary
tracks and the derived mass and age, instead favoring a somewhat younger BD
than the host star’s inferred age. This tension can potentially be resolved
using atmosphere and evolution models that consider both clouds and
disequilibrium chemistry simultaneously, or by additional future spectra at
higher resolution or in other band passes. Photometric measurements alone
remain consistent with the luminosity predicted by evolutionary tracks. Our
work highlights the importance of considering complexities like clouds,
disequilibrium chemistry, and composition when comparing spectral models to
evolutionary tracks.

Brown dwarfs (185), Direct imaging (387), Infrared spectroscopy (2285),
Astrometry (80)

††journal: ApJ††facilities: Subaru Telescope, National Observatory of Japan
(NAOJ)††software: Numpy (Harris et al., 2020), Pandas(pandas development team,
2020), Scipy (Virtanen et al., 2020), Astropy (Astropy Collaboration et al.,
2022)

##  1 Introduction

Table 1: HD 33632 System Properties

HD 33632 A
---
Mass ($M_{\odot}$) | $1.1\pm 0.1$ | (Brandt et al., 2021a)
$\log_{10}(L/L_{\odot})$ | $-4.62^{+0.04}_{-0.08}$ | (Currie et al., 2020)
Distance (pc) | $26.39\pm 0.18$ | (Gaia Collaboration, 2020)
Age (Gyr) | $1.7\pm 0.4$ | (Brandt et al., 2021a)
[Fe/H] (dex) | $-0.22\pm 0.03$ | (Soubiran et al., 2016)
$J_{\text{2MASS}}$ | $5.430\pm 0.020$ | (Cutri et al., 2003)
$H_{\text{2MASS}}$ | $5.193\pm 0.015$ | (Cutri et al., 2003)
$K_{\text{S,2MASS}}$ | $5.166\pm 0.020$ | (Cutri et al., 2003)
HD 33632 B
Mass ($M_{\odot}$) | $0.215\pm 0.029$ | (Brandt et al., 2021a)
Separation (′′) | $33.99086\pm 0.00003$ | (Gaia Collaboration, 2020)
$J_{\text{2MASS}}$ | $10.38\pm 0.03$ | (Cutri et al., 2003)
$H_{\text{2MASS}}$ | $9.863\pm 0.019$ | (Cutri et al., 2003)
HD 33632 Ab
| (Currie et al., 2020) | (Brandt et al., 2021a) | This Work
$J_{\text{MKO}}$ | $16.91\pm 0.11$ | … | …
$H_{\text{MKO}}$ | $16.00\pm 0.09$ | … | $15.96\pm 0.08$
$K_{\text{S,MKO}}$ | $15.37\pm 0.09$ | … | $15.54\pm 0.12$
$L_{\text{P,MKO}}$ | $13.67\pm 0.15$ | … | …
| Fitted Parameters
Mass ($M_{J}$) | $46.4^{+8.1}_{-7.5}$ | $50.0^{+5.5}_{-5.0}$ | $51.6^{+5.4}_{-4.8}$
a (au) | $21.2^{+4.1}_{-4.5}$ | $23.6^{+3.2}_{-4.5}$ | $23.7^{+3.1}_{-3.8}$
e | $0.18^{+0.20}_{-0.14}$ | $0.12^{+0.18}_{-0.09}$ | $0.143^{+0.13}_{-0.089}$
i (∘) | $39.4^{+8.0}_{-20}$ | $45.2^{+4.7}_{-11}$ | $42.2^{+5.5}_{-11.0}$
$\omega$ (∘) | $151^{+155}_{-131}$ | $0^{+86}_{-140}$ | $282^{+60}_{-252}$
$T_{0}$ (JD) | $2440900^{+4500}_{-13100}$ | $2468815^{+18000}_{-5800}$ | $2466542^{+3096}_{-5100}$
PA ascending (∘) | $38.2^{+7.2}_{-7.0}$ | $39.3^{+5.7}_{-6.5}$ | $38.0^{+5.3}_{-5.3}$
$\lambda_{\text{ref}}$ (∘) | $21^{+14}_{-10}$ | $202^{+14}_{-9.5}$ | $208^{+10}_{-8}$
Period (yrs) | $91^{+27}_{-27}$ | $107^{+21}_{-28}$ | $107^{+22}_{-23}$
Parallax (mas) | $37.647\pm 0.071$ | $37.8952\pm 0.00547$ | $37.8953\pm 0.0060$
Barycenter $\mu_{\alpha}$ (mas/yr) | $-144.90^{+0.13}_{-0.13}$ | $-144.934^{+0.074}_{-0.067}$ | $-144.923^{+0.073}_{-0.070}$
Barycenter $\mu_{\delta}$ (mas/yr) | $-135.21^{+0.32}_{-0.32}$ | $-134.99^{+0.26}_{-0.24}$ | $-134.86^{+0.27}_{-0.23}$
RV Zero Point (m/s) | $-129^{+31}_{-68}$ | $1^{+1.6}_{-0.74}$ | $-168^{+61}_{-38}$
RV Jitter (m/s) | $2.53^{+0.19}_{-0.18}$ | $2.52^{+0.19}_{-0.17}$ | $17.9^{+4.3}_{-3.3}$

Note. — Reference time for $\lambda_{\textrm{ref}}$ is 2010.0.

Figure 1: HD 33632 Ab imaged with Subaru/CHARIS in 2022. Both images are post-
processed with pyklip ADI+SDI mode and have identical color scaling. The
region outside of 1.1 arcseconds has been masked for clarity. An additional
source is present at the very edge of the full CHARIS field (not shown), but
has been previously identified as a background star in Currie et al. (2020).

Independent measurements of the mass, age, effective temperature, and surface
gravity of brown dwarfs (BDs) are exceptionally valuable for testing models of
substellar formation and evolution. While temperature and surface gravity can
be estimated in BDs from photometry and spectra alone (e.g. Knapp et al. 2004;
Rice et al. 2009; Martin et al. 2017), those properties cannot be translated
into a precise mass determination. Furthermore, mass and age are degenerate
for BDs of a given temperature and luminosity (Burrows et al., 1989, 2001;
Baraffe et al., 2003), meaning it is only possible to make model independent
estimates of both mass and age for BDs that occur in binaries or as companions
to host stars. In systems with a host star, mass can be measured via
gravitational dynamics and age can be bounded by stellar age indicators. These
so-called benchmark BDs (Pinfield et al., 2006) with independent dynamical
masses, plus stellar age and metallicity estimates, make a up a relatively
small, but valuable fraction of the BD population. On the order of a few dozen
BDs have dynamical mass measurements among thousands of identified BDs.

Directly imaged BDs that are simultaneously detected with astrometric and/or
radial velocity accelerations of a host star are an important source of
benchmark BDs. Astrometric acceleration and radial velocity measurements over
long baselines constrain the mass of the BD, while its emission spectrum and
luminosity can be directly measured with high-contrast imaging spectrometers
(e.g. Gemini Planet Imager, Macintosh et al. 2014, VLT/SPHERE, Beuzit et al.
2019, Subaru/CHARIS, Groff et al. 2016). Furthermore, long-term imaging of the
BD can quantify orbital eccentricity with implications for formation history
(e.g. Forgan et al. 2018; Bowler et al. 2020). The interest of BDs with
concurrent direct imaging and spectroscopy, and measured dynamical masses, is
elevated by a small, but increasing number of BDs in this category with masses
or ages inconsistent with evolutionary models (e.g. Gl 229 B, Brandt et al.
2020, 2021a, HD 4113 C, Cheetham et al. 2018, HD 176535 Ab Li et al. 2023,
Lewis et al. in prep.). These objects are interesting probes of the complexity
of substellar evolution since explaining them may require modifications or
additions to current BD models. Unresolved binaries and factors like clouds
have been suggested (Li et al., 2023) as possible explanations for the tension
between models and observations.

HD 33632 Ab is a previously discovered (Currie et al., 2020) BD companion to a
F8V host star. The dynamical mass has been estimated using absolute astrometry
from the Hipparcos-Gaia Catalog of Accelerations (HGCA) updated for Gaia early
Data Release 3 (Brandt, 2018, 2021) and relative astrometry from two epochs of
imaging from Keck/NIRC2 and Subaru/CHARIS published at discovery (Currie et
al., 2020; Brandt et al., 2021a). A low-resolution ($R\sim 19$) JHK spectrum
from Subaru/CHARIS has also been matched to observational template spectra of
field BDs in Currie et al. (2020), classifying HD 33632 Ab as a $1.7\pm
0.4\,$Gyr, $50^{+5.6}_{-5.0},\text{M}_{\text{J}}$ object near the L-T spectral
transition. Concurrent with this study, high-resolution K-band spectra have
been obtained with Keck/KPIC to measure the rotational and radial velocity,
and chemical composition of HD 33632 Ab (Hsu et al., 2024).

In this paper, we provide further analysis of the mass, orbit, and spectra of
HD 33632 Ab utilizing a new epoch of data from the Coronagraphic High Angular
Resolution Imaging Spectrograph (CHARIS), an integral field spectrograph (IFS)
at the Subaru telescope, in H- and K- band high-resolution mode ($\sim
4\times$ greater spectral resolution than the broadband mode). In Section 2 we
describe our new observations and the data reduction methodology, followed by
relative astrometry measurements in Section 3.1, and orbit and mass fitting in
3.2. We then extract spectra in 4.1 and fit them to templates and models in
Section 4.2, including fits to new generations of disequilibrium and cloudy BD
spectral models and evolutionary tracks not explored in past work. Finally, we
compare our results for the derived mass, effective temperature, and surface
gravity from our analysis to previous measurements of age, and to BD evolution
models in Section 5.

##  2 Observations and Data Reduction

We observed HD 33632 A with Subaru/CHARIS (Groff et al., 2016) behind the
Subaru Coronagraphic Extreme Adaptive Optics (SCExAO) system (Guyon et al.,
2010; Vogt et al., 2010) on UT 2022-10-07. A summary of the properties of the
HD 33632 system are compiled in Table 1. The M dwarf companion HD 33632 B is
separated by $\sim 34^{\prime\prime}$ ($\sim 900\,$au projected distance), and
is therefore not observed. We acquired 58 exposures with integration times of
40 seconds (38 minutes total) in the H band filter in CHARIS’s high-resolution
mode (R$\sim 75$, 20 spectral samples from $1.44746$ to $1.79091\,\mu$m)
between 13:56 and 14:40 UTC, followed immediately by 69 K band exposures (17
spectral samples from $2.01514$ to $2.36803\,\mu$m) until 15:54 UTC. The first
three K band exposures have 40-second integration times, while the remaining
were increased to 60 seconds, all of which are used in data reduction ($68$
minutes total). Transit of HD 33632 occurred at 14:31 UT, such that the
observations approximately bookend the time of maximum field rotation. The
observations were taken with the 139 mas inner working angle (IWA) Lyot
coronagraph and a 25 nm “astrogrid” (Sahoo et al., 2020) added by fast
deformable mirror (DM) modulation.

Initial data reduction steps including reads-to-ramp conversion and $\chi^{2}$
microspectra extraction to create the 3D IFS data cubes are performed with the
CHARIS pipeline described in Brandt et al. (2017). The resulting cubes are 135
by 135 pixels with a $2^{\prime\prime}.2\times 2^{\prime\prime}.2$ field-of-
view ($16.15\,$ mas per lenslet). Post-processing steps including image
registration using satellite spots and modelling and subtraction of the
stellar PSF are performed with the CHARIS arm of the open-source pyKLIP
package (Wang et al., 2015; Chen et al., 2023). For initial data reduction and
astrometry, we utilize Karhunen Loéve Image Processing (KLIP, Soummer et al.
2012; Pueyo 2016) implementations of both angular differential imaging (ADI,
Marois et al. 2006) and spectral differential imaging (SDI, Racine et al.
1999) PSF subtraction techniques. We find that signal-to-noise (S/N) of the
companion is maximized around 11 KL modes and an angular and spectral movement
exclusion criteria of 2 pixels. S/N is measured with the SNR map feature of
pyKLIP, which calculates the noise in concentric annuli from the center with a
small masking region to minimize bias from the real point source. Post-
processed images are shown in Figure 1.

Table 2: HD 33632 Ab Measured Relative Astrometry  Date | $\Delta\alpha$ | $\Delta\delta$ | Sep. | PA
---|---|---|---|---
(UT) | (mas) | (mas) | (mas) | (∘)
2018-10-18 | $-761\pm 5$ | $-176\pm 4$ | $781\pm 5$ | $257.0\pm 0.4$
2018-11-01 | $-753\pm 5$ | $-178\pm 5$ | $774\pm 5$ | $256.7\pm 0.4$
2020-08-31 | $-740\pm 5$ | $-95\pm 3$ | $746\pm 5$ | $262.8\pm 0.4$
2022-10-07 | $-708.1\pm 4.5$ | $16.4\pm 5.0$ | $708.3\pm 4.5$ | $271.3\pm 0.4$
Figure 2: Posterior distributions of orbital elements and primary and
secondary mass of HD 33632 Ab. These results closely match those of Brandt et
al. (2021a). A visualization of the orbit is shown in the upper right for
varying secondary masses. The reference time for $\lambda_{\text{ref}}$ is
2010.0. Figure 3: Derived orbits with varying secondary masses compared to
astrometric and radial velocity data. From top left to bottom right plots show
separation and position angle between primary and secondary from Subaru/CHARIS
and Keck/NIRC2 imaging (Currie et al. 2020; this work), barycentric proper
motions from Hipparcos and Gaia (Perryman et al., 1997; Gaia Collaboration et
al., 2023a), and radial velocities from the Lick (Fischer et al., 2014).

##  3 Astrometry and Orbit Fitting

###  3.1 Relative Astrometry

We measure new relative astrometry of HD 33632 Ab using the Bayesian KLIP
forward-model (KLIP-FM) astrometry feature of pyKLIP (Wang et al., 2015).
Details of the algorithm are described in Pueyo (2016) and Wang et al. (2016),
and will only be summarized here.

The KLIP-FM astrometry methodology is, at its core, PSF fitting with the
critical addition of forward-modelling to account for systematic over or self-
subtraction that ADI, SDI, and KLIP impart on the companion PSF which may
distort the astrometry. For CHARIS, a companion PSF model is produced using
the satellite spots prior to the subtraction of the stellar PSF and a best
guess of the companion’s spectrum. The PSF model is created as an average of
all four satellite spots across all exposures on a per wavelength basis. We
visually verify that the satellite spots have stable PSFs over the entire
sequence. The relative brightness of the PSFs are then modulated by the
approximate observed companion spectrum. The approximate companion spectrum
does not need to be highly precise, with $20\%$ spectral errors causing a
factor of a few less variation in the resulting astrometry than speckle noise.

This companion PSF model is used as input for Bayesian Markov chain Monte
Carlo (MCMC, implemented with the emcee Python package, Foreman-Mackey et al.
2013) simulations to estimate the most likely centroid position of the
companion PSF and its associated uncertainty, which can be translated to
astrometric separation and position angle (PA). An additional complication of
the forward-model is that stellar speckle noise (the dominant noise source) is
correlated. pyKLIP considers correlated noise on the spatial extent of single
speckles ($\sim\lambda/D$) using Gaussian process methods to model covariance,
which is an additional fitted parameter in the Bayesian framework and is
reflected in the final astrometric uncertainties.

We perform astrometry on the H\- and K-band image cubes separately and compare
the results as an assessment of bias and uncertainty. For the companion
spectrum, we try a linear interpolation of the spectrum observed by Currie et
al. (2020), as well as our own retrieved spectra recursively, described in the
next section, with negligible difference in the results. No astrometric
calibrator (typically a well-characterized binary or compact star cluster) was
observed along with HD 33632, therefore we use the results of Chen et al.
(2023) for instrumental calibration of PA and platescale. We adopt the PA zero
point nearest to our observations ($-2.03^{\circ}\pm 0.44$), which
coincidentally is the same as the average PA zero point from all tabulated
measurements. Chen et al. (2023) demonstrate that the PA and platescale are
stable to $0.02$ mas/lenslet and $0.08^{\circ}$ respectively. Our relative
astrometry is presented along with previous relative astrometry in Table 2.

Initial astrometry of the K-band data resulted in an unexpected measured
separation $1.5\sigma$ greater than the separation measured in H-band
astrometry, while the measured PA is well within $1\sigma$ between H and K
image cubes. To search for a source of the discrepancy, we varied all tuneable
parameters of the PSF subtraction and astrometric fit, including KLIP
aggressiveness, the input spectra, contrast, and pixel position of the
companion, correlated noise parameters, and MCMC chain parameters, all with
inconsequential results. As a next step, we performed a leave-one-out test
where we performed all post-processing and KLIP-FM astrometry steps again
while sequentially leaving single data cubes out of the analysis. H-band data
proved to be exceptionally stable through this analysis with less than a
single milliarcsecond variation in separation from the removal of any data
cube, while K-band data showed much greater sensitivity to individual data
cubes, with as much as 5 milliarcseconds in variability ($1\sigma$ of the
reported uncertainty in separation from speckle noise). We then compared the
KLIP-FM astrometry to more simplistic centroid fitting with a 2D guassian to
the companion PSF in the post-processed data cubes. Centroid fitting of both
H\- and K-band data gives astrometric separation consistent with the
separation reported by KLIP-FM for H band. We therefore conclude from these
tests that there may be an image centering issue for the K-band astrometry
from KLIP-FM and choose to adopt the H-band KLIP-FM relative astrometry as our
data point for this epoch. We should note that adopting either the H\- or
K-band KLIP-FM astrometry from this step results in derived orbital
parameters, described in the next section, that are consistent within
$1\sigma$ of each other such that this issue is insignificant for final
results.

We also check the relative astrometry of the background source previously
identified by Currie et al. (2020). Due to proximity to the edge of the field,
we are unable to perform a KLIP-FM fit of its location, however, a simple
centroid fit places it at [E,N]′′
$=[-0.^{\prime\prime}68,-1.^{\prime\prime}32]$, only $\sim 30\,$ mas from its
expected location predicted solely from the proper motion of HD 33632 A.

###  3.2 Orbit Fitting and Mass Estimate

We fit the orbit and mass of HD 33632 Ab using the orvara python package
(Brandt et al., 2021b), which jointly fits radial velocity (RV) observations,
relative astrometry, and absolute astrometry using a Bayesian MCMC framework.
This updates the work of Currie et al. (2020) and Brandt et al. (2021a) for HD
33632 with the addition of our new epoch of relative astrometry for the BD
companion, but with otherwise the same methodology. This new epoch occurs
approximately two years after the last epoch, covering a small, but measurable
fraction of the estimated century long orbital period. We include the same
data previously used by Brandt et al. (2021a) for HD 33632 Ab, namely the
Hipparcos-Gaia Catalog of Accelerations updated for Gaia Data Release 3
(Brandt, 2021; Gaia Collaboration et al., 2023a), Lick archival radial
velocity data (Fischer et al., 2014), and previous relative astrometry from
Subaru/CHARIS and Keck/NIRC2 observations by Currie et al. (2020). Brandt et
al. (2021a) showed that the inclusion of the outer companion HD 33632 B had
insignificant impact on the orbital fit for HD 33632 Ab (and otherwise
increased the complexity of the MCMC sampling); therefore, we ignore the outer
companion in our orbit and mass fit. Our fitting priors are log-flat for semi-
major axis, companion mass, and radial velocity jitter, gaussian for the
primary mass, $\sin{i}$ for inclination, and flat for all other parameters,
identical to those of Brandt et al. (2021a). We also test uniform priors for
companion mass and semi-major axis, with nearly identical results to log-flat
priors.

A corner plot of the posterior distributions for the orbital elements and
primary and companion masses is shown in Figure 2. Orbit plots compared to the
observed relative astrometry, proper motions, and radial velocities are shown
in Figure 3. The $\chi^{2}$ residual between the best fit orbit and observed
separation and PA are $1.03$ and $3.19$ respectively. Overall, analysis with
addition of our new epoch of relative astrometry produces very similar results
to those from Currie et al. (2020) and Brandt et al. (2021a). Argument of
periastron, which is poorly constrained in all fits, is the only parameter not
within $1\sigma$ of prior results. A summary of the derived orbital properties
and companion mass is presented in Table 1. Zero eccentricity cannot be ruled
out, but is $<0.35$ with $95\%$ confidence.

Contemporaneous to this study, Hsu et al. (2024) have published a refined
orbit and dynamical mass for HD 33632 Ab which does not include our new epoch
of relative astrometry, but does consider a new direct radial velocity
measurement of HD 33632 Ab. Using the orbitize! orbit fitting code (Blunt et
al., 2020), they report a dynamical mass of $37\pm 7\,\text{M}_{\text{J}}$,
$\sim 2\sigma$ lower than our analysis. Their derived orbital parameters
otherwise mostly differ by approximately $1\sigma$ from our results. We do not
consider this radial velocity measurement in this paper, and a future paper
should therefore re-analyze the orbit and mass using all available data. We do
discuss the results of our spectral fitting in the next section in the context
of both of these derived masses.

Figure 4: Top: H\- and K-band extracted spectra from 2022 compared to the
broadband spectrum from 2020 presented in Currie et al. (2020), also from
Subaru/CHARIS. Both datasets come from Subaru/CHARIS. H and K are separately
scaled to best match the broadband data. Bottom: Comparison of observed
spectra (ADI only) to spectra of SpeX Prism optical and infrared standards
(Burgasser & Splat Development Team, 2017; Burgasser, 2007; Looper et al.,
2007; Burgasser et al., 2004). Spectral type T0 is the best fit to the data,
but L8 and T1 are also plausible fits. Note the top plot is in $F_{\nu}$ units
while the bottom plots are in $F_{\lambda}$ units.

##  4 Spectral Extraction and Analysis

###  4.1 Spectral Extraction and Spectrophotometry

H\- and K-band spectra are extracted using the KLIP-FM functionality of
pyKLIP. The process is similar to that used for relative astrometry, but in
this case the brightness of the same wavelength dependent forward model PSF is
adjusted to match the companion PSF after KLIP post-processing. More details
can be found in Greenbaum et al. (2018). We extract spectra with both the
ADI+SDI and ADI-only modes to check for notable systematic differences.
ADI+SDI and ADI-only spectra are consistent, although uncertainties in the
ADI+SDI spectra are greater. Uncertainties on spectral data points are
calculated using an injection and recovery technique. 11 copies of the forward
model PSF with the same spectrum as the extracted spectrum of the real
companion are injected into the data cubes before post-processing at the same
separation with an even angular spacing ($30^{\circ}$ between injections). We
extract the spectrum of each with the same methodology as the real companion.
Our error bars are equal to the full range of flux values retrieved from the
injected PSFs for each wavelength (i.e. they are more conservative than a
$1\sigma$ limit of the variation). Our H and K extracted spectra for both ADI
and ADI+SDI are shown in the top panel of Figure 4 in comparison to the
broadband JHK spectrum from Currie et al. (2020), also from Subaru/CHARIS.

We perform spectrophotometry to obtain passband magnitudes for comparison to
Currie et al. (2020). The H- and K- band $f_{\lambda}$ spectra are multiplied
by the response function of the $H_{\text{MKO}}$ and $K_{\text{S,MKO}}$
filters respectively, then integrated and converted to a magnitude using the
$f_{\lambda}$ zero points (Tokunaga et al., 2002; Tokunaga & Vacca, 2005;
Rodrigo et al., 2012; Rodrigo & Solano, 2020). We measure
$H_{\text{MKO}}=15.92\pm 0.08$ and $K_{\text{S,MKO}}=15.54\pm 0.12$.
$H_{\text{MKO}}$ is within $1\sigma$ of that reported by Currie et al. (2020),
while $K_{\text{S,MKO}}$ is $\sim 1.5\sigma$ fainter. The $K-$ band CHARIS
spectrum does not fully cover the $K_{\text{S,MKO}}$ bandpass by a small
margin, which likely explains the slight discrepancy. In any case, the
discrepancy is not large enough to impact interpretation of HD 33632 Ab’s
spectral type.

Figure 5: $\chi_{\nu}^{2}$ residuals from fits of ATMO-2020 grid models
(Phillips et al., 2020) (top) and Sonora-Bobcat grid models (Marley et al.,
2021) (bottom) to extracted H- and K-band spectral data. H-band residuals are
blue, K-band are red. The size of the circle is inversely proportional to the
residual squared — larger symbols indicate better fits. Columns represent
models of equilibrium chemistry (CEQ), weak, and strong non-equilibrium
chemistry (NEQ, vertical mixing) for ATMO-2020, and varying metallicity for
Bobcat models. Contours of constant mass and age are overplotted as blue
shaded and gray dotted lines respectively. Regions of the preferred mass and
age from dynamical constraints and stellar age indicators are shaded, with
double shading in the temperature and gravity region that agrees with both the
independently-constrained mass and age. The largest circles show the preferred
effective temperature and gravity from spectral fitting, which can be compared
with the grey shaded regions indicating previous constraints. Fits to both
ATMO-2020 and Sonora-Bobcat models prefer a hotter, younger BD than the age
suggests, although they may be nearly consistent within the grid spacing.
Solar metallicity Bobcat models are the best fitting and most consistent
models with mass and age of the sets shown here.

###  4.2 Spectral Fitting

We fit the extracted spectra from §4.1 to BD spectral templates and
atmospheric models to estimate the effective temperature, surface gravity,
cloudiness, strength of disequilibrium chemistry, and C/O ratio in the
atmosphere of HD 33632 Ab. When compared to evolutionary tracks in the
following subsections, these variables provide a consistency check on the
independent dynamical mass and age estimates for this BD, in addition to a
basic characterization of atmospheric conditions.

####  4.2.1 H and K Spectral Fits

To fit each template or model to the observed spectra, we arbitrarily
normalize the spectral data then multiply the model by a factor such that the
reduced $\chi^{2}$ residual ($\chi_{\nu}^{2}$) is minimized. Each set of
template and model spectra come with varying native resolutions, all of which
are higher than our observed spectral resolution. Each one is therefore
convolved with a Gaussian and resampled to match our observed spectral
resolution and sampling before it is compared to the observed data. To avoid
additional spectral uncertainty from the relative flux calibration of H and K
bands, we fit H and K independently, calculating a separate $\chi_{\nu}^{2}$
value for both bands for each model.

SpeX Prism Library of Spectra — The SpeX Prism library is a large collection
of low-resolution ($R=75-200$) observed BD and dwarf star spectra between
$0.65$ and $2.55$ microns (Burgasser & Splat Development Team, 2017). We use
our own code to compare our spectral data to the L and T dwarf optical and
infrared standards which include 15 spectra sequentially from L0 to T8
spectral types (Burgasser & McElwain, 2006; Burgasser et al., 2004, 2007;
Burgasser, 2007; Looper et al., 2007; Burgasser et al., 2006). We confirm the
results of Currie et al. (2020), finding that the best fit to both H and K
bands is the T0 template, with T1 and L8 also consistent with our spectra. A
comparison of our spectra to these spectral templates is shown in the bottom
panels of Figure 4.

ATMO-2020 Model Grids — ATMO-2020 is a set of cloudless 1D radiative-
convective models for T and Y dwarf objects with self-consistent evolution
tracks (Phillips et al., 2020). In addition to a grid with equilibrium
chemistry, it also includes two grids with “weak” and “strong” disequilibrium
chemistry that results from atmospheric vertical mixing. This disequilbrium
chemistry is described by a gravity dependent $\log K_{zz}$ parameter that is
greater for lower gravity objects, but is otherwise constant throughout a
given atmosphere. Plots showing the $\chi_{\nu}^{2}$ parameter fit with
varying surface gravity and effective temperature models are shown in the top
half of Figure 5. Spectral models that better fit the data (have lower
$\chi_{\nu}^{2}$) are represented by larger circles, with the best fitting
models having the largest circles. Fits to H- and K- band data are shown as
blue and red circles respectively. Our K- band spectrum is somewhat
degenerate, with multiple surface gravity and effective temperature models
having $\chi_{\nu}^{2}<1$, while no model has $\chi_{\nu}^{2}<1$ with our
H-band spectrum. The best fitting models are $T_{\mathrm{eff}}\sim 1700\,$K
and $\log g\sim 5.0$. Overplotted in the figure are contours of mass and age
from the ATMO-2020 evolutionary tracks, with the estimated mass and age limits
of HD 33632 Ab shaded in grey. Based on these evolutionary models, our best
spectral fits could almost be consistent with the derived dynamical mass and
age for HD 33632 Ab within the grid spacing, although they are slightly ($\sim
100-200\,K$) hotter than expected. Increased vertical mixing / disequilibrium
chemistry slightly favors cooler models, but the results are only subtly
different from equilibrium chemistry. We also fit to models with differing
metallicities, the results of which are not shown for ATMO-2020, but are
consistent with the metallicity results from comparison with Sonora Bobcat
models described next.

Sonora Bobcat Model Grids — Sonora Bobcat is another independent set of
cloudless 1D radiative-convective models with self-consistent evolution, and
C/O ratio as an additional model variable (Marley et al., 2021). We plot the
result of our fits to Bobcat models in the bottom half of Figure 5, which is
analogous to the top half for ATMO-2020 models. While the fitting results are
very similar to those of ATMO-2020, fits to Bobcat prefer a slightly cooler
BD. $H-$ band data are a much better fit to Sonora-Bobcat than ATMO-2020, and
prefer a BD that is $\sim 50\%$ younger than stellar age estimates, although
once again, it may be nearly consistent within the model grid spacing. A metal
rich atmosphere produces overall much poorer fits, while a metal poor
atmosphere suggests younger, lower gravity models. While there is
observational evidence that HD 33632 A is metal poor (Soubiran et al., 2016;
Chen et al., 2000; Mishenina et al., 2004; Valenti & Fischer, 2005; Gaia
Collaboration et al., 2023a, b), the best fit models in that regime are
incompatible with the estimated stellar age and dynamical mass. Only solar C/O
ratios are tested for Bobcat fits, but it is a variable for fits to Sonora Elf
Owl models, where low vertical mixing results with varying C/O ratios would be
expected to give very similar results to that of Sonora Bobcat. High-
resolution spectra from Hsu et al. (2024) support an approximately solar C/O
ratio for HD 33632 Ab.

Sonora Elf Owl Model Grids — Sonora Elf Owl is a recent development in the
Sonora family of BD model grids, which adds disequilibrium chemistry
(atmospheric vertical mixing) as a variable (Mukherjee et al., 2024). Compared
to ATMO-2020 (Phillips et al., 2020), $\log K_{zz}$ in Elf Owl is not set by
gravity and has a range of $\log K_{zz}$ for each atmosphere of a given
temperature, gravity, metallicity, and C/O ratio. Our fitting results with Elf
Owl are shown in Figure 6. The results are similar to those of Sonora Bobcat.
Increased vertical mixing (greater $\log K_{zz}$) notably favors cooler models
the stronger it becomes when compared to our K-band spectrum. This is the same
finding as with ATMO-2020. This could suggest that the discrepancy in age
between our favored spectral models and the stellar age could be explained by
increasing vertical mixing. Indeed, retrievals on high-resolution spectra of
HD 33632 Ab in Hsu et al. (2024) greatly prefer strong disequilibrium
chemistry. However, goodness-of-fit to our H-band spectrum degrades with
increasing $\log K_{zz}$. Therefore, we cannot claim that a high vertical
mixing constant is definitively preferred to a low constant based on our data
alone. We also fit against models varying in C/O ratio. Subsolar C/O results
in a preference for cooler models when compared against our K-band spectra,
but these models are simultaneously poor fits to H-band. Supersolar C/O
conversely leads to a preference for hotter models, which are fit well in H
and K, but would be less consistent with mass and age, assuming that BD
evolutionary tracks for supersolar C/O are very similar to that of solar C/O.
These trends are consistent across metallicity. Given this pattern, solar or
subsolar C/O ratios are preferred on the basis of better agreement with
dynamical mass and stellar age, which is also consistent with the finding of a
roughly solar C/O for HD 33632 Ab from Hsu et al. (2024). Note that the impact
of C/O ratio, metallicity, gravity, disequilibrium chemistry, and clouds on
low-resolution spectra like these can be partially degenerate, and firm
constraints cannot be made without flexible forward models that can tune each
parameter independently. No self-consistent evolutionary tracks are yet
available to plot for Elf Owl.

Sonora Diamondback Model Grids — Sonora Diamondback adds clouds to the Sonora
family of models, parameterized by sedimentation efficiency $f_{\mathrm{sed}}$
(Morley et al., 2024). Thicker clouds have a lower value of
$f_{\mathrm{sed}}$. Results of fitting to the Diamondback grid are shown in
Figure 7. The results clearly favor models that are either cloud-free or have
thin clouds. The best models with thin clouds are slightly cooler than those
of cloud-free, however, they also have slightly lower preferred surface
gravities, and are therefore not necessarily better fits with the dynamical
mass and stellar age. Cloud-free is a better fit to H and K spectra
individually, but the proximity of HD 33632 Ab to the L-T transition (Currie
et al., 2020) indicates that it probably has at least some clouds (e.g.
Burrows et al. 2006; Charnay et al. 2018).

####  4.2.2 Broadband JHK Spectral Fits

As a check of the results from the previous subsection, we also perform the
same spectral fitting analysis using the broadband spectrum ($R\sim 19$) from
(Currie et al., 2020) that includes J-band and will briefly summarize the
results here. When the broadband spectrum is compared against Sonora Bobcat
models, the fits are very poor across the grid range, with no model spectra
able to replicate the relative fluxes of molecular band peaks in J, H and K.
The relative flux in J-band is much higher in all cloudless models than in the
observed broadband spectrum, seen in Figure 4. The addition of clouds in the
Sonora-Diamondback models reduces the relative flux of $J-$band, bringing the
model spectrum into much better agreement with the data. This is strong
evidence for the presence of clouds on HD 33632 Ab, which cannot be well-
constrained by spectra in individual NIR bands, but has greater impact in
broadband. Comparison of the broadband spectrum to Sonora Diamondback prefers
an $f_{\mathrm{sed}}$ parameter in the range $3.0-8.0$, with thicker clouds
and cloudless models being significantly worse fits. The best-fit temperature
and gravity for $f_{\mathrm{sed}}=8.0$ is within the model grid spacing of the
estimated mass and age of HD 33632 Ab, and may be fully consistent with
evolutionary tracks with that cloud parameter.

####  4.2.3 Possibility of Lower Mass for HD 33632 Ab

Hsu et al. (2024) have recently published a direct measurement of the radial
velocity of HD 33632 Ab, which they use to derive a $2\sigma$ lower dynamical
mass than we have presented ($37\,\text{M}_{\text{J}}$ compared to
$51.6\,\text{M}_{\text{J}}$). Here we briefly re-assess the best-fit spectral
models in the context of the lower mass from Hsu et al. (2024). A dynamical
mass of $37\,\text{M}_{\text{J}}$ corresponds to a roughly $-0.25$ shift in
the $log{g}$ gravity of the shaded preferred mass region of Figures 5 and 7.
For ATMO-2020, Sonora Bobcat, and Sonora Diamondback models, this would
improve the overlap between the best-fitting models and the dynamical mass for
solar metallicity and cloud-free models. There is otherwise no significant
change in the preference for models with differing metallicity, C/O ratio,
cloud sedimentation parameter, or strength of disequilibrium chemistry.
However, this lower mass would simultaneously increase the disagreement
between the age and temperature of the BD, requiring HD 33632 Ab to be either
$\sim 500\,$K hotter than expected for its age, or significantly younger than
the estimated age of the host star. In short, a lower mass would be more in
line with best-fit spectral models, but would be more difficult to explain
from an evolutionary standpoint.

##  5 Summary

We obtained a new epoch of Subaru/CHARIS imaging and spectroscopy of HD 33632
Ab at high angular resolutions H- and K-bands. We extract relative photometry
and update the orbit and dynamical mass of the BD, finding a mass
($51.6^{+5.4}_{-4.8}\,\text{M}_{\text{J}}$) and orbit
($a=23.7^{+3.1}_{-3.8}\,\text{au}$ $e=0.143^{+0.13}_{-0.089}$) consistent with
Currie et al. (2020) and Brandt et al. (2021a), but a dynamical mass $\sim
2\sigma$ higher than Hsu et al. (2024).

Spectra in H and K are independently compared to new generations of BD
atmosphere and evolutionary models, which include disequilibrium chemistry via
vertical mixing, clouds, and variable metallicity and C/O ratio. While our
spectra are consistent with the broadband spectrum from (Currie et al., 2020)
and a spectral type near the L-T transition, no spectral model is
simultaneously a best fit to the observed spectra and consistent with
evolutionary tracks predicted from the derived age and preferred mass range
from this study or Hsu et al. (2024). Instead, the best spectral fits tend to
prefer a $\sim 50\%$ younger object than the stellar age based on a hotter
derived temperature and/or lower gravity than expected. The best fit spectral
models have effective temperatures of $1600-1800\,$K, whereas the stellar age
($\sim 1-2.5\,$Gyr) and dynamical mass range predicts temperatures of $\sim
900-1500\,$K, with lower predicted temperature for lower masses. We find that
disequilibrium chemistry and light clouds independently reduce the discrepancy
by favoring a lower effective temperature, though the inferred surface
gravities remain slightly too low. More complex forward modeling of the
spectrum which considers both disequilibrium chemistry and clouds
simultaneously, and predicts evolution, may be able to resolve the
disagreement. Indeed, some plausible spectral models (i.e. Sonora Diamondback)
are consistent with our derived mass and age within the grid spacing, and
photometry alone is consistent with the expected luminosity of the BD from
evolutionary tracks (Currie et al., 2020; Brandt et al., 2021a). We do not
consider if HD 33632 Ab being an unresolved binary could be consistent with
the age, mass, and spectra, however, this explanation seems presently
unnecessary given the plausibility of atmospheric only explanations.

Future imaging and data releases from Gaia will refine the orbit and dynamical
mass. Continued radial velocity monitoring should also be pursued to improve
dynamical mass measurements, with precise radial velocity monitoring being
more powerful in the near term than imaging alone. Foremost, our study
highlights the necessity of considering complexities like clouds and
disequilibrium chemistry when comparing benchmark BDs to evolutionary models,
especially for BDs near the L-T transition.

Figure 6: Same as Figure 5, but for fits to Sonora Elf Owl disequilibrium
models. Increasing vertical mixing, parameterized by $K_{zz}$, prefers cooler
models for fits to K data only, while fits to H data are less sensitive. No
mass or age contours are presently available to plot for Elf Owl. Columns
represent models of differing metallicity, while rows have changing vertical
mixing constants. Plotted models have a solar C/O ratio. Figure 7: Same as
Figure 5, but for fits to Sonora Diamondback cloudy models. Cloud free or thin
clouds are preferred to the thickest clouds for fits to H and K data alone. As
noted in Section 4.2.2, fitting Sonora Diamondback models to broadband data
strongly prefers some clouds to cloud free. Lower $f_{sed}$ represents thicker
clouds. Columns represent models of differing metallicity, while rows have
changing $f_{sed}$, which changes cloud thickness. Plotted models have a solar
C/O ratio. Only the cloudless and $f_{sed}=2.0$ models have mass and age
contours supplied. The plotted contours are the same mass and age values as
those plotted in Figure 5.

This research is based in part on data collected at the Subaru Telescope,
which is operated by the National Astronomical Observatory of Japan. The
authors wish to recognize and acknowledge the very significant cultural role
and reverence that the summit of Maunakea has always had within the indigenous
Hawaiian community. We are most fortunate to have the opportunity to conduct
observations from this mountain. This research has benefited from the SpeX
Prism Spectral Libraries, maintained by Adam Burgasser at
https://cass.ucsd.edu/ajb/browndwarfs/spexprism/. This work has made use of
data from the European Space Agency (ESA) mission Gaia
(https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and
Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular the institutions
participating in the Gaia Multilateral Agreement. This material is based upon
work supported by the National Science Foundation Graduate Research Fellowship
under Grant No. 2021-25 DGE-2034835 for author BLL. Any opinions, findings,
and conclusions or recommendations expressed in this material are those of the
authors(s) and do not necessarily reflect the views of the National Science
Foundation.

## References

  * Astropy Collaboration et al. (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167, doi: 10.3847/1538-4357/ac7c74
  * Baraffe et al. (2003) Baraffe, I., Chabrier, G., Barman, T. S., Allard, F., & Hauschildt, P. H. 2003, A&A, 402, 701, doi: 10.1051/0004-6361:20030252
  * Beuzit et al. (2019) Beuzit, J. L., Vigan, A., Mouillet, D., et al. 2019, A&A, 631, A155, doi: 10.1051/0004-6361/201935251
  * Blunt et al. (2020) Blunt, S., Wang, J. J., Angelo, I., et al. 2020, AJ, 159, 89, doi: 10.3847/1538-3881/ab6663
  * Bowler et al. (2020) Bowler, B. P., Blunt, S. C., & Nielsen, E. L. 2020, The Astronomical Journal, 159, 63, doi: 10.3847/1538-3881/ab5b11
  * Brandt et al. (2021a) Brandt, G. M., Dupuy, T. J., Li, Y., et al. 2021a, AJ, 162, 301, doi: 10.3847/1538-3881/ac273e
  * Brandt (2018) Brandt, T. D. 2018, ApJS, 239, 31, doi: 10.3847/1538-4365/aaec06
  * Brandt (2021) —. 2021, ApJS, 254, 42, doi: 10.3847/1538-4365/abf93c
  * Brandt et al. (2020) Brandt, T. D., Dupuy, T. J., Bowler, B. P., et al. 2020, The Astronomical Journal, 160, 196, doi: 10.3847/1538-3881/abb45e
  * Brandt et al. (2021b) Brandt, T. D., Dupuy, T. J., Li, Y., et al. 2021b, AJ, 162, 186, doi: 10.3847/1538-3881/ac042e
  * Brandt et al. (2017) Brandt, T. D., Rizzo, M., Groff, T., et al. 2017, Journal of Astronomical Telescopes, Instruments, and Systems, 3, 048002, doi: 10.1117/1.JATIS.3.4.048002
  * Burgasser (2007) Burgasser, A. J. 2007, ApJ, 659, 655, doi: 10.1086/511027
  * Burgasser et al. (2006) Burgasser, A. J., Burrows, A., & Kirkpatrick, J. D. 2006, ApJ, 639, 1095, doi: 10.1086/499344
  * Burgasser et al. (2007) Burgasser, A. J., Looper, D. L., Kirkpatrick, J. D., & Liu, M. C. 2007, ApJ, 658, 557, doi: 10.1086/511518
  * Burgasser & McElwain (2006) Burgasser, A. J., & McElwain, M. W. 2006, AJ, 131, 1007, doi: 10.1086/499042
  * Burgasser et al. (2004) Burgasser, A. J., McElwain, M. W., Kirkpatrick, J. D., et al. 2004, AJ, 127, 2856, doi: 10.1086/383549
  * Burgasser & Splat Development Team (2017) Burgasser, A. J., & Splat Development Team. 2017, in Astronomical Society of India Conference Series, Vol. 14, Astronomical Society of India Conference Series, 7–12, doi: 10.48550/arXiv.1707.00062
  * Burrows et al. (1989) Burrows, A., Hubbard, W. B., & Lunine, J. I. 1989, ApJ, 345, 939, doi: 10.1086/167964
  * Burrows et al. (2001) Burrows, A., Hubbard, W. B., Lunine, J. I., & Liebert, J. 2001, Reviews of Modern Physics, 73, 719, doi: 10.1103/RevModPhys.73.719
  * Burrows et al. (2006) Burrows, A., Sudarsky, D., & Hubeny, I. 2006, ApJ, 640, 1063, doi: 10.1086/500293
  * Charnay et al. (2018) Charnay, B., Bézard, B., Baudino, J.-L., et al. 2018, The Astrophysical Journal, 854, 172, doi: 10.3847/1538-4357/aaac7d
  * Cheetham et al. (2018) Cheetham, A., Ségransan, D., Peretti, S., et al. 2018, A&A, 614, A16, doi: 10.1051/0004-6361/201630136
  * Chen et al. (2023) Chen, M., Wang, J. J., Brandt, T. D., et al. 2023, RAS Techniques and Instruments, 2, 620, doi: 10.1093/rasti/rzad039
  * Chen et al. (2000) Chen, Y. Q., Nissen, P. E., Zhao, G., Zhang, H. W., & Benoni, T. 2000, A&AS, 141, 491, doi: 10.1051/aas:2000124
  * Currie et al. (2020) Currie, T., Brandt, T. D., Kuzuhara, M., et al. 2020, ApJ, 904, L25, doi: 10.3847/2041-8213/abc631
  * Cutri et al. (2003) Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2003, VizieR Online Data Catalog: 2MASS All-Sky Catalog of Point Sources (Cutri+ 2003), VizieR On-line Data Catalog: II/246. Originally published in: 2003yCat.2246….0C
  * Fischer et al. (2014) Fischer, D. A., Marcy, G. W., & Spronck, J. F. P. 2014, ApJS, 210, 5, doi: 10.1088/0067-0049/210/1/5
  * Foreman-Mackey et al. (2013) Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306, doi: 10.1086/670067
  * Forgan et al. (2018) Forgan, D. H., Hall, C., Meru, F., & Rice, W. K. M. 2018, MNRAS, 474, 5036, doi: 10.1093/mnras/stx2870
  * Gaia Collaboration (2020) Gaia Collaboration. 2020, VizieR Online Data Catalog: Gaia EDR3 (Gaia Collaboration, 2020), VizieR On-line Data Catalog: I/350. Originally published in: 2021A&A…649A…1G; doi:10.5270/esa-1ug, doi: 10.26093/cds/vizier.1350
  * Gaia Collaboration et al. (2023a) Gaia Collaboration, Vallenari, A., Brown, A. G. A., et al. 2023a, A&A, 674, A1, doi: 10.1051/0004-6361/202243940
  * Gaia Collaboration et al. (2023b) Gaia Collaboration, Creevey, O. L., Sarro, L. M., et al. 2023b, A&A, 674, A39, doi: 10.1051/0004-6361/202243800
  * Greenbaum et al. (2018) Greenbaum, A. Z., Pueyo, L., Ruffio, J.-B., et al. 2018, AJ, 155, 226, doi: 10.3847/1538-3881/aabcb8
  * Groff et al. (2016) Groff, T. D., Chilcote, J., Kasdin, N. J., et al. 2016, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9908, Ground-based and Airborne Instrumentation for Astronomy VI, ed. C. J. Evans, L. Simard, & H. Takami, 99080O, doi: 10.1117/12.2233447
  * Guyon et al. (2010) Guyon, O., Martinache, F., Garrel, V., et al. 2010, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 7736, Adaptive Optics Systems II, ed. B. L. Ellerbroek, M. Hart, N. Hubin, & P. L. Wizinowich, 773624, doi: 10.1117/12.857878
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357, doi: 10.1038/s41586-020-2649-2
  * Hsu et al. (2024) Hsu, C.-C., Wang, J. J., Xuan, J. W., et al. 2024, arXiv e-prints, arXiv:2405.08312, doi: 10.48550/arXiv.2405.08312
  * Knapp et al. (2004) Knapp, G. R., Leggett, S. K., Fan, X., et al. 2004, AJ, 127, 3553, doi: 10.1086/420707
  * Li et al. (2023) Li, Y., Brandt, T. D., Brandt, G. M., et al. 2023, MNRAS, 522, 5622, doi: 10.1093/mnras/stad1315
  * Looper et al. (2007) Looper, D. L., Kirkpatrick, J. D., & Burgasser, A. J. 2007, AJ, 134, 1162, doi: 10.1086/520645
  * Macintosh et al. (2014) Macintosh, B., Graham, J. R., Ingraham, P., et al. 2014, Proceedings of the National Academy of Science, 111, 12661, doi: 10.1073/pnas.1304215111
  * Marley et al. (2021) Marley, M. S., Saumon, D., Visscher, C., et al. 2021, ApJ, 920, 85, doi: 10.3847/1538-4357/ac141d
  * Marois et al. (2006) Marois, C., Lafrenière, D., Doyon, R., Macintosh, B., & Nadeau, D. 2006, The Astrophysical Journal, 641, 556, doi: 10.1086/500401
  * Martin et al. (2017) Martin, E. C., Mace, G. N., McLean, I. S., et al. 2017, ApJ, 838, 73, doi: 10.3847/1538-4357/aa6338
  * Mishenina et al. (2004) Mishenina, T. V., Soubiran, C., Kovtyukh, V. V., & Korotin, S. A. 2004, A&A, 418, 551, doi: 10.1051/0004-6361:20034454
  * Morley et al. (2024) Morley, C. V., Mukherjee, S., Marley, M. S., et al. 2024, arXiv e-prints, arXiv:2402.00758, doi: 10.48550/arXiv.2402.00758
  * Mukherjee et al. (2024) Mukherjee, S., Fortney, J. J., Morley, C. V., et al. 2024, arXiv e-prints, arXiv:2402.00756, doi: 10.48550/arXiv.2402.00756
  * pandas development team (2020) pandas development team, T. 2020, pandas-dev/pandas: Pandas, latest, Zenodo, doi: 10.5281/zenodo.3509134
  * Perryman et al. (1997) Perryman, M. A. C., Lindegren, L., Kovalevsky, J., et al. 1997, A&A, 323, L49
  * Phillips et al. (2020) Phillips, M. W., Tremblin, P., Baraffe, I., et al. 2020, A&A, 637, A38, doi: 10.1051/0004-6361/201937381
  * Pinfield et al. (2006) Pinfield, D. J., Jones, H. R. A., Lucas, P. W., et al. 2006, MNRAS, 368, 1281, doi: 10.1111/j.1365-2966.2006.10213.x
  * Pueyo (2016) Pueyo, L. 2016, ApJ, 824, 117, doi: 10.3847/0004-637X/824/2/117
  * Racine et al. (1999) Racine, R., Walker, G. A. H., Nadeau, D., Doyon, R., & Marois, C. 1999, PASP, 111, 587, doi: 10.1086/316367
  * Rice et al. (2009) Rice, E. L., Barman, T., Mclean, I. S., Prato, L., & Kirkpatrick, J. D. 2009, The Astrophysical Journal Supplement Series, 186, 63, doi: 10.1088/0067-0049/186/1/63
  * Rodrigo & Solano (2020) Rodrigo, C., & Solano, E. 2020, in XIV.0 Scientific Meeting (virtual) of the Spanish Astronomical Society, 182
  * Rodrigo et al. (2012) Rodrigo, C., Solano, E., & Bayo, A. 2012, SVO Filter Profile Service Version 1.0, IVOA Working Draft 15 October 2012, doi: 10.5479/ADS/bib/2012ivoa.rept.1015R
  * Sahoo et al. (2020) Sahoo, A., Guyon, O., Lozi, J., et al. 2020, AJ, 159, 250, doi: 10.3847/1538-3881/ab88cd
  * Soubiran et al. (2016) Soubiran, C., Le Campion, J.-F., Brouillet, N., & Chemin, L. 2016, A&A, 591, A118, doi: 10.1051/0004-6361/201628497
  * Soummer et al. (2012) Soummer, R., Pueyo, L., & Larkin, J. 2012, ApJ, 755, L28, doi: 10.1088/2041-8205/755/2/L28
  * Tokunaga et al. (2002) Tokunaga, A. T., Simons, D. A., & Vacca, W. D. 2002, Publications of the Astronomical Society of the Pacific, 114, 180, doi: 10.1086/338545
  * Tokunaga & Vacca (2005) Tokunaga, A. T., & Vacca, W. D. 2005, PASP, 117, 421, doi: 10.1086/429382
  * Valenti & Fischer (2005) Valenti, J. A., & Fischer, D. A. 2005, ApJS, 159, 141, doi: 10.1086/430500
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261, doi: 10.1038/s41592-019-0686-2
  * Vogt et al. (2010) Vogt, F., Martinache, F., Guyon, O., et al. 2010, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 7736, Adaptive Optics Systems II, ed. B. L. Ellerbroek, M. Hart, N. Hubin, & P. L. Wizinowich, 773612, doi: 10.1117/12.857915
  * Wang et al. (2015) Wang, J. J., Ruffio, J.-B., De Rosa, R. J., et al. 2015, pyKLIP: PSF Subtraction for Exoplanets and Disks, Astrophysics Source Code Library, record ascl:1506.001.  http://ascl.net/1506.001
  * Wang et al. (2016) Wang, J. J., Graham, J. R., Pueyo, L., et al. 2016, AJ, 152, 97, doi: 10.3847/0004-6256/152/4/97