Evaluating the effects of columnar NO₂ on the accuracy of aerosol optical properties retrievals

We aim to evaluate the NO₂ absorption effect in aerosol columnar properties, namely the aerosol optical depth (AOD), Ångström exponent (AE), and single scattering albedo (SSA), derived from sun–sky radiometers in addition to the possible retrieval algorithm improvements by using more accurate characterization of NO₂ optical depth from co-located or satellite-based real-time measurements. For this purpose, we employ multiannual (2017–2022) records of AOD, AE, and SSA collected by sun photometers at an urban and a suburban site in the Rome area (Italy) in the framework of both the Aerosol Robotic Network (AERONET) and SKYNET networks. The uncertainties introduced in the aerosol retrievals by the NO₂ absorption are investigated using high-frequency observations of total NO₂ derived from co-located Pandora spectroradiometer systems in addition to spaceborne NO₂ products from the Tropospheric Monitoring Instrument (TROPOMI). For both AERONET and SKYNET, the standard network products were found to systematically overestimate AOD and AE. The average AOD bias found for Rome is relatively low for AERONET (∼ 0.002 at 440 nm and ∼ 0.003 at 380 nm) compared to the retrieval uncertainties but quite a bit higher for SKYNET (∼ 0.007). On average, an AE bias of ∼ 0.02 and ∼ 0.05 was estimated for AERONET and SKYNET, respectively. In general, the correction seems to be low for areas with low columnar NO₂ concentrations, but it is still useful for low AODs (< 0.3), where the majority of observations are found, especially under high NO₂ pollution events. For the cases of relatively high NO₂ levels (> 0.7 DU), the mean AOD bias was found within the range 0.009–0.012 for AERONET, depending on wavelength and location, and about 0.018 for SKYNET. The analysis does not reveal any significant impact of the NO₂ correction on the derived aerosol temporal trends for the very limited data sets used in this study. However, the effect is expected to become more evident for trends derived from larger data sets and in the case of an important NO₂ trend. In addition, the comparisons of the NO₂-modified ground-based AOD data with satellite retrievals from the Deep Blue (DB) algorithm of the NASA Moderate Resolution Imaging Spectroradiometer (MODIS) resulted in a slight improvement in the agreement of about 0.003 and 0.006 for AERONET and SKYNET, respectively. Finally, the uncertainty in assumptions on NO₂ seems to have a non-negligible impact on the retrieved values of SSA at 440 nm leading to an average positive bias of about 0.02 (2 %) in both locations for high NO₂ loadings (> 0.7 DU).

1 Introduction

Atmospheric particles have both direct and indirect effects on Earth's radiation budget and climate (IPCC, 2021). Direct radiative forcing arises from the interaction of aerosols with solar radiation through absorption and scattering processes (Hobbs, 1993). As an indirect impact, aerosols play an important role in cloud formation and properties by acting as cloud condensation nuclei on which water vapor condenses and by influencing the cloud albedo and lifetime (Rosenfeld et al., 2014). Moreover, heterogeneous chemical reactions can take place on the surfaces of atmospheric particles, thus having a crucial effect on atmospheric chemistry and composition. Examples of such aerosol-driven reactions are those that lead to stratospheric ozone depletion in the polar regions (Solomon et al., 1986). In addition to their footprint on radiative forcing and climate, aerosols adversely affect human health and have been associated with a wide variety of health issues such as respiratory and neurological diseases, cancer, diabetes, cardiovascular diseases, and hypertension (e.g., Lelieveld et al., 2015; Molina et al., 2020, and references therein).

The above effects of airborne particulate matter on Earth's climate and human health strongly depend on the intra-annual variations in its loading and properties. The most widely used variable for the estimation of columnar aerosol concentration in the atmosphere is the multiwavelength aerosol optical depth (AOD). Aerosol optical properties are monitored globally by satellite, e.g., the Moderate Resolution Imaging Spectroradiometer (MODIS) and ground-based networks of sun photometers like the Aerosol Robotic Network (AERONET; Holben et al., 1998), SKYNET (Nakajima et al., 2020), or the Global Atmosphere Watch Precision Filter Radiometer (GAW-PFR) network (Kazadzis et al., 2018a). Ground-based remote sensing allows accurate AOD retrievals, i.e., of the order of 0.01–0.02, depending on the AOD wavelength (Kazadzis et al., 2018b), which are in fact widely used as a validation reference for satellite- or model-based AOD products (e.g., Chu et al., 2002; Remer et al., 2005; Green et al., 2009; Levy et al., 2010; Li et al., 2015; Sherman et al., 2016; Gkikas et al., 2021; Di Tomaso et al., 2022) and used as input for various modeling initiatives (e.g., Benedetti et al., 2018).

However, AOD retrieval from sun photometers includes some assumptions in order to take into account all the non-aerosol effects in the retrieval spectral range. In particular, AOD retrievals are sensitive to the assumptions on the concentration of atmospheric trace gases absorbing in the instrument spectral bands considered, among which are ozone (O₃) and nitrogen dioxide (NO₂). The exact effect of trace gases in the retrieval at a particular bandwidth depends also on their absorption cross section. For the case of NO₂, as filter radiometers retrieve the AOD in certain wavelength bands based on their filter responsivity, such retrievals, especially in the standard wavelengths of 380 and 440 nm (AERONET), have to be corrected for the NO₂ optical depth. Currently, some AOD retrievals do not take NO₂ optical depth into consideration when deriving AOD (e.g., SKYNET; Nakajima et al., 2020; GAW-PFR; Kazadzis et al., 2018a), while others use satellite-based climatological NO₂ data sets for estimating it (e.g., AERONET; Giles et al., 2019). In the case of the GAW-PFR network, the error introduced in AOD retrievals by NO₂ absorption can be assumed to be negligible due to the low NO₂ concentrations observed in the GAW remote stations (the annual mean values of NO₂ optical depth are in general < 0.001; Kazadzis et al., 2018a). However, especially over polluted areas, NO₂ is characterized by a rather short lifetime and high spatiotemporal variations, due to inhomogeneous local emission patterns and photochemical destruction (e.g., Richter et al., 2005; Boersma et al., 2008; Tzortziou et al., 2014, 2015; Drosoglou et al., 2017; Fan et al., 2021). Although the stratospheric component of NO₂ is quite stable spatially, the tropospheric NO₂ is highly variable in space and time and can bias the calculation of AOD if neglected (Arola and Koskela, 2004; Boersma et al., 2004). Hence, areas with high tropospheric NO₂ emission will tend to have greater proclivity for deviating from climatological mean values, which might not be representative of the actual NO₂ loading and spatial distribution in the atmosphere, introducing potential errors in AOD calculations in those spectral regions with a significant NO₂ absorption footprint.

Satellite observations with improved spatial and temporal resolution, e.g., the Sentinel-5 Precursor TROPOspheric Monitoring Instrument (S5P/TROPOMI), models, or co-location with surface-based Pandora instruments from the Pandonia Global Network (PGN) spectroradiometers (Cede et al., 2020) measuring the total column of NO₂ may assist in reducing the uncertainty in the NO₂ optical depth contribution in later versions of AOD retrieval algorithms. In the present study, we aim at evaluating if and how much AOD, in addition to its spectral variability, i.e., the associated Ångström Exponent (AE), and single scattering albedo (SSA) retrievals could be improved by applying a specific correction using synchronous and co-located measurements of the total NO₂ column from the Pandonia network spectroradiometers. To this end, we exploit the unique configuration of twin observational sites in the Rome area (Italy), where multiannual (2017–2022) records of both multispectral AOD observations and columnar NO₂ measurements are available both in the city center and in a suburban location. High-frequency measurements of total NO₂ performed by co-located Pandora spectroradiometer systems were used to evaluate the current uncertainty in the retrievals of aerosol properties. Aerosol retrieval modifications based on Pandora NO₂ measurements are proposed for both AERONET and SKYNET. In addition, relatively high spatially resolved NO₂ observations from the S5P/TROPOMI satellite sensor were used to demonstrate the possibility of applying the corrections globally. A first attempt to investigate the impact of those corrections on AOD and AE annual trends is also conducted.

2 Instrumentation, data, and methodology

2.1 The target area and relevant observational sites

Rome is the capital and the most populous city of Italy with almost 3 million inhabitants and one of the most densely populated cities in the European Union (ISTAT, 2021). It is located about 24 km east of the Tyrrhenian Sea, surrounded by an extensive undulating plain, and crossed by the Tiber and Aniene rivers. The city is part of the Lazio administrative region in the central part of the Italian Peninsula. The economic activities in the metropolitan area are characterized by the absence of heavy industrial facilities and are related mainly to the services and high-technology sectors, as well as commercial activities and tourism. The city air quality is strongly affected by local emission sources, such as transportation and domestic heating, but it is also markedly affected by local circulation and mid- to long-range transport events of sea salt, wildfires, and Saharan dust (e.g., Ciardini et al., 2012; Gobbi et al., 2013; Barnaba et al., 2017; Valentini et al., 2020; Di Bernardino et al., 2021).

