The operational methane retrieval algorithm for TROPOMI

Abstract

This work presents the operational methane retrieval algorithm for the Sentinel 5 Precursor (S5P) satellite and its performance tested on realistic ensembles of simulated measurements. The target product is the column-averaged dry air volume mixing ratio of methane (XCH${}_{\mathrm{4}}$), which will be retrieved simultaneously with scattering properties of the atmosphere. The algorithm attempts to fit spectra observed by the shortwave and near-infrared channels of the TROPOspheric Monitoring Instrument (TROPOMI) spectrometer aboard S5P.

The sensitivity of the retrieval performance to atmospheric scattering properties, atmospheric input data and instrument calibration errors is evaluated. In addition, we investigate the effect of inhomogeneous slit illumination on the instrument spectral response function. Finally, we discuss the cloud filters to be used operationally and as backup.

We show that the required accuracy and precision of $<$ 1 % for theXCH${}_{\mathrm{4}}$ product are met for clear-sky measurements over land surfaces and after appropriate filtering of difficult scenes. The algorithm is very stable, having a convergence rate of 99 %. The forward model error is less than 1 % for about 95 % of the valid retrievals. Model errors in the input profile of water do not influence the retrieval outcome noticeably. The methane product is expected to meet the requirements if errors in input profiles of pressure and temperature remain below 0.3 % and 2 K, respectively. We further find that, of all instrument calibration errors investigated here, our retrievals are the most sensitive to an error in the instrument spectral response function of the shortwave infrared channel.

Introduction

Methane (CH${}_{\mathrm{4}}$) is the most important anthropogenic greenhouse gas after carbon dioxide (CO${}_{\mathrm{2}}$). While it occurs in smaller concentrations, it has a higher global warming potential per molecule than CO${}_{\mathrm{2}}$. An accurate understanding of CH${}_{\mathrm{4}}$ sources and sinks is essential for a reliable prediction of climate change. Space-based measurements can provide continuous and global monitoring of CH${}_{\mathrm{4}}$, leading to much-needed improved constraints on the surface fluxes.

Several current and future satellite missions measure CH${}_{\mathrm{4}}$ abundances in the Earth's atmosphere. The SCanning Imaging Absorption spectroMeter for Atmospheric CartograpHY (SCIAMACHY) aboard ENVISAT () was the first space-based instrument to measure atmospheric CH${}_{\mathrm{4}}$ with sensitivity down to the Earth's surface. were the first to use these measurements to constrain CH${}_{\mathrm{4}}$ surface fluxes. After loss of contact with ENVISAT in 2012, the Greenhouse gases Observing SATellite (GOSAT) is currently the only satellite measuring atmospheric CH${}_{\mathrm{4}}$ . While GOSAT has a higher sensitivity and spatial resolution than SCIAMACHY, it has a fairly low spatial sampling. In 2017, the TROPOspheric Monitoring Instrument (TROPOMI) will be launched aboard the Sentinel 5 Precursor (S5P) satellite, and it will provide CH${}_{\mathrm{4}}$ measurements as one of its key products with unprecedented high precision, spatial resolution and global daily coverage.

The common goal of the above-mentioned missions is to provide atmospheric CH${}_{\mathrm{4}}$ concentrations with sufficient accuracy and spatiotemporal coverage to allow the assessment of CH${}_{\mathrm{4}}$ sources through inverse modelling. The observation strategy relies on measuring spectra of sunlight, backscattered by the Earth's surface and atmosphere, in the shortwave infrared (SWIR) spectral range. Absorption features of CH${}_{\mathrm{4}}$ molecules allow for retrieval of its atmospheric concentration with high sensitivity down to the Earth's surface where the main CH${}_{\mathrm{4}}$ sources are located. The applicability of such measurements for estimating source strengths, however, strongly depends on the precision and accuracy achieved. Residual systematic biases must be well below 1 % to facilitate inverse modelling .

Scattering by aerosols and cirrus is one of the major challenges for retrievals of CH${}_{\mathrm{4}}$ from space-based SWIR observations. While contamination by optically thick clouds can be filtered out reliably, optically thin scatterers are much harder to detect, yet they still modify the light path of the observed backscattered sunlight. This can lead to underestimation or overestimation of the true CH${}_{\mathrm{4}}$ column if not appropriately accounted for. The net light path effect strongly depends on the amount, size and height distribution of the scatterers as well as on the reflectance of the underlying surface . Therefore, retrieval strategies rely on inferring the target gas concentration either simultaneously with atmospheric scattering properties or with a light path proxy.

introduced the “proxy” approach for CH${}_{\mathrm{4}}$ retrieval from SCIAMACHY measurements around 1600 nm, by using the simultaneously retrieved CO${}_{\mathrm{2}}$ column as a light path proxy. The proxy approach relies on the assumptions that scattering effects cancel in the ratio of the CH${}_{\mathrm{4}}$ column and the CO${}_{\mathrm{2}}$ column, and that a prior estimate of the CO${}_{\mathrm{2}}$ column is sufficiently accurate to recalculate the CH${}_{\mathrm{4}}$ column from the measured CH${}_{\mathrm{4}}/$ CO${}_{\mathrm{2}}$ ratio. In this case scattering is ignored in the forward modelling. Further applications of the proxy approach for CH${}_{\mathrm{4}}$ retrieval from SCIAMACHY are described by and . For GOSAT, the proxy approach has been successfully applied by and .

Alternatively, scattering-induced light path modification can be taken into account by simultaneously inferring the atmospheric CH${}_{\mathrm{4}}$ concentration and physical scattering properties of the atmosphere. Such “physics-based” methods have been developed for space-based CO${}_{\mathrm{2}}$ and CH${}_{\mathrm{4}}$ measurements from SCIAMACHY, GOSAT and the Orbiting Carbon Observatory (OCO); see e.g. . The physics-based methods make use of the oxygen (O${}_{\mathrm{2}}$) A-band in the near-infrared (NIR) around 760 nm and absorption bands of the target absorber in the SWIR spectral range. The advantage of physics-based methods for CH${}_{\mathrm{4}}$ retrieval compared to proxy methods is that they do not depend on prior information on the CO${}_{\mathrm{2}}$ column. On the other hand, the physics-based algorithms are more complex and computationally expensive. In addition, they may be limited by the information content of the measurement with respect to aerosol properties and related forward model errors in the description of aerosols. A detailed comparison between the two methods for GOSAT is provided by .

TROPOMI has four spectral channels in the ultraviolet (UV), visible (VIS), near-infrared (NIR) and shortwave infrared (SWIR), with spectral ranges of 270–320, 310–495, 675–775 and 2305–2385 nm, respectively. We use TROPOMI measurements in the NIR and SWIR for CH${}_{\mathrm{4}}$ retrievals. This spectral range does not allow for a light-path-proxy approach, and thus the effect of aerosols and cirrus needs to be accounted for by using a physics-based method as described above. The goal of this paper is to present the CH${}_{\mathrm{4}}$ retrieval algorithm for TROPOMI and investigate its sensitivity to algorithm assumptions, atmospheric input data and instrument calibration errors and filtering criteria. To this end, we simulated realistic TROPOMI measurements for aerosol- and cirrus-loaded atmospheres under clear-sky and cloudy conditions.

The outline of this paper is as follows. We start with the methodology in Sect. , giving an overview of the instrument, the retrieval algorithm and the methane data product. Then, a detailed error sensitivity study is presented in Sect. , based on methane retrievals on a clear-sky global ensemble of simulated spectra. In Sect. , a study of cloud filtering is performed on a TROPOMI orbit of simulated spectra that cover a realistic range of cloud parameters. Conclusions are presented in Sect. .

