No thick carbon dioxide atmosphere on the rocky exoplanet TRAPPIST-1 c

Seven rocky planets orbit the nearby dwarf star TRAPPIST-1, providing a unique opportunity to search for atmospheres on small planets outside the Solar System1. Thanks to the recent launch of the James Webb Space Telescope (JWST), possible atmospheric constituents such as carbon dioxide (CO2) are now detectable2,3. Recent JWST observations of the innermost planet TRAPPIST-1 b showed that it is most probably a bare rock without any CO2 in its atmosphere4. Here we report the detection of thermal emission from the dayside of TRAPPIST-1 c with the Mid-Infrared Instrument (MIRI) on JWST at 15 µm. We measure a planet-to-star flux ratio of fp/f⁎ = 421 ± 94 parts per million (ppm), which corresponds to an inferred dayside brightness temperature of 380 ± 31 K. This high dayside temperature disfavours a thick, CO2-rich atmosphere on the planet. The data rule out cloud-free O2/CO2 mixtures with surface pressures ranging from 10 bar (with 10 ppm CO2) to 0.1 bar (pure CO2). A Venus-analogue atmosphere with sulfuric acid clouds is also disfavoured at 2.6σ confidence. Thinner atmospheres or bare-rock surfaces are consistent with our measured planet-to-star flux ratio. The absence of a thick, CO2-rich atmosphere on TRAPPIST-1 c suggests a relatively volatile-poor formation history, with less than \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${9.5}_{-2.3}^{+7.5}$$\end{document}9.5−2.3+7.5 Earth oceans of water. If all planets in the system formed in the same way, this would indicate a limited reservoir of volatiles for the potentially habitable planets in the system.

