iCOR Atmospheric Correction on Sentinel-3/OLCI over Land: Intercomparison with AERONET, RadCalNet, and SYN Level-2

: To validate the iCOR atmospheric correction algorithm applied to the Sentinel-3 Ocean and Land Color Instrument (OLCI), Top-of-Atmosphere (TOA) observations over land, globally retrieved Aerosol Optical Thickness ( AOT ), Top-of-Canopy (TOC) reﬂectance, and Vegetation Indices (VIs) were intercompared with (i) AERONET AOT and AERONET-based TOC reﬂectance simulations, (ii) RadCalNet surface reﬂectance observations, and (iii) SYN Level 2 (L2) AOT , TOC reﬂectance, and VIs. The results reveal that, overall, iCOR’s statistical and temporal consistency is high. iCOR AOT retrievals overestimate relative to AERONET, but less than SYN L2. iCOR and SYN L2 TOC reﬂectances exhibit a negative bias of ~ − 0.01 and − 0.02, respectively, in the Blue bands compared to the simulations. This diminishes for RED and NIR, except for a +0.02 bias for SYN L2 in the NIR. The intercomparison with RadCalNet shows relative differences < ± 6%, except for bands Oa02 (Blue) and Oa21 (NIR), which is likely related to the reported OLCI “excess of brightness”. The intercomparison between iCOR and SYN L2 showed R 2 = 0.80–0.93 and R 2 = 0.92–0.96 for TOC reﬂectance and VIs, respectively. iCOR’s higher temporal smoothness compared to SYN L2 does not propagate into a signiﬁcantly higher smoothness for TOC reﬂectance and VIs. Altogether, we conclude that iCOR is well suitable to retrieve statistically and temporally consistent AOT , TOC reﬂectance, and VIs over land surfaces from Sentinel-3/OLCI observations.


Introduction
The European Commission's Copernicus Programme Sentinel-3 mission [1] was designed to ensure the long-term collection and operational delivery of high-quality measurements to Copernicus ocean, land, and atmospheric services. The mission's main objective is to measure ocean and land surface color, sea surface topography, and sea and land surface temperatures with high accuracy and reliability to support ocean forecasting systems and environmental and climate monitoring. The mission's definition is driven by the need for continuity in provision of the European Remote Sensing Satellite (ERS, [2]), the Environmental Satellite (ENVISAT, [3]) and Système Probatoire Observation de la Terre-Végétation (SPOT-VGT, [4]) data, with improvements in instrument performance and coverage.
Sentinel-3 currently comprises two satellites (Sentinel-3A and Sentinel-3B, hereafter referred to as S3A and S3B), launched on 16 February 2016 and 25 April 2018, respectively. Sentinel-3C and -3D are currently planned for launch no earlier than 2023 and 2025, In (i), a globally representative picture of the iCOR TOC reflectance retrieval accuracy is obtained using ground-truth TOC reflectances obtained from OLCI TOA reflectance observations and Radiative Transfer Model (RTM) simulations. AERONET AOT and Total Column Water Vapor (TCWV) retrievals are used as auxiliary input for the TOA to TOC inversion, assuming that the AERONET retrievals provide the most accurate characterization of the atmospheric state. Similar approaches have been applied within ACIX-I and -II. In (ii), the iCOR OLCI TOC reflectance retrievals are compared with observed Bottom-of-Atmosphere (equivalent to TOC) reflectances at RadCalNet sites that are convolved on the OLCI Spectral Response Functions (SRFs). Finally, in (iii), iCOR OLCI TOC reflectances are intercompared to SYN L2 TOC reflectances on their statistical and temporal consistency over five globally distributed Regions of Interest (ROIs).
The paper's outline is as follows. Section 2 introduces the iCOR and SYN atmospheric correction algorithms, while Section 3 provides details on the validation data and methods, respectively. Validation results are described and discussed in Sections 4 and 5, respectively, and final conclusions are drawn in Section 6.

iCOR
iCOR uses the Moderate-Resolution Atmospheric Radiance and Transmittance Model-5 (MODTRAN5, [16]) for its radiative transfer calculations, works with Look-Up Tables (LUTs) to speed up the retrievals, and is a surface adaptive atmospheric correction algorithm. It identifies whether a pixel is water or land and subsequently applies a dedicated atmospheric correction. iCOR runs without user interaction and derives the required input parameters from the image. In August 2017, iCOR was released to the broader user community as a Sentinel Application Platform (SNAP) toolbox plug-in for Landsat-8/OLI and Sentinel-2/MSI. iCOR4S3 was released in 2019 and is also available as a SNAP plug-in. An iCOR4S3 processing service is available at ESA's Grid Processing on Demand for Earth Observation Application (G-POD, https://gpod.eo.esa.int/).
The iCOR image based AOT retrieval divides a TOA reflectance image into macropixels of about 15 × 15 km 2 , being large enough to include high spectral variation and small enough to assume spatial atmospheric homogeneity. In the first retrieval step, the lowest TOA radiance values in the different visible bands are searched to determine an upper AOT boundary value for each macro-pixel. Subsequently, this AOT value is refined based on the spectral variation within the macro-pixel, using a multi-parameter end-member inversion technique. Five pixels with high spectral contrast (selected on TOA NDVI values) are represented by a linear combination of three pre-defined default vegetation spectra and a soil spectrum.
If a macro-pixel contains too many clouds and/or if the internal spectral TOA reflectance variability is too small (i.e., no sufficient variability in the TOA NDVI values), no AOT retrieval is performed and the AOT is obtained through a bilinear interpolation from surrounding macro-pixels for which the retrieval was successful. Currently, for iCOR4S3 no maximum distance for this macro-pixel interpolation is set. However, in the version available on G-POD users can set a scene-specific AOT for cases without successful macropixel AOT retrieval. Further, currently no AOT fallback solution has been implemented for iCOR4S3. For Sentinel-2/MSI, the operational iCOR workflow contains a fallback to Copernicus Atmosphere Monitoring Service (CAMS) AOT data and it is foreseen that a similar approach will be implemented in a future iCOR4S3 version. Finally, to up-scale the macro-pixel values to the nominal pixel resolution, a cubic convolution method is used. iCOR applies a rural aerosol model in the TOA to TOC reflectance conversion.