Methodology

The S5P satellite has a designed 7-year lifetime and will fly in a sun-synchronous orbit at 824 km altitude. Its single payload, TROPOMI, is a push-broom imaging spectrometer with a wide swath of 2600 km and a ground pixel of 7 $\times$ 7 km${}^{\mathrm{2}}$ in exact nadir. Approximately, TROPOMI observes a full swath per second, which results in $\sim$ 216 spectra per second. The instrument comprises two spectrometer modules, the first containing the UV, VIS and NIR spectral channels, and the second dedicated to the SWIR channel. We use the NIR and SWIR channels with spectral resolutions of 0.38 and 0.25 nm and spectral sampling ratios of 2.8 and 2.5, respectively . Since the NIR and SWIR detectors are incorporated in different instrument modules, the NIR spectra will be co-registered with the SWIR spectra before performing CH${}_{\mathrm{4}}$ retrievals. Examples of simulated TROPOMI NIR and SWIR spectra are shown in Fig. .

Simulated TROPOMI spectra in the NIR and SWIR channels for a scenario from the global ensemble; see Sect. . The scenario, located in a rural area in Utah (USA), is observed at nadir and SZA $={50}^{\circ }$.

[Figure omitted. See PDF]

The retrieval algorithm

The S5P operational CH${}_{\mathrm{4}}$ retrieval algorithm is based on RemoTeC, which was orginally developed for CO${}_{\mathrm{2}}$ and CH${}_{\mathrm{4}}$ retrievals from OCO and GOSAT obervations . A first performance study on TROPOMI measurements has been done by , where a detailed description of the algorithm can be found. Here, we summarise the essentials and elaborate on the modifications made since then.

The forward model

The algorithm aims to infer the atmospheric state vector $\mathit{x}$ from spectral measurements $\mathit{y}$ in the NIR (757–774 nm) and SWIR (2305–2385 nm) ranges. This requires a forward model $\mathbf{F}$ that can accurately compute the measurement given the atmospheric state: $$\mathit{y}=\mathbf{F}\left(\mathit{x}\right)+{\mathit{e}}_y+{\mathit{e}}_{\mathrm{F}},$$ where ${\mathit{e}}_y$ is the measurement noise error, and ${\mathit{e}}_{\mathrm{F}}$ is the forward model error. The forward model incorporates the linearised radiative transfer model LINTRAN . LINTRAN simulates the radiance at the top of the atmosphere $I^{\mathrm{TOA}}$ on a fine internal spectral grid. The model needs a solar irradiance spectrum on the internal spectral grid as input, which is inferred from the daily solar measurements of TROPOMI using the deconvolution approach by and . The simulated radiance measurement is obtained by a spectral convolution of $I^{\mathrm{TOA}}$ with the instrument spectral response function (ISRF).

In addition to accuracy, an important requirement of the algorithm is computational speed. To be able to process TROPOMI's huge number of measurements, the CPU time per retrieval has to be on the order of seconds. Our algorithm achieves this by, among other things, using the linear k-method for multiple scattering . Single scattering is calculated line by line as it is computationally less expensive. To our knowledge, our algorithm is among the fastest in use for CH${}_{\mathrm{4}}$ full-physics retrievals, with an average CPU time of 7–10 s per retrieval .

For a given model atmosphere, the forward model simulates spectra of backscattered sunlight, taking absorption and scattering by molecules and particles into account. The model atmosphere is defined for 36 pressure-equidistant vertical layers, with the top at 0.1 hPa and bottom at the surface pressure $p_{\mathrm{surf}}$. We calculate $p_{\mathrm{surf}}$ by interpolating the meteorological surface pressure, from the European Centre for Medium-Range Weather Forecasts (ECMWF), on the surface elevation, from a digital elevation map (). The absorbing trace gases of interest are O${}_{\mathrm{2}}$ in the NIR band, and CH${}_{\mathrm{4}}$, H${}_{\mathrm{2}}$O and CO in the SWIR band. The first guess layer subcolumns of these gases are calculated from input profiles of CH${}_{\mathrm{4}}$ and CO (from the global chemical transport model TM5; ), and temperature, humidity and pressure profiles (from ECMWF forecast data). Molecular absorption features are calculated using appropriate spectroscopic databases . Here, we evaluate the absorption cross sections on a 72-layer pressure-equidistant grid to account for the strong temperature dependence and pressure dependence of the cross sections. Molecular scattering properties are given by Rayleigh theory. Particulate absorption and scattering are computed with Mie theory using tabulated aerosol properties by .

In our algorithm, the aerosol type is characterised by the refractive index and size distribution. The complex refractive index is fixed at $n=1.4-0.01i$ in the O${}_{\mathrm{2}}$ A-band and $n=1.47-0.008i$ in the SWIR band, which are averaged values derived from the ECHAM5-HAM model . The size distribution is described by a power-law function with size parameter $\mathit{\alpha }$ (e.g. ):