Rome's air quality is monitored on a regular basis by standard in situ instrumentation. These measurements are complemented by multiplatform, long-term observations of aerosol and trace gases performed by a variety of ground-based remote sensing instruments such as sun–sky radiometers, Raman and elastic lidars, automated lidar ceilometers, Pandora, Brewer, and differential optical absorption spectroscopy (DOAS) spectrophotometers (e.g., Di Ianni et al., 2018; Iannarelli et al., 2021; Diémoz et al., 2021). In this study, we used remote sensing measurements of columnar NO₂ and aerosol properties performed in two stations located in the greater area of Rome. More specifically, observations were obtained from an urban station (APL-SAP hereafter) located at the Atmospheric Physics Laboratory of the Physics Department of the Sapienza University of Rome in the city center (41.90^∘ N, 12.52^∘ E; altitude 75 m a.s.l. – above sea level) and a suburban site at the southern east edge of the city in the National Research Council (CNR) Institute of Atmospheric Sciences and Climate (ISAC) Rome Atmospheric Supersite (CIRAS) in Tor Vergata, Rome (41.84^∘ N, 12.65^∘ E; altitude 117 m a.s.l.). These two observational sites, along with the rural station of the CNR Institute of Atmospheric pollution Research (IIA) in Montelibretti, contribute to the Boundary-layer Air Quality-analysis Using Network of Instruments (BAQUNIN) supersite (Iannarelli et al., 2021) and to several national and international observing networks.

2.2 Aerosol data sets

2.2.1 AERONET

