Cloud detection from multi-angular polarimetric satellite measurements using a neural network ensemble approach

This paper describes a neural network cloud masking scheme from PARASOL (Polarization and Anisotropy of Reflectances for Atmospheric Science coupled with Observations from a Lidar) multi-angle polarimetric measurements. The algorithm has been trained on synthetic measurements and has been applied to the processing of 1 year of PARASOL data. Comparisons of the retrieved cloud fraction with MODIS (Moderate Resolution Imaging Spectroradiometer) products show overall agreement in spatial and temporal patterns, but the PARASOL neural network (PARASOL-NN) retrieves lower cloud fractions. Comparisons with a goodness-of-fit mask from aerosol retrievals suggest that the NN cloud mask flags fewer clear pixels as cloudy than MODIS (∼ 3 % of the clear pixels versus ∼ 15 % by MODIS). On the other hand the NN classifies more pixels incorrectly as clear than MODIS (∼ 20 % by NN, versus ∼ 15 % by MODIS). Additionally, the NN and MODIS cloud mask have been applied to the aerosol retrievals from PARASOL using the Remote Sensing of Trace Gas and Aerosol Products (RemoTAP) algorithm. Validation with AERONET shows that the NN cloud mask performs comparably with MODIS in screening residual cloud contamination in retrieved aerosol properties. Our study demonstrates that cloud masking from multi-angle polarimeter (MAP) aerosol retrievals can be performed based on the MAP measurements themselves, making the retrievals independent of the availability of a cloud imager.

1 Introduction

Aerosols and clouds, acting as key components of the climate system, highlight a significant challenge to understand, quantify, and predict climate change (Bellouin et al., 2020). Anthropogenic aerosols have the potential to cause a radiative forcing that is comparable in magnitude to greenhouse gases but with an opposite sign. While the climate impact of greenhouse gases is relatively well understood, the cooling effect (negative forcing) caused by aerosols is the largest source of uncertainty in the latest Intergovernmental Panel on Climate Change (IPCC) assessment, as well as in previous assessments (Arias et al., 2021).

Aerosols affect Earth's climate by scattering and absorbing radiation (aerosol–radiation interactions) and acting as condensation nuclei for cloud droplets and as ice nucleating particles to promote ice formation. Substantial changes in anthropogenic aerosol emissions occurred in the industrial era, as the result of fossil fuel burning in transport, industry, and energy production. The most uncertain aspects in aerosol radiative forcing are related to aerosol–cloud interactions (Bellouin et al., 2020). By acting as cloud condensation nuclei (CCN), aerosols affect the cloud droplet number concentration (N_d) and, consequently, the cloud albedo, causing the radiative forcing due to aerosol–cloud interaction (RFaci) (Twomey, 1974). Subsequently, rapid adjustments take place in, e.g., cloud fraction (CF) and liquid water path (LWP) that result from an initial change in N_d. The combination of RFaci and adjustments results in effective radiative forcing due to aerosol–cloud interaction (ERFaci). Apart from their interactions with aerosols, clouds alone also represent an important uncertainty in future climate predictions. Specifically, the impact of changing temperatures on clouds, known as cloud feedback, remains poorly understood (Zelinka et al., 2022).

Satellite-based remote sensing provides essential information for the understanding and quantification of aerosol–radiation interactions (Myhre et al., 2009; Lacagnina et al., 2015, 2017; Chen et al., 2022) and aerosol–cloud interactions (Gryspeerdt et al., 2017; Hasekamp et al., 2019b; Quaas et al., 2020; Gryspeerdt et al., 2023). For aerosols, the richest set of information on aerosol properties can be obtained from instruments that measure both intensity and polarization of scattered sunlight at multiple wavelengths and multiple viewing angles for one ground pixel (Mishchenko and Travis, 1997; Hasekamp and Landgraf, 2007; Dubovik et al., 2019). In this paper, we refer to such an instrument as a multi-angle polarimeter (MAP). The POLarization and Directionality of Earth Reflectances (POLDER) instrument has flown in three different incarnations and has pioneered space-based MAP remote sensing. Until now, the only instrument of this type that has provided a multi-year data record is the POLDER-3 instrument on the PARASOL (Polarization and Anisotropy of Reflectances for Atmospheric Science coupled with Observations from a Lidar) micro-satellite, in orbit between 2004 and 2013. The successor mission of POLDER is the 3MI (Fougnie et al., 2018) on ESA/EUMETSAT MetOp SG-A satellite. The instrument suite of the NASA Phytoplankton, Aerosol, Cloud, and ocean Ecosystem (PACE) mission (Werdell et al., 2019), launched in February 2024, has the capability to improve substantially in both aerosol and cloud retrievals by performing multi-angle polarimetric measurements at higher accuracy, more wavelengths (hyper-spectral), and more viewing angles (hyper-angular) than existing instruments. PACE carries two polarimeters. (i) The SPEXone instrument (Hasekamp et al., 2019a), contributed by the Netherlands, provides new and improved measurements of aerosol microphysical properties (size distribution, refractive index, shape) and optical properties (aerosol optical depth, AOD, and single scattering albedo, SSA). (ii) The Hyper-Angular Rainbow Polarimeter-2 (HARP-2) has the capability to provide improved cloud properties (N_d, droplet size distribution, cloud phase) from the cloud bow in polarization (Grosvenor et al., 2018). The PACE mission is the first to provide polarimetric retrievals of aerosol and cloud properties in more than a decade.

