Open Access
Issue
A&A
Volume 711, July 2026
Article Number A87
Number of page(s) 7
Section Planets, planetary systems, and small bodies
DOI https://doi.org/10.1051/0004-6361/202660276
Published online 09 July 2026

© The Authors 2026

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. This email address is being protected from spambots. You need JavaScript enabled to view it. to support open access publication.

1 Introduction

Beyond elemental abundances, isotopic ratios trace the chemical environment of the protoplanetary disc and the relative roles of gas and solids during accretion. A familiar Solar System example is the deuterium-to-hydrogen ratio (D/H): Jupiter and Saturn have D/H consistent with the protosolar nebula, whereas Uranus and Neptune are enriched, which is commonly attributed to accretion of D-rich ices beyond the water snow line (Feuchtgruber et al. 2013). Carbon isotopes provide an analogous diagnostic for exoplanets. In discs, fractionation can shift 12C/13C away from the local interstellar medium (ISM) value of approximately 68 (Milam et al. 2005), for instance through isotopic ion-exchange reactions (Langer & Penzias 1993) and isotope-selective photodissociation of CO by ultraviolet radiation (Visser et al. 2009). These processes can enrich 13C in ices relative to the gas, so a planet's atmospheric 12C/13C depends on where it formed and on the balance between gas and solid carbon accretion (Zhang et al. 2021a; Bergin et al. 2024). Highresolution (ℜ ≳ 30 000) near-infrared spectroscopy can measure 12C/13C in substellar atmospheres. The CO overtone bands in the K band around 2.3 μm include resolved lines from both 12CO and 13CO, allowing retrievals to disentangle their relative abundances (Mollière & Snellen 2019; Zhang et al. 2021a). The first exoplanet measurement targeted the young super-Jupiter YSES 1 b and initially reported 12C/13C = 3110+17Mathematical equation: $^{12}{\rm{C}}{/^{13}}{\rm{C = 31}}_{ - 10}^{ + 17}$ (Zhang et al. 2021a); follow-up work revised this to 12C/13C = 88 ± 13 (Zhang et al. 2024), highlighting the need for high resolution and adequate signal-to-noise ratio (S/N) when constraining isotopic ratios.

Improved observations and retrieval methods have made 12C/13C measurements feasible for an increasing range of substellar objects and cool stars (Crossfield et al. 2019; Zhang et al. 2021b; Gandhi et al. 2023; Xuan et al. 2024a; Grasser et al. 2025; González Picos et al. 2025). Access to the fundamental CO band at 4.6 μm with the JWST has enabled precise measurements of 12C/13C in a number of directly-imaged planets and brown dwarfs, including three of the four HR 8799 planets (Ruffio et al. 2026), providing a direct view of how 12C/13C varies with orbital radius within a single system. Several isotopic studies of young and field brown dwarfs have reported 12C/13C values ranging from ISM-like to solar and above (Xuan et al. 2024b; de Regt et al. 2024; Grasser et al. 2025; de Regt et al. 2026). Observations of M dwarfs have also provided valuable context for the interpretation of isotopic ratios, suggesting that the wide spread of 12C/13C in nearby M dwarfs is due to long-term Galactic chemical evolution (González Picos et al. 2025a), which can produce high 12C/13C at old ages (12C/13C ≥ 200 at ≈ 10 Gyr) and low 12C/13C ≈ 60 for young, metal-rich objects in the solar neighbourhood (Romano 2022).

2 Observations

We observed β Pic b with CRIRES+ at the VLT over 11 nights between October 2024 and May 2025, integrating one hour per night at ℜ ≈ 100 000 in the K band. The data were taken under programme 114.27DX.001 (PI: M. Kenworthy; Kenworthy et al. 2026) with a fixed instrument setup across all epochs. Unlike previous CRIRES+ observations (Landman et al. 2024), we placed the star outside the slit (45 deg slit angle) to improve the star-to-planet flux ratio and limit detector noise, allowing 120 s exposures instead of 10 s. Under favourable conditions (low airmass and seeing), the planet-to-star contrast is up to ~4 times higher than in that analysis, which we attribute to the reduced detector noise and improved suppression of scattered starlight. Combined with our novel spectral extraction (see below), we obtain S/N per pixel at native resolution between 2 and 9 across nights, well above the ≲ 1 reported previously (Landman et al. 2024), making this the highest-S/N high-resolution spectrum of β Pic b to date.

We reduced the data with the excalibuhr pipeline (Zhang et al. 2024) (see also de Regt et al. 2024; González Picos et al. 2025b). The reduction included dark and flat corrections, spectral-order tracing, lamp-based dispersion calibration applied to the science frames, nod-subtraction sky removal, and stacking of the frames acquired at each nod position. We then constructed a one-dimensional spectrum at the planet location using a novel spectral-extraction scheme tailored to faint sources on a strong slit background, such as close-in companions. Spatial rows with significant planetary signal were co-added with weights derived from the cross-correlation strength along the slit (Fig. 1). We calculated the cross-correlation map by cross-correlating a planetary template (Teff=1700K, log g=4.0, [M/H]=0.0; Morley et al. 2024) with the data at each row position after subtracting the fitted stellar model. The stellar model was jointly fitted with the planetary spectrum at every slit position using the linear model outlined in Sect. 3.2.

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

