1 INTRODUCTION

Estimation of regional and global ocean primary productivity is crucial to the understanding of the role of oceans in both Earth's carbon cycles and global climate change. The proof-of-concept Coastal Zone Color Scanner (CZCS), which operated from 1978 to 1986, aptly demonstrated that general views of the marine biosphere obtained by spaceborne ocean color instruments could provide valuable data at spatial and temporal scales incomparable with shipboard or moored instrumentation. Since the cessation of the CZCS operation, global operational monitoring of the marine biosphere has been provided by several satellite ocean color missions, e.g., the Sea-viewing Wide Field-of-view Sensor (SeaWiFS) aboard the Orbview-2 satellite operating from 1997 to 2011(Hooker et al., 1992 ; Hooker and Esaias, 1993), Moderate Resolution Imaging Spectroradiometer (MODIS) aboard the EOS Terra and Aqua satellites launched in 1999 and 2002, respectively (Salomonson et al., 1989 ; Esaias et al., 1998), and Medium Resolution Imaging Spectrometer (MERIS) aboard the Envisat satellite launched in 2002 and retired in 2012(Rast et al., 1999). The deployment of these sensors marked a considerable advance on the CZCS. Atmospheric correction (AC) algorithms have been developed for these ocean color sensors for data retrieved over water (Gordon, 1978 ; Gordon and Castaño, 1987, 1989 ; Gordon and Wang, 1994a, b ; Antoine and Morel, 1998, 1999 ; Nobileau and Antoine, 2005). These AC algorithms remove the contribution of the atmosphere from the spectral radiance measured directly by the spaceborne sensor at the top of the atmosphere to derive simultaneously both the normalized water-leaving radiance L_wn (λ)(or remote sensing reflectance R_rs (λ)) and the aerosol optical thickness (AOT). It is well known that appropriate AC is an important procedure in ocean color remote sensing. The accuracy and effectiveness of the AC directly affect remotely retrieved ocean bio-optical products, such as the concentration of chlorophyll a, concentration of suspended particulate matter, concentration of colored dissolved organic matter, inherent optical properties of seawater (absorption coefficient a, and backscattering coefficient b_b), and these bio-optical products are retrieved from the R_rs (λ) or L_wn (λ). According to the NASA's technical report series on SeaWiFS, the product uncertainty of derived waterleaving radiances should be within 5%(Mueller and Austin, 1992). This value represents an objective stated prior to the launch of SeaWiFS and the actual uncertainties are larger than this for many marine regions. Thus, assessing the uncertainty of the parameters of the retrieved R_rs (or L_wn) and AOT and evaluating the AC algorithms are major requirements for the development of time series of ocean color products suitable for marine biogeochemistry and climate investigations (Zibordi et al., 2006 ; Antoine et al., 2008).

Many studies conducted in various ocean areas have evaluated the R_rs (or L_wn) and AOT of different ocean color satellite missions. For open oceans, Bailey and Werdell (2006) demonstrated that the median absolute percentage difference (MPD) for L_wn (λ) ranges between a low of 10.78% for the 412- nm band and a high of 16.64% for the 555-nm, when comparing SeaWiFS data with in situ measurements from the SeaWiFS Bio-optical Archive and Storage System (Werdell et al., 2003). Mélin et al.(2010) showed that the MPD for the AOT increases from 26% at 443 nm to 42% at 865 nm for SeaWiFS and from 39% to 51% at the same wavelengths for MODIS, using field observations from the Aerosol Robotic Network (AERONET)(Holben et al., 1998). It should be noted that although satellite-derived AOTs show high uncertainty, the accuracy of satellitederived L_wn (λ) in open oceans is generally considered satisfactory.

Coastal regions are generally complex in terms of both water constituents and atmospheric conditions. These can vary significantly with time and location and cause greater uncertainties for ocean color products in coastal waters than over the open ocean. Thus, considerable local validation has been performed in different coastal areas (Bulgarelli et al., 2003 ; Mélin et al., 2003, 2007, 2010 ; Pinkerton et al., 2003 ; Wang and Shi, 2005 ; Wang et al., 2005, 2007, 2009 ; Zibordi et al., 2006 ; Aiken et al., 2007 ; Ohde et al., 2007 ; Sørensen et al., 2007 ; Sun et al., 2009, 2010 ; Cui et al., 2010 ; He et al., 2011 ; Jamet et al., 2011). Because of the influence of the Yellow and Yangtze rivers and human activities, the highly turbid water and highly mixed atmosphere in the region of the East China Seas (ECS) have made ocean color retrievals very challenging (He et al., 2000 ; Li et al. 2003; Quinn et al., 2004). Currently, three Chinese ocean color sensors (the Chinese Ocean Color and Temperature Scanner aboard HY-1B and Medium Resolution Spectral Imagers aboard FY-3A and FY- 3B) are in operation; however, high-quality operational data products are still unavailable. Thus, the development of a suitable ocean color retrieval algorithm for the ECS is required urgently, and a systematic validation of the current operational ocean color data products in the ECS is a necessary part of this process.