Currently, aerosol properties from these instruments can only be retrieved in cloud-free areas or for areas where an aerosol is located above a cloud (Waquet et al., 2009, 2013; Knobelspiesse et al., 2015). Therefore, it is important to be able to identify cloud-free pixels (cloud screening) for which an aerosol retrieval algorithm can provide meaningful output. Stap et al. (2015b) have shown that the goodness of fit of aerosol retrieval from MAP measurements can be used for cloud screening. The idea behind this is that for a measurement that is contaminated by a cloud, the forward model of the aerosol retrieval cannot find a good fit with the measurement because the forward model of a clear-sky aerosol retrieval cannot fit the characteristic spectral and angular signals caused by scattering on cloud particles (e.g., the cloud bow in polarization, as shown in Fig. 1), as demonstrated by Stap et al. (2015a, b, 2016). However, a disadvantage of this approach is that it requires large computation time because the aerosol retrieval has to be applied on all pixels (cloudy and cloud-free), and the cloud-free pixels can only be identified after the retrieval. Here it is important to note that the aerosol retrieval procedure itself is computationally expensive because it needs online radiative transfer calculations in order to provide a sufficiently accurate result. Cloudy pixels would even require larger computation time because of the large optical thickness. Given that the results for the cloudy pixels (∼ 80 % of all pixels; Krijger et al., 2007) are discarded, processing of these cloudy pixels can be considered a waste of computation time. To avoid this redundant computation, cloud masking should be applied before the retrieval.

Figure 1 True reflectance from MODIS (a) and polarized radiance (b) at 865 nm measured by PARASOL as a function of scattering angle on 2 px (one covered by dust and another by cloud), 1 July 2008. The “peak” of polarized radiance near scattering angle 140° is the “cloud bow”, a distinctive angular polarization feature of liquid clouds.

The goal of this work is to develop a cloud screening procedure for aerosol retrievals from MAP instruments. Here, we focus on the POLDER-3/PARASOL instrument (hereafter simply referred to as PARASOL) as it is the only MAP that has provided a multi-year data set. Currently, the most accurate cloud mask for aerosol retrievals from PARASOL is provided by MODIS (Moderate-Resolution Imaging Spectroradiometer). The MODIS cloud mask (Ackerman et al., 1998), which provides measurements within 3 min from PARASOL (for the period PARASOL was a part of the NASA A-Train), is based on input signals from visible and infrared bands, which detect the high, spectrally flat reflectance and low brightness temperature feature of clouds. Although the cloud mask from MODIS has been shown to be useful for performing cloud screening for PARASOL aerosol retrievals, there are two important reasons to develop a cloud screening algorithm based on PARASOL alone: first of all, not all PARASOL pixels have a corresponding MODIS cloud retrieval. Only for the period 2005–2009 were PARASOL and MODIS measurements co-located in time and location because both instruments flew in the same orbit (NASA A-Train). So, for the period 2010–2013, no MODIS measurements are available to perform cloud screening for PARASOL. Also for 2005–2009, there are sometimes orbits missing in the MODIS data because of instrument switch-offs or processing problems. Another motivation for developing a PARASOL-only cloud mask is that the MAP measurements of PARASOL contain unique sensitivity to clouds, and this may result in better capability to screen for clouds than MODIS, e.g., in the case of dust aerosols (Wang et al., 2016).

Cloud fraction retrievals from PARASOL were performed operationally at 18 km × 18 km (Zeng et al., 2011). These cloud fraction retrievals are based on thresholds on measured radiances, polarized radiances, and apparent oxygen pressure at 6 km×6 km. To avoid the choice of a priori thresholds which define whether a pixel is cloud contaminated or not, we develop a retrieval of cloud fraction (CF) and use the retrieved cloud fraction to define one or more cloud masks (e.g., a strict and loose one). The cloud mask is at the native spatial resolution of PARASOL, which is ∼ 6 km × 6 km. Given that the MAP measurements are also affected by aerosol properties and surface reflection properties, a classical retrieval algorithm (by radiative transfer calculations and inversion) for the retrieval of CF would be very complex because it would need, in addition to CF, to include many fit parameters related to aerosol (e.g., amount, size, refractive index, layer height) and surface bidirectional reflection distribution function (BRDF) parameters. Apart from the fact that such a retrieval algorithm requires a prohibitive amount of computing power, it may also have the risk of not finding an optimal solution (e.g., by ending at a local minimum of the inversion cost function). The use of neural networks (NNs) is a promising method because of the efficiency in computation and the possibility to define a reduced state vector (with just CF), while still taking into account the effect of aerosol and surface reflection. NNs have been used successfully in the polarimetric remote sensing of aerosols by, e.g., Di Noia et al. (2017), Gao et al. (2021a), and Gao et al. (2021b), as well as for the polarimetric remote sensing of cloud microphysical properties by Di Noia et al. (2019).

The paper is organized as follows: Sect. 2 introduces the data used in the study; Sect. 3 describes the training methodology of the neural network with a consistency check based on a synthetic evaluation set; and Sect. 4 shows the data processing of 1 year (2008) PARASOL measurements and a comparison with MODIS, a goodness-of-fit mask from RemoTAP aerosol retrieval, and AERONET data. Finally, Sect. 5 summarizes and concludes this study.