SYN L2
The Sentinel 3 operational SYNERGY processor is divided into two modules. First, a SYN Level-1 (L1) module is dedicated to the co-registration of OLCI and SLSTR acquisitions and to the production of internal datasets collecting all OLCI/SLSTR radiometric Remote Sens. 2021, 13, 654 4 of 24 measurements and contextual parameters (e.g., meteorological input data, as well as geolocation and angular data) on the same reference grid. Subsequently, these internal datasets are considered as inputs to a SYN Level-2 (L2) module, including the aerosol retrieval and aerosol correction processing steps.
In the atmospheric correction algorithm, the AOT is retrieved using collocated OLCI and SLSTR TOA radiances. These acquisitions are filtered for cloud contamination and averaged onto a coarser grid than the original images (15 × 15 pixels, i.e.,~4.5 × 4.5 km 2 ). The AOT retrieval is performed following an approach combining Dark Object and multiple viewing angle methods derived from existing work for MERIS and AATSR [17]. For a given atmospheric aerosol model (currently a continental model is used) and a parameterized AOT computed at 0.55 µm, a set of surface reflectances is derived from OLCI/SLSTR observations. The derived error metric is then minimized with respect to angular and spectral constraints, using the Brent and Powell methods for univariate and multivariate minimization, respectively. The 300 m AOT is obtained from the relative contributions of the three surrounding super-pixel AOTs on a 15 × 15-pixel area centered at a given 300 m pixel. This interpolated 300 m AOT is then used in the atmospheric correction to obtain the SYNERGY Surface Directional Reflectance (SDR, equivalent to TOC reflectance).
The resulting AOT and SDR associated with all relevant SLSTR and OLCI channels are projected on a regular 300 m resolution grid to be included in the SYN L2 product.

Regions of Interest and iCOR Retrievals
Five regions of interests (ROIs, see Table 1) were defined, covering different types of land cover and atmospheric conditions: Western Europe (WEUR), Northern Brazil (NBRA), Central Africa (CAFR), Eastern China (CHIN), and Northern Australia (NAUS). For January-September 2019, OLCI Level-1 Product Dissemination Units (PDUs, 3-min time frames) were processed by the iCOR4S3 implementation deployed on ESA's G-POD facility, using the standard iCOR4S3 settings [18]. This processing yielded AOT and TOC reflectance retrievals for 16 of the 21 OLCI bands (see Figure 1), omitting the spectral bands affected by oxygen (Oa13, Oa14, and Oa15) and water vapor absorption (Oa19 and Oa20). For January-September 2019, OLCI Level-1 Product Dissemination Units (PDUs, 3min time frames) were processed by the iCOR4S3 implementation deployed on ESA's G-POD facility, using the standard iCOR4S3 settings [18]. This processing yielded AOT and TOC reflectance retrievals for 16 of the 21 OLCI bands (see Figure 1), omitting the spectral bands affected by oxygen (Oa13, Oa14, and Oa15) and water vapor absorption (Oa19 and Oa20). To restrict the number of PDUs retrieved for the ROIs, the maximum cloud cover percentage was arbitrarily set to 20% for CAFR, NBRA, and CHIN, 50% for WEUR, and 5% for NAUS (see Table 1). These percentages were chosen as a trade-off between having sufficient PDUs for the analyses and reducing processing time. Cloud coverage percentages vary between ROIs, as cloud occurrence differs between the ROIs. Therefore, to ensure the automatic selection of a comparable number of PDUs for all ROIs, different thresholds were applied.

AERONET
AERONET is a network of ground-based Sun photometers that measure atmospheric aerosol properties. The measurement system comprises a solar-powered CIMEL Electronique 318A spectral radiometer that measures Sun and sky radiances at a number of fixed wavelengths within the visible and near-infrared (VNIR) spectrum. AERONET provides continuous cloud-screened observations of spectral AOT, TCWV, and inversion aerosol products in various aerosol regimes. For more than 25 years, AERONET has provided long-term, continuous, and readily accessible public domain databases of aerosol optical, microphysical, and radiative properties for aerosol research. The applications are among others the characterization and validation of satellite retrievals, as well as providing synergy with other aerosol databases. The network requires the standardization of instruments and calibration procedures, as well as unified processing and data distribution. AERONET provides globally distributed spectral AOT observations, aerosol inversion To restrict the number of PDUs retrieved for the ROIs, the maximum cloud cover percentage was arbitrarily set to 20% for CAFR, NBRA, and CHIN, 50% for WEUR, and 5% for NAUS (see Table 1). These percentages were chosen as a trade-off between having sufficient PDUs for the analyses and reducing processing time. Cloud coverage percentages vary between ROIs, as cloud occurrence differs between the ROIs. Therefore, to ensure the automatic selection of a comparable number of PDUs for all ROIs, different thresholds were applied.

AERONET
AERONET is a network of ground-based Sun photometers that measure atmospheric aerosol properties. The measurement system comprises a solar-powered CIMEL Electronique 318A spectral radiometer that measures Sun and sky radiances at a number of fixed wavelengths within the visible and near-infrared (VNIR) spectrum. AERONET provides continuous cloud-screened observations of spectral AOT, TCWV, and inversion aerosol products in various aerosol regimes. For more than 25 years, AERONET has provided long-term, continuous, and readily accessible public domain databases of aerosol optical, microphysical, and radiative properties for aerosol research. The applications are among others the characterization and validation of satellite retrievals, as well as providing synergy with other aerosol databases. The network requires the standardization of instruments and calibration procedures, as well as unified processing and data distribution. AERONET provides globally distributed spectral AOT observations, aerosol inversion (optical properties) products, and TCWV data for various aerosol regimes. The Version 3 AOT data are computed at three data quality levels: Level 1.0 (unscreened), Level 1.5 (cloud-screened and quality controlled), and Level 2.0 (quality-assured). In this study, Level 1.5 data were used. Table 1 provides a list of selected AERONET stations, the corresponding ROIs, the maximum set cloud cover, and the amount of processed OLCI frames.