Detection of β Pic b. Panel a: calibrated detector image of one echelle order with central wavelength 2330 nm. The scattered starlight is visible as a bright stripe close to the centre of the slit. Panel b: cross-correlation of the data and a planet template after starlight subtraction across the slit. Panel c: cross-correlation function at locations with significant planet signal. Line transparency follows the S/N per row in Panel b. The dashed line shows the cross-correlation on the opposite side of the slit, indicating the typical noise level.

3 Methods

High-resolution spectroscopy resolves individual molecular lines and is sensitive to the pressure-temperature (P–T) structure. K-band observations probe H2O and CO well but constrain clouds poorly: the narrow bandpass and lack of continuum information make it difficult to separate cloud opacity from metallicity and P–T effects. For β Pic b, this yields a strong metallicity–cloud degeneracy in spectroscopic analyses (e.g. Landman et al. 2024). We therefore fit the CRIRES+ spectra jointly with archival 1.0–5.0 μm broadband photometry (Morzinski et al. 2015). The Y, J, and H bands respond to Teff and cloud opacity, whereas the L and M bands are less cloud-sensitive and anchor the bolometric luminosity. Combining the line-resolved spectrum with the photometric spectral energy distribution (SED) tightens constraints on the atmospheric structure and clouds compared with CRIRES+ alone, as shown for other directly-imaged planets whose continuum is lost in reduction (e.g. Ruffio et al. 2026). A limitation is that the photometry is single-epoch and not contemporaneous with the spectroscopy, which could bias the retrieval if the atmosphere is highly variable. β Pic b is nevertheless expected to be relatively stable (≲5%) as it is an early-L object, based on L-dwarf variability studies (Vos et al. 2020, 2022), and cloud properties are not expected to significantly affect carbon isotope measurements (Xuan et al. 2024a; de Regt et al. 2026).

3.1 Forward model

The observed K-band spectrum contains contributions from the planet, the partially scattered stellar point spread function, and telluric absorption by Earth's atmosphere. At the 0.48″ projected separation of ß Pic b, the raw starlight-to-planet flux ratio is approximately 100, so the star overwhelms the planet signal. Rather than subtracting it empirically, we included all three components simultaneously in a forward model, following an approach similar to that of González Picos et al. (2025b) and inspired by Landman et al. (2024).

We generated planetary spectra with petitRADTRANS, using P–T profiles drawn1 from the Sonora Diamondback grid of model atmospheres (Morley et al. 2024). This choice couples the thermal structure to the cloud properties and chemical equilibrium in a physically motivated way. The grid spans effective temperature Teff, surface gravity log g, metallicity [M/H], the sedimentation parameter fsed, which controls how efficiently cloud particles rain out of the atmosphere, and the eddy diffusion coefficient Kzz, which regulates vertical mixing of the cloud material. We obtained the atmospheric composition assuming chemical equilibrium, including condensation and rainout (Kitzmann et al. 2024). Molecular line opacities for 12CO, 13CO, H2O, CH4, and HF were included alongside collision-induced absorption from H2–H2 and H2–He pairs (see González Picos et al. 2025b for linelists and references). Condensate opacity was included using the parametrisation of Mollière et al. 2020, with the cloud particle size set by the eddy diffusion coefficient (Ackerman & Marley 2001) and a log-normal size distribution set by the standard deviation σg. We considered Fe and Mg2SiO4 condensates as possible cloud-forming species (Gao et al. 2020; Vos et al. 2023). We included three independent free parameters for each species: cloud-base pressure Pbase,s, cloud-base mass fraction Xbase,s, and sedimentation parameter fSed,s- The model spectrum was Doppler-shifted, rotationally broadened by direct integration over the planetary surface (Carvalho & Johns-Krull 2023) using the projected rotational velocity υ sin i, and convolved with a Gaussian kernel matching the instrumental resolution.

We extracted the telluric transmission empirically from the data. Because β Pictoris is an early-type star (Crifo et al. 1997) with a spectrum that closely resembles a blackbody in the Κ band, dividing the observed stellar spectrum by a Planck function at the stellar effective temperature yields a high-quality telluric absorption spectrum for each night (similar to Landman et al. 2024). We accounted for residual wavelength-dependent variation in the stellar contamination across the slit with a linear-spline scaling function fitted simultaneously with the stellar model, following González Picos et al. (2025b). This differs from the principal component analysis (PCA) approach used by Landman et al. (2024).

Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

Detection of 13CO in the atmosphere of β Pic b. Panel a: combined observed spectrum from the six nights with S/N > 5 in the 13CO overtone bands, shifted to the planet rest-frame and binned to 10 km s–1 for display. The two models are averages of the individual best-fit models with and without 13CO. Panel b: residuals between the observed spectrum and the model without 13CO, with a 13CO-only model overplotted. Panel c: cross-correlation of the residuals, averaged over the full wavelength range and the six nights, with a 13CO template. The dashed line shows the template auto-correlation function, and the faint grey lines show the individual-night cross-correlations.

3.2 Bayesian retrieval