The Aerosol Robotic Network (AERONET) is a ground-based passive remote sensing aerosol monitoring network initiated by NASA and expanded by several national and international networks and collaborators (Holben et al., 1998). For more than 2 decades, AERONET has been delivering continuous, long-term data sets of aerosol optical, microphysical, and radiative properties to support aerosol studies and the validation of spaceborne retrievals. The network uses the Cimel CE318-T Sun Sky Lunar multispectral photometers and provides the standardization of instrument calibration and data acquisition, in addition to centralized data processing and distribution. The AERONET public domain database provides retrievals of spectral AOD, inversion products, and precipitable water at a global scale (https://aeronet.gsfc.nasa.gov/, last access: 21 October 2022).

Figure 1 Time series of monthly averaged AOD (a) and AE (b) measurements over APL-SAP (AERONET and SKYNET) and CNR-ISAC (AERONET). Note that AERONET AOD and AE correspond to the wavelength channels of 440 and 440–870 nm, respectively, whereas SKYNET AOD and AE refer to 400 and 400–1020 nm, respectively. The shaded areas correspond to the monthly 1σ standard deviation.

In this study, we employed level 1.5 quality-assured retrievals of AOD at 380, 440, 500, 675, and 870 nm, along with AE at 440–870 nm from the version 3 processing algorithm (Giles et al., 2019; Sinyuk et al., 2020). Level 1.5 data are cloud screened and quality assured, but final calibration has not been applied to them. However, they represent a good tradeoff between quality and readiness, considering that our approach aims to perform a near-real-time improvement on aerosol products. In the standard AERONET AOD retrieval, the NO₂ optical depth is estimated from monthly climatological values of total NO₂ from the Ozone Monitoring Instrument (OMI/Aura) Level-3 retrievals during the 2004–2013 period at 0.25^∘ by 0.25^∘ spatial resolution and the NO₂ absorption coefficients from Burrows et al. (1998). The observations over the CNR-ISAC station used in this work cover the period from March 2017 to mid-August 2022, in which synchronous data from the co-located Pandora instrument are also available. The respective period for APL-SAP is from April 2017 through early September 2022. The aerosol data sets for both locations are presented in Fig. 1. The average AE is 1.23 ± 0.4 and 1.31 ± 0.5 at APL-SAP and CNR-ISAC, respectively, while the average AOD is about 0.18 ± 0.1 at both stations. AOD has a quite marked yearly cycle, with higher AOD values recorded during summer months, i.e., about 0.22 ± 0.1 and 0.21 ± 0.1 at APL-SAP and CNR-ISAC, respectively. AE is also higher during summer, with a mean value of 1.26 ± 0.4 for APL-SAP and 1.38 ± 0.5 for CNR-ISAC.

2.2.2 SKYNET

The SKYNET network, established at the beginning of the 2000s, is a ground-based radiation observation network dedicated to aerosol, cloud, and solar radiation interaction research using the Prede POM sun sky radiometers (Takamura and Nakajima, 2004; Nakajima et al., 2020). It is based on the collaboration and maintenance by several universities and research institutes around the world. This network imposes the standardization of instrument calibration, data acquisition, and data processing and implements two data analysis flows (SR-CEReS and ESR-MRI), mainly based on the SKYRAD.pack, a software package implemented for the POM sky radiometer (e.g., Nakajima et al., 1996; https://www.skynet-isdc.org/methodology.php, last access: 21 October 2022). In contrast to AERONET AOD retrieval methodologies, no correction for NO₂ optical depth is applied in the calculation of SKYNET AOD (e.g., Campanelli et al., 2004; Estellés et al., 2012). Here, we used the ESR-MRI/SUNRAD processor version 0.9 level 2 AOD at 400, 500, 675, 870, and 1020 nm and AE at 400–1020 nm data sets over APL-SAP from late September 2017 to May 2022, which are open-access and available online (https://www.skynet-isdc.org/data.php, last access: 9 June 2023). The SKYNET time series used in our analysis is also illustrated in Fig. 1. The calculated mean AOD and AE are 0.18 ± 0.1 and 1.23 ± 0.4, respectively. These values are similar to the AERONET APL-SAP averages mentioned in Sect. 2.2.1, though they correspond to slightly different wavelengths. SKYNET also reports higher values on average during summer, i.e., 0.22 ± 0.1 and 1.38 ± 0.5 for AOD and AE, respectively.

2.2.3 MODIS Deep Blue data

The Moderate Resolution Imaging Spectroradiometer (MODIS) is a key sensor on board the NASA Terra and Aqua satellites flying, respectively, since 2000 and 2002. Terra MODIS (descending node; about 10:30 UTC) and Aqua MODIS (ascending node; about 13:30 UTC) are observing the entire Earth's surface every 1 to 2 d, acquiring data in 36 spectral bands ranging in wavelength from 0.4 to 14.4 µm, with a spatial resolution of 1 km at nadir (except for a few bands with higher spatial resolution).

Inversion of MODIS observations allows retrievals of several geophysical quantities. Here, we used the aerosol AOD products retrieved using the MODIS Deep Blue (DB) algorithm (Hsu et al., 2004, 2006, 2013). The basic principle of DB algorithms is to utilize the precalculated land surface reflectance database in deep blue bands (0.412 µm), where surface reflectance is relatively lower than those in longer bands. In particular, we used the collection 6.1 DB AOD products for both Aqua and Terra satellites. More details about the DB algorithm are in Hsu et al. (2013) and references therein. The spatial resolution of this product is 10 km. Wei et al. (2019) highlighted that the DB algorithm is relatively more stable and less affected by changes in atmospheric and surface conditions with respect to the Dark Target algorithm (Levy et al., 2013), showing better performances in urban areas for slightly polluted cases, such as the area of Rome. They also highlighted that collection 6.1 AOD products perform better than the previous collections, especially in Europe and North America. The MODIS DB products used in this study are available at the Level-1 and Atmosphere Archive and Distribution System Distributed Active Archive Center (LAADS DAAC; http://ladsweb.nascom.nasa.gov, last access: 21 October 2022).

2.3 Total NO₂ observations

2.3.1 Pandora spectroradiometers

Pandora instruments are compact spectrometers that perform spectral measurements, with a high temporal resolution, of direct solar irradiance and scattered radiance for the retrieval of total and tropospheric column densities of atmospheric trace gases (e.g., NO₂, O₃, and HCHO) that affect air quality, in addition to their near-surface concentrations and vertical profiles (e.g., Herman et al., 2009; Tzortziou et al., 2012, 2015). The total NO₂ vertical column data sets used in the present study were obtained from the Pandora spectrometer no. 115 that has been operating at CNR-ISAC since March 2017 and the Pandora systems (no. 117 and no. 138) that have both been deployed at APL-SAP since April 2016 and within the period August 2019–October 2020, respectively. The above time series have been affected by the COVID-19 lockdown period during February–May 2020 (Campanelli et al., 2021). The monthly averaged values from both stations are presented in Fig. 2 and intercompared in the scatterplot of Fig. 3. On average, the Pandora total NO₂ column over APL-SAP is about 0.07 % higher compared to the CNR-ISAC NO₂.

Figure 2 Time series of monthly NO₂ total column from Pandora instruments over APL-SAP (blue line) and CNR-ISAC (yellow line). The shaded areas correspond to the 1σ standard deviation of the monthly averaged values. The NO₂ concentration is clearly affected by the COVID-19 lockdown that took place during February–May 2020.

Figure 3 Monthly NO₂ total column from Pandora over CNR-ISAC against synchronous APL-SAP observations. The gray shaded area corresponds to the 95 % confidence interval of the linear regression fit (red line).

Pandora total NO₂ column product is derived from the direct-sun measurements in the UV-VIS spectral range 280–530 nm, with an average resolution of 0.6 nm by means of the Blick software and the algorithm implemented therein, as described by Cede (2021). The data sets employed for this work were obtained with the direct-sun retrieval code “nvs3” and the Blick processor version 1.8. Pandora instruments are part of the Pandonia Global Network (PGN; Cede et al., 2020) and have been fully characterized, following the calibration procedures presented by Müller et al. (2020). The recorded raw spectrally resolved radiation measurements are centrally processed for the retrieval of atmospheric trace gas products, which are all publicly available online (https://www.pandonia-global-network.org/, last access: 21 October 2022). In the current study, high- (flags 0 and 10) and medium-quality (flags 1 and 11) data are employed. Information on the quality control of Pandora products can be found in Cede (2021). Pandora NO₂ retrievals have been compared and validated with other ground-based and spaceborne observations during several field campaigns (e.g., Flynn et al., 2014; Martins et al., 2016; Lamsal et al., 2017; Herman et al., 2018; Kreher et al., 2020). Total NO₂ data from the Pandora instrument no. 117 located at APL-SAP have been compared with NO₂ observations retrieved by the co-located MkIV Brewer spectrophotometer (with serial no. 067), revealing a correlation coefficient above 0.96 and a negligible absolute median bias of 0.002 DU (Diémoz et al., 2021). According to Herman et al. (2009), the Pandora direct-sun total NO₂ has a clear-sky precision of 0.01 DU in the slant column and a nominal estimated accuracy of 0.1 DU in the vertical column. In the same study, a systematic difference of less than 1 % was found between the relative slant columns of Pandora and a MultiFunction Differential Optical Absorption Spectroscopy (MFDOAS) instrument.

As already mentioned in Sect. 2.2.1, AERONET uses climatological values from OMI L3 products for the estimation of NO₂ optical depth in AOD retrievals. The corresponding OMI total NO₂ ranges between about 0.2 and 0.3 DU, with an average value of 0.26 ± 0.02 DU. The time series of the Pandora columnar NO₂ differences from the AERONET climatological values for both urban (APL-SAP) and suburban (CNR-ISAC) locations is illustrated in the upper panel of Fig. 4. Pandora NO₂ data are time-interpolated to AERONET measurements. The percentage frequency distributions of the absolute Pandora-OMI deviation for both locations are also presented (Fig. 4; lower panel). About 89 % of the APL-SAP and 87 % of the CNR-ISAC data pairs show an OMI climatology systematic underestimation of NO₂ (positive deviations in Fig. 4). AERONET aerosol retrievals seem to significantly underestimate the NO₂ abundance over urban and suburban locations, with an average absolute difference between the actual Pandora measurements and the estimations from satellite climatology of about 0.15 ± 0.19 DU (61.5 ± 71.5 %) and 0.16 ± 0.18 DU (61.5 ± 67.2 %) for APL-SAP and CNR-ISAC, respectively. This underestimation of the NO₂ levels over urban locations, characterized by strong spatial gradients, can be attributed to the fact that OMI climatology cannot capture the temporal and spatial NO₂ variability within an urban context (e.g., Drosoglou et al., 2017; Herman et al., 2019). Thus, the derived differences in total NO₂ are highly correlated to the Pandora measurements. The majority of PGN-OMI biases lie within 0–0.5 DU, corresponding to Pandora values lower than 1 DU. More specifically, 90 % of the PGN NO₂ data over APL-SAP differ within −0.14 DU (−50 %) and 0.44 DU (150 %) from OMI climatology, while the respective deviation ranges between −0.14 and 0.51 DU (−50 %–170 %) for CNR-ISAC. However, there are quite a few cases (∼ 9.5 % and ∼ 8.8 % for APL-SAP and CNR-ISAC, respectively) of higher PGN values (< 2 DU), leading to larger deviations (up to ∼ 1.6 DU for APL-SAP and ∼ 1.5 DU for CNR-ISAC).

Figure 4 (a) Time series of the Pandora total NO₂ deviation from AERONET NO₂ climatological values (OMI) for both APL-SAP and CNR-ISAC. (b) The corresponding relative frequency distributions of Pandora-OMI deviation for both locations.

2.3.2 TROPOMI

The Tropospheric Monitoring Instrument (TROPOMI) is a nadir-viewing spectrometer on board the Sentinel-5 Precursor (S5P) satellite, which was launched on 13 October 2017. Since August 2019, TROPOMI has a pixel size of 5.5 km × 3.5 km (the initial resolution was 7 km × 3.5 km). NO₂ columns are retrieved using the backscatter solar radiation detected in the spectral window of 405–465 nm (van Geffen et al., 2015) by applying the DOAS technique (Platt, 1994; Platt and Stutz, 2008). The operational TROPOMI NO₂ products are generated using the algorithm described by van Geffen et al. (2022), which is an improvement of the NO₂ DOMINO algorithm (Boersma et al., 2011) developed by the Royal Netherlands Meteorological Institute (KNMI) for the OMI satellite sensor measurements. Both near-real-time (NRTI) and offline (OFFL) NO₂ data sets are retrieved using the KNMI standard algorithm (Eskes et al., 2022; Eskes and Eichmann, 2022). NRTI data files are available within 3 h from the measurement, whereas the OFFL data are processed in offline mode, and the respective files are generated a few days after the sensing time (van Geffen et al., 2022).

In this study, the OFFL NO₂ retrievals are employed, which are the main S5P/TROPOMI product. The extracted NO₂ data set covers the period October 2018–August 2022 and includes observations obtained from several processor versions, beginning with version 01.02.00 before March 2019 and going up to version 02.04.00 after July 2022. The total NO₂ column was calculated from the sum of the tropospheric and stratospheric components, which is preferred over the TROPOMI total NO₂ product for comparisons with ground-based data because the latter suffers from retrieval uncertainties due to its significant dependence on the ratio of the a priori tropospheric and stratospheric columnar data (van Geffen et al., 2022). Additionally, the satellite pixels have been filtered to keep only those with a QA (quality assurance) value > 0.75, corresponding to cloud radiance fraction < 0.5 (Eskes and Eichmann, 2022). The S5P/TROPOMI NO₂ products have been downloaded from the Sentinel-5P Pre-Operations Data Hub of the Copernicus Open Access Hub (https://scihub.copernicus.eu/, last access: 21 October 2022).

For visualization purposes, the averages of the summed NO₂ column re-gridded on a 500 m grid are plotted for the greater Rome area (Fig. 5). The data used in Fig. 5 cover the period from 2018 to 2021, excluding the COVID-19 lockdown period (February–May 2020) in order to prevent the average NO₂ values from being affected by the low values observed during that period.

Figure 5 S5P/TROPOMI summed the total NO₂ column averaged for the period 2018–2021, excluding the COVID-19 lockdown period. The data are gridded on a 500 m grid. The locations of the two observational sites used in this study are also reported for reference.

2.4 AOD and AE corrections for NO₂ absorption

2.4.1 AOD retrievals

The methodology to derive AOD (also referred to as τ) from photometric measurements is based on the Lambert–Beer law (Eq. 1), which describes light attenuation by atmospheric components. I₀(λ) is the intensity of the incident light and I(λ) denotes the radiation intensity after traversing through the atmosphere at a specific wavelength λ.

$$\begin{array}{}\text{(1)}& I\left(\mathit{\lambda }\right)=I_{\mathrm{0}}\left(\mathit{\lambda }\right)\cdot e^{-\left(m_{\mathit{\tau }}\left(\mathit{\lambda }\right)\mathit{\tau }\left(\mathit{\lambda }\right)+m_{\mathrm{R}}\left(\mathit{\lambda }\right){\mathit{\tau }}_{\mathrm{R}}\left(\mathit{\lambda }\right)+{\sum }_j{m}_j\left(\mathit{\lambda }\right){\mathit{\tau }}_j\left(\mathit{\lambda }\right)\right)},\end{array}$$

$$\begin{array}{}\text{(2)}& \begin{array}{rl}{\displaystyle \frac{\ln I\left(\mathit{\lambda }\right)}{\ln I_{\mathrm{0}}\left(\mathit{\lambda }\right)}}=& -(m_{\mathit{\tau }}\left(\mathit{\lambda }\right)\mathit{\tau }\left(\mathit{\lambda }\right)+m_{\mathrm{R}}\left(\mathit{\lambda }\right){\mathit{\tau }}_{\mathrm{R}}\left(\mathit{\lambda }\right)\\ & +\sum _j{m}_j\left(\mathit{\lambda }\right){\mathit{\tau }}_j\left(\mathit{\lambda }\right))\end{array}\end{array}$$

The quantities τ and τ_R describe the optical depth of radiation extinction due to aerosols (Mie scattering) and atmospheric molecules (Rayleigh scattering), whereas m_τ and m_R are the respective air mass factors. ∑_jm_jτ_j represents the sum of the extinction due to absorption from atmospheric gases (Eq. 3), with this depending on the wavelength.

$$\begin{array}{}\text{(3)}& \begin{array}{rl}\sum _j{m}_j\left(\mathit{\lambda }\right){\mathit{\tau }}_j\left(\mathit{\lambda }\right)=& m_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right){\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)+m_{{\mathrm{O}}_{\mathrm{3}}}\left(\mathit{\lambda }\right){\mathit{\tau }}_{{\mathrm{O}}_{\mathrm{3}}}\left(\mathit{\lambda }\right)\\ & +m_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}\left(\mathit{\lambda }\right){\mathit{\tau }}_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}\left(\mathit{\lambda }\right)+\mathrm{\dots }\end{array}\end{array}$$

In our study, we investigate the effects of using an independent, direct measurement of ${\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)$ rather than the climatological value used in the AERONET inversion in determining the AOD (τ). Thus, by combining Eq. (2) with Eq. (3), assuming that the air mass factor in direct-sun measurements is equal to sec(θ) for both aerosol and NO₂, where θ is the solar zenith angle, and absorption from all the other gaseous components stays the same, the difference in AOD due to the different estimation of NO₂ optical depth is obtained by Eq. (4):

$$\begin{array}{}\text{(4)}& \mathrm{\Delta }\mathit{\tau }\left(\mathit{\lambda }\right)={\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}\left(\mathit{\lambda }\right)-{\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{AER}}\left(\mathit{\lambda }\right),\end{array}$$

where ${\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{AER}}$ is the NO₂ absorption optical depth climatology used by AERONET, and ${\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}$ is the optical depth calculated from Pandora NO₂ measurements. The latter is derived using Eq. (5):

$$\begin{array}{}\text{(5)}& {\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}\left(\mathit{\lambda }\right)={\mathit{\sigma }}_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)\cdot c_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}.\end{array}$$