*Corresponding author(s).E-mail(s): zieba@mpia.de; Seven rocky planets orbit the nearby dwarf star TRAPPIST-1, providing a unique opportunity to search for atmospheres on small planets outside the Solar System [1].Thanks to the recent launch of the James Webb Space Telescope (JWST), possible atmospheric constituents such as carbon dioxide (CO 2 ) are now detectable [2,3].Recent JWST observations of the innermost planet TRAPPIST-1 b showed that it is most probably a bare rock without any CO 2 in its atmosphere [4].Here we report the detection of thermal emission from the dayside of TRAPPIST-1 c with the Mid-Infrared Instrument (MIRI) on JWST at 15 µm.We measure a planet-to-star flux ratio of f p /f * = 421 ± 94 parts per million (ppm) which corresponds to an inferred dayside brightness temperature of 380 ± 31 K.This high dayside temperature disfavours a thick, CO 2 -rich atmosphere on the planet.The data rule out cloud-free O 2 /CO 2 mixtures with surface pressures ranging from 10 bar (with 10 ppm CO 2 ) to 0.1 bar (pure CO 2 ).A Venus-analogue atmosphere with sulfuric acid clouds is also disfavoured at 2.6σ confidence.Thinner atmospheres or bare-rock surfaces are consistent with our measured planet-tostar flux ratio.The absence of a thick, CO 2 -rich atmosphere on TRAPPIST-1 c suggests a relatively volatile-poor formation history, with less than 9.5 +7. 5  −2.3 Earth oceans of water.If all planets in the system formed in the same way, this would indicate a limited reservoir of volatiles for the potentially habitable planets in the system.
Little is known about the compositions of terrestrial exoplanet atmospheres, or even whether atmospheres are present at all.The atmospheric composition depends on many unknown factors, including the initial inventory of volatiles, outgassing resulting from volcanism, and possible atmospheric escape and collapse [see e.g. , 5].Atmospheric escape may also depend on the spectral type of the host star: planets around M dwarfs may be particularly vulnerable to atmospheric loss during the long pre-main sequence phase [6].The only way to robustly determine whether a terrestrial exoplanet has an atmosphere is to study it directly, through its thermal emission, reflected light, or transmission spectrum.The tightest constraints on atmospheric properties so far have come from observations of the thermal emission of LHS 3844 b, GJ 1252 b, and TRAPPIST-1 b.The measurements revealed dayside temperatures consistent with no redistribution of heat on the planet and no atmospheric absorption from carbon dioxide [4,7,8].These results motivate observations of cooler planets, which may be more likely to retain atmospheres.
We observed four eclipses of TRAPPIST-1 c with MIRI on JWST in imaging mode.The observations took place on 27 October, 30 October, 6 November, and 30 November 2022 as part of General Observer programme 2304.Each visit had a duration of approximately 192 minutes, covering the 42-minute eclipse duration of TRAPPIST-1 c as well as out-of-eclipse baseline to correct for instrumental systematic noise.The observations used the MIRI F1500W filter, a 3 µm-wide bandpass centred at 15 µm, which covers a strong absorption feature from CO 2 .Across the four visits, we collected 1,190 integrations in total using the FULL subarray.See Methods for further details on the design of the observations.We performed four independent reductions of the data using the publicly available Eureka! code [9] as well as several custom software pipelines.Each reduction extracted the light curve of TRAPPIST-1 using aperture photometry (see Methods and Table 2).We then fitted the light curves with an eclipse model and a range of different parameterizations for the instrumental systematics, including a polynomial in time, exponential ramps, and decorrelation against the position and width of the point spread function (PSF).For the different analyses, the scatter of the residuals in the fitted light curves had a root mean square (rms) variability ranging from 938 -1,079 ppm, within 1.06 -1.22 times the predicted photon noise limit when using a corrected gain value [10].We estimated the eclipse depths using Markov chain Monte Carlo The phase-folded secondary eclipse light curve of TRAPPIST-1 c, measured with the JWST/MIRI imager at 15 µm.The eclipse is centred at orbital phase 0.5 and has a measured depth of fp/f * = 421 ± 94 ppm.The light curve includes four visits (that is, four eclipses), each spanning approximately 3.2 hours.To make the eclipse more easily visible, we binned the individual integrations (grey points) into 28 orbital phase bins (black points with 1σ error bars).The light curve was normalized and divided by the best-fit instrument systematic model.The best-fit eclipse model is shown with the solid red line.The data and fit presented in this figure are based on the SZ reduction, one of the four independent reductions we performed in this work.
(MCMC) fits to the data, which marginalized over all the free parameters in the analysis.The resulting eclipse depths from the four data analyses are consistent and agree to well within 1σ (see Table 4).The phase-folded light curve from one of the reductions can be seen in Figure 1.To determine the final eclipse depth, we took the mean value and uncertainty from the different reductions.To account for systematic error owing to differences in data reduction and modelling choices, we also added an additional 6 ppm to the uncertainty in quadrature, which corresponds to the standard deviation in the eclipse depth between the four analyses.The resulting eclipse depth is f p /f * = 421 ± 94 ppm.
From the measured eclipse depth, we derive a brightness temperature of 380 ± 31 K for TRAPPIST-1 c.The innermost planet in the system, TRAPPIST-1 b, was found to have a brightness temperature of 503 +26 −27 K [4].Compared with previous detections of thermal emission from small (R p < 2 R ⊕ ) rocky planets (see Fig. 5) these temperatures are more than 500 K cooler (the previous lowest measured brightness temperature was 1, 040 ± 40 K for LHS 3844 b [7]).TRAPPIST-1 c is the first exoplanet with measured thermal emission that is comparable with the inner planets of the Solar System; Mercury and Venus have equilibrium temperatures of 440 K and 227 K, respectively, assuming uniform heat redistribution and taking the measured Bond albedo values (A B,Mecury = 0.068, A B,Venus = 0.76) from [11,12].Our measured temperature for TRAPPIST-1 c is intermediate between the two limiting cases for the atmospheric circulation for a zero-albedo planet: zero heat redistribution (430 K; expected for a fully absorptive bare rock), versus global heat redistribution (340 K; expected for a thick atmosphere).This intermediate value hints at either a moderate amount of heat redistribution by an atmosphere (ε = 0.66 +0.26 −0.33 ) or a non-zero Bond albedo for a rocky surface (A B = 0.57 +0.12 −0.15 ) (following the parameterization described in [13]).
To further explore which possible atmospheres are consistent with the data, we compared the dayside flux with a grid of cloud-free, O 2 -dominated models with a range of surface pressures (0.01 bar -100.0 bar) and CO 2 contents (1 ppm -10,000 ppm).Also, we generated cloud-free, pure CO 2 atmospheres using the same surface pressures.The models account for both atmospheric heat redistribution and absorption by constituent gasses [2,7,14] and assume a Bond albedo of 0.1 (see Methods).O 2 /CO 2 mixtures are expected for hot rocky planets orbiting late M-type stars as the planet's H 2 O photodissociates and escapes over time, leaving a desiccated atmosphere dominated by O 2 [6,15,16].Substantial CO 2 (up to about 100 bar) is expected to accumulate from outgassing and does not escape as easily as H 2 O [17,18].For these mixtures, the predicted eclipse depth decreases with increasing surface pressure and with increasing CO 2 abundance, owing to the strong CO 2 absorption feature centred at 15 µm.Strong inversions for a planet in this parameter space are not expected [19].With our measured eclipse depth, we rule out all thick atmospheres with surface pressures P surf ≥ 100 bar (see Fig. 2).For the conservative assumption that the CO 2 content is at least 10 ppm, we rule out P surf ≥ 10 bar.For cloud-free, pure CO 2 atmospheres we can rule out surface pressures P surf ≥ 0.1 bar.As the TRAPPIST-1 planets have precisely measured densities, interior-structure models can give constraints on the atmospheric surface pressures, that is, higher surface pressures would decrease the observed bulk density of the planet.Our findings here agree with these models, which put an upper limit of 160 bar (80 bar) on the surface pressure at a 3σ (1σ) level [20].
We also compared the measured dayside brightness with several physically motivated forward models inspired by Venus.The insolation of TRAPPIST-1 c is just 8% greater than that of Venus [21], so it is possible that the two planets could have similar atmospheric chemistry.We used a coupled   climate-photochemistry model to simulate an exact Venus-analogue composition (96.5%, CO 2 3.5% N 2 , and Venus lower atmospheric trace gases), both with and without H 2 SO 4 aerosols [3] (see Methods).The assumed surface pressure was 10 bar, which would produce similar results to a true 93 bar Venus-analogue, because for both cases, the emitting layer and cloud deck lie at similar pressures.We find that these cloudy and cloud-free Venus-like atmospheres are disfavoured at 2.6σ and 3.0σ, respectively (see Fig. 3 for the 10 bar cloudy Venus spectrum).The cloudy case is marginally more consistent with the data because the SO 2 aerosols locally warm the atmosphere, providing a warmer emission temperature within the core of the 15 µm band, and therefore a larger secondary eclipse depth.
Finally, we compared the measured flux with bare-rock models with a variety of surface compositions, including basaltic, feldspathic, Fe-oxidized (50% nanophase haematite, 50% basalt), granitoid, metal-rich (FeS 2 ), and ultramafic compositions [22].We also considered space weathering for these models, as TRAPPIST-1 c should have been substantially weathered owing  The CO 2 feature overlaps directly with the MIRI F1500W filter used for these observations.The two limiting cases for the atmospheric circulation for a zero-albedo planet (zero heat redistribution, that is, instant reradiation of incoming flux and global heat redistribution) are marked with dashed lines.Two cloud-free, O 2 /CO 2 mixture atmospheres are shown with purple and red solid lines.They show decreased emission at 15 µm owing to CO 2 absorption.A bare-rock model assuming an unweathered ultramafic surface of the planet with a Bond albedo of 0.5 is shown by the solid black line (see text for more information on weathering, including a full comparison of our measurement to a suite of surfaces in Fig. 9).The cloudy Venus forward model with a surface pressure of 10 bars is shown with a solid yellow line.
to its proximity to the host star.On the Moon and Mercury, space weathering darkens the surface by means of the formation of iron nanoparticles [23].On TRAPPIST-1 c, this process would similarly darken the surface and therefore increase the eclipse depth.We find that all bare-rock surfaces are consistent with the data (see Fig. 3 for an unweathered ultramafic surface and Fig. 9 for all surfaces that we considered).Overall, fresh low-albedo surfaces (for example basalt) or weathered surfaces are all compatible with the data, comparable with the likely bare-rock exoplanet LHS 3844 b [7].The highest albedo models, unweathered feldspathic and granitoid surfaces, are a marginally worse fit (consistent at the 2σ level).
To put our results into context with the formation history of the planet, we ran a grid of atmospheric evolution models over a range of initial water inventories (0.1 -100 Earth oceans) and extreme ultraviolet (XUV) saturation fractions for the host star (10 −4 -10 −2 ) (see Fig. 4).The model incorporates outgassing, escape of water vapour and oxygen, and reaction of oxygen with the magma ocean [15].For an XUV saturation fraction of 10 −3 being a typical value for a low-mass star [24], we find that the final surface pressure of oxygen could range over several orders of magnitude (0.1 -100 bar), depending on the initial water inventory (see Fig. 4).Our measured eclipse depth disfavours surface pressures at the high end of this range (greater than 100 bar) for conservative CO 2 abundances, implying that TRAPPIST-1 c most likely formed with a relatively low initial water abundance of less than 9.5 +7.5 −2.3 Earth oceans.For higher CO 2 abundances (> 10 ppm), we rule out surface pressures greater than 10 bars, implying that the planet formed with less than 4 +1.3 −0.8 Earth oceans.Our result suggests that rocky planets around M-dwarf stars may form with a smaller volatile inventory or experience more atmospheric loss than their counterparts around Sun-like stars.This finding motivates further study of the other planets in the TRAPPIST-1 system to assess whether a low volatile abundance is a typical outcome, particularly for the planets in the habitable zone.Fig. 4 Final oxygen atmospheric pressure for TRAPPIST-1 c after 7.5 Gyr of energy-limited escape.We explore different initial planetary water abundances and the amount of XUV the planet receives during the star's saturated activity period [6], described as a fraction of its total bolometric luminosity.The vertical lines represent the nominal XUV saturation fraction of log 10 (L XU V /L bol ) = −3.03+0. 23 −0.12 as estimated by [25].We assume an escape efficiency of 0.1.The white numbers are the contour values for the logarithm of the atmospheric pressure in bars.Our upper limit on surface pressure of 10 -100 bars implies an initial water abundance of approximately 4 -10 Earth oceans.