We sampled the posterior distribution with nested sampling as implemented in PyMultiNest (Buchner 2016). We ran each night independently with 800 live points and a constant efficiency of 5%, following common settings in the literature (Zhang et al. 2021b; de Regt et al. 2024).

The free parameters comprised the atmospheric structure (Teff, log g, [M/H], fsed, log Kzz), cloud properties (cloud-base pressure, cloud-base mass fraction, and opacity from Fe and Mg2SiO4 condensates; Vos et al. 2023; Morley et al. 2024), planetary radius, radial velocity, ν sin i, C/O ratio, and 12C/13C. A noise inflation parameter b accounted for underestimated uncertainties or model–data mismatches. It was defined as σscaled2=σdata2+10bMathematical equation: $\sigma _{{\rm{scaled}}}^2 = \sigma _{{\rm{data}}}^2 + {10^b}$, ensuring correct error propagation into the posteriors (Line et al. 2015; Mollière et al. 2025).

We adopted the likelihood of Landman et al. (2024) (following Ruffio et al. 2019) and evaluated it separately for each night and nod position, as in González Picos et al. (2025b). For a given night and nod, log  =12[ N log(2π)+ log det C+(dm)C1(dm) ],Mathematical equation: ${\rm{log }}{\cal L} {\rm{ = }} - {1 \over 2}\left[ {N{\rm{ log}}\left( {2\pi } \right){\rm{ + log det }}{\bf{C}} + {{\left( {{\bf{d}} - {\bf{m}}} \right)}^}{{\bf{C}}^{ - 1}}\left( {{\bf{d}} - {\bf{m}}} \right)} \right],$(1) m=Mϕ,ϕ=argminϕ0 dMϕ 22,Mathematical equation: $\matrix{ {{\bf{m}} = {\bf{M}}\phi ,} & {\phi = \arg \mathop {\min }\limits_{\phi \ge 0} \left\| {{\bf{d}} - {\bf{M}}\phi } \right\|_2^2,} \cr }$(2)

where d is the data vector and C is the covariance built from the scaled uncertainties. For each order–detector spectrum, the design matrix M contains the planet atmospheric spectrum with telluric transmission (m0) and nstar shifted spline components of the on-axis stellar spectrum (González Picos et al. 2025b). We solve for ϕ at each likelihood evaluation with scipy.optimize.nnls (Lawson & Hanson 1995), so stellar contamination and low-frequency systematics are constrained without extra non-linear parameters. The planet amplitude ϕ0 gives the wavelength-dependent planet-to-star flux ratio, which we use to track night-to-night variations driven by observing conditions and adaptive optics (AO) performance.

4 Results

We present the combined spectrum and average best-fit model in Fig. 2, together with the cross-correlation detection of 13CO. The data quality varies substantially between nights (S/N per pixel between 2 and 9 after error scaling), which we attribute to changing observing conditions. Individual-night constraints on Teff, [M/H], and C/O are shown in Fig. 3, and those on 12C/13C are shown in Fig. 4. Nights with S/N > 5 agree within the uncertainties, while lower-S/N nights show systematic shifts. We therefore combined the posterior distributions of the six best nights (S/N > 5) with equal weights to obtain the final constraints and incorporate the night-to-night scatter into the uncertainty budget. Alternative weighting schemes, such as Bayesian model averaging (Nixon et al. 2024) or S/N-weighted combination, assign negligible weight to all but the highest-quality nights. We therefore adopted equal weights and set a threshold of S/N > 5 to retain the scatter as a conservative estimate of systematic uncertainties for nights with comparable S/N.

Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Retrieved atmospheric parameters. We show each individual night and the combined posterior of the six nights with S/N>5. The Teff converges to a value consistent with the spectral type of β Pic b. The metallicity shows scatter and a bimodal tendency, which we attribute to the cloud–metallicity degeneracy. The C/O ratio is less affected by this degeneracy but still shows substantial epoch-to-epoch scatter.

Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Measurement of 12C/13C in the atmosphere of β Pic b. Top: combined posterior from the six nights with S/N per pixel > 5. Bottom: per-night posteriors. The S/N uncertainties are estimated from the difference between the median S/N in the two nod positions. Colours encode the nested-sampling evidence of each run relative to the lowest-evidence case; higher values indicate stronger support.

4.1 Atmospheric parameters

4.1.1 Temperature structure

From the joint posterior, we obtain Teff=162928+30KMathematical equation: ${T_{{\rm{eff}}}} = 1629_{ - 28}^{ + 30}{\rm{K}}$. This is in line with the spectral type of ß Pic b and with earlier SED and low-resolution spectroscopic analyses (Bonnefoy et al. 2013; Baudino et al. 2015; Worthen et al. 2024; Ravet et al. 2025). The retrieved pressure-temperature profiles systematically reach the lower end of the grid at fsed ≤ 1 (Fig. A.1), so the best-fitting P–T structures match models that include thick clouds. For nights with S/N > 5, the retrieved Teff values agree well with the adopted value; lower-S/N epochs show larger scatter but remain consistent with it within 1σ. Our Teff is inferred directly from the P–T profile fit. It is therefore model-dependent and sensitive to the cloud-opacity assumptions built into the Sonora Diamond-back grid (Morley et al. 2024). Estimates of Teff from a full spectral energy distribution should be more reliable, because they provide an empirical reading of the bolometric flux.