$$n\left(r\right)=\left\{\begin{array}{ll}A,& \text{if}r\le r_{\mathrm{1}}.\\ A(r/r_{\mathrm{1}}{)}^{-\mathit{\alpha }},& \text{if}{r}_{\mathrm{1}}<r\le r_{\mathrm{2}}.\\ \mathrm{0},& \text{if}r>r_{\mathrm{2}},\end{array}\right.$$ where $r_{\mathrm{1}}=0.1$ $\mathrm{\mu }$m, $r_{\mathrm{2}}=10$ $\mathrm{\mu }$m, $r$ is the particle radius and $A$ is a normalisation constant. The amount of aerosol and its vertical distribution is provided by the vertically integrated column number density $N_{\mathrm{aer}}$ and a normalised Gaussian height distribution $h\left(z_k\right)$ with $z_k$ the height of layer $k$: $$h\left(z_k\right)=B\exp \left(-,{\displaystyle \frac{\mathrm{4}\ln \mathrm{2}(z_k-z_{\mathrm{aer}}{)}^{\mathrm{2}}}{w^{\mathrm{2}}}}\right),$$ where $w=2000$ m, $z_{\mathrm{aer}}$ is the central height and $B$ is a normalisation constant. Thus, in layer $k$ with thickness $\mathrm{\Delta }{z}_k$, the layer subcolumn $n_k$ is given by $$n_k=N_{\mathrm{aer}}h\left(z_k\right)\mathrm{\Delta }{z}_k.$$ suggest that the exact choice of the fixed-value parameters, such as the refractive indices or the width of the height distribution, does not affect the retrievals significantly. To support this claim, we have studied the effect of different values for the fixed-value aerosol parameters on our retrievals, and we arrive at the same conclusion; see Appendix .

We have modified the retrieval forward model with respect to to account for chlorophyll fluorescence emission in a simplified manner as described by . In short, scattering of the fluorescence emission is ignored, and solely the spectral shape and absorption features by oxygen (O${}_{\mathrm{2}}$) are modelled. This allows fluorescence to be treated as a simple additive term to the radiance before convolution with the ISRF. The surface emission at the top of the atmosphere (TOA) is then modelled as $$F_{\text{s}}(\mathit{\lambda }{)}^{\text{TOA}}=F_{\text{s},755}^{\text{surf}}(\mathrm{1}-s(\mathit{\lambda }-755))e^{-{\mathit{\tau }}_{{\mathrm{O}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)/\mathit{\mu }},$$ where ${\mathit{\tau }}_{{\mathrm{O}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)$ is the vertical optical thickness of oxygen, $\mathit{\mu }$ is the cosine of the viewing zenith angle and $\mathit{\lambda }$ is the wavelength in nanometres. $F_{\text{s},755}^{\text{surf}}$ and $s$ represent the chlorophyll emission at 755 nm and its spectral slope, respectively.

The state vector

The state vector $\mathit{x}$ consists of 25 elements; see Table . The target is a 12-layer pressure-equidistant vertical profile of CH${}_{\mathrm{4}}$ partial column number densities. The 12 retrieval layers are related to 36-layer model atmosphere by the shared interfaces; i.e. each retrieval layer is divided in three sublayers for the forward model calculations. Retrieved ancillary parameters are the total column number densities of H${}_{\mathrm{2}}$O and CO, three scattering parameters $N_{\mathrm{aer}}$, $\mathit{\alpha }$ and $z_{\mathrm{aer}}$ (related to amount, size and height), the surface albedo (constant-term and first-order spectral dependence

Note that the algorithm has the option to fit higher-order spectral dependence of the surface albedo. It will be investigated whether this is needed with real data.

and in the NIR and SWIR band), two terms to account for a spectral shift of the measurement (NIR and SWIR band) and two terms to account for chlorophyll fluorescence emission, $F_{\text{s},755}^{\text{surf}}$ and $s$ (see Eq. ). An overview of the state vector elements and their a priori values is given in Table . We note that for the operational methane retrievals, the a priori value for $F_{\text{s},755}^{\text{surf}}$ will come from a dedicated fluorescence retrieval scheme based on the Fraunhofer lines. In this study however, we simply took an a priori $F_{\text{s},755}^{\text{surf}}=\mathrm{0}$ photons s${}^{-\mathrm{1}}$ m${}^{-\mathrm{2}}$ nm${}^{-\mathrm{1}}$. Further, we would like to mention that our retrieved CO column should be regarded as a diagnostic for the CH${}_{\mathrm{4}}$ retrieval. The official TROPOMI CO product has a dedicated retrieval algorithm described by .

State vector elements and their a priori values or references to source.

State vector element	A priori value (units)
CH${}_{\mathrm{4}}$ subcolumns in 12 vertical layers	TM5 (cm${}^{-\mathrm{2}}$)
CO total column	TM5 (cm${}^{-\mathrm{2}}$)
H${}_{\mathrm{2}}$O total column	derived from ECMWF specific humidity (cm${}^{-\mathrm{2}}$)
Aerosol column $N_{\mathrm{aer}}$	derived from aerosol optical thickness $=$ 0.1 at 760 nm (cm${}^{-\mathrm{2}}$)
Aerosol size parameter $\mathit{\alpha }$	3.5 (–)
Aerosol height parameter $z_{\mathrm{aer}}$	5000 (m)
Lambertian surface albedo in the NIR band	maximum of measured reflectance in the NIR band (–)
First-order spectral dependence surface albedo in the NIR band	0 (–)
Lambertian surface albedo in the SWIR band	maximum of measured reflectance in the SWIR band (–)
First-order spectral dependence surface albedo in the SWIR band	0 (–)
Spectral shift NIR	0 (nm)
Spectral shift SWIR	0 (nm)
Fluorescence emission at 755 nm $F_{\text{s},755}^{\text{surf}}$	0 (photons s${}^{-\mathrm{1}}$ m${}^{-\mathrm{2}}$ nm${}^{-\mathrm{1}}$)
Fluorescence spectral slope $s$	0 (–)

The inversion procedure

The state vector is found by inverting Eq. (), where the inverse method is based on a Philips–Tikhonov regularisation scheme . Regularisation is required because the inverse problem is ill-posed; i.e. the measurements $\mathit{y}$ typically contain insufficient information to retrieve all state vector elements independently. Philips–Tikhonov regularisation aims to reduce contributions from measurement noise to the retrieved state vector while retaining valuable information. Because the forward model $\mathbf{F}\left(\mathit{x}\right)$ is non-linear in $\mathit{x}$, the inversion is performed iteratively by a step-size-controlled Gauss–Newton scheme, where at each iteration step the forward model is linearised.

The inverse algorithm finds $\mathit{x}$ by minimising the cost function that is the sum of the least-squares cost function and a side constraint weighted by the regularisation parameter $\mathit{\gamma }$ according to $$\widehat{\mathit{x}}=\underset{x}{min}\left(\Vert ,{{\mathbf{S}}_y}^{-\mathrm{1}/\mathrm{2}},\left(\mathbf{F},\left(,\mathit{x},\right),-,\mathit{y}\right),|,{|}^{\mathrm{2}},+,\mathit{\gamma },\Vert ,\mathbf{W},(,\mathit{x},-,{\mathit{x}}_{\text{a}},\left),\right|,{|}^{\mathrm{2}}\right),$$ where ${\mathbf{S}}_y$ is the diagonal measurement error covariance matrix, which contains the noise estimate. ${\mathit{x}}_{\text{a}}$ is an a priori state vector, and $\mathbf{W}$ is a diagonal weighting matrix that renders the side constraint dimensionless and ensures that only the CH${}_{\mathrm{4}}$ parameters and the scattering parameters contribute to its norm: $W_{jj}=\mathrm{1}/x_{\mathrm{a},j}$ for the CH${}_{\mathrm{4}}$ column number densities and the three aerosol parameters, and $W_{jj}={10}^{-\mathrm{7}}/x_{\mathrm{a},j}$ for all other state vector elements. The latter are thus retrieved in a least-squares sense. For determining $\mathit{\gamma }$, the L-curve criterion is applied in the baseline algorithm. However, $\mathit{\gamma }$ may be fine-tuned later using real observations.

The CH${}_{\mathrm{4}}$ data product

Although we retrieve methane in $n=12$ sublayers, there is virtually no profile information in the measurement. The degree of freedom of signal of the retrieved methane profile is about 1. Therefore, the methane data product is given as a column-averaged dry air mixing ratio XCH${}_{\mathrm{4}}$. This quantity is obtained from the methane entries of the retrieved state vector $x_i$ through $${\mathrm{XCH}}_{\mathrm{4}}=\sum _{i=\mathrm{1}}^n{x}_i/V_{\mathrm{air},\mathrm{dry}},$$ where $V_{\mathrm{air},\mathrm{dry}}$ is the dry air column (calculated from meteorological input surface pressure and the water vapour profile).

To interpret the retrieved XCH${}_{\mathrm{4}}$ correctly, one also needs the column averaging kernel ${\mathbf{A}}_{\mathrm{col}}$ that describes the sensitivity of the retrieved CH${}_{\mathrm{4}}$ column to changes to the true methane profile (see , for details): $$A_{\mathrm{col},i}={\displaystyle \frac{\partial \sum _{i=\mathrm{1}}^n{x}_i}{\partial x_{\mathrm{true},i}}}.$$ As an example, in Fig.  we show the column averaging kernel corresponding to the methane retrieval performed on the simulated spectra of Fig. . We see that the column averaging kernel is around 1 in the lower atmosphere. From Eq. () it is clear that the closer ${\mathbf{A}}_{\mathrm{col}}$ is to 1, the more the retrieved column represents the true column (see also Eq. ). This illustrates that the retrieval of methane columns from the SWIR has a nearly ideal sensitivity to methane in the troposphere and the tropospheric boundary layer.

For validation and interpretation purposes it is important to realise that the retrieved XCH${}_{\mathrm{4}}$ is related to the true methane profile ${\mathit{x}}_{\mathrm{true}}$ and the a priori profile ${\mathit{x}}_{\mathrm{a}}$ as $$\begin{array}{ccc}{\mathrm{XCH}}_{\mathrm{4}}& =& \sum _{i=\mathrm{1}}^n\left(A_{\mathrm{col},i},x_{\mathrm{true},i},+,(,\mathrm{1},-,A_{\mathrm{col},i},),x_{\mathrm{a},i}\right)/V_{\mathrm{air},\mathrm{dry}}\\ & & +& \mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},\text{F}}+\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},y},\end{array}$$ where $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},\text{F}}$ is the bias caused by forward model errors, and $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},y}$ is the retrieval noise due to measurement noise. The standard deviation of the retrieval noise, i.e. the precision ${\mathit{\sigma }}_{{\mathrm{XCH}}_{\mathrm{4}}}$, follows from the error covariance matrix ${\mathbf{S}}_x$, which describes the effect of measurement noise on the retrieved parameters (see e.g. , Sect. 4.3): $${\mathit{\sigma }}_{{\mathrm{XCH}}_{\mathrm{4}}}={\displaystyle \frac{\sqrt{{\sum }_{i=\mathrm{1}}^n{\sum }_{j=\mathrm{1}}^n{S}_{x,i,j}}}{V_{\mathrm{air},\mathrm{dry}}}}.$$