Methods
Fig. 5 Comparison of small exoplanets with measured infrared emission.Following [8], we show the normalized dayside brightness temperature for super-Earths (Rp < 2R ⊕ ) with measured thermal emission, as a function of planet size (a) and maximum equilibrium temperature (b).The brightness temperatures are normalized relative to predictions for a bare rock with zero albedo and zero heat redistribution, Teq,max.The thermal emission of TRAPPIST-1 c has been detected in this work at 15 µm.The other planets are TRAPPIST-1 b (T1b in plot; also at 15 µm) and planets that have been observed with Spitzer's IRAC Channel 2 at 4.5 µm.The uncertainties on the radius for the planets in the TRAPPIST-1 system are smaller than the marker symbol.Error bars show 1σ uncertainties.

JWST MIRI Observations
As part of JWST General Observer (GO) program 2304 (principal investigator (PI): L. Kreidberg) [26], we observed four eclipses of the planet TRAPPIST-1 c (see Table 1).They were taken on 25 October, 27 October, 30 October, and 6 November 2022 with JWST's MIRI instrument using the F1500W filter.The observations used the FULL subarray with FASTR1 readout and 13 groups per integration.Each visit had a duration of approximately 3.2 hours.We did not perform target acquisition for any of the visits because it was not enabled for MIRI imaging observations during cycle 1.However, the blind pointing precision of JWST was perfectly sufficient to place the target well centered on the field of view of the full array (74" x 113").Figure 6 shows one of the integrations with the FULL array.

Data Reduction
We performed four different reductions of the data collected for JWST program GO 2304.The assumptions made by the reductions are listed in Table 2.In the following, we describe the individual reductions.

Data Reduction SZ
For our primary data reduction and data analysis we used the open-source Python package Eureka![9] which is an end-to-end pipeline for time series observations performed with JWST or the Hubble space telescope (HST).We started our reduction with the raw uncalibrated ("uncal") FITS files which we downloaded from the Mikulski Archive for Space Telescopes (MAST) and followed the multi-stage approach of Eureka! to generate a light curve for TRAPPIST-1 c. Eureka! has been previously successfully used to reduce and analyse the first JWST observations of exoplanets [27][28][29][30][31].
Stages 1 and 2 of Eureka! serve as a wrapper of the JWST Calibration Pipeline [32] (version 1.8.2.).Stage 1 converts groups to slopes and applies basic detector-level corrections.We used the default settings for all steps in this stage but determined a custom ramp-jump detection threshold for each visit by minimizing the median absolute deviation (MAD) of the final light curves.This step detects jumps in the up-the-ramp signal for each pixel by looking for outliers in each integration that might be caused by events such as cosmic rays.We determined a best jump detection threshold of 7σ, 6σ, 7σ and 5σ for visits 1, 2, 3, and 4, respectively, compared with the default value of 4σ set in the JWST pipeline.In stage 2, we only skipped the photom step to leave the data in units of DN/s and not convert into absolute fluxes.In stage 3 of Eureka!, we first masked pixels in each visit that were flagged with an "DO NOT USE" data quality entry, indicating bad pixels identified by the JWST pipeline.Next, we determined the centroid position of the star by fitting a 2D Gaussian to the source.JWST remained very stable during our observations of TRAPPIST-1 c and our target stayed well within a 0.01-pixel area (see Table 1 and Fig. 7).We recorded the centroid position in x and y and the width of the 2D Gaussian in x and y over time to be used in the fitting stage.Next, we determined the best target and background apertures by minimizing the rms of the final light curve.We therefore determined a target aperture of 4 pixels and a background annulus from 25 to 41 pixels from the centroid for each visit.
The light curves show a ramp-like trend at the beginning of the observations, which has already been observed in previous JWST MIRI observations and is most likely caused by charge trapping (see, for example, [10]).We decided to remove the first 10 integrations from each visit, corresponding to approximately 6 minutes or 3% of the data per visit, so that we do not have to also model this initial ramp.Finally, we checked for significant outliers in the final light curves by performing an iterative 5σ outlier clipping procedure.However,

Data Reduction ED
For the second data reduction, we also used the Eureka!pipeline [9] for stages 1 to 5. We also started from the uncal.fitsfiles and used the default jwst pipeline settings with the exception of the ramp-fitting weighting parameters in stage 1 that we set to uniform instead of default, as it slightly improved the rms of our residuals.This improvement can be explained by the fact that the default ramp-fitting algorithm uses a weighting of the ramp that gives additional weight to the first and last groups of the ramp, which can be problematic when the number of groups is small, such as for TRAPPIST-1 c (only 13 groups).Indeed, the first and last groups can be affected by effects such as the reset switch charge decay or saturation.Thus, to ensure that we fit the ramp correctly, we used an unweighted algorithm that applies the same weight to all groups.Furthermore, in stage 2 we turned off the photom step.Then, in stage 3, we defined a subarray region ([632, 752],[450, 570]), masked the pixels flagged in the DQ array, interpolated bad pixels and performed aperture photometry on the star with an aperture size that minimized the rms of the residuals for each visits.For each integration, we recorded the centre and width of the PSF in the x and y directions after fitting a 2D Gaussian.We computed the background on an annulus of 20 to 35 pixel (centred on the target) and subtracted it.We note that the choice of the background annulus has little impact on the light curve.We did not remove any integrations a priori but, in stage 4 we sigma clipped 4σ outliers compared with the median flux calculated using a 10-integrations-width boxcar filter.Then, for each visit for aperture photometry, we chose the aperture radius that led to the smaller rms.These radii were 3.7, 4.0, 3.6, and 3.8 pixels, respectively (see Table 2).

Data Reduction MG
We reduced the data using the following methodology.Starting from the uncal.fitsfiles, we calibrated them using the two first stages of the Eureka!pipeline [9].We performed a systematic exploration of all the combinations of all Eureka!stage 1 options, and we selected the combination resulting in the most precise light curves.Our selected combination corresponds to the default jwst pipeline settings, except for (1) the ramp-fitting weighting parameter set to uniform, and (2) the deactivation of the jump correction.The rest of the reduction was done using a pipeline coded in IRAF and Fortran 2003.
It included for each calibrated image (1) a change of unit from MJy/sr to recorded electrons, (2) the fit of a 2D Gaussian function on the profile of the star to measure the subpixel position of its centroid and its full width at half maximum (FWHM) in both directions, and (3) the measurement of the stellar and background fluxes using circular and annular apertures, respectively, with IRAF/DAOPHOT [33].Finally, the resulting light curves were normalized and outliers were discarded from them using a 5σ clipping with 20-min moving median algorithm.For each visit, the radius of the circular aperture used to measure the stellar flux was optimized by minimizing the standard deviation of the residuals.For each stellar flux measurement, the corresponding error was computed taking into account the star and background photon noise, the readout noise, and the dark noise, and assuming a value of 3.1 el/ADU for the gain (E.Ducrot, private communication).See Table 2 for more details.