4.1.2 Chemical composition

We recover [ M/H ]=0.200.12+0.16Mathematical equation: $\left[ {{\rm{M/H}}} \right] = 0.20_{ - 0.12}^{ + 0.16}$ and C/O = 0.52 ± 0.03. Both are compatible with recent work that favours roughly solar to mildly super-solar metallicity and near-solar C/O (Ravet et al. 2025), but they disagree with the metallicity and chemistry inferred from GRAVITY low-resolution spectroscopy (GRAVITY Collaboration 2020). Our C/O exceeds the value reported from CRIRES+ by Landman et al. (2024) (C/O = 0.41 ± 0.04). Their limited constraints on metallicity, with subsolar [M/H] in all models except the one that adopted the GRAVITY P–T structure, likely propagate into C/O, which cautions against a direct comparison. The C/O ratios from individual nights approach the solar value (≈0.59; Asplund et al. 2021) as the S/N increases (right panel of Fig. 3). Recent C/O measurements for young brown dwarfs in the β Pictoris moving group (Liu et al. 2026) suggest near-solar compositions, consistent with our findings and with other young brown dwarfs (González Picos et al. 2024; Grasser et al. 2025).

4.1.3 Cloud properties

The spectra show tentative evidence for condensate opacity: the retrieval allows non-zero contributions from Fe and Mg2SiO4 clouds, with posteriors that differ between epochs (Fig. A.1). A tentative SiO detection in the M band by Parker et al. 2024 suggests that not all silicon is locked in silicates, which would favour Fe as the dominant cloud-forming species. Unlike Landman et al. (2024), who found cloud parameters to be largely unconstrained and argued for cloud-free solutions because their P–T structures did not cross any of the considered condensation curves, we treated the cloud-base pressure and mass fraction of each species as free parameters. We therefore did not tie condensation strictly to optical constants and the intersection of the P–T profile. Several high-S/N nights favour cloud opacity near the bottom of the photosphere (Pbase ≈ 1 bar) with elevated base mass fractions (log Xbase ≈ −4 to −1), whereas lower-S/N nights yield weak or absent constraints. A full treatment of cloud composition is left to future work, which would benefit from a broader wavelength range, in particular at bluer wavelengths where the spectrum is more sensitive to cloud properties (e.g. Zhang et al. 2023).

4.1.4 Spin

The retrieved υ sin i values are in the range of 18-19 km s–1 and agree broadly with recent high-resolution estimates (19.9 ± 1.0 km s−1; Landman et al. 2024), but they are systematically lower than some earlier values (e.g. Snellen et al. 2014). We attribute the difference mainly to the rotational broadening implementation, which integrates the visible hemisphere and includes angle-dependent intensities rather than a singletemperature slab (Mollière et al. 2019). A convolution kernel with a limb-darkening correction may overestimate υ sin i. Applying the rotational broadening to intensity spectra at different incident angles gives a υ sin i value that naturally includes limb darkening. Small differences in ν sin i between nights may reflect surface heterogeneity, which can alter the line profiles, but we do not explore this further in the present work.

4.2 Carbon isotope ratio

From the combined posterior, we measure 12C/13C = 5815+18Mathematical equation: $^{12}{\rm{C}}{/^{13}}{\rm{C = 58}}_{ - 15}^{ + 18}$ in the atmosphere of ß Pic b (Fig. 4). This is consistent with the present-day local ISM (12C/13C = 68 ± 14; Milam et al. 2005) and below the solar photospheric value (12C/13C = 91.4 ± 1.3; Lyons et al. 2018). An independent GRAVITY+ analysis (von Stauffenberg et al. 2026; companion paper) finds a somewhat higher ratio, around the solar value and compatible with our measurement within 1σ. We emphasise the importance of multi-epoch validation of 12C/13C to ensure that the measurement is not biased by low-S/N data or model-data mismatches. Measurements from individual nights show discrepancies up to 3σ (see Fig. 4), which can lead to incorrect conclusions about the planet's carbon isotopic composition if assessed in isolation. This is the case for two of our six nights with S/N > 5, which may be affected by systematic errors in the retrieval or data reduction (see Fig. A.1). The remaining four higher-S/N epochs show good agreement within 1σ. By combining the posteriors from several epochs of comparable quality, we obtain a more robust measurement, with uncertainties that reflect the night-to-night scatter. Our final adopted value is compatible with all individual night measurements within 1σ.

Ravet et al. (2025) tentatively inferred much lower 12C/13C for ß Pic b (down to ≈11) from a cross-correlation detection of 13CO in medium-resolution GRAVITY data ( ~ 4000). Such values would imply strong 13C enrichment and are neither found in our retrieval nor supported by the GRAVITY+ study (von Stauffenberg et al. 2026). Low-S/N detections can bias isotopic ratios towards low 12C/13C, as Zhang et al. (2024) illustrated for YSES 1 b by revising a SINFONI-based estimate of 31 to 88 with higher-S/N, higher-resolution CRIRES+ data. We therefore stress the need to confirm 12C/13C at adequate resolving power and S/N.