Together with XCH${}_{\mathrm{4}}$ and ${\mathbf{A}}_{\mathrm{col}}$, the precision ${\mathit{\sigma }}_{{\mathrm{XCH}}_{\mathrm{4}}}$ is given in the main data product. Note that for simulations, the bias $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},\mathrm{F}}$ can be calculated from Eq. (), since all other terms are either known from the simulations or calculated by the retrieval algorithm. In Sect. , we evaluate $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},\mathrm{F}}$, $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},y}$ and their sensitivity to input errors.

Column averaging kernel for a typical methane retrieval. The retrieval was performed using the simulated NIR and SWIR spectra in Fig. .

[Figure omitted. See PDF]

Data filtering

Our algorithm has been designed to be efficient and accurate under certain assumptions; e.g. the atmosphere is cloud-free and plane-parallel with optically thin scatterers. Therefore, it is essential to filter out those scenes for which these assumptions break down to ensure data quality and reach the required accuracy. Moreover, if filtering can be done a priori, for example using cloud detection, we can considerably reduce the number of performed retrievals and hence computation time.

A priori filtering

For operational cloud filtering, measurements from the Visible Infrared Imaging Radiometer Suite (VIIRS) aboard the Suomi-NPP satellite will be used. S5P is foreseen to operate in loose formation with Suomi-NPP, meaning that both missions observe the same ground scene with a time delay of about 5 min. In case VIIRS data are not available, we have developed a backup cloud filter using H${}_{\mathrm{2}}$O retrievals from the weak and strong absorption band, assuming a non-scattering atmosphere. Here we make use of the fact that clouds and aerosols will modify the optical path length in the two bands differently due to the different absorption strengths . Thus, there will be a difference between the H${}_{\mathrm{2}}$O column retrieved from the strong band and the H${}_{\mathrm{2}}$O column retrieved from the weak band in the presence of clouds.

Based on tests with simulations, we chose 2329–2334 nm as the weak absorption band and 2367–2377 nm as the strong absorption band (see also ). In Fig. , we show the absorption features of H${}_{\mathrm{2}}$O in the SWIR range and highlighted the proposed weak and strong absorption band. Because the weak-band and the strong-band retrievals are differently affected by scattering, the ratio $\frac{\left|{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}_{\mathrm{weak}},-,{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}_{\mathrm{strong}}\right|}{{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}_{\mathrm{strong}}}$ is strongly correlated with cloud contamination. Scenes in which this ratio exceeds a certain threshold can be flagged as cloudy and filtered out. Based on a realistic ensemble of synthetic measurements (see Sect. ), we find that a threshold of 0.08 filters out most cloudy scenes and keeps most of the clear-sky scenes. However, a direct cloud filter based on VIIRS data is shown to be superior to this indirect two-band cloud approach. Therefore, the two-band cloud filter will only be used as backup when VIIRS data are unavailable.

Simulated reflectance spectrum showing the absorption features of H${}_{\mathrm{2}}$O in the SWIR spectral range for the same scenario as in Fig. . The red shaded area represents the weak absorption band, while the blue area represents the strong absorption band of H${}_{\mathrm{2}}$O. The difference between the water column retrieved from the weak and strong band is used for the backup cloud filter.

[Figure omitted. See PDF]

Further, before performing any CH${}_{\mathrm{4}}$ retrievals, we filter out cases with solar zenith angle (SZA) larger than ${70}^{\circ }$, and viewing zenith angle (VZA) larger than ${50}^{\circ }$. These thresholds have been derived from simulations and shall be fine-tuned after launch using real observations.

A posteriori filtering

identified an a posteriori filter based on retrieved scattering parameters: SWIR aerosol optical thickness ${\mathit{\tau }}_{\mathrm{swir}}$ at 2350 nm, size parameter $\mathit{\alpha }$ and height parameter $z_{\mathrm{aer}}$. Following this work, we filter out retrievals with $${\displaystyle \frac{{\mathit{\tau }}_{\mathrm{swir}}\cdot z_{\mathrm{aer}}}{\mathit{\alpha }}}>120\text{m}.$$ This indicates that the retrieval algorithm has difficulties with scenes that have optically thick scattering layers of large particles at high altitude. Additionally, we filter out cases with retrieved albedo in the SWIR band smaller than 0.02. More a posteriori filters (e.g. based on goodness of fit) will be determined after launch using real observations.

Sensitivity studies

Here, we evaluate the sensitivity of the retrieved XCH${}_{\mathrm{4}}$ forward model error to instrument errors and auxiliary input to the retrieval algorithm (e.g. meteorological data). To put the errors in perspective, the S5P product accuracy of XCH${}_{\mathrm{4}}$ was envisioned to be within 2 % and the product precision to be within 0.6 % . Accuracy is defined as the mean deviation from the truth, and precision is the variation due to random processes such as instrument noise. More recently, the requirements have been slightly reformulated as 1 % bias and 1 % precision . From the 1 % bias, 0.8 % is reserved for forward model errors and 0.6 % for instrument-related errors.

Synthetic measurements: the global ensemble

We performed detailed sensitivity studies of the CH${}_{\mathrm{4}}$ algorithm on a global ensemble of simulated spectra consisting of land-only, clear-sky scenes on a ${2.79}^{\circ }\times {2.8125}^{\circ }$ latitude $\times$ longitude grid. This ensemble is, to a large extent, identical to the one used by . It contains realistic aerosol- and cirrus-loaded scenes for 4 days, one per season. The treatment of aerosols and cirrus in the simulations is far more complex than in the retrieval forward model, where only one effective aerosol type is considered; see Sect. . For the simulations, the aerosol physical properties and vertical distributions are derived from the global aerosol model ECHAM5-HAM for five different chemical species and on a superposition of seven log-normal size distributions. The aerosol optical thickness is derived from MODIS observations . Furthermore, the simulations contains cirrus with optical thickness and vertical distribution based on CALIOP measurements . Finally, the surface albedo in the NIR is taken from the MODIS land albedo product in the 841–876 nm channel. For the albedo in the SWIR, the SCIAMACHY surface albedo product at 2350 nm is used . For simplicity, we assume a constant surface albedo within the NIR and SWIR band in our simulations. We expect no difficulties to fit a more realistic spectral dependent surface albedo based on our experience with real GOSAT data. We refer to Figs. 2 and 3 of for the geographical distribution of the total optical thickness and surface albedo, respectively, used in the simulations. The measurements are simulated for the nadir viewing direction and a solar zenith angle that is representative of TROPOMI with an overpass time of 13:30 local time. We have 8633 simulated measurements in total.