Data Reduction PT
We performed an additional analysis using the level 2 (flux-calibrated) "calints" science products as processed by the Space Telescope Science Institute and hosted on the MAST archive.We determined centroid positions and average seeing FWHM values in the x and y dimensions with a 2D Gaussian fit to the star.We performed fixed-aperture photometry with circular apertures centred on the source centroids, with radii ranging from 3.2-5.0pixels in 0.1-pixel increments.We also performed variable-aperture photometry using circular apertures with radii set to c times a smoothed time series of the measured FWHM values, where c ranged from 0.75-1.25 in increments of 0.05.We smoothed the FWHM values using a 1D Gaussian kernel with a standard deviation of 2. For both fixed-aperture and variable-aperture photometry, we measured the background using a circular annulus with an inner radius of 30 pixels and an outer radius of 45 pixels.We subtracted the sigma-clipped mean of the pixel values within this annulus from the source counts in each frame, using a clipping level of 4σ.Finally, we recorded the values of a grid of background-subtracted pixels interior to the average photometric aperture size surrounding the source centroid in each frame.We used normalized time series of these pixel values to test whether pixel-level decorrelation (PLD) methods developed for minimizing intrapixel effects in Spitzer Space Telescope data [34] are warranted in the analysis of JWST/MIRI time-series data.
We excluded the first integration of each visit from our analysis as the measured source flux in this exposure was found to be significantly lower than the remainder of the time series for each of the four visits.We checked for outliers in each visit by performing sigma clipping with a threshold of 4σ, but no exposures were flagged with this step.We then selected the aperture size and method (fixed or variable) that minimized the out-of-eclipse scatter for each visit for use in our analysis.We found that fixed-aperture photometry provided the best performance in each case, with optimal radii of 4.4, 4.1, 3.9, and 3.5 pixels for the four visits, respectively.

Data Analysis
We fitted each of the reductions to extract an eclipse depth value.The different assumptions for the four global fits are listed in Table 4.

Data Analysis SZ
We fitted the eclipse light curve using the open-source python MCMC sampling routine emcee [35].Our full fitting model, F (t), was the product of a batman [36] eclipse model, F eclipse (t) and a systematic model, F sys (t).We fit the systematics of JWST with a model of the following form: where F polynom is a polynomial in time and F x (t), F y (t), F σx (t), and F σy (t) detrend the light curve against a time series of the centroid in x and y and the width of the PSF in x and y, respectively.Before fitting the full light curve consistent out of the four visits, we first determined the best systematic model for each visit by minimizing the Bayesian Information Criterion (BIC) [37][38][39].We tried a range of polynomials ranging from zeroth order to third order and detrended for the shift in x-and y-pixel positions or for the change in the width of the PSF in time.The best final combination of polynomials and detrending parameters for each visit are listed in Table 4.Our eclipse model used the predicted transit times from [40] which accounts for the transit-timing variations (TTVs) in the system and we allowed for a non-zero eccentricity.
We also accounted for the light travel time, which is approximately 16 seconds for TRAPPIST-1 c, that is, its semi-major-axis is about 8 light-seconds.We fixed the other parameters of the planet and system, such as the ratio of the semi-major axis to stellar radius a/R * , the ratio of the planetary radius to stellar radius R p /R * , and the inclination i, to the values reported in Ref [40].
We decided to also supersample the light curve by a factor of 5 in our fitting routine because the sampling of the data (≈ every 40 seconds) is comparable with the ingress/egress duration of 200 seconds [40].Our global fit consisted of 32 free parameters: 6 physical (the eccentricity, the argument of periastron, and an eclipse depth for each visit), 22 parameters to fit for the systematics, and 4 free parameters that inflated the uncertainties in the flux for each visit.These four free parameters are necessary because the current gain value on the Calibration References Data System (CRDS) has been empirically shown to be wrong for MIRI data [10].For our global MCMC, we used 128 walkers (= 4 times the number of free parameters), 150,000 steps, and discarded the first 20% of steps (= 30,000 steps) as burn-in.This corresponds to approximately 80 times the autocorrelation length.After calculating a weighted average of the four eclipse depths, we get an eclipse depth of f p /f * = 431 +97 −96 ppm for this reduction.Figure 8 shows the Allan deviation plots of the residuals for each of the visits and the global fit.The rms of the residuals as a function of bin size follows the inverse square root law, which is expected for Gaussian noise.

Data Analysis ED
Once we obtained the light curve for each visit from stage 4 of the Eureka!pipeline we used the Fortran code trafit which is an updated version of the adaptive MCMC code described in [41][42][43].It uses the eclipse model of [44] as a photometric time series, multiplied by a baseline model to represent the other astrophysical and instrumental systematics that could produce photometric variations.First, we fit all visits individually.We tested a large range of baseline models to account for different types of external sources of flux variations/modulations (instrumental and stellar effects).This includes polynomials of variable orders in time, background, PSF position on the detector (x, y) and PSF width (in x and y).Once the baseline was chosen, we ran a preliminary analysis with one Markov chain of 50,000 steps to evaluate the need for rescaling the photometric errors through the consideration of a potential underestimation or overestimation of the white noise of each measurement and the presence of time-correlated (red) noise in the light curve.After rescaling the photometric errors, we ran two Markov chains of 100,000 steps each to sample the probability density functions of the parameters of the model and the physical parameters of the system, and assessed the convergence of the MCMC analysis with the Gelman & Rubin statistical test [45].For each individual analysis, we used the following jump parameters with normal distributions: M ⋆ , R ⋆ , T eff,⋆ , [Fe/H], t 0 , b all priors were taken from [46] except for the transit timings, which were derived from the dynamical model predictions by [40].We fixed P , i and e to the literature values given in [40,46].The eclipse depth that we computed for each visit individually were 445 ± 193 ppm, 418 ± 173 ppm, 474 ± 158 ppm, and 459 ± 185 ppm, respectively.We then performed a global analysis with all four visits, using the baseline models derived from our individual fits for each light curve.Again, we performed a preliminary run of one chain of 50,000 steps to estimate the correction factors that we then apply to the photometric error bars and then a second run with two chains of 100,000 steps.The jump parameters were the same as for the individual fits except for the fact that we fixed t 0 and allowed for transit timing variations (TTV) to happen for each visit (each transit TTV has an unconstrained uniform prior).We used the Gelman & Rubin statistic to assess the convergence of the fit.We measure an eclipse depth of 423 +97 −95 ppm from this joint fit.

