An improved air mass factor calculation for nitrogen dioxide measurements from the Global Ozone Monitoring Experiment-2 (GOME-2)

An improved tropospheric nitrogen dioxide (NO2) retrieval algorithm from the Global Ozone Monitoring Experiment-2 (GOME-2) instrument based on air mass factor (AMF) calculations performed with more realistic model parameters is presented. The viewing angle dependency of surface albedo is taken into account by improving the GOME-2 Lambertian-equivalent reflectivity (LER) climatology with a directionally dependent LER (DLER) dataset over land and an ocean surface albedo parameterisation over water. A priori NO2 profiles with higher spatial and temporal resolutions are obtained from the IFS (CB05BASCOE) chemistry transport model based on recent emission inventories. A more realistic cloud treatment is provided by a cloudsas-layers (CAL) approach, which treats the clouds as uniform layers of water droplets, instead of the current cloudsas-reflecting-boundaries (CRB) model, which assumes that the clouds are Lambertian reflectors. On average, improvements in the AMF calculation affect the tropospheric NO2 columns by ±15 % in winter and ±5 % in summer over largely polluted regions. In addition, the impact of aerosols on our tropospheric NO2 retrieval is investigated by comparing the concurrent retrievals based on ground-based aerosol measurements (explicit aerosol correction) and the aerosolinduced cloud parameters (implicit aerosol correction). Compared with the implicit aerosol correction utilising the CRB cloud parameters, the use of the CAL approach reduces the AMF errors by more than 10 %. Finally, to evaluate the improved GOME-2 tropospheric NO2 columns, a validation is performed using ground-based multi-axis differential optical absorption spectroscopy (MAXDOAS) measurements at different BIRA-IASB stations. At the suburban Xianghe station, the improved tropospheric NO2 dataset shows better agreement with coincident ground-based measurements with a correlation coefficient of 0.94.

Abstract. An improved tropospheric nitrogen dioxide (NO 2 ) retrieval algorithm from the Global Ozone Monitoring Experiment-2 (GOME-2) instrument based on air mass factor (AMF) calculations performed with more realistic model parameters is presented. The viewing angle dependency of surface albedo is taken into account by improving the GOME-2 Lambertian-equivalent reflectivity (LER) climatology with a directionally dependent LER (DLER) dataset over land and an ocean surface albedo parameterisation over water. A priori NO 2 profiles with higher spatial and temporal resolutions are obtained from the IFS (CB05BASCOE) chemistry transport model based on recent emission inventories. A more realistic cloud treatment is provided by a cloudsas-layers (CAL) approach, which treats the clouds as uniform layers of water droplets, instead of the current cloudsas-reflecting-boundaries (CRB) model, which assumes that the clouds are Lambertian reflectors. On average, improvements in the AMF calculation affect the tropospheric NO 2 columns by ±15 % in winter and ±5 % in summer over largely polluted regions. In addition, the impact of aerosols on our tropospheric NO 2 retrieval is investigated by comparing the concurrent retrievals based on ground-based aerosol measurements (explicit aerosol correction) and the aerosolinduced cloud parameters (implicit aerosol correction). Compared with the implicit aerosol correction utilising the CRB cloud parameters, the use of the CAL approach reduces the AMF errors by more than 10 %. Finally, to evaluate the improved GOME-2 tropospheric NO 2 columns, a validation is performed using ground-based multi-axis differential optical absorption spectroscopy (MAXDOAS) measurements at different BIRA-IASB stations. At the suburban Xianghe station, the improved tropospheric NO 2 dataset shows better agreement with coincident ground-based measurements with a correlation coefficient of 0.94.