While only investigated the scattering-induced error, we increased the inconsistency between the simulation forward model and the retrieval forward model model. We have attempted to include the most important contributions to the forward model error, except for errors due to the underlying spectroscopic database, which have been investigated elsewhere; see and . The simulations were computed using a line-by-line radiative transfer model, whereas the retrieval method uses the linear k-method. In addition, the simulations have a higher vertical and spectral resolution than used in the retrieval. Furthermore, we have added chlorophyll fluorescence emission. In the simulations, fluorescence is modelled to have a double Gaussian spectral shape , which is different from the linear spectral shape assumed in the retrieval forward model (see Eq. ). As in the retrieval scheme, we neglect scattering of fluorescence emission for simplicity. The fluorescence at the TOA then becomes $$F_{\text{s}}\left(\mathit{\lambda }\right)=F_{\text{s},755}\left(\sum _{i=\mathrm{1},\mathrm{2}},A_i,e^{\frac{-(\mathit{\lambda }-{\mathit{\lambda }}_i{)}^{\mathrm{2}}}{{\mathit{\sigma }}_i^{\mathrm{2}}}}\right)e^{-{\mathit{\tau }}_{{\mathrm{O}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)/\mathit{\mu }}.$$ For the parameters $A_{\mathrm{1}}$, $A_{\mathrm{2}}$, ${\mathit{\lambda }}_{\mathrm{1}}$, ${\mathit{\lambda }}_{\mathrm{2}}$, ${\mathit{\sigma }}_{\mathrm{1}}$ and ${\mathit{\sigma }}_{\mathrm{2}}$, we use the same values as in . We only included fluorescence emission for scenes in the global ensemble with albedo${}_{\mathrm{NIR}}/$ albedo${}_{\mathrm{SWIR}}>\mathrm{5}$ as a rough selection criterion for regions with vegetation.

After convolving the simulated TOA radiance with the ISRF, the spectra are superimposed with instrument noise from the TROPOMI noise model . For the NIR, the noise consists solely of shot noise, while for the SWIR, the noise is composed of both shot noise and a signal-independent term. The corresponding continuum signal-to-noise ratios are 500 in the NIR and 100 in the SWIR for a reference scene with surface albedo $A_{\text{s}}=0.05$, viewing zenith angle $\mathrm{VZA}={\mathrm{0}}^{\circ }$ and solar zenith angle $\mathrm{SZA}={70}^{\circ }$.

The baseline performance of the operational CH${}_{\mathrm{4}}$ algorithm is tested on the simulated global ensemble described above. In Fig. , we show a world map of the bias $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},\mathrm{F}}$ before and after applying the a posteriori filters based on retrieved scattering parameters and albedo. Note that the error due to measurement noise $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},y}$ has been subtracted from the total XCH${}_{\mathrm{4}}$ error. This is a random error and is evaluated separately in Sect.  In Fig. , the cumulative probability distribution of the absolute XCH${}_{\mathrm{4}}$ retrieval error is shown (blue line). We get a convergence rate of 99 %, and 53 % of the converged retrievals pass the filters. Finally, 94 % of the valid retrievals have an absolute error $<$ 1 %. In Fig. , we also plotted the cumulative probability distribution in case fluorescence emission is not fitted (red line). Then, we have a convergence rate of 95 and 92 % of the valid retrievals have an error $<$ 1 %. Thus our retrieval results are improved by fitting fluorescence, which is why we have included this in the baseline.

We note that most of the rejected retrievals are found in the desert and dust regions in northern Africa and central Asia in April and July where our simulations capture dust storm events. For real TROPOMI measurements though, we will less likely filter out large regions because TROPOMI has higher spatial resolution than our simulations. Furthermore, our simulations only represent 1 day per season with a mean atmospheric state. Therefore our results merely indicate that certain regions will have fewer days with successful retrievals in certain seasons, not that entire regions will be rejected.

Residual forward model XCH${}_{\mathrm{4}}$ errors for the baseline retrieval method for the clear-sky global TROPOMI measurement ensemble before (a) and after (b) a posteriori filtering. The ensemble covers scenes for a day in January (JAN), April (APR), July (JUL) and October (OCT) as described in Sect. .

[Figure omitted. See PDF]

Cumulative probability distribution of the absolute XCH${}_{\mathrm{4}}$ forward model error for the baseline retrieval (blue line) and for the case that fluorescence emission is not fitted (red line).

[Figure omitted. See PDF]

Sensitivity to atmospheric input data

We investigated the effect of imprecise atmospheric input from TM5 and ECMWF on the CH${}_{\mathrm{4}}$ retrievals. The results are summarised in Fig. .

A priori CH${}_{\mathrm{4}}$ profile

For the baseline performance test in Sect.  we used the true profile of methane ${\mathit{x}}_{\mathrm{true}}$ as the a priori profile for the retrievals. The bias $\mathrm{\Delta }{\mathrm{XCH}}_{\mathrm{4},\mathrm{F}}$ as defined in Eq. () should not depend on the choice of the a priori profile because this effect is accounted for by the averaging kernel. To test this, we take the zonal mean as the a priori profile, where we averaged over all longitudes within a ${2.79}^{\circ }$ latitude bin, with a column deviation up to $\pm$2 %. To illustrate the effect on the global ensemble, we evaluate the root mean square (rms) of the XCH${}_{\mathrm{4}}$ bias (i.e. the total XCH${}_{\mathrm{4}}$ error minus the contribution due to noise) of all retrievals that pass the a posteriori filters. Figure  (first panel, blue line) shows that the a priori CH${}_{\mathrm{4}}$ does not influence the retrieval accuracy, in terms of the rms of the XCH${}_{\mathrm{4}}$ bias, nor the stability, in terms of the convergence rate or number of valid retrievals.

A priori H${}_{\mathrm{2}}$O profile

To investigate the sensitivity to errors on the assumed H${}_{\mathrm{2}}$O profile, we follow the same procedure as in Sect. . The error on the prior H${}_{\mathrm{2}}$O profile is established in the same way as for CH${}_{\mathrm{4}}$, i.e. by taking a normalised zonal mean profile per latitude bin. Note that for H${}_{\mathrm{2}}$O, there is an additional (minor) influence on the retrieval of XCH${}_{\mathrm{4}}$ through the dry air column. The H${}_{\mathrm{2}}$O column error is varied up to $\pm$10 %. Figure  (left panel, red line) shows that the prior H${}_{\mathrm{2}}$O profile has negligible influence on the rms of the XCH${}_{\mathrm{4}}$ bias, increasing it with $<$ 0.01 %. There is a small effect on the convergence rate, reducing it from 99 to 97 % and leading to a reduction of valid retrievals from 53 to 50 %. We note that taking a zonal mean profile for H${}_{\mathrm{2}}$O represents a worst case in terms of accuracy for the specific humidity of ECMWF.

Pressure