6SV Simulations Using AERONET Observations
The earlier introduced TOC reflectance simulation approach was adapted from the Atmospheric Correction Intercomparison Exercise [12], using the Second Simulation of a Satellite Signal in the Solar Spectrum [19] Vector (6SV) V2.1 RTM, downloaded from http://6s.ltdri.org/, accessed 15 April 2020. The 6SV AERONET input data were the AOT, converted to 0.55 µm using the Ångström parameter between 0.44 and 0.675 µm, as well as TCWV and ozone concentrations. To obtain an indication on aerosol model impact to the iCOR retrievals, 6SV was run for its four standard aerosol models: continental, maritime, urban, and desert. See [20] for an overview of the 6SV aerosol models' optical properties and relative contributions of dust-like, water-soluble, oceanic, and soot components. The simulations were performed over AERONET stations with different land surface types (e.g., homogeneous and heterogeneous vegetation, bare soils, urban areas, desert areas, etc.). In total, 11 AERONET sites from the existing list of the ACIX-II AERONET land sites were selected (see Table 1) as a trade-off between global distribution, land surface type, and OLCI versus AERONET matchup availability.
Starting from the list of iCOR G-POD runs (see Table 1 for their distribution over the five ROIs), for each AERONET station and OLCI frame for which iCOR was run, the corresponding Sentinel-3A or -3B OLCI Level-1 TOA reflectance, geolocation, angular files, and SYN flags file were collected. Subsequently, when all files were available and the AERONET station was confirmed to be located in the given OLCI frame, the AERONET station's Level 1.5 data were read and checked for available AOT, TCWV, and ozone observations within +/−30 min of the OLCI frame time stamp.
The AERONET AOT at 0.44 µm and Ångström coefficient derived over the 0.44-0.675 µm spectral range were used to calculate the AOT at 0.55 µm: With λ and λ o being 0.55 µm and 0.44 µm, respectively, and α the Ångström coefficient. The median AOT 0.55 value obtained in the 60 min around the OLCI frame time stamp was taken as 6SV input. For ozone and TCWV, the arithmetic means of the available observations were used.
6SV input files were composed for OLCI frames for which the AERONET input data were available. The following variables are required by 6SV: and atmospheric observations to derive TOA reflectance that is generally used for postlaunch radiometric calibration and validation of optical imaging sensor data [15]. In this study, only TOC reflectance data were considered and intercompared with iCOR TOC reflectances. The TOC reflectance from RadCalNet is provided at 10 nm spectral intervals from 0.4-2.5 µm (depending on the site's instrumentation). It is acquired every 30 min between 09:00 and 15:00 Local Standard Time (LST) and only in cloud-free conditions ( Table 2). The network currently includes four sites located in the USA, France, China, and Namibia. A new sandy site in China called Baotou has recently been added. In this study, the choice of sites suitable for the intercomparison analysis is based on the assessment of the "surface reflectance variability" at 500 × 500 m 2 . Considering that a Sentinel-3 OLCI FR pixel is approximately 300 m, a surface reflectance variability less than 3% is considered acceptable for the analysis, resulting in only the Gobabeb (GONA, Lat: 23.600 • S, Lon: 15.119 • E) and Railroad Valley (RVUS, Lat: 38.497 • N, Lon: 115.690 • W) sites considered suitable for this research. GONA and RVUS are desert sites, characterized by sand and gravel with some widely scattered dry grass, and by a high-desert playa surrounded by mountains, respectively. The current uncertainty estimates in the in-situ surface reflectance measurements are 3.5-5.3% for RVUS and 3.5-5.0% for GONA. Figure 2 shows the average and standard deviation in TOC reflectance measured at GONA (black) and RVUS (red) over the considered time over the two sites. The plots demonstrate a high reflectance stability over time, with standard deviations of~0.02 at GONA and~0.03 at RVUS.

Sentinel-3 SYN L2
Sentinel-3 Synergy L2 (SYN L2) products are provided in Non-Time Critical timeliness (NTC, available 30-48 h after acquisition) and include daytime and land-only datasets. Surface reflectance and aerosol parameters are provided on a 300 m regular orthogeolocated grid, similar to the OLCI L2 products. Some sub-sampled annotations (solar and viewing azimuth/zenith angles) and meteorological data (mean sea-level pressure, More information on RadCalNet, its instrumentation, and analysis procedures can be found in [15].

Sentinel-3 SYN L2
Sentinel-3 Synergy L2 (SYN L2) products are provided in Non-Time Critical timeliness (NTC, available 30-48 h after acquisition) and include daytime and land-only datasets. Surface reflectance and aerosol parameters are provided on a 300 m regular ortho-geolocated grid, similar to the OLCI L2 products. Some sub-sampled annotations (solar and viewing azimuth/zenith angles) and meteorological data (mean sea-level pressure, total column ozone, and water vapor, all derived from ECMWF data) are also provided in the OLCI tie points data grids.
TOC reflectances and their associated error estimates are provided for each solar SLSTR and OLCI channel, except for the oxygen absorption bands Oa13, Oa14, and Oa15, the water vapor absorption bands Oa19 and Oa20, and the dedicated SLSTR S4 cloud detection band. Aerosol parameters, such as AOT at 0.550 µm (referenced as T550) and Ångström coefficient around this wavelength (referenced as A550) are included in separate files and are also provided on a 300 m resolution grid. The datasets are associated with contextual parameters, such as geographical position, time stamps, quality flags concerning the SYN L2 processing and OLCI and SLSTR input quality, and classification flags.

Sampling Strategy
For the intercomparison of iCOR and SYN L2 TOC retrievals with corresponding 6SV TOC reflectance simulations, for each AERONET station and OLCI band only dates at which iCOR, SYN L2, and the AERONET data (to calculate the 6SV reference TOC reflectance) were available were selected. Subsequently, for the three datasets the arithmetic mean and standard deviation were computed for those pixels in the 3 × 3-pixel area centered at the station's geolocation for which both SYN L2 and iCOR had valid values.
In the pairwise data extraction for the statistical consistency analysis at frame level, systematic subsampling was performed to reduce the number of data points, extracting every 10th pixel in the X and Y directions. For the statistical analysis over all ROIs and PDUs, an additional random sampling of 2% was applied.
The CLOUD_flags layer as included in the SYN L2 product (CLOUD, CLOUD_AMBIG UOUS, CLOUD_MARGIN, and SNOW_ICE) is used to exclude observations perturbed by clouds or snow/ice. Because currently no list of recommended flags when interpreting SYN L2 data is available, no additional flagging based on SYN L2 flag bands was applied. However, SYN L2 observations with negative AOT at 0.55 µm (~0.002% of the pixels) were discarded.
By applying the ratios in these vegetation indices (VIs), any 'noise' resulting from the observation geometry (directionality) is expected to be reduced, because these directional effects are similar in the visible and near-infrared bands [21], thereby providing a better indication on the temporal variations.

RadCalNet Intercomparison
The methodology adopted for comparing iCOR OLCI TOC reflectances with in-situ RadCalNet observations followed an approach similar to [22]. First, OLCI Level-1 data were selected considering cloud-free conditions over the RadCalNet sites and satellite in -situ data availability. OLCI-L1 FR selected data were then processed using iCOR4S3 deployed on ESA's G-POD. Afterwards, the iCOR OLCI TOC reflectance and satellite acquisition time from the pixel over the RadCalNet sites were extracted. Finally, the RadCalNet TOC reflectance data were convolved using the OLCI mean SRF and interpolated in time to the satellite overpass time, before comparing the two datasets. A description of the workflow is provided below, while the workflow is depicted in Figure 3.
-situ data availability. OLCI-L1 FR selected data were then processed using iCOR4S3 deployed on ESA's G-POD. Afterwards, the iCOR OLCI TOC reflectance and satellite acquisition time from the pixel over the RadCalNet sites were extracted. Finally, the RadCalNet TOC reflectance data were convolved using the OLCI mean SRF and interpolated in time to the satellite overpass time, before comparing the two datasets. A description of the workflow is provided below, while the workflow is depicted in Figure 3.
1. RadCalNet TOC reflectances were cubically interpolated from 30 min to 1 sec for every wavelength; 2. TOC reflectance values at different wavelengths were extracted at the sensor overpass time; 3. TOC reflectances from point 2 were cubically interpolated from 10 nm to 0.1 nm ( ); 4. TOC reflectances from point 3 were convolved with the OLCI mean SRF ( ), using Equation (5)
TOC reflectance values at different wavelengths were extracted at the sensor overpass time; 3.

Validation Metrics
The validation metrics used in the intercomparison are listed in Table 3. To quantify the deviation between iCOR, SYN L2 TOC reflectances and the 6SV reference simulations, and between iCOR or SYN L2 AOT and the AERONET AOT, the Accuracy, Precision, and Uncertainty (often abbreviated as APU) statistical metrics [23] were calculated. For these cases, X refers to the reference (6SV, RadCalNet or AERONET), and Y to iCOR or SYN L2.
To identify the relationship between iCOR and SYN L2 TOC reflectances and AOT, the geometric mean (GM) regression model is used. Such an orthogonal (model II) regression is appropriate, because-unlike when comparing to an absolute reference-both datasets are subject to noise. By applying an eigen decomposition to the X (SYN L2) and Y (iCOR) covariance metrics, two eigenvectors are obtained that describe the principal axes of the point cloud [24], i.e., the regression line. The GM regression intercept and slope value are added as quantitative information related to the scatterplots. The Root Mean Squared Difference (RMSD) expresses the difference magnitude between two datasets from 0 and is an expression of the overall difference, similar to the Uncertainty. Table 3. List of validation metrics, with n the number of valid samples used for the comparison, σ(X) and σ(Y) the standard deviation of X and Y, σ(X, Y) the co-variation of X and Y, andX andŶ estimated using the GM regression model. P(d i ), P(d i+1 ) and P(d i+2 ) are three consecutive observations on dates d i , d i+1 , and d i+2 .

Validation Metric Formula
Accuracy (Acc) or mean bias Root of the unsystematic mean product difference (RMPDu) Root of the systematic mean product difference (RMPDs) The GM regression model is also used to differentiate between systematic and random differences, providing additional information on the difference's nature between two datasets. In order to be comparable to the RMSD's magnitude, the root of the systematic and unsystematic mean product difference is used (RMPD s and RMPD u ). The coefficient of determination (R 2 ) indicates the agreement or covariation between two datasets with respect to the linear regression model, summarizing the total explained variance by this model.
The temporal smoothness δ [21,25] of AOT, TOC reflectance, and the VI time series is a measure for its consistency throughout time. Temporal smoothness is evaluated by taking three consecutive observations and by computing the absolute difference value between the center P(d n+1 ) and the corresponding linear interpolation between the two extremes P(d n ) and P(d n+2 ) The Time Series Smoothness Index (TSI) is an estimate of the time series noise [26], and is defined such that a lower value indicates less temporal noise and thus a smoother time series profile. Figure 4 shows the AOT and TOC reflectance intercomparison between iCOR and SYN (Y axis) with AERONET and 6SV, respectively (X axis), for S3A (left column) and S3B (right column) and for OLCI bands Oa03, Oa09, Oa17, and Oa21. Table 4 presents the number of selected data points included for the intercomparison, separated for S3A and S3B. Note that for CAFR (Bujumbura) no collocations were found for both S3A and S3B, resulting from both quasi-persistent cloud coverage and irregular AERONET observations in general. In the TOC reflectance scatterplots, the error bars in the X direction span the standard deviation of the 6SV-simulated TOC reflectances over the four aerosol models (continental, maritime, urban, and desert). , (e,f) Oa09, (g,h) Oa17, and (i,j) Oa21. SYN and iCOR values are indicated by red and blue dots, respectively. The black dashed lines indicate the 1:1 line, while the red and blue solid lines denote the GM linear regressions for SYN and iCOR, respectively. The error bars in the X direction span the standard deviation of the 6SV-simulated TOC reflectances for the four aerosol models: continental, maritime, urban, and desert. The legend text presents the APU values, as well as the GM linear regression equations and corresponding R 2 values. Note the different X and Y axis ranges for Oa03, Oa09 and Oa17, Oa21. The scatterplots show that both iCOR and SYN AOT retrievals overestimate relative to AERONET for both S3A and S3B, but that the SYN AOT overestimation is substantially larger. Further, the SYN AOT retrievals for S3B have a larger overestimation than for S3A. Another observation is that for S3B, iCOR AOT > SYN AOT for very low AERONET AOT (<~0.1). The larger overestimation for SYN compared to iCOR also reflects in higher (less favorable) APU values than iCOR. For the Blue band (Oa03) TOC reflectances, both SYN and iCOR generally exhibit an underestimation, with the SYN underestimation being larger, especially for S3B. This underestimation gradually reduces with increasing wavelength, due to decreasing aerosol contribution to atmospheric scattering and absorption. Further, the sometimes large error bars for the 6SV TOC reflectance indicate a large spread in the obtained 6SV TOC reference reflectances when different aerosol models are applied. This aerosol model susceptibility indicates that additional research on dynamic aerosol model application for iCOR (which currently uses the MODTRAN rural aerosol model) could be beneficial to further optimize the iCOR TOC reflectance retrieval performance. For Oa09 (Red), the SYN and iCOR TOC reflectances are generally in good agreement with their 6SV counterparts, with a marginal positive bias (Acc = 0-0.002) for SYN, while iCOR shows a slight negative bias (Acc =~−0.01). For Oa17 and Oa21 (NIR), agreement between SYN, iCOR, and 6SV has increased, which can be expected due to the small aerosol contributions in these bands. iCOR results show virtually no bias relative to 6SV, while for SYN a slight positive bias is seen (Acc = +0.017 and Acc = +0.021 for S3A and S3B, respectively).

Intercomparison with 6SV Simulations Using AERONET Input
The APU and R 2 spectral signatures are shown in Figure 5. The results are presented separately for S3A and S3B in the left and right columns, respectively. For all statistical metrics, a similar signature can be seen for S3A and S3B. For Accuracy, SYN has lower negative values (larger negative bias) than iCOR up to~0.5 µm. Further, the SYN negative bias changes into positive from~0.55 µm onwards, while for iCOR a small negative bias remains until 0.8 µm, becoming slightly positive beyond this wavelength. The Precision evolution for iCOR and SYN is similar, especially for S3A. However, for S3B differences in Precision are larger, with SYN > iCOR up to~0.65 µm and iCOR > SYN beyond this wavelength, with differences further increasing (i.e., more scattered TOC reflectance retrievals) for λ > 0.75 µm. Uncertainty is larger for SYN compared to iCOR in the BLUE and NIR channels, which indicates a larger spread for SYN. The R 2 increases towards 1 for both iCOR and SYN with increasing wavelength, with SYN having lower values than iCOR at λ < 0.55 µm, with this difference being more prominent for S3B.

Intercomparison with RadCalNet Observations
iCOR OLCI TOC reflectances were intercompared with convolved RadCalNet TOC reflectances at the GONA and RVUS sites, following the methodology described in Section 3.3.4. Initially, to reduce any Bi-directional Reflectance Distribution Function (BRDF) impact, only satellite images over the RadCalNet sites with VZA < 10° were considered. However, this constraint provided only a few matchups. Thanks to an ad-hoc processing by the Centre National d'Études Spatiales (CNES) that provided BRDF-corrected data in OLCI acquisition geometry [27], all possible VZAs could be considered, which allows for inclusion of more matchups in the analysis. The former approach was applied to the RVUS site, the latter to the GONA site. Table 5 shows the datasets used in the analysis for the two sites.

Intercomparison with RadCalNet Observations
iCOR OLCI TOC reflectances were intercompared with convolved RadCalNet TOC reflectances at the GONA and RVUS sites, following the methodology described in Section 3.3.4. Initially, to reduce any Bi-directional Reflectance Distribution Function (BRDF) impact, only satellite images over the RadCalNet sites with VZA < 10 • were considered. However, this constraint provided only a few matchups. Thanks to an ad-hoc processing by the Centre National d'Études Spatiales (CNES) that provided BRDF-corrected data in OLCI acquisition geometry [27], all possible VZAs could be considered, which allows for inclusion of more matchups in the analysis. The former approach was applied to the RVUS site, the latter to the GONA site. Table 5 shows the datasets used in the analysis for the two sites.
Given the OLCI FR data 300 m resolution and the surface reflectance variability at the RadCalNet sites, the iCOR TOC (surface) reflectance values were extracted for 1 pixel at GONA and for an average of 3 × 3 pixels centered at RVUS.  Figure 6 shows the average relative difference and standard deviation at (a) GONA and (b) RVUS, together with (c,d) the corresponding APU values for GONA and RVUS, respectively. Band Oa01 (0.4 µm) was excluded from the analysis, because the RadCalNet acquisitions start from 0.400 µm and the OLCI SRF for Oa01 starts from 0.388 µm. -38 * -15/12/2018-08/10/2019 Given the OLCI FR data 300 m resolution and the surface reflectance variability at the RadCalNet sites, the iCOR TOC (surface) reflectance values were extracted for 1 pixel at GONA and for an average of 3 × 3 pixels centered at RVUS. Figure 6 shows the average relative difference and standard deviation at (a) GONA and (b) RVUS, together with (c,d) the corresponding APU values for GONA and RVUS, respectively. Band Oa01 (0.4 μm) was excluded from the analysis, because the RadCalNet acquisitions start from 0.400 μm and the OLCI SRF for Oa01 starts from 0.388 μm. In general, at GONA the average relative differences are within ±6% for all bands, except for Band Oa02 (0.412 μm) and Oa21 (1.02 μm). The larger relative difference of these two bands relate to OLCI S3A having an "excess of brightness", mainly in the Blue bands and the Oa21 (NIR) band (see https://sentinel.esa.int/web/sentinel/technicalguides/sentinel-3-olci/mission-performance, accessed 18 November 2020). It is emphasized that, because in the comparison with the 6SV TOC reflectance simulations (see Section 4.1) both iCOR and 6SV start from the OLCI Level 1B data, the excess of brightness is present in both datasets and is thus cancellated in the intercomparison results. iCOR on OLCI S3A and OLCI S3B show similar values, also when the dataset contains a different number of products. Figure 6b shows that the Accuracy and Uncertainty for bands Oa17, Oa18, and Oa21 are higher (implying less agreement) than for the other bands, reaching In general, at GONA the average relative differences are within ±6% for all bands, except for Band Oa02 (0.412 µm) and Oa21 (1.02 µm). The larger relative difference of these two bands relate to OLCI S3A having an "excess of brightness", mainly in the Blue bands and the Oa21 (NIR) band (see https://sentinel.esa.int/web/sentinel/technical-guides/ sentinel-3-olci/mission-performance, accessed 18 November 2020). It is emphasized that, because in the comparison with the 6SV TOC reflectance simulations (see Section 4.1) both iCOR and 6SV start from the OLCI Level 1B data, the excess of brightness is present in both datasets and is thus cancellated in the intercomparison results. iCOR on OLCI S3A and OLCI S3B show similar values, also when the dataset contains a different number of products. Figure 6b shows that the Accuracy and Uncertainty for bands Oa17, Oa18, and Oa21 are higher (implying less agreement) than for the other bands, reaching Oa21 Accuracy (Uncertainty) peak values of 0.039 (0.041) and 0.044 (0.046) for iCOR on S3A and S3B, respectively.
At the RVUS site (Figure 6c,d), the average relative differences are within ±6% for all bands, similar to GONA, except for Oa21 (+8%). Band Oa02 shows a negative difference (−3.5%) in comparison with the GONA TOC reflectances. This means that for Oa02 iCOR surface reflectance overestimates at GONA and underestimates at RVUS. Concerning the APU, Accuracy and Uncertainty show a gradual increase with wavelength. Accuracy  The boxplot shows a slightly lower median AOT for iCOR and a larger AOT interquartile range for SYN L2. In addition, the SYN L2 AOT maximum value is much higher, with outliers up to AOT = 3.2. Similar results were observed over all ROIs, except for NAUS, where AOT is generally low, interquartile ranges very small, and iCOR having slightly higher AOT than SYN L2. The GM regression density plot indicates relatively good AOT consistency for values around 0.2. iCOR exhibits higher values than SYN L2 for AOT <~0.2, while the opposite occurs at larger AOT. As for the boxplots, similar results are observed over all ROIs, except for NAUS. The AOT bias frequency histogram shows that the histogram is skewed towards positive values, indicating for a relatively large proportion of pixels with SYN L2 AOT > iCOR AOT. Table 6 summarizes the validation metrics overall and per ROI. For all ROIs, intercepts are above 0 and slopes are below 1. Bias is lowest for NAUS and highest for CHIN. The systematic differences (RMPDs) are large compared to the RMSD, i.e., the difference between the regression line and the 1:1 line is large.  (Figure 6c,d), the average relative differences are within ±6% for all bands, similar to GONA, except for Oa21 (+8%). Band Oa02 shows a negative difference (−3.5%) in comparison with the GONA TOC reflectances. This means that for Oa02 iCOR surface reflectance overestimates at GONA and underestimates at RVUS. Concerning the APU, Accuracy and Uncertainty show a gradual increase with wavelength. Accuracy (Uncertainty) maximum values are 0.029 (0.032) for band Oa21. The Precision shows values <0.01 for all bands.

Intercomparison with Sentinel-3 SYN L2
4.3.1. Statistical Consistency iCOR Versus SYN L2 Figure 7 shows (a) the iCOR versus SYN L2 AOT boxplot, (b) GM regression density plot, and (c) the bias frequency histogram. The boxplot shows a slightly lower median AOT for iCOR and a larger AOT interquartile range for SYN L2. In addition, the SYN L2 AOT maximum value is much higher, with outliers up to AOT = 3.2. Similar results were observed over all ROIs, except for NAUS, where AOT is generally low, interquartile ranges very small, and iCOR having slightly higher AOT than SYN L2. The GM regression density plot indicates relatively good AOT consistency for values around 0.2. iCOR exhibits higher values than SYN L2 for AOT <~0.2, while the opposite occurs at larger AOT. As for the boxplots, similar results are observed over all ROIs, except for NAUS. The AOT bias frequency histogram shows that the histogram is skewed towards positive values, indicating for a relatively large proportion of pixels with SYN L2 AOT > iCOR AOT. Table  6 summarizes the validation metrics overall and per ROI. For all ROIs, intercepts are above 0 and slopes are below 1. Bias is lowest for NAUS and highest for CHIN. The systematic differences (RMPDs) are large compared to the RMSD, i.e., the difference between the regression line and the 1:1 line is large.   (Figure 8) show that for all ROIs the iCOR interquartile range is larger than SYN L2 in bands Oa01-Oa05. Except for NAUS, the iCOR  TOC reflectance boxplots per band (Figure 8) show that for all ROIs the iCOR interquartile range is larger than SYN L2 in bands Oa01-Oa05. Except for NAUS, the iCOR TOC reflectances are slightly higher in bands Oa01-Oa03. NAUS is a particular case because the AOT in general is very low and the iCOR AOT values are higher than SYN L2 AOT (see Figure 7). For bands Oa06-Oa21, iCOR has a lower value range, and the TOC reflectances are lower. These findings are in line with the intercomparison results over the AERONET sites (Figure 4).
Remote Sens. 2021, 13, x FOR PEER REVIEW 17 of 26 TOC reflectances are slightly higher in bands Oa01-Oa03. NAUS is a particular case because the AOT in general is very low and the iCOR AOT values are higher than SYN L2 AOT (see Figure 7). For bands Oa06-Oa21, iCOR has a lower value range, and the TOC reflectances are lower. These findings are in line with the intercomparison results over the AERONET sites (Figure 4). Figure 8. Boxplots of TOC reflectance comparing SYN L2 (Blue) and iCOR (green) over (a-e) the 5 defined ROIs (see Table  1 for their definition).  Table 1 for their definition).
The GM regression density plots over all ROIs (Figure 9) show intercepts close to 0, indicating that high consistency is achieved for low TOC reflectance values. Slopes are generally below 1, i.e., with increasing reflectance values SYN L2 becomes increasingly higher than iCOR. These differences at higher reflectances are larger at the shorter wavelengths. The GM regression density plots over all ROIs (Figure 9) show intercepts close to 0, indicating that high consistency is achieved for low TOC reflectance values. Slopes are generally below 1, i.e., with increasing reflectance values SYN L2 becomes increasingly higher than iCOR. These differences at higher reflectances are larger at the shorter wavelengths. The bias frequency histograms ( Figure 10) show that the bias between SYN L2 and iCOR reflectances is centered around 0, except for Oa06 and Oa09, indicating a systematic bias for these intermediate wavelengths. This is confirmed in Figure 11, which shows the overall systematic (RMPDs) and unsystematic difference (RMPDu) as functions of wavelength. The systematic bias (the difference between the regression line and the 1:1 line) is relatively low for the Blue bands and beyond ~0.75 μm and shows a peak at intermediate wavelengths, whereas the unsystematic bias (i.e., scatter) is virtually unrelated to wavelength. For both RMPDs and RMPDu, examination of the separate ROIs (plots are not shown) indicate that the spectral evolutions are similar, except for NAUS, where the systematic bias remains low at longer wavelengths. The bias frequency histograms ( Figure 10) show that the bias between SYN L2 and iCOR reflectances is centered around 0, except for Oa06 and Oa09, indicating a systematic bias for these intermediate wavelengths. This is confirmed in Figure 11, which shows the overall systematic (RMPDs) and unsystematic difference (RMPDu) as functions of wavelength. The systematic bias (the difference between the regression line and the 1:1 line) is relatively low for the Blue bands and beyond~0.75 µm and shows a peak at intermediate wavelengths, whereas the unsystematic bias (i.e., scatter) is virtually unrelated to wavelength. For both RMPDs and RMPDu, examination of the separate ROIs (plots are not shown) indicate that the spectral evolutions are similar, except for NAUS, where the systematic bias remains low at longer wavelengths.
The observed bias results from the band combinations. In both NDVI and EVI, band Oa09 is used, for which the systematic bias iCOR-SYN L2 is smaller compared to the other bands. For the NDWI, only bands Oa17 and Oa21 are used, which both show a similar bias between iCOR and SYN L2. In this case, the band combination thus neutralizes the bias between the TOC reflectances.   The NDVI, NDWI, and EVI bias frequency histograms SYN L2-iCOR ( Figure 12) indicate that primarily iCOR NDVI is higher than SYN L2 NDVI. This bias is also visible, although less obvious, for EVI. For NDWI, a slight positive bias (i.e., iCOR NDWI < SYN L2 NDWI) is observed. The mean (standard deviation) differences are −0.019 (0.054), 0.003 (0.008), and −0.001 (0.052) for NDVI, NDWI, and EVI, respectively.
The observed bias results from the band combinations. In both NDVI and EVI, band Oa09 is used, for which the systematic bias iCOR-SYN L2 is smaller compared to the other bands. For the NDWI, only bands Oa17 and Oa21 are used, which both show a similar bias between iCOR and SYN L2. In this case, the band combination thus neutralizes the bias between the TOC reflectances.  The NDVI, NDWI, and EVI bias frequency histograms SYN L2-iCOR ( Figure 12) indicate that primarily iCOR NDVI is higher than SYN L2 NDVI. This bias is also visible, although less obvious, for EVI. For NDWI, a slight positive bias (i.e., iCOR NDWI < SYN L2 NDWI) is observed. The mean (standard deviation) differences are −0.019 (0.054), 0.003 (0.008), and −0.001 (0.052) for NDVI, NDWI, and EVI, respectively.
The observed bias results from the band combinations. In both NDVI and EVI, band Oa09 is used, for which the systematic bias iCOR-SYN L2 is smaller compared to the other bands. For the NDWI, only bands Oa17 and Oa21 are used, which both show a similar bias between iCOR and SYN L2. In this case, the band combination thus neutralizes the bias between the TOC reflectances. A similar conclusion can be drawn from the Geometric Mean Regression (GMR) density plots (Figure 13). Over all ROIs, all VIs show a very high coefficient of determination (R 2 = 0.92-0.96). The systematic bias between iCOR and SYN L2 NDVI is relatively large (0.01-0.03), resulting from an intercept well above 0 (except for NAUS) and slopes close to 1. The systematic bias is lower for EVI (0.00 to 0.02) and very low for NDWI (below 0.005). The unsystematic bias (i.e., scatter around the regression line) is also larger for NDVI (0.02-0.07) and EVI (0.01-0.07), compared to NDWI (within 0.01). A similar conclusion can be drawn from the Geometric Mean Regression (GMR) density plots (Figure 13). Over all ROIs, all VIs show a very high coefficient of determination (R 2 = 0.92-0.96). The systematic bias between iCOR and SYN L2 NDVI is relatively large (0.01-0.03), resulting from an intercept well above 0 (except for NAUS) and slopes close to 1. The systematic bias is lower for EVI (0.00 to 0.02) and very low for NDWI (below 0.005). The unsystematic bias (i.e., scatter around the regression line) is also larger for NDVI (0.02-0.07) and EVI (0.01-0.07), compared to NDWI (within 0.01). Figure 12. NDVI, NDWI, and EVI bias frequency histograms (SYN L2-iCOR).
A similar conclusion can be drawn from the Geometric Mean Regression (GMR) density plots ( Figure 13). Over all ROIs, all VIs show a very high coefficient of determination (R 2 = 0.92-0.96). The systematic bias between iCOR and SYN L2 NDVI is relatively large (0.01-0.03), resulting from an intercept well above 0 (except for NAUS) and slopes close to 1. The systematic bias is lower for EVI (0.00 to 0.02) and very low for NDWI (below 0.005). The unsystematic bias (i.e., scatter around the regression line) is also larger for NDVI (0.02-0.07) and EVI (0.01-0.07), compared to NDWI (within 0.01).

Temporal Consistency iCOR vs. SYN L2
Based on the temporal evolution plots, the TSI was computed over the common valid observations for the 11 AERONET sites. Results are presented in Figure 14 and show that iCOR values are considerably lower than those of SYN L2, indicating a higher AOT temporal smoothness. However, this higher iCOR AOT temporal smoothness does not propagate into a significantly higher temporal smoothness in the TOC reflectance and VI.

Temporal Consistency iCOR vs. SYN L2
Based on the temporal evolution plots, the TSI was computed over the common valid observations for the 11 AERONET sites. Results are presented in Figure 14 and show that iCOR values are considerably lower than those of SYN L2, indicating a higher AOT temporal smoothness. However, this higher iCOR AOT temporal smoothness does not propagate into a significantly higher temporal smoothness in the TOC reflectance and VI. An AOT, TOC reflectance, and VIs timeseries example is given in Figure 15 for Aubière LAMP (France, WEUR). For each AERONET station, the TSI was computed over the common valid observations. Figure 15 shows the earlier discussed persistent overestimation of SYN L2 AOT compared to AERONET. In addition, iCOR tends to overestimate AOT, but the difference compared to AERONET AOT is smaller than for SYN L2. This confirms the results presented in Figure 4. For the TOC reflectances, similar temporal patterns for iCOR and SYN L2 are seen. Some outliers are visible, which are possibly related to undetected thin clouds or cloud shadows. In addition, the temporal plots for the VIs show a high degree of similarity. Timeseries over AERONET stations in other ROIs (not shown) reveal a consistent behavior. An AOT, TOC reflectance, and VIs timeseries example is given in Figure 15 for Aubière LAMP (France, WEUR). For each AERONET station, the TSI was computed over the common valid observations. Figure 15 shows the earlier discussed persistent overestimation of SYN L2 AOT compared to AERONET. In addition, iCOR tends to overestimate AOT, but the difference compared to AERONET AOT is smaller than for SYN L2. This confirms the results presented in Figure 4. For the TOC reflectances, similar temporal patterns for iCOR and SYN L2 are seen. Some outliers are visible, which are possibly related to undetected thin clouds or cloud shadows. In addition, the temporal plots for the VIs show a high degree of similarity. Timeseries over AERONET stations in other ROIs (not shown) reveal a consistent behavior.

Discussion
The various intercomparisons generally show a high agreement between iCOR AOT, TOC reflectance, and VIs with the corresponding reference datasets. For AOT, iCOR has a higher agreement with AERONET observations than SYN L2. However, although iCOR has a slightly better agreement with 6SV than SYN L2 for the TOC reflectances, the differences are considerably smaller than for AOT. Despite AOT being important in the atmospheric correction, its impact appears rather limited on the differences between iCOR and SYN L2. iCOR and SYN L2 apply different aerosol models, rural and continental, respectively.
As was shown in the 6SV intercomparison results in Figure 4, applying a different aerosol model can substantially affect the obtained TOC reflectance. The rural aerosol model used by iCOR comprises 70% water-soluble and 30% dust-like particles [28], whereas the SYN L2 continental aerosol model from 6S is composed of 29% water-soluble, 70% dust-like, and 1% soot particles [20]. The resulting single scattering albedos are 0.94 and 0.89 for iCOR and SYN L2, respectively, which implies that more of the incoming and outgoing reflected solar irradiance is absorbed for SYN L2, especially in case of multiple interactions within the land surface-atmosphere system. However, it is not straightforward to assess the final effect of this difference, as it depends on, among others, the underlying surface, the AOT (higher probability of multiple scattering events), and atmospheric absorption becoming dominant over scattering with increasing wavelength (see e.g., [29]). The RadCalNet intercomparison shows that over bare desert surfaces iCOR provides consistent results for all bands, except for larger differences in the Blue bands and Oa21 (NIR). However, these larger differences for the mentioned bands can be attributed to the reported "excess of brightness" for OLCI on S3A. As this excess of brightness was observed on various intercomparisons of OLCI-A and OLCI-B TOA observations, any aerosol influence in the Blue would impact both sensors. More information on OLCI-A and -B performance is available in [30]. In the iCOR versus SYN L2 intercomparison, primarily for NDVI a systematic bias was shown, related to the relatively high systematic bias for band Oa17 and low bias for Oa09. In the case of NDWI, the combination using bands Oa17 and Oa21, which have a similar bias, results in a neutralization, thereby resulting in good consistency.
Because iCOR and SYN L2 are currently the only known atmospheric correction algorithms that are applied to Sentinel-3/OLCI observations, any conclusions that arise from this study are preliminary. Nevertheless, as iCOR has participated in ACIX-I [12] and its successor ACIX-II, some additional performance indication for usage on Sentinel-2/MSI is available. In ACIX-I, an AOT overestimation for AERONET AOT <~0.5 was found, which decreased for larger AOT values, resulting in an overall mean AOT difference of 0.15 ± 0.15. For surface reflectance, its retrievals exhibited a small bias compared to AERONET-based 6S TOC reflectance <~0.2, which increased for brighter surfaces, and resulted in an overall Accuracy (bias) of 0.000 (band 2, 0.492 µm)-0.032 (band 8, 0.833 µm).

Conclusions
This paper presented the validation of the iCOR atmospheric correction algorithm applied to Sentinel-3A and -B Ocean and Land Colour Instrument (OLCI) observations over land. iCOR-retrieved Aerosol Optical Thickness (AOT), Top-Of-Canopy (TOC) reflectance, and resulting Vegetation Indices (VIs) were compared with Aerosol Robotic Network (AERONET) AOT, AERONET-based Second Simulation of a Satellite Signal in the Solar Spectrum Vector (6SV)-simulated TOC reflectance, RadCalNet TOC reflectance observations, and corresponding SYN Level 2 (L2) retrievals. The comparisons were performed over 11 AERONET stations and 389 OLCI frames over 5 ROIs: Western Europe (WEUR), Northern Brazil (NBRA), Central Africa (CAFR), Eastern China (CHIN), and Northern Australia (NAUS). The comparisons focused on the statistical and temporal consistency of the retrieved reflectances and VIs. Further, the OLCI iCOR TOC reflectances were compared with collocated and convolved RadCalNet TOC reflectances at the Gobabeb (GONA, Namibia) and Railroad Valley (RVUS, USA) sites.
The intercomparison with observed AOT and 6SV-simulated reference TOC reflectances over AERONET stations reveals that for AOT both iCOR and SYN L2 overestimate relative to AERONET, with the SYN L2 overestimation being larger than for iCOR. For the TOC reflectances in band Oa03 (Blue), both SYN and iCOR generally underestimate, which gradually reduces with increasing wavelength due to the decreasing aerosol contribution on atmospheric scattering and absorption. For Oa09, SYN, iCOR, and 6SV TOC reflectances are generally in good agreement, with only marginal positive and negative bias for SYN and iCOR, respectively. For the NIR bands (Oa17 and Oa21), the iCOR retrievals show virtually no bias relative to 6SV, while SYN shows a slight positive bias.
The intercomparison of iCOR TOC reflectances with RadCalNet convolved TOC reflectance observations showed that there is good consistency for all bands for both S3A and S3B, except for the Oa02 (Blue) and Oa21 (NIR) bands. For the NIR bands, there is an overestimation of iCOR compared to RadCalNet. In Oa02, an over-and an underestimation are seen at GONA and RVUS, respectively.
The statistical consistency analysis between iCOR and SYN L2 confirms the higher AOT for SYN L2 and a larger AOT dynamic range, except for NAUS, where AOT is generally very low. In general, a high consistency between iCOR and SYN L2 is observed at low TOC reflectance values. At higher reflectance values, SYN L2 reflectances are larger than iCOR. Systematic bias is relatively low for the Blue bands and around 0.7 µm. Except for NAUS, the iCOR TOC reflectances are slightly higher in bands Oa01-Oa03. For bands Oa06-Oa21, iCOR has a lower dynamic range and TOC reflectances are lower than SYN L2. Temporal AOT evolution plots from iCOR4S3, SYN L2, and AERONET show a consistent overestimation of SYN L2 AOT, and-to a lesser extent-iCOR relative to AERONET. For the TOC reflectances and VIs, temporal patterns for iCOR and SYN L2 are similar. The iCOR AOT TSI is considerably lower (indicative of less noise) than SYN L2 AOT and is comparable to that of AERONET, but the smoother temporal iCOR AOT pattern does not result in a significantly smoother iCOR TOC reflectance and VIs time series than those of SYN L2.
In conclusion, after the successful validation of iCOR over inland waters, this validation study shows that iCOR applied to Sentinel-3/OLCI also performs well over land surfaces. Overall, a good consistency with the SYN L2 TOC reflectance products is observed. In addition, in comparison to SYN L2, iCOR AOT and TOC reflectance values are shown to be in better agreement with the AERONET AOT and simulated 6SV TOC reflectance reference, respectively. Future iCOR developments will involve: (i) the possibility to include additional aerosol models, and (ii) include Copernicus Atmospheric Monitoring Service (CAMS) AOT data for cases in which the AOT retrieval does not converge to a solution, using a similar approach as for iCOR on Sentinel-2/MSI. iCOR's previous operational applications, continuous development, and agility to different processing infrastructures ensure a stable, efficient, and high-quality processing performance, thereby making iCOR suitable for various additional operational implementations.