The quantity ${\mathit{\sigma }}_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)$ in Eq. (5) refers to the absorption cross section of NO₂ at wavelength λ (Burrows et al., 1998), and $c_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}$ is the total NO₂ column from Pandora instrument. The modified AOD values (τ_{AER_mod}) are obtained from the standard AERONET AOD (τ_AER) by applying the following equation:

$$\begin{array}{}\text{(6)}& \begin{array}{rl}{\mathit{\tau }}_{\mathrm{AER}\mathrm{\_}\mathrm{mod}}\left(\mathit{\lambda }\right)=& {\mathit{\tau }}_{\mathrm{AER}}\left(\mathit{\lambda }\right)-(\left({\mathit{\sigma }}_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)\cdot c_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}\right)\\ & -{\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{AER}}\left(\mathit{\lambda }\right)).\end{array}\end{array}$$

The same approach was also applied to the SKYNET AOD data. However, since the SKYNET retrievals assume ${\mathit{\tau }}_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{SKYNET}}=\mathrm{0}$, Eqs. (4) and (6) are modified as follows:

$$\begin{array}{}\text{(7)}& {\displaystyle }\mathrm{\Delta }\mathit{\tau }\left(\mathit{\lambda }\right)={\mathit{\sigma }}_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)\cdot c_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}},\text{(8)}& {\displaystyle }{\mathit{\tau }}_{\mathrm{SKYNET}\mathrm{\_}\mathrm{mod}}\left(\mathit{\lambda }\right)={\mathit{\tau }}_{\mathrm{SKYNET}}\left(\mathit{\lambda }\right)-\left({\mathit{\sigma }}_{{\mathrm{NO}}_{\mathrm{2}}}\left(\mathit{\lambda }\right)\cdot c_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}\right),\end{array}$$

where τ_SKYNET(λ) denotes the standard SKYNET AOD at spectral channel λ, and τ_{SKYNET_mod}(λ) is the modified AOD at wavelength λ.

2.4.2 AE retrievals

The spectral variability in AOD is generally expressed as follows:

$$\begin{array}{}\text{(9)}& {\displaystyle }\mathit{\tau }=\mathit{\beta }\cdot {\mathit{\lambda }}^{-\mathit{\alpha }},\text{(10)}& {\displaystyle }\ln \mathit{\tau }=\ln \mathit{\beta }-\mathit{\alpha }\cdot \ln \mathit{\lambda },\end{array}$$

where α stands for the Ångström exponent (AE).

The AERONET AE product (Eck et al., 1999) is calculated by applying a least squares regression fit on Eq. (10), using the AOD and wavelength logarithms for each non-polarized wavelength channels in different spectral ranges (i.e., 340–440, 380–500, 440–675, 440–870, and 500–870 nm). The negative slope of this linear fit is the Ångström exponent α (Eq. 11).

$$\begin{array}{}\text{(11)}& \mathit{\alpha }=-{\displaystyle \frac{N\sum \ln {\mathit{\tau }}_i\ln {\mathit{\lambda }}_i-\sum \ln {\mathit{\lambda }}_i\sum \ln {\mathit{\tau }}_i}{N\sum {\left(\ln {\mathit{\lambda }}_i\right)}^{\mathrm{2}}-{\left(\sum \ln {\mathit{\lambda }}_i\right)}^{\mathrm{2}}}}\end{array}$$

Here, we also investigate the impact of using synchronous Pandora total NO₂ data in an AOD algorithm (as described in Sect. 2.4.1) on AE retrievals. To do this, the AERONET AE product in the range 440–870 nm was used along with the AOD of non-polarized channels included in this range, i.e., 440, 500, 675 and 870 nm. AE was recalculated based on Eq. (11), using the modified AOD at wavelengths 440 and 500 nm obtained from Eq. (6). For the other channels (675 and 870 nm) in which NO₂ absorption is negligible, the standard AOD data from AERONET were employed.

For SKYNET, AE is calculated by applying a least squares regression fit on Eq. (10), using the AOD and wavelength logarithms at all wavelengths (400, 500, 675, 870, and 1020 nm). Again, AOD was recalculated using Eq. (8) only at wavelengths 400 and 500 nm, where the impact of the NO₂ absorption is significant.

The difference in AE due to the different estimation of NO₂ optical depth in AOD retrievals is expressed as follows:

$$\begin{array}{}\text{(12)}& \mathrm{\Delta }\mathit{\alpha }\left(\mathit{\lambda }\right)=\mathit{\alpha }\left(\mathit{\lambda }\right)-{\mathit{\alpha }}_{\mathrm{mod}}\left(\mathit{\lambda }\right),\end{array}$$

where α_mod(λ) represents the modified AE data, and α(λ) denotes the AE standard product from the AERONET or SKYNET network.

2.5 Trend calculations

In this study, we also evaluate the impact of modified AOD and AE retrievals, as described in Sect. 2.4.1 and 2.4.2, on aerosol temporal trends. This is only a first attempt to investigate the possible effect of NO₂ absorption on the AOD and AE trends, since the data sets used here are quite short for statistically meaningful calculations.

The annual trends in AOD and AE were estimated by applying the weighted least squares fitting technique introduced by Weatherhead et al. (1998) and previously adopted in several aerosol trend analysis studies from space and the ground (e.g., Zhang and Reid, 2010; Yoon et al., 2012; Logothetis et al., 2021). The applied linear trend model is based on the following formula:

$$\begin{array}{}\text{(13)}& Y_t=\mathit{\mu }+\mathit{\omega }{X}_t+{\mathit{\epsilon }}_t,t=\mathrm{1},\mathrm{\dots },T,\end{array}$$

where Y_t is the monthly average aerosol property of interest, μ is a constant term representing the linear fit offset at the start of the time series, ω stands for the magnitude of the trend per year, and ε_t is the monthly average noise not represented by the linear fit. $X_t=t/\mathrm{12}$ is the decimal number of years since the first month of the time series, t is the month index, T denotes the total number of months, and $T/\mathrm{12}$ is the total number of years in the time series.

In order to account for data variability due to severe aerosol events and cloud disturbance, we introduced a monthly weighting factor w_t into the linear fitting procedure (Eq. 14; Yoon et al., 2012). This weighting factor is defined as the square root of the number of observations available each month n_t divided by the monthly standard deviation σ_t (Eq. 15).

$$\begin{array}{}\text{(14)}& {\displaystyle }{\mathit{\chi }}^{\mathrm{2}}\left(\mathit{\mu },\mathit{\omega }\right)=\sum _{t=\mathrm{1}}^T{\left(w_t\cdot \left(Y_t-\mathit{\mu }-\mathit{\omega }{X}_t\right)\right)}^{\mathrm{2}},\text{(15)}& {\displaystyle }{w}_t={\displaystyle \frac{\sqrt{n_t}}{{\mathit{\sigma }}_t}}\end{array}$$

In order to derive statistically significant monthly mean values, a minimum number of 10 observations in a daily basis was ensured. In addition, qualified monthly averages require the availability of measurements from at least 10 d per month. Data were filtered based on the above criteria, and days and/or months that did not fulfill them were excluded from the data sample used in the trend calculations. It should be noted that the data sets employed in this study are quite short for statistically meaningful aerosol trend analysis. However, this is a first attempt to investigate the impact of modified AOD and AE calculations on the derived temporal trends.