An ISM-like 12C/13C indicates that the carbon isotopic composition was not substantially altered during planet formation and evolution. At a separation of approximately 10 au from its A-type host star (Lacour et al. 2021), β Pic b lies interior to the CO snow line (Öberg et al. 2011). If the planet formed near its current location, the available carbon budget was predominantly in the gas phase and was therefore probably not subject to low-temperature fractionation processes that can enrich 13C in ices (Woods & Willacy 2009).

Measurements of carbon isotope ratios in young giant planets so far span 12C/13C ≈ 30–120 (see Fig. 5). The literature shows no clear dependence of 12C/13C on atmospheric properties or orbital distance. Additional measurements with improved precision are needed to assess whether carbon isotope ratios are useful formation tracers. The current sample is small, with around a dozen measurements (see Fig. 5), and the uncertainties are often too large to draw firm conclusions. We expect this sample to grow substantially with next-generation facilities such as the Extremely Large Telescope (ELT)ZMETIS (Brandl et al. 2021), which will access the fundamental CO band around 4.6 μm at high resolution and with much higher S/N.

Future constraints on stellar or disc 12C/13C would help contextualise our atmospheric measurement. High-resolution spectroscopy of the A-type host remains impractical, but (sub)millimetre interferometry can in principle constrain carbon isotopes in circumstellar gas (Cataldi et al. 2018; Yoshida et al. 2022). For β Pictoris, the debris disc is collisionally replenished by exocometary material, so millimetre-line 12C/13C measurements would trace reprocessed solids and released volatiles rather than a pristine nebular reservoir. Isotopic measurements of objects in the same association may therefore provide a useful baseline: a young brown dwarf in the β Pictoris moving group, 2MASSI J0443+0002, has a reported value of 12C/13C = 70 ± 5 (Liu et al. 2026), comparable to the ISM value and consistent with our result within the uncertainties.

Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Carbon isotope ratio as a function of the semi-major axis for young, directly-imaged companions. Measurements in the HR 8799 system for planets e, b, and c (Ruffio et al. 2026) are shown with squares and connected with a dotted line. The other measurements are for HD 984 b (Costes et al. 2024), GQ Lup b (González Picos et al. 2025b), VHS 1256 b (Gandhi et al. 2023), YSES 1 b (Zhang et al. 2024), AB Pic b (Gandhi et al. 2025), and DH Tau b and HIP 79098 b (Xuan et al. 2024a). The present-day local ISM (12C/13C = 68 ± 14; Milam et al. 2005) and solar photospheric value (12C/13C = 91.4 ± 1.3; Lyons et al. 2018) are overlaid. The horizontal line spans the typical CO snowline location for a solar-type host star (22–36 au). For stars hotter than the Sun, the CO snow line is expected at distances greater than 30 au, so β Pic b is located interior to the CO snow line.

5 Conclusion

We presented eleven epochs of K-band CRIRES+ spectroscopy of β Picb at ≈ 100 000, reduced with excalibuhr and extracted with a cross-correlation-weighted scheme suited to faint companions on a starlight-dominated background. A forward model combining the planet, scattered starlight, and telluric absorption, fitted jointly to the high-resolution spectrum and broadband photometry, yields a 13CO detection and 12C/13C = 5815+18Mathematical equation: $^{12}{\rm{C}}{/^{13}}{\rm{C = 58}}_{ - 15}^{ + 18}$, consistent with the local ISM and below the solar photo-spheric value. The retrieval simultaneously constrains the thermal and chemical structure (Teff=162928+30KMathematical equation: ${T_{{\rm{eff}}}} = 1629_{ - 28}^{ + 30}{\rm{K}}$, [ M/H ]=0.20012+0.16Mathematical equation: $\left[ {{\rm{M/H}}} \right] = 0.20_{ - 012}^{ + 0.16}$, C/O = 0 52 ± 0 03, and tentative condensate opacity). These results broadly align with recent work while differing from several earlier CRIRES+ and GRAVITY-based inferences. The added photometry eases metallicity–cloud degeneracies inherent to high-resolution spectra alone and places the metallicity of β Picb in the solar to mildly super-solar range. Observations of younger, gas-rich protoplanetary discs may provide useful context for linking atmospheric 12C/13C to formation chemistry when stellar constraints are unavailable. Extending high-resolution isotopic measurements to additional companions in the same system, notably the inner super-Jupiter β Pic c, will eventually test whether in situ 12C/13C tracks orbital distance together with the population-level benchmarks for young, directly-imaged planets.

Acknowledgements

We thank the anonymous referee for their helpful comments and constructive suggestions, which significantly improved the manuscript. D.G.P. and I.S. acknowledge NWO grant OCENW.M.21.010. Support for this work was provided by the NL-NWO Spinoza (SPI.2022.004). Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere under ESO programme(s) 114.27DX.001. This work used the Dutch national e-infrastructure with the support of the SURF Cooperative using grant no. EINF-4556. JLB acknowledges funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme under grant agreement No 805445. Software: NumPy (Harris et al. 2020), SciPy (Virtanen et al. 2020), Matplotlib (Hunter 2007), petitRADTRANS (Mollière et al. 2019), species (Stolker et al. 2020), fastchem (Kitzmann et al. 2024), PyMultiNest (Buchner 2016), Astropy (Astropy Collaboration 2022), corner (Foreman-Mackey 2016), ExoMol (Tennyson et al. 2024), HITEMP (Rothman et al. 2010), pyROX (de Regt et al. 2025).