To date, because of the shortage of in situ data, only limited validation activities have been conducted in the ECS region. For example, He et al.(2007) used in situ ECS data, collected between 2003 and 2007, to provide a preliminary evaluation of the performance of MERIS-derived R_rs. They found that the operational MERIS Level 2.0 R_rs data were generally overestimated with an average MPD of 71% over the bands of 412 to 665 nm. Based on field measurements acquired in the ECS in the spring and autumn of 2003, Sun et al.(2009) revealed that the MODIS-derived R_rs data were mostly underestimated with an averaged MPD of 25% over the bands of 412 to 667 nm; the poorest performance occurred at 412 nm (MPD=39%). Using in situ AOTs at 870 nm measured in the ECS from April 2006 to November 2007, Li and Chen (2010) showed that accuracy of MODIS-derived AOTs varied with locations and seasons. All these retrieval uncertainties have been preliminarily attributed to two potential reasons: unsuitable aerosol models and the near-infrared (NIR) water-leaving reflectance correction method (He et al., 2007 ; Sun et al., 2009 ; Li and Chen, 2010).

To explore how best to improve the AC algorithm for the Chinese in-orbit and scheduled ocean color sensors over the ECS, the focus of the present study was on the validation and comparison of MERIS R_rs and AOT data in the ECS, which were derived using two existing operational AC algorithms: the ESA MERIS processor MEGS 7.4.1 and NASA SeaDAS 6.1(S61). The MERIS data were used because the MERIS sensor has some outstanding features, such as high spectral and radiometric resolutions, dual spatial resolution, and programmable spectral bands. Furthermore, and its time-series data spans an entire decade. The in situ observations were obtained from field cruises, and Aerosol Robotic Network (AERONET) and Skyradiometer Network (SKYNET) coastal and island stations. Moreover, both algorithms are operational methods commonly used around the world to obtain ocean color data products. Different aerosol models (traditional and new) and NIR correction methods (considering suspended particulate matter scattering and not) were implemented in both algorithms and the results compared to provide better understanding of how the current algorithms could be improved for use in the ECS. The validation results of the present study provide important information not only for algorithm developers but also for a wide range of users of ocean color data.

2 DATA AND ATMOSPHERIC CORRECTION ALGORITHMS 2.1 In situ data

The in situ observations used in this study comprised water radiometric measurements and atmospheric aerosol optical properties. Table 1 provides information relating to these data and to the corresponding responsible principal investigators. The water radiometric measurements of R_rs (λ) were obtained from the Ocean Remote Sensing Institute ocean optics data base (ORSIO² DB) collected during 11 field experiments (No. 18-28) in the ECS (including the Yellow and the East Seas) between 2002 and 2008(covering four seasons). The aerosol optical properties of the AOT (denoted as τ_a (λ)), Ångström exponent (denoted as α), and single scattering albedo (SSA; denoted as ω_a) were obtained from 6 field cruises (No. 22-27) between 2006 and 2008(ORSIO² DB) and 13 AERONET (No. 1-13) and 4 SKYNET (No. 14-17) stations.

The R_rs (λ) data were collected using a hyperspectral tethered spectral radiometric buoy (HyperTSRB; Satlantic, Inc., Halifax, Canada), which was deployed at least 30 m away from the ship for about 10-30 min at each field site, and the data were subsequently processed using Satlantic's analytical software (ProSoft 6.3). The HyperTSRB measured upwelling radiance 0.65 m below the sea surface, L_u (λ, 0.65), and downwelling irradiance just above the sea surface, E_d (λ, 0⁺) at 256 channels between the wavelengths of 400 and 800 nm (about 3.3-nm spectral resolution). The value of L_u (λ, 0.65) was extrapolated to just below the sea surface and just above the sea surface using the diffuse attenuation coefficient for upwelling radiance K_L (λ) and the n-squared law for radiance:

where n ≈1.345(Mobley, 1994) is the real index of refraction of seawater and ρ ≈0.021(Preisendorfer and Mobley, 1986 ; Mobley, 1994) is the Fresnel reflectance of the surface. The value of K_L (λ) was computed based on the empirical algorithms proposed by Austin and Petzold (1984) and Morel (1988). Then, the value of R_rs (λ) was computed using L_w (λ, 0⁺) and E_d (λ, 0⁺):

The AOT τ_a (λ) is an index used to represent the integral of light extinction by aerosols within an atmospheric column. Ångström (1929) proposed that the optical thickness of aerosols in the atmosphere could be determined using the following formula:

where β is the atmospheric turbidity parameter (τ_a at λ = 1.0 μm) and α is the Ångström exponent, which is inversely related to the average size of the particles of the aerosol; the smaller the particles, the larger the exponent. According to Eq.3, the relationship of the AOT at the two different wavelengths of λ₁ and λ₂ can be described as:

The aerosol SSA ω_a is a measure of the effectiveness of scattering relative to extinction (the sum of scattering and absorption) for light encountering the atmospheric aerosol particles. In radiative transfer studies, the parameter ω_a is the ratio of scattering aerosol optical depth to the total aerosol optical depth of the atmosphere. It is a dimensionless quantity (range 0-1) for which a value of unity implies that all aerosol particle extinction is due to scattering and a value of zero implies that all extinction is due to absorption. The atmospheric aerosol optical properties used in the present study were obtained from different sources, including AERONET and SKYNET stations and cruises. A brief introduction to the data collection and processing methods is given in the following.