2.6 GRASP algorithm

The Generalized Retrieval of Atmosphere and Surface Properties (GRASP; Dubovik et al., 2021) is a state-of-the-art inversion algorithm based on a statistically optimized multiterm least squares method (LSM) proposed by Dubovik (2004). GRASP has been applied to numerous applications covering a vast variety of instruments and, interestingly, to very different combinations between them. Among the different applications of GRASP, it is possible to find GRASP/POLDER-3 (Chen et al., 2020), GRASP/AOD (Torres et al., 2017), OLCI/GRASP (Chen et al., 2022), the combination of active lidar measurements and ground-based radiometry (Lopatin et al., 2013, 2021; Román et al., 2018; Herreras et al., 2019), the retrieval of all-sky cameras (Román et al., 2017, 2022), or, for example, applications to in situ measurements including polar nephelometers (Espinosa et al., 2017, 2019; Schuster et al., 2019).

The GRASP scientific core was borne from the heritage of the AERONET inversion algorithm (Dubovik and King, 2000; Dubovik et al., 2000; Dubovik, 2004; King and Dubovik, 2013). At the same time, as discussed above and by Dubovik et al. (2011, 2021), the possibilities of GRASP have been extended due to the totally generalized nature of the inversion module and the continuous developments of the forward model.

For this study, GRASP has been used to mimic AERONET standard retrieval in order to understand the effects of the NO₂ concentration on the retrieved SSA at 440 nm. In this case, two different approaches were followed for the GRASP algorithm. First of all, GRASP has been used as close as possible to the standard AERONET retrieval, which means that the input measurements of the algorithm are the total optical depth (TOD) and the almucantar sky measurement routine at 440, 675, 870, and 1020 nm. In the first approach (GRASP/AERONET NO₂ hereafter), the NO₂ absorption is taken into account, using OMI climatology, exactly as done for AERONET. On the other hand, GRASP flexibility allows the use of different assumptions of the gaseous properties. Therefore, in addition to the standard approach, the aerosol retrieval has also been done using the total columnar NO₂ concentrations provided by the Pandora spectrometers co-located with AERONET instruments at the two stations selected for this study. This methodology will hereafter be referred to as GRASP/Pandora NO₂. Thus, in addition to the standard AERONET retrieval products, GRASP has provided aerosol retrieval using these more accurate NO₂ concentrations. The NO₂ absorption features were calculated more precisely from those concentrations by using a k-distribution approach or the “kbin” code (Doppler et al., 2014a, b) to speed up the calculations.

3 Results and discussion

3.1 Differences in AOD and AE retrievals using Pandora NO₂ data

The differences in AOD (Δτ) at 440 nm and, thus, of its spectral variability through the AE (Δα at 440–870 nm) correcting for measured NO₂ effects with respect to the standard AERONET retrievals are illustrated in Fig. 6 for both the Rome CNR-ISAC and APL-SAP stations. The frequency distributions of AOD, Δτ, and Δα are also included in Fig. 6. Δτ is defined as the standard minus the modified AOD (τ_AER–τ_{AER_mod}; see Eqs. 4–6). Similarly, Δα is defined as α_AER–α_{AER_mod} (Eq. 12). The derived values are presented versus the AOD at 440 nm and are color coded with respect to the Pandora NO₂ retrievals. The dependency of Δτ on NO₂ is quite clear. As expected, higher Δτ absolute values are obtained for higher NO₂ concentrations, regardless of the initial measured AOD. Also, the absolute percentage of Δτ with respect to the AOD is higher for lower aerosol loadings, which means that the impact of the NO₂ correction is more significant on lower AODs. This fact is also clear from Δα, which is higher not only for higher NO₂ but also for lower AOD values as well. Interestingly, based on Fig. 6, the highest Pandora NO₂ retrievals (reddish colors) are not associated with the highest AOD values, indicating that in Rome the high AOD loadings are not strictly associated with high NO₂ pollution events. In fact, high AODs are frequently related to the long-range transport of elevated layers of desert dust, fire plumes, or a combination of both (e.g., Barnaba et al., 2011; Gobbi et al., 2019; Campanelli et al., 2021; Andrés Hernandez et al., 2022). Hence, it might be worth modifying the aerosol retrievals for high NO₂ in those pollution-related events with low to medium AOD levels. More about AOD and aerosol type climatology for the Rome area can be found in Di Ianni et al. (2018) and in Campanelli et al. (2022).

Figure 6 The differences in the modified AERONET AOD at 440 nm (a, b) and AE at 440–870 nm (c, d) over CNR-ISAC and APL-SAP from the standard products illustrated with respect to the standard AERONET AOD measurements at 440 nm and the actual NO₂ observed by Pandora (color scale). The corresponding distributions of all variables are also included.

In general, considering the climatological value chosen for Rome in AERONET retrievals, the use of actual, coincident NO₂ measurements on the calculations of aerosol properties still seems to be useful for AOD < 0.3, while being quite low (less than 10 %) for AOD > 0.5 and almost negligible for AOD > 0.8. In most cases, AERONET retrievals seem to overestimate AOD and AE. However, there are cases of underestimation, especially in AE retrievals, which seems to be higher for lower AODs. Those underestimations correspond to overestimation of NO₂ from satellite monthly climatological values used in AERONET retrievals. The estimated AOD and AE deviations are below 0.01 and 0.1, respectively, for the majority of observations, i.e., about 96 %–98 % of occurrences for both CNR-ISAC and APL-SAP (see also distributions in Fig. 6). The average AOD bias is between 0.002 ± 0.003 and 0.003 ± 0.003 (with the higher values observed at 380 nm), while the average AE bias is ∼ 0.02 ± 0.03. Overall, the mean AOD bias is low compared to the estimated uncertainties for the standard AERONET product, i.e., 0.01–0.02 (with the higher errors observed in the UV; Sinyuk et al., 2020). However, the mean AOD bias for the cases of high NO₂ levels (> ∼ 0.7 DU) is ∼ 0.011 ± 0.003 at 440 nm and ∼ 0.012 ± 0.003 at 380 nm for APL-SAP and ∼ 0.009 ± 0.003 at 440 nm and ∼ 0.010 ± 0.003 at 380 nm for CNR-ISAC, which is comparable to the AERONET reported uncertainties. The estimated mean bias of AE retrievals for the cases with high NO₂ (> ∼ 0.7 DU) is ∼ 0.08 ± 0.04 for both Rome sites. The threshold for NO₂ has been selected as being the average Pandora NO₂ (∼ 0.4) calculated from the whole data set plus 2 times the standard deviation.

The results for SKYNET observations are similar (Fig. 7), but only positive Δτ and Δα values are derived, indicating the overestimation of the aerosol properties, since the NO₂ optical depth is not considered in the standard retrieval processes (see Eqs. 7–8). Δτ is defined as τ_SKYNET–τ_{SKYNET_mod} (see Eqs. 7–8), and Δα stands for α_SKYNET–α_{SKYNET_mod} (Eq. 12). In addition, the derived deviations in aerosol properties reach higher values compared to AERONET. Especially AE differences extend up to a value of about 0.7, which is more than double compared to AERONET results. Interestingly, these quite large Δα values (> 0.3) correspond to relatively low NO₂ loadings (< 1.2 DU). The differences observed between the two networks can be partly attributed to the different wavelength channels used for AOD and AE retrievals. Similar to AERONET, the derived AOD and AE biases for SKYNET are below 0.01 and 0.1, respectively, for the majority of observations (i.e., about 85 % of occurrences for AOD and about 90 % for AE; see also distributions in Fig. 7). The overall average AOD bias is ∼ 0.007 ± 0.003, which can be assumed to be low, considering that Nakajima et al. (2020) have estimated a root mean square difference (RMSD) of about 0.03 for wavelengths < 500 nm in city areas in AOD comparisons with other networks. However, the mean AOD bias for the cases with high NO₂ levels (> ∼ 0.7 DU) is found to be about 0.018 ± 0.003, which is comparable to the RMSD value reported by Nakajima et al. (2020). The overall average AE bias calculated in this study is ∼ 0.05 ± 0.04, whereas the AE bias averaged over the high NO₂ cases is about 0.10 ± 0.05.

Figure 7 The differences in the modified SKYNET AOD at 400 nm (a) and AE at 400–1020 nm (b) over APL-SAP from the standard products illustrated with respect to the standard SKYNET AOD measurements at 400 nm and the actual NO₂ observed by Pandora (color scale). The corresponding distributions of all variables are also included. Note that the spectral channels for the retrievals and the axis scales are different compared to AERONET.