2 Data description

2.1 PARASOL

PARASOL was operational from 2004 to 2013, being part of the NASA A-Train satellite constellation in synthesis with MODIS/AQUA (multi-spectral imager), CALIPSO (lidar), and CloudSat until 2009. PARASOL has provided intensity measurements in nine spectral bands (443, 490, 565, 670, 763, 765, 865, 910, 1020 nm) and linear polarization (Stokes parameters Q and U) in three spectral bands (490, 670, 865 nm), viewing a ground scene under (up to) 16 different viewing geometries (Fougnie et al., 2007). The level-1 (non-cloud-screened) observations are available on a sinusoidally projected grid of $\sim \mathrm{6}\mathrm{km}\times \mathrm{6}\mathrm{km}$ pixels, named the full-resolution (FR) grid. In this work, we only use PARASOL measurements within the latitude ranges from 60° S to 60° N that have at least 14 viewing angles available.

2.2 PARASOL RemoTAP aerosol retrievals

The Remote sensing of Trace gas and Aerosol Products (RemoTAP) algorithm (Hasekamp et al., 2011; Fu and Hasekamp, 2018; Fu et al., 2020; Lu et al., 2022) is a flexible algorithm that can be used for the retrieval of optical and microphysical aerosol properties from multi-angle polarimeter (MAP) measurements, retrieval of aerosol properties and trace gas columns (e.g., carbon dioxide (CO₂) and methane (CH₄)) from spectroscopic measurements, or for joint retrievals using multiple instruments (e.g., MAP and spectrometer together). RemoTAP is based on iteratively fitting a linearized radiative transfer model (Hasekamp and Landgraf, 2002, 2005; Schepers et al., 2014) to the measurements of intensity and polarization of light reflected by the Earth's atmosphere and surface. It has large flexibility in the definition of parameters to be retrieved and allows retrievals over land, ocean, and clouds.

For the setup of the aerosol retrieval, we utilize a parametric three-mode aerosol description, where the size distribution is characterized by three log-normal modes (i.e., N_modes=3), comprising one fine mode and two coarse modes (dust and soluble). The state vector for the fine mode includes parameters such as effective radius r_eff, effective variance v_eff, aerosol column number N_aer, spherical fraction f_sph, and refractive index coefficients c_k, which correspond to the standard refractive index spectra of inorganic aerosol (real part), black carbon (imaginary part), organic carbon (imaginary part), and water. The dust mode consists of non-spherical dust, and its state vector includes r_eff, N_aer, and a coefficient for the imaginary part of the dust refractive index. The fixed parameters are f_sph=0, v_eff=0.6, and c_k=1 for the real part of the dust refractive index. The third mode represents a coarse soluble mode, and its state vector includes r_eff, N_aer, and a coefficient c_k of the inorganic refractive index spectrum. The fixed parameters for this mode are f_sph=1 and v_eff=0.6. A detailed description can be found in Lu et al. (2022).

2.3 Cloud fraction from MODIS-Aqua cloud mask product

The MODIS Cloud Mask product used in this work is a level-2 product generated at 1 km (at nadir) spatial resolutions from MODIS-Aqua (MYD35_L2). The algorithm, as described in Ackerman et al. (1998), employs a series of visible and infrared threshold and consistency tests to specify the confidence that an unobstructed view of the Earth's surface is observed. An indication of shadows affecting the scene is also provided. There are four cloud flag categories in the product: confidently cloudy, uncertain cloudy, probably clear, and definitely clear. The cloud fraction (referred to as MODIS cloud fraction hereafter) is calculated as the fraction of confidently and uncertain cloudy flagged 1 km resolution MODIS pixels within a 6 km × 6 km PARASOL grid.

2.4 AERONET data

The RemoTAP aerosol retrievals are validated by data from the Aerosol Robotic Network (AERONET; Holben et al., 2001). The aerosol optical depth (AOD) at 550 nm is from AERONET level-2 data (version 3.0). The single scattering albedo (SSA) at 550 nm is from AERONET level-2 Almucantar retrieval inversion products (Dubovik and King, 2000; Dubovik et al., 2000, 2002).

3 Neural network design

3.1 Training set generation

Our neural network approach for the retrieval of cloud fraction from PARASOL data was decomposed into two subtasks, with different neural networks dedicated to pixels over land and ocean, respectively.

The training of the NNs in this work uses synthetic measurements created with the radiative transfer model (RTM): LINTRAN (Schepers et al., 2014). In our setup, LINTRAN computes the top-of-atmosphere intensity vector I with Stokes parameters I, Q, and U as components as a function of wavelength and viewing-solar geometry. In the simulation process, the ocean reflection properties are parameterized with wind speed as in Cox and Munk (1954) and the chlorophyll-a concentration as in Chowdhary et al. (2006) and Fan et al. (2019). For the simulations over land, the surface bidirectional reflectance distribution function (BRDF) is parameterized using the Ross–Li model (Wanner et al., 1995), and the bidirectional polarization distribution function (BPDF) is parameterized as in Maignan et al. (2009). Liquid clouds are described by spherical particles with a gamma size distribution, with the refractive index of water taken from Hess et al. (1998). For ice clouds, hexagonal crystals with varying aspect ratios and surface distortions are used as proxies for ice crystals with variable complex shapes (van Diedenhoven et al., 2020). Uniform distributions for liquid/ice cloud effective radius are used in order to make sure that the NN has equal sensitivity to clouds with both small and large particles. The aerosol size distribution is based on three log-normal modes as in Lu et al. (2022), where each mode is described by the effective radius, effective variance, complex refractive index (wavelength dependent), aerosol optical depth (AOD) at 550 nm, fraction of spherical particles, and aerosol layer height. The aerosols are considered homogeneously distributed between 0 and 1 km in our simulations.