References

  1. Ackerman, A. S., & Marley, M. S. 2001, ApJ, 556, 872 [Google Scholar]
  2. Asplund, M., Amarsi, A. M., & Grevesse, N. 2021, A&A, 653, A141 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Astropy Collaboration (Price-Whelan, A. M., et al.) 2022, ApJ, 935, 167 [NASA ADS] [CrossRef] [Google Scholar]
  4. Baudino, J.-L., Bézard, B., Boccaletti, A., et al. 2015, A&A, 582, A83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Bergin, E. A., Bosman, A., Teague, R., et al. 2024, ApJ, 965, 147 [NASA ADS] [CrossRef] [Google Scholar]
  6. Bonnefoy, M., Boccaletti, A., Lagrange, A.-M., et al. 2013, A&A, 555, A107 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Brandl, B., Bettonvil, F., van Boekel, R., et al. 2021, The Messenger, 182, 22 [NASA ADS] [Google Scholar]
  8. Buchner, J. 2016, PyMultiNest: Python interface for MultiNest, Astrophysics Source Code Library [record ascl:1606.005] [Google Scholar]
  9. Carvalho, A., & Johns-Krull, C. M. 2023, RNAAS, 7, 91 [NASA ADS] [Google Scholar]
  10. Cataldi, G., Brandeker, A., Wu, Y., et al. 2018, ApJ, 861, 72 [NASA ADS] [CrossRef] [Google Scholar]
  11. Costes, J. C., Xuan, J. W., Vigan, A., et al. 2024, A&A, 686, A294 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Crifo, F., Vidal-Madjar, A., Lallement, R., Ferlet, R., & Gerbaldi, M. 1997, A&A, 320, L29 [NASA ADS] [Google Scholar]
  13. Crossfield, I. J. M., Lothringer, J. D., Flores, B., et al. 2019, ApJ, 871, L3 [NASA ADS] [CrossRef] [Google Scholar]
  14. de Regt, S., Gandhi, S., Snellen, I. A. G., et al. 2024, A&A, 688, A116 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. de Regt, S., Gandhi, S., Siebenaler, L., & González Picos, D. 2025, arXiv e-prints [arXiv:2510.20870] [Google Scholar]
  16. de Regt, S., Snellen, I. A. G., González Picos, D., et al. 2026, A&A, 707, A210 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Feuchtgruber, H., Lellouch, E., Orton, G., et al. 2013, A&A, 551, A126 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Foreman-Mackey, D. 2016, J. Open Source Softw., 1, 24 [Google Scholar]
  19. Gandhi, S., de Regt, S., Snellen, I., et al. 2023, ApJ, 957, L36 [NASA ADS] [CrossRef] [Google Scholar]
  20. Gandhi, S., de Regt, S., Snellen, I., et al. 2025, MNRAS, 537, 134 [Google Scholar]
  21. Gao, P., Thorngren, D. P., Lee, E. K. H., et al. 2020, Nat. Astron., 4, 951 [NASA ADS] [CrossRef] [Google Scholar]
  22. González Picos, D., Snellen, I. A. G., de Regt, S., et al. 2024, A&A, 689, A212 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. González Picos, D., Snellen, I., & de Regt, S. 2025a, Nat. Astron., 9, 1692 [Google Scholar]
  24. González Picos, D., Snellen, I. A. G., de Regt, S., et al. 2025b, A&A, 693, A298 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Grasser, N., Snellen, I. A. G., de Regt, S., et al. 2025, A&A, 698, A252 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. GRAVITY Collaboration (Nowak, M., et al.) 2020, A&A, 633, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  28. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  29. Kenworthy, M. A., Landman, R., Vanderburg, A., et al. 2026, MNRAS, 549, stag1060 [Google Scholar]
  30. Kitzmann, D., Stock, J. W., & Patzer, A. B. C. 2024, MNRAS, 527, 7263 [Google Scholar]
  31. Lacour, S., Wang, J. J., Rodet, L., et al. 2021, A&A, 654, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Landman, R., Stolker, T., Snellen, I. A. G., et al. 2024, A&A, 682, A48 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Langer, W. D., & Penzias, A. A. 1993, ApJ, 408, 539 [NASA ADS] [CrossRef] [Google Scholar]
  34. Lawson, C. L., & Hanson, R. J. 1995, Solving Least Squares Problems (SIAM) [Google Scholar]
  35. Line, M. R., Teske, J., Burningham, B., Fortney, J. J., & Marley, M. S. 2015, ApJ, 807, 183 [NASA ADS] [CrossRef] [Google Scholar]
  36. Liu, Y., Zhang, Y., Xuan, J. W., et al. 2026, ApJ, 1004, 172 [Google Scholar]
  37. Lyons, J. R., Gharib-Nezhad, E., & Ayres, T. R. 2018, Nat. Commun., 9, 908 [Google Scholar]
  38. Milam, S. N., Savage, C., Brewster, M. A., Ziurys, L. M., & Wyckoff, S. 2005, ApJ, 634, 1126 [Google Scholar]
  39. Mollière, P., & Snellen, I. A. G. 2019, A&A, 622, A139 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Mollière, P., Wardenier, J. P., van Boekel, R., et al. 2019, A&A, 627, A67 [Google Scholar]
  41. Mollière, P., Stolker, T., Lacour, S., et al. 2020, A&A, 640, A131 [Google Scholar]
  42. Mollière, P., Kühnle, H., Matthews, E. C., et al. 2025, A&A, 703, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Morley, C. V., Mukherjee, S., Marley, M. S., et al. 2024, ApJ, 975, 59 [NASA ADS] [CrossRef] [Google Scholar]
  44. Morzinski, K. M., Males, J. R., Skemer, A. J., et al. 2015, ApJ, 815, 108 [NASA ADS] [CrossRef] [Google Scholar]
  45. Nixon, M. C., Welbanks, L., McGill, P., & Kempton, E. M.-R. 2024, ApJ, 966, 156 [Google Scholar]
  46. Öberg, K. I., Murray-Clay, R., & Bergin, E. A. 2011, ApJ, 743, L16 [Google Scholar]
  47. Parker, L. T., Birkby, J. L., Landman, R., et al. 2024, MNRAS, 531, 2356 [Google Scholar]
  48. Ravet, M., Bonnefoy, M., Chauvin, G., et al. 2025, A&A, 704, A325 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Romano, D. 2022, A&A Rev., 30, 7 [Google Scholar]
  50. Rothman, L. S., Gordon, I. E., Barber, R. J., et al. 2010, J. Quant. Spec. Radiat. Transf., 111, 2139 [NASA ADS] [CrossRef] [Google Scholar]
  51. Ruffio, J.-B., Macintosh, B., Konopacky, Q. M., et al. 2019, AJ, 158, 200 [Google Scholar]
  52. Ruffio, J.-B., Xuan, J. W., Chachan, Y., et al. 2026, Nat. Astron., 10, 511 [Google Scholar]
  53. Snellen, I. A. G., Brandl, B. R., de Kok, R. J., et al. 2014, Nature, 509, 63 [Google Scholar]
  54. Stolker, T., Quanz, S. P., Todorov, K. O., et al. 2020, A&A, 635, A182 [EDP Sciences] [Google Scholar]
  55. Tennyson, J., Yurchenko, S. N., Zhang, J., et al. 2024, J. Quant. Spec. Radiat. Transf., 326, 109083 [NASA ADS] [CrossRef] [Google Scholar]
  56. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Med., 17, 261 [Google Scholar]
  57. Visser, R., van Dishoeck, E. F., & Black, J. H. 2009, A&A, 503, 323 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  58. von Stauffenberg, A., Sauter, J., Mollière, P., et al. 2026, A&A, 711, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Vos, J. M., Biller, B. A., Allers, K. N., et al. 2020, AJ, 160, 38 [NASA ADS] [CrossRef] [Google Scholar]
  60. Vos, J. M., Faherty, J. K., Gagné, J., et al. 2022, ApJ, 924, 68 [NASA ADS] [CrossRef] [Google Scholar]
  61. Vos, J. M., Burningham, B., Faherty, J. K., et al. 2023, ApJ, 944, 138 [NASA ADS] [CrossRef] [Google Scholar]
  62. Woods, P. M., & Willacy, K. 2009, ApJ, 693, 1360 [Google Scholar]
  63. Worthen, K., Chen, C. H., Law, D. R., et al. 2024, ApJ, 964, 168 [NASA ADS] [CrossRef] [Google Scholar]
  64. Xuan, J. W., Hsu, C.-C., Finnerty, L., et al. 2024a, ApJ, 970, 71 [CrossRef] [Google Scholar]
  65. Xuan, J. W., Wang, J., Finnerty, L., et al. 2024b, ApJ, 962, 10 [NASA ADS] [CrossRef] [Google Scholar]
  66. Yoshida, T. C., Nomura, H., Furuya, K., Tsukagoshi, T., & Lee, S. 2022, ApJ, 932, 126 [NASA ADS] [CrossRef] [Google Scholar]
  67. Zhang, Y., Snellen, I. A. G., Bohn, A. J., et al. 2021a, Nature, 595, 370 [NASA ADS] [CrossRef] [Google Scholar]
  68. Zhang, Y., Snellen, I. A. G., & Mollière, P. 2021b, A&A, 656, A76 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  69. Zhang, Z., Mollière, P., Hawkins, K., et al. 2023, AJ, 166, 198 [NASA ADS] [CrossRef] [Google Scholar]
  70. Zhang, Y., González Picos, D., de Regt, S., et al. 2024, AJ, 168, 246 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A Posterior distributions