The World Meteorological Organization (WMO, 2005) states that, when comparing AOD retrieved from sun photometers, 95 % of the AOD differences should lie within $\pm (\mathrm{0.005}+\mathrm{0.01}/m)$ of AOD, where m is the optical air mass. The first term of the equation (0.005) represents the maximum tolerance for the uncertainty due to the atmospheric parameters used for the AOD calculation (additional atmospheric trace gas corrections, i.e., ozone and NO₂ and Rayleigh scattering), while the second term ($\mathrm{0.01}/m$) describes the calibration-related relative uncertainties, for which the WMO recommends an upper limit of 1 % (e.g., Cuevas et al., 2019; Kazadzis et al., 2018a). Based on the above, although the average deviations found in this study are low compared to the retrieval uncertainties, they cannot be considered negligible, especially the average systematic underestimation of AOD of about 0.007 from SKYNET, also bearing in mind that there are locations with much higher average NO₂ compared to the city of Rome.

Table 1 Deviation of Pandora total NO₂ column from satellite climatology used for AERONET retrievals and differences in the modified AERONET and SKYNET AOD and AE from the standard products over CNR-ISAC and APL-SAP calculated using actual Pandora total NO₂ observations, as well as daily and monthly averaged values of NO₂. Note that the spectral channels used in AERONET retrievals are 380 and 440 nm for AOD and 440–870 nm for AE, whereas for SKYNET the wavelength channels are 400 and 400–1020 nm for AOD and AE, respectively.

The statistics showing mean differences in AOD and AE AERONET and SKYNET retrievals using actual, coincident NO₂ measurements are presented in Table 1. AERONET AOD retrievals at 380 nm are also included in the table. In addition, deviations of AOD and AE using daily or monthly averages of NO₂ in AERONET and SKYNET observations are also investigated. Table 1 shows that the average deviations of AOD and AE values do not change significantly, regardless of whether the actual Pandora NO₂ measurements or the daily or monthly mean values are used for the retrievals. The percentage differences for AOD lie within the range 1.2 %–1.9 % for AERONET, while they are more than doubled (5.3 %–5.7 %) for SKYNET. For the standard aerosol products of the latter, NO₂ optical depth is not considered. The estimated percentage differences for AE are within 1.2 %–1.7 % and 2.6 %–2.9 % for AERONET CNR-ISAC and APL-SAP, respectively, and between 7 %–7.9 % for SKYNET APL-SAP. It should be noted that the spectral channels used in AERONET retrievals are 380 and 440 nm for AOD and 440–870 nm for AE, whereas SKYNET data refer to 400 and 400–1020 nm for AOD and AE, respectively.

3.2 AOD and AE retrievals based on TROPOMI NO₂ data

Satellite sensors perform measurements globally and provide information on the air quality, even over regions that lack ground-based observations. However, as already mentioned for OMI in Sect. 2.3.1, the spatial resolution of the satellite retrievals is limited by the pixel size. Co-located S5P/TROPOMI observations, characterized by an improved spatial and temporal resolution compared to previous satellite missions (e.g., OMI), were also employed to investigate whether the ground-based retrievals of aerosol properties could be improved on a global scale. Again, the approach described in Sect. 2.4.1 and 2.4.2 was applied by replacing the Pandora total NO₂ ($c_{{\mathrm{NO}}_{\mathrm{2}}\mathrm{PGN}}$) with the corresponding columnar retrievals from TROPOMI. Based on the current satellite footprint (5.5 km × 3.5 km), a radius of 5 km around each ground-based station was selected for the spatial co-location. The TROPOMI NO₂ data were time-interpolated to AERONET and SKYNET measurements. Despite the improved spatial resolution of TROPOMI, the NO₂ corrections using TROPOMI data are expected to be less accurate than those performed with the Pandora product. For example, Lambert et al. (2021) showed a bias between TROPOMI and Pandora total NO₂ column ranging from −23 % over polluted stations to +4.1 % over clean areas, with a median bias of −7.1 %, in the frame of the standard validation process of TROPOMI Level 2 NO₂ products. Other studies have concluded similar results. For example, Zhao et al. (2020) showed a negative bias for the standard TROPOMI total NO₂ product in the range 23 %–28 % over urban and suburban environments and a positive bias of 8 %–11 % at a rural site, while Park et al. (2022) showed 26 %–29 % negative bias and R² within 0.73–0.76 over the Seoul metropolitan area in South Korea.

Table 2 Similar to Table 1, here we are using TROPOMI measurements instead of Pandora total NO₂ for the estimation of NO₂ abundance in AERONET and SKYNET aerosol retrievals.

The statistical metrics of the averaged deviations of the modified AERONET and SKYNET AOD and AE retrievals using actual, co-located TROPOMI NO₂ measurements from the network standard products are presented in Table 2. Similar to Sect. 3.1 and Table 1, the deviations of AOD and AE retrievals derived by employing daily or monthly mean TROPOMI total NO₂ were also investigated. The average deviations of AOD and AE values do not change significantly, regardless of whether the actual TROPOMI NO₂ measurements or the daily mean values are used for the retrievals. This behavior is expected, considering that TROPOMI overpasses occur once or twice per day, and hence, they do not capture daily variations in NO₂. In the case of the monthly averaged TROPOMI NO₂ data, the estimated differences between the standard and modified aerosol products drop notably for AERONET. However, there are still differences compared to OMI NO₂ climatology due to the improved spatial resolution of the TROPOMI pixel. The average AOD bias is ∼ 0.001 ± 0.001 (with the higher values observed at 380 nm), while the average AE bias is ∼ 0.01 ± 0.01 for both AERONET stations. For the cases of high NO₂ levels (> ∼ 0.7 DU), the mean AOD bias is ∼ 0.004 ± 0.001 at 440 nm and ∼ 0.005 ± 0.002 at 380 nm for APL-SAP and ∼ 0.003 ± 0.001 at both 440 and 380 nm for CNR-ISAC. The estimated mean bias of AE retrievals for the cases with high NO₂ (> ∼ 0.7 DU) is ∼ 0.05 ± 0.04 and ∼ 0.02 ± 0.01 for APL-SAP and CNR-ISAC, respectively. In the case of SKYNET, the overall average AOD bias is ∼ 0.005 ± 0.002 for AOD and ∼ 0.04 ± 0.03 for AE. For the high NO₂ cases, a mean AOD bias of about 0.011 ± 0.002 and an average AE bias of ∼ 0.07 ± 0.04 were calculated. Interestingly, the deviations of SKYNET retrievals using monthly TROPOMI data are very similar to those derived using the actual overpasses or daily averaged TROPOMI NO₂, probably due to the fact that the NO₂ optical depth is not included in the standard network AOD retrieval processes.

The percentage differences for AOD lie within the range 0.2 %–0.9 % for AERONET and are about 3.8 %–3.9 % for SKYNET, which are much lower compared to those derived using Pandora NO₂ (see Table 1). The estimated percentage differences for AE are ∼ 0.8 %–0.9 % and ∼ 1.6 %–1.7 % for AERONET CNR-ISAC and APL-SAP, respectively, and about 4 % for SKYNET APL-SAP using actual or daily TROPOMI data. It should be noted again that the spectral channels used in AERONET retrievals are 380 and 440 nm for AOD and 440–870 nm for AE, whereas SKYNET data refer to 400 and 400–1020 nm for AOD and AE, respectively.

3.3 Case study: impact of high Pandora NO₂ on low AOD

In order to investigate further the impact of high NO₂ during pollution events on the retrieval of relatively low levels of AOD, we used measurements performed at APL-SAP on 25 June 2020, the morning of which there was a high NO₂ event. In the upper panels of Fig. 8, the total NO₂ measured from Pandora during that day is illustrated. For AERONET (left panels of Fig. 8), the satellite climatological values used in the retrieval of standard AOD product and their deviations from Pandora NO₂ are also displayed. The standard and NO₂-modified AOD and AE data from both AERONET and SKYNET (see also Sect. 2.4.1 and 2.4.2), in addition to the magnitude of the respective differences (Δτ and Δα), are presented in the middle and lower panels of Fig. 8.

Figure 8 Case study over APL-SAP on 25 June 2020 for both AERONET and SKYNET. (a, b) Pandora total NO₂ column and its deviation from climatology. (c, d) AOD (solid blue line), its improvement using Pandora NO₂ (dashed blue line), and the magnitude of improvement (light orange line and right y axis). (e, f) Similar to panels (c) and (d) but for AE retrievals. Note that the spectral channels for the retrievals are different for the two networks.

The differences in the AOD and AE retrievals from both networks are significant only within a time span of about 3 h around the high NO₂ event (∼ 07:00–10:00 UT) and can be assumed to be negligible for the rest of the day when the NO₂ levels remain quite low. The median AOD bias for AERONET is about 0.003, with a maximum of about 0.02 at the peak of the event. The median and maximum AE biases are 0.014 and 0.11, respectively. It can be also noted that, in the case of SKYNET, both AOD (median value of ∼ 0.008 with a maximum of ∼ 0.03) and AE deviations (median and maximum values of ∼ 0.03 and 0.10, respectively) are a bit higher compared to the respective AERONET deviations of synchronous data. This can be mainly attributed to the fact that SKYNET standard AOD retrieval processes do not account for the NO₂ absorption and can be partly explained by the different channels used in the detectors of the two networks.

3.4 Impact on AOD and AE trends