An erroneous pressure affects the retrieval of XCH${}_{\mathrm{4}}$ in two ways: first of all, through the pressure dependence of the cross sections and, secondly, through the retrieved air column that is used to convert the CH${}_{\mathrm{4}}$ total column to the dry air mixing ratio, XCH${}_{\mathrm{4}}$. The latter will introduce a retrieval error of the same magnitude as the pressure error. To evaluate the net effect of a pressure error, the prior pressure profile is perturbed with a scaling factor up to $\pm$0.3 %, corresponding to $\left|\mathrm{\Delta }{P}_{\mathrm{surf}}\right|\approx$ 3 hPa. We expect a better accuracy from the ECMWF surface pressure together with the digital elevation map . Figure  (right panel, blue line) shows that the increase in the rms of the XCH${}_{\mathrm{4}}$ bias is $<$ 0.15 %. There is no visible effect on the stability of the algorithm.

Temperature

An error in the temperature will propagate to the XCH${}_{\mathrm{4}}$ retrievals though the temperature dependence of the cross sections. To investigate this effect, the temperature profile is offset up to $\pm$2 K. Figure  (right panel, blue line) shows that the increase in the rms of the XCH$\mathrm{4}$ bias is $<$ 0.15 %. There is a small effect on the stability of the algorithm, reducing the convergence rate to 97 % and the number of valid retrievals to 47 %.

Influence of errors in atmospheric input on the accuracy (panels 1 and 2) and stability (panels 3 and 4) of XCH${}_{\mathrm{4}}$ retrievals. The upper and lower $x$ axes refer to perturbations. Accuracy refers to the XCH${}_{\mathrm{4}}$ bias, and stability refers to the fraction of converged and valid retrievals. The profiles of methane and water have been perturbed in panels 1 and 3; the pressure and temperature profiles have been perturbed in panels 2 and 4.

[Figure omitted. See PDF]

Sensitivity to instrument errors

We investigated the effect of different possible instrument and calibration errors on the CH${}_{\mathrm{4}}$ retrievals. The results are summarised in Table . Below we discuss each effect separately.

Signal-to-noise ratio

The simulated spectra include instrument noise as described in Sect. . The precision is given by the standard deviation of the retrieval noise ${\mathit{\sigma }}_{{\mathrm{XCH}}_{\mathrm{4}}}$; see Eq. (). The world map in Fig.  shows the precision relative to the retrieved XCH${}_{\mathrm{4}}$. Typically the precision is better than the accuracy. The signal-to-noise ratio only becomes a limiting factor for scenarios with snow-covered ground and large SZA, which is why we filter for SZA $<$ 70${}^{\circ }$ and albedo $>$ 0.02 to keep this error relatively small.

Relative precision of XCH${}_{\mathrm{4}}$ due to the instrument noise for the a posteriori filtered data set for the clear-sky ensemble in Fig. .

[Figure omitted. See PDF]

Instrument spectral response function

The synthetic measurements were created by convolving the underlying line-by-line spectra with a Gaussian ISRF with a full width half maximum (FWHM) of 0.38 nm in the NIR band and 0.25 nm in the SWIR band. The mission requirements state that the ISRF shall be known within 1 % of its maximum . This is achieved approximately by varying the FWHM by 1 %. Table  give the results for retrievals with an assumed error in the FWHM. We find that the retrievals are mostly sensitive to the accuracy of the ISRF in the SWIR band, leading to an increase of the rms XCH${}_{\mathrm{4}}$ bias from 0.56 % (baseline) to 0.67 %. We note that of all the instrument errors investigated here, the ISRF gives the largest error contribution to the methane retrievals. Thus, our study indicates that accurate calibration of the ISRF should have high priority.

Spectral calibration

According to the instrument requirements, the centre wavelengths of spectral channels are known within 2 pm (picometres). For our global ensemble, a spectral shift of 2 pm has negligible effect on the error characteristics. Here, we evaluate the effect of a wavelength shift of $\mathrm{1}/20$ of the spectral sampling distance, i.e. 10 and 5 pm for the NIR and SWIR band, respectively. The reference retrieval fits a spectral shift. To test this fitting option, the synthetic spectra were shifted with a constant wavelength: $${\mathit{\lambda }}_{k,\mathrm{meas}}={\mathit{\lambda }}_k+\mathrm{\Delta }\mathit{\lambda },$$ where ${\mathit{\lambda }}_k$ is the real wavelength, ${\mathit{\lambda }}_{k,\mathrm{meas}}$ is the measured wavelength at pixel $k$ and $\mathrm{\Delta }\mathit{\lambda }$ the spectral shift. In Table , results are given for the case that a spectral shift is not fitted and is fitted (between brackets). When fitted, the performance is as good as for the reference retrieval, i.e. simulations without an error in the spectral position. This is as expected and indicates that the spectral shift fitting is robust.

Optionally a wavelength-dependent shift, $\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{squeeze}}$, can be fitted. To test this fitting option, the synthetic wavelength grid was “squeezed”: $${\mathit{\lambda }}_{k,\mathrm{meas}}={\mathit{\lambda }}_k+\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{squeeze}}\times ({\mathit{\lambda }}_k-{\mathit{\lambda }}_{\mathrm{mid}})/({\mathit{\lambda }}_{\mathrm{end}}-{\mathit{\lambda }}_{\mathrm{mid}}),$$ where ${\mathit{\lambda }}_{\mathrm{end}}$ and ${\mathit{\lambda }}_{\mathrm{mid}}$ are the wavelengths at the end and middle of the band, respectively. Table  shows the performance of retrievals with this assumed error on the measured wavelength grid. Since a spectral squeeze is a second-order effect and has a much smaller impact on the XCH${}_{\mathrm{4}}$ retrievals than a spectral shift, it is not fitted in the baseline. However, our results show that, if needed, the option to fit a spectral squeeze can be used reliably.

Effect of instrument calibration errors on convergence rate, fraction of valid retrievals after filtering and rms values of XCH${}_{\mathrm{4}}$ bias and precision for the global ensemble. Note that all sensitivities include the baseline forward model error, caused mainly by aerosol and cirrus scattering. The terms between brackets are for the cases where the relevant quantity is also retrieved. For each instrument calibration error, multiple simulation runs were performed with all combinations of errors in NIR and SWIR channels. The results shown here correspond to the runs with poorest performance in terms of the rms error.

	Convergence	Valid retrievals	rms of XCH${}_{\mathrm{4}}$ bias	rms of XCH${}_{\mathrm{4}}$ precision
