本文来自作者[lejiaoyi]投稿,不代表言希号立场,如若转载,请注明出处:https://lejiaoyi.cn/zixun/202506-1052.html
我们使用无名,Nirhiss56,SupremeSpoon57和TranspectRosproscopy58管道对Niriss/Soss/Soss/Soss/Soss 13,54,55的Eclipse观察到了四个单独的减少,用于单个还原步骤59。所有管道均围绕官方的STSCI JWST减少管道60构建,并增加了用于系统的定制校正步骤,例如1/F噪声,黄道十二折叠背景和宇宙射线。从原始未校准的数据(无名,Nirhiss和最高勺)或第1阶段(转镜镜)产物进行降低,直到提取分光光度计光曲线。
我们使用无名管道59通过光谱提取来减少未校准数据产品的WASP-18 B观测值。我们使用了JWST管道版本1.6.0,CRD(校准参考数据系统)版本11.16.5和CRDS上下文JWST_0977.PMAP减少。首先,除了dark_current步骤外,我们经历了JWST管道阶段1的所有步骤。我们跳过了黑暗的减法步骤,以避免由于缺乏高保真参考文件而引入额外的噪音。在坡道拟合步骤之后,我们浏览了JWST Pipeline阶段2的分配_WC,SRCTYPE和FLAT_FIELD步骤;我们跳过背景步骤,并应用自己的自定义程序来处理背景。由于不需要绝对的通量校准,因此我们跳过了路径和光子步骤。我们通过计算所有帧59的列中的第二个衍生物的乘积来执行异常检测。我们将框架分为4×4像素的窗口,然后为此计算第二个衍生物的本地中值和标准偏差。我们标记所有距窗口中位数≥4σ的像素。此外,我们标志着显示无效或负数的像素。所有标记的像素设置为等于其窗户的当地中间。我们使用以下例程更正背景系统。首先,我们将第(x [5,400],y [0,20])的部分确定为基质96子阵列的区域,其中光谱顺序对计数的贡献很小。接下来,我们计算上述区域内STSCI JDOX用户文档中提供的中位数和模型背景之间的缩放系数。我们将分布的第16个百分点视为缩放值,并从所有集成中减去缩放背景。我们密切关注这些观察结果的1/F校正,因为光谱痕量变化的大小高度依赖于次生食物中的波长。所以, 当缩放迹线计算1/F噪声时,我们会独立考虑所有列。此外,我们将这种校正分为两个部分,因为我们观察到倾斜事件20,21在第1,336个整合左右(扩展数据图2),这可能是由于主要镜像的一个细分市场突然运动,导致痕迹形态的变化,突然显示出横向变化。首先,我们在倾斜事件之前和之后计算中位列;我们使用集成300–900和1,350–1,900。然后,我们在集成i处定义一个给定的列j和行k,作为缩放中间列和1/f噪声的总和。使用JWST管道返回的错误εi,j,k,我们解决了Mi,j和ni,J的值,从而最小化了观察到的和缩放柱之间的卡方。通过施加∂χ2/∂mi,j,∂χ2/∂ni,j = 0来获得缩放Mi,j和1/f噪声Ni,J的这些值
和
然后,我们从所有列和集成中减去Ni,J的测量值。我们将所有具有JWST管道返回的非零数据质量标志的像素的误差设置为∞(∞在此定义为正无穷大的IEEE 754浮点数表示),以使它们在1/f噪声的拟合中不被考虑。我们还将误差设置为区域x [76,96],y [530,1,350]的检测器,其中第二个光谱顺序的一部分是可见的。这是适当的处理方法,因为由于波长覆盖率的差异,1/F噪声在每个顺序上独立缩放。校正1/f噪声后,我们通过计算所有列的高斯滤波器的痕迹的最大值来追踪检测器上订单1的位置。我们使用样条函数进一步平滑痕量质心的位置。最后,我们使用转谱法进行一阶的盒谱提取。光谱镜检查。孔径谱的例程,孔径为30像素。
我们使用Nirhiss Python开源数据还原管道,如参考文献中所述。59。总而言之,该管道使用Eureka!62从阶段0 JWST输出到阶段2校准,该阶段校准了,该校准应用了检测器级校正,可产生计数率图像并校准单个接触。从第2阶段的“ Calints”适合文件,我们使用Nirhiss纠正背景源,1/F噪声,宇宙射线删除,痕量识别和光谱提取。我们遵循两个步骤进行跟踪识别。首先,我们使用“ jwst_niriss_spectrace_0022.fits”的参考文件(y值偏移约25个像素)从“ jwst_niriss_spectrace_0022.fits”中的跟踪位置(x,y)位置来识别订单2。我们掩盖了此区域,以至于它不会污染跟踪识别或后来的途径。接下来,我们使用nirhiss.tracing.mask_method_edges函数。该技术使用Scikit-Image(一种开源图像处理软件包63)使用Canny Edge-tection方法来标识订单1的边缘。该方法使用高斯函数的导数来识别具有最大梯度的区域。从此步骤中,沿最大值的电势边缘缩小到1像素曲线。这将导致列出顺序1的轮廓的图像。我们从订单1的顶部和底部边缘确定沿圆柱的中值位置,并通过拟合四阶多项式来平滑迹线。我们发现,四阶多项式最佳符合这两个订单的整体形状,而六阶多项式过度拟合,而二阶多项式则属于订单轮廓。我们使用迹线掩盖了踏入Nirhiss背景例程时订单1的位置。对于背景处理,我们遵循参考文献中提出的类似方法。59,即,我们在STSCI JDOX用户文档上的模型背景上确定了一个没有光谱迹线的区域,并将该区域扩展到同一区域。我们使用区域x [4,250],y [0,30],找到平均缩放系数为0.6007。我们将此缩放因子应用于模型背景,并从所有集成中减去它。接下来,我们以类似的方式去除1/f噪声,并将此1/F噪声处理扩展到放弃外部的整合(0-1,250和1,900–2,500)。我们使用L.A.Cosmic Technique删除宇宙射线64。最后,我们使用最佳提取程序提取光谱,这是一种强大的手段,可以同时删除进一步的不良像素/宇宙射线事件,同时将每个像素在每个像素上放置不均匀的权重以消除由空间配置文件产生的失真65。我们使用归一化图像来最好地捕获独特的NIRIS/SOSS空间轮廓。
我们遵循参考文献中提出的最高勺子的类似方法。59。我们从原始未校准的数据文件开始,我们从Mikulski档案库(MAST)档案中下载了该文件,并通过最高吊杆1阶段1对其进行处理,该阶段1执行探测器级校准,包括超级次数减法,饱和度标志,跳转检测,跳高检测和RAMPAMP FITT。与以前的管道一样,我们不执行任何暗电流减法。SupremeSpoon进一步处理小组级别的1/F噪声。这是通过减去所有放生内集成的中间堆栈来完成的,并通过白光曲线缩放到每个单个集成的磁通水平,以创建一个差异图像,以揭示特征1/f条纹。然后,从相应的集成中减去nth差异图像的列的中位数。如参考文献中所述,Supreme Spoon还可以通过将STSCI提供的Sudstrip96 SOSS背景模型缩放到每个集成的通量水平,直接在1/F校正步骤之前删除黄道背景信号。59。然后,我们通过Supreme Spoon Stage 2传递了阶段1加工的文件,该文件执行了进一步的校准,例如Flat Fielding和Hot像素插值。我们使用一个简单的盒子孔提取,在订单1迹线上以30像素为中心的宽度为30像素,我们在像素水平上提取出色的光谱,因为与二阶重叠所产生的稀释量已显示为可忽略不计66,67。迹线的Y像素位置是通过EdgetRigger算法67确定的。我们发现提取的痕量位置与在调试过程中测量并包含在Spectrace参考文件中的位置,因此我们使用默认的JWST波长解决方案。
我们遵循参考文献中讨论的转镜管道采用的类似方法。59。此减少从_rateInts.ts开始,适用于STSCI的JWST Pipeline处理的文件。我们将STSCI JDOX用户文档上提供的黄道背景模型缩放到了Pixels X [10,250],Y [10,30]的框中观察到的2D光谱。然后从每个集成中减去缩放背景。总而言之,删除1/F噪声的过程如下:我们将所有集成的中值从每个集成中减去,这主要是1/F噪声。然后,我们将这种残余噪声的列中位数逐列中位数,考虑到距痕迹中心20像素的像素,并假定它代表了图像的1/F噪声的结构。然后从每列中减去这些值。对于频谱提取,我们使用了镜头镜检查。光谱镜检查。孔径谱的例程,光圈宽度为30像素。我们通过取下所有光谱的合并中值和标记与此中位光谱偏离5σ的异常值的合并中位数来删除由宇宙射线引起的提取光谱的异常值或偏离像素。然后,通过取下相邻垃圾箱的平均值来纠正标记的波长箱。
从上述恒星光谱提取管道中,我们创建并拟合模型到分光光度光度曲线F(T,λ)。光曲线由三个不同的信号组成:整个部分相位曲线和日食FP的行星通量,恒星通量f*和Systematics s,我们通过等式将其模拟为时间t和波长λ的函数(3)。
由于关注的主要科学量是行星信号fp(t,λ,θ),其中θ代表行星轨道参数,我们的目标是正确地表征和纠正恒星通量和系统。曲线拟合分为两个单独的步骤:(1)我们拟合白光曲线,(2)我们对每个分光光度计箱进行单独拟合。从白光曲线获得的轨道参数的值固定在分光光度法曲线中。以下各节描述了我们对行星通量FP(t,λ,θ),恒星变异性F*(t,λ)和系统S(t,λ)的处理。
尽管这些观察结果的主要目标是黄蜂-18 B的次要食物,但我们还捕获了放弃前和放置后基线期间其相位曲线的一部分。在观测过程中,行星旋转107°,揭示了有关其大气空间分布的大量信息。为了对行星通量进行建模,我们考虑了二阶谐波功能
其中0≤f(t,θ)≤1是轨道参数θ和P的函数的行星投影磁盘的时间依赖性,可见的分数是轨道周期。谐波由一个描述行星通量变化的半值得fn的术语以及谐波达到最大值的时间TN。使用蝙蝠侠Python软件包68模型的归一化二级光曲线计算可见分数。二阶谐波函数为JWST69的噪声层提供了足够的精确度。我们适合决定日食的形状和持续时间的轨道参数:上等连接的时间TSEC(U [2459802.78,2459802.98]),影响参数B(n [0.360,0.0262]),以及半大轴对stellar gaus axis axelar gaus garius a/radius a/r a/r r*(n [3.496,0.0292])。考虑到半高轴和冲击参数的正常先验使用从苔丝获得的约束的中值和1σ不确定性(请参阅扩展数据表1)作为先验的中心和标准偏差。使用这些先验,因为从几个苔丝的转运中获得了精确的约束,并且与我们的次级折叠光曲线相反,它们没有与潜在的不均匀日期的相关性。我们选择保持A/R*和B的自由,而不是将它们固定到苔丝值中,以确保从白色光曲线拟合中检索到的其他参数在苔丝的不确定性上被边缘化。我们假设轨道是圆形的,这是通过TESS分析证明的,因为它发现强贝叶斯信息标准(BIC)和Akaike Information Criterion71(AIC)对非式轨道的偏好。考虑到WASP-18B与其宿主星的紧密距离,强烈的潮汐相互作用导致行星及其宿主的椭圆形变形。过去的研究表明,WASP-18B的这种变形约为2.5×10-3 rp,导致阶统一PPM的通量变化, 因此16,72忽略不计。最后,在第一阶(0.85μm)的下端附近,也对WASP-18B反射的恒星光中观察到的通量也有贡献。但是,几何反照率的上限(ag< 0.048 (ref. 73) and Ag = 0.025 ± 0.027 (ref. 74)) obtained from TESS correspond to a reflected-light contribution of <35 ppm near 0.8 μm. We therefore do not consider a specific term for the reflected light and instead assume that this will be fit by the second-order harmonic function.
We consider three phenomena that can lead to temporal changes in the observed stellar flux: stellar activity, A, ellipsoidal variations, E, and Doppler boosting, D. Stellar activity, generally caused by the presence and movement of star spots on the stellar hemisphere visible to the observer, leads to variations in the observed stellar spectrum on a timescale that is on the order of the stellar rotation period. Past observations of WASP-18 have constrained this period to be Prot ≈ 5.5 days (ref. 75). Despite this relatively short period with respect to stars of similar physical properties (for example, the effective temperature, Teff, and luminosity, L)76, the star shows abnormally low activity in the ultraviolet and X-ray domains, possibly because of tidal interactions with WASP-18b disrupting its dynamo77,78,79. As we expect the rotational modulation to be on a relatively long timescale compared with our observations and its amplitude to be low, we do not directly fit this term and instead assume it to be handled by the systematics model. For the ellipsoidal variation and Doppler boosting, they are both caused by the influence of WASP-18b on its host star. Although the ellipsoidal deformation of WASP-18b leads to a negligible impact on the phase curve, the same is not true for its host. The stellar ellipsoidal effect, with maxima fixed at quadrature when the projected area is at its highest, is found to be of semiamplitude 172.2 ppm from the TESS analysis. Following previous analyses of HST observations16, we consider ellipsoidal variation to be achromatic and fix its amplitude to the TESS value for the full NIRISS/SOSS wavelength range. The Doppler boosting effect is a result of the Doppler shift of the stellar spectral energy distribution as its radial velocity varies throughout its orbit. We fix this amplitude to 21.8 ppm, with the maximum at phase 0.25, as carried out in ref. 16. The observed stellar flux is therefore described as the sum of the ellipsoidal variation and Doppler boosting to the mean stellar flux in time .
The white and spectrophotometric light curves show two distinct important systematics: a sudden drop in flux caused by a tilt event shortly after the beginning of full eclipse, as well as high-frequency variations in the flux throughout the observations caused by small changes in the morphology of the trace. We track the trend of these systematics throughout the observations by performing incremental principal components analysis (PCA) with the open-source scikit-learn80 package on the processed detector images (Extended Data Fig. 2). The first principal component is the tilt event, which we use to determine the exact integration in which it occurs. We handle the tilt event in the white and spectrophotometric light curves by fitting for an offset in flux for the data after the 1,336th integration. We also observe two principal components analogous to the y position and full width at half maximum (FWHM) of the trace in time. We find that these two components are correlated to short-frequency variations in the light curves and therefore detrend linearly against them at the fitting stage. We find that, despite having a lower variance than the y position, the variation of the FWHM has a larger effect on the light curves when using box extraction. Finally, we fit for a linear trend in time to account for any further potential stellar activity and instrumental effect.
We note that considering a second-order polynomial trend or higher in time for the systematics results in notable correlation with the partial-phase-curve information. However, a linear-trend systematics model has been found sufficient to fit previous NIRISS/SOSS observations59,81. Furthermore, we find that the curvature around the secondary eclipse increases monotonically with wavelength, which is expected from the planetary signal and inconsistent with stellar activity as well as instrumental effects.
Light-curve fitting is performed on the extracted spectrophotometric observations using the ExoTEP framework22. With the orbital parameter values constrained from the white-light curve (see Extended Data Table 1), we solely fit for the planetary flux and systematics for each spectrophotometric light curve. The retrieved values of a/R* and b are within 1σ of the TESS values and such deviations do not affect the retrieved thermal emission spectrum as it is insensitive to the orbital parameters compared with transmission spectroscopy. We further explore the impact of our uncertainties on the retrieved map in the ‘Secondary-eclipse mapping’ section and find it to be robust against variations of the orbital parameters. We chose a resolution of 5 pixels per bin for our spectrum, corresponding to 408 spectrophotometric bins, to mitigate potential correlation between wavelengths in the atmospheric retrievals, as pixels in the spectral direction are not independent. All fits for the 408 bins are then done independently to obtain the planetary flux at secondary eclipse for all bins. Light-curve fits are performed using the affine-invariant Markov chain Monte Carlo (MCMC) ensemble sampler emcee82, using 20,000 steps and four walkers per free parameter. The first 12,000 steps, 60% of the total amount, are discarded as burn-in to ensure that the samples are taken after the walkers have converged. The samples from the white and spectrophotometric light curves are used to produce two science products: the detrended white-light curve and dayside thermal emission spectrum. The detrended white-light curve is obtained by dividing out the best-fit systematics model and subtracting the stellar variability from the light curve to isolate the planetary signal. For the dayside thermal emission, the median values and uncertainties of Fp(Tsec, λ) are computed from the samples of the parameters of equation (4).
During the TESS Primary Mission, the WASP-18 system was observed in sectors 2 and 3 (22 August to 18 October 2018). The full-orbit phase curve was analysed in several previous publications73,83, which reported a robust detection of the planet’s secondary eclipse and high signal-to-noise measurements of the planet’s phase-curve variation and signals corresponding to the ellipsoidal distortion and Doppler boosting of the host star. During the continuing TESS Extended Mission, the spacecraft reobserved WASP-18 in sectors 29 and 30 (26 August to 21 October 2020). We carried out a follow-up phase-curve analysis of all four sectors’ worth of TESS data, following the same methods used previously.
The data from the TESS observations were processed by the Science Processing Operations Center (SPOC) pipeline, which yielded near-continuous light curves at a 2-min cadence. As well as the raw extracted flux measurements contained in the simple aperture photometry light curves, the SPOC pipeline also produced the pre-search conditioning light curves, which have been corrected for common-mode instrumental systematics trends that are shared among all sources on the corresponding detector. We used these pre-search conditioning light curves for our phase-curve analysis. After dividing the light curves into individual segments that are separated by the scheduled momentum dumps of the spacecraft, we fit each segment to a combined phase curve and systematics model. The astrophysical phase-curve model consists of two components describing the planetary and stellar fluxes:
The Doppler boosting D and ellipsoidal distortion E semiamplitudes are defined as before. Here the orbital phase of the planet is defined relative to the mid-transit time T0: ϕ = 2π(t − T0). The phase-curve contribution of the planet has a single mode with a semiamplitude of Fatm and oscillates around the average relative planetary flux . The transit and secondary-eclipse light curves ϕt and ϕe were modelled using batman with quadratic limb darkening. In this parametrization, the secondary-eclipse depth (that is, dayside flux) and nightside flux are and , respectively. For the systematics model, we used polynomials in time and chose the optimal polynomial order for each segment individually that minimized the BIC.
We used ExoTEP to calculate the posterior distributions of the free parameters through a joint MCMC fit of all light-curve segments. To reduce the number of free parameters in the joint fit, we first carried out fits to smaller groups of light-curve segments corresponding to each TESS sector and then divided the light-curve segments by the best-fit systematics model. In the final joint fit of the systematics-corrected TESS light curve, no further systematics parameters were included. We accounted for time-correlated noise (that is, red noise) by fixing the uncertainty of all data points within each sector to the standard deviation of the residuals, multiplied by the fractional enhancement of the average binned residual scatter from the expected Poisson noise scaling across bin sizes ranging from 30 min to 8 h (ref. 83). As well as the phase-curve parameters described above (, Fatm, D and E), we allowed the mid-transit time T0, orbital period P, relative planetary radius Rp/R*, scaled orbital semi-major axis a/R*, impact parameter b and quadratic limb-darkening coefficients u1 and u2 to vary freely.
The results of our TESS phase-curve fit are listed in Extended Data Table 1. The revised values for the orbital ephemeris and transit shape parameters are statistically consistent with the published results from the previous TESS phase-curve analyses73,83, while being substantially more precise. We used the median and 1σ uncertainties of a/R* and b as Gaussian priors in the NIRISS/SOSS white-light-curve fit. We also used the median values of P, Rp/R*, D and E as fixed parameters for the NIRISS/SOSS white and spectrophotometric light curves fit. We obtained a secondary-eclipse depth of 357 ± 14 ppm and a nightside flux that is consistent with zero. All three phase-curve amplitudes were measured at high signal-to-noise ratio: Fatm = 177.5 ± 5.7 ppm, D = 21.8 ± 5.2 ppm and E = 172.2 ± 5.6 ppm.
To perform eclipse mapping, we use both ThERESA40 and the methods in ref. 39, which are separate implementations of the same process introduced in ref. 38 when applied to white-light curves, to cross-check our results. First, we generate a basis set of light curves from spherical harmonic maps84,85 with degree lmax and then transform these light curves to a new, orthogonal basis set of ‘eigencurves’ using PCA. Each eigencurve corresponds to an ‘eigenmap’, the planetary flux map that, when integrated over the visible hemisphere at each exposure time, generates the corresponding eigencurve.
We then fit the white-light curve with a linear combination of a uniform-map light curve, the N most informative (largest eigenvalue) eigencurves and a constant offset term to adjust for the fact that the observed planetary flux during eclipse (when the planet is entirely blocked by the star) must be equal to 0 and to allow for adjustments to light-curve normalization. Because the eigenmaps represent differences from a uniform map, it is possible to recover a fit that contains regions of negative planet emission. This is physically impossible, so we impose a positivity constraint on the total flux map. Although the eigenmaps are mathematically defined across the entire planetary sphere, our observations only constrain the portion of the planet that is visible during the observation, so we only enforce the flux-positivity condition in the visible region of the planet. Although this positivity condition could introduce a Lucy–Sweeney86 bias near zero flux, we note that our fitted maps (and GCM predictions) are far from negative, and increasing this boundary to, for example, 300 K leads to no change in the results. We test all combinations of lmax ≤ 6 and N ≤ 8 using a least-squares minimization and select the optimal values by minimizing the BIC.
We find that the fit with the lowest BIC to the broadband light curve has lmax = 5 and N = 5. However, the fit with lmax = 2 and N = 5 was only slightly less preferred, so here we explore the inferred brightness distribution from both solutions. Extended Data Fig. 8 shows the resulting light curve for the lmax = 5, N = 5 solution after sequential subtraction of each eigencurve. The preference for a fit with five eigencurves is driven by the residuals in ingress and egress, which can be seen by eye and are not sufficiently corrected for with a uniform map. Including the uniform-map light curve and the constant term, the fit thus contained a total of seven free parameters. We used a MCMC procedure to estimate parameter uncertainties. For the analysis following ref. 39, we test for convergence of the MCMC by ensuring that the chain length is 50 times the autocorrelation timescale, whereas for the analysis using ThERESA40, we use the Gelman–Rubin convergence test87 and achieve values ≤1.00006.
The resulting weights of each eigencurve are then applied to the corresponding eigenmaps to generate a flux map of the planet. We convert the star-normalized flux map to brightness temperature by assuming that the planet is a blackbody and the star emits as a PHOENIX41 spectrum calculated with PyMSG88, both integrated over the NIRISS/SOSS throughput. We estimate temperature-map uncertainties by computing a subsample of maps from the MCMC posterior distribution and calculating 68.3%, 95.5% and 99.7% quantiles at each location, including the effects of uncertainties in planetary radius, stellar temperature and stellar log g.
Figure 4 shows the resulting broadband brightness temperature map for the lmax = 5, N = 5 case and longitudinal brightness temperature profiles for both the lmax = 5, N = 5 and lmax = 2, N = 5 cases, calculated by averaging meridian flux at each longitude weighted by cos(latitude)2. Furthermore, we compare the equatorial slices to predictions from several GCMs (see ‘GCMs’ section). We note that not all structures on a planetary map will leave an observable signature in a secondary-eclipse light curve. When comparing GCMs to secondary-eclipse maps, it is important to only compare the components of GCM maps that can be physically accessed with eclipse mapping. Therefore, we use the methods in ref. 89 to separate each GCM map into the ‘null space’, or components that are inaccessible to eclipse-mapping observations, and the ‘preimage’, or components that are accessible through mapping. Figure 4 compares the longitudinal temperature trends from only the preimage of each GCM to the observed map. We find that both map solutions agree on a steep gradient in temperature near the limbs, which is well matched by GCM predictions. Furthermore, both maps show a temperature distribution roughly symmetrical in longitude about the substellar point. However, the two maps disagree on the exact shape of the brightness distribution. The lmax = 5, N = 5 map shows an extended hot-plateau region of roughly constant temperature from −40° to +40° in longitude, whereas the lmax = 2, N = 5 map shows a more concentrated hotspot with a steady decrease in temperature away from the substellar point. As these maps both fit the data with similar BIC, the present data do not give us the necessary precision to determine which solution represents the true temperature distribution of WASP-18b. Future observations at higher precision may distinguish between these two modes of solutions.
To test our ability to constrain latitudinal structures, we performed two eclipse-mapping fits: the eigenmapping fit presented above and a fit in which the initial basis set of maps is a longitudinal Fourier series that is constant with latitude (see Extended Data Fig. 9). Both methods retrieve similar longitudinal temperature structures, with steep gradients near the limb and an extended hot plateau. However, the constant-with-latitude map is also able to fit the data well, indicating a lack of constraints on latitudinal features within the uncertainties on the data. This is not unexpected, as the relatively low impact parameter (b = 0.34) of WASP-18b results in a lower amount of latitudinal information contained in the secondary-eclipse signal. Further observations of WASP-18b, or of planets with higher impact parameter, may enable us to pull latitudinal signals out of the noise.
Our eclipse mapping assumes that the orbital parameters of the system are precisely known relative to data uncertainties, a safe assumption with Spitzer data38. With the JWST, data quality may be high enough that uncertainties on orbital parameters impart substantial uncertainty on the mapping results. To test this, we ran analyses with impact parameter, orbital semi-major axis and eclipse time fixed to values ±1σ. In some cases, this led to a ‘hotspot’ model such as the red one in Fig. 4 being preferred over a ‘plateau’ model such as the blue one, which is unsurprising given their similar statistical preference. However, all resulting maps were well within the uncertainties of one of the two models presented. We note that, although the eccentricity is kept fixed to zero throughout this analysis, as justified by the preference for a circular orbit in the TESS analysis, considering an eccentric orbit would allow to the first order for variations in mid-eclipse time and eclipse duration90. Past photometric and radial-velocity observations of WASP-18b have found a small but non-zero eccentricity for WASP-18b on the order of e = 0.008 (refs. 91,92,93), corresponding to an offset of the time of mid-eclipse of 9 s, as well as a difference of 120 s between the transit and eclipse durations. These differences in eclipse timing and duration are of the same magnitude as those induced when varying Tsec, a/R* and b (8 s for the mid-eclipse time and 90 s for the eclipse duration). Therefore, performing the mapping considering a circular orbit while varying Tsec, a/R* and b is analogous to effects that could be expected from an eccentric orbit.
We perform 1D atmospheric retrievals on the NAMELESS reduction at a resolution of 5 pixels per bin (408 bins) using four techniques with varying levels of physical assumptions: a self-consistent radiative–convective–thermochemical equilibrium grid retrieval (ScCHIMERA), a chemical equilibrium with free temperature–pressure profile retrieval (SCARLET), a free-chemistry retrieval with thermal dissociation (HyDRA) and a free-chemistry retrieval with abundances assumed constant with altitude (POSEIDON). None of the retrievals used here consider the presence of clouds, as the dayside of giant exoplanets are expected to be cloudless at the equilibrium temperature Teq = 2,429 K (ref. 24) of WASP-18b (ref. 94). All retrieval methods considered the same PHOENIX stellar spectrum23, produced using previously published parameters for WASP-18 (that is, Teff = 6,435 K, log g = 4.35 and [Fe/H] = 0.1 (refs. 24,95)), to convert from model planet flux spectra to Fp/Fs values. We chose to use a model stellar spectrum instead of the extracted spectrum to avoid the possible introduction of systematic errors through the process of absolute flux calibration.
We use a 1D radiative–convective–thermochemical equilibrium (1D-RCTE) grid-retrieval-based method, ScCHIMERA7,96, with the opacity sources described in refs. 2,97. These 1D-RCTE models assume cloud-free 1D-RCTE using the methods described in ref. 98 to solve for the net flux divergence across each layer of the atmosphere and the Newton–Raphson iteration scheme in ref. 99 to march towards an equilibrium vertical temperature structure. The NASA Chemical Equilibrium with Applications 2 (CEA2) routine100 is used to compute the thermochemical equilibrium gas and condensate mole fractions for hundreds of relevant species. Opacities are computed with the correlated-k random-overlap-resort-rebin method101. Input elemental abundances from ref. 102 are scaled to a given metallicity ([M/H]) and carbon-to-oxygen (C/O) ratio.
Using the planetary and stellar parameters of WASP-18b, we produced a grid of 2,730 1D-RCTE models and resulting top-of-atmosphere thermal emission spectra spanning the atmospheric C/O (0.01–2.0), [M/H] (−2.0–2.0, for which 0 is solar, 1 is 10× etc.) and heat redistribution (f, 1.0–2.8, in which 1 = full, 2 = dayside and 2.67 is the maximum value, as defined in ref. 16). We also include a scale factor, AHS (allowed to vary from 0.5 to 2.0), that multiplies the planetary flux by a constant to account for a hotspot area fraction emitting most of the observed flux51. The PyMultiNest103 routine is used to sample the 1D-RCTE spectra through interpolation (and subsequent binning to the data wavelength bins) to obtain posterior probability constraints on the above parameters. We have made public our grid models, including temperature–pressure profiles, molecular abundances and emission spectra, as well as extra figures showing the posteriors of retrieved parameters and the impact of each parameter on the spectrum. This can be found on Zenodo at https://doi.org/10.5281/zenodo.7332105.
We use the SCARLET atmospheric-retrieval framework22,104,105,106 to perform a chemical-equilibrium retrieval with a free temperature–pressure profile on our retrieved dayside thermal emission spectrum. The SCARLET forward model computes the emergent disk-integrated thermal emission for a given set of molecular abundances, temperature structure and cloud properties. The forward model is then coupled to the affine-invariant MCMC ensemble sampler emcee82 to constrain the atmospheric properties. Owing to the high dayside temperature and large pressures examined through thermal emission spectroscopy, we assume that the atmosphere is in thermochemical equilibrium. For the equilibrium chemistry, we consider the following species: H2, H, H− (refs. 107,108), He, Na, K (refs. 109,110), Fe, H2O (ref. 111), OH (ref. 112), CO (ref. 112), CO2 (ref. 112), CH4 (ref. 113), NH3 (ref. 114), HCN (ref. 115), TiO (ref. 116), VO (ref. 117) and FeH (ref. 118). The abundances of these species are interpolated in temperature and pressure using a grid of chemical-equilibrium abundances from FastChem 2 (ref. 119), which includes the effects of thermal dissociation for all the species included in the model. These abundances also vary with the atmospheric metallicity, [M/H] (U[−3, 3]), and carbon-to-oxygen ratio, C/O (U[0, 1]), which are considered as free parameters in the retrieval. For the temperature structure, we use a free parametrization120 that here fits for N = 10 temperature points (U[100, 4,500] K) with fixed spacing in log-pressure (P = 102–10−6 bar). Although this parametrization is free, it is regularized by a prior punishing for the second derivative of the profile using a physical hyperparameter, σs, with units of kelvin per pressure decade squared (K dex−2). This prior is implemented to prevent overfitting and non-physical temperature oscillations at short pressure-scale lengths. For this work, we use a hyperparameter value of σs = 1,000 K dex−2, corresponding to a low punishment against second derivatives, as we want the retrieval to explore freely the temperature–pressure profile parameter space. We note that further lowering this punishment does not affect the retrieved temperature–pressure profile. Finally, we fit for an area fraction AHS (U[0, 1]) that is multiplied directly with the thermal emission spectrum, for a total of 14 free parameters. This factor is used to compensate for the presence of a hotspot that, although taking up only a portion of the planetary disk, contributes almost completely to the observed emission51. For the retrieval, we use four walkers per free parameter and consider the standard chi-square likelihood for the spectrum fit. We run the retrieval for 25,000 steps and discard the first 15,000 steps, 60% of the total amount, to ensure that the samples are taken after the walkers have converged. Spectra are initially computed using opacity sampling at a resolving power of R = 15,625, which is sufficient to simulate JWST observations121, convolved to the instrument resolution and subsequently binned to the retrieved wavelength bins.
Owing to the large mass of WASP-18b of Mp = 10.4 MJ, an important amount of rocky and icy material can be accreted without markedly changing the overall metallicity of the planet. As a zeroth-order estimate, we assume that the planet formed with exactly solar metallicity. Then, the mass of metals accreted needed to increase the overall metallicity to N times solar is given by ZMp(N − 1), in which Z = 0.0134 (ref. 25) is the solar metal mass fraction. From this relation, and assuming that the envelope of the planet is well mixed, we relate our retrieved metallicity probability posterior to the mass of metals accreted.
We quantify the impact of the stellar spectrum considered for the analysis on the retrieved atmospheric properties by running the same retrieval while varying the stellar spectrum. First, we explore the impact of the PHOENIX stellar-model parameters by varying them within their 1σ uncertainties (50 K for Teff, 0.05 for log g and 0.1 for [M/H]). We find these variations to have minimal impact on the retrieved metallicity with the measured median varying at most by 0.04 dex (about 0.15σ). The same is true for the C/O upper limit, with all retrieved upper limits being within 0.06 of the retrieval with the standard stellar parameters, with the exception of the Teff = 6,385 K case, which retrieves C/O < 0.64 at 3σ but does not affect our conclusions on the formation and migration history of WASP-18b. Second, we test the impact of the type of stellar model considered by also running the retrieval with an ATLAS9 stellar model122 (Teff = 6,435 K, log g = 4.35, [M/H] = 0.1). We find that the effect on the results are similar to those observed when varying the stellar parameters within their 1σ uncertainty, with a retrieved metallicity measurement of and C/O 3σ upper limit of 0.603. Finally, we test the effect of using the flux-calibrated spectrum on the retrieved atmospheric properties. The use of a flux-calibrated spectrum, measured directly from the NIRISS/SOSS observations, was avoided owing to some slight issues found in the CRDS reference files available at present used in the photom step of the jwst pipeline stage 2. The most recent reference file, photom_0034, is able to reproduce accurately the continuum from the PHOENIX model considered in the main retrieval but shows substantial noise in the observed spectrum. We also looked at reference file photom_0037, which was produced from ground data and does not account for the larger-than-expected throughput that was observed on sky20. Despite this, we perform a retrieval on the flux-calibrated stellar spectrum obtained by smoothing the response curve of reference file photom_0034 with a median filter of width 100. The retrieval ran on the flux-calibrated stellar spectrum retrieves a metallicity [M/H] of and a C/O 3σ upper limit of 0.739.
We also quantify the impact of the choice of reduction on the retrieved atmospheric properties by performing the same retrieval on the four spectra shown in Extended Data Fig. 3. We find that all reductions retrieve metallicities that are within 1σ of the NAMELESS reduction, with [M/H] values of , and for the nirHiss, transitspectroscopy and supreme-SPOON reductions, respectively. We also retrieve C/O 3σ upper limits of 0.749, 0.602 and 0.627 in that same order. We note that the slightly higher metallicity retrieved from the supreme-SPOON reduction is most probably because of the downward slope longwards of 2 μm, possibly caused by dilution of the signal through the process of background subtraction or 1/f correction, which is not observed in the nirHiss and transitspectroscopy reductions. We also find that the nirHiss retrieves larger uncertainties on the measured [M/H] and C/O, caused by the larger scatter at short wavelengths, which is possibly introduced through the optimal extraction process, as this effect is not seen in the reductions using box extraction.
We use two independent atmospheric-retrieval codes to perform free-chemistry retrievals on the dayside thermal emission spectrum of WASP-18b: y11,26,123,124, including the effects of thermal dissociation, and 125,126, which here assumes constant-with-depth abundances for all chemical species.
y consists of a parametric forward atmospheric model coupled to a Python implementation of the MultiNest127 nested-sampling Bayesian parameter-estimation algorithm128, PyMultiNest103. The inputs to the parametric model are the deep-atmosphere mixing ratios for each of the chemical species included, the temperature–pressure profile parameters (six free parameters) and a dilution parameter (area fraction) to account for 3D effects on the dayside51. Given the high dayside temperatures of WASP-18b, we consider high-temperature opacity sources and the effects of thermal dissociation, as described in ref. 11. We include opacity resulting from the molecular, atomic and ionic species with spectral features in the 0.8–2.8 μm range, which are expected in a high-temperature, H2-rich atmosphere: collision-induced absorption owing to H2–H2 and H2–He (ref. 129), H2O (ref. 112), CO (ref. 112), CO2 (ref. 112), HCN (ref. 130), OH (ref. 112), TiO (ref. 116), VO (ref. 117), FeH (ref. 131), Na (ref. 110), K (ref. 110) and H− (refs. 107,108). The line-by-line absorption cross-sections for these species are calculated following the methods in ref. 132, using data from the references listed. The H− free-free and bound-free opacity is calculated using the methods in refs. 107,108, respectively. y includes the effects of thermal dissociation of H2O, TiO, VO and H− as a function of pressure and temperature, following the method in ref. 9. In particular, the abundance profiles of these species are calculated following equations (1) and (2) in ref. 9, in which the deep abundance of each species (A0) is a parameter in the retrieval and the α, β and γ parameters are those given in Table 1 of that same work. The abundance profiles of the remaining chemical species are assumed to be constant with depth.
y uses the parametric temperature–pressure profile in ref. 133, which has been used extensively in exoplanet atmosphere retrievals, including ultra-hot Jupiters such as WASP-18b (ref. 11). The temperature parametrization is able to capture thermally inverted, non-inverted and isothermal profiles, spanning the range of possible thermal structures for ultra-hot Jupiters. The y retrievals also include an area fraction parameter, AHS, which multiplies the emergent emission spectrum by a constant factor to account for the dominant contribution of the hotspot51. Given the input chemical abundances, temperature–pressure profile parameters and area fraction, the model thermal emission spectrum is calculated at a resolving power of R ≈ 15,000, convolved to the instrument resolution, binned to the data resolution and compared with the data to calculate the likelihood of the model instance. Detection significances are calculated for specific chemical species by comparing the Bayesian evidences of retrievals, which include/exclude the species in question105. These detection significances factor in the ability of the retrieval to fit the observations with a different temperature profile and/or other chemical species, when the species in question is not included. Because thermal emission spectra are very sensitive to the atmospheric temperature profile, changes in the temperature structure can, in some cases, slightly compensate for the absence of a particular chemical species, contributing to a lower detection significance.
We also use y to test the sensitivity of the free-chemistry retrievals to the limits of the log-normal prior distributions assumed for the chemical abundances. For the y retrieval including thermal-dissociation effects, we test two scenarios: wide, uninformative priors for all 11 species included (log mixing ratio ranging from −12 to −1) and slightly more restricted priors for the refractory species included (log mixing ratio ranging from −12 to −4 for TiO, VO and FeH, −12 to −2 for Na and K and −12 to −1 for the remaining species). The more restricted prior limits for the refractory species are motivated by the relatively lower abundances expected for these species compared with the volatile species, across a range of metallicities and C/O ratios134,135,136.
We find that the atmospheric properties retrieved with y are consistent within 1σ for the two choices of prior limits. With the wide priors, the y retrieval infers a H2O log mixing ratio of , with double-peaked posterior probability distributions for the H2O and TiO abundances. Although the dominant posterior peaks correspond to approximately solar H2O and TiO abundances, the second, lower-likelihood peaks correspond to about 30 times and about 104 times supersolar H2O and TiO abundances, respectively. Such an extreme TiO abundance warrants scepticism and may, for example, be a result of the well-known degeneracy between chemical abundances and the atmospheric temperature gradient (see also ref. 5). The retrieved H2O abundance in this case is consistent with the chemical-equilibrium retrievals and self-consistent 1D radiative–convective models described above, although with a larger uncertainty owing to the double-peaked posterior distribution. When the restricted priors are used, the low-likelihood, high-abundance peaks in the H2O and TiO posterior distributions are no longer present, and the retrieved H2O abundance is , in excellent agreement with the chemical-equilibrium retrievals and self-consistent 1D radiative–convective models. We note that the retrieved temperature–pressure profiles and abundances of CO, CO2, HCN, OH, H− and FeH are unaffected by the choice of prior limits discussed above. The abundances of Na and K are unconstrained in both cases. Although the two choices of prior limits give consistent results, the expectation of chemical equilibrium in the atmosphere of WASP-18b, as well as the unlikelihood of an approximately 104 times supersolar TiO abundance, motivate the use of the restricted priors on the refractory chemical abundances.
We note that, for either choice of prior, the retrieved deep-atmosphere abundance of VO is substantially higher than that inferred by the chemical-equilibrium retrieval (Extended Data Fig. 5). This is because of a difference in the thermal-dissociation prescriptions; in the y retrieval, thermal dissociation results in a markedly depleted VO abundance at the photosphere, whereas in the chemical-equilibrium retrieval, thermal dissociation begins at lower pressures. Furthermore, the posterior distribution for the VO abundance peaks at highly supersolar values (about 104 times solar). Such a high abundance is physically unlikely and may indicate the presence of further sources of optical opacity not included in the retrieval.
We also conduct free-chemistry retrievals, without the inclusion of thermal dissociation, using . is an atmospheric-retrieval code originally designed for the interpretation of exoplanet transmission spectra125. We have recently extended to include secondary-eclipse emission spectra modelling and retrieval capabilities. For an ultra-hot Jupiter such as WASP-18b, in which the dayside can be assumed clear, the emission forward model of calculates the emergent flux by means of a standard single-stream prescription without scattering
in which Ilayer bot and Ilayer top are, respectively, the upwards specific intensity incident on the lower layer boundary and the intensity leaving the upper layer boundary, dτlayer is the differential vertical optical depth across the layer, μ = cosθ specifies the ray direction and B(Tlayer, λ) is the blackbody spectral radiance at the layer temperature. Using the boundary condition Ideep(μ, λ) = B(Tlayer, λ), propagates equation (8) upwards to determine the emergent intensity at the top of the atmosphere. The emergent planetary flux is determined through
in which W are the Gaussian quadrature weights corresponding to each μ (taken here as second-order quadrature over the interval μ [0, 1]). Finally, the planet–star flux ratio seen at Earth is given by
in which Rp,phot is the effective photosphere radius137,138 at wavelength λ (evaluated at τ(λ) = 2/3). Because the calculation of Rp,phot requires a reference radius boundary condition to solve the equation of hydrostatic equilibrium, we prescribe r(P = 10 mbar) as a free parameter.
For the WASP-18b retrieval analysis, we calculated emission spectra through opacity sampling and explored the parameter space using MultiNest through its Python wrapper PyMultiNest103,127. solves the radiative transfer on an intermediate-resolution wavelength grid (here, R = 20,000 from 0.8 to 3.0 μm), onto which high-spectral-resolution (R ≈ 106), pre-computed cross-sections139 are downsampled. For WASP-18b, we consider the following opacity sources: H2–H2 (ref. 140) and H2–He (ref. 140) collision-induced absorption, H2O (ref. 111), CO (ref. 141), CO2 (ref. 142), HCN (ref. 115), H− (ref. 108), OH (ref. 143), FeH (ref. 118), TiO (ref. 116), VO (ref. 117), Na (ref. 144) and K (ref. 144). We prescribed uniform-in-altitude mixing ratios, defined by a single free parameter for each of the chemical species included. The PyMultiNest retrievals with use 2,000 live points and the six-parameter temperature–pressure profile133 outlined above in the description of y. accounts for the dominant contribution of the hotspot by prescribing the 10 millibar radius as a free parameter, which is subsequently converted into an equivalent AHS posterior by comparison with the white-light planet radius.
We note that both y and yield consistent retrieval results when thermal dissociation is not considered in the y retrievals. However, the inclusion of thermal dissociation results in substantially different retrieved H2O abundances (see Extended Data Fig. 5) and temperature–pressure profiles, which are in agreement with the chemical-equilibrium retrievals and self-consistent 1D radiative–convective models described above.
To further justify the use of chemical-equilibrium models in our analysis of WASP-18b, we produce a grid of disequilibrium chemistry forward models to assess whether disequilibrium effects might strongly shape our observations. We begin by calculating the atmospheric temperature–pressure structure of WASP-18b under radiative–convective equilibrium with the HELIOS145,146 radiative-transfer code. Next, we calculate altitude-dependent mixing ratios of chemical species under this temperature–pressure structure with the VULCAN147 1D chemical-kinetics code, using an N–C–H–O reaction network that includes ionization and recombination of Fe, Mg, Ca, Na, K, H and He. We use the current version of VULCAN (VULCAN2 (ref. 148)), which includes optional photochemistry and parametrizes the transport flux of chemical species with eddy diffusion, molecular diffusion, thermal diffusion and vertical advection. We revise this code to include the effect of photoionization. Finally, we generate emission spectra with the PLATON radiative-transfer code149 at the resolution and wavelength range of NIRISS/SOSS. Our PLATON emission spectrum calculations use the code branch that allows varying chemical mixing ratios as a function of altitude150. We modify PLATON to calculate bound-free and free-free H− opacity as a function of altitude; this alteration is necessary to assess whether disequilibrium abundance H− opacity could mute spectral features more strongly than predicted by equilibrium chemistry. Furthermore, we modify PLATON to accept higher-temperature (T >3,000 K)我们使用Helios-K Code151,152计算的不透明度文件。
对于我们的一组模型,我们将涡流扩散系数KZZ从107 cm -2 s -1到1013 cm -2 s -1改变,对于给定的模拟,将其保持恒定。我们在许多数量级的KZZ上进行此扫描,以了解化学可能对观察到的发射光谱产生的最大影响。尽管KZZ是垂直混合的有限描述符,预计会随着高度的函数而变化(例如,参考文献153),但我们假设我们的正向模型支架是该行星的预期垂直混合行为。如果我们近似垂直风速度V和大气尺度的高度h(参考文献154,155),则该声明是由我们的GCM进一步激励的。我们的运动学MHD GCM的几天平均kzz(请参阅“ GCMS”部分)约为108 cm-2 s-1,最大的天数kzz约为109 cm-2 s-1,良好,远布在我们的藤栅范围内。我们的模型网格还可以切换分子扩散和光化学的包含。作为包括光化学的输入,我们使用了适用于参考文献的WASP-18的出色频谱。156。光谱是通过从参考文献中加入合成光谱来构建的。157和紫外线HD 222368的紫外线测量值是300 nm的国际紫外线探险家。
我们发现,我们包含了不平衡化学效应,例如光化学,分子扩散,热扩散和涡流扩散 - 产生的光谱与假设化学平衡的光谱并不强烈差异。实际上,在化学平衡和化学不平衡光谱下产生的光谱产生的光谱之间的所有差异均小于10 ppm。预计这一协议是因为在低分辨率发射光谱检查的压力水平上的化学不被预测通过光化学或混合(例如,参考文献148)会强烈修改。此外,黄蜂18b的高温意味着大气中的化学时间尺度很短,使化学反应比相关的不平衡时间刻度更快。
总体而言,我们的化学不平衡前向模型网格表明,此处考虑的不平衡化学效应不会强烈影响Niriss/Soss波段中黄蜂18b的发射光谱。
将GCM的套件与检索到的日期频谱和日期温度图进行了比较。
我们使用SPARC/MITGCM53对黄蜂18B的3D大气结构进行建模。该模型使用MITGCM158和辐射转移方程在球形几何形状中求解了原始方程,并使用当前的1D辐射转移模型159求解。我们使用相关的K框架基于逐线不良的不良范围生成不相差160。我们的模型假设元素丰度161和化学平衡气相组成162。我们的模型自然考虑了热解离的效果9。我们使用了25 s的时间步骤,并在过去100天内进行了大约300个接地天的模拟。
我们包括通过瑞利(Rayleigh)拖动参数化的进一步的阻力来源,每个模型单个恒定时间尺度,该限制确定流动降低的效率。在整个行星气氛中,阻力是恒定的。我们在模型之间从τdrag= 103到106 s(有效的拖动)和一个具有τdrag=∞的无拖网模型之间的时间尺度变化。我们的阻力范围涵盖了从无阻力的,风循环箱到拖曳主导的循环的过渡。参考文献中更详细地描述了我们使用的特定WASP-18B模拟。16。
我们使用的第二个模型是具有修订的纠察辐射转移方案163的运动学MHD GCM(参考文献47中详细描述)。由于该行星的高重力,我们选择在T31的水平分辨率(对应于赤道的大约3°分辨率)的水平分辨率下,将行星从100 bar到10-4 bar建模。我们使用参考文献中描述的运动学MHD处方。44,已用于HD Jupiters HD 209458B和HD 189733B(参考文献164)的模型,以及超热木星WASP-76B(参考文献47,52)。这种阻力处方假设由内部发电机产生的全局偶极磁场。由于这种几何形状,我们的阻力时间尺度仅作为瑞利阻力术语,仅在东 - 最新动量方程(垂直于磁场线的影响流),并计算为
其中B是所选的全局磁场强度,ϕ是纬度,ρ是密度,η是磁性电阻率。这个时间尺度是在本地和经常计算的,使得时间表在整个模型中的数量级有多个数量级以上,并响应大气温度的变化。在τmag小于行星轨道的1/20的位置中,将最小的时间尺度截止(大约103 s)应用于数值稳定性。我们选择建模一系列磁场强度(0、5、10和20 g),因为其真实值尚不清楚。
赞 (3)
评论列表(3条)
我是言希号的签约作者“lejiaoyi”
本文概览: 我们使用无名,Nirhiss56,SupremeSpoon57和TranspectRosproscopy58管道对Niriss/Soss/Soss/Soss/Soss 13...
文章不错《超热木星黄蜂18B的宽带热发射光谱》内容很有帮助