JWST-TST High Contrast: JWST/NIRCam Observations of the Young Giant Planet β Pic b

This paper is part of a series by the James Webb Space Telescope (JWST) Telescope Scientist Team, ¹² which uses Guaranteed Time Observer (GTO) time awarded by NASA in 2003 (PI: M. Mountain) for studies in three different subject areas: (a) Transiting Exoplanet Spectroscopy (lead: N. Lewis); (b) Exoplanet and Debris Disk High-contrast Imaging (lead: M. Perrin); and (c) Local Group Proper-motion Science (lead: R. van der Marel). Here, we present new results in the second area; previously reported results across these areas include Grant et al. (2023) and A. Gressier et al. (2024, in preparation), Ruffio et al. (2023) and Rebollido et al. (2024), and Libralato et al. (2023), respectively. A common theme of these investigations is the desire to pursue and demonstrate science for the astronomical community at the limits of what is made possible by the exquisite optics and stability of JWST.

In this paper, we continue to present the first coronagraphic observations of the β Pic system with the JWST (Gardner et al. 2006, 2023) Near-InfraRed Camera (NIRCam; Rieke et al. 2003, 2023), spanning a wavelength range from ∼1.7 to 5 μm. We previously reported in Rebollido et al. (2024) on the edge-on debris disk around β Pic as seen in these same NIRCam data, along with coronagraphic observations at 15.5 and 23 μm obtained with the Mid-InfraRed Instrument (MIRI; Rieke et al. 2015; Wright et al. 2015, 2023), which led to the discovery of a dramatic curved "tail" of dust that appears to extend sharply away from the plane of the disk and may be unbound debris from a recent major collision of large planetesimals.

In this work, we turn our attention to the planets within the β Pic system. The outer giant planet β Pic b is detected in all six observed NIRCam filters close to its greatest orbital elongation. The inner giant planet β Pic c remains undetected behind the coronagraphic masks (COMs) of NIRCam. Both known planets β Pic b and β Pic c remain undetected in the MIRI coronagraphy data, due to the tremendously bright debris disk in the mid-IR, so we choose to not further discuss those MIRI data here, except for one brief usage in Section 3.3 for evaluating the nature of additional candidate companions detected in the NIRCam data. However, there are also MIRI Medium Resolution Spectroscopy (MRS; Wells et al. 2015; Argyriou et al. 2023) observations of the β Pic system (program 1294—PI: C. Chen) that detect the outer planet β Pic b at longer wavelengths (∼5–7 μm; Worthen et al. 2024); we make use of their spectrum in joint analyses below.

The remainder of this paper is organized as follows. Section 2 presents the observations and data reduction, including two independent approaches to point-spread function (PSF) subtraction. Section 3 presents the measured photometry and astrometry of β Pic b, and detection limits on lower-mass planets in the outer system. In Section 4, we discuss the implications for β Pic b's atmospheric and bulk properties and its evolutionary status, using the NIRCam measurements, prior literature values, and exoplanetary model grids. We summarize our findings and conclude in Section 5.

The NIRCam observations targeted β Pic in six different filters from ∼1.7 to 5 μm, two of which were observed using NIRCam's short-wavelength (SW) channel with a pixel scale of ∼31 mas (F182M and F210M) and four of which were observed using its long-wavelength (LW) channel with a pixel scale of ∼63 mas (F250M, F300M, F335M, and F444W). The F444W wideband filter was included because it is the most sensitive to additional, as-yet-undiscovered planets in the system (e.g., Carter et al. 2023). The SW observations were taken using NIRCam's 210R round COM; these were the first ever science observations taken with SW channel coronagraphy. The LW observations were taken using the 335R mask. These masks were chosen because they are the smallest available round COMs, providing the best performance at small angular separations while retaining access to the entire 360° field of view. We note that the observations were taken sequentially in the two channels, as parallel operation of the SW and LW channels was not yet enabled at the time.

The observations were conducted using the standard high-contrast imaging strategy of two rolls on the science target plus a PSF reference star. The science target β Pic was observed at two roll angles offset by a position angle (PA) of ∼10° from one another. This enables PSF subtraction using angular differential imaging (ADI) techniques (Liu 2004; Marois et al. 2006), which are frequently used to detect faint companions (e.g., Bowler 2016; Currie et al. 2023a). In addition, a PSF reference star (α Pic) was observed in sequence with the science target using the 5-POINT-BOX ¹⁵ small-grid dither (SGD) pattern. This enables reference differential imaging (RDI), including classical PSF subtraction methods, which are typically less prone to self- or oversubtraction of extended circumstellar structure such as the debris disk around β Pic (Lafrenière et al. 2007; Soummer et al. 2012). The SGD can be used to enhance the quality of the PSF subtraction (i.e., to suppress residual stellar speckles in the PSF-subtracted images) by sampling changes in the shape of the PSF as a function of the relative misalignment between the star and the COM (Soummer et al. 2014). This is necessary because the target acquisition (TA) procedure for NIRCam coronagraphic imaging can have uncertainties of typically up to ∼20 mas as measured in flight (Girard et al. 2022; Rigby et al. 2023). Table 1 shows a summary of the NIRCam coronagraphic observations.

Table 1. Summary of the JWST/NIRCam Observations of the Young Giant Planet β Pic b Presented in This Paper (JWST Program ID 1411)

Target	Filter	Readout Pattern	Dither Pattern	Nints/Ngroups/Nframes	${T}_{\exp }$	PA
					(s)	(deg)
β Pic	F182M	RAPID	NONE	90/4/1	1506.9	84.34
β Pic	F182M	RAPID	NONE	90/4/1	1506.9	94.34
α Pic	F182M	RAPID	5-POINT-BOX	10/4/1	837.2	76.46
β Pic	F210M	RAPID	NONE	90/4/1	1506.9	84.34
β Pic	F210M	RAPID	NONE	90/4/1	1506.9	94.34
α Pic	F210M	RAPID	5-POINT-BOX	10/4/1	837.2	76.46
β Pic	F250M	BRIGHT2	NONE	80/10/2	1710.5	84.57
β Pic	F250M	BRIGHT2	NONE	80/10/2	1710.5	94.57
α Pic	F250M	BRIGHT2	5-POINT-BOX	8/10/2	855.2	76.39
β Pic	F300M	BRIGHT2	NONE	80/10/2	1710.5	84.57
β Pic	F300M	BRIGHT2	NONE	80/10/2	1710.5	94.57
α Pic	F300M	BRIGHT2	5-POINT-BOX	8/10/2	855.2	76.39
β Pic	F335M	SHALLOW4	NONE	35/10/4	1833.4	84.57
β Pic	F335M	SHALLOW4	NONE	35/10/4	1833.4	94.57
α Pic	F335M	SHALLOW4	5-POINT-BOX	4/10/4	1047.7	76.39
β Pic	F444W	SHALLOW4	NONE	35/10/4	1833.4	84.57
β Pic	F444W	SHALLOW4	NONE	35/10/4	1833.4	94.57
α Pic	F444W	SHALLOW4	5-POINT-BOX	4/10/4	1047.7	76.39