In this section, a first attempt is made to investigate the effect of the modified AOD and AE retrievals based on the Pandora total NO₂ observations on the annual trends of those aerosol properties. The annual trends of AERONET/SKYNET AOD and AE over both APL-SAP and CNR-ISAC sites, calculated by applying the approach described in Sect. 2.5, and their uncertainties (standard errors in the regression slope) are presented in Table 3.

Table 3 AOD and AE trends and their uncertainties for both standard and modified AERONET and SKYNET products over CNR-ISAC and APL-SAP. Note that the spectral channels used in AERONET retrievals are 440 nm for AOD and 440–870 nm for AE, whereas those for SKYNET are 400 and 400–1020 nm for AOD and AE, respectively. The trend uncertainties refer to the standard error in the regression slope. The differences are calculated on the absolute trend values.

It should be noted here that the aerosol data sets from the two networks correspond to slightly different time periods. In addition, there are significant gaps in the time series from CNR-ISAC due to instrument problems, and the COVID-19 lockdown period (February–May 2020) has been excluded from the data analysis. Therefore, the results in Table 3 are mainly intended to highlight how a different NO₂ correction may affect the aerosol trends and should be interpreted separately for each individual site. Interpretation of the trend significance for the Rome area is not possible using only this short period of time (∼ 5.5 years), considering that the estimated trends are quite small and the uncertainties introduced by linear regression are relatively high.

One aspect shown here is that the difference in the AOD and AE trends for the two data sets (original and modified NO₂) is comparable with the calculated trends. As expected, AE trends with and without NO₂ correction show relatively higher differences, as AE is much more sensitive to spectral AOD changes. However, the linear fitting uncertainty of AE is also high. NO₂ effects on AOD trends would be more obvious in the case of a significant NO₂ trend during a certain period. A thorough long-term trend analysis is out of the scope of this work and could be the topic for a future study.

3.5 Impact on the intercomparison of ground-based and satellite AOD data

In this section, we have analyzed a potential effect of considered NO₂ corrections on the agreement of AERONET and SKYNET AOD products with relevant satellite data. Indeed, it is well known that most satellite retrievals are validated against ground-based measurements of AOD that are considered to be a ground truth. Moreover, most satellite retrieval algorithms are substantially tuned to closely match AERONET observations. For example, all MODIS algorithms, including DB, rely, in one way or another, on AERONET dynamic aerosol models and climatologies of AERONET retrievals. Nonetheless, since MODIS retrievals fundamentally rely on MODIS radiances that are fully independent of AERONET data, some inaccuracies in the assumptions, such as those regarding the NO₂ amount, can cause some additional biases between AERONET and MODIS AOD results.

To evaluate the effects of the proposed correction, we have compared AERONET and SKYNET AOD products against MODIS DB AOD products at 470 nm for the 2017–2022 period. In the intercomparison, we considered only MODIS DB AOD products for which the distance between the center of the pixel and the AERONET site location (APL-SAP or CNR-ISAC) does not exceed 5 km. Furthermore, we considered all the AERONET (or SKYNET) AOD data within ±30 min from the MODIS satellite overpasses. In order to guarantee the quality of the data, we used MODIS DB AOD with a QA index ≥ 2, which corresponds to good and very good products (Wei et al., 2019).

The intercomparison has been performed using MODIS DB AOD at 470 nm. Consequently, we computed the AERONET and SKYNET AOD at 470 nm, thus exploiting the AE. The AERONET AOD at 470 nm was calculated using the standard AERONET AOD at 440 nm and AE at 440–870 nm. Similarly, the SKYNET AOD at 470 nm was computed using the standard SKYNET AOD at 400 nm and AE at 400–1020 nm. The NO₂-modified AERONET and SKYNET AOD values at 470 nm were also computed with the same approach, and the AOD and AE retrievals have been modified using the Pandora NO₂ data.

Figure 9 Intercomparison of MODIS DB with standard (a–c) and modified (d–f) ground-based AOD at 470 nm for both CNR-ISAC (a, d) and APL-SAP (b, c, e, f) sites against AERONET AOD (a, b, d, e) and SKYNET AOD (c, f). The y = x lines and MODIS DB EE envelopes ±(0.05 + 20 %) are plotted as dashed lines. The intercomparison was performed considering a maximum distance between the center of the MODIS DB pixel and the site location of 5 km and Δt_max​​​​​​​ (time between MODIS and AERONET/SKYNET observations) of ±30 min.

We observe a generally satisfactory agreement between the ground-based (both AERONET and SKYNET) and MODIS DB AOD data, with a Pearson correlation (r) higher than 0.7. In general, MODIS DB AOD slightly overestimates the AOD observed by the sun photometers. The bias (calculated as satellite minus sun photometer AOD) between MODIS DB and the different ground-based data sets before the correction (upper panels of Fig. 9) varies from −0.009 for SKYNET APL-SAP data (−0.008, considering AERONET) to 0.027 for AERONET CNR-ISAC. AERONET data, available for both sites, highlight a lower agreement for the CNR-ISAC site, with a bias about 3 times larger with respect to the APL-SAP site. The correction introduces a slight change of about 0.003 in the agreement between MODIS DB and AERONET AOD products and of 0.006 between MODIS and SKYNET data (lower panels of Fig. 9). Figure 9 also shows an improvement in the percentage of MODIS AOD data falling within the expected error (EE) of ±(0.05 + 20 %; Hsu et al., 2013) for APL-SAP by adopting the correction for both AERONET and SKYNET.

Figure 10 (a–c) Absolute correction as a function of the corresponding MODIS DB AOD data and PGN NO₂ data (color scale) for both CNR-ISAC and APL-SAP sites using AERONET and SKYNET AOD. The analysis was performed considering a maximum distance between the center of the MODIS DB pixel and the site location of 5 km and Δt_max of ±30 min. (d, e) Absolute correction of the MODIS DB AOD data for both CNR-ISAC and APL-SAP, using AERONET AOD as a function of the corresponding MODIS DB AOD data and the absolute difference between PGN and OMI climatological NO₂ (color scale).

In Fig. 10, we show the absolute correction (computed as the difference between original AERONET/SKYNET AOD data at 470 nm and modified ones) as a function of the MODIS DB AOD and the NO₂ column retrieved by the Pandora instruments located at APL-SAP and CNR-ISAC sites (upper panels). As already highlighted, we observe that the correction only depends on the NO₂ amount and not on the AOD. Figure 10 also highlights that, although the improvement is relatively low on average, the correction can be larger than 10 %/15 % in many cases.

This intercomparison exercise demonstrated that the proposed correction slightly improves the agreement between MODIS DB AOD data and AERONET and SKYNET AOD products, even if, on average, it is not statistically significant. Nevertheless, as shown in Fig. 10, the improvement becomes significant when the differences between the NO₂ values observed by Pandora and the OMI NO₂ climatology are also significant (lower panels in Fig. 10). Furthermore, since the proposed correction depends on the amount of NO₂, the improvement is more evident in the correspondence of high values of NO₂ (upper panels in Fig. 10), typical of highly polluted areas such as the urban area of Rome (APL-SAP). Also, a slight improvement is also achieved in the suburban area of Rome (CNR-ISAC). Finally, in the case of SKYNET AOD products, the systematic overestimation, due to neglected NO₂ extinction in the official retrieval chain, is eliminated.

3.6 Impact on SSA

One of the main impacts of accurate characterization of the columnar NO₂ concentration is certainly expected on the retrieved values of SSA in spectral ranges coinciding with NO₂ absorption. In order to quantify this effect, the sensitivity of the AERONET retrieval of SSA at 440 nm has been tested. As previously explained (Sect. 2.6), two different GRASP approaches have been applied to this purpose, namely the GRASP/AERONET NO₂ and the GRASP/Pandora NO₂. Despite the close methodological basis between GRASP and AERONET retrievals, the divergence in the development of both algorithms has led to some differences in the retrieved products. Thus, in order to ensure that the difference in the retrieved SSA at 440 nm is produced exclusively by the changes in the description of NO₂ absorption and to avoid the inclusion of any other sources of discrepancy, the GRASP code has been used in both approaches instead of the standard AERONET SSA product.

Figure 11 Comparisons of SSA at 440 nm obtained with GRASP, following the standard AERONET procedure (x axis), and a similar approach but precisely accounting for NO₂ concentration (y axis) from the co-located Pandora instruments in two different stations, namely APL-SAP (a) from March 2017 to November 2020 and CNR-ISAC (b) from April 2017 to September 2021. The data have been filtered to show retrievals corresponding to NO₂ concentration higher than 0.7 DU. The color of the circles is an indicator of the density of points; i.e., colors closer to red indicate a higher number of points close together. The absolute mean bias error (MBE; percent in parentheses), the root mean square error (RMSE) and the correlation coefficient of the linear fit are also shown in the figure. The probability density functions of the difference between both methodologies (GRASP/Pandora NO₂–GRASP/AERONET NO₂) can be found in the lower panels, correspondingly, for each station. The probability density functions for SSA values higher or lower than 0.9 are also included.