The AERONET program (Holben et al., 1998) is a federation of ground-based remote sensing aerosol networks, established by NASA and PHOTONS (PHOtométrie pour le Traitement Opérationnel de Normalisation Satellitaire), which provides qualityassured aerosol optical properties with which to assess and validate satellite retrievals. The AERONET data consist of optical thickness and microphysical parameters of aerosols in the atmosphere, derived from measurements by CE-318 sun photometers (Cimel Electronique, Paris, France; Dubovik et al., 2000, 2002, 2006 ; Dubovik and King, 2000 ; Smirnov et al., 2000 ;O'Neill et al., 2003). The AOT data are derived from measurements of direct solar irradianceby the sun photometers, based on the Beer-Bouguer law, in eight spectral bands centered at 340, 380, 440, 500, 670, 870, 940, and 1 020 nm. The microphysical properties, such as SSA, are derived from diffuse sky radiance measurements, generally in four bands at 440, 670, 870, and 1 020 nm. Three data quality levels of AOT data are provided by AERONET: Level 1.0(unscreened), Level 1.5(cloud-screened), and Level 2.0(cloud-screened and quality-assured). The microphysical properties are derived from these levels with additional quality checks. Level 2.0 data were used in this work.

SKYNET is an observational network established in the last decade for the purposes of understanding the interaction between aerosols, cloud, and radiation in the atmosphere of East Asia. The observation sites, which are equipped with several automatic instruments including sky radiometers, pyranometers, and pyrgeometers, are located mainly in East Asia. The levels of diffuse and direct radiation are both measured either by POM-01 sky radiometers in seven channels centered at 315, 400, 500, 675, 870, 940, and 1 020 nm and/or POM-02 sky radiometers in 11 channels centered at 315, 340, 380, 400, 500, 670, 870, 940, 1 020, 1 600, and 2 200 nm (PREDE Co. Ltd., Japan) at each site. Atmospheric parameters such as the AOT, SSA, volume size distribution, and refractive index of the aerosol particles are retrieved using SKYRAD. pack version 4.2(Nakajima et al., 1996). To separate cloud-free data from cloud-affected data, a cloudscreening algorithm (Khatri and Takamura, 2009) is applied, which involves three tests based on algorithms proposed by Long and Ackerman (2000), Kaufman et al.(2006), and Dubovik et al.(2000). As no further quality assurance criteria are applied to SKYNET data, only the cloud-screened SKYNET AOT data were used in this study and the SSA data were rejected.

The aerosol optical parameters from six cruises of the ORSIO² DB were all measured using a shipborne sky radiometer (POM-01 MKII, Prede Co. Ltd., Tokyo, Japan). The instrument took measurements of the direct solar irradiance and the diffuse sky radiance distribution at seven (POM-01) wavelengths. Analysis of these data was performed using SKYRAD.pack version 4.2, from which 813 records of 57 days were obtained. To eliminate the effects of cloud, a cloudscreening algorithm was implemented, which included the spectral variability test proposed by Kaufman et al.(2006) and some statistical analysis tests of Smirnov et al.(2000). After the cloudscreening procedure, 620 records were considered as clear-sky data. The microphysical properties of SSA were further screened with additional quality checks (solar zenith angle ＞50°, scattering angles >80°, sky error <10%, τ_a (440) no less than 0.4), which were determined according to the AERONET screening criteria. Only 65 SSA records passed the screening.

2.2 MERIS data

The MERIS sensor has medium spectral and high radiometric resolutions, and apart from the high spectral imagers, it is recognized as one of the optimal ocean color sensors and therefore, MERIS data were used in this study.

The MERIS Level 1.0b and Level 2.0 data of Reduced Resolution (MER_RR_1P and MER_ RR_2P) with a ground sampling distance of 1.2 km at nadir (May 2002 to October 2009) used in this study were obtained from the ESA. For the validation analysis performed in this work, the MERIS R_rs and aerosol optical properties processed by the ESA MEGS 7.4.1 and NASA S61 were considered.

The spectral normalized surface reflectance ρ_wn (λ), τ_a (865), and Ångström exponent α were routinely provided by the ESA MERIS Level 2.0 data products, which were processed using the atmospheric correction procedure over water of the ESA MEGS 7.4.1. The ρ_wn (λ) was converted to R_rs (λ) by

and with the τ_a (865) and Ångström exponent α, the τ_a (λ) was computed using the power law relationship in Eq.4. The R_rs (λ) and τ_a (λ) were also derived from MERIS Level 1.0b data through the atmospheric correction algorithms of the NASA S61 software. The Ångström exponent α was computed from the logtransformed τ_a (λ) and associated wavelengths by means of a least squares fit to Eq.3.

2.3 Atmospheric correction algorithms

Both existing operational AC algorithms (ESA MERIS processor MEGS 7.4.1 and NASA S61) were developed initially based on the scheme of Gordon and Wang (1994a, b) for Case I water. For the oceanatmosphere system, the radiance L_t (λ) received by an ocean color sensor at the top of the atmosphere can be written as (Gordon and Wang, 1994a):

where t is diffuse transmittance and T is direct transmittance. The goal of the AC is to retrieve an accurate estimate of the spectral water-leaving radiance L_w (or R_rs) by removing the whitecap radiance L_wc, sun glint radiance L_g, Rayleigh radiance L_r, combined radiance of aerosols, and the multiple interaction term between molecules and aerosols (referred to as the aerosol contribution) L_a + L_ra from the radiance at the satellite altitude L_t. In brief, the general approach is to use two NIR bands for which L_w is known to make an assessment of the aerosol contribution and to extrapolate these into the visible spectrum (Gordon and Wang, 1994a).

2.3.1 NASA operational AC algorithm

In the scheme of Gordon and Wang (1994a), 12 candidate aerosol models from Shettle and Fenn (1979) were used to build lookup tables for computing the aerosol contributions. The aerosol spectral optical properties, as well as the coefficients describing the relationship between multiple and single scattering reflectance for different aerosol models were precalculated and stored in these LUTs. As the single scattering reflectance ρ_as is related directly to the aerosol phase function, SSA, and optical thickness, the value of ρ_a + ρ_ra can be computed for each model using the LUTs.

In the AC procedure of NASA's S61(version 6.1) used in this work, a suite of 80 newly developed candidate aerosol models were used, while 12 aerosol models from Shettle and Fenn (1979) were used in earlier versions before SeaDAS 6.x. The new models were developed based on observed size distributions and SSAs from maritime AERONET sites, including both open-ocean (eight sites) and coastal (three sites in the Chesapeake Bay region) environments. Table 2 shows the optical properties of the new models. In Table 2, the Ångström exponent is the mean value for a given aerosol model and the SSA is at the wavelengths of 440 and 865 nm.

In optically complex coastal waters, where the water-leaving reflectance in the NIR is non-negligible because of backscattering by phytoplankton, detritus, and suspended sediments, the black assumption used for open-ocean water no longer holds. Thus, some methods that compensate for the NIR water-leaving radiance have been developed (Moore et al., 1999 ; Hu et al., 2000 ; Ruddick et al., 2000 ; Siegel et al., 2000 ; Stumpf et al., 2003 ; Lavender et al., 2005 ; Kuchinke et al., 2009 ; Bailey et al., 2010). In S61, the NIR water-leaving reflectance is estimated using an NIR iterative-correction scheme (Stumpf et al., 2003 ; Bailey et al., 2010) whenever the black assumption is not valid. An empirical spectral dependence of the backscattering (Lee et al., 2010) is used to derive the NIR water-leaving contribution.

2.3.2 ESA operational AC algorithm

Details of the ESA operational AC scheme used in the MERIS MEGS 7.4.1 are given in Antoine and Morel (1998, 1999). The main differences between this scheme and that of Gordon and Wang (1994a, b), as given by Antoine and Morel (1998), are as follows:(1) in this scheme, there is no need to rely on or to produce the single scattering reflectance ρ_as with a view to identify the aerosol type;(2) the transfer from one wavelength to another one is directly performed by means of the changes in τ_a (λ), in contrast with the ρ_as (λ) in the method of Gordon and Wang (1994a). In addition, besides the maritime, coastal, and rural aerosol models from Shettle and Fenn (1979), which were similar to those used by Gordon and Wang (1994a), the aerosol assemblages also included aerosol models of the blue, dust, continental, and H₂ SO₄ from Santer et al.(1999), Moulin et al.(2001), and WCRP (1986). The optical properties of the aerosol models used in the ESA AC scheme are shown in Table 3.

In turbid water, the bright-water correction algorithm proposed by Moore et al.(1999) was implemented in the ESA operational AC scheme. In the algorithm, the water-leaving reflectances at the three bands of 705, 775, and 865 nm were used by a model of the absorption and backscattering properties of sediment-dominated waters, with the suspended sediment concentration as a free parameter.

3 PROCEDURE 3.1 Match-up procedures

The match-up analyses were performed separately for the MERIS data obtained by the two processors (ESA MEGS 7.4.1 and NASA S61), and the exclusion criteria in the match-up procedure were formulated based on the work of Bailey and Werdell (2006), Mélin et al.(2010), Bailey and Wang (2001), Wang et al.(2000), and Bailey et al.(2000), which removed invalid or redundant data from further consideration. Figure 1 shows the flowchart of the match-up procedure that comprised six steps.

The MERIS remote sensing reflectance and AOT from the ESA Level 2.0 products and field observations were first screened using this match-up procedure. Then, the MERIS Level 1.0b data acquired on days with successful match-ups were downloaded and processed by S61. The same match-up scheme was performed again between the MERIS data processed by S61 and the in situ data. Note that different temporal and spatial windows were adopted in the match-up schemes of R_rs and τ_a. For R_rs, a temporal window of ±3 h and a spatial window of 3×3 pixels were used, whereas for τ_a, a temporal window of ±1 h and a spatial window of 25×25 pixels were selected.

In step 3 of the match-up procedure shown in Fig. 1, the flags “CLOUD”, “WATER”(not land), “PCD_19”(uncertain aerosol type and optical thickness, corresponding to a failure of the AC, were was used in the match-up procedure for τ_a), “PCD_1_13”(uncertain normalized surface reflectance, corresponding to a failure of the AC, was used in the match-up procedure for R_rs), and “LOW_SUN”(high solar zenith angle) were considered in the match-up procedure for the ESA MERIS Level 2.0 data. Pixels flagged as “HIGH_GLINT”(glint), “CLOUD”(cloud), “CLDICE”(probable cloud or ice contamination), “LAND”(land), “ATMFAIL”(AC failure), “HISATZEN”(sensor view zenith angle>60°), and “HISOLZEN”(solar zenith angle>70°) were excluded from the MERIS data processing by S61.

In step 6, the CoV is the spatial coefficient of variation, which is defined as the ratio of the standard deviation and average computed over the selected pixels. In the match-up procedures for R_rs, only the screening criterion of MERIS CoV of τ_a (865) was checked, while for τ_a, the field τ_a (870) and MERIS τ_a (865) were both considered.

Finally, 9 and 8 R_rs match-ups were obtained and 143 and 77 with τ_a match-ups were available for the ESA MEGS 7.4.1 and NASA S61 processors, respectively. The numbers of match-ups for each cruise or site are given in Tables 4 and 5, and the distribution of the match-up sites is displayed in Fig. 2.

3.2 Statistical indicators

The statistical indicators including the mean ratio (MR), mean percentage difference (MPD), root mean square difference (RMS), and coefficient of determination (R²) were used in the present study, and these are defined as follows:

where x_i^sat is the i^th satellite-derived value, x _i ^mea is the i^th field observed value, and N is the number of match-up data. As the ranges of R_rs and τ_a vary considerably at different bands and various locations in the ECS, the relative RMS, defined in Eq.9, was used in this work. A least squares fit was also adjusted within the match-up data with slope and intercept.

4 RESULT

The performances of the ocean color operational AC algorithms used by NASA and ESA were evaluated separately by validating the spectral remote sensing refl ectance R_rs (λ), aerosol optical thickness τ_a (λ), Ångström exponent α, and SSA retrieved from the MERIS data with in situ measurements.

4.1 Evaluation of MERIS R_rs

Figure 3 shows the comparison of spectral remote sensing reflectance at the match-up sites, and the MERIS data derived by the ESA MEGS 7.4.1 processor and NASA S61 software are respectively denoted as ESA and S61. At Station 9, only the matchup for ESA was derived and it clearly shows an overestimation of R_rs (λ). This overestimation could have been caused by cloud-affected pixels because the match-up for S61 failed the screening procedure outlined in Fig. 1. Thus, Station 9 was considered an outlier in the following analysis.

Scatter plots of R_rs match-up points (including Station 9) for ESA and S61 are shown in Fig. 4, and the statistical parameters (excluding Station 9) are given in Table 6. Band ratio algorithms are commonly used in operational bio-optical algorithms and thus, the comparison of ratios R_rs (443)/ R_rs (560), R_rs (490)/ R_rs (560), and R_rs (510)/ R_rs (560), as well as the maxima of these three ratios, are also presented in Fig. 4 and Table 6.

The values of MERIS R_rs (λ) are overestimated at short wavelengths (MR ＞1) and underestimated at long wavelengths (MR<1) by both the ESA and the S61 algorithms, as shown by Fig. 4 and the mean ratios in Table 6. The largest overestimation occurs at 413 nm with MRs of 1.59 and 1.56 and the largest underestimation occurs at 620 nm with MRs of 0.73 and 0.66 for ESA and S61, respectively. Accordingly, the MPDs and RMSs of MERIS R_rs (λ) derived by ESA and S61 vary with wavelength: the largest error occurs in the blue bands of 413 and 443 nm (particularly the 413-nm band), followed by the red bands of 620 and 665 nm, and finally, the greenyellow bands from 490 to 560 nm. The maximum MPDs reach 86.83% and 71.96 at 413 nm and the RMSs are 1.131 2 and 0.835 8 for ESA and S61, respectively. The coefficient of determination (R²) in the blue bands (413 and 443 nm) is also lower than in the green and red bands. Thus, the MERIS R_rs (λ) derived with both ESA and S61 shows higher accuracy at bands from 490 to 560 nm than at other bands, which is consistent with the results of Antoine et al.(2008), Cui et al.(2010), and Zibordi et al.(2006).

As the MRs in Table 6 decrease from 443 to 560 nm for both ESA and S61, the three individual band ratios, as well as the maximum band ratio, are all overestimated, and the errors at the three individual bands decrease with increasing wavelength, as shown in Fig. 4(h -k) and Table 6. The MPDs of R_rs (443)/ R_rs (560) reach 32.84% and 35.12% for ESA and S61, respectively; however, the MPDs of the maximum band ratio are 12.69% and 12.97%. The reason for this is that the value of R_rs (443) is usually smaller than R_rs (490) and R_rs (510) at individual matchup sites in the present study, as shown in Fig. 3, and the maximum band ratio is dominated by R_rs (490)/ R_rs (560) and R_rs (510)/ R_rs (560), not by R_rs (443)/ R_rs (560).

It should be noted that 8 common R_rs match-ups were obtained using 291 in situ data under the strict match-up procedures. Moreover, they are located in different water areas, including very turbid water, turbid water, and clear water, as can be seen in Figs. 3 -5. The error statistics calculated from only N =8 match-ups may lack sufficient confidence and therefore, the match-ups are used to demonstrate the different cases and the error statistics are given simply for reference.

The water type (turbid or clear) can be reflected by the spectral shapes and amplitudes of R_rs. Taking R_rs (560) as the turbidity index, the relationship of the maximum band ratio with turbidity is given in Fig. 5. Although R_rs (560) is a simple turbidity index, the match-up stations are arranged reasonably according to increasing turbidity (R_rs (560)), except for Station 8. Actually, the water at Station 8 is clearer than at the other stations. It can be seen from Fig. 5 that the max band ratios decrease with increasing turbidity. Thus, it can also be concluded from Fig. 4k that the overestimation errors of the maximum band ratios for ESA and S61 decrease with increasing turbidity. This suggests that the MERIS chlorophyll concentrations derived in clear water can have lower accuracy than in turbid waters in the ECS region.

4.2 Evaluation of MERIS τ_a α

Using the Ångström exponent α and AOT at 870 nm as the indicator of aerosol type and atmospheric turbidity, data of 143 and 77 match-ups for ESA and S61, respectively, are plotted in Fig. 6. To show the correlation between the Ångström exponent and the AOT at 870 nm, these match-up points are arranged in ascending order of in situ τ_a (870) and α.

Observations can be made in relation to Fig. 6a, b:(1) the ESA MERIS τ_a (870) roughly follows the ascending in situ τ_a (870) in the large range 0-0.4, and it tends to be overestimated compared with the in situ data;(2) few match-ups are available in the range 0.2-0.4 for S61 MERIS data; however, a rough trend of S61 MERIS τ_a (870) following the ascending in situ data can also be observed in the range 0-0.2. The S61 MERIS τ_a (870) are obviously overestimated at low τ_a (870) values, but slightly underestimated at high τ_a (870) values;(3) the ESA MERIS Ångström exponent data are generally underestimated compared with in situ data. The underestimation error is large at low τ_a (870) values and it decreases with the increase of τ_a (870);(4) the S61 MERIS Ångström exponent data are also underestimated at low τ_a (870) values, but overestimated at high τ_a (870) values. This is just the opposite to the error analysis for S61 MERIS τ_a (870).

From Fig. 6c and d, it can be seen that the overestimation error of ESA MERIS Ångström exponent increases with the increase of in situ α, while the S61 MERIS α is first overestimated at low α values, and then underestimated at high α values. In addition, a slight trend can be observed from the in situ data (blue circles and triangles) in Fig. 6, whereby small τ_a (870) values occur with large α values, and vice versa. Furthermore, it must be highlighted that this conclusion might be limited by the locations of the match-up points (i.e., they are all over the ECS).

Different match-up points are obtained for ESA (142) and S61(77) MERIS data. In order to compare the performances of these two operational AC algorithms in a fair and quantitative way, 77 common match-up data are used. Scatter plots of MERIS α and τ_a at four bands (440, 500, 675, and 870 nm) versus in situ data are displayed in Fig. 7. Note that three outliers of ESA MERIS data occur in Fig. 7(a -d), which could be related to the effect of clouds. Considering that outliers can have considerable effect on statistical results and thus, lead to misleading interpretations, the 3 outliers were removed and only the 74 common match-up points were used in the final statistical analysis, the results of which are summarized in Table 7.

Figure 7e and Table 7 show that the ESA MERIS α is obviously underestimated with an MR of 0.63, whereas the S61 MERIS α agrees well with the in situ data with a high MR of 0.99. However, the retrievals of MERIS α obtained with the ESA and the S61 algorithms both exhibit low accuracy, with MPDs of 35.15% and 40.56% and RMSs of 0.865 8 and 0.484 9, respectively.

The MERIS τ_a (λ) is overestimated by both the ESA and the S61 algorithms with the 4-band averaged MR of 1.39 and 1.38, respectively. For the ESA algorithm, the MR, MPD, and RMS decrease with decreasing wavelength. For the S61 algorithm, the MR, MPD, and RMS also decrease from 870 to 500 nm, but then they increase from 500 to 440 nm. In the blue-green bands (such as 440 and 500 nm), the ESA algorithm provides a more accurate estimation of τ_a (λ) than S61, but in the red-NIR bands (such as 675 and 870 nm), the S61 algorithm is the more accurate of the two. This indicates that for the ESA algorithm, the underestimated α and overestimated τ_a (870) just give a more accurate estimation of τ_a (λ) at blue-green wavelengths by extrapolation from the NIR spectral region. For the S61 algorithm, the MERIS τ_a (λ) at the blue-green wavelengths are apparently overestimated because of the lower bias of α and the overestimated τ_a (870). An overestimation of τ_a (λ) might lead to an underestimation of the water-leaving radiance in the AC process.

The Ångström exponent α can be used to indicate whether an appropriate aerosol model set, in terms of the particle size distribution, is being used for the MERIS data processing. As shown in Table 7, the S61 MERIS α data show better agreement with the in situ data (MR=0.99) than the ESA MERIS data (MR=0.63). Thus, the aerosol size distributions defined for the aerosol models used in the S61 processing are more reasonable than in the ESA algorithms. However, both the ESA and the S61 MERIS α data were estimated with low accuracy (i.e., high MPD and RMS). This could be partly due to the lack of specific aerosol models over the ECS, which can be seen from Fig. 6. In Fig. 6b and d, the disagreement between the S61 MERIS α and in situ measurements occurs at low τ_a (870) and high α, or at high τ_a (870) and low α. The low τ_a (870) and high α are likely associated with aerosols dominated by fine particles and thus, they could be affected by urban or biomass burning aerosols. High τ_a (870) and low α are likely associated with dust aerosols. In Fig. 6a and c, the difference between the ESA MERIS α and in situ measurements increases with the increase of in situ α and decrease of in situ τ_a (870), and thus aerosols dominated by the lack of fine particles could be the main problem in this instance.

There are several factors that can cause the estimation error of MERIS τ_a (870), e.g., the aerosol models, NIR ocean contribution correction methods, and radiative transfer modeling (LUTS). The overestimation error of MERIS τ_a (870) derived by the S61 software is lower than that derived by the ESA processor. This could be caused by the different aerosol models used, i.e., as defined by Ahmad et al.(2010) for S61 and by Shettle and Fenn (1979) for ESA. It has already been highlighted by Wang et al.(2005) and Ahmad et al.(2010) that the Shettle and Fenn (1979) models overestimate the AOT over oceans. The radiative transfer modeling used to generate the aerosol LUTS can also introduce some errors into the estimation of MERIS τ_a (870). The uncertainty resulting from the aerosol models and radiative transfer modeling is usually within 15% for open-ocean waters, according to results presented by the International Ocean Color Coordinating Group (IOCCG, 2010). However, in the present study, the ESA and S61 MERIS τ_a (870) is obviously overestimated with MR>1.42, MPD>50%, and RMS>0.80. This large overestimation is very likely because the ESA and S61 NIR correction methods are unsuitable for the ECS. As is known, the optical properties of the particulate backscattering are complex in the ECS, and very significant ocean NIR contributions are evident in the China east coastal region throughout the year (Wang et al., 2007 ; Zhang et al., 2010). It must be stated that an incorrect estimation of the contribution of the NIR waterleaving reflectance can also lead to the selection of inappropriate aerosol types in the AC procedure, which could reduce the retrieval accuracy of the Ångström exponent α. Note that even though the AOTs are overestimated by the AC algorithms, the ocean color products (remote sensing reflectance) can still be derived accurately (IOCCG, 2010). Similar results can be found in this study. The MERIS τ_a (λ) are significantly overestimated from 500 to 675 nm by the ESA and the S61 algorithms (MR>1.25, MPD>34%, RMS>0.69), while the MERIS R_rs (λ) show good agreement with the in situ data from 490 to 560 nm (MR about 1.07, MPD within 14%, RMS within 0.17).

It was also found in this study that the MERIS τ_a (λ) are overestimated less by the ESA and S61 algorithms at the blue band (440 nm) with MRs of 1.19 and 1.36, MPDs of 32% and 39%, and RMSs of 0.69 and 0.71, respectively, compared with the longer wavelengths from 500 to 870 nm. However, the ESA and S61 MERIS R_rs (λ) both show larger overestimation errors in the blue bands (especially 413 nm) than in the longer wavelengths, with MRs of up to 1.59 and 1.56, MPDs of up to 87% and 72%, and RMSs of up to 1.13 and 0.84, respectively. Similar overestimation results have also been reported for ESA MERIS R_rs (λ) data by Zibordi et al.(2006) for a coastal site in the northern Adriatic Sea, by Antoine et al.(2008) for a deep clear ocean site in the northwestern Mediterranean Sea, and by Cui et al.(2010) for turbid coastal waters in the Chinese Bo Sea. This overestimation is most likely due to the lack of on-orbit vicarious calibration for MERIS, something that has been employed for both SeaWiFS and MODIS (IOCCG, 2010).

4.3 Evaluation of MERIS SSA

SSAs are not operational products for ocean color missions. Thus, the SSA properties of the aerosol models used by the ESA and NASA S61 algorithms are compared with in situ data. The Ångström exponent α and SSA are two important parameters for aerosols and thus, we evaluated the SSA properties using α-SSA density scatter plots.

Suites of aerosol models, which are assumed to represent the realistic atmospheric aerosol conditions over oceans, are usually used in current operational AC algorithms, as mentioned in Section 2.3. The aerosol model is defined by its particle size distribution and refractive indices (both real and imaginary parts). The Ångström exponent α is related to the particle size distribution, and the SSA is dominated by the refractive indices. Therefore, the Ångström exponent and SSA from ground-based measurements can be used to indicate whether the aerosol models used in the AC algorithms are suitable. Satellite-field matchup results of the Ångström exponent α have been given in the previous section. As the SSA parameters are not provided in the ESA MERIS Level 2.0 data, the SSA properties of aerosol models used by the ESA and S61 algorithms are compared directly with the field observations in this section, and the match-up procedure of Section 2.4 is not performed for the SSA parameters.

The colored background in (a, b, d, and e) is the density scatter plot of in situ data. Black triangles indicate the ESA aerosol models and the red crosses represent the S61 aerosol models. Histograms of in situ SSA data are also given in (c and f).

Figure 8 provides comparisons of the SSA and Ångström exponent α from the ESA and S61 aerosol models with those of ground-based measurements. The ground-based measurements include the AERONET and ship-borne data from April 2006 to July 2008. It must be highlighted that the range of in situ α in Fig. 8 might not represent the natural variety of α data in the ECS because these data were acquired under limited temporal and spatial ranges. As shown by the in situ data processing of Section 2.1, additional quality checks were applied on the SSA data and thus, the numbers of in situ quality-assured SSA data are reduced in comparison with the Ångström exponent data. In Fig. 8, only the in situ α corresponding to the SSAs are used. Considering that the SSA changes with wavelength and that it is related to the Ångström exponent α, the density scatter plots of in situ α versus SSA (440) and SSA (870) are used in Fig. 8. The α and SSA properties of the ESA and S61 aerosol models are overplotted on the density scatter plots. The histograms of AERONET and ship-borne SSA data are also shown separately in Fig. 8c and f.

It is apparent from the density plots of Fig. 8a, b, d, and e that the SSA ranges for both the ESA and the S61 aerosol models are not large enough to cover the entire range of in situ measurements. Except for the dust aerosol models in the ESA algorithms (α<0.5, SSA at 440 nm<0.93), the aerosol models of the ESA and S61 are all non- or weakly absorbing, and most of their SSA values are concentrated within the range 0.95-1.0. However, strongly absorbing aerosols play an important role in the ECS, especially in the coastal seas, which can be seen from the density scatter plots and histograms of the AERONET and ship-borne measurements in Fig. 8a, b, d, and e. Note that the in situ aerosol properties collected at the AERONET coastal and island sites shown in Fig. 8a -c are different from the ship-borne data in Fig. 8d -f. As the former were greatly affected by land activities, the SSA values of the former are smaller than the latter and thus, the absorption of the former is stronger than the latter.

As for the Ångström exponent α, the properties of the ESA and S61 aerosol models can cover the in situ data; however, the α distribution of the S61 aerosol models (red crosses in Fig. 8a, b, d, and e is more reasonable than that of the ESA (black triangles). The S61 aerosol model values are distributed evenly over the range of [0, 2.2], while most of the ESA aerosol model values are located in the range of [0, 1.5].

As shown in Fig. 8, absorbing and strongly absorbing aerosols (with low SSA at 440 nm<0.95) exist over the ECS, especially in the coastal regions (most AERONET sites used in this study are located near the coast). Therefore, it is necessary to consider how best to detect and predict the contribution of these absorbing aerosols over the ECS, although no quantitative analysis of the effects of the absorbing aerosols on the R_rs estimations is provided in this study. Failure to correct for their effects could introduce an artifact in the retrieved water-leaving reflectance, especially at short wavelengths (Fukushima and Toratani, 1997 ; Gordon et al., 1997 ; Moulin et al., 2001 ; Li et al., 2003 ; Ransibrahmanakul and Stumpf, 2006).

5 CONCLUSION

The objective of this study was to explore how best to improve the AC scheme to acquire high-quality operational data products for Chinese in-orbit and scheduled ocean color sensors. This was achieved using MERIS data products of spectral remote sensing reflectance (R_rs), aerosol optical thickness (τ_a), and Ångström exponent (α) derived by two existing operational AC algorithms (ESA MEGS 7.4.1 and NASA S61) over the ECS between May 2002 and October 2009, which were compared with in situ measurements. The in situ R_rs data comprised shipborne measurements from the ECS, and the in situ aerosol optical properties data included measurements from coastal and island stations of AERONET and SKYNET around the ECS, and ship-borne data collected over the ECS. The validation exercise yielded 8 R_rs and 74 τ_a / α common match-ups for both the ESA and the NASA S61 algorithms through a strict match-up procedure.

The statistical analysis of satellite-derived and in situ R_rs match-ups exhibited comparable accuracy for the ESA and S61 algorithms with MPDs varying from 9% to 13% in the 490-560 nm spectral range. The values of satellite-derived R_rs (λ) were found to be overestimated at short wavelengths and underestimated at long wavelengths. Significantly high values of 87% and 72% were found at 413 nm for the ESA and NASA operational algorithms, respectively.

The results from the analysis of satellite-derived τ_a and α showed low accuracy for both the ESA and S61 algorithms. The MERIS τ_a (λ) was overestimated by both the ESA and S61 algorithms with averaged MPDs of 46% and 43% over four bands. The estimation errors of τ_a from 675 to 870 nm were higher than from 440 to 500 nm. The ESA algorithm provided the more accurate estimation of τ_a (λ) at blue-green bands (440 and 500 nm), while the S61 algorithm was better at red-NIR bands (675 and 870 nm). The ESA MERIS α was obviously underestimated with an MR of 0.63, whereas the S61 MERIS α agreed well with the in situ data with a high MR of 0.99. This suggests that compared with the ESA algorithm, the aerosol model set, in terms of the particle size distribution, used in the S61 atmospheric correction is more appropriate to the ECS region.

To clarify whether the NASA operational algorithm performs better than the ESA algorithm for the retrieval of α, scatter plots of α-SSA density were prepared. A qualitative comparison of the SSAs from the ESA and NASA S61 aerosol models with ground measurements was performed. This revealed that strongly absorbing aerosols exist over the ECS, but that the SSA ranges of both the ESA and the S61 aerosol models are not large enough to cover the in situ measurements.

In summary, to establish an appropriate operational AC algorithm over the ECS, in terms of the aerosol models, the aerosol particulate size distribution of fine and coarse particles should be emphasized and the proportion of strongly absorbing aerosols should be increased. In terms of the NIR correction of the water-leaving radiance, a suitable spectral model for the particulate backscattering probability in the ECS is required.

The AC algorithm, as a key and complex component in ocean color remote sensing, has remained an important issue to be resolved since the 1970s. In the special case of the ECS, the aerosols and water constituents have a complex composition that apparently changes with season and sea area. Thus, a good understanding of the aerosol optical properties is a prerequisite for the development of an accurate AC algorithm for the ECS. However, high-quality field aerosol properties over the sea are usually even more difficult to acquire than bio-optical and radiometric quantities. This lack of in situ observations has slowed the process of obtaining a satisfactory AC algorithm and thus, has prevented the acquisition of high-quality operational data products over the ECS. Therefore, a comprehensive in situ dataset covering all seasons and sea conditions is required urgently if an effective ocean color operational algorithm for the ECS is to be developed (He et al., 2014).

6 ACKNOWLEDGMENT

Thanks to the principal investigators and their staffs for establishing and maintaining the Ocean Remote Sensing Institute Ocean Optics Data Base (ORSIO² DB), and the AERONET and SKYNET stations used in this work. Thanks to the ESA Dragon project ID2566 and ESA EO project ID6760 for providing MERIS data.