Note. The science target is β Pic and the PSF reference target is α Pic. The observations were taken on 2023 March 18. The total exposure time ${T}_{\exp }$ sums up all five dither positions for the reference target.

Download table as:  ASCIITypeset image

All data were reduced with the spaceKLIP ¹⁶ community pipeline (Kammerer et al. 2022b; Carter et al. 2023). This pipeline was developed during JWST commissioning and Early Release Science and has since been updated and improved steadily. The data reduction process is briefly summarized in the following, with the individual steps described in more detail in the subsequent subsections.

• 1.  
  Reduce the raw ("uncal") files with the jwst ¹⁷ stage 1 and 2 pipelines (Bushouse et al. 2023), first to "rateints" and then to flux-calibrated "calints" files, with several custom adaptions specific to NIRCam coronagraphy data that are implemented in spaceKLIP.
• 2.  
  Use the spaceKLIP image processing library to clean bad pixels and align the PSFs in preparation for PSF subtractions.
• 3.  
  Perform PSF subtractions. Several complementary methods were employed, as described below. In particular, we used the pyKLIP ¹⁸ community pipeline (Wang et al. 2015) invoked through spaceKLIP to perform PSF subtractions using projections on Karhunen–Loève eigenimages (Karhunen–Loève Image Processing or KLIP; Soummer et al. 2012). We also applied model-constrained reference differential imaging (MCRDI; Lawson et al. 2022) as an alternative approach to help separate disk and planetary fluxes.
• 4.  
  For the KLIP subtractions, we use the pyKLIP forward-modeling routines (Pueyo 2016) with NIRCam coronagraphy PSF models generated through spaceKLIP to extract companion astrometry and photometry and compute contrast curves.

2.2.1. Calibration of Individual Images

In the first stage, we downloaded the "uncal" files from the MAST archive and processed them with the Coron1Pipeline and Coron2Pipeline within spaceKLIP. These two pipelines are customized implementations of the jwst Detector1Pipeline and Image2Pipeline. We used version 1.12.1 of the jwst pipeline and calibration context jwst_1174.pmap for reducing the NIRCam data. As in Carter et al. (2023), we disabled the flagging of pixels that are diagonal neighbors to saturated pixels, skipped the pipeline's dark subtraction step, flagged a 4 pixel wide border around the edge of the subarrays to be used as "pseudo" reference pixels, and set the jump rejection threshold to 4. In addition, and as already discussed in Rebollido et al. (2024), we also employed a custom 1/f stripe noise correction step that is now implemented in spaceKLIP and yields a mostly cosmetic improvement. The outputs of this stage are flux-calibrated "calints" files in units of megajanskys per steradian. This reduction made use of the first in-flight photometric calibrations for NIRCam coronagraphy delivered in fall 2023. ¹⁹ That on-sky flux calibration accounts for the reduced throughput of the NIRCam coronagraphic Lyot stops (pupil plane masks), but not for the attenuation of the focal plane masks for sources located at small angular separation from the observed target. ²⁰ It also accounts for the wavelength-dependent transmission of the COM substrate (Krist et al. 2010). The attenuation of the coronagraphic (focal plane) masks is accounted for later in the companion fitting step, when the model PSF of the companion is computed (see Section 3.1).

In the second stage, we prepared the reduced images for processing with pyKLIP. First, we median-subtracted each image to remove the sky background flux, then cleaned bad pixels. Due to the extended debris disk visible in all NIRCam images, the automatic bad pixel identification routine in spaceKLIP did not work satisfyingly well and we used a custom bad pixel map in addition to the DO_NOT_USE pixels flagged by the jwst pipeline. This custom map was obtained by inspecting the cleaned images by eye and manually flagging bad pixels that the jwst pipeline had missed. While we found no additional bad pixels in the SW images, we identified and manually flagged 28 bad pixels in the LW images. The routines that were used to clean the flagged bad pixels are the same as the ones in Carter et al. (2023).

We include additional steps to mitigate the flux from the debris disk and to deal with spatial undersampling. The shortest filters in each channel are spatially undersampled, which can lead to numerical artifacts when shifting or rotating images. Before recentering and aligning the images, we first applied a Gaussian high-pass filter to subtract spatially extended and smooth flux from the debris disk. We explored a variety of filter sizes, ranging from standard deviations of 1–9 pixels, as well as no high-pass filtering at all. To address the undersampling, we then blurred all images with a Gaussian kernel with a FWHM of

$\begin{eqnarray}&&{\mathrm{FWHM}}_{\mathrm{filter}}=\sqrt{{\mathrm{FWHM}}_{\mathrm{desired}}^{2}-{\mathrm{FWHM}}_{\mathrm{current}}^{2}}\end{eqnarray} \tag{ 1 }$

$\begin{eqnarray}&&=\sqrt{{(2.3\cdot 1.5)}^{2}-{\left(\displaystyle \frac{{\lambda }_{\min }}{D\cdot s}\right)}^{2}}.\end{eqnarray} \tag{ 2 }$

Here, ${\lambda }_{\min }$ is the minimum wavelength of the filter bandpass, D is the effective telescope diameter (5.2 m for NIRCam coronagraphy, due to the undersized Lyot stops), s is the detector pixel scale, and the desired FWHM is 2.3 · 1.5 = 3.45 pixels. ²¹ The blurring helps to avoid "Fourier ripples" when numerically shifting or rotating undersampled images (due to the Gibbs phenomenon; e.g., Gottlieb & Shu 1997).

2.2.2. Image Alignment and TA Performance