Thumbnail: Fig. A.1 Refer to the following caption and surrounding text. Fig. A.1

Posterior distributions for three representative nights: January 26 (S/N≈3), February 26 (S/N≈7), and March 25 (S/N≈8), 2025. The two higher-S/N nights were deliberately chosen from the S/N > 5 subset to illustrate that even high-S/N single-epoch retrievals can show significant discrepancies in some atmospheric parameters, underscoring the precision limitations of single-epoch analysis with moderate-S/N data. Additional scatter may arise from degeneracies between metallicity, surface gravity, and cloud properties. The sedimentation parameter fsed used for interpolating the Sonora Diamondback P–T profiles consistently reaches the grid lower bound (fsed ≤ 1), corresponding to thick clouds. Cloud-opacity sedimentation is fitted separately for each cloud species, but remains weakly constrained and varies between nights. The upper panel shows the 1, 2, and 3-σ temperature envelopes for each night, with Fe and Mg2SiO4 condensation curves overlaid (dotted). Shaded regions show the integrated emission contribution functions, indicating the vertical extent of the photosphere.

All Figures

Thumbnail: Fig. 1 Refer to the following caption and surrounding text. Fig. 1

Detection of β Pic b. Panel a: calibrated detector image of one echelle order with central wavelength 2330 nm. The scattered starlight is visible as a bright stripe close to the centre of the slit. Panel b: cross-correlation of the data and a planet template after starlight subtraction across the slit. Panel c: cross-correlation function at locations with significant planet signal. Line transparency follows the S/N per row in Panel b. The dashed line shows the cross-correlation on the opposite side of the slit, indicating the typical noise level.