Data Analysis MG
Our data-analysis methodology was the same as that used by ED, that is, we used the Fortran 2003 code trafit to perform a global analysis of the four light curves, adopting the Metropolis-Hasting MCMC algorithm to sample posterior probability distributions of the system's parameters.Here too, we tested for each light curve a large range of baseline models, and we adopted the ones minimizing the BIC.They were (1) a linear polynomial of time for the first visit, (2) a cubic polynomial of time and a linear polynomial of the y position for the second visit, (3) a linear polynomial of time and of the x position for the third visit, and (4) a cubic polynomial of time and of the y position for the fourth visit.We also performed a preliminary analysis (composed of one Markov Chain of 10,000 steps) to assess the need to rescale the photometric errors for white and red noise.We then performed two chains of 500,000 steps each (with the first 20% as burn-in).The convergence of the analysis was checked using the Gelman & Rubin statistical test [45].The jump parameters of the analysis, that is, the parameters perturbed at each step of the MCMC chains, were (1) for the star, the logarithm of the mass, the logarithm of the density, the effective temperature, and the metallicity, and (2) for the planet, the planet-to-star radius ratio, the occultation depth, the cosinus of the orbital inclination, the orbital parameters √ e cos ω and √ e sin ω (with e the orbital eccentricity and ω the argument of pericentre), and the timings of the transits adjacent to each visit.We assumed normal prior distributions for the following parameters based on the results from reference [40]: M * = 0.0898 ± 0.023, R * = 0.1192 ± 0.0013, T eff = 2566 ± 26 K, and [Fe/H] = 0.05 ± 0.09 for the star; (R p /R * ) 2 = 7123 ± 65 ppm, b = 0.11 ± 0.06, e = 0 + 0.003 (semi-gaussian distribution) for the planet.We also tested the assumption of a circular orbit and obtained similar results.For each visit, we considered for the timings of the two adjacent transits normal prior distributions based on the predictions of the dynamical model of reference [40].At each step of the MCMC, the orbital position of the planet could then be computed for each time of observation from the timings of the two adjacent transits and from e and ω, and taking into account the approximately 16s of light-travel time between occultation and transit.This analysis led to the value of 414 ± 91 ppm for the occultation, and to an orbital eccentricity of 0.0016 +0.0015 −0.0008 consistent with a circular orbit.Under the assumption of a circular orbit, our analysis led to an occultation depth of 397 ± 92 ppm, in excellent agreement with the result of the analysis assuming an eccentric orbit.
We also performed a similar global analysis, but allowing for different occultation depths for each visit.The resulting depths were 400± 163 ppm, 374 ± 184 ppm, 421 ± 187 ppm, and 403 ± 202 ppm, i.e. they were consistent with a stable thermal emission of the planet's dayside (at this level of precision).Similar to data reduction SZ, we also did create Allan deviation plots for this particular data reduction.The best-fit residuals as a function of bin size from each visit do generally follow the inverse square root law (see Fig. 8 for the Allan deviation plots of data reduction SZ).
Finally, we computed the brightness temperature of the planet at 15 µm from our measured occultation depth using the following methodology.We measured the absolute flux density of the star in all the calibrated images, using an aperture of 25 pixels large enough to encompass the wings of its PSF.We converted these flux densities from MJy/sr to mJy, and computed the mean value of 2.559 mJy and the standard deviation of 0.016 mJy.We added quadratically to this error of around 0.6% a systematic error of 3%, which corresponds to the estimated absolute photometric precision of MIRI (P.-O.Lagage, private communication).It resulted in a total error of 0.079 mJy.Multiplying the measured flux density by our measured occultation depth led to a planetary flux density of 1.06 ± 0.23 µJy.Multiplying again this result by the square of the ratio of the distance of the system and the planet' radius, and dividing by π, led to the mean surface brightness of the planet's dayside.Applying Planck's law, we then computed the brightness temperature of the planet, while its error was obtained from a classical error propagation.Our result, for this specific reduction, was 379 ± 30 K, to be compared with an equilibrium temperature of 433 K computed for a null-albedo planet with no heat distribution to the nightside.
It is also worth mentioning that applying the same computation on the star itself led to a brightness temperature of 1867 ± 55 K, which is significantly lower than its effective temperature.

Data Analysis PT
We began our analysis by determining which time-series regressors (if any) should be included for fitting systematics in the photometry on the basis of the BIC.Our total model is the product of a batman eclipse model (F eclipse ) and a systematics model (F syst ) to the data, which has a general form of Here, F polynom is a polynomial in time, F x (t) and F y (t) are time series of the target centroids in x and y, F FWHM (t) is the time series of average FWHM values for the source determined with a 2D Gaussian fit, and F ramp is an exponential function that accounts for ramp-up effects.F PLD (n, t) is the linear combination of n basis pixel time series, and it has a form of Here, Pi (t) is the normalized intensity (from 0-1) of pixel i at time t and C i is the coefficient of pixel i determined in the fit.PLD was developed to mitigate systematic intrapixel effects in Spitzer/Infrared Array Camera (IRAC) data [34], in which the combination of source PSF motion and intrapixel gain variations introduced percent-level correlated noise in time-series data [e.g., 47].
In our analysis, we tested forms of F polynom ranging from degree 0-3 and different sets of PLD basis pixels including the brightest 1, 4, 9, 16, 25, and 36 pixels.For each visit, we explored grids of every possible combination of the components of F syst (t).For each combination, we first initialized the coefficients of each component using linear regression.We then used emcee to perform an MCMC fit of the total eclipse and systematic model to the visit data.We ran 2ν+1 walkers for 10,000 steps in each fit, in which ν represents the number of free parameters in the total model.The first 1,000 steps of these chains were discarded as burn-in.We fit for seven physical parameters in our calculation of F eclipse , these being the orbital period, a/R * , orbital inclination, eccentricity, longitude of periastron, eclipse depth, and time of secondary eclipse.Gaussian priors were assigned to these parameters with means and standard deviations set by their measurements reported in [40].We also placed Gaussian priors on the coefficients of the components of F syst , with means set by the linear regression fit and standard deviations set to the absolute value of the square root of those values.
We calculated the BIC of the best-fitting model that resulted from the MCMC analysis, and then selected the form of F syst that minimized the BIC.The form of F syst that we determined for each visit with this approach consisted of only an F polynom component.The first visit was best fit by a linear polynomial, whereas the remaining three were best fit by a quadratic polynomial.
With the form of F syst (t) determined for each visit, we then performed a joint fit of all four eclipses.This fit included 18 total free parameters: 7 physical and 11 for fitting systematics (see Table 4).We ran this fit with 64 chains for 50,000 steps, discarding the first 5,000 steps for burn-in.We measured a resulting eclipse depth of 418 +90 −91 ppm from this fit.

Brightness Temperature Calculation
The following analysis was based on stage 0 (.uncal) data products preprocessed by the JWST data processing software version number 2022 3b, and calibrated with Eureka! as described above in the "Data Reduction MG" section.We computed the brightness temperature of the planet at 15 µm from our measured occultation depth using the following methodology.We measured the absolute flux density of the star in all the calibrated images, using an aperture of 25 pixels large enough to encompass the wings of its PSF.We converted these flux densities from MJy/sr to mJy, and computed the mean value of 2.559 mJy and the standard deviation of 0.016 mJy.We added quadratically to this error of about 0.6% a systematic error of 3%, which corresponds to the estimated absolute photometric precision of MIRI (P.-O.Lagage, private communication).It resulted in a total error of 0.079 mJy.Multiplying the measured flux density by our measured occultation depth led to a planetary flux density of 1.06 ± 0.23 µJy.Multiplying again this result by the square of the ratio of the distance of the system and the planet' radius, and dividing by π, led to the mean surface brightness of the planet's dayside.Applying Planck's law, we then computed the brightness temperature of the planet, whereas its error was obtained from a classical error propagation.Our result, for the MG reduction, was 379 ± 30 K, to be compared with an equilibrium temperature of 433 K computed for a null-albedo planet with no heat distribution to the night side.It is also worth mentioning that applying the same computation on the star itself led to a brightness temperature of 1867 ± 55 K, which is significantly lower than its effective temperature.