Image alignment to precisely register the science images and reference PSFs is critical for achieving high-quality PSF subtractions. For each filter independently, the first science target image of each observation was recentered on the position of the target star. This is a difficult task, because the star is located behind the COM so that its center cannot be easily determined. We followed the same procedure that has already used by Greenbaum et al. (2023) to center the star based on its speckle pattern with the help of a simulated PSF computed from a contemporaneous wave-front map of the telescope using WebbPSF_ext ²² (Perrin et al. 2014; Girard et al. 2022). Greenbaum et al. (2023) reported a centroiding error of ∼7 mas for this procedure. However, given the bright debris disk around β Pic, the centroiding error is likely larger in our case. Once the first image of each observation was recentered, we then used the same image registration routine as in Kammerer et al. (2022b) to align all subsequent science and PSF reference target images to the first one. As discussed in Rigby et al. (2023), the pointing stability of JWST is ∼1 mas, and we are able to recover a ∼1 mas rms jitter between individual images. We are also able to recover the injected 5-POINT-BOX SGD for the PSF reference target and to measure the relative TA offset between the science and the PSF reference target.

Figure 1 shows this alignment for one of the SW and one of the LW filters, revealing good TA performance with an offset of only ∼4 mas for the SW channel, but poor performance with an offset of ∼50 mas for the LW channel PSF reference observation. This is atypical; in fact, we believe it to be one of the largest observed TA offsets in the ensemble of NIRCam coronagraphy data yet obtained. We subsequently determined that the reason for this outlier TA performance is that the SIMBAD coordinates and proper motion of α Pic are not sufficiently accurate, ²³ so that after the initial telescope slew and fine guidance sensor acquisition, α Pic was offset by ∼07 from the expected position. As a result, its PSF did not entirely fit within the NIRCam LW TA subarray. This resulted in an inaccurate centroid measurement. We note that the same systematic position offset of α Pic was also seen in the NIRCam SW and MIRI TA images, but for those modes the relative size of the PSF with respect to the TA subarray is smaller, so that the PSFs did still fit within the TA subarrays and the centroid measurements were accurate. We expect the poor NIRCam LW TA performance to impact the quality of the RDI PSF subtraction at small separations from the target for the LW data sets, specifically to reduce the contrast performance by a factor of ∼3–4 at a separation of 1'' with respect to what would have been achievable with a good TA performance (see Figure 16 in Girard et al. 2022).

Zoom In Zoom Out Reset image size

After image registration, we then averaged all aligned integrations in an individual exposure (each dither position is an individual exposure) to speed up the subsequent processing with pyKLIP. This is beneficial for rerunning the reduction several times with different parameters later on in the paper. We note that averaging the images has a negligible impact on the KLIP reduction because the line-of-sight pointing and PSF of JWST are extremely stable. For comparison, we also ran the pyKLIP processing with our baseline parameters once without averaging the images and found that ∼99.9% of the signal is contained in the first KLIP mode of each individual exposure, which represents the noise-weighted average of the individual images in that exposure. In other words, there is negligible information loss in averaging together all integrations within an exposure.

2.2.3. PSF Subtractions (KLIP)

We tested and compared several methods for PSF subtractions, including RDI (Lafrenière et al. 2007), ADI (Marois et al. 2006), and MCRDI (Lawson et al. 2022). The former two were employed using the implementation of the KLIP algorithm (Soummer et al. 2012) in pyKLIP (Wang et al. 2015) through the spaceKLIP community pipeline and are described in this section, and MCRDI is described in the following section.

To remove post-coronagraph residual stellar speckles from the NIRCam images of β Pic, we employed ADI, RDI, and ADI+RDI techniques through pyKLIP. The objective of the KLIP algorithm is to project each science target image onto a covariance-weighted orthogonal basis of eigenimages of the PSF reference library, which consists of the PSF reference images of all five dither positions in the case of RDI, or the science target images from the other telescope roll in the case of ADI. The reference library constructed from the SGD observations of the PSF reference star will hence capture changes in the PSF shape as a function of the TA offset in order to provide a good representation of the science target image, which is similarly affected by a TA uncertainty of typically up to ∼20 mas (Girard et al. 2022). Then the projection of the science target image onto the reference library is subtracted from the science target image itself, ideally leaving behind flux from off-axis companions and circumstellar structure.

For β Pic, RDI subtraction alone does not work well, due to the extended debris disk contributing significant flux in all NIRCam images and additionally the poor PSF reference star TA performance in the LW channel. The KLIP routine tends to oversubtract the residual stellar speckles because it tries to minimize the total residual flux, including the flux of the debris disk. ADI alone yields a much better speckle subtraction, but its performance at small angular separations is limited, due to the small telescope roll angle of only 10°.

Combining ADI and RDI yields the visually best performance. The KLIP-subtracted NIRCam images are shown in Figure 2. We did not split the images into multiple annuli or subsections because we did not find any improvement in fake companion injection and recovery tests from doing that (see Appendix A). We note that these PSF subtractions are optimized for measurements of the planet β Pic b; different choices are optimal for the study of the debris disk, as discussed in Rebollido et al. (2024).

Zoom In Zoom Out Reset image size

2.2.4. PSF Subtractions (MCRDI)

In MCRDI, the effects of RDI oversubtraction (e.g., Pueyo 2016) are avoided by explicitly assuming the presence of circumstellar flux while the stellar PSF model is constructed (Lawson et al. 2022). For this purpose, a synthetic model of the circumstellar scene is optimized using standard forward-modeling techniques for an initial unconstrained RDI reduction. A final MCRDI reduction is then carried out using the resulting best-fit model (Lawson et al. 2022).

The MCRDI procedure predominantly follows the approach outlined in Rebollido et al. (2024), but uses slightly different input data resulting from the pre-processing details described in Sections 2.2.1 and 2.2.2 (e.g., the larger blurring kernel to mitigate the effects of the Gibbs phenomenon). The circumstellar model prescription effectively follows the one in Rebollido et al. (2024)—assuming a circumstellar scene that is the superposition of a simple ring-like disk and a point source. We point out that this simplified disk model is only used for computing the PSF subtraction coefficients in the RDI framework; the PSF-subtracted scene images do however still comprise the disk in its full complexity. For convolution, we use a more finely spatially sampled grid of synthetic PSFs from WebbPSF (Perrin et al. 2014), now sampling the origin along with 12 logarithmically spaced radial positions at each of eight linearly spaced azimuthal positions (for 97 spatial samples in total). The finer spatial sampling has no discernible impact on the final MCRDI image, but does slightly change the disk-model-subtracted residuals within the Inner Working Angle (IWA). Otherwise, the MCRDI procedure is as described in Rebollido et al. (2024).

The MCRDI-reduced NIRCam images are shown in Figure 3, zoomed in to show just the area centered on β Pic b, after the subtraction of the disk model for each filter. Fairly good PSF subtractions are achieved for the two SW filters, but a small blob can be seen in the residuals of the F250M data toward the northeast of the host star position. We note that a similarly shaped blob can also be seen in other of the LW filters if the color stretch is optimized for the residuals instead of the planet β Pic b and the blob's brightness scales with the brightness of the other residual stellar speckles. We hence attribute this blob to the subpar TA performance for the LW observations.