Introduction
Tropospheric nitrogen dioxide (NO 2 ) is an important air pollutant that harms the human respiratory system, even over short exposure periods (Gamble et al., 1987;Kampa and Castanas, 2008), and contributes to the formation of tropospheric ozone, urban haze, and acid rain (Charlson and Ahlquist, 1969;Crutzen, 1970;McCormick, 2013). Besides natural sources of nitrogen dioxide such as soil emissions and lightning, the combustion-related emission sources from anthropogenic activities like fossil fuel consumption, car traffic, and biomass burning produce substantial amounts of nitrogen oxides (NO x =NO 2 + NO).
Satellite measurements from the Global Ozone Monitoring Experiment (GOME) , the SCanning Imaging Absorption SpectroMeter for Atmospheric CHartographY (SCIAMACHY) (Bovensmann et al., 1999), the Ozone Monitoring Instrument (OMI) (Levelt et al., 2006), and the Global Ozone Monitoring Experiment-2 (GOME-2) (Callies et al., 2000;Munro et al., 2016) have produced global NO 2 measurements on long timescales. New-generation instruments like the TROPOspheric Monitoring Instrument (TROPOMI) (Veefkind et al., 2012) aboard the Sentinel-5 Precursor satellite and geostationary missions such as the Sentinel-4  will continue this record and deliver NO 2 datasets with a high spatial resolution and short revisit time.
The GOME-2 instrument, which is the main focus of this study, is included on a series of MetOp satellites as part of the EUMETSAT Polar System (EPS). The first GOME-2 was launched in October 2006 aboard the MetOp-A satellite, and a second GOME-2 was launched in September 2012 aboard MetOp-B. This dataset will be further extended by a third GOME-2 aboard the MetOp-C satellite, which was launched in December 2018. GOME-2 is a nadir-scanning UV-Vis spectrometer that measures the solar irradiance and Earth's backscattered radiance with four main optical channels: this covers the spectral range between 240 and 790 nm with a spectral resolution between 0.26 and 0.51 nm. The default swath width of GOME-2 is 1920 km, enabling a global coverage in ∼ 1.5 d. The along-track dimension of the instantaneous field of view is ∼ 40 km, whereas the acrosstrack dimension depends on the integration time used for each channel. For the 1920 km swath, the maximum ground pixel size is 80 km × 40 km in the forward scan, which remains almost constant over the full swath width. In a tandem operation between MetOp-A and MetOp-B from July 2013 onwards, a decreased swath of 960 km and an increased spatial resolution of 40 km × 40 km are employed by GOME-2/MetOp-A (Munro et al., 2016). GOME-2 provides morning observations of NO 2 at about 09:30 LT (local time), which complements early afternoon measurements e.g. from OMI or TROPOMI. The GOME-2 NO 2 measurements have been widely used in trend studies, satellite dataset intercomparisons, and NO 2 emission estimations (e.g. Mijling et al., 2013;Hilboll et al., 2013Hilboll et al., , 2017Krotkov et al., 2017;Irie et al., 2012;Gu et al., 2014;Miyazaki et al., 2017;Ding et al., 2017).
The NO 2 retrieval algorithm for the GOME-2 instrument contains three steps: (1) the spectral fitting of the slant column (concentration along the effective light path) using the differential optical absorption spectroscopy (DOAS) method (Platt and Stutz, 2008) from the measured GOME-2 (ir)radiances, (2) the separation of stratospheric and tropospheric contributions using a modified reference sector method, (3) and the conversion of the tropospheric slant column to a vertical column using a tropospheric air mass factor (AMF) calculation. The quality of GOME-2 NO 2 mea-surements is strongly related to the calculation of the AMF, which is determined by a radiative transfer model, depending on a set of model parameters, such as viewing geometry, surface albedo, vertical distribution of NO 2 , cloud, and aerosol. The model parameters, generally taken from external databases, contribute substantially to the overall AMF uncertainty, which is estimated to be in the range of 30 %-40 % (Lorente et al., 2017).
The surface is normally assumed to be Lambertian with an isotropic diffuse reflection independent of the viewing and illumination geometry in the NO 2 retrieval (e.g. Boersma et al., 2011;van Geffen et al., 2019;Liu et al., 2019b). However, due to the occurrences of retroreflection and shading effects (mainly over rough surfaces such as vegetation) and specular reflection (mainly over smooth surfaces like water), the Lambertian assumption is not always fulfilled. To account for the geometry-dependent surface scattering characteristics, the surface bidirectional reflectance distribution function (BRDF) (Nicodemus et al., 1992) has been considered in previous studies (e.g. Zhou et al., 2010;Lin et al., 2014Lin et al., , 2015Noguchi et al., 2014;Vasilkov et al., 2017;Lorente et al., 2018;Laughner et al., 2018;Qin et al., 2019), mainly based on measurements from the MODerate resolution Imaging Spectroradiometer (MODIS) over land. However, due to the use of different instruments, biases are possibly introduced in the NO 2 retrieval. In addition, owing to the generally unavailable full surface BRDF under all conditions and the complexity of accounting for the BRDF, most of the current NO 2 and cloud retrievals still rely on Lambertian surface reflection (e.g. Krotkov et al., 2017;Boersma et al., 2018;van Geffen et al., 2019;Loyola et al., 2018;Desmons et al., 2019).
To account for the varying sensitivity of the satellite to NO 2 at different altitudes, a priori vertical profiles of NO 2 are required, and they are generally prescribed using a chemistry transport model. The importance of the a priori NO 2 profiles used in the retrieval has already been recognised and motivated the use of model data with a high spatial resolution and/or high temporal resolution (e.g. Valin et al., 2011;Heckel et al., 2011;Russell et al., 2011;McLinden et al., 2014;Yamaji et al., 2014;Kuhlmann et al., 2015;Lin et al., 2014;Boersma et al., 2016;Laughner et al., 2016). Within the Monitoring Atmospheric Composition and Climate (MACC) European project, a global data assimilation system for atmospheric composition forecasts and analyses has been developed and is running operationally in Copernicus Atmosphere Monitoring Service (CAMS, http: //atmosphere.copernicus.eu, last access: 12 February 2020). The CAMS system relies on a combination of satellite observations with state-of-the-art atmospheric modelling (Flemming et al., 2017); therefore, the European Centre for Medium Range Weather Forecasts (ECMWF) numerical weather prediction Integrated Forecast System (IFS) was extended to include modules to describe the atmospheric composition Inness et al., 2015;Mor-crette et al., 2009;Benedetti et al., 2009;Engelen et al., 2009;Agustí-Panareda et al., 2016). Profile forecasts from CAMS are planned to be applied in the operational NO 2 retrieval algorithm for the Sentinel-4 (Sanders et al., 2018) and Sentinel-5 (van Geffen et al., 2018) missions and will provide the advantages of operational implementation and high resolution. Lately, an advanced IFS system, referred to as IFS (CB05BASCOE) (Huijnen et al., 2016(Huijnen et al., , 2019 or IFS (CBA) for short, operates at improved horizontal, vertical, and temporal resolutions based on recent emission inventories, providing an improved profile "representativeness".
Clouds influence the NO 2 retrieval due to their increased reflectivity, their shielding effect on the NO 2 column below the cloud, and multiple scattering that enhances absorption inside the cloud (Martin et al., 2002;Liu et al., 2004;Stammes et al., 2008;Kokhanovsky and Rozanov, 2008). The presence of clouds is taken into account in the NO 2 AMF calculation using cloud parameters based on the Optical Cloud Recognition Algorithm (OCRA) and the Retrieval Of Cloud Information using Neural Networks (ROCINN) (Loyola et al., 2007. OCRA/ROCINN has been applied in the operational retrieval of trace gases from GOME (Van Roozendael et al., 2006), GOME-2 Hao et al., 2014;Liu et al., 2019b), and TROPOMI Theys et al., 2017;Loyola, 2020). The latest version of OCRA/ROCINN Loyola et al., 2018) provides two sets of cloud products: one treats clouds as ideal Lambertian reflectors in a "clouds-asreflecting-boundaries" (CRB) model, and the second treats clouds as uniform layers of water droplets in a "clouds-aslayers" (CAL) model. The CAL model, which allows for the penetration of photons through the cloud, is more realistic than the CRB model, which screens the atmosphere below the cloud (Rozanov and Kokhanovsky, 2004;Richter et al., 2015).
Aerosol scattering and absorption influence the top-ofatmosphere radiances and the light path distribution. The radiative effect of scattering aerosols and clouds is comparable (i.e. the albedo effect, the shielding effect, and multiple scattering), whereas the presence of absorbing aerosols generally reduces the sensitivity to NO 2 within and below the aerosol layer by decreasing the number of photons scattered back from this region to the satellite (Leitão et al., 2010). Because cloud retrieval does not distinguish between clouds and aerosols, the effect of aerosol on the AMF is normally corrected using an "implicit aerosol correction" by assuming that the effective clouds retrieved as Lambertian reflectors (i.e. using the CRB model); this accounts for the effect of aerosols on the light path (Boersma et al., 2004. Previous works have also applied an "explicit aerosol correction" for OMI pixels by considering additional aerosol parameters (e.g. Lin et al., 2014Lin et al., , 2015Kuhlmann et al., 2015;Castellanos et al., 2015;Liu et al., 2019a;Chimot et al., 2019) and have reported large biases related to the implicit aerosol correction for polluted cases, which is likely due to the fact that the simple CRB model can not fully describe the effects inherent to aerosol particles (Chimot et al., 2019).
The operational GOME-2 NO 2 products are generated using the GOME Data Processor (GDP) algorithm and are provided by DLR in the framework of EUMETSAT's Satellite Application Facility on Atmospheric Composition Monitoring (AC-SAF). The retrieval algorithm of total and tropospheric NO 2 from GOME-2 has been introduced by Valks et al. (2011Valks et al. ( , 2017 as implemented in the current operational GDP version 4.8. An updated slant column retrieval and stratosphere-troposphere separation have been presented by Liu et al. (2019b), and an improved AMF calculation is described in this paper, which will be implemented in the next version of GDP.
In the context of AC-SAF (Hassinen et al., 2016), the NO 2 data derived from the GOME-2 GDP algorithm are being validated at BIRA-IASB by comparison with correlative observations from ground-based multi-axis differential optical absorption spectroscopy (MAXDOAS) (Pinardi et al., , 2015Pinardi, 2020). The MAXDOAS instrument collects scattered sky light in a series of line-of-sight angular directions extending from the horizon to the zenith. High sensitivity towards absorbers near the surface is obtained for the smallest elevation angles, whereas measurements at higher elevations provide information on the rest of the column. This technique allows for the determination of vertically resolved abundances of atmospheric trace species in the lowermost troposphere (Hönninger et al., 2004;Wagner et al., 2004;Wittrock et al., 2004;Heckel et al., 2005).
This work is organised as follows. In Sect. 2, we briefly introduce the reference retrieval algorithm for GOME-2 NO 2 measurements, which was described in detail in Liu et al. (2019b). In Sect. 3, we improve the AMF calculation in the reference retrieval algorithm by accounting for the dependency of surface albedo on direction over land and over water, applying the advanced IFS (CBA) a priori NO 2 profiles with higher model resolution, and implementing the more realistic CAL cloud model. We investigate the properties of the implicit aerosol correction for aerosol-dominated scenes by comparing it to the explicit aerosol correction in Sect. 4. Finally, we show a validation of the GOME-2 tropospheric NO 2 columns using MAXDOAS datasets in Sect. 5.
2 Reference retrieval for GOME-2 NO 2 measurements As described in Liu et al. (2019b), the NO 2 slant column retrieval applies an extended 425-497 nm wavelength fitting window  to include more NO 2 structures and an improved slit function treatment to compensate for the long-term and in-orbit drifts of the GOME-2 slit function. The uncertainty in the NO 2 slant columns is ∼ 4.4 × 10 14 molec. cm −2 and is calculated from the average slant column error using a statistical method (Valks et al., 2011, see Sect. 6.1 therein). To determine the strato-spheric NO 2 components, the STRatospheric Estimation Algorithm from Mainz (STREAM) method  with an improved treatment of polluted and cloudy pixels is adopted. The uncertainty in the GOME-2 stratospheric columns is ∼ 4-5 × 10 14 molec cm −2 for polluted conditions based on the daily synthetic GOME-2 data and is ∼ 1-2 × 10 14 molec cm −2 for monthly averages.
Mainly focusing on the third retrieval step, we apply the tropospheric AMF M conversion (Palmer et al., 2001;Boersma et al., 2004) to account for the average light path through the atmosphere: where m l represents the box-air-mass factors (box-AMFs) in layer l, x l represents the partial columns from the a priori NO 2 profiles, and c l is a correction coefficient to account for the temperature dependency of the NO 2 cross section (Boersma et al., 2004;Nüß et al., 2006;Bucsela et al., 2013).
The box-AMF m l values are derived using the multilayered multiple scattering LIDORT  radiative transfer model and are stored in a look-up table (LUT) as a function of several model inputs (b), including GOME-2 viewing geometry, surface pressure, and surface albedo. Table 1 summarises the ancillary parameters used in the AMF calculation. The surface albedo is described by a monthly Lambertianequivalent reflectivity (LER) database (Tilstra et al., 2017), derived from GOME-2 measurements for the years 2007-2013 with a spatial resolution of 1.0 • × 1.0 • (latitude, longitude) for standard grid cells and 0.25 • × 0.25 • for coastlines . The LER is retrieved by matching the simulated reflectances to the Earth reflectance measurements for cloud-free scenes found using a statistical method (Koelemeijer et al., 2003;Kleipool et al., 2008;Tilstra et al., 2017).
The daily a priori NO 2 profiles are obtained from the TM5-MP three-dimensional chemistry transport model (Williams et al., 2017) with a horizontal resolution of 1 • ×1 • for 34 vertical layers, as summarised in Table 2. The model is driven by ECMWF ERA-Interim meteorological reanalysis (Dee et al., 2011) and updated every 3 h with the interpolation of fields for the intermediate time periods. Compared with previous versions of the TM model (e.g. Williams et al., 2009;Huijnen et al., 2010), which have been commonly used in tropospheric NO 2 retrieval studies (e.g. Boersma et al., 2011;Chimot et al., 2016;Lorente et al., 2017), the main advantages of TM5-MP are the better spatial resolution (1 • ×1 • ), the updated NO x emissions (year-specific MACCity emission inventory, Granier et al., 2011), and the improved chemistry scheme (an expanded version of the modified CB05 chemistry scheme, Williams et al., 2013). The a priori profiles are determined for the GOME-2 overpass time (09:30 LT) and interpolated to the centre of the GOME-2 pixel based on four nearest-neighbour TM5-MP cell centres.
In the presence of clouds, the AMF is derived based on the independent pixel approximation (Cahalan et al., 1994), which assumes that the AMF is a linear combination of a cloudy-sky AMF M cl and a clear-sky AMF M cr : where ω is the cloud radiance fraction. M cl is determined using Eq. (1) with the cloud surface regarded as a Lambertian reflector and with m l = 0 for layers below the cloud top pressure c p . ω is derived from the GOME-2 cloud fraction c f as follows: where I cr is the backscattered radiance for a clear scene derived using LIDORT, and I cl is for a cloudy scene. Note that the cloud fraction c f is a radiometric cloud fraction instead of a geometric one. The GOME-2 cloud properties are derived by the OCRA and the ROCINN algorithms (Loyola et al., 2007Loyola et al., 2018;Lutz et al., 2016). As clouds generally have a higher reflectivity than the ground, OCRA calculates the radiometric cloud fractions by comparing the measured reflectances in three broadband wavelength regions across the UV-Vis-NIR region with corresponding cloudfree background composite maps using a RGB colour space approach. The monthly cloud-free background map is calculated from GOME-2/MetOp-A measurements for the years 2008-2013, accounting for instrumental degradation and dependencies on the viewing zenith angle (VZA), latitude, and season. With the radiometric cloud fractions from OCRA as inputs, ROCINN retrieves the cloud top pressures (cloud top heights) and cloud albedo (cloud optical depth) by comparing the simulated and measured satellite radiances in the O 2 A-band around 760 nm using regularisation theory. Based on the independent pixel approximation and the CRB cloud model, the ROCINN algorithm treats the clouds as Lambertian surfaces.
The tropospheric NO 2 column calculation is complicated in the case of cloudy conditions. For many measurements over cloudy scenes, the cloud top is well above the NO 2 pollution in the boundary layer, and the enhanced tropospheric NO 2 concentrations cannot be detected by GOME-2 if the clouds are optically thick. Therefore, the tropospheric NO 2 column is only calculated for GOME-2 observations with a cloud radiance fraction ω less than 0.5. Note that the "belowcloud amount" (i.e. the amount of NO 2 below the cloud top) for these partly cloudy conditions is implicitly accounted for via the cloudy-sky AMF M cl (in which m l = 0 for layers below the cloud top). Table 1. Ancillary parameters used to derive GOME-2 tropospheric NO 2 columns. See Table 2 for details regarding the chemistry transport models used to obtain the a priori NO 2 profiles.

Surface albedo
The dependency of surface reflection on the incoming and outgoing directions is mathematically described by the BRDF (Nicodemus et al., 1992), which shows a "hot spot" of increased reflectivity in the backward scattering directions over rough surfaces like vegetation and a strong forward scattering peak near "sun glint" geometries over smooth surfaces such as water. In this study, we account for the direction dependency of surface albedo for the GOME-2 LER climatology by applying a directionally dependent LER (DLER) dataset over land surfaces (see Sect. 3.1.1) and by implementing an ocean surface albedo parameterisation over water surfaces (see Sect. 3.1.2).

Over land
To account for the surface BRDF in our NO 2 AMF calculation over land, the surface reflectivity is described by a GOME-2 DLER dataset ) that captures the VZA dependency. Compared with the traditional GOME-2 LER climatology (Tilstra et al., 2017), which is derived from a range of viewing angles (∼ 115 • for GOME-2 measurements covering the directions from east to west), the GOME-2 DLER dataset is derived by dividing the range of viewing angles into five segments and applying the same retrieval method as in the traditional GOME-2 LER determination for each segment with a parabolic fit to parameterise the VZA dependency. The main idea of this VZA-dependency parameterisation is to use the VZA as a proxy for observation geometry over land, as the solar zenith angle (SZA) and relative azimuth angle (RAA) are nearly constant at a given latitude and, thus, have been captured in the original GOME-2 LER dataset.
For each GOME-2 measurement, the surface DLER α DLER is calculated as follows: where VZA θ is positive on the western side of the orbit swath and negative on the eastern side of the orbit swath. c 0 , c 1 , and c 2 are parabolic fitting coefficients depending on latitude, longitude, month, and wavelength. The nondirectional LER α LER is taken from the traditional GOME-2 LER climatology. Note that no directionality is provided by the DLER dataset over water (without sea ice cover), which is mainly due to the dependency on parameters such as wind speed and chlorophyll concentration that can not be easily cast into climatology. Additionally, due to the strong solar and viewing angle dependency of specular reflection, changes in the solar position in the course of a month influence the albedo over water bodies much more than over land, and this influence is modelled and described in Sect. 3.1.2. Figure 1a-c show the traditional GOME-2 LER climatology, the GOME-2 DLER dataset over land, and their differences on 3 February and 5 August 2010. The DLER data show a stronger increase -by ∼ 0.02 -for the western view- ing direction over vegetation, ∼ 0.05 over desert, and ∼ 0.2 over snow and ice, due to the increasing BRDF in the backward scattering direction. A slight change of up to 0.01 is found over vegetation and desert with an enhancement in the central part of the orbit swath and a reduction on the eastern side of the orbit swath. This effect is larger over snow and ice, due to the forward scattering peak or double scattering peak in the BRDF pattern for snow (Dumont et al., 2010). The difference in surface albedo is generally larger in winter due to the change in surface conditions and/or sun elevation, with the exception of desert areas.   Figure 3 shows the differences in tropospheric NO 2 columns retrieved using the surface LER and DLER datasets for a given day and for the monthly average in February and August 2010. The daily differences in the tropospheric NO 2 columns are consistent with Fig. 1c, with a larger impact found over polluted regions. Taking Spain on 3 February 2010 as an example, the approximately 0.005 smaller surface DLER in the central part of the orbit swath results in a lower sensitivity to tropospheric NO 2 columns in the AMF calculation; therefore, the AMF decreases and the tropospheric NO 2 columns increases by ∼ 1×10 14 molec. cm −2 (3 %). In contrast, the surface DLER is approximately 0.02 higher on the western side of the orbit swath over eastern China on the same day; thus, the tropospheric NO 2 column is approximately 3 × 10 15 molec. cm −2 (11 %) lower. The monthly differences in the tropospheric NO 2 columns show a larger reduction in winter: by more than 5×10 14 molec. cm −2 over regions such as central Europe, South Africa, India, and eastern China, and by ∼ 1 × 10 14 molec. cm −2 over areas such as the eastern US, Southeast Asia, and Mexico.
Accounting for the BRDF effect, the surface DLER captures the cross-track dependency of surface albedo, such as the increased reflectivity in the backward scattering viewing geometries, which is in agreement with studies that have applied the BRDF product from MODIS to describe the dependency of land surface reflectance on illumination and viewing geometry (e.g. Zhou et al., 2010;Noguchi et al., 2014;Vasilkov et al., 2017;Lorente et al., 2018;Laughner et al., 2018;Qin et al., 2019). With a good agreement with the established MODIS BRDF product , both in the absolute sense and in the directional/angular dependency, the GOME-2 DLER dataset is derived from measurements of the instrument itself, which is consistent with the GOME-2 NO 2 observations, considering the illumination conditions, observation geometry, and instrumental characteristics; therefore, the use of GOME-2 DLER introduces no additional bias caused by the instrumental differences.