In the text
Thumbnail: Fig. 2 Refer to the following caption and surrounding text. Fig. 2

Detection of 13CO in the atmosphere of β Pic b. Panel a: combined observed spectrum from the six nights with S/N > 5 in the 13CO overtone bands, shifted to the planet rest-frame and binned to 10 km s–1 for display. The two models are averages of the individual best-fit models with and without 13CO. Panel b: residuals between the observed spectrum and the model without 13CO, with a 13CO-only model overplotted. Panel c: cross-correlation of the residuals, averaged over the full wavelength range and the six nights, with a 13CO template. The dashed line shows the template auto-correlation function, and the faint grey lines show the individual-night cross-correlations.

In the text
Thumbnail: Fig. 3 Refer to the following caption and surrounding text. Fig. 3

Retrieved atmospheric parameters. We show each individual night and the combined posterior of the six nights with S/N>5. The Teff converges to a value consistent with the spectral type of β Pic b. The metallicity shows scatter and a bimodal tendency, which we attribute to the cloud–metallicity degeneracy. The C/O ratio is less affected by this degeneracy but still shows substantial epoch-to-epoch scatter.

In the text
Thumbnail: Fig. 4 Refer to the following caption and surrounding text. Fig. 4

Measurement of 12C/13C in the atmosphere of β Pic b. Top: combined posterior from the six nights with S/N per pixel > 5. Bottom: per-night posteriors. The S/N uncertainties are estimated from the difference between the median S/N in the two nod positions. Colours encode the nested-sampling evidence of each run relative to the lowest-evidence case; higher values indicate stronger support.

In the text
Thumbnail: Fig. 5 Refer to the following caption and surrounding text. Fig. 5

Carbon isotope ratio as a function of the semi-major axis for young, directly-imaged companions. Measurements in the HR 8799 system for planets e, b, and c (Ruffio et al. 2026) are shown with squares and connected with a dotted line. The other measurements are for HD 984 b (Costes et al. 2024), GQ Lup b (González Picos et al. 2025b), VHS 1256 b (Gandhi et al. 2023), YSES 1 b (Zhang et al. 2024), AB Pic b (Gandhi et al. 2025), and DH Tau b and HIP 79098 b (Xuan et al. 2024a). The present-day local ISM (12C/13C = 68 ± 14; Milam et al. 2005) and solar photospheric value (12C/13C = 91.4 ± 1.3; Lyons et al. 2018) are overlaid. The horizontal line spans the typical CO snowline location for a solar-type host star (22–36 au). For stars hotter than the Sun, the CO snow line is expected at distances greater than 30 au, so β Pic b is located interior to the CO snow line.

In the text
Thumbnail: Fig. A.1 Refer to the following caption and surrounding text. Fig. A.1

Posterior distributions for three representative nights: January 26 (S/N≈3), February 26 (S/N≈7), and March 25 (S/N≈8), 2025. The two higher-S/N nights were deliberately chosen from the S/N > 5 subset to illustrate that even high-S/N single-epoch retrievals can show significant discrepancies in some atmospheric parameters, underscoring the precision limitations of single-epoch analysis with moderate-S/N data. Additional scatter may arise from degeneracies between metallicity, surface gravity, and cloud properties. The sedimentation parameter fsed used for interpolating the Sonora Diamondback P–T profiles consistently reaches the grid lower bound (fsed ≤ 1), corresponding to thick clouds. Cloud-opacity sedimentation is fitted separately for each cloud species, but remains weakly constrained and varies between nights. The upper panel shows the 1, 2, and 3-σ temperature envelopes for each night, with Fe and Mg2SiO4 condensation curves overlaid (dotted). Shaded regions show the integrated emission contribution functions, indicating the vertical extent of the photosphere.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.