Zoom In Zoom Out Reset image size

The debris disk around β Pic is clearly detected in all JWST/NIRCam images and thus affects the extracted photometry of the giant planet β Pic b. The impact of this disk on the extracted planet flux in the KLIP subtractions can be complex. Residual disk flux may be added on top of the planet flux, leading to an overestimation of the planet flux, or the disk may be oversubtracted and part of the planet flux may be removed, leading to an underestimation of the planet flux. Self-subtraction artifacts from running ADI on extended sources may further introduce biases in the planet flux that differ from these expectations. An important parameter that controls the impact of the disk on the extracted planet flux is the size of the high-pass filter that is used to remove (part of) the disk flux from the NIRCam images (see Section 2.2.1). Fake companion injection and recovery tests (see Appendix A) reveal that the best performance in the two SW filters is achieved using ADI alone (with a high-pass filter size of 7 pixels), whereas the best performance in the four LW filters is achieved by combining ADI and RDI (with a high-pass filter size of 3 pixels), so that we adopt these as the baseline approach for extracting planet photometry in the respective cases. We note that even with these optimized parameters, the retrieved flux in the F250M data set is still much higher than the injected one in these tests, showing that the coarse detector sampling at 2.5 μm together with the extended debris disk prohibits an accurate companion flux measurement using KLIP techniques.

In the MCRDI approach, rather than having parameters for spatial filtering, instead there are additional free parameters for the disk model, which are fit as part of the MCRDI process. In effect, these disk model parameters become nuisance parameters that must be solved for simultaneously, increasing the overall dimensionality of the fitting process, but providing a loosely physical approach for subtracting the disk's light around the location of the planet. As in Rebollido et al. (2024), we note that the disk model fit in MCRDI is intentionally simplified and not intended as a physically correct model of the β Pic debris disk, but it does suffice to largely remove the disk light while avoiding biases from oversubtraction (Lawson et al. 2022).

Given the substantial impact of β Pic's bright disk on the PSF subtractions, much effort and iteration went into optimizing these methods and their parameters until the two independent PSF subtraction approaches yielded consistent measurements for the photometry of β Pic b. This cross-validation was essential for establishing confidence in the removal of subtraction systematics. Potential sources of systematic errors were taken into account by inflating the error bars on the planet photometry, as described in the third paragraph of Section 3.1.

We extracted new photometry for the giant planet β Pic b in all six observed filters, spanning a wavelength range of ∼1.7–5 μm. This new photometry is especially interesting, since L- and M-band observations from the ground are typically affected by a high thermal background from the Earth's atmosphere, leading to large uncertainties in the measured photometry (e.g., Stolker et al. 2020). Our new measurements from JWST are not affected by this issue. Furthermore, they also open up bands such as F182M and F300M, which are inaccessible from the ground due to telluric absorption from water vapor in the Earth's atmosphere.

The photometry was obtained by forward modeling the planet PSF in all of the six observed NIRCam bands and fitting it to the KLIP-subtracted images using forward-model matched filtering (Pueyo 2016). The methods used follow the ones in Kammerer et al. (2022b) and Carter et al. (2023), but now using new in-flight photometric calibrations for NIRCam coronagraphy instead of the preflight ones. Briefly, a model PSF of β Pic b was computed at its expected distance from the COM center with the WebbPSF_ext package, using the closest-in-time available JWST wave-front measurement from MAST. This PSF model accounts for the attenuation by the COM. Due to different TA and therefore COM offsets, this needed to be done separately for each of the two telescope rolls. The model PSFs were then fed into the pyKLIP BKA routine (Wang et al. 2016), which forward models the KLIP-processed PSF to account for companion self-subtraction (due to ADI) and KLIP algorithm throughput losses. Finally, the forward-modeled PSFs were fitted to the KLIP-subtracted images using Markov Chain Monte Carlo (MCMC) sampling with emcee (Foreman-Mackey et al. 2013). The planet photometry (and astrometry: see the next section), including uncertainties, were obtained from the converged MCMC chains. An advantage of JWST over ground-based imagers is that the planet photometry can directly be measured in the images, since they are flux-calibrated by the jwst pipeline, except for the COM throughput, which is accounted for separately in our PSF forward models. We note that this entire procedure is integrated within a single user-friendly function within spaceKLIP. Figure 4 shows our best-fit forward-model PSFs compared to the data.

Zoom In Zoom Out Reset image size

The new β Pic b NIRCam photometry is shown in Table 2 and Figure 5. The KLIP photometry was obtained using the baseline extraction parameters defined in Section 2.3. Besides the statistical uncertainties from the MCMC fit of the forward-modeled PSF to the data, we estimate systematic uncertainties of ∼2% from the uncertainty of the absolute flux calibration of NIRCam, ²⁴ ∼1% from numerical inaccuracies in the forward-modeled PSFs from WebbPSF, and ∼2% from the uncertainty on the COM throughput (due to the uncertainty on the companion/mask position), which need to be added in quadrature to the statistical uncertainties. The companion flux uncertainty from the uncertainty on the companion/mask position was derived by shifting the companion and mask closer together and farther away from each other by their respective position uncertainties and evaluating the resulting change in flux using WebbPSF. The uncertainty in the WebbPSF models is motivated by Weisz et al. (2024), who find less than 1% photometric errors for WebbPSF models in regular NIRCam imaging, and we budget 1% for this effect to be slightly conservative for the more complex optical situation of coronagraphy. In addition, and based on the results from the fake companion injection and recovery tests (Appendix A), we added another systematic error of 5% to the KLIP photometry, which accounts for the typical offset that we find between the injected and recovered fake companion flux.

Zoom In Zoom Out Reset image size

Figure 5 shows the KLIP and MCRDI photometry together with existing photometry and spectra and a DRIFT-PHOENIX model atmosphere of β Pic b from the literature. The model atmosphere was fitted to the literature data only (and not to the new NIRCam data), using the species ²⁵ toolkit (Stolker et al. 2020). The purpose of this figure is to show whether the new NIRCam photometry agrees with the existing literature data and qualitatively assess whether the NIRCam data points at wavelengths inaccessible from the ground are broadly consistent with models used to interpret the ground-based data; model atmosphere fitting including the new NIRCam photometry will be more rigorously conducted in Section 4.1. Figure 5 shows that the new NIRCam photometry agrees well with a DRIFT-PHOENIX model atmosphere fitted to previously existing literature data of β Pic b. Both the KLIP and MCRDI data points fall within 3σ from the prediction of the model atmosphere.