Baseline	99 %	53 %	0.56 %	0.43 %
$\mathrm{\Delta }{\mathrm{FWHM}}_{\mathrm{NIR}}=-\mathrm{1}$ %	99 %	53 %	$0.56$ %	0.43 %
$\mathrm{\Delta }{\mathrm{FWHM}}_{\mathrm{SWIR}}=-\mathrm{1}$ %	99 %	49 %	0.67 %	0.44 %
$\mathrm{\Delta }{\mathrm{FWHM}}_{\mathrm{NIR}/\mathrm{SWIR}}=-\mathrm{1}$ %	99 %	49 %	0.67 %	0.44 %
$\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{shift},\mathrm{NIR}}=-10$ pm	97 %	48 %	0.57 %	0.42 %
$\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{shift},\mathrm{SWIR}}=\mathrm{5}$ pm	99 %	50 %	0.99 %	0.42 %
$\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{shift},\mathrm{NIR}/\mathrm{SWIR}}=-10/\mathrm{5}$ pm	99 % (99 %)	46 % (53 %)	1.02 % (0.56 %)	0.41 % (0.43 %)
$\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{squeeze},\mathrm{NIR}}=10$ pm	99 %	53 %	0.56 %	0.42 %
$\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{squeeze},\mathrm{SWIR}}=\mathrm{5}$ pm	99 %	52 %	0.62 %	0.42 %
$\mathrm{\Delta }{\mathit{\lambda }}_{\mathrm{squeeze},\mathrm{NIR}/\mathrm{SWIR}}=10/\mathrm{5}$ pm	96 % (99 %)	53 % (53 %)	0.63 % (0.56 %)	0.42 % (0.43 %)
$I_{\mathrm{offset},\mathrm{NIR}}=-0.1$ %	99 %	53 %	0.57 %	0.43 %
$I_{\mathrm{offset},\mathrm{SWIR}}=0.1$ %	99 %	50 %	0.58 %	0.43 %
$I_{\mathrm{offset},\mathrm{NIR}/\mathrm{SWIR}}=-0.1/0.1$ %	99 %	50 %	0.59 %	0.43 %
$G_{\mathrm{NIR}}=1.02$	99 %	53 %	0.56 %	0.43 %
$G_{\mathrm{SWIR}}=1.02$	99 %	53 %	0.58 %	0.42 %
$G_{\mathrm{NIR}}=G_{\mathrm{SWIR}}=1.02$	99 %	52 %	0.58 %	0.42 %

Methane bias due to heterogenous slit illumination for spatially varying surface reflection over a marsh scene at Siberia close to the river Ob at latitude 62.8${}^{\circ }$ N and longitude 72.1${}^{\circ }$ E. Measurement simulations are performed with the instrument model by for an instantaneous field of view of 3.4 km across the slit and 7.0 km along the slit.

[Figure omitted. See PDF]

Radiometric offset: additive factor

The effect of an unknown systematic offset in the Earth radiance is investigated. The offsets in the NIR and SWIR bands are independently varied with $\pm$0.1 % of the continuum. Table  shows the effect of a radiometric offset to the XCH${}_{\mathrm{4}}$ retrievals. We note that a radiometric offset in the SWIR band causes a larger XCH${}_{\mathrm{4}}$ error than an offset in the NIR band. The latter is partly compensated by the retrieved fluorescence.

Radiometric gain: multiplicative factor

The absolute radiometric accuracy of the measurement of the Earth spectral radiance shall be better than 2 % according to the system requirements. To investigate the effect of such an error, the synthetic spectra were multiplied by a scaling factor $G$. Table  shows that there is negligible effect of an error of 2 % in radiometric gain. This error is largely compensated by the retrieved surface albedo. It follows that interaction between surface albedo and aerosols has a negligible impact for gain errors $<$ 2 %.

Heterogenous slit illumination

For the two-dimensional TROPOMI push-broom spectrometer, light in across-slit dimension is dispersed by the instrument grating in order to spectrally resolve the received signal. The along-slit direction is aligned across flight direction to achieve the desired spatial resolution. For a homogeneous illumination of the instrument slit, the spectral instrument response is characterised extensively during the pre-flight calibration of the TROPOMI instrument, and it is used as baseline to simulate the TROPOMI radiometric measurement in our retrieval. In space, however, the instrument slit will be illuminated inhomogeneously due to ground scene heterogeneities on scales smaller than the instrument's field of view. Inhomogeneous illumination across the slit leads to a distortion of the ISRF as described by , and , and it can affect the retrieval accuracy of the TROPOMI methane product. Small-scale heterogeneities of the ground scene are generally caused by spatial variations of surface reflection and by broken clouds. Because of our strict cloud filtering, spatial variations in surface reflection are the only cause of methane retrieval biases due to inhomogeneous slit illumination. To evaluate the effect of surface scene heterogeneity on our XCH${}_{\mathrm{4}}$ product, we employ the instrument model described by for both the NIR and SWIR bands of our retrieval. Furthermore, we use a high spatial resolution MODIS albedo map for the 50 $\times$ 50 km${}^{\mathrm{2}}$ marsh region in central Siberia with structures in the surface reflection due to ponds, shown in Fig. . Depending on the scene heterogeneity in the flight direction, the XCH${}_{\mathrm{4}}$ error shows an oscillation structure with a maximum amplitude $\le$ 0.4 %, a standard deviation of 0.12 % and a mean error of $-0.01$ %. For a particular temporal and spatial sampling of the scene, a pseudo-random scatter is introduced to the XCH${}_{\mathrm{4}}$ product. This means that overall the effect can be considered small. One may consider this error as a limitation when interpreting very localised sources in surroundings of heterogenous surface reflection, but for most applications some averaging either in time or space will be done, which reduces this error.

Cloud filtering

The global ensemble from as used in Sect.  cannot be used to evaluate the cloud filters, because it consists purely of cloud-free scenes. Therefore, the cloud filters are tested using a new ensemble of synthetic measurements, which also include scenes with water clouds.

Synthetic measurements: the TROPOMI test orbit

Simulations of a TROPOMI orbit. Panel (a) shows MODIS cloud fraction resampled on the orbit's ground pixels. Panel (b) gives the XCH${}_{\mathrm{4}}$ bias of all processed pixels that converged. Panel (c) and (d) give XCH${}_{\mathrm{4}}$ bias of the valid retrievals after cloud filtering with MODIS data and the backup cloud filter, respectively. Here, we also applied the a posteriori filter based on retrieved scattering parameters and albedo.

[Figure omitted. See PDF]

Error statistics of XCH${}_{\mathrm{4}}$ retrievals from the L1B orbit using different cloud filters in comparison to the global clear-sky orbit. The performance of the MODIS filter is expected to be comparable with the operational cloud filter using VIIRS data.

	Retrievals with
	$\left|{\mathrm{XCH}}_{\mathrm{4}},,\text{error}\right|$ $<$ 1 %	$\left|{\mathrm{XCH}}_{\mathrm{4}},,\text{error}\right|$ $<$ 0.5 %	rms of XCH${}_{\mathrm{4}}$ bias	rms of XCH${}_{\mathrm{4}}$ precision
Global clear-sky ensemble	94 %	78 %	0.56 %	0.43 %
L1B orbit with MODIS filter	96 %	80 %	0.56 %	0.31 %
L1B orbit with backup cloud filter	94 %	79 %	0.71 %	0.27 %

To test the performance of the proposed backup cloud filter, we have constructed synthetic TROPOMI L1B radiance spectra for an entire orbit that passes over Africa. We used realistic viewing geometries from a TROPOMI orbit simulator provided by the Royal Netherlands Meteorological Institute (KNMI). The meteorological data are from ECMWF. We used CO and CH${}_{\mathrm{4}}$ model profiles from TM5 . The surface elevation comes from a digital elevation map constructed by KNMI based on USGS and NASA data . The aerosol and cirrus properties and surface albedos are taken from the global ensemble, which means that these are the same for all TROPOMI ground pixels within a $\sim$ 3${}^{\circ }$ $\times$ 3${}^{\circ }$ latitude $\times$ longitude box. While the global ensemble used for the sensitivity studies is land-only and clear-sky, the test orbit contains cloudy scenes over land and ocean. Cloud fraction, cloud optical thickness and cloud top pressure are obtained from MODIS Aqua measurements at $\mathrm{5}\times \mathrm{5}$ km${}^{\mathrm{2}}$. Here, the cloud top height is derived from the cloud top pressure and the surface pressure, also provided by MODIS. The cloud properties are collocated in time and space to the TROPOMI orbit and used in the measurement simulation. For fractional clouds, the independent-pixel approximation is used to combine the cloudy and clear-sky parts of the scene. For reference, the cloud fraction used in the simulations is shown in Fig. , upper left. We note that 28 % of the TROPOMI ground pixels are fully cloud-free for this test orbit. Globally, one would expect, on average, 20 % cloud-free pixels . Resampling of the auxiliary data on the TROPOMI ground pixels has been performed using the Multi-Instrument Preprocessor (MIPrep) developed at SRON.