The comparisons of the SSA at 440 nm obtained with both methodologies for the two stations for the complete data set (not shown) do not show a clear influence of the change in the NO₂ concentration. High correlations (R > 0.98) and a mean bias error (MBE < 0.002) very close to zero are obtained. The mean NO₂ column concentration for the retrievals presented here is 0.4 DU. Thus, in general, the analyzed improvements are not expected to produce an important change in the retrieved parameters at 440 nm in conditions with relatively low NO₂ absorption. However, in the cases where NO₂ concentration is elevated compared to the climatologically expected range, significant changes in the SSA at 440 nm retrievals can be appreciated. Figure 11 shows the comparisons of the SSA at 440 nm obtained with GRASP, following an AERONET-like approach (x axis) and the approach with the new NO₂ concentrations provided by Pandora (y axis) and filtered for NO₂ concentrations higher than 0.7 DU, which corresponds to the average NO₂ plus 2 times the standard deviation. The two stations are correspondingly represented in the left and right panels. As it can be noted, for both stations in conditions of high NO₂ concentrations, there is a consistent positive bias of ∼ 0.02 (∼ 2 %). However, a high correlation (R > 0.96) and root mean square errors (RMSE < 0.03) are also observed. Previous studies found SSA retrieval uncertainties in the range of 0.02–0.03 (Eck et al., 2003; Corr et al., 2009; Jethva et al., 2014; Kazadzis et al., 2016), whereas the correction, when high NO₂ is recorded, is usually higher. Thus, it is clear that in conditions of high NO₂ concentrations an accurate characterization of this gas is necessary in order to avoid noticeable bias in the affected AERONET channel around 440 nm.

4 Summary and conclusions

The retrievals of aerosol properties from sun photometers may be affected by NO₂ absorption in the observed spectral range, and thus, accurate assumptions on NO₂ concentrations are highly desirable. Currently, some ground-based aerosol networks, such as SKYNET, do not take NO₂ optical depth into consideration in AOD retrieval processes, while others (e.g., AERONET) use satellite-based NO₂ climatology for estimating it. However, significant errors could be introduced in the AOD retrievals, especially over urban areas, where NO₂ variability can be high and also the occurrence of high NO₂ events is more frequent. Such errors may occur only in the cases where NO₂ is not taken into account or the used NO₂ climatology underestimates such high NO₂ events.

Actual co-located surface-based NO₂ measurements (e.g., from Pandora instruments) or spaceborne observations with improved spatial and temporal resolution (e.g., S5P/TROPOMI) may be helpful for reducing the uncertainty in the NO₂ optical depth contribution in later versions of the AOD retrieval algorithms. In this study, we evaluated the possible improvements of AOD and AE retrievals by applying a specific correction using synchronous and co-located measurements of the total NO₂ column from Pandora spectroradiometers and the TROPOMI satellite sensor. For this purpose, we used multiannual (2017–2022) observations from both AERONET and SKYNET multispectral AOD observations co-located with Pandora instruments and collected over two locations in Rome (Italy) with different anthropic pressure (one in the city center and the other in a suburban area).

The deviations of the NO₂-modified AOD retrievals from the network standard products were investigated. AERONET-used NO₂ climatology was found to systematically underestimate Pandora-measured NO₂ over both sites. The impact of the correction is higher in the case of SKYNET, since the NO₂ optical depth is not considered at all in the standard retrieval processes of that network. At the same time, the observed differences in the results between the two networks can also be partly explained by the different channels used for the retrievals. For both AERONET and SKYNET, a low but systematic AOD overestimation was found. Although in most of the cases the differences are lower than 0.01 for AOD and lower than 0.1 for AE retrievals, the correction can still be useful for lower AODs (< 0.3) where the majority of observations are found, especially under high NO₂ pollution events. The mean AOD bias derived for the high NO₂ cases (> ∼ 0.7 DU) is ∼ 0.011 ± 0.003 at 440 nm and ∼ 0.012 ± 0.003 at 380 nm for AERONET APL-SAP and ∼ 0.009 ± 0.003 at 440 nm and ∼ 0.010 ± 0.003 at 380 nm for AERONET CNR-ISAC. The mean AE bias for the high NO₂ is ∼ 0.08 ± 0.04 for both Rome AERONET sites. In the case of SKYNET, the mean bias for the cases with high NO₂ levels (> ∼ 0.7 DU) is ∼ 0.018 ± 0.003 and ∼ 0.10 ± 0.05 for AOD and AE, respectively. Overall, the average biases in AOD retrievals are systematic but within the reported AOD uncertainties. However, they are important enough to be reported here, as AOD retrieval uncertainties not linked with instrument calibration (e.g., Rayleigh, ozone, and NO₂-related optical depths) are considered to have an upper limit of 0.005 as a goal for sun photometers, according to WMO (2005). As expected, the effect of improved NO₂ assumption in the retrievals is more evident in both AOD and AE when the actual synchronous ground-based Pandora NO₂ measurements are employed, compared to the situations when the used correction was based on daily or monthly averaged Pandora data or TROPOMI NO₂ retrievals. The use of TROPOMI NO₂ data is a demonstration of the possibility for corrections on a global scale. However, the underestimation of NO₂ concentrations by TROPOMI compared to Pandora NO₂ data for Rome leads to lower AOD corrections.

In addition, a first attempt to evaluate the impact of those corrections on AOD and AE annual trends was conducted. However, the aerosol data sets employed in this trend analysis are quite short for a robust trend analysis. Here only quantitative comparisons are performed for each individual data set, i.e., corresponding to specific instrument and site, before and after the NO₂-based correction. Although the effect of NO₂ on the derived trends seems to be insignificant, and the linear fit trend calculations introduce uncertainties similar to or higher than the NO₂ effects on AOD, the more pronounced impact may be expected for trends derived from larger data sets and in the case of a significant NO₂ trend.

We also investigated the possible effects of the proposed NO₂ optical depth correction on the agreement between ground-based and spaceborne AOD retrievals. In particular, we compared MODIS DB AOD retrievals at 470 nm with AERONET and SKYNET AOD products. In general, the agreement between ground-based (both AERONET and SKYNET) and MODIS DB AOD is quite good, revealing a correlation coefficient (r) higher than 0.7. The use of Pandora NO₂ in the sun photometer retrievals introduces a slight improvement in the absolute values of ∼ 0.003 in the agreement between MODIS DB and AERONET AOD and an improvement of ∼ 0.006 between MODIS and SKYNET observations. Although the impact on the comparisons between spaceborne and ground-based observations of AOD is quite small, it can be quite useful for eliminating or decreasing possible biases in the intercomparisons of satellite and ground-based data in situations with NO₂ concentrations typical of highly polluted areas.

Finally, we investigated the impact of using a precise characterization of the total NO₂ concentration on the SSA retrieval at 440 nm from AERONET measurements. For this, the GRASP algorithm was used to evaluate the effect of NO₂ correction on AERONET aerosol retrievals obtained by inverting TOD and almucantar radiances at 440, 675, 870, and 1020 nm. GRASP aerosol retrieval, using the actual total NO₂ concentration provided by the co-located Pandora over both stations selected for this study, were compared with GRASP retrievals mimicking AERONET operational retrievals. The results showed that, in general, the effect in the retrieved parameters at 440 nm under low NO₂ absorption conditions was not significant. At the same time, for the cases with high NO₂ loadings (> 0.7 DU), important changes in the retrieved SSA were observed, with an average positive bias of 0.02 (2 %) for both locations.

In general, the effect of NO₂ absorption can be relatively important in the retrievals of aerosol properties, especially AE, AOD, and SSA at 440 and 380 nm, when NO₂ is not included in the retrieval algorithms or in cases where NO₂ absorption is significantly higher than the NO₂ climatology used. If NO₂ absorption is taken from climatological data, then the accuracy of such approach may not be sufficient at locations where NO₂ has high diurnal variability during high NO₂ concentration episodes that cannot be captured by the satellite climatology. In such situations, the use of accurate co-located NO₂ observations, e.g., by Pandora instruments, is highly desirable. Thus, based on the results of this study, the effect of NO₂ correction could be considered relatively small for a large fraction of the observations; nonetheless, the correction has certainly contributed towards lowering the uncertainty in AOD and, especially, aerosol SSA provided by sun photometers.

In future studies, the effect of NO₂ correction on the absorption Ångström exponent (AAE) could be explored. AAE is an aerosol optical property that describes the absorption variation with respect to wavelength and is significantly influenced by particle size, shape, and chemical composition used for aerosol characterization and apportionment studies (e.g., Schuster et al., 2006). Since AAE is a function of spectral AOD and SSA, the NO₂ correction for certain AOD wavelengths and SSAs, shown in this study, is expected to impact the AAE calculations towards lower values (as the NO₂-corrected AOD is systematically lower, and the corrected SSA is higher).

Finally, the improved technology including real-time NO₂ monitoring (e.g., the Pandonia network), real-time satellite-based products at high spatial resolution (e.g., TROPOMI), and the more precise NO₂ products foreseen (e.g., from Sentinel 4) tend to positively contribute towards improving retrieved aerosol properties in the spectral range (∼ 380–440 nm) affected by NO₂ absorption.