Notably, all NIRCam data points fall below the prediction of the model atmosphere. In the L and M bands, this makes the new NIRCam photometry consistent with the existing VLT/NACO photometry, which also falls slightly below the prediction of the model atmosphere. In the K band, however, there is a ∼10%–20% discrepancy between the new NIRCam photometry and the existing VLTI/GRAVITY spectrum. We note that the blurring, although being accounted for in the forward modeling, might be responsible for blending some of the planet flux with the disk flux and leading to an underestimation of the planet flux. The disk being brighter at shorter wavelengths makes this especially problematic in the K band. We further note that the fake companion injection and recovery tests (Figure 14) show that the retrieved planet flux declines with increasing high-pass filter size. In these tests, a rather large FWHM was ideal for retrieving the injected flux for the fake companion, but it could well be that for the real planet β Pic b, a smaller FWHM might be ideal, so that we are underestimating its flux. This is extremely difficult to quantify, though, since the PSF- and disk-subtraction residuals are not exactly the same on the northeast and southwest side of the debris disk. Despite these uncertainties, it is worth mentioning that our NIRCam photometry is in good agreement with the NACO data point at ∼2.1 μm, though. Besides issues with the data reduction, astrophysical variability might also play a role for the measured discrepancies. While variability in the planet itself is expected to be below ∼5%, based on the observed variability of low-mass brown dwarfs and planetary mass objects (Metchev et al. 2015; Vos et al. 2020, 2022), the findings from Rebollido et al. (2024) suggest that the debris disk around β Pic is a highly dynamic environment. Moreover, K. Worthen et al. (2024, in preparation) find significant variability between their MIRI/MRS and archival Spitzer spectra (Lu et al. 2022) of the debris disk supporting this hypothesis. Nevertheless, the impact that this disk variability would have on the flux of the planet β Pic b is expected to be only a fraction of the disk variability itself and hence smaller than the observed discrepancy between GRAVITY and NIRCam, NACO, and GPI (see also Section 4.1), so that systematics in the flux calibration remain the more likely hypothesis.

Due to the complications of the debris disk affecting the planet flux measured in the KLIP-subtracted images and the advantages of using MCRDI in such cases, as discussed in Section 2.3, we decided to use the planet flux measurements from the MCRDI reductions as our baseline for the atmospheric characterization of β Pic b in Section 4.

From the KLIP-subtracted JWST/NIRCam images of the β Pic system, we also extracted relative astrometry of the well-known outer giant planet β Pic b. This astrometry is shown in Table 3 and Figure 6. It was obtained using KLIP ADI+RDI and a high-pass filter size of 3 pixels for all observed NIRCam bands. ²⁶ While the photometry of β Pic b was taken from the MCRDI reduction, we found that the planet astrometry from MCRDI is rather sensitive to the exact choice of disk model parameters (e.g., a one-component vs two-component disk model). On the other hand, the KLIP astrometry is much more consistent for different high-pass filter size choices, so we decided to use this more robust measurement.

Zoom In Zoom Out Reset image size

In the NIRCam images, β Pic b is detected close to its maximum orbital elongation on the northeast side of the disk. The astrometric measurements from NIRCam are not yet competitive with high-contrast imagers from the ground, mainly due to the difficulties in finding the position of the star behind the COM and the missing distortion correction in spaceKLIP. Together, they lead to systematic uncertainties on the order of ∼10 mas in the NIRCam astrometry. We further expect that the numerous residual stellar speckles caused by the poor reference star TA performance in the LW channel negatively affect the precision of the NIRCam astrometry, so that its systematic offset with respect to the GRAVITY orbit shown in Figure 6 is not a surprise. This expectation is supported by the LW channel astrometry generally agreeing worse with the GRAVITY orbit than the SW channel astrometry. Given these systematic errors and the precision of the existing GRAVITY astrometry of only ∼0.1–0.3 mas (Lacour et al. 2021), we did not attempt any orbital fits for the β Pic system with our new NIRCam measurements. Instead, we only show the orbits of both planets inferred by Lacour et al. (2021) from the GRAVITY data and overplot the new β Pic b NIRCam astrometry to demonstrate that we indeed detected the planet at its expected location. Improving the NIRCam astrometry and making it competitive with ground-based imagers is left for future work.

JWST/NIRCam coronagraphy has demonstrated exquisite sensitivity to new companions at wide separations from the host star, down to the sub-Jupiter-mass regime in the wideband filters (e.g., Carter et al. 2023). Here, we use the new NIRCam data of β Pic to search for previously undiscovered companions.

By visually inspecting the KLIP-subtracted images, we identify five other sources at separations >5'' from β Pic in addition to the known outer giant planet β Pic b (see Table 4 and Figure 7). We note again that NIRCam coronagraphy does not achieve the IWA required to detect the inner known giant planet β Pic c at a separation of ∼100 mas at the time of our observations. With the exception of one source (BG4), the additional sources are only detected in the NIRCam LW channel and are visible at a low number of KL modes (≤4) in all four LW filters. BG4 is also visible in the MIRI F1550C four-quadrant phase-mask data from Rebollido et al. (2024). Three of the five sources (BG1, BG2, and BG3) are located to the east of the edge-on debris disk and the other two sources (BG4 and BG5) are located to the west of the disk. The two sources to the west are clearly resolved and appear diffuse, so that we classify them as background galaxies. Two of the three sources to the east (BG2 and BG3) also appear resolved, so that we count four background galaxies in total. Only one of the sources to the east (BG1) appears point-like and unresolved. Its approximate position with respect to β Pic is ΔRA ≈ 6'' and ΔDec ≈ 8''. The approximate astrometry of all five sources can be found in Table 4.

Zoom In Zoom Out Reset image size