Emission modelling for TRAPPIST-1 c
We generated various emission spectra for TRAPPIST-1 c to compare them to our measured eclipse depth at 15 µm.These models include (1) bare-rock spectra, (2) O 2 /CO 2 mixture atmospheres and pure CO 2 atmospheres, and (3) coupled climate-photochemical forward models motivated by the composition of Venus.In the following, we describe each of these models.

Bare Rock
Our bare-rock model is a spatially resolved radiative transfer model and computes scattering and thermal emission for a variety of surface compositions.For each composition, the surface's radiative equilibrium temperature is computed on a 45x90 latitude-longitude grid, assuming TRAPPIST-1 c is tidally locked.Surface reflectance and emissivity data are from [22], which were derived from reflectance spectra of rock powders or minerals measured in the laboratory combined with an analytical radiative-transfer model [48].These data have previously been used to model surface albedos and emission spectra of barerock exoplanets [7,22,49].Here we consider six compositions as well as a blackbody: basaltic, feldspathic, Fe-oxidized (50% nanophase haematite, 50% basalt), granitoid, metal-rich (FeS 2 ), and ultramafic (see Fig. 9).Given the uncertainty in the measured eclipse depth, we assume a Lambertian surface with isotropic scattering and emission, and neglect the angular dependency of the surface reflectance and emissivity, which would depend on the surface roughness and regolith particle size [22].Sensitivity tests show that these surface-model assumptions are indistinguishable within the current precision of the TRAPPIST-1 c measurements (not shown).Furthermore, albedos and spectra of bare rocks in the Solar System are modified by space weathering, so we also consider the impact of space weathering on TRAPPIST-1 c.The timescale for lunar space weathering through exposure to the solar wind has been estimated to range from ∼ 10 5 years to ∼ 10 7 years [50,51].We extrapolate from the lunar value to TRAPPIST-1 c using scaling relations from a stellar-wind model [52].We find that the space-weathering timescale for TRAPPIST-1 c is significantly shorter than the lunar value, about 10 2 − 10 3 years, largely because of the planet's small semi-major axis.An exposed surface on TRAPPIST-1 c should therefore have been substantially weathered.To simulate the impact of space weathering on unweathered surfaces, we incorporate the same approach as that in [23,53].The surface composition is modelled as a mixture of a fresh host material (described above) and nanophase metallic iron using Maxwell-Garnett effective medium theory.The refractive index of metallic iron is taken from [54].

Simple 1+D O 2 /CO 2 Mixtures
We construct a grid of O 2 -dominated model atmospheres with a range of surface pressures and mixing ratios of CO 2 .These are broadly representative of a plausible outcome of planetary atmosphere evolution, in which water in the atmospheres of terrestrial planets orbiting late-type M dwarfs is photolysed and the H is lost, leaving a large O 2 reservoir [15,55].The atmosphere models we construct are 1D models following the approach presented in [2], with adiabatic pressure-temperature profiles in the deep atmosphere and isothermal pressure-temperature profiles above 0.1 bar (for thicker atmospheres, P >0.1 bar) or the skin temperature (for thinner atmospheres).This approach uses DISORT [56,57] to calculate radiative transfer in 1D through the atmosphere to generate emission spectra.
We do consider how an atmosphere can transport heat to the nightside.To include heat transport to the nightside, we implement the analytic approach of [14]; we use the redistribution factor f calculated in equation (3) of that work for each of the models in the grid.We assume that both the surface Bond albedo and the top-of-the-atmosphere Bond albedo are 0.1.We construct a grid of O 2 -dominated model atmospheres with surface pressures from 0.01 to 100 bar (in 1-dex steps) and CO 2 mixing ratios from 1 ppm to 10,000 ppm (in 1-dex steps).We also generate pure CO 2 atmospheres with the same surface pressures.For the thicker atmospheres (P surf ≥ 1 bar) we set the thermopause (in which the atmosphere transitions from adiabatic to isothermal) to 0.01 bar.

Coupled Climate-Photochemical Venus-like Atmospheres
We use a 1.5D coupled climate-photochemical forward model [3,58,59, VPL Climate] that explicitly models day and night hemispheres with layer-by-layer, day-night advective heat transport driven by simplified versions of the 3D primitive equations for atmospheric transport, to simulate plausible atmospheric states for TRAPPIST-1 c for cloudy Venus-like scenarios.VPL Climate uses SMART [60] with DISORT [56,57] for spectrum-resolving radiative transfer for accuracy and versatility for both the climate modelling and the generation of the resulting planetary spectra.The model has been validated for Earth [61] and Venus [62], but is capable of modelling a range of atmospheric states.
Due to the early luminosity evolution of the star, TRAPPIST-1 c would have been subjected to very high levels of radiation [63], and so we would anticipate evolved atmospheres that had undergone atmospheric and possibly ocean loss [6].We start with the self-consistently coupled climate-photochemical Venus-like atmospheres generated for an evolved TRAPPIST-1 c from [3], with 96.5% CO 2 and 3.5% N 2 and assume Venus lower atmospheric trace gases and self-consistent generated sulfuric acid aerosols.We use these atmospheres as a starting point for 1.5D clear-sky Venus-like atmospheres (0.1, 1, 10 bars) and 1.5D cloudy Venus-like atmospheres (10 bars) with sulfuric acid haze.Note that 10 bar Venus-like atmospheres will produce similar results to a 93 bar Venus-like atmosphere due to the emitting layer being above or at the cloud deck, which is at a similar pressure for the 10 bar and 93 bar cases.All the modelled clear-sky Venus atmospheres produce 15 µm CO 2 features with depths spanning 134-143 ppm, with the cloudy 10 bar Venus centered at 181 ppm.Because H 2 SO 4 aerosols are likely to condense in the atmosphere of a Venus-like planet at TRAPPIST-1 c's orbital distance [3], we show the dayside spectrum for the 10 bar cloudy Venus for comparison with the data in Figure 3.The emitting layer (cumulative optical depth 1) for the cloud aerosols occurs at 7 mbars in this atmosphere, although the 15 µm CO 2 absorption is sufficiently strong that it emits from a comparable pressure level in the core of the band.The observations rule out a self-consistent Venus-like atmosphere for TRAPPIST-1 c to 2.6 sigma.

Atmospheric Escape models
We use energy-limited atmospheric escape models [6,15] from a steam atmosphere to explore the amount of atmospheric escape that TRAPPIST-1 c may have experienced over its lifetime.The model assumes that escape occurs in the stoichiometric ratios of H/O in water vapour, allows for escape of oxygen, and reaction of with the magma ocean.The model transitions from magma ocean to passive stagnant-lid outgassing when surface temperatures drop below the silicate melting point.Escape continues throughout all tectonic stages.In Figure 4, we show the final amount of O 2 gas left in the atmosphere after 7.5 Gyr of evolution for a range of planetary water abundances and XUV saturation fractions.For typical saturation fractions of 10 −3 [25,64], our observations suggest that the planet likely had a relatively low starting volatile abundance.We note that these models are likely upper limits on thermal escape and more detailed models of escape, especially incorporating other gases such as CO 2 and N 2 , are needed in the future to confirm these results.We also estimate total ion-driven escape fluxes due to stellar wind interactions of a minimum of 1-3 bars over the planet's lifetime, assuming constant stellar wind over time [65].We also considered the extended pre-main sequence for a star like TRAPPIST-1.We used the stellar evolution models of [63] for a 0.09 M ⊙ star to approximate the pre-main sequence evolution of the star.