The distribution of the different input parameters for the training set is described in Table 1. The cloud properties are taken from a random distribution. To roughly represent the true cloud fraction distribution and also to emphasize small cloud fractions (which matters for cloud masking), 40 % of cloud fractions are generated between 0.2 and 1, 20 % are between 0 and 0.2, and 20 % completely clear and 20 % completely cloudy pixels are also included in the data set. The aerosol and surface (land and ocean) properties are from randomly picked pixels of RemoTAP global retrieval for the year 2008. The combination of solar zenith angles, viewing zenith angles, and relative azimuth angles is randomly picked from PARASOL level-1 measurements for the year 2008. Here it should be noted that only the measurements performed at a minimum of 14 angles are considered for training the neural networks. The choice was made in order to avoid dealing with the difficulty of having an input vector of variable size or, as an alternative, of passing input vectors with missing data to the neural networks.

Table 1 Details of the statistical distributions of the aerosol and cloud parameters used to generate the training data sets. Distribution of “RemoTAP” means properties are randomly taken from the 2008 PARASOL-RemoTAP L2 database.

The NN was designed to produce cloud fraction (CF) as output, so the partially cloudy measurements and “true” cloud fraction are needed in the training process. To model the intensity vector I for a partially cloudy pixel with cloud fraction f, the independent pixel approximation (IPA) was used:

$$\begin{array}{}\text{(1)}& {\mathit{I}}_{\text{ipa}}=f{\mathit{I}}_{\text{cloudy}}+(\mathrm{1}-f){\mathit{I}}_{\text{clear}},\end{array}$$

where the vectors I_cloudy and I_clear correspond to simulations for cloudy conditions (either fully covered by liquid cloud or ice cloud) and clear conditions. The cloud screening task for aerosol retrieval pays more attention to small cloud fractions, where a larger sensitivity is required. Therefore, instead of feeding the cloud fraction directly into the NN, we apply the inversion of the sigmoid function:

$$\begin{array}{}\text{(2)}& T=\ln \left({\displaystyle \frac{f}{\mathrm{1}-f}}\right),\end{array}$$

where T is the training target fed in the NN, and f is the original cloud fraction. Here it should also be noted that f is clipped between 10⁻⁵ and ($\mathrm{1}-{\mathrm{10}}^{-\mathrm{5}}$), as the function T has no definition at f=0 or f=1. We have also investigated other functions for T (logarithmic, linear, step functions), but the inverse sigmoid provides the best results for cloud screening.

Based on the same set of aerosol and surface properties, synthetic forward calculations for three scenes should be generated: fully clear, fully covered by liquid cloud, and fully covered by ice cloud (simply referred to as clear scene, liquid-cloudy scene, and ice-cloudy scene hereafter). We do not consider situations that are partly covered by both ice and liquid clouds as experiments indicate that this does not improve the NN performance. For each clear–cloudy combination, 20 randomly generated cloud fractions (10 for liquid clouds and 10 for ice clouds) were applied using IPA. Moreover, considering the situations where a cloud can be observed in a viewing angle but missed at another angle, because of the parallax effect or cloud movement during the multi-angle data acquisition (Stap et al., 2016), 20 % of pixels are assigned perturbed cloud fractions in different viewing angles during the IPA procedure. The cloud fraction angular perturbation range is restricted within both 0.2 absolute value and 100 % of the cloud fraction itself, which corresponds best to the angular patterns seen in the PARASOL measurements. The full training set is randomly divided into two parts: about 7.5 million samples in the training set and 1.6 million in the test set.

Training NNs on the whole 7.5 million sample training set requires considerably large computational time and memory. Therefore, instead of training one NN on the whole training set, we choose the “neural network ensemble” approach (Hansen and Salamon, 1990), in which the whole training set is equally and randomly separated into several parts (here we divide the set into 16 parts), and on each part of the training set, an NN is individually trained. The final output is defined as the average of outputs from all the NN ensembles:

$$\begin{array}{}\text{(3)}& f_{\text{final}}=S\left({\displaystyle \frac{\mathrm{1}}{N_{\text{ens}}}}\sum _i^{N_{\text{ens}}}{T}_i\right),\end{array}$$

where S is the sigmoid function, N_ens is the number of NN ensembles, and T_i is the direct output of the ith NN ensemble (defined in Eq. 2). As is already explained and tested in Di Noia et al. (2019), such an approach can achieve similar performance as training one NN on the full training set.

The measurement noise is modeled as a Gaussian random number with a zero mean and a standard deviation of 1 %–3 % relative noise for intensity and 0.012 absolute noise for degree of linear polarization (DoLP). Such noise is added to the training set as a form of regularization, as explained in Bishop (1995).