To determine the nature of the point-like source (BG1), we extract its photometry by fitting forward-modeled PSFs to the KLIP-subtracted images, similar to how we measured the photometry of β Pic b. Figure 8 shows color–magnitude diagrams with the point-like source indicated by an orange star and an evolutionary track for an 18.5 Myr old bound companion in blue. For reference, β Pic b is also shown by a red star. The evolutionary track was obtained by combining ATMO equilibrium chemistry models from Phillips et al. (2020) with Bern EXoplanet cooling curves (BEX) hot- and cold-start models from Linder et al. (2019). This was necessary to cover the entire planet mass range from ∼0.1 to 13 M_Jup to which JWST is sensitive. We first interpolated the evolutionary models at an age of 18.5 Myr, and then stitched together the ATMO and BEX models by interpolating linearly between them in log(mass) space. The transition between the two models happens between 0.5 and 1 M_Jup, where both models are defined. We focused on chemical equilibrium models without clouds to maintain consistency among the different companion mass ranges. We further assumed a distance of 19.44 pc (van Leeuwen 2007) for β Pic. The point-like source is far too blue in order to be a young and bound planetary mass companion. Instead, its colors are more consistent with a background star. We note that we did also inspect archival Hubble Space Telescope (HST)/Space Telescope Imaging Spectrograph and WFC3 data (program IDs 12551—PI: D. Apai and 11150—PI: J. Graham, respectively), but did not identify any background star whose position and proper motion would agree with the point-like source seen in the NIRCam data. Given the position of the source near the edge-on debris disk and the proper motion of β Pic, it seems possible that the source might have been right behind the disk in the archival HST images. In any case, its colors indicate it is almost certainly a star.

Zoom In Zoom Out Reset image size

Due to their unparalleled sensitivity, the NIRCam data also allow us to put deep constraints on the presence of additional companions in the β Pic system. First, contrast curves were computed from the KLIP-subtracted images, assuming a noise distribution following a student's t-distribution and correcting for small sample statistics following Mawet et al. (2014). To avoid being impacted by residual flux from the debris disk when computing the contrast curves, we here consider the high-pass-filtered data (with a high-pass filter size of 3 pixels). Nevertheless, since some amount of residual noise from the debris disk remains in the high-pass-filtered images, we always show two kinds of detection limits: (1) those in a region far away from the disk; and (2) those in a small region around the disk midplane. Then, using the ATMO + BEX hot-start evolutionary models introduced in the previous paragraph, the companion mass sensitivity of our observations was derived from the contrast curves. Figure 9 shows the contrast and companion mass limits obtained for the different filters as a function of the angular separation from the host star β Pic. The contrast with respect to β Pic was obtained by comparing the flux measured in the NIRCam images to the flux of β Pic according to a BOSZ A6V-type stellar atmosphere model (Bohlin et al. 2017) fitted to the archival 1–5 μm photometry of β Pic obtained from VizieR. ²⁷ The contrast limits reach down to better than 1e − 7 beyond ∼4'' (∼80 au). In the background-limited regime, this is about 10 times deeper in the L band and 30 times deeper in the M band if compared to the NACO detection limits reported in Stolker et al. (2019). The ADI+RDI and ADI reductions perform fairly similar at close-in and wide separations and there is an intermediate regime (∼05–2'') where ADI+RDI performs better than ADI alone. We note that this regime is where the coronagraphic optics create a strong quasi-static speckle field in the images. The RDI-only reduction performs the worst. As discussed in Section 2.2.2, this is due to the poor reference star subtraction that resulted from the suboptimal TA performance for the reference star observations of α Pic in the LW channel.

Zoom In Zoom Out Reset image size

Looking at the companion mass limits, the F444W wideband filter outperforms all other filters by a factor of ≳10 at wide separations and reaches a 5σ detection limit of ∼0.05 M_Jup (approximately a Neptune/Uranus mass) beyond ∼4'' (∼80 au) from the disk midplane. The far deeper sensitivity of the F444W filter has multiple reasons. First, the wideband filters have a significantly increased throughput compared to the medium-band filters, resulting in higher sensitivity. Second, the contrast between a young giant planet and an A-type host star is smaller at longer wavelengths (e.g., Beiler et al. 2023; Leggett & Tremblin 2023), which benefits the F444W filter the most. We can hence rule out the presence of companions away from the disk midplane with masses above ∼0.05 M_Jup between ∼4'' and 10'' (∼80–200 au) from β Pic. In the disk midplane, the companion mass limits are worse, due to residual photon noise and flux from the disk. We reach a 5σ sensitivity of ∼1 M_Jup beyond 2'' (∼40 au) from β Pic.

Figure 10 shows the position of β Pic b in a NIRCam color–magnitude and color–color diagram together with AMES-Dusty (Allard et al. 2001) and Sonora Bobcat (Marley et al. 2021) isochrones. For reference, two other substellar companions that have already been observed with JWST/NIRSpec spectroscopy (TWA 27 B—Luhman et al. 2023; Manjavacas et al. 2024; and VHS 1256 B—Miles et al. 2023) are also shown. These diagrams can inform the evolution of substellar companions from the correlation of their luminosity and spectral type and from correlations between the spectral slopes of two different wavelength regimes, respectively. Both diagrams were made with the species toolkit (Stolker et al. 2020). Given the dynamical mass of ∼9 ± 2 M_Jup of β Pic b (Nowak et al. 2020; Brandt et al. 2021), the young giant planet agrees better with the AMES-Dusty isochrone, which is not surprising, given that previous works have already identified β Pic b to have a cloudy atmosphere (e.g., Bonnefoy et al. 2013; Currie et al. 2013; Morzinski et al. 2015; Chilcote et al. 2017). Both isochrones assume equilibrium chemistry, which plays a role for the depth of the methane feature discussed below. While Sonora Bobcat describes a cloudless atmosphere with only rainout chemistry included, AMES-Dusty completely neglects the gravitational settling of the grains, so that layers of clouds automatically build up through condensation. The dust clouds absorb light at shorter wavelengths and re-emit it in the near-IR, leading to a redder color and thus an increased brightness in the near-IR if compared to a cloudless atmosphere at the same surface gravity (i.e., planet mass). In the color–color diagram, both isochrones predict that objects will become bluer in F300M–F444W with increasing mass, because the effective temperature also increases, moving the blackbody peak to shorter wavelengths. These colors are hence probing mostly the pseudo-continuum of the atmosphere. The combination of the F300M and the F335M filters provides a tracer of the methane absorption feature at ∼3.3 μm. The F300M−F335M color first becomes redder at small masses, because the depth of the methane absorption feature in the F335M filter decreases. Then, at a mass of about ∼7 M_Jup for AMES-Dusty and ∼14 M_Jup for Sonora Bobcat, the methane feature completely disappears and objects become bluer with increasing mass, as expected from a pure blackbody. It can be seen that the giant planet β Pic b is already in the regime where the methane feature has disappeared. We note, however, that the evolutionary models used here do not include disequilibrium chemistry, so that methane absorption already appears at higher effective temperatures compared to what has been observed (Yamamura et al. 2010; Zahnle & Marley 2014; Moses et al. 2016; Stolker et al. 2020). For a low-gravity L-type object such as β Pic b, methane absorption is not expected.