Interior structure model
We use an interior-structure model to perform an MCMC retrieval on the planetary mass and radius of TRAPPIST-1 c, and the possible stellar Fe/Si of TRAPPIST-1.The estimated Fe/Si mole ratio of TRAPPIST-1 is 0.76±0.12[66], which is lower than the Solar value, Fe/Si = 0.97 [67].Our interiorstructure model solves a set of differential equations to compute the density, pressure, temperature, and gravity as a function of radius in a one-dimensional grid [68,69].The interior model presents two distinct layers: a silicate-rich mantle and an Fe-rich core.On top of the mantle, we couple the interior model with an atmospheric model to compute the emission and the Bond albedo.These two quantities enable us to solve for radiative-convective equilibrium, find the corresponding surface temperature at the bottom of the atmosphere, and find the total atmospheric thickness from the surface up to a transit pressure of 20 mbar [20,70].We consider a H 2 O-dominated atmosphere, with 99% H 2 O and 1% CO 2 .Our 1D, k-correlated atmospheric model prescribes a pressure-temperature profile comprised of a near-surface convective layer and an isothermal region on top.In the regions of the atmosphere where the temperature is low enough for water to condense and form clouds, we compute the contribution of these to the optical depth and their reflection properties as described in [71,72].
The posterior distribution function of the surface pressure retrieved by our MCMC indicates a 1σ confidence interval of 40±40 bar for TRAPPIST-1 c. Surface pressures between 0 and 120 bar would be compatible with our probability density function within 2σ [73].Oxygen is more dense than H 2 O. Consequently, for a similar surface pressure, an O 2 -rich atmosphere would be less extended than the H 2 O-dominated envelope we consider in our coupled interior-atmosphere model.This means that the density of TRAPPIST-1 c could be reproduced with an oxygen-rich atmosphere with a lower surface pressure as low as our H 2 O upper limit, 80 bar.

Stellar Properties
The stellar properties of TRAPPIST-1 have been constrained with observations of the total luminosity of the star, L * = 4πD 2 * F bol (based on broadband photometry to obtain the bolometric flux of the star, F bol , and a distance measured with GAIA, D * ), a mass-luminosity relation [74] to obtain the stellar mass, M * (L * ), with uncertainty, as well as a precise stellar density, ρ * , thanks to modelling of the seven transiting planets [40,46,75].These combine to give the stellar radius and effective temperature, R * and T eff, * , The planets' properties have also been measured precisely in relation to the star using the depths of transit, yielding R p /R * , and transit-timing variations, yielding M p /M * [40].To convert the secondary eclipse depth, δ = F p /F * , into a brightness temperature of the planet requires an estimate of the brightness temperature of the star: or in which I b, * , I b,p are the mean surface brightness of the star and planet in the MIRI band at full phase (that is, secondary eclipse), respectively.The ratio R p /R * is well constrained from the transit depth, whereas the brightness temperature of the star can be measured with an absolute calibration of the stellar flux in the MIRI band, F * (e.g.[46]).The stellar intensity may then be computed as: in which Ω * is the solid angle of the star.Because our estimate of R * is proportional to M .For a given value of R * , this surface brightness can be translated into a brightness temperature, T b, * , and with the equation 7 above, we can compute the intensity and therefore the surface brightness of the planet, T b,p to be 380 ± 31 K using the eclipse depth and the stellar flux density.We also estimate the stellar brightness temperature in the MIRI band with an atmospheric model for the star relating T b, * in the MIRI band to the T eff , as α = T b, * /T eff, * .We have accomplished this with the state-of-the-art SPHINX model for low-temperature stars [76] and assumed T eff = 2566 K [40], yielding α = 0.72 at 14.87 µm.We also compute the α from JWST spectrophotometric observations with a flux of 2.599 ± 0.079 mJy at 14.87 µm, yielding α = 0.71 ± 0.02.The MIRI images are flux-calibrated (with an internal error of 3%).We measure the stellar flux in all images within an aperture large enough to encompass the whole PSF, and then compute the mean and the standard deviation.We compute the total error on the measurement to be the 3% larger than this standard deviation.As the unit of flux in MIRI images is given in Jy/str, we multiply the measured fluxes by the angular area covered by a pixel in str to yield units of Jy.
The stellar brightness temperature scales linearly with effective temperature and metallicity in the MIRI wavelength range, and scales inversely with surface gravity of the star.The effective temperature, however, scales as T eff ∝ M ), with stellar mass (or radius) relative to the estimate based on the measured flux.The estimate of α, therefore, may have a significant imprecision given the possible heterogeneity of the stellar atmosphere, as well as the inherent uncertainties involved in modelling late-type stellar atmospheres.Both the synthetically derived α and those from observations match within 2σ uncertainty, lending credence to empirical mass-luminosity relations and synthetic atmosphere-model-derived stellar brightness temperatures.Note, however, that the mass-luminosity relation is only calibrated with a handful of low-mass stars in binaries [74], and hence its applicability to TRAPPIST-1 may be tenuous; this may thus be the weakest link in determining the stellar parameters.Assumption-driven deviations between synthetic models for late-type stars and empirically calibrated methods both still remain a significant challenge in truly understanding these hosts.

Eclipse Timing Variations
Dynamical modelling of the TRAPPIST-1 system [40] gives a precise forecast of the times of transit and eclipse for all seven planets.These have been used in the planning of the observations, and can also be compared with the measured times.
The times of eclipse can be offset from the mid-point between the times of transit due to four different effects: 1) the light-travel time across the system  [77], 2) non-zero eccentricity [78], 3) non-uniform emission from an exoplanet [79] 1 , and 4) eclipse-timing variations due to perturbations by other planets in the system.Of these three effects, the second effect is typically the largest, which can be used to constrain one component of the eccentricity vector of the transiting planet [78].
In Table 6 we list the measured eclipse times from the four different reductions and in Figure 10 we compare them with the forecast from [40].To make the forecast, we used the posterior probability of the timing model to compute the times of transit and eclipse, and then we calculate the time of eclipse minus the mean of the two adjacent transits of planet c to derive an "Eclipse timing offset".This offset should be zero for a circular, unperturbed orbit with negligible light-travel time (which is about 16 seconds, or 1.8×10 −4 days for TRAPPIST-1 c).The dynamical modelling is constrained by the times of transit, which place some constraint on the eccentricity of the orbit of planet c (in particular, the mean or free eccentricity could be non-zero).The uncertainty on the eccentricity leads to uncertainty on the times of secondary eclipse.Our forecast models for the eclipse timing offset have a 1σ uncertainty of ∼3.5 minutes at the measured times of eclipse (approximately 0.0024 days).
The measured times were taken from four analyses (by SZ, PT, ED, and MG) in which a broad prior was placed on the times of transit, whereas the duration and depth were constrained to the measured values of the four eclipses.The times of each eclipse were then free to vary, and the posterior times of transit were inferred using MCMC (ED/MG/PT) or nested sampling (SZ).The four analyses give good agreement on the values, but have significant differences between the uncertainties.
Overall, the forecast eclipse timing offsets agree well with the measured times, within 1-2σ offsets.The uncertainties on the measured times are comparable with the forecast uncertainties, and so, in future work, we hope to use these measured eclipse times to further constrain the eccentricity vector of the orbit of planet c.This may help to constrain tidal damping models of planet c, but it may also constrain tidal damping of all of the planets as the free eccentricity vector of planet c is tightly correlated with those of the other planets due to the "eccentricity-eccentricity" degeneracy present in transiting planet systems [80].a data reduction and data analysis of the four visits for this work and contributed an eclipse depth value.CM, DK, XL, RH, AL, and VM ran theoretical models for the planet's atmosphere and surface.LA ran models on the planet's interior structure.AI and EA modeled the stellar spectrum.LS modeled the atmospheric escape for the planet.LK, MG, VM, DK, RH, CM, LS, EA, FS EB, and AM contributed to the observing proposal.