The input variables for the NNs are the intensity, DoLP at 14 viewing angles, and the corresponding viewing geometries (solar zenith angle, viewing zenith angle, relative azimuth angle, and scattering angle). The intensity and DoLP, as a function of wavelength and viewing angle, are compressed using a principal component analysis (PCA) before the training. A total of 25 principal components are retained for radiance and 33 for DoLP, in order to compress the measurements, which are the same numbers as in Di Noia et al. (2015).

The NNs are implemented using PyTorch (version 1.11.0, https://pytorch.org/, last access: 11 October 2021) and designed as multi-layer perceptrons (MLPs). The training of the NNs utilizes the backpropagation (BP) algorithm (Rumelhart et al., 1986) and batch training with a batch size of 12 000. The Adam optimizer (Kingma and Ba, 2014) is employed to minimize the RMSE (root mean square error) loss function. The neural network architecture consists of three hidden layers with 80 neurons each layer for over-ocean scenes and 40 neurons for over-land scenes. These setups are chosen based on comparison with a goodness-of-fit cloud mask applied to real PARASOL measurements.

3.2 Synthetic results

Before applying the NNs to real PARASOL measurements, it was tested on the 1.6 million sample test set outside the training set. In Table 2, the results of the neural network cloud fraction retrieval on the test set for both the land and ocean schemes are presented, showcasing the performance metrics including bias, mean absolute error (MAE), correlation (corr), and root mean square error (RMSE).

Table 2 Bias, MAE, correlation, and RMSE of cloud fraction on synthetic test set.

To derive a cloud mask from the cloud fraction predictions, we need to set thresholds for the cloud fraction to determine whether a pixel is cloudy or not. Therefore, to find a suitable threshold, the performance of the cloud mask should be assessed, and we define the following three criteria:

1. Information loss. This is calculated as the fraction of clear pixels that are wrongly flagged as cloudy, relative to all clear pixels. A low information loss is desired for not discarding too many clear pixels where aerosol information can be retrieved.

2. Effectiveness. This is calculated as the fraction of correctly flagged cloudy pixels with respect to all cloudy pixels. When the “effectiveness” is higher, fewer actually cloudy pixels are wrongly flagged as clear, so there are fewer cases in which aerosol retrievals are applied on actually cloudy pixels.

3. Overall agreement. This is defined as the fraction of correctly flagged pixels.

As an example, the performance of the corresponding cloud masks over land on the synthetic test set is summarized in Fig. 2. Here we evaluate two cases, one with a very strict definition of truly cloud-free pixels ($\text{CF}<{\mathrm{10}}^{-\mathrm{4}}$) and one with a somewhat looser definition (CF<0.01). Then we evaluate the NN cloud mask for different thresholds of the retrieved CF. The overall agreement, effectiveness, and information loss decrease when thresholds for the retrieved cloud fraction (CF) are increased. This is because a looser cloud mask (a larger threshold) is less likely to incorrectly discard clear pixels; however, it also has a higher chance of including cloudy pixels. Comparing the results for the very strict definition of cloud-free ($\text{CF}<{\mathrm{10}}^{-\mathrm{4}}$) and the looser one (CF<0.01), we see that the NN has difficulty in identifying scenes with $\text{CF}<{\mathrm{10}}^{-\mathrm{4}}$, as the overall agreement drops by 15 %–20 % points when using the stricter definition of cloud-free. The same conclusions also hold for the retrievals over ocean. Furthermore, we also obtain very similar performance on the synthetic measurement data set created using an aerosol layer height that varies between 1 and 8 km. This demonstrates that the NN cloud fraction retrieval is not hampered by the fact that it is trained with the assumption that aerosols are located in the lowest 1 km of the atmosphere.

Figure 2 Performance of cloud mask over land on the synthetic test set for different thresholds. (a) “true” CF<0.01 as truly clear and (b) $\text{CF}<{\mathrm{10}}^{-\mathrm{4}}$ as truly clear. The red bar stands for overall agreement. The blue bar represents the effectiveness of masking out cloud, and the green bar is the information loss.

For applying a cloud mask before further aerosol retrievals, an important aspect is that the cloud mask has a low information loss (keeps more clear pixels) instead of a very high effectiveness (fewer redundant aerosol retrievals on cloudy pixels). Therefore, we suggest that a threshold of 0.05 for the retrieved cloud fraction should be good for aerosol retrieval.

4 Application to PARASOL data

The test results on the synthetic data set provide only an initial consistency check to conclude the algorithm works theoretically. To further evaluate the performance of the NNs, we apply it to PARASOL L1 data in 2008. The results are compared with MODIS for the full year and with the RemoTAP goodness-of-fit cloud mask for 2 d. In addition, the MODIS and NN cloud mask are applied to RemoTAP aerosol retrievals for the year 2008 and the remaining aerosol properties are validated with AERONET.

4.1 Global cloud fraction distribution in 2008

Figure 3 shows maps of the annual mean cloud fraction for both the PARASOL neural network (PARASOL-NN) retrieval and MODIS. Both MODIS and PARASOL-NN see high levels of annual cloud fraction average in regions around the Equator and around 60° S. In the desert areas (Sahara, Kalahari, Mojave, Atacama, and Great Victoria desert), MODIS and PARASOL-NN observe quite a low annual cloud fraction average. The latitudinal cloud fraction distribution also agrees with the trend between PARASOL-NN and MODIS. However, PARASOL-NN retrieves a cloud fraction that is significantly lower (∼ 0.15 on average) than MODIS.

Figure 3 Whole-year cloud fraction average (on each $\sim \mathrm{6}\mathrm{km}\times \mathrm{6}\mathrm{km}$ PARASOL grid, shown in the colored world map) and latitudinal distribution (cloud fraction average along each longitude as a function of latitude, shown in the histogram). (a) NN annual cloud fraction average and (b) MODIS annual cloud fraction average.

Figure 4 shows the cloud fraction comparison between PARASOL-NN and MODIS for spring (March–May), summer (June–August), autumn (September–November), and winter (December–February). Similar to the full-year average, the NN retrieves a lower cloud fraction than MODIS for each season, but they both capture the same seasonal dependence of clouds in the tropics, subtropics, and also mid-latitudes. In tropical regions, near the Equator, cloudiness tends to remain relatively high throughout the year. Subtropical regions (20–30° latitude), on the contrary, exhibit lower cloud fraction compared to the tropics. Seasonal variation of cloud fraction can be observed in mid-latitudes (30–60°), but more significant variation appears in the monsoon-dominated regions. As an example, the Indian peninsula experiences distinct seasonal differences of a large cloud fraction in summer (due to the influence of the monsoon) and a low cloud fraction in winter. Both MODIS and PARASOL-NN agree well with those seasonal trends, but PARASOL-NN produces a lower cloud fraction.

Figure 4 Seasonal cloud fraction average on each ∼ 6 km × 6 km PARASOL grid, shown in the colored world map. From the top to the bottom row, spring (March–May, a, b), summer (June–August, c, d), autumn (September–November, e, f), and winter (December–February, g, h) and latitudinal distribution (cloud fraction average along each longitude as a function of latitude, shown in the histogram). Panels (a), (c), (e), and (g) show the result from the NN cloud fraction and panels (b), (d), (f), and (h) is from the MODIS cloud fraction.

The main reason that PARASOL-NN retrieves a lower cloud fraction than MODIS is that PARASOL-NN retrieves a cloud fraction close to 1 (0.95–1) for significantly fewer pixels (∼ 50 %) than MODIS. As shown in Fig. 5, the number of pixels where PARASOL-NN gives a cloud fraction larger than 0.95 is about half of that of MODIS. For these pixels, the NN often retrieves cloud fraction in the range 0.70–0.95. Fortunately, this underestimation will not impact the performance of the cloud mask for aerosol retrieval, which requires the sensitivity mainly to small cloud fraction. Another piece of information from Fig. 5 is that PARASOL-NN predicts more cloud fraction below 0.05 than MODIS does, which also contributes to the underestimation in cloud fraction compared to MODIS. It should be noted that MODIS is expected to overestimate cloud fraction at the PARASOL ∼ 6 km × 7 km resolution because it only provides cloud fraction 0 (cloud-free) or 1 (fully cloudy) at the native MODIS ∼ 1 km × 1 km resolution.

Figure 5 Frequency histogram of annual cloud fraction from MODIS and PARASOL-NN. The NN retrieves much fewer cases where the CF is very close to 1 (0.95–1.0). Instead, it retrieves many more cloud fractions in the range 0.70–0.95. This feature leads to an overall underestimation compared to the MODIS cloud fraction.

Part of the reason that PARASOL-NN retrieves more small cloud fractions is that on the one hand, PARASOL-NN seems to “miss” small sub-pixel clouds, but on the other hand, MODIS seems to sometimes misinterpret dust aerosol as a cloud. We illustrate this in Fig. 6 with two examples where the cloud fractions from PARASOL-NN and MODIS are over-plotted on a background-corrected true reflectance image from MODIS (https://worldview.earthdata.nasa.gov/, last access: 26 January 2023). Over the Red Sea (panels a–c), MODIS interprets pixels with desert dust as “confidently cloudy”, while PARASOL-NN can correctly see it as clear without missing the clouds in the area. This difference could be attributed to the contribution of the polarization information in the NN input, which is sensitive to particle size and shape, allowing it to effectively distinguish between clouds and thick aerosols (Hasekamp, 2010; Stap et al., 2015b). Panels (d)–(f) show a partly cloudy area with small clouds (the ocean area to the south of Somalia). MODIS seems to see too large a cloud fraction for parts of this scene, while the NN predicts too low a cloud fraction. A possible reason for the cloud fraction underestimation by the NN is the movement of the clouds during the time of acquisition of all angles, leading to different cloud fractions sampled per viewing angle. Although we try to mitigate this effect by adding samples with different cloud fractions to each viewing angle in the NN training set, some small moving clouds are still missed.

Figure 6 Cloud fraction around the Red Sea (a–c) and to the south of Somalia (d–f) on 1 July 2008. Panels (a) and (d) are the MODIS corrected true reflectance from the NASA Worldview website. Panels (b) and (e) show the NN cloud fraction and panels (c) and (f) is MODIS cloud fraction.

Figure 7 shows a comparison between a strict MODIS (threshold 0.01) cloud mask and PARASOL-NN cloud mask for different cloud fraction thresholds. The overall agreement is about 80 %. The effectiveness decreases somewhat towards the larger cloud fraction threshold, which means that the percentage of filtered “true” (MODIS) cloudy pixels decreases with the threshold increases, whereas the information loss also decreases towards the larger cloud fraction threshold, which means that the percentage of wrongly flagged “true” (MODIS) clear pixels decreases with the increase of the threshold.

Figure 7 NN cloud mask compared with MODIS in the full-year processing. The red bar stands for overall agreement. The blue bar represents the effectiveness of masking out cloud, and the green bar is the information loss.

In conclusion, underestimation is observed, especially in areas where the cloud fraction is typically large (e.g., tropical regions along the Equator), which might indicate that NN is not suitable for some applications in cloudiness research. Nevertheless, the performance of PARASOL-NN at pixels with a large cloud fraction does not affect the suitability for cloud screening before aerosol retrievals, where the focus is on pixels with a small cloud fraction.

4.2 Comparison to RemoTAP goodness-of-fit mask

As an additional assessment of the capability of the PARASOL-NN cloud mask, we perform RemoTAP aerosol retrievals on all PARASOL pixels, both clear and cloudy, on 2 d (1 January and 1 July 2008) to acquire a goodness-of-fit-based cloud mask. The high computational cost prevents us from a whole year retrieval without cloud screening, but this 2 d period contains both situations in summer and in winter, and thus should be reasonably representative.

The goodness of fit from a MAP aerosol retrieval is considered a good identifier of clear pixels, although it is important to realize that the RemoTAP goodness-of-fit filter is not fully independent of the PARASOL-NN cloud mask because they both rely on the same PARASOL level-1 measurements. The thought behind the goodness-of-fit cloud mask is that the forward model in the aerosol retrieval scheme cannot reproduce cloud signals, which results in a large χ². Based on previous works (Stap et al., 2015a, b), χ²<5 is a good threshold to identify cloud-free pixels.

The NN cloud fraction, MODIS cloud fraction, and RemoTAP goodness-of-fit mask (blue area for χ²<5) on 1 July 2008 are shown in Fig. 8. Both the NN and MODIS cloud fractions reveal a similar spatial pattern, while NN predicts more pixels with a cloud fraction less than 0.1 (plotted in blue) than MODIS, consistent with the finding that the NN cloud fraction is overall smaller than for MODIS, as expected from the annual comparison. In addition, both of them usually lead to cloud fractions either larger than 0.9 (mostly cloudy) or less than 0.1 (mostly clear). The goodness-of-fit mask has a smaller cloud-free area than both NN and MODIS. Here, it should be noted that χ²>5 can also happen due to other reasons than cloud contamination, i.e. pixels with large measurement errors, pixels for which the inversion ends in a local minimum, or pixels where the aerosol and/or surface model in the RemoTAP forward model is not suited to describe the actual atmospheric and/or surface conditions.

Figure 8 NN cloud fraction (a), MODIS cloud fraction (b), and goodness-of-fit mask (c, blue area for χ²<5) on 1 July 2008.

Figure 9 shows a comparison of the NN and MODIS cloud masks with the goodness-of-fit mask for four cloud fraction thresholds. From this comparison, it follows that the information loss of NN is considerably less than that of MODIS for most thresholds. This indicates that NN keeps more cloud-free pixels (as identified by the goodness-of-fit mask), and this feature benefits a larger aerosol data coverage. The effectiveness of NN, however, is lower than that of MODIS, leading to more redundant aerosol retrievals that result in large χ². Nevertheless, the effectiveness decrease is relatively small. As a balance of effectiveness and information loss, the threshold of 0.05 is recommended, and it will be used to test the PARASOL-NN cloud mask effect on retrieved aerosol properties in the next section.

Figure 9 Comparison of the cloud masks with goodness-of-fit mask on 1 January and 1 July 2008. (a) The comparison on 1 January. (b) The comparison on 1 July. The comparison between the NN cloud mask and goodness-of-fit mask is shown by shaded bars, while the comparison between the MODIS and goodness-of-fit mask is in unshaded bars. The red bars stand for overall agreement. The blue bars represent the effectiveness of masking out cloud, and the green bars are the information loss.

4.3 Effect on retrieved aerosol properties

Undetected clouds can cause substantial biases in retrieved aerosol properties, such as an overestimation in AOD. However, too strict a cloud mask leads to reduced data coverage, especially in areas important in the study of aerosol–cloud interactions. In this section, we assessed the quality of the retrieved aerosol properties after applying the PARASOL-NN or MODIS cloud mask.

We perform aerosol retrievals using RemoTAP for all AERONET-collocated PARASOL pixels in 2008 without applying any prior cloud filter. Then we evaluate the retrieved aerosol properties against AERONET after applying different cloud masks. We focus on the validation with AERONET on three aerosol properties:

1. the aerosol optical depth (AOD; τ) at 550 nm;

2. the single scattering albedo (SSA) at 550 nm;

3. the Ångström exponent (AE) for the wavelength pair 440–865 nm, defined as

   $$\begin{array}{}\text{(4)}& \mathit{\alpha }=-{\displaystyle \frac{\log ({\mathit{\tau }}_{\mathrm{440}}/{\mathit{\tau }}_{\mathrm{865}})}{\log (\mathrm{440}/\mathrm{865})}}.\end{array}$$

We evaluate the performance against the requirements formulated by the Global Climate Observing System (GCOS). For AOD the GCOS requirement is that the AOD error should be smaller than 0.03 or 10 % (whichever is greater). For AERONET validation, this requirement has been modified in the Aerosol_cci study (Popp et al., 2016) to 0.04 or 10 % to also take into account the uncertainty in AERONET AOD. For SSA the GCOS requirement is that the error should be smaller than 0.03. This requirement is not modified for AERONET evaluation given that the 0.03 requirement is considered already loose (Popp et al., 2016). For AE we use a requirement of 0.2 in line with what was used in the ESA HARPOL project (https://www.sron.nl/harpol, last access: 8 January 2023). It should be noted that the error for AOD, AE, and SSA is defined as the absolute value of the substraction between truth (AERONET data) and retrieval (PARASOL RemoTAP retrieval).

Figure 10 shows an overview of the percentage of pixels within the accuracy requirements, the RMSE, and the bias corresponding to the three cloud masks. For both AOD and AE, applying a cloud mask (either MODIS or PARASOL-NN) can slightly improve the percentage of retrievals within the requirements, which suggests that either the goodness-of-fit-filtered data still contain some cloud-contaminated pixels or that the cloud filters screen out some other difficult pixels. Moreover, it is interesting to see the bias of AE is negative when applying the goodness-of-fit mask only and the PARASOL-NN cloud mask, while that of MODIS is positive instead. Nevertheless, for AE the biases for the three different cloud filters are also small. The overall performance of the three different cloud masks is similar, but PARASOL-NN keeps more pixels than MODIS (NN keeps 163 905 px, MODIS keeps 160 791 px, and goodness of fit keeps 213 578 px), which is also consistent with the results in the 2 d (1 January and 1 July 2008) comparison.

Figure 10 Effect of different cloud masks on retrieved aerosol properties. (a) Percentage of retrievals within requirements (GCOS/CCI requirement for AOD and SSA. HARPOL requirements for AE). (b) RMSE. (c) Bias of the retrievals after applying additional cloud masks: NN cloud mask (red bar), MODIS cloud mask (blue bar), and only goodness-of-fit mask (green bar). Note that NN cloud mask and MODIS cloud mask are applied on top of goodness-of-fit mask.

5 Conclusion

This paper presents an algorithm to filter clouds for aerosol retrievals from multi-angle, multi-wavelength polarimetric measurements. The proposed approach is based on neural networks trained on synthetic measurements, where variations in aerosol and surface properties, as well as cloud macro- and microphysical properties, were taken into account. Separate NNs have been trained for scenes over ocean and over land. With the NN we retrieve a cloud fraction at the native PARASOL level-1B spatial resolution of ∼ 6 km × 6 km, and this cloud fraction is subsequently used for cloud masking.

The neural network algorithm has been applied to process the entire PARASOL data set for the year 2008, and a comparison to the MODIS cloud fraction has been performed for seasonally and yearly averaged data. The comparison shows an overall similarity in spatial patterns between the PARASOL-NN and MODIS cloud fraction, while the PARASOL-NN cloud fraction is around 0.15 smaller than MODIS. The lower cloud fraction is mainly caused by the fact that NN retrieves a cloud fraction close to 1 (0.95–1) for much fewer pixels than MODIS. Also, NN retrieves more pixels with a small cloud fraction (< 0.05) because on the one hand NN sometimes misses small sub-pixel clouds, and on the other hand MODIS sometimes misinterprets dust aerosol as clouds.

The comparison of the PARASOL-NN and MODIS cloud masks with the RemoTAP goodness-of-fit mask for 2 d (1 January and 1 July) shows the NN keeps more clear pixels for aerosol retrieval (low information loss), while the effectiveness in filtering clouds is lower than MODIS. Nevertheless, the decrease in effectiveness is acceptable due to the substantially lower information loss. Taking the balance of effectiveness and information loss, the NN cloud mask with a CF threshold of 0.05 is recommended for cloud screening before aerosol retrievals.

We also studied the effect of the cloud masks on retrieved aerosol properties. There is no significant difference in the ability to remove the cloud-contaminated aerosol retrievals between PARASOL-NN and MODIS, while NN keeps more pixels than MODIS.

Based on the experiments above, it can be concluded that in the aerosol retrieval process, the NN cloud mask based on MAP measurements alone can replace the prior cloud mask provided by a cloud imager like MODIS. However, a goodness-of-fit mask is still needed to be applied after aerosol retrieval finished, because aerosol retrievals can also fail for reasons other than cloud contamination, i.e. pixels with large measurement errors, pixels for which the inversion ends in a local minimum, or pixels where the aerosol and/or surface model in the RemoTAP forward model is not suited to describe the actual atmospheric and/or surface conditions.

The proposed algorithm could have a potential future application in analyzing data from recently developed multi-angle polarimeters. By adjusting the instrument-specific factors, such as the number of viewing angles and spectral channels and tuning noise configurations during the neural network's training process, the algorithm can be adapted to work with instruments like the 3MI (Fougnie et al., 2018) on ESA/EUMETSAT MetOp SG-A satellite; the multi-angle polarimeter (MAP) instrument on the Copernicus CO2M mission (Spilling and Thales, 2021); and SPEXone (Hasekamp et al., 2019a) and the Hyper-Angular Rainbow Polarimeter (HARP-2), both of which were launched successfully in February 2024 on the NASA Plankton Aerosol, Cloud, ocean Ecosystem (PACE) mission (Werdell et al., 2019). Application of the NN cloud screening approach to these new instruments will provide insight into the sensitivity of the approach to measurement uncertainty, number of viewing angles, and number of spectral bands.