Zoom In Zoom Out Reset image size

To derive a set of physical parameters of β Pic b based on our JWST/NIRCam photometry, we fitted model atmospheres to its observed spectrophotometry using the species toolkit (Stolker et al. 2020). Besides the new JWST/NIRCam photometry, we included the Gemini/GPI spectra from Chilcote et al. (2017), the VLTI/GRAVITY spectrum from GRAVITY Collaboration et al. (2020), the recent JWST/MIRI spectrum from Worthen et al. (2024), the Magellan/VisAO and Gemini/NICI photometry from Males et al. (2014), and the VLT/NACO photometry from Bonnefoy et al. (2011), Currie et al. (2013), Stolker et al. (2019), and Stolker et al. (2020). We explored three different model atmosphere grids: the DRIFT-PHOENIX (Helling et al. 2008), the Exo-REM (Baudino et al. 2015; Charnay et al. 2018), and the BT-Settl (Allard et al. 2012) grids. All of them assume radiative-convective equilibrium to calculate the temperature structure of the atmosphere and include photospheric absorption by dust clouds. Cloudy model atmospheres have been found to provide good fits to the observed spectrophotometry of β Pic b by many previous works (e.g., Bonnefoy et al. 2013; Currie et al. 2013; GRAVITY Collaboration et al. 2020; Stolker et al. 2020). A cloudy atmosphere is also expected for β Pic b based on its effective temperature and the condensation temperatures of common refractory species in its atmosphere (e.g., Marley & Robinson 2015). To be able to fit the different model grids to the data, we first bin them to the spectral resolution of the respective instrument and calculate synthetic photometry for all considered filters. Then we interpolate the binned spectra and synthetic photometry linearly between the different grid points and infer the posterior distribution of the grid parameters using nested sampling with the PyMultiNest ²⁸ package (Buchner et al. 2014). Uniform distributions were chosen for the priors in the parameter ranges shown in Table 5.

Table 5. Prior Ranges Adopted for the Considered Model Atmosphere Grids

Model	Teff	$\mathrm{log}g$	(Fe/H)	C/O	R
	(K)	(dex)			(RJup)
BT-Settl	1000–3000	3.5–5.5	⋯	⋯	0.5–2.5
DRIFT-PHOENIX	1000–3000	3.0–5.5	–0.6–0.3	⋯	0.5–2.5
Exo-REM	1000–2000	3.0–5.0	–0.5–1.0	0.1–0.8	0.5–2.5

Note. The priors were distributed uniformly within the quoted ranges.

Download table as:  ASCIITypeset image

To investigate the characterization power of NIRCam by itself, we perform two types of fits: (1) using only the new NIRCam photometry; and (2) including the existing literature data in addition to the new NIRCam photometry. Given that β Pic b is such a well-studied object, these tests are ideally suited to check which atmospheric parameters our chosen set of NIRCam filters can constrain by themselves and which not. This knowledge will be especially interesting for future NIRCam observations of previously unknown exoplanetary companions. The best-fit parameter values for the three considered model atmosphere grids and the two types of fits are shown in Table 6. We note that fits using only the existing literature data without the new NIRCam photometry do not significantly differ from the ones including the new NIRCam photometry and hence are not shown or further discussed here. There are multiple reasons for this. First, the GPI and GRAVITY spectra have so many data points that the weight of the NIRCam photometry is very small. Second, there is already NACO photometry in the L and M bands and the error bars of the NIRCam photometry are too large to further constrain the best-fit models in this wavelength regime.

Table 6. Best-fit Bulk Parameters Inferred for the Young Giant Planet β Pic b Using Different Model Atmosphere Grids

Model	Teff	$\mathrm{log}g$	(Fe/H)	C/O	R	$\mathrm{log}L/{L}_{\odot }$
	(K)	(dex)			(RJup)	(dex)
BT-Settl (All data)	${1640}_{-4}^{+3}$	${3.51}_{-0.01}^{+0.01}$	⋯	⋯	${1.68}_{-0.01}^{+0.01}$	$-{3.71}_{-0.00}^{+0.00}$
BT-Settl (NIRCam only)	${1567}_{-71}^{+96}$	${3.84}_{-0.10}^{+0.17}$	⋯	⋯	${1.64}_{-0.17}^{+0.11}$	$-{3.82}_{-0.03}^{+0.02}$
DRIFT-PHOENIX (All data)	${1723}_{-4}^{+5}$	${4.00}_{-0.04}^{+0.05}$	${0.14}_{-0.01}^{+0.02}$	⋯	${1.44}_{-0.01}^{+0.01}$	$-{3.76}_{-0.00}^{+0.00}$
DRIFT-PHOENIX (NIRCam only)	${1671}_{-63}^{+60}$	${4.04}_{-0.10}^{+0.09}$	$-{0.20}_{-0.30}^{+0.34}$	⋯	${1.40}_{-0.08}^{+0.09}$	$-{3.84}_{-0.03}^{+0.03}$
Exo-REM (All data)	${1452}_{-3}^{+4}$	${3.76}_{-0.04}^{+0.03}$	${0.80}_{-0.04}^{+0.05}$	${0.37}_{-0.01}^{+0.01}$	${2.03}_{-0.02}^{+0.02}$	$-{3.76}_{-0.01}^{+0.01}$
Exo-REM (NIRCam only)	${1450}_{-62}^{+75}$	${3.79}_{-0.10}^{+0.09}$	${0.34}_{-0.52}^{+0.46}$	${0.53}_{-0.22}^{+0.14}$	${1.88}_{-0.12}^{+0.12}$	$-{3.82}_{-0.04}^{+0.04}$

Note. For each model grid, we separately report the results of fits to all available data and to the NIRCam data only.

Download table as:  ASCIITypeset image

Figure 11 shows the best-fit model atmospheres for β Pic b for all six considered cases. The fits that consider only the new NIRCam photometry favor an object that emits less flux over the shown spectral range and has a smaller bolometric luminosity. Notably, the NACO photometry and the GPI H-band spectrum appear to fit better to the model atmosphere derived from the NIRCam photometry alone. This means that the NACO and GPI H-band data seem to be more consistent with the new NIRCam photometry than with the existing GRAVITY spectrum, which pulls the entire near-IR spectrum of β Pic b up to higher fluxes in the fits that include the literature data. We note that these fits use an equal weight for each data point, so that the GRAVITY spectrum, although correlations between spectral channels are taken into account, dominates the fit. Investigating the origin of the systematic differences between the GRAVITY spectrum and the other instruments is beyond the scope of this paper, but we note in Section 3.1 that variability in the scattered light from the debris disk could play a role.

