VLT/FORS2 comparative transmission spectroscopy II: confirmation of a cloud-deck and Rayleigh scattering in WASP-31b, but no potassium?
We present transmission spectroscopy of the hot-Jupiter WASP-31b using FORS2 on the VLT during two primary transits. The observations cover a wavelength range of 400–840 nm. The light curves are corrupted by significant systematics, but these were to first order invariant with wavelength and could be removed using a common-mode correction derived from the white light curves. We reach a precision in the transit depth of 140 ppm in 15 nm bins, although the precision varies significantly over the wavelength range. Our FORS2 observations confirm the cloud-deck previously inferred using HST/STIS. We also re-analyse the HST/STIS data using a Gaussian process model, finding excellent agreement with earlier measurements. We reproduce the Rayleigh scattering signature at short wavelengths ( Å) and the cloud-deck at longer wavelengths. However, our FORS2 observations appear to rule out the large potassium feature previously detected using STIS, yet it is recovered from the HST/STIS data, although with reduced amplitude and significance ( ). The discrepancy between our results and the earlier STIS detection of potassium ( ) is either a result of telluric contamination of the ground-based observations, or an underestimate of the uncertainties for narrow-band features in HST/STIS when using linear basis models to account for the systematics. Our results further demonstrate the use of ground-based multi-object spectrographs for the study of exoplanet atmospheres, and highlight the need for caution in our interpretation of narrow-band features in low-resolution spectra of hot-Jupiters.
keywords:methods: data analysis, stars: individual (WASP-31), planetary systems, techniques: spectroscopic, techniques: Gaussian processes
Transit and radial velocity surveys have made remarkable progress in understanding the population of extrasolar planets in our Galaxy. Obtaining spectroscopy of them is the next step in understanding the chemistry and physical processes in their atmospheres. Transiting planets enable such measurements by temporally resolving planets from their host stars, rather than spatially separating light from the planet and star. One such technique, transmission spectroscopy, measures the apparent radius of a planet during primary transit as a function of wavelength. This is the altitude at which the planet becomes opaque to starlight, which depends on the opacity, and therefore the composition and physical structure of the planet’s atmosphere (Seager & Sasselov, 2000; Brown, 2001).
Space-based observations have led the way in our understanding of exoplanet atmospheres (e.g. Charbonneau et al., 2002; Pont et al., 2008; Huitson et al., 2012; Berta et al., 2012; Pont et al., 2013; Kreidberg et al., 2014; Nikolov et al., 2015; Sing et al., 2016), however ground-based observations are rapidly increasing in importance, with the adoption of multi-object spectrographs (MOS) to perform differential spectro-photometry (e.g. Bean et al., 2010; Bean et al., 2011; Crossfield et al., 2013; Gibson et al., 2013a, b; Jordán et al., 2013; Kirk et al., 2016; Lendl et al., 2016; Mallonn & Strassmeier, 2016; Stevenson et al., 2014). The FOcal Reducer and low dispersion Spectrograph (FORS2; Appenzeller et al., 1998) on the European Southern Observatory’s (ESO) Very Large Telescope (VLT) was the first to demonstrate this technique successfully (Bean et al., 2010), but its impact has been limited by significant systematic effects related to the non-uniformity of the Linear Atmospheric Dispersion Corrector (LADC). This component has been recently upgraded (Boffin et al., 2015, 2016), and Sedaghati et al. (2015) has since demonstrated the improved performance of the instrument for exoplanet spectroscopy.
Here we report on the use of FORS2 to observe the transmission spectrum of the hot-Jupiter, WASP-31b. This is part of a small survey to re-observe targets that we have already observed using the Hubble Space Telescope (HST) (Sing et al., 2016), and our first results for WASP-39b are already reported in Nikolov et al. (2016). We selected targets for which their spectra already show evidence for alkali absorption and scattering by aerosols that should be easily detectable from the ground. Our survey has two general aims: to test the performance of FORS2 by observing objects with known spectroscopic features; and to confirm the signals detected with HST. This is particularly important, as the field of exoplanet spectroscopy has been limited by our ability to model instrumental time-series systematics (e.g. Gibson et al., 2011; Gibson et al., 2012a), and observing planets with multiple telescopes gives us the opportunity to verify and refine our current results and methodology.
WASP-31b is an inflated hot-Jupiter discovered by Anderson et al. (2011), with a mass and radius of and , respectively. It orbits a late F-type dwarf () with a period of 3.4 days. The optical and near-infrared transmission spectra have already been observed by Sing et al. (2015) using the Space Telescope Imaging Spectrograph (STIS) and Wide Field Camera 3 (WFC3), revealing strong Rayleigh- and Mie-scattering by aerosols, and a strong K feature (but no Na), that perhaps indicates a sub-stellar Na/K abundance ratio, and has interesting implications for the formation and/or evolution of the planet. These features make this an excellent target for our FORS2 observations, which should be able to recover and confirm these features. This paper is structured as follows; in Sect. 2 we describe the observations and data reduction and in Sect. 3 we present our light curve analysis and extraction of the transmission spectra. Finally in Sects. 4 and 5 we present our results and conclusions.
2 FORS2 Observations
Two transits of WASP-31b were observed using the 8.2-m ‘Antu’ telescope (Unit Telescope 1) of the VLT with FORS2: a general purpose imager, (multi-object) spectrograph and polarimeter. Observations were taken on the nights of 2016 February 15 and 2016 March 3, as part of program 096.C-0765 (PI: Nikolov), and followed a similar observing strategy to that described in Nikolov et al. (2016). The first transit was observed using the grism GRIS600B (hereafter 600B). Science exposures covered 5.2 hours, with 266 exposures of 40s (except the first 2). The second transit was observed using the grism GRIS600RI (hereafter 600RI) in combination with the GG435 order blocking filter, and science exposures lasted 5 hours, with 325 exposures typically ranging from 25–30s, and were adjusted through the night due to varying conditions.
FORS2 has an imaging field-of-view of 6.86.8 arcminutes squared, and consists of two 2k 4k CCDs arranged with a small (few arcsecond) gap in the cross-dispersion axis. We observed the target (V11.7) and 5 comparison stars simultaneously in multi-object (MXU) mode, using a custom mask designed from FORS2 pre-imaging to accurately place the slits. Observations were taken in binning mode, giving a pixel scale of 0.25/pixel. For our analysis, only the brightest comparison star was used (V11.4), as the other stars were significantly fainter and did not affect the final light curves. An acquisition image is shown in Fig. 1, showing the arrangement of the CCDs, the target and comparison star, and the approximate positions and shape of the slits. The width of all slits was 22, and the lengths 50 and 65 for the target and comparison star, respectively. This resulted in seeing-limited resolution of for the first night (FWHM varied from 3–7 pixels, at ), and for the second night (FWHM2–3 pixels, at ). We also constructed a calibration mask matching the science mask, but with 1 slit widths. This was used to obtain flat fields that would more closely match the stellar PSF than the wide slit, and arcs with narrower features for more precise wavelength calibration.
Data were bias subtracted and flat fielded using the FORS2 pipeline; however, these procedures made little difference to the final results, and we decided to extract spectra directly from the raw frames. This is not surprising, given that stable spatial variations in the detector sensitivity will cancel out in the differential photometry, and there were no high spatial frequency variations in the flat field structure. Furthermore, it is impossible to properly flat-field the data using a single flat field, as this cannot account for the difference in the slit transmission (for the sky background) and the stellar PSF, nor changes in the shape and position of the stellar PSFs that arise due to pointing errors and seeing variation. In the end the FORS2 pipeline was only used for the wavelength calibration, which used arc lamp exposures taken with the calibration mask.
The spectra for the target and comparison stars were extracted using a custom pipeline written using IRAF
The cross-dispersion PSF was fitted for each pixel channel to obtain a time-series of the full-width at half-maximum (FWHM) and the -position of the stars (averaged over wavelength). We also determined movement in the dispersion direction () by cross-correlating the spectra using the H line for the blue grism and the O A band telluric feature for the red grism, after fitting a low-order polynomial to normalise the continua. The time-series spectra for each star were interpolated to the same wavelength scale using the -shifts determined from the cross-correlation, and also from cross-correlation between the target and comparison star, again using the same spectral features.
We then proceeded to construct differential light curves from the time-series spectra. First, a ‘white’ light curve was constructed for each transit by integrating the target and comparison stars’ fluxes over a broad wavelength range, as indicated in Fig. 2, and dividing the time-series of the target star’s flux by the comparison star’s flux. We avoided using the edges of the spectra, and where the signal-to-noise was lower, although this has negligible influence on our results. We also extracted ‘spectral’ light curves by integrating over narrower wavelength channels, again shown in Fig. 2, and discussed further in Sect. 3.2. The white light curves and spectral light curves are shown in Figs. 3, 4 and 5. Finally, for the spectral response functions we extract limb darkening parameters using the PyLDTk toolkit (Parviainen & Aigrain, 2015), which calculates limb darkening coefficients and uncertainties for a range of limb darkening laws using the spectral libraries of Husser et al. (2013). We used the stellar parameters and uncertainties for WASP-31b derived in the discovery paper (Anderson et al., 2011).
3.1 White Light Curve Analysis
The white light curves and the spectral light curves exhibit large scale (1 percent level) systematics for both transit datasets. Fortunately, the systematics are mainly ‘grey’; in other words to first order they do not vary with wavelength. We use this fact to correct the spectral light curves using a “common-mode correction” to high accuracy and therefore extract precise transmission spectra of the planet. A Gaussian process (GP) is used to model the systematics simultaneously with the deterministic transit model to untangle the signals. This is a similar procedure to that used in Gibson et al. (2013a) and Gibson et al. (2013b). The primary difference is that the systematics are larger, and we therefore use prior information on the transit parameters in order to determine the systematics signal to higher accuracy. A similar process is used for our FORS2 observations described in Nikolov et al. (2016); however, the systematics were significantly smaller in this case, and a parametric model was used to construct the common-mode correction.
The exact cause of the systematics is unknown, but they are most likely due to remaining spatial inhomogeneities in the instrument/telescope throughput early in the optical path, as the underlying systematic variations tend to be largest when the telescope is rotating quickly, and they are strongly common-mode. It is worth noting however, that the target and comparison stars are located on different detectors. The most-likely cause is the LADC, which is unlikely to have a completely uniform throughput following removal of the anti-reflection coating. This is the only component of the instrument that rotates relative to the field. Other possibilities include variable vignetting of the field-of-view as the field rotates relative to the telescope structure, or scattered light, although we see no relationship between the size of the systematic variations and the phase/position of the moon. We also found no obvious correlations between the systematics and the auxiliary measurements (i.e. postition/width of the spectral trace, rotator angle/speed), and therefore fit the systematics using a time-correlated model.
We proceed by fitting each of the white transit light curves using our GP model. GPs were introduced to model stochastic signals in exoplanet time-series in Gibson et al. (2012a), where the reader is referred for details. See Gibson et al. (2013a, b) for the application of GPs to similar datasets, and Gibson (2014) for a comparison to other commonly-used techniques for systematics analysis.
Our GP is simply the multivariate Gaussian probability distribution defined as:
Here, and are vectors containing the time and flux measurements, respectively. is the multivariate Gaussian, where the mean vector is described by a transit function with parameter vector . The covariance matrix is defined by a covariance function or kernel, with parameters , and each element defined as
where is a hyperparameter that specifies the maximum covariance (and therefore the amplitude of the systematics), is the time difference, is the inverse length scale (defining the typical length scale of stochastic variations in the time-series), is the Kronecker delta, and is the white noise component, assumed to be the same for all points in the time-series.
Our light curve model (i.e. the mean function) is calculated using the equations of Mandel & Agol (2002) using quadratic limb darkening (with parameters and ). We fitted for the central transit time (), the system scale (), the planet-to-star radius ratio (), the impact parameter (), and two parameters describing a linear function of time for the baseline (, ). The period was held fixed, and Gaussian priors were placed on , , and , using the values reported in Sing et al. (2015), or directly derived from them. We also placed Gaussian priors on the planet-to-star radius ratio, , by determining the weighted mean reported over the wavelength ranges covered by the two white light curves. The limb darkening parameters were constrained using Gaussian priors with the best fit values and uncertainties from PyLDTk. These are summarised in Tab. 1. This naturally constrains the white light curve parameters to previously derived values, and enables a more accurate recovery of the common-mode systematics. The downside is that the overall level of the transmission spectrum is not derived from our data; however, this will not affect the conclusions of this study.
The kernel hyperparameters (, and ) were variable in the fitting process; however, we fitted for and , which is equivalent to imposing priors of the form , and is the natural choice of parameterisation for scale parameters. For example, fitting for the log length scale () or log inverse length scale () are mathematically equivalent, i.e. this implies the same prior information, but this is not true when fitting for or directly. The length scale was also constrained to be no finer than the time-sampling, and no longer than the twice the full time-span using a truncated uniform prior in log space. This was to aid convergence of the fitting algorithm.
We fit for both of the white light curves simultaneously. To obtain the best fitting model we first optimised the posterior with respect to the transit and kernel hyperparameters using a differential evolution algorithm
We justified our use of the quadratic limb darkening law by comparing the light curves to those generated using the non-linear limb darkening law, using the limb darkening parameters generated from PyLDTk. We found that for the white light curves and spectroscopic channels, the difference in transit depth was typically , and conclude that our choice of limb darkening law has little impact on our study.
Once we determined the best-fitting parameters, we extracted the mean of the Gaussian process conditioned on the observed data, and separated the systematic component from the transit model in order to generate a systematics-only model. This includes the linear baseline model. The white light curves and their best fit models are shown in Fig. 3, along with their corresponding systematics models and residuals. We then proceed to model the spectroscopic light curves, using the systematics models derived here.
|Both||3.405886 days (fixed)|
3.2 Spectroscopic Light Curve Analysis
We first created low-resolution spectroscopic light curves, using uniform bins of width 150 Å, with the edges smoothed using a Tukey window of width 5 Å– i.e. a cosine lobe. This was to mitigate the effect of sharp edges on the spectral channels. These correspond to the blue and red response curves shown in Fig. 2, which are multiplied by the target spectrum to show the spectral response for each channel. We extracted 16 of these low-resolution spectral channels for the 600B grism, and 22 for the 600RI grism. The resulting spectroscopic light curves are shown in Figs. 4 and 5. Clearly, they exhibit strong systematics, the main component of which is invariant in wavelength, and similar in shape to the white light curve systematics. We proceed by first removing these common-mode systematics by dividing each spectroscopic light curve by the systematics models derived from the respective white light curve fits. We also remove high-frequency systematics by subtracting the residuals from the white light curves and their best-fitting models. Such high-frequency systematics are often observed in such datasets, e.g. Gibson et al. (2013a, b), and likely arise due to spatial variations in the atmospheric throughput (e.g. uneven cloud cover), or additional instrumental effects; however, it made little difference for these datasets. The results of this correction are also shown in Figs. 4 and 5.
The use of a common-mode correction enables us to recover higher-precision measurements of the relative planet-to-star radius ratio, and hence the transmission spectrum; however, it is important to note that it comes at the cost of (potentially) applying an offset to the overall depth of the transmission spectrum, as it assumes a single, best-fit systematics model. In other words our final posterior probability distributions are conditioned on a single common-mode correction, which has an associated uncertainty on the transit depth. The crucial (and reasonable) assumption is that the shape of the transmission spectrum does not change with different random draws of the systematics model. This is verified by finding consistent transmission spectra using different methods to find the common-mode correction, including excluding the priors from the white light curve fits. A full analysis would require marginalisation over transmission spectra generated using different common-mode corrections drawn from the posterior distribution of the white light curve fits, after correcting for offsets in the spectra; however, this is beyond the scope of this paper and is unlikely to impact our conclusions.
The uncertainty in the offset of the transmission spectrum is related to the uncertainty in the white light curve fits. In this case, we impose the overall depth of the transmission spectrum using our white light curve priors derived from HST/STIS data (Sing et al., 2015), as discussed in Sect. 3.1. The absolute transit depth could still be derived in principle if the common-mode correction is not applied, but we would be marginalising the transit parameters over the same systematics signal, which would result in inflated error bars, and a highly correlated transmission spectrum.
We proceed to fit the common-mode corrected light curves with exactly the same model used for the white light curves; however, we fix the central transit time, the system scale, and the impact parameter at the value determined from the white light curve fits. To be conservative, rather than use the uncertainties from PyLDTk as Gaussian priors, we arbitrarily increased these to have a standard deviation of 0.1. This did not have a major impact on our results. Again, we find the global maximum using a differential evolution algorithm. We then clip data points that are more than 4 from the best-fit GP model. This typically resulted in clipping only 1–2 points per transit, and the arbitrary nature of this procedure does not change our results. We then explore the posterior distribution using a Markov-Chain Monte-Carlo (MCMC) analysis, using the implementation described in Gibson (2014) and references therein. For each light curve we used 4 chains each of length 60,000, and discarded the first half of the chains. We checked for convergence using the Gelman and Rubin statistic (Gelman & Rubin, 1992). The best-fit models are shown in Figs. 4 and 5, and the derived values for the planet-to-star radius ratio are shown in Tab. 2 and plotted in Fig. 6. Finally, we compared the fitted white noise values to the uncertainties calculated from the photon noise (excluding scintillation), finding these were typically 5-20% larger than the theoretical calculations. However, this simple comparison neglects the contribution of the systematics to the error budget, which in this case depends on the covariance terms and .
|Wavelength||Radius ratio||Limb Darkening|
|Centre [Range] (Å)||c1||c2|
3.3 Sodium and Potassium
In addition to the low-resolution light curves described in the previous section, we also search for signs of Na and K absorption by constructing high-resolution light curves. In particular this was motivated by trying to confirm the prominent K signature observed by Sing et al. (2015). For Na, we used 30 Å bins centred at 5892.9 Å, and for K we used 75 Å bins centred at 7681.5 Å. This corresponds to the same bins used in Sing et al. (2015). We also extracted light curves in the neighbouring continuum regions using the same bin widths, two at each side of each feature, in order to ensure the systematics were of similar amplitude. This resulted in sets of 5 higher-resolution light curves centred around the Na and K features. Both Na and K are covered by the 600RI grism, and only Na is covered by the 600B grism. These narrow spectral bins are shown in Fig. 2, and the corresponding light curves in Fig. 7. The systematics follow the same functional form as for the wider channels, and we fitted the light curves using the same procedure as before, using the common-mode correction.
The central Na (binned for both grisms) and K features are plotted in Fig. 6. No evidence for Na or K is found in our dataset. This is surprising as Sing et al. (2015) reported significant K absorption, finding a planet-to-star radius ratio of 0.13340.0020, which should be easily detectable using our FORS2 data. We therefore decided to explore this further by constructing differential light curves around the K feature. First, we simply took the light curves for the central K band, and divided it through by the sum of the other 4 continuum light curves (i.e. those plotted in Fig. 7). This is shown in the left of Fig. 8. We calculated the differential light curve depth from the Sing et al. (2015) transmission spectrum, using the central K channel and the neighbouring 4 channels, finding . The dashed line in Fig. 8 shows this signal over-plotted on the differential light curve, indicating that the signal should easily be detectable by eye in our light curves. The baseline (solid line) was determined using a GP fit to the light curve, using the same kernel as before, and a quadratic baseline in time as the mean function. The quadratic baseline was preferred by both the Bayesian Information Criterion (BIC; Schwarz, 1978) and the Akaike Information Criterion (AIC; Akaike, 1974). To set an upper limit on the absorption, we also fitted for a differential transit as a mean function to the GP (we neglected the effects of differential limb darkening in this analysis). We found an eclipse depth of , indicating no detection of K, and we can place a 3 upper limit of on the K absorption above the continuum. In addition, the AIC and BIC both favoured the model without an eclipse.
We also constructed a differential light curve ignoring the comparison star, i.e. by using only the raw flux from WASP-31 centred on the K feature, and the 4 neighbouring raw channels as reference to correct for the Earth’s atmosphere. The light curve is shown on the right of Fig. 8, again alongside a GP fit and the projected differential eclipse assuming the K signal from Sing et al. (2015). There is clearly a larger variation in the baseline, due to the differential throughput of the Earth’s atmosphere, particularly prominent as the K channel spans an O telluric feature. Nonetheless, a visual inspection appears to rule out excess K absorption from WASP-31b. We performed a similar fit to determine an eclipse depth, finding , and setting a 3 upper limit on the K absorption of . Again, the AIC and BIC favoured the model without an eclipse.
These results appear to rule out the large K absorption inferred from HST/STIS. However, it is important to emphasise that the K feature falls on a strong O telluric feature that complicates detection of K from the ground. We note that other studies using FORS2 have reported problems extracting information around this telluric feature (Sedaghati et al., 2016). In addition, our light curves suffer from significant instrumental systematics across the full spectral range, and we are assuming that the systematics are consistent for different narrow-band spectroscopic channels, and therefore cancel in the differential light curves. The interpretation of this result is discussed further in Sect. 4.
|Wavelength||Radius ratio||Limb Darkening|
|Centre [Range] (Å)||c1||c2|
3.4 Re-analysis of HST/STIS data
In order to resolve this contradiction in the K feature, and to produce a joint transmission spectrum, we re-analysed the HST/STIS data of WASP-31b that was first presented in Sing et al. (2015). These data were analysed using linear basis models to account for the systematics, using model selection via the BIC to choose the ‘correct’ model. While this has proven to be a reliable technique for STIS analyses (e.g. Sing et al., 2016; Nikolov et al., 2016), this method has also been shown to be the wrong statistical approach for determining transit parameters when faced with multiple systematics models, and consequently may underestimate the contribution of the choice of systematics models to the final transmission spectrum (Gibson, 2014). We therefore reanalyse the data using a Gaussian process model to account for this additional source of uncertainty in the transmission spectrum, to test if this is a possible explanation for the strong K feature detected from HST/STIS, yet absent in our FORS2 data.
The STIS data consist of three transits of WASP-31b: two observed with the G430L grating, and one with the G750L grating. We refer the reader to Sing et al. (2015) for further details, and comment only on the differences in the reduction applied here.
We started with the flt images, which were processed using the latest version of the STIS pipeline (calstis v3.4). This performed bias and flat-field corrections. The G750L data show significant fringing at long wavelengths, but we chose not to implement a fringing correction, as this is stable in time and therefore should not impact the light curves as long as the pointing is stable. Given the long exposure times (278 seconds), it was important to remove cosmic rays from the images. This was preformed by treating each pixel in the STIS array as a time-series. Each time-series was first divided through by the time-series averaged along each row (dispersion direction) of the array in order to remove the transit signal, and normalised. Outliers were identified by those more than 5 standard deviations from the median, and replaced by the value of a fitted interpolating function,
Spectral extraction was performed using a similar procedure to the FORS2 data, using the same aperture of 13 pixels as used by Sing et al. (2015), after background subtraction. We extracted the same diagnostic information as for the FORS2 data. However, we chose not to align the spectra to correct for movement in the dispersion direction, as the shift was sub-pixel. We checked that this didn’t influence our final light curves. The white light curves were extracted as before, using a wide pass-band spanning the majority of the wavelength regions observed with STIS. The spectral light curves were extracted using the same spectral bins as Sing et al. (2015) to enable easy comparison with these results, and also using the same bins as the FORS2 data to enable a combined analysis.
Rather than modelling the systematics purely as time-correlated noise, we also considered auxiliary inputs that describe the transient state of the observational setup. This is standard procedure for space based datasets (e.g. Brown et al., 2001; Gilliland & Arribas, 2003; Pont et al., 2007; Sing et al., 2011; Gibson et al., 2011; Evans et al., 2013; Nikolov et al., 2014; Sing et al., 2016). This is easy to perform using GPs, and we simply replace the time-dependent kernel used for the FORS2 data with the multidimensional kernel previously used in Gibson et al. (2012a); Gibson et al. (2012b):
and again specify the maximum covariance and noise terms, respectively. is the th element of the th auxiliary input, and are inverse length scales, one for each of the input dimensions.
As inputs to the GP kernel, we used the orbital phase of HST, and the - and -positions of the spectral trace on the detector. These were included as they were the same inputs used by Sing et al. (2015). We also tested including the width and angle of the spectral trace, and found that they made negligible difference to the final spectra. We again chose to fit for and . Fits were performed in the same way as before, first optimising the (hyper-)parameters using a global optimiser, and then using an MCMC to explore the joint posterior distribution. We again used four chains of length 60,000, and discarded the first half of the chains.
We first fitted the white light curve, allowing the kernel hyperparameters (, , ) to vary along with , and . The values of , and were fixed at those determined by Sing et al. (2015). The systematics component for each of the white light curves was isolated from the transit model, and used for the common-mode correction.
The spectroscopic light curves were first divided through by the common-mode correction. We then fitted each of them using the same model, determining the planet-to-star radius ratio for each. The results are plotted in Fig. 9, using the same bins as the Sing et al. (2015) analysis, and in Fig. 10, using the same bins as the FORS2 data. Tab. 4 lists the results. The fitted white noise values for the spectroscopic light curves were on average higher than the calculated photon noise for the two G430L transits, and higher for the G750L transit, consistent with that reported in Sing et al. (2015). Again, we note that this neglects the contribution of systematics to the uncertainties.
|Wavelength||Radius ratio||Limb Darkening|
|Centre [Range] (Å)||c1||c2|
WASP-31b is an inflated hot-Jupiter with mass and radius of and , respectively, a surface gravity 4.56 m/s, and a zero-albedo equilibrium temperature of 1580 K (Anderson et al., 2011). This results in a scale height of 1220 km, or 0.0014 . Sing et al. (2015) used nightly photometry over 3 years to monitor the activity of WASP-31, concluding that the effects of stellar variability are negligible in the interpretation of WASP-31b’s transmission spectrum.
Using the HST/STIS data described previously, Sing et al. (2015) detected a Rayleigh scattering signature at short wavelengths ( nm), a cloud-deck (flat spectrum) at longer wavelengths, and a prominent K feature at 768 nm. Following Sing et al. (2015), we interpret our transmission spectrum by fitting with a two component model including scattering at short wavelengths, and a flat spectrum at longer wavelengths. At short wavelengths we assume the scattering is dominated by a species with a power-law absorption cross-section with index : . Lecavelier Des Etangs et al. (2008) showed that this results in a constant gradient with respect to , given by
where is the temperature of the atmosphere, is the mean molecular weight (assumed to be 2.3 times the proton mass), is the surface gravity, is the Boltzmann constant, and is the planet’s radius. At longer wavelengths we assume a flat spectrum, defined by the planetary radius given by the scattering law at a transition wavelength, . Sing et al. (2015) found K, and a transition wavelength nm.
4.1 FORS2 Transmission Spectrum
Our FORS2 transmission spectrum of WASP-31b is shown in Fig. 6. The broad-band spectrum is consistent with the results presented in Sing et al. (2015). We discuss the non-detection of K in the Sect. 4.4. The uncertainties range from a factor of 2 lower than our STIS results at the centre of the FORS2 grisms (corresponding to ), to a factor of 2 larger at the edges of the grisms (). This reflects the importance of the common-mode correction to our FORS2 results, which is most accurate in the central wavelength bands, as well as the fact that the sensitivity falls off at the edges of the wavelength coverage, in particular at the bluest wavelengths covered by the 600B grism.
We use a Levenburg-Marquardt algorithm
We also fit the spectrum using a linear fit, and a flat spectrum. The model with the scattering signature and cloud-deck is preferred by the AIC (AIC = 0.1), whereas the flat model (cloud-deck only) is preferred by the BIC ( BIC = 1.2). Our FORS2 data therefore do not distinguish clearly between these two models. We also performed the same fits fixing the transition wavelength to 510 nm, i.e. that determined by Sing et al. (2015). In this case both the AIC and BIC marginally prefer the two-component model (AIC = 1.1, BIC = 0.4), although again this is not a significant distinction.
The Rayleigh scattering slope is therefore only tentatively detected by the FORS2 data alone. This is unsurprising, as the sensitivity of FORS2 at the bluest wavelengths (450 nm) drops off sharply
4.2 STIS Transmission Spectrum
Our results from the re-analysis of the STIS data using GP models are shown in Fig. 9. They closely match those reported in Sing et al. (2015). The uncertainties are typically 10% larger for the G430L grism, and 30% larger for the G750L grism. This is explained by the extra flexibility allowed by the GP model, and the fact that it effectively marginalises over a large range of systematics models, rather than select a single one. Furthermore, our GP model assumes a joint systematics model of the inputs rather than independent signals summed together (as is the case for linear-basis models). Our model would be conceptually similar to adding cross-terms (e.g. , , etc) to the linear basis models, and marginalising over a large range of systematics models, although we note that our GP effectively marginalises over an infinite number of basis models (Gibson et al., 2012a; Gibson, 2014).
Qualitatively, the Rayleigh scattering signature at wavelengths Å, and the cloud deck signature at longer wavelengths are reproduced. The models from Sing et al. (2015) and Barstow et al. (2017) are also plotted in Fig. 9, and our broad-band results remain consistent with these models. We did not reproduce the strong 4.3 K feature, but still detect an excess absorption in the narrow K-band. Using the neighbouring 4 points as reference, we determine the excess absorption in the narrow K feature to be (2.5), indicating tentative evidence for K, but not a conclusive detection. This is likely due to the different systematics treatment presented here.
4.3 Combined Broad-band Transmission Spectrum
The final transmission spectrum of WASP-31b is shown in Fig. 10, after combining the FORS2 and STIS datasets using a weighted average where the bins overlap. The values are provided in Tab. 5. We again fitted the two-component scattering plus cloud-deck model to the spectrum, finding K, and a transition wavelength nm. These values are consistent with Sing et al. (2015), although we measure a slightly shallower gradient at short wavelengths. The two-component model is strongly preferred over the flat model (AIC = 6.8, BIC = 6.2), indicating that our results confirm the interpretation of Sing et al. (2015), that of a scattering aerosol dominating the opacity at short wavelengths, with a cloud-deck dominating at longer wavelengths. The model provides a good fit to the data, giving with 33 degrees of freedom. Note that the narrow Na and K bins were not included in the fit.
Assuming the equilibrium temperature of 1580 K, we find , consistent with a Rayleigh scattering signature (). Alternatively, assuming Rayleigh scattering, we find a temperature of K, consistent with the equilibrium temperature, and also with temperatures where condensates might form in the atmosphere ( K, e.g. Fortney, 2005).
4.4 The Potassium Feature
We also searched for Na and K absorption in the FORS2 data, using narrow bands centred on the respective doublets, as discussed in Sect. 3.3. We did not find evidence for Na or K detection in our FORS2 dataset, which is inconsistent with the 4.3 detection of K reported in Sing et al. (2015). The differential light curves discussed in Sect. 3.3 and presented in Fig. 8 provide further evidence that the K absorption is not present in our FORS2 data at the level measured by Sing et al. (2015). The interpretation of this is complicated by the K feature falling on a large O telluric feature. In principal the comparison star should account for varying levels of telluric absorption in the atmosphere as the telescope tracks the target; however, it is often the case that systematics are larger both in telluric features, and where there are sharp boundaries in the target or comparison stars’ spectra. It is nonetheless difficult to imagine a scenario where unaccounted-for systematics hide such a prominent K feature in the differential light curves, but do not appear elsewhere in the light curve – in other words, the systematics would have to be closely correlated to the transit shape, and approximately match the amplitude of the K feature.
The O feature also modifies the shape of the K bin, but again this is unlikely to dilute the signal enough to hide the K. For broad K absorption from the planet this is apparent in the shape of the narrow spectral channels in Fig. 2. However, it is still possible that the K feature is concentrated in the narrow cores of each line of the doublet, that is not resolved in our spectral bins. To test this hypothesis we used the ESO SkyCalc tool (Noll et al., 2014), to compute the telluric absorption spectrum at high resolution for the extremes of airmass covered by our observations. It is possible that one of the lines of the doublet (but not both) overlaps with a deep telluric feature with absorption , depending on the exact barycentric correction. However, the resulting telluric features are extremely narrow, with FWHM Å. Significant dilution of the K signal would require the K lines to exactly match up with telluric features, and be concentrated in similarly narrow cores. This scenario seems unlikely, as the STIS K feature is detected using a 75 Å bin, and therefore such narrow features would be diluted by a large factor when integrating over a wide spectral bin, requiring a much larger planet-to-star radius ratio than reported at low-resolution. Even if this was the case for the one line in the doublet, this would dilute the signal by which would still be detectable in our data. Furthermore, assuming the planet is tidally locked, the cores are likely to be wider than Å from rotational broadening alone.
We also re-analysed the HST/STIS data revealing an excess K absorption feature of ( ). This is of lower amplitude and with a marginally larger uncertainty than Sing et al. (2015), who reported (). Given the lower statistical significance, the K feature from the STIS data could be re-interpreted as a statistical fluctuation in the data, and therefore consistent with the lack of detection found with FORS2. The original detection could therefore be dismissed as an underestimate of the uncertainties for narrow band features. Nonetheless, our re-analysis does still provide tentative evidence for excess K absorption, and it remains suspicious that the one place our ground-based data disagree with the original STIS analysis is close to a deep telluric signal. Without a more detailed understanding of the telluric absorption (which is dependent on the unknown shape of the planet’s K feature) we cannot completely dismiss it conspiring to hide the K, either from additional time-series systematics, dilution of the signal due to the narrow K cores overlapping with telluric lines, or indeed a combination of both.
Our combined transmission spectrum results in a K detection of , therefore does not confirm the interpretation of excess K absorption found in Sing et al. (2015), although neither can it rule it out. The interpretation of our results largely depends on the choice between the Linear Basis Model or GP approach, and on the effect of the O telluric feature on our ability to detect K using our FORS2 data, which we cannot definitively quantify. We therefore cannot completely rule out K absorption at the level reported by Sing et al. (2015), although our FORS2 results and STIS re-analysis do raise questions about this interpretation.
|Wavelength||Radius ratio||Limb Darkening|
|Centre [Range] (Å)||c1||c2|
We have presented a ground-based transmission spectrum of the hot-Jupiter WASP-31b using two transits to cover the wavelength range nm, using the recently upgraded FORS2 instrument on the VLT. This is the second paper in a series to re-observe targets already studied in detail with HST, both to confirm the signals detected with HST, and to verify the use of FORS2 for transmission spectroscopy. The results for our first target, WASP-39b, were presented in Nikolov et al. (2016). We also presented a re-analysis of the HST/STIS transmission spectrum of WASP-31b presented in Sing et al. (2015).
Our FORS2 observations suffer from significant systematic effects, showing percent level variations in both the broad-band and spectroscopic light curves. The exact cause(s) of the systematics is unknown, but the most likely culprit is the LADC, or some other spatial variation in the throughput early in the optical path. Luckily, the observed systematics are predominantly common-mode, i.e. they are invariant with wavelength to first order. We use this fact to extract a systematics model from the broad-band (‘white’) light curves, and use this to correct the spectroscopic light curves prior to extracting the transmission spectrum.
The FORS2 observations are consistent with the HST/STIS results and confirm the presence of a cloud-deck in the atmosphere of WASP-31b, owing to the flat spectrum across most of the optical range, and the lack of broad Na and/or K pressure broadened wings. We do not unambiguously detect the presence of Rayleigh scattering from the FORS2 observations alone. We also cannot reproduce the large K absorption reported in Sing et al. (2015), either using high-resolution channels in the transmission spectrum, or from differential light curves. As the K feature falls within a strong telluric O band, the FORS2 data do not conclusively rule out a strong K core in WASP-31b.
We also re-analysed the HST/STIS data using a Gaussian process model, finding consistent results with Sing et al. (2015), only with slightly larger uncertainties. Our combined analysis confirms the overall picture from Sing et al. (2015), that of a Rayleigh scattering haze dominating the opacity at short wavelengths, and a cloud-deck dominating at longer wavelengths. Our GP model detects excess absorption around the K feature, although with a larger uncertainty and slightly lower amplitude, resulting in a reduced detection significance of 2.5 , and 2.2 when combined with the FORS2 data.
The combined FORS2 and STIS transmission spectrum therefore only shows tentative evidence for K absorption, and the low-significance detection from the STIS data could be interpreted as a statistical fluctuation. We therefore cannot confirm the sub-solar Na/K abundance inferred by Sing et al. (2015). However, the interpretation of our narrow-band results depends on the choice of systematics models applied to the STIS data (Linear-basis models vs Gaussian processes), and on the impact of the O telluric absorption on the K feature, which is difficult to definitively quantify. Our results therefore raise some doubts on the interpretation of narrow-band features in low-resolution spectra of hot-Jupiters. Nonetheless, the overall agreement of the broad-band FORS2 and STIS spectra, combined with our previous results from Nikolov et al. (2016), is a powerful demonstration of the stability of both instruments, and in our current interpretation of scattering signatures in hot-Jupiters (e.g. Sing et al., 2016). Resolving the discrepancy of the K feature will require further observations, preferably at higher-resolution.
This work is based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere under ESO programme 096.C-0765. N. P. G. gratefully acknowledges support from the Royal Society in the form of a University Research Fellowship. N. N, D. K. S, and T. M. E. acknowledge funding from the European Research Council under the European Unions Seventh Framework Programme (FP7/2007-2013) / ERC grant agreement no. 336792. J. K. B. is supported by a Royal Astronomical Society Research Fellowship. P.A.W. acknowledges the support of the French Agence Nationale de la Recherche (ANR), under program ANR-12-BS05-0012 ‘Exo-Atmos’. We are grateful to the developers of the NumPy, SciPy, Matplotlib, iPython and Astropy packages, which were used extensively in this work (Jones et al., 01; Hunter, 2007; Pérez & Granger, 2007; Astropy Collaboration et al., 2013).
- pubyear: 2016
- pagerange: VLT/FORS2 comparative transmission spectroscopy II: confirmation of a cloud-deck and Rayleigh scattering in WASP-31b, but no potassium?–VLT/FORS2 comparative transmission spectroscopy II: confirmation of a cloud-deck and Rayleigh scattering in WASP-31b, but no potassium?
- IRAF is distributed by the National Optical Astronomy Observatory, which is operated by the Association of Universities for Research in Astronomy (AURA) under cooperative agreement with the National Science Foundation
- PyRAF is a product of the Space Telescope Science Institute, which is operated by AURA for NASA
- For a textbook level introduction see Rasmussen & Williams (2006), or Bishop (2006) for a more general introduction in the context of Bayesian inference.
- as implemented in the SciPy package, based on Storn & Price (1997)
- In this case we used a GP with a squared exponential kernel, but the exact nature of the interpolating function has little impact.
- as implemented in SciPy
- This could be improved by using the blue-sensitive CCD, but this is only available in visitor mode.
- Akaike H., 1974, IEEE Transactions on Automatic Control, 19, 716
- Anderson D. R., et al., 2011, \aap, 531, A60
- Appenzeller I., et al., 1998, The Messenger, 94, 1
- Astropy Collaboration et al., 2013, \aap, 558, A33
- Barstow J. K., Aigrain S., Irwin P. G. J., Sing D. K., 2017, \apj, 834, 50
- Bean J. L., Miller-Ricci Kempton E., Homeier D., 2010, \nat, 468, 669
- Bean J. L., et al., 2011, \apj, 743, 92
- Berta Z. K., et al., 2012, \apj, 747, 35
- Bishop C. M., 2006, Pattern Recognition and Machine Learning. Springer
- Boffin H., et al., 2015, The Messenger, 159, 6
- Boffin H. M. J., et al., 2016, \procspie, 9908, 99082B
- Brown T. M., 2001, \apj, 553, 1006
- Brown T. M., Charbonneau D., Gilliland R. L., Noyes R. W., Burrows A., 2001, \apj, 552, 699
- Charbonneau D., Brown T. M., Noyes R. W., Gilliland R. L., 2002, \apj, 568, 377
- Crossfield I. J. M., Barman T., Hansen B. M. S., Howard A. W., 2013, \aap, 559, A33
- Evans T. M., et al., 2013, \apjl, 772, L16
- Fortney J. J., 2005, \mnras, 364, 649
- Fortney J. J., Shabram M., Showman A. P., Lian Y., Freedman R. S., Marley M. S., Lewis N. K., 2010, \apj, 709, 1396
- Gelman A., Rubin D. B., 1992, Stat. Sci., 7, 457
- Gibson N. P., 2014, \mnras, 445, 3401
- Gibson N. P., Pont F., Aigrain S., 2011, \mnras, 411, 2199
- Gibson N. P., Aigrain S., Roberts S., Evans T. M., Osborne M., Pont F., 2012a, \mnras, 419, 2683
- Gibson N. P., et al., 2012b, \mnras, 422, 753
- Gibson N. P., Aigrain S., Barstow J. K., Evans T. M., Fletcher L. N., Irwin P. G. J., 2013a, \mnras, 428, 3680
- Gibson N. P., Aigrain S., Barstow J. K., Evans T. M., Fletcher L. N., Irwin P. G. J., 2013b, \mnras, 436, 2974
- Gilliland R. L., Arribas S., 2003, Instrument Science Report NICMOS 2003-001
- Huitson C. M., Sing D. K., Vidal-Madjar A., Ballester G. E., Lecavelier des Etangs A., Désert J.-M., Pont F., 2012, \mnras, 422, 2477
- Hunter J. D., 2007, Computing In Science & Engineering, 9, 90
- Husser T.-O., Wende-von Berg S., Dreizler S., Homeier D., Reiners A., Barman T., Hauschildt P. H., 2013, A&A, 553, A6
- Jones E., Oliphant T., Peterson P., 2001–, SciPy: Open source scientific tools for Python, http://www.scipy.org/
- Jordán A., et al., 2013, \apj, 778, 184
- Kirk J., Wheatley P. J., Louden T., Doyle A. P., Skillen I., McCormac J., Irwin P. G. J., Karjalainen R., 2016, preprint, (arXiv:1611.06916)
- Kreidberg L., et al., 2014, \nat, 505, 69
- Lecavelier Des Etangs A., Pont F., Vidal-Madjar A., Sing D., 2008, \aap, 481, L83
- Lendl M., et al., 2016, \aap, 587, A67
- Mallonn M., Strassmeier K. G., 2016, \aap, 590, A100
- Mandel K., Agol E., 2002, \apjl, 580, L171
- Nikolov N., et al., 2014, \mnras, 437, 46
- Nikolov N., et al., 2015, \mnras, 447, 463
- Nikolov N., Sing D. K., Gibson N. P., Fortney J. J., Evans T. M., Barstow J. K., Kataria T., Wilson P. A., 2016, \apj, 832, 191
- Noll S., Kausch W., Kimeswenger S., Barden M., Jones A. M., Modigliani A., Szyszka C., Taylor J., 2014, \aap, 567, A25
- Nortmann L., Pallé E., Murgas F., Dreizler S., Iro N., Cabrera-Lavers A., 2016, \aap, 594, A65
- Parviainen H., Aigrain S., 2015, MNRAS, 453, 3821
- Pérez F., Granger B. E., 2007, Computing in Science and Engineering, 9, 21
- Pont F., et al., 2007, \aap, 476, 1347
- Pont F., Knutson H., Gilliland R. L., Moutou C., Charbonneau D., 2008, \mnras, 385, 109
- Pont F., Sing D. K., Gibson N. P., Aigrain S., Henry G., Husnoo N., 2013, \mnras, 432, 2917
- Rasmussen C. E., Williams K. I., 2006, Gaussian Processes for Machine Learning. The MIT Press
- Schwarz G., 1978, Ann. Statist., 6, 461
- Seager S., Sasselov D. D., 2000, \apj, 537, 916
- Sedaghati E., Boffin H. M. J., Csizmadia S., Gibson N., Kabath P., Mallonn M., Van den Ancker M. E., 2015, \aap, 576, L11
- Sedaghati E., et al., 2016, \aap, 596, A47
- Sing D. K., et al., 2011, \mnras, 416, 1443
- Sing D. K., et al., 2015, \mnras, 446, 2428
- Sing D. K., et al., 2016, \nat, 529, 59
- Stevenson K. B., Bean J. L., Seifahrt A., Désert J.-M., Madhusudhan N., Bergmann M., Kreidberg L., Homeier D., 2014, \aj, 147, 161
- Storn R., Price K., 1997, Journal of Global Optimization, 11, 341