Fig. 1
Fig. 1 Eclipse light curve of TRAPPIST-1 c taken with MIRI F1500W.The phase-folded secondary eclipse light curve of TRAPPIST-1 c, measured with the JWST/MIRI imager at 15 µm.The eclipse is centred at orbital phase 0.5 and has a measured depth of fp/f * = 421 ± 94 ppm.The light curve includes four visits (that is, four eclipses), each spanning approximately 3.2 hours.To make the eclipse more easily visible, we binned the individual integrations (grey points) into 28 orbital phase bins (black points with 1σ error bars).The light curve was normalized and divided by the best-fit instrument systematic model.The best-fit eclipse model is shown with the solid red line.The data and fit presented in this figure are based on the SZ reduction, one of the four independent reductions we performed in this work.

Fig. 2
Fig. 2 Grid plot comparing a suite of atmospheric models to the measured eclipse depth.Comparison between the measured eclipse depth and a suite of different O 2 /CO 2 , cloud-free atmospheres for TRAPPIST-1 c with varying surface pressures and compositions.Darker grid cells indicate that we more significantly rule out this specific atmospheric scenario.The number in each cell is the absolute difference between each model and the observations in units of sigma.The lower the modelled atmosphere is in the grid, the higher its surface pressure.The rightmost column shows pure CO 2 atmospheres.The other columns are O 2 -dominated atmospheres with different amounts of CO 2 ranging from 1 ppm (= 0.0001%) to 10,000 ppm (= 1%).

Fig. 3
Fig.3Observed flux of TRAPPIST-1 c and various emission models.Simulated emission spectra compared with the measured eclipse depth of TRAPPIST-1 c (red diamond, with the vertical error bar representing the 1σ uncertainty on the measured eclipse depth).The CO 2 feature overlaps directly with the MIRI F1500W filter used for these observations.The two limiting cases for the atmospheric circulation for a zero-albedo planet (zero heat redistribution, that is, instant reradiation of incoming flux and global heat redistribution) are marked with dashed lines.Two cloud-free, O 2 /CO 2 mixture atmospheres are shown with purple and red solid lines.They show decreased emission at 15 µm owing to CO 2 absorption.A bare-rock model assuming an unweathered ultramafic surface of the planet with a Bond albedo of 0.5 is shown by the solid black line (see text for more information on weathering, including a full comparison of our measurement to a suite of surfaces in Fig.9).The cloudy Venus forward model with a surface pressure of 10 bars is shown with a solid yellow line.

Fig. 6
Fig. 6 Example of a MIRI integration using the FULL array.An integration taken during our observations showing the MIRI imager focal plane.The majority of the FULL array is taken up by the imager field of view on the right side.TRAPPIST-1 is centred on the imager highlighted by the red arrow.The left side of the imager was not used in our analysis and consists out of the Lyot coronagraph (top left) and the three 4-quadrant phase masks coronagraphs (lower left).

Fig. 7
Fig. 7 Diagnostic plot of all four visits taken during JWST program GO 2304 based on the SZ reduction.Every column corresponds to a visit.a and b.The top and second rows show the raw and background flux in units of electrons per integration per pixel, respectively.The raw flux is referring to the flux level within the target aperture before the subtraction of the background flux.c -f.The following rows are depicting the properties of the centroid over time.We fitted a 2D Gaussian distribution to the target at every integration to determine its x and y positions on the detector.∆σx and ∆σy describe change in the width of the 2D Gaussian with time.The integrations were taken approximately every 40 seconds.The lower four rows were additionally binned to 5 minutes (= 8 integrations) shown with the solid black lines.Due to stronger systematics, we excluded the first 10 integrations in the SZ reduction shown by the grey region at the beginning of each visit.

Fig. 8
Fig. 8 Allan deviation plots.a -d.Allan deviation plots of the individual visits: rootmean-square (rms) of the best-fit residuals from data reduction SZ as a function of the number of data points per bin shown in black.e.The same, but for the combined dataset.A bin size value of one corresponds to no binning.The red line shows the expected behaviour if the residuals are dominated by Gaussian noise.The absolute slope of this line is 1/ √ bin size, following the inverse square root.The rms of our residuals closely follow this line, showing that our residuals are consistent with uncorrelated photon shot noise.

Fig. 9
Fig.9Measured eclipse depth compared to a suite of simulated bare-rock emission spectra.a. Secondary eclipse spectra for various fresh surface compositions, assuming that TRAPPIST-1 c is a bare rock.High-albedo feldspathic and granitoid surfaces are cool and fit the data moderately poorly (2σ), as does a low-albedo and hot blackbody surface (1.7σ).b.Space weathering via formation of iron nanoparticles (npFe) lowers the albedo at short wavelengths, thereby increasing the surface's temperature and its secondary eclipse depth.An intermediate-albedo fresh ultramafic surface would fit the data well, but the fit becomes marginal after taking into account the influence of strong space weathering (1.6σ, or about 90% confidence).The vertical error bar on our 15 µm measurement represents the 1σ uncertainty on the observed eclipse depth.

Fig. 10
Fig. 10 Measured eclipse times compared with the predicted eclipse times.The points show the measured eclipse timing offsets (defined as the time of eclipse minus the mean of the two adjacent transit times of planet c) from four different analyses.The error bars correspond to the 16th and 84th percentiles of the eclipse time posterior.The dark (light) green shaded region shows the 1-(2-)sigma confidence intervals forecast from the transit-timing analysis from [40].

Table 1
Summary of the observations in JWST program GO 2304.

Table 2
Details of the four different data reductions.
no integrations were removed during this process, leaving us with 288, 287, 287, and 288 integrations for the four visits, respectively.

Table 4
Details of the four data analyses.See Methods for more details on the individual fits.The uncertainties on the eclipse depth fp/f * are 1σ.

Table 6
Measured eclipse times.Eclipse times (in BJD TDB ) determined by the four different reductions by fitting an eclipse model to each visit.An offset of 2459880.0 days was subtracted from each of the values in the table.