Zoom In Zoom Out Reset image size

Comparing, for each model, the fits with the NIRCam data alone to the ones including existing literature data, Figure 12 shows that the retrieved posteriors for the effective temperature, the surface gravity, and the radius are mostly consistent with one another within 1σ. As observed frequently in the literature, there are however systematic differences between the different model grids. In general, the retrieved effective temperature ranges from ∼1450–1750 K, which is consistent with previous studies (e.g., Bonnefoy et al. 2014b; GRAVITY Collaboration et al. 2020). However, it becomes evident that at the current photometric precision, our specific set of NIRCam filters is not suitable for constraining the planet's metallicity or C/O abundance ratio. Figure 12 illustrates that the posterior distributions of these quantities remain unconstrained regardless of the chosen model atmosphere grid. When including literature data in the fits, the metallicity and C/O abundance ratio posterior distributions converge sharply. With the Exo-REM grid, we retrieve a C/O abundance ratio of ${0.37}_{-0.01}^{+0.01}$, which is consistent with the value retrieved in GRAVITY Collaboration et al. (2020) with the GRAVITY K-band spectrum alone. The retrieved metallicity is strongly model-dependent. We note that we can obtain much tighter constraints on the metallicity and C/O abundance ratio from our NIRCam observations alone when artificially increasing the photometric precision to ∼1%. When doing so, the uncertainty on the retrieved metallicity drops to ∼0.05 dex and the uncertainty on the retrieved C/O ratio drops to ∼0.1. We suspect that the six photometric points from NIRCam alone can constrain these quantities quite tightly, because they cover such a broad wavelength range (∼1.7–5 μm). However, similar to the above, there remain systematic differences between the different model grids, so that the retrieved parameter values need to be treated with caution.

Zoom In Zoom Out Reset image size

For all three model grids, though, we find that the luminosity retrieved from the NIRCam photometry alone is significantly smaller than the one retrieved from the fit including literature data. This is not completely surprising, given that the new NIRCam photometry falls consistently below the best-fit DRIFT-PHOENIX model atmosphere fitted to literature data shown in Figure 5. Our specific set of NIRCam filters covers both the water bands at ∼2 and ∼2.5 μm and the CO/CO₂ bands in the ∼4–5 μm regime, which are the main carbon- and oxygen-bearing species in the atmosphere of β Pic b. However, we are lacking medium-band observations between 4 and 5 μm, which could provide a better handle on the relative strength of the CO and CO₂ features.

The formation history of β Pic b has been investigated in several studies, and most of them concluded that either its luminosity–age relationship (e.g., Bonnefoy et al. 2014b; Marleau & Cumming 2014; Nowak et al. 2020) or its atmospheric composition (e.g., GRAVITY Collaboration et al. 2020) are most consistent with a moderate formation entropy (≳9.75 k_B baryon^–1), suggesting a formation by warm-start core accretion (Pollack et al. 1996; Spiegel & Burrows 2012; Mordasini 2013) followed by strong planetesimal enrichment in the inner regions of the β Pic system (i.e., within the CO₂ ice line, e.g., Öberg et al. 2011; Eistrup et al. 2018). Based on the new JWST/NIRCam photometry, which improves the constraints on β Pic b's thermal emission especially in the ∼3–5 μm regime, we obtain slightly updated constraints on the planet's bolometric luminosity and radius. We note that while the ∼5–7 μm MIRI/MRS spectrum of β Pic b is included in the fits, it does not contribute significantly to the estimated bolometric luminosity of the planet, because we let the spectrum scale up or down freely in these fits, following the same procedure as in Worthen et al. (2024). Moreover, as in Stolker et al. (2020), we keep the GRAVITY and GPI Y–J-band spectra fixed while also fitting a scale factor for the GPI H-band spectrum, which otherwise does not agree with the GRAVITY spectrum. A similar discrepancy between the GPI H-band spectrum and the GRAVITY and SPHERE spectra of the brown dwarf companion HD 206893 B has been reported by Kammerer et al. (2021).

We compare these updated values to evolutionary model predictions from the New Generation Planetary Population Synthesis (NGPPS) simulation (Emsenhuber et al. 2021) in Figure 13. The NGPPS simulation uses a global end-to-end planetary formation and evolution model based on the core accretion process to predict all common observational quantities of exoplanets and exoplanetary systems. It can be seen that the dynamical mass from Nowak et al. (2020) and Brandt et al. (2021) is in tension with the luminosity measured from the spectral energy distribution (SED) fits, which suggests that β Pic b has a mass closer to the deuterium-burning limit. This has already been noted by Brandt et al. (2021). The same conclusion can be drawn from the inferred planet radius, which (regardless of the used model atmosphere grid) is also in tension with the one inferred from the NGPPS and the dynamical mass of β Pic b. Investigating the reason for this discrepancy is beyond the scope of this paper. However, we note that Lacour et al. (2021) find a more consistent dynamical mass of 11.9 ± 3 M_Jup for β Pic b when only considering the radial velocities and the astrometry, ignoring the Gaia–Hipparcos proper-motion anomaly.

Zoom In Zoom Out Reset image size

Our inferred bolometric luminosity is consistent with the one reported in Morzinski et al. (2015) and Chilcote et al. (2017) to within ±0.05 dex, but the comparison of dynamical mass and NGPPS predicts a bolometric luminosity that is ∼0.5 dex lower. We further note that the radius inferred by the Exo-REM grid is very large, almost unphysical if compared to the NGGPS population. However, it is roughly consistent with the radius of 1.9 R_Jup inferred by GRAVITY Collaboration et al. (2020), also using the Exo-REM grid. With the free retrieval code petitRADTRANS, the same authors obtain a radius of 1.36 ± 0.02 R_Jup, which is much closer to the one that we infer using the DRIFT-PHOENIX model atmosphere grid. The large discrepancies between the inferred radii can originate from the different cloud prescriptions in the different model atmosphere grids and retrieval codes. Clouds can have a similar effect on a planet's SED as a change in effective temperature and radius and are therefore often correlated with these quantities. In the literature, the DRIFT-PHOENIX grid and the cloud prescription therein has been found to provide good fits to young and dusty low-gravity objects (Patience et al. 2012; Bonnefoy et al. 2014a; Lachapelle et al. 2015) such as β Pic b (Stolker et al. 2020), HD 95086 b (De Rosa et al. 2016), and HD 206893 B (Kammerer et al. 2021).