Performance of cloud filters

First, we show the performance of the XCH${}_{\mathrm{4}}$ retrievals on the test orbit when no a posteriori data filtering is applied, only a priori filtering of ocean pixels and SZA $>$ 70${}^{\circ }$ and VZA $>$ 50${}^{\circ }$; see Fig. , upper right. It is clear that cloud-contaminated measurements lead to large XCH${}_{\mathrm{4}}$ errors ($>$ 2 %).

Assuming that we have cloud data from VIIRS, we would then be able to filter out the cloudy pixels almost perfectly. To illustrate the effect, we filtered out pixels with cloud fraction $>$ 0.02 in Fig. , lower left. Note that in this case we have also applied a posteriori filtering based on retrieved scattering parameters and albedo. One is then left with valid retrievals of $\sim$ 3 % of all simulations in the test orbit.

In comparison, the performance of the backup cloud filter based on the difference between the H${}_{\mathrm{2}}$O column retrieved from strong and weak bands (see Sect. ) is shown in Fig. , lower right. The backup cloud filter removes most cloudy pixels, but some remain. In Table  and Fig. , the statistics of XCH${}_{\mathrm{4}}$ retrievals on the orbit are summarised. After cloud filtering with MODIS data (representative of operational VIIRS data), the results for the test orbit are comparable to the clear-sky global ensemble. However, the backup cloud filter is less effective. The rms of the XCH${}_{\mathrm{4}}$ bias is then 0.71 % instead of the 0.56 % that is expected for the operational VIIRS cloud mask.

Cumulative probability distribution of the absolute XCH${}_{\mathrm{4}}$ bias for the simulated level 1b orbit. Cloud-contaminated measurements are filtered using MODIS data (blue line) or the backup cloud filter (red line).

[Figure omitted. See PDF]

Conclusions

This paper describes the algorithm baseline of the operational methane retrievals from the S5P measurements. The level 2 product includes the column-averaged dry air mixing ratio XCH${}_{\mathrm{4}}$, the column averaging kernel and the noise standard deviation. In order to account for the effect of aerosols and cirrus, the algorithm developed retrieves the methane column simultaneously with effective scattering parameters related to particle amount, size and height distribution. The choice of scattering parameters reflects the information content of the measurements as closely as possible. The retrieval algorithm uses the radiance and irradiance measurements in the SWIR 2305–2385 nm band and additionally in the NIR band between 757 and 774 nm (O${}_{\mathrm{2}}$ A-band). The forward model of the retrieval algorithm uses online radiative transfer calculations, fully including multiple scattering in an efficient manner. Absorption cross sections of the relevant atmospheric trace gases and optical properties of aerosols are calculated from lookup tables. The inversion is performed using Phillips–Tikhonov regularisation in combination with a reduced-step-size Gauss–Newton iteration scheme.

To test the developed algorithm we generated two ensembles of simulated measurements that cover the range of scenes that will likely be encountered by the S5P instrument: one clear-sky global ensemble and one test orbit containing cloud-contaminated measurements. Overall, the algorithm developed performs well in correcting for the effect of aerosols and cirrus clouds on the retrieved XCH${}_{\mathrm{4}}$. For both ensembles, $\sim$ 80 % of the cases have an XCH${}_{\mathrm{4}}$ error $<$ 0.5 % and $\sim$ 95 % have an error $<$ 1 %. To achieve this, it is necessary to apply a priori filtering of cloud-contaminated scenes and a posteriori filtering based on retrieved parameters. We illustrated the performance of the proposed backup cloud filter based on retrievals of H${}_{\mathrm{2}}$O from weak and strong absorption bands in the SWIR under the assumption of a non-scattering atmosphere. It should be noted that the cloud filter based on S5P measurements itself is less efficient than the VIIRS cloud mask for water clouds.

Apart from forward model errors induced by aerosols, we also studied effects of model errors in temperature, pressure and water vapour profiles. We expect to stay within product requirements for errors in input profiles of water, pressure and temperature below 10 %, 0.3 % and 2 K, respectively. Another relevant source of errors to the CH${}_{\mathrm{4}}$ data product could be spectroscopic errors. This has been studied in detail by and . Note that a study is ongoing to improve the spectroscopic data for TROPOMI's SWIR spectral range . Concerning instrument errors, we found that the most critical error source is an error in the ISRF in the SWIR band. To conclude, we have shown that for a compliant instrument, our algorithm provides a methane product that meets the requirements.

Data availability

The underlying research data for the simulated atmospheric profiles and TROPOMI measurement simulation and retrievals are available upon request from Haili Hu ([email protected]). HITRAN spectroscopic line parameters are available through HITRAN online (http://hitran.org), and the line parameters from are available in their Supplement.

Sensitivity of retrievals to fixed-value aerosol parameters

Effect of fixed-value aerosol parameters on convergence rate, fraction of valid retrievals after filtering and rms values of XCH${}_{\mathrm{4}}$ bias and precision for the global ensemble. See text for values used in the baseline run.

	Convergence	Valid retrievals	rms of XCH${}_{\mathrm{4}}$ bias	rms of XCH${}_{\mathrm{4}}$ precision
Baseline	99 %	53 %	0.56 %	0.43 %
$n=1.37-0.0005i$	99 %	47 %	0.60 %	0.40 %
$n=1.55-0.02i$	99 %	53 %	0.57 %	0.43 %
$w=\mathrm{1}$ km	99 %	53 %	0.56 %	0.42 %
$w=\mathrm{3}$ km	99 %	53 %	0.56 %	0.41 %

Our retrieval algorithm retrieves three effective aerosol parameters, namely total column, size parameter and central height. Other aerosol parameters are assumed to be fixed, such as the refractive indices $n$ and width $w$ of the aerosol layer (i.e. FWHM of Gaussian height distribution). The values in the baseline retrieval are fixed at $n_{\mathrm{NIR}}=1.4-0.01i$, $n_{\mathrm{SWIR}}=1.47-0.008i$ and $w=\mathrm{2}$ km, respectively. To justify this assumption, we studied the effect of the fixed-value aerosol parameters on the methane retrievals for our synthetic global ensemble as described Sect. . We have performed four retrieval test runs with refractive indices of $1.37-0.0005i$ and $1.55-0.02i$, representative of the majority of aerosol types , and with widths of 1 and 3 km. For each test run the other values are kept to the baseline values. For simplicity, we used the same refractive index for NIR and SWIR in the test runs. The results are summarised in Table  from which we conclude that, within a reasonable range, the exact value of the fixed aerosol parameters does not affect the methane retrievals significantly. We do see a slight deterioration of our results when a relatively low refractive index is used.

Acknowledgements

We thanks M. Sneep (The Royal Netherlands Meteorological Institute, KNMI) and S. Houweling (SRON) for providing us with auxiliary data needed to simulate TROPOMI measurements. This research has been funded in part by the TROPOMI national programme from the Netherlands Space Office (NSO). André Butz is supported by the Deutsche Forschungsgemeinschaft through the Emmy Noether Programme, grant BU2599/1-1 (RemoteC). Edited by: J. Kim Reviewed by: two anonymous referees