Over water
The surface reflectivity over water is described with improved GOME-2 LER data using an ocean surface albedo parameterisation (Jin et al., 2004(Jin et al., , 2011 to account for the direction dependency. Based on atmospheric radiation measurements and the Coupled Ocean and Atmosphere Radiative Transfer (COART) model (Jin et al., 2006), the parameterisation developed by Jin et al. (2011) derives the surface reflectivity for the direct (sun) and diffuse (sky) incident radiation separately and further divides each of them into contributions from surface and water respectively. This parameterisation has been used to derive ocean surface albedo (e.g. Séférian et al., 2018) and to generate satellite NO 2 products (e.g. Laughner et al., 2018).
Following Jin et al. (2011), the ocean surface albedo α total is defined as follows: where α s dir and α s dif are the respective direct and diffuse contribution of the surface reflection, and α w dir and α w dif are the respective direct and diffuse contribution of the volume scattering of water below the surface. The respective direct and diffuse fraction of downward surface flux, f dir and f dif (f dif = 1 − f dir ), are calculated using the online COART model (https://satcorps.larc.nasa.gov/jin/coart.html, last access: 12 February 2020). The direct surface albedo α s dir , which is one main component of the total ocean surface Figure 3. Differences in the GOME-2 tropospheric NO 2 columns retrieved using the GOME-2 LER and DLER datasets for a given day and for the monthly average in February and August 2010. Only measurements with a cloud radiance fraction less than 0.5 are included. albedo, describes the contribution of Fresnel reflection depending on the incident angle, the refractive index of seawater (1.343 at 460 nm), and the slope distribution of the ocean surface (defined by Cox and Munk (1954) and related to a wind speed of 5 m s −1 from the climatological mean). The diffuse surface albedo α s dif is difficult to formulate analytically due to its variation with atmospheric conditions; thus, it is parameterised practically to be 0.06 for an assumed wind speed of 5 m s −1 . The direct water volume albedo α w dir is considered for Case 1 waters (open ocean waters dominated by phytoplankton and associated products that comprise 99 % of the ocean) and is primarily affected by the chlorophyll concentration (0.2 mg m −3 from the global ocean average). The diffuse water volume albedo α w dif is defined by the α w dir at an effective incident direction, i.e. arccos(0.676), and is calculated to be 0.0145. The direct fraction of downward surface flux f dir is calculated using a radiative transfer simulation with a mid-latitude summer atmosphere, a marine aerosol optical depth of 1 (at 550 nm), and a 100 m ocean depth with the average Petzold phase function for ocean particle scattering. Figure 4 shows the parameterised ocean surface albedo for a non-glint situation as well as its four albedo components as a function of incident angle. The overall shape of the total ocean surface albedo α total is dependent on the incident angle with a peak near 70 • , which is similar to Jin et al. (2004) and Laughner et al. (2018). The surface component (α s dir + α s dif ) is larger than the water volume component (α w dir + α w dif ), particularly for larger incident angles. The direct component (α s dir + α w dir ) increases with incident angle and has lower values than the diffuse component (α s dif + α w dif ) for smaller inci- . Parameterised ocean surface albedo for a non-glint condition and its albedo components due to direct and diffuse surface reflection and direct and diffuse water volume scattering as a function of incident angle. dent angles (below 55 • ) and higher values for larger incident angles. The relative contribution of the diffuse component to the total ocean surface albedo f dif increases from ∼ 0.65 to ∼ 1 with incident angle. It is worth noting that the four albedo components are independent of each other and are, thus, flexible to update or replace. Based on measurements over a long period (2007-2018 for version 3.1), the GOME-2 LER climatology mainly provides the diffuse component (α s dif + α w dif ) over water bodies with minimised impact of direct contribution. Therefore, we replace the simplified expression of α s dif + α w dif in Jin et al. (2011) with values taken from the GOME-2 LER climatology. This scheme enables the consideration of the direction dependency for the GOME-2 LER climatology over water with minimal bias introduced. In addition, most of the ocean surface albedo studies (e.g. Ohlmann, 2003;Jin et al., 2004Jin et al., , 2011Li et al., 2006;Laughner et al., 2018) employ a straightforward assumption that the SZA is the only directional parameter involved in the parameterisation, specifically the incident angle is assumed to be equivalent to the SZA in the Fresnel reflection calculation. In this work, we apply the full equation to derive the local incident angle that also considers the dependencies on VZA and RAA, and we additionally implement the Cox-Munk sun glitter model over glint-contaminated regions. See Cox and Munk (1954) and Gordon (1997) for more details on the configuration and derivation. Figure 1b and d present the calculated ocean surface albedo and the differences with respect to values taken from GOME-2 LER climatology on 3 February and 5 August 2010. Consistent with Vasilkov et al. (2017), the improved ocean surface albedo shows up to 0.015 higher values at larger SZAs and VZAs, where the higher incident angles result in stronger Fresnel reflections, and up to 0.025 higher values over areas affected by sun glint, which is typically the eastern swath of GOME-2 orbits. Figure 5 shows the impact of using updated ocean surface albedo on our GOME-2 NO 2 retrieval for a given day and for the monthly average in February and August 2010. The tropospheric NO 2 columns are mainly reduced over the polluted coastal regions that have large NO 2 concentrations and for large SZAs and VZAs. For instance, the ocean surface albedo around Spain increases by ∼ 0.01 on 3 February 2010, leading to a decrease in the tropospheric NO 2 columns by up to 8 × 10 14 molec. cm −2 (9 %). The monthly average of the tropospheric NO 2 columns decreases in winter by more than 3 × 10 14 molec. cm −2 near the coastal area, e.g. around the US, eastern China, and Brazil, and by up to 1 × 10 14 molec. cm −2 along the shipping lanes, e.g. in the Mediterranean Sea, the Red Sea, and maritime Southeast Asia.

A priori NO 2 profile
In regions with strong gradients with respect to NO x emissions in space and time, the significant variation of surface NO 2 can only be captured in a model with sufficient horizontal, vertical, and temporal resolution. The advanced IFS (CBA) (Huijnen et al., 2016(Huijnen et al., , 2019 global chemistry forecast and analysis system combines the stratospheric chemistry scheme developed for the Belgian Assimilation System for Chemical ObsErvations (BASCOE, Skachko et al., 2016) and the modified CB05 tropospheric chemistry scheme (Williams et al., 2013). As summarised in Table 2, the spatial resolution of IFS (CBA) is a reduced Gaussian grid at a spectral truncation of T255, which is equivalent to a grid spacing of ∼ 80 km globally (∼ 0.7 • ). The model is run with the standard 137 hybrid sigma-pressure layers as also used operationally in the forecast and reanalysis model from the ECMWF. From this, we select a vertical discretisation based on 69 vertical levels up to 0.1 hPa with ∼ 12 layers in the boundary layer for further processing. An essential difference compared with TM5-MP is that the chemistry is an integral part of the meteorological forecast model in IFS (CBA).
Here we use the forecast model from cycle 45r2, which is initialised daily using ERA5 meteorology. Additionally, anthropogenic emissions are based on the recently prepared CAMS_GLOB_ANT v2.1 emission inventory (Granier et al., 2019), as illustrated in Appendix A, and day-specific biomass burning emissions are taken from GFASv2.1 (Kaiser et al., 2012). The NO 2 data are available on an hourly basis, with profiles at the satellite measurement time obtained by linear interpolation. Figure 6 shows an intercomparison of area-averaged and monthly averaged profiles from TM5-MP and IFS (CBA) at the GOME-2 overpass time (09:30 LT) over western Europe (44-53 • N, 0-7 • E) and eastern China (21-41 • N, 110-122 • E) in February and August 2010. Generally, TM5-MP and IFS (CBA) show similar mean profile shapes over the two regions. In February, IFS (CBA) shows a larger boundary layer concentration and a sharper transition to the free troposphere over western Europe and larger NO 2 gradients in the free troposphere over eastern China. In August, the NO 2 concentrations in the free troposphere are lower than in February for both models due to the reduced emissions and the reduced lifetime of NO 2 , and a larger surface layer NO 2 gradient is found for the IFS (CBA) model for both regions. Figure 7 shows the daily TM5-MP and IFS (CBA) a priori NO 2 profiles over the Netherlands (52.8 • N, 4.7 • E) and China (39.1 • N, 118.0 • E) on 3 February 2010 as examples. IFS (CBA) shows a higher surface layer NO 2 concentration (a steeper profile shape) and yields a tropospheric AMF that is reduced by 0.21 over the Netherlands, which will enhance the retrieved tropospheric NO 2 column. In contrast, the tropospheric AMF increases by 0.06 over China due to the larger NO 2 gradients in the free troposphere (a less steep profile shape) modelled by IFS (CBA). Figure 8 shows the differences in the tropospheric NO 2 columns retrieved using TM5-MP and IFS (CBA) a priori NO 2 profiles for a given day and for the monthly average in February and August 2010. The differences are consistent with the changes in the profile shapes in Figs. 6 and 7. The use of IFS (CBA) generally increases the tropospheric NO 2 columns over polluted regions, e.g. over western Europe, the eastern US, and Argentina, by up to 2 × 10 15 molec. cm −2 and decreases the values over areas such as central Africa, South Africa, and Brazil, by up to 1 × 10 15 molec. cm −2 . In February, however, a strong enhancement of ∼ 7 × 10 15 molec. cm −2 is found over northern Germany and Poland, and a strong reduction of ∼ 4 × 10 15 molec. cm −2 is found over the North China Plain. The Figure 5. Differences in the GOME-2 tropospheric NO 2 columns retrieved using the original GOME-2 LER climatology and the GOME-2 LER data improved by the ocean surface albedo parameterisation for a given day and for the monthly average in February and August 2010. Only measurements with a cloud radiance fraction less than 0.5 are included. To quantify the effect of model resolution, a more detailed analysis was implemented for IFS (CBA) that used 1 • grids for the horizontal resolution, 34 layers for the vertical resolution, and a 2 h time step for the temporal resolution. These values are of the same order of magnitude as the model resolution of TM5-MP and other chemistry transport models currently employed in the satellite retrieval of global NO 2 (e.g. van Geffen et al., 2019;Lorente et al., 2017;Boersma et al., 2018;Liu et al., 2019b). AMFs differ by more than 0.02 for both example locations due to differences in the horizontal and vertical resolutions. The current 2 h temporal sampling and subsequent linear interpolation between the sampling points is sufficient for the retrieval of tropospheric NO 2 columns (not shown). When   Fig. 7. Only measurements with a cloud radiance fraction less than 0.5 are included. a coarser spatial resolution is used, the "domain-averaged" profiles generally show an increased surface NO 2 concentration for the unpolluted domain and the opposite for the emission source. Consequently, the AMF is underestimated for unpolluted areas and overestimated for polluted areas. When the number of layers is reduced, the coarser sampling points can not accurately capture the large NO 2 gradients at low altitudes, particularly for polluted regions where the measurement sensitivity of the satellite decreases significantly towards the surface. Figure 9 shows the absolute and relative differences in tropospheric NO 2 columns retrieved by altering the model resolutions for IFS (CBA) a priori NO 2 profiles in February 2010. In Fig. 9a, the increase in the spatial resolution (1 • vs. 0.7 • grid) changes the tropospheric NO 2 columns by up to 7 × 10 14 molec. cm −2 or 20 % for polluted regions, confirming the importance of applying a priori NO 2 profiles with better spatial resolution (Heckel et al., 2011;Lin et al., 2014;Kuhlmann et al., 2015). Larger relative differences are observed over cities surrounded by rural ar-eas, coastal regions, isolated islands, and shipping lanes, where the use of high spatial resolutions more accurately captures the NO 2 emissions and chemistry for a priori profiles. In Fig. 9b, the improvement in the vertical resolution (34 vs. 69 layers) enhances the tropospheric NO 2 columns by up to 5 × 10 14 molec. cm −2 or 15 %. Increasing the number of layers generally better resolves the NO 2 vertical variation, especially for the lowest model layers where the box-AMF decreases significantly in the polluted cases. Consequently, the tropospheric AMFs are lower and the tropospheric NO 2 columns are higher for polluted regions (Heckel et al., 2011;Russell et al., 2011;Valin et al., 2011). For unpolluted regions, the differences are generally small (within ±2 × 10 14 molec. cm −2 or ±3 %). In addition, the use of different temporal resolution (a 2 h vs. a 1 h time step) generally has a negligible impact on the tropospheric NO 2 columns (less than 2 × 10 14 molec. cm −2 or 3 %, not shown).

Cloud correction
For cloudy scenarios, the retrieval of tropospheric NO 2 is affected by the cloud parameters due to the variation in scene albedo and the photon path redistribution in the troposphere. As described in Sect. 2, the cloudy-sky AMFs are calculated by the independent pixel approximation using GOME-2 cloud parameters: radiometric cloud fraction from OCRA as well as cloud top pressure (cloud top height) and cloud albedo (cloud optical depth) from ROCINN. To improve the cloud correction in our NO 2 retrieval, the CAL model from the ROCINN cloud algorithm , where the clouds are treated as optically uniform layers of light-scattering water droplets, is applied. The CAL model is more representative of the real situation than the CRB model (where the clouds are idealised as Lambertian reflectors with zero transmittance), as it allows the penetration of photons through the cloud layer. This treatment takes the multiple scattering of light inside the cloud and the contribution of the atmospheric layer between the cloud bottom and the ground into account. Figure 10 shows the differences in cloud top heights obtained with the CRB and CAL models from GOME-2 measurements in February and August 2010. Consistent with Loyola et al. (2018) (see Figs. 3 and 13 therein), the cloud top heights from CAL are generally higher by ∼ 0.9 km on average. Stronger increases (up to 2 km) are found over regions with thick and high clouds, such as the Intertropical Convergence Zone and the Western Hemisphere Warm Pool, which is very similar to the findings of Lelli et al. (2012;see Fig. 12 therein). In general, the CRB-based cloud retrieval underestimates the cloud top height due to the neglect of oxygen absorption throughout a cloud layer (Vasilkov et al., 2008) and, thus, the misinterpretation of the smaller top-of-atmosphere reflectance as a lower cloud layer (Saiedy et al., 1967). The retrieved cloud height is normally close to the middle, i.e. the optical centroid of clouds (Ferlay et al., 2010;Richter et al., 2015), which can be considered to be a reflectance-weighted height located inside a cloud. Additionally, as the enhanced multiple scattering is not fully taken into account in the CRBbased cloud retrieval, the retrieved cloud height is normally close to the altitude of the middle.
In Fig. 10, higher cloud top heights are mainly found over land surfaces characterised by the presence of a large amount of absorbing aerosols when using CRB, i.e. over regions with strong desert dust emissions, such as the Sahara, the Arabian Desert, and the deserts in Australia, as well as regions with strong biomass burning emissions, such as South America, South Africa, and Southeast Asia. Over these areas, ROCINN likely retrieves an effective aerosol height close to the top of aerosol layer, depending on the type of absorbing aerosols and on the aerosol optical depth. The presence of strongly absorbing aerosols, which typically have large aerosol optical depth and/or are located at high altitudes (up to ∼ 8 km), reduces the fraction of photons reaching the lowest part of the atmosphere. In order to approximate this shortened light path, the CRB-based cloud retrieval has to put the Lambertian reflector at a higher altitude Chimot et al., 2016). This effect is larger for aerosol layers at higher altitudes and is also dependent on geometry parameters like SZA, on surface properties such as surface albedo, and on the accuracy of radiometric cloud fractions from OCRA. See Sect. 4 for further discussion on the CRB and CAL cloud models in the presence of aerosols.
To apply the CAL cloud model in our NO 2 AMF calculation, a single scattering albedo of 1 and an asymmetry parameter of 0.85 for water clouds are assumed for the radiative transfer calculation; these values are consistent with those used in the cloud retrieval . Cloud observations with a fitting root mean square (rms) greater than 1 × 10 −4 or a number of iterations greater than 20 are filtered out. The NO 2 box-AMFs are derived through the pixelspecific radiative transfer calculation instead of the interpolation from a LUT with fixed reference points. This requires no projection from the layer coordinate of the NO 2 profile to the coordinate assumed in the LUT and also requires no linear interpolation based on the model parameters. Figure 11 shows an example of the box-AMFs derived for clear and cloudy sky using the CRB and CAL models over Italy (45.3 • N, 11.2 • E) on 1 February 2010. The cloud information and the calculated tropospheric AMFs are also reported. Compared with the clear-sky box-AMFs, the CALbased cloudy-sky box-AMFs increase above the cloud layer (albedo effect) and decrease below the cloud layer (shielding effect). Compared with the CRB model, the use of the CAL model considers the sensitivities inside and below the cloud layer and increases the cloudy-sky AMF by 0.3, which consequently decreases the retrieved tropospheric NO 2 column by 2.5×10 15 molec. cm −2 (12 %), based on the polluted NO 2 profile with most of the NO 2 concentration located near the surface. Only measurements with a cloud radiance fraction less than 0.5 are included. Figure 10. Differences in cloud top heights retrieved using the ROCINN_CAL and ROCINN_CRB cloud models in February and August 2010. Only measurements with a cloud fraction less than 0.3 are included. Observations with a fitting rms greater than 1 × 10 −4 or a number of iterations greater than 20 are filtered out. Figure 12 shows the differences in the tropospheric NO 2 columns retrieved using the CRB and CAL models in February and August 2010. The use of the CAL model decreases the tropospheric NO 2 columns by more than 1 × 10 14 molec. cm −2 over polluted regions. Larger values are found in winter (up to 3 × 10 15 molec. cm −2 ), when most of the NO 2 concentrations are located at the surface and the cloud fractions are generally larger due to the seasonal variation of clouds. Figure 13 shows the tropospheric NO 2 columns retrieved using the improved AMF calculation and the differences from the reference data in February and August 2010. Larger differences are found in winter over the polluted regions. For instance, the tropospheric NO 2 columns are reduced by more than 1 × 10 15 molec. cm −2 over China and India in February as well as Brazil and South Africa in August. Increased val- Figure 11. The box-AMFs for clear and cloudy sky using the ROCINN_CAL and ROCINN_CRB cloud models over Italy (45.3 • N, 11.2 • E) on 1 February 2010. The tropospheric AMF is given next to each label. The ROCINN_CRB cloud top pressure is shown as a green horizontal dotted line, and the ROCINN_CAL cloud top and base pressures are shown as brown horizontal dotted lines. The cloud radiance fraction is 0.47, the cloud optical depth is 6.85, the SZA is 69 • , the VZA is 3 • , and the RAA is 42 • . ues are found over regions such as Mexico, Argentina, and Russia. Table 3 summarises the individual changes and the combined effect of our improved AMF calculation on the retrieved tropospheric NO 2 columns over western Europe, eastern China, the eastern US, and central Africa. Increases in the GOME-2 surface albedo reduce the tropospheric NO 2 columns by 2 %-6 %. The use of IFS (CBA) a priori NO 2 profiles affects (mostly increases) the tropospheric NO 2 columns by up to 21 %, and the use of CAL cloud parameters decreases the values by up to 14 %. The combined effect of individual improvements yields to an average change in the tropospheric NO 2 columns of within ±15 % in winter and ±5 % in summer over polluted regions.

Combined impact
The uncertainty in the improved AMF calculation is likely reduced relative to the reference retrieval if the improved surface albedo dataset, a priori NO 2 profiles, and cloud parameters are considered, as they are the main causes of AMF uncertainty (Lorente et al., 2017). The uncertainty in the AMF calculation for polluted conditions is estimated to have improved from 10 % to 45 % for the reference retrieval to between 10 % and 35 % in this work.

Implicit aerosol correction
Aerosols can increase or decrease the sensitivity to tropospheric NO 2 , depending on the NO 2 and aerosol vertical distribution as well as the optical and physical properties of the particles (Martin et al., 2003;Leitão et al., 2010). As the OCRA/ROCINN cloud retrieval does not distinguish between clouds and aerosols, the aerosol effect is assumed to be corrected implicitly in the AMF calculation via the effective cloud parameters (i.e. aerosols are treated as clouds). Figures 14 and 15 show the land surface RGB image with active fire locations from the Moderate Resolution Imaging Spectroradiometer (MODIS) on board the Terra (10:30 LT) and the OCRA/ROCINN cloud products measured by GOME-2 (09:30 LT) over eastern China and central Africa on a given day respectively. The MODIS dataset (https://worldview.earthdata.nasa.gov/, last access: 12 February 2020) describes the cloud or aerosol amount and fire locations (for central Africa). For both regions, large aerosol loads are found in the RGB image for cloud-free areas, such as the Beijing-Tianjin-Hebei economic region in eastern China and burning regions across central Africa, where the aerosol loads are identified as thin clouds (a cloud optical depth of ∼ 5) near the surface (a cloud top height of ∼ 3 km) with cloud fractions of up to 0.18. Therefore, we assume that the thin clouds near the surface in the OCRA/ROCINN cloud products are attributed to aerosol loads for measurements with a cloud radiance fraction less than 0.5 or a cloud fraction less than 0.2, and we evaluate the accuracy of the implicit aerosol correction by comparing it with the explicit aerosol correction. For this purpose, the explicit correction for aerosols is implemented using ground-based aerosol observations at Xianghe (39.75 • N, 116.96 • E), which is a suburban site surrounded by heavily industrialised areas in northeastern China (Clémer et al., 2010;Hendrick et al., 2014;Vlemmix et al., 2015), and at Bujumbura (3.38 • S, 29.38 • E), which is located in the central African country of Burundi and surrounded by intensive biomass burning activities (Gielen et al., 2017), as indicated in Figs. 14 and 15 respectively. Our analysis is further limited to satellite measurements with a cloud optical depth less than 5 and a cloud top height less than 3 km to reduce the cloud contamination. With this selection, the aerosol concentrations are generally low or moderate (an aerosol optical depth less than 1).
The explicit modelling of aerosol scattering and absorption for the AMF calculation is implemented by introducing the aerosol optical properties (i.e. single scattering albedo and phase function) and vertical distributions (i.e. extinction vertical profiles) into the radiative transfer calculation. The single scattering albedo describes the fraction of the aerosol light scattering over the extinction, and the phase function describes the angular distribution of the scattered light intensity. In this study, we apply the Henyey-Greenstein phase function with an asymmetry parameter (the first moment of phase function) describing the asymmetry between forward and backward scattering. The long-term statistics of the single scattering albedo and the asymmetry parameter at 440 nm are derived for Xianghe and Bujumbura using the Version 3 Level 2.0 inversion products from the sun photometer radiance measurements from AERONET (Holben et al., 1998;Giles et al., 2019) (http://aeronet.gsfc.nasa.gov/, last access: 12 February 2020). Monthly mean parameters are calcu- Figure 12. Differences in the GOME-2 tropospheric NO 2 columns retrieved using the ROCINN_CRB and ROCINN_CAL cloud models in February and August 2010. Only measurements with a cloud radiance fraction less than 0.5 are included. Cloud observations with a fitting rms greater than 1 × 10 −4 or a number of iterations greater than 20 are filtered out. Figure 13. GOME-2 tropospheric NO 2 columns retrieved using the improved algorithm, and the differences from the reference data in February and August 2010. Only measurements with a cloud radiance fraction less than 0.5 are included. Cloud observations with a fitting rms greater than 1 × 10 −4 or a number of iterations greater than 20 are filtered out. lated based on up to 7 years of observations (2010-2016 for Xianghe and 2013-2016 for Bujumbura) that are available within ±1 h of the GOME-2 overpass time (09:30 LT) for each month.
Xianghe is located ∼ 60 km southeast of Beijing and belongs to the highly urbanised Beijing-Tianjin-Hebei economic region on the North China Plain; this area experiences heavy anthropogenic aerosol emissions, especially in winter due to enhanced domestic heating. Mixtures of desert dust and urban-industrial aerosols mainly affect the region in spring (March-May). Based on the monthly climatology of AERONET measurements, the single scattering albedo at Xianghe is 0.91 on average with a maximum in July (0.96) and low values in winter (∼ 0.87), which are mostly related to the black carbon emissions (Yang et al., 2011). The asymmetry parameter ranges between 0.7 and 0.75, which is consistent  with values from the urban aerosol models in East Asia (Lee and Kim, 2010). Bujumbura is located in tropical Africa; this is an area that is typically affected by biomass burning emissions (mainly during the local dry seasons -June to August and January to February) and, to a lesser extent, by anthropogenic emissions (throughout the year with negligible seasonal variations). The single scattering albedo at Bujumbura is higher between March and May (∼ 0.9), which is related to the major wet season, and lower between July and August and be- tween December and January (0.83-0.87), which coincides with intensive agricultural activities and the transport of forest fire emissions in the surrounding regions (Gielen et al., 2017). The asymmetry parameter is 0.69 on average, which is in agreement with values in biomass burning aerosol models (Torres et al., 2013).
The colocated aerosol extinction vertical profiles at 477 nm are taken from the MAXDOAS measurements at Xianghe from March 2010 to December 2016 (Clémer et al., 2010) and at Bujumbura from December 2013 to December 2015 (Gielen et al., 2017). The MAXDOAS data are used to derive aerosol information based on the oxygen collision complexes (O 4 ) absorption, as the O 4 vertical profile is generally constant and, thus, capable of describing the influence of aerosol scattering and absorption on the photon path (Wagner et al., 2004;Frieß et al., 2006). The MAXDOAS technique can reliably determine the aerosol extinction profiles in the lower troposphere (Frieß et al., 2016), where most aerosols are located over Xianghe and Bujumbura. We colocate the space-and ground-based measurements by selecting GOME-2 pixels within 50 km of the stations and averaging the MAXDOAS aerosol profiles within ±1 h of the GOME-2 overpass time (09:30 LT). Figure 16 shows typical NO 2 box-AMFs, the simulated TM5-MP NO 2 profile, and the MAXDOAS aerosol extinction profile for Xianghe on 21 November 2013. The MAX-DOAS aerosol profile follows an exponentially decreasing shape with a peak of aerosol loads close to the ground (950 hPa or 0.4 km). The NO 2 follows the same profile shape and is well mixed with aerosol. Depending on the seasonal variation, local emissions, and transport process, aerosol profiles with a peak at elevated heights (up to 900 hPa or 1 km) are also observed (not shown) due to a long residence time. The discontinuity of box-AMFs corrected using the CRB cloud model is introduced by the effective clouds (see Eq. 2), below which the cloudy box-AMFs are zero. Due to the overestimated cloud altitudes from the CRB-based cloud retrieval (see Sect. 3.3), the CRB-based implicit aerosol correction underestimates the tropospheric AMF by 14 %; this bias is largely reduced by applying the CAL cloud model (6 %), which causes a gradual reduction in box-AMFs towards the surface and shows better agreement with the shape from the explicit aerosol correction. Figure 17 shows the same data as Fig. 16 but for Bujumbura on 9 September 2015. Compared with the data at Xianghe, the aerosol profile at Bujumbura shows a lower aerosol amount but an uplifted layer of aerosol loads (820 hPa or 1.8 km), whereas NO 2 continues to peak at the surface. The difference in AMF between the implicit and explicit aerosol corrections decreases from 15 % using CRB to 5 % using CAL. Figure 18 presents the relative biases in the tropospheric NO 2 columns retrieved assuming no aerosol correction (i.e. applying the clear-sky AMFs) and assuming implicit aerosol correction via the CRB and CAL cloud models for Xianghe from March 2010 to December 2016. Only measurements with a cloud radiance fraction less than 0.5, a cloud optical depth less than 5, and a cloud top height less than 3 km are included. The relative biases introduced by assuming no aerosol correction vary between −30 % and 31 % with an average of 7 % for GOME-2 pixels, which is in agreement with previous studies focusing on the industrialised part of eastern China (Ma et al., 2013;Lin et al., 2014Lin et al., , 2015Kuhlmann et al., 2015;Chimot et al., 2016;Wang et al., 2017). Due to the overestimated shielding effect, the tropospheric NO 2 columns retrieved using the CRB-based implicit aerosol correction are 33 % larger on average than those retrieved using the explicit aerosol correction, and the differences are largely reduced by applying the CAL cloud model (9 %). Enhanced differences are found for larger cloud radiance fraction values, which is probably due to the increased pollution level (NO 2 columns) compared with clear sky conditions , as the cloud (radiance) fraction is highly correlated with the MAXDOAS aerosol optical depth (correlation coefficient of 0.7 and regression slope of 0.17, not shown). Figure 19 shows the same data but for Bujumbura from December 2013 to December 2015. The explicit aerosol correction yields tropospheric NO 2 columns that are 6 % smaller than the clear-sky tropospheric NO 2 columns on average, which is consistent with Martin et al. (2003) and Castellanos et al. (2015). The average difference between the tropospheric NO 2 columns from the implicit and explicit aerosol corrections decreases from 15 % using the CRB model to 5 % using the CAL model.
In Figs. 18 and 19, the relative biases introduced by the CAL-based implicit aerosol correction are close to the values assuming no aerosol correction, addressing the complexities related to the tropospheric NO 2 measurements in the presence of aerosols. In combination with the cloud model error, errors related to the implicit aerosol correlation can result from the different radiative effects of scattering clouds and absorbing aerosols as well as the different characteristic sizes and phase functions of clouds and aerosols in general. The errors may be additionally enhanced in the presence of actual clouds. Therefore, future works include the further quantitative interpretation of OCRA/ROCINN cloud parameters for aerosol-dominated scenes and the impact on the NO 2 retrieval algorithm.

Tropospheric NO 2 validation
A validation of our improved GOME-2 tropospheric NO 2 columns is based on BIRA-IASB ground-based MAXDOAS measurements at the Xianghe, Uccle, Bujumbura, and Observatoire de Haute Provence (OHP) stations, as introduced in Appendix B. Xianghe station, owing to its polluted suburban nature, is the best site for validation (Liu et al., 2019b). The satellite-based NO 2 data are commonly underestimated for urban polluted stations like Uccle, due to the averaging of a local source over a pixel size (80 km×40 km/40 km×40 km for GOME-2) larger than the horizontal sensitivity of the ground-based measurements which is approximately a few kilometres to tens of kilometres (Irie et al., 2011;Wagner et al., 2011;Ortega et al., 2015). Difficulties may arise when small local sources are present in remote locations, such as the Bujumbura (urban station in Burundi) or OHP (background station in southern France) sites (Pinardi et al., 2015;Gielen et al., 2017). For the validation of GOME-2 measurements, the satellite data are filtered for clouds (a cloud radiance fraction less than 0.5), and the closest valid pixel within 50 km of the stations is compared to the ground-based MAX-DOAS data, which are linearly interpolated to the GOME-2 overpass time (09:30 LT), if original data exist within ±1 h. Figure 20 shows the time series and scatter plot of the comparison of the daily and monthly means between the improved GOME-2 tropospheric NO 2 columns and the groundbased MAXDOAS measurements at Xianghe, including the statistical information on the correlation coefficient, slope, and intercept of orthogonal regression analysis. The monthly mean values from the GOME-2 and MAXDOAS measurements indicate good agreement with similar seasonal variations in the tropospheric NO 2 column. A correlation coefficient of 0.94, a regression slope of 0.69, and an intercept of 0.41 × 10 15 molec. cm −2 are derived when comparing the monthly mean values. These results are qualitatively similar to previous validation exercises at other sites, from other satellites, and using other NO 2 products (Celarier et al., 2008;Kramer et al., 2008;Chen et al., 2009;Irie et al., 2012;Ma et al., 2013;Wu et al., 2013;Kanaya et al., 2014;Lamsal et al., 2014;Wang et al., 2017;Drosoglou et al., 2017Drosoglou et al., , 2018. Similar figures for previous GDP products can be found in Liu et al. (2019b) and on the AC-SAF validation website (http://cdop.aeronomie.be/ validation/valid-results, last access: 12 February 2020). Figure 21 presents the daily and monthly mean absolute and relative differences of GOME-2 and MAXDOAS measurements at Xianghe. The differences are within ±1 × 10 16 molec. cm −2 on average with a mean difference of −3.8×10 15 molec. cm −2 . The NO 2 levels are underestimated by 9.9 % by GOME-2 with a standard deviation of ±21 %, which is mostly explained by the relatively low sensitivity of spaceborne measurements near the surface, the gradientsmoothing effect, and the aerosol shielding effect. These effects are often inherent to the different measurement types or Atmos. Meas. Tech., 13, 755-787, 2020 www.atmos-meas-tech.net/13/755/2020/  the specific conditions of the validation sites and also to the remaining impact of structural uncertainties (Boersma et al., 2016), such as the impact of the choice of the a priori NO 2 profiles and/or the albedo database assumed for the satellite AMF calculations. Similar figures to Figs. 20 and 21 for the reference retrieval can be found in Liu et al. (2019b), and figures for previous GDP products can be found on the AC-SAF validation website (http://cdop.aeronomie.be/validation/valid-results, last access: 12 February 2020). Tables B1-B4 in Appendix B summarise the statistics for the improved retrieval, the reference retrieval, and the operational GDP 4.8 product at Xianghe, Uccle, Bujumbura, and OHP respectively. As discussed in Pinardi et al. (2015), for remote (Bujumbura) and background (OHP) stations, the mean bias is considered to be the best indicator of the validation results due to the relatively small variability in the measured NO 2 . In urban (Uc-cle) and suburban (Xianghe) situations, the NO 2 variability is large enough, and the correlation coefficient provides a good indication of the coherence of the satellite and ground-based datasets, although larger differences in terms of slope and mean bias can be expected for the urban case because satellite measurements smooth out the local NO 2 hot spots. Compared with the reference retrieval, the improvement of the algorithm leads to an increase in the correlation coefficient from 0.9 to 0.94 under suburban (Xianghe) conditions and a decrease in the mean relative bias from −43 % to −37 % under urban (Uccle) conditions; small impacts on the mean bias are found for Bujumbura and OHP. Compared with the operational GDP 4.8 product, however, both the reference retrieval and the improved retrieval show a significant improvement for all stations.
As Xianghe is located in a highly polluted region, the presence of large aerosol loads increases the uncertainty on the Figure 18. Scatter plot of the relative biases in the GOME-2 tropospheric NO 2 columns assuming no aerosol correction (a) and assuming an implicit aerosol correction using the ROCINN_CRB and ROCINN_CAL cloud models (b) with respect to the cloud radiance fraction at Xianghe station from March 2010 to December 2016. Only measurements with a cloud radiance fraction less than 0.5, a cloud optical depth less than 5, and a cloud top height less than 3 km are included. Cloud observations with a fitting rms greater than 1 × 10 −4 or a number of iterations greater than 20 are filtered out. The mean value and standard deviation are given next to each label in the legend.  Table 4. The mean difference (MD) and median difference (AD) (SAT-GB in 10 15 molec. cm −2 ), standard deviation (SD, in 10 15 molec. cm −2 ), standard error of the mean (SE, in 10 15 molec. cm −2 ), and correlation coefficient R, slope S, and intercept I (in 10 15 molec. cm −2 ) of the orthogonal regression for the daily GOME-2 tropospheric NO 2 product compared to MAXDOAS data at Xianghe. Intermediate results for different surface albedo and a priori NO 2 profiles are reported for completely clear-sky conditions (a cloud radiance fraction of 0) for a total of 73 GOME-2 pixels.  Table 5. Similar to Table 4 but for different aerosol corrections for aerosol-dominated conditions (a cloud radiance fraction less than 0.5, a cloud optical depth less than 5, and a cloud top height of less than 3 km) for a total of 146 GOME-2 pixels. Results are calculated using DLER surface albedo and IFS (CBA) a priori profiles.   and increases the complexity of validation. As the cloud retrieval can hardly distinguish between clouds and aerosols (see Sect. 4), the "typical" validation results in Table B1 (a cloud radiance fraction less than 0.5) do not show improvements e.g. in the mean bias  . Daily (grey dots) and monthly mean (black dots) absolute and relative GOME-2 (SAT) and MAXDOAS (GB) time series differences for Xianghe station. The histogram of the daily differences is also given, showing the mean and median difference. The total time series of absolute and relative monthly differences are given in the bottom-right panel. and the slope of the regression line when compared to the reference retrieval. To separate this effect, the validation at Xianghe is differentiated under completely clear-sky conditions in Table 4 (a cloud radiance fraction of 0) and aerosoldominated cases in Table 5 (a cloud radiance fraction less than 0.5, a cloud optical depth less than 5, and a cloud top height less than 3 km).
To summarise the improvements of each of the changes discussed in previous sections, Table 4 reports the statistical results including the biases and regression analyses for the use of different surface albedo and a priori NO 2 profiles at Xianghe station for completely clear-sky conditions (a cloud radiance fraction of 0). Compared with the reference retrieval (based on GOME-2 surface LER climatology and the TM5-MP a priori profile), better results are obtained with the improved algorithm (based on the surface DLER dataset and the IFS (CBA) a priori profile) with a median difference of −1.0 × 10 15 molec. cm −2 , which will be used to further test for aerosol correction type below. Table 5 presents the statistical results for the retrievals with no aerosol correction and the CAL-based implicit aerosol correction at Xianghe station for aerosol-dominated cases (a cloud radiance fraction less than 0.5, a cloud optical depth less than 5, and a cloud top height less than 3 km). Consistent with Fig. 18, the GOME-2 NO 2 columns retrieved using the CAL-based implicit aerosol correction are higher than the results assuming no aerosol correction, which improve the biases relative to the MAXDOAS measurements (e.g. median difference from −2.8 × 10 15 molec. cm −2 to −2.3 × 10 15 molec. cm −2 ), as well as the standard deviation, correlation coefficient, and regression parameters.

Conclusions
The operational GOME-2 NO 2 dataset, generated using the GDP algorithm at DLR, was introduced in detail by Valks et al. (2011Valks et al. ( , 2017 and has been successfully applied in many studies. An improved AMF calculation with a more accurate knowledge of surface albedo, the a priori NO 2 profile, as well as cloud and aerosol correction is described in this paper and is expected to be implemented in an upcoming version of GDP in combination with Liu et al. (2019b).
The viewing angle dependency of surface albedo is taken into account by improving the currently used GOME-2 surface LER climatology (Tilstra et al., 2017). Over land, the surface albedo is described by a GOME-2 DLER dataset , which is determined by dividing the GOME-2 orbit swath into five segments and retrieving the traditional surface LER for each segment based on the data from the respective part of orbit swath. Compared with the nondirectional GOME-2 LER climatology, the use of the DLER dataset improves the underestimation of the surface albedo on the western side of the GOME-2 orbit (backscattering geometry) and increases the AMFs by up to 15 % for polluted regions. Over water, the surface albedo is improved with an ocean surface albedo parameterisation (Jin et al., 2011), in which the albedo is parameterised for the direct and diffuse incident radiation separately. We update the simplified expression of the diffuse contribution with more realistic values taken from the GOME-2 LER data, and we improve the description of the dependency on the viewing direction for the parameterisation. The resulting surface albedo increases over sun glint areas and polluted coastal regions with large SZAs and VZAs, for which the tropospheric NO 2 columns are reduced by up to 10 %.
Higher-resolution a priori profiles, obtained from the IFS (CBA) chemistry transport model with recent emission inventories, provide a better description of the spatial and temporal variability in the NO 2 fields. Compared with the currently used TM5-MP profiles, the application of IFS (CBA) profiles affects the tropospheric NO 2 columns by up to 7 × 10 15 molec. cm −2 for polluted regions, which is mainly due to the differences in the model specifications and model resolutions. To quantify the influence of model resolutions, we implement an analysis in which we alter the horizontal, vertical, and temporal resolutions for IFS (CBA). Changing the horizontal resolution from 1 to 0.7 • affects the tropospheric NO 2 columns by up to 20 %, with enhanced values for emission sources and the opposite for their unpolluted surroundings. When the vertical resolution changes from 34 to 69 layers, the tropospheric NO 2 columns increase by up to 15 % due to the capture of small box-AMFs at low altitudes. Small differences (< 3 %) are found for a temporal resolution (time step) of 1 or 2 h.
The CAL model from the ROCINN cloud algorithm, which treats clouds as uniform layers of water droplets, allows for the penetration of photons through the clouds and provides more realistic cloud parameters than the current CRB model, which treats the clouds as idealised Lambertian reflectors. The application of CAL cloud parameters in the AMF calculation considers the sensitivities inside and below the cloud layers and reduces the tropospheric NO 2 columns by up to 3 × 10 15 molec. cm −2 for polluted regions.
As the cloud retrieval does not distinguish between clouds and aerosols, the aerosol correction is implicitly implemented in the AMF calculation using the cloud parameters. To evaluate the accuracy of the implicit aerosol correction through a cloud model, we explicitly account for the aerosol effect using ground-based aerosol measurements for aerosoldominated conditions. For Xianghe, a suburban site in China with primarily anthropogenic aerosol emissions, and Bujumbura, a remote site in tropical Africa that is typically affected by biomass burning aerosols, aerosol optical properties from AERONET measurements and extinction vertical profiles from correlative MAXDOAS measurements are applied. Assuming the explicit aerosol correction as a reference, the use of an implicit aerosol correction via the CAL cloud model yields a bias that is 24 % smaller than the CRB cloud model for Xianghe and 10 % smaller for Bujumbura.
A validation of the improved NO 2 measurements is performed by comparing the GOME-2 tropospheric NO 2 dataset with BIRA-IASB ground-based MAXDOAS measurements. At the suburban Xianghe station, the GOME-2 NO 2 measurements show similar seasonal variation to the MAX-DOAS dataset with a monthly averaged difference of −9.9 % (−3.8 × 10 15 molec. cm −2 in absolute) and a correlation coefficient of 0.94, which indicates good agreement. The application of the new surface albedo, the a priori NO 2 profile, and the cloud correction in the AMF calculation improves the biases, the correlation coefficients, and the regression parameters for Xianghe.
In the future, further studies focusing on the cloud correction will be implemented due to its importance in the AMF calculation. The BRDF effect on cloud parameters will be considered by implementing the GOME-2 DLER dataset in the cloud retrieval from ROCINN, providing a consistent treatment of surface albedo for both NO 2 and cloud retrieval. Note that the BRDF effect is not discussed for OCRA, because no surface albedo climatology is directly needed, and the correction for VZA dependency has been applied in the cloud fraction retrieval as a proxy for BRDF constellation (see Lutz et al., 2016, Sect. 2.2.2 therein). In addition, the interpretation of the OCRA/ROCINN cloud product for aerosol-dominated scenes and the impact on the NO 2 retrieval algorithm will be further investigated in future studies. Furthermore, the NO 2 algorithm will be adapted to measurements from the TROPOMI instrument with a spatial resolution as high as 7 km × 3.5 km (5.5 km × 3.5 km from August 2019 onwards).
Appendix A: Comparison of anthropogenic NO x emissions Figure A1 shows the emission maps of the MACCity (Granier et al., 2011) and CAMS_GLOB_ANT v2.1 (Granier et al., 2019) inventories for Europe in August 2010. The annual total anthropogenic emissions are 70.8 Tg NO yr −1 for MACCity and 73.7 Tg NO yr −1 for CAMS_GLOB_ANT v2.1. Significant changes in local emission patterns are clear in Fig. A1.

Appendix B: Validation for different MAXDOAS stations
The validation of our improved GOME-2 tropospheric NO 2 columns is based on ground-based MAXDOAS NO 2 measurements at Xianghe station ( Tables B1-B4 summarise the validation of the improved retrieval, the reference retrieval, and the current operational GDP 4.8 product using MAXDOAS measurements at the four BIRA-IASB stations. Xianghe is a typical suburban station in northeastern China, Uccle is an urban station in Belgium, Bujumbura is a small city in a remote region in tropical Africa, and OHP (southern France) is largely rural but occasionally influenced by polluted air masses transported from neighbouring cities.  Table B1. Mean difference (MD, SAT-GB in 10 15 molec. cm −2 ), standard deviation (SD, in 10 15 molec. cm −2 ), standard error of the mean (SE, in 10 15 molec. cm −2 ), relative difference (RD, (SAT-GB)/GB in percent), and correlation coefficient R, slope S, and intercept I (in 10 15 molec. cm −2 ) of the orthogonal regression for the monthly mean GOME-2 tropospheric NO 2 product when compared to MAXDOAS data at the suburban Xianghe station. Values for the improved retrieval (this study) are given, and values for the reference retrieval and the current operational GDP 4.8 product are reported. Results are reported for a cloud radiance fraction less than 0.5 for a total of 77 GOME-2 monthly points. Data availability. The current operational GDP version 4.8 NO 2 data from GOME-2 (Valks et al., , 2017 can be ordered via the FTP server and the EUMETSAT Data Centre (https://acsaf.org/, last access: 12 February 2020). The improved dataset is currently available upon request.
Author contributions. SL and PV developed the retrieval framework. SL processed the data, analysed the results, and contributed to the MAXDOAS validation. GP analysed the MAXDOAS validation results. JX processed the cloud data and provided expertise regarding the aerosol correction. AA and RL provided expertise regarding the OCRA/ROCINN-based cloud correction.
LGT provided the GOME-2 surface albedo dataset and expertise regarding the surface albedo. VH provided the IFS (CB05BASCOE) model data and expertise regarding the a priori profiles. FH provided the MAXDOAS measurements. MVR contributed to the MAXDOAS validation. SL, PV, and GP wrote the paper with comments from all authors.
Competing interests. The authors declare that they have no conflict of interest.
Acknowledgements. This work is funded by the DLR-DAAD Research Fellowships 2015 programme (grant no. 57186656; reference no. 91585186) and is undertaken in the framework of the EUMETSAT AC-SAF project. We acknowledge the Belgian Science Policy Office (BELSPO) for supporting part of this work through the PRODEX programme B-ACSAF project. We thank EUMETSAT for the GOME-2 ground segment interfacing work and for the provision of GOME-2 Level 1 products. We thank the UPAS team for the developmental work on the "Universal Processor for UV/Vis Atmospheric Spectrometers" (UPAS) system at DLR. We are thankful to Henk Eskes (KNMI) for creating the TM5-MP a priori NO 2 profiles. We acknowledge the free use of the GOME-2 surface LER database that is created by KNMI and provided through the AC-SAF of EUMETSAT. We also acknowledge the free use of the online COART radiative transfer model established at the NASA Langley Research Center. We are grateful for the use of imagery from the NASA Worldview application (https://worldview.earthdata.nasa.gov/, last access: 12 February 2020) as part of the NASA Earth Observing System Data and Information System (EOSDIS). We are also grateful to the PIs of the AERONET sites used in this study for maintaining their instruments and providing their data to the community.
Financial support. This research has been supported by the DLR-DAAD Research Fellowships 2015 programme organised by the German Academic Exchange Service (DAAD) and the German Aerospace Center (DLR).
The article processing charges for this open-access publication were covered by a Research Centre of the Helmholtz Association.
Review statement. This paper was edited by Lok Lamsal and reviewed by two anonymous referees.