1. Introduction

In Ocean Color (OC) applications for deep and shallow waters the accuracy of retrievals of chlorophyll-a concentrations, water optical properties and depth depends on the quality of the estimated remote sensing reflectance,_{$R_{rs}(\lambda )$}. One of the significant uncertainties in this estimation is associated with the characterization of the ocean surface and removal of the sky component reflected from it, especially in windy conditions.

For satellite applications, the reflectance coefficient of skylight from the ocean surface is included in the atmospheric correction algorithm, specifically in the calculation of the Rayleigh and aerosol components [1–3]. This coefficient is the one of the main parameters in the data processing for the NASA Aerosol Robotic Network for Ocean Color (AERONET-OC) which determines normalized water leaving radiance and _{$R_{rs}\left(\lambda \right)$}using SeaPRISM instruments installed on ocean platforms. For shipborne measurements, the reflectance coefficient is pre-calculated for the specific geometries, which are recommended for such observations [4,5]. In all these cases the reflectance coefficient for wind roughened surfaces are estimated based on Cox-Munk [6] statistics, where characteristics of wave slopes were measured and the reflected radiance was parameterized as a function of the wind speed. This approach is included in multiple scalar [7] and vector [1,8–12] radiative transfer models. It was emphasized that polarization effects should be taken into account for the accurate estimation of the impact of windy surface and reflectance coefficient [1,10,13]. Recent studies show that the reflectance coefficient can depend on the wavelength [14] and on aerosol optical thickness (AOT) [15–17], and such dependencies are related to the wind speed and can possibly affect the quality of the near surface _{$R_{rs}$} measurements and atmospheric correction [18]. Using vector radiative transfer simulations (VRT) Zhang et al [17] analyzed the contribution of the Sun glint to the total reflectance from the surface emphasizing its spectral features. Hyperspectral characteristics of the ocean surface are of great interest due to the development of novel hyperspectral atmospheric correction techniques for the Hyperspectral Imager for Coastal Ocean (HICO) [19] in the preparation for the NASA Plankton, Aerosol, Cloud, ocean Ecosystem mission (PACE) [20]. The goal of this paper is to extend the scope of VRT simulations for a broad range of observational geometries and conditions, compare them with available experimental data and evaluate the possible impact on shipborne measurements, operation of AERONET-OC and satellite observations.

2. Theory and modeling

In the field of OC, the main parameter used for the retrieval of water properties is the remote sensing reflectance_{$R_{rs}(\lambda )$}, defined as the ratio of the water leaving radiance _{$L_w({\theta }_v,{\phi }_v,\lambda )$}to the down-welling irradiance _{$E_d(\lambda )$} [21]:

$R_{rs}({\theta }_v,{\phi }_v,\lambda )=L_w({\theta }_v,{\phi }_v,\lambda )/E_d(\lambda ),$ (1)where _{${\theta }_v,{\phi }_v$}are the viewing zenith (VZA) and azimuth angles respectively and_$\lambda$is the wavelength. By convention, the azimuth angle _{${\phi }_v$} is equal to 0° when the Sun and the sensor are in opposition. In above surface ocean observations, assuming that Sun glint is avoided and there are no whitecaps, _{$L_w({\theta }_v,{\phi }_v,\lambda )$} is determined as:

$L_w({\theta }_v,{\phi }_v,\lambda )=L_t({\theta }_v,{\phi }_v,\lambda )-L_{surf}({\theta }_v,{\phi }_v,\lambda ),$ (2)where _{$L_t$}is the total radiance measured by the sensor and _{$L_{surf}$}is the radiance component of the sea surface corresponding to the skylight reflected by the air-sea interface toward the sensor. _{$L_{surf}$} depends on the polarization state of skylight and can be expressed with the Stokes vector formalism, as follows [13]:

$I_{surf}\left({\mu }_v,{\phi }_v\right)={\mu }_i^{}{R}_{surf}\left({\mu }_i,{\phi }_i;{\mu }_v,{\phi }_v\right)F_0{e}^{-{\tau }_{atm}/{\mu }_0}+{\displaystyle \underset{0}{\overset{2\pi }{\int }}{\displaystyle \underset{0}{\overset{1}{\int }}{\mu }^{\prime }}}{R}_{surf}\left({\mu }^{\prime },\phi \text{'};{\mu }_v,{\phi }_v\right)I_{sky}\left({\mu }^{\prime },\phi \text{'}\right)d{\mu }^{\prime }d\phi \text{'},$ (3)where _{$I={\left(L,Q,U,V\right)}^T$} is the Stokes vector of light with L the radiance and Q, U, V the polarization terms. I_surf and I_sky are the surface and skylight downward component Stokes vectors and R_surf is the reflection matrix of the rough sea surface projected in the reference plane, _{$F_0(\lambda )$} is the extraterrestrial solar irradiance, _{${\tau }_{atm}$}is the total optical thickness of the atmosphere including molecular and aerosol components, µ_i and µ_v are the cosines of the Sun and viewing zenith angles, φ_i and φ_v are the Sun and viewing azimuths, respectively, _{$\mu \text{'}$} and _{$\phi \text{'}$} are respectively the cosine of the zenith angle and the azimuth angle of the incident skylight.

Note that L_surf corresponds to the first element of the I_surf vector and L_sky is the sky radiance, which corresponds to the first element of the I_sky vector. From Eq. (3), it can be seen that accurate computation of L_surf relies on a full knowledge of the skylight directional distribution and its state of polarization in each celestial direction, which are both dependent on aerosol load and type in real world conditions. This kind of computation remains unpractical for current instrument setups and, therefore, some approximations have to be used. For this purpose, the reflectance coefficient was introduced as:

$\rho ({\theta }_v,{\phi }_v,\lambda )=L_{surf}({\theta }_v,{\phi }_v,\lambda )/L_s(\pi -{\theta }_v,{\phi }_v,\lambda ).$ (4)

Due to this approximation, if the composition of the atmosphere changes for a given wavelength or there is a change of observation wavelength, the resulting redistribution of light can cause a change in the _$\rho$ coefficient. The Eq. (2) then can be written as

$L_w({\theta }_v,{\phi }_v,\lambda )=L_t({\theta }_v,{\phi }_v,\lambda )-\rho L_s(\pi -{\theta }_v,{\phi }_v,\lambda )$ (5)

Recommended azimuth angles for near surface observations are _{${\varphi }_v=90°$}or _{${\varphi }_v=135°$} [4]. For a flat surface and neglecting polarization of light, the _$\rho$ coefficient corresponds to the Fresnel coefficient defined by the viewing angle [22] and the indices of refraction of the air and water, but in the presence of ocean waves it is also a function of the illumination conditions and wind speed. For the illumination angle_{${\theta }_i=40°$}${\varphi }_v=90°$,_{${\theta }_v=40°$}, the flat surface coefficient is_{${\rho }_f=0.0256$}, while for a wind speed of W = 5m/s,_{$\rho =0.028$} according to [4]. This last value of the coefficient is usually used in the shipborne above water observations and AERONET-OC processing and is assumed to be independent of the wavelength for the whole range of 400-900 nm which are of interest for OC. It should be noted that a small dependence of _$\rho$ on the wavelength (about 5-10%) exists due to the spectral variation of the sea water refractive index [23, 24]. Also [17], showed that small amount of the Sun glint radiance _{$L_g$}normalized by the sky radiance _{$L_s$}, _{${\rho }_g=L_g/L_s$}, sometimes can be difficult to separate from the sky reflectance_$\rho$, so_{${\rho }_g$}increases towards the red wavelengths since _{$L_g$}has weak spectral dependence and sky radiance decreases from blue to red.

For satellite observations, the total radiance at the TOA in the assumption of no whitecaps is [2, 25]:

$L_{TOA}(\lambda )=L_R(\lambda )+L_a(\lambda )+T(\lambda )L_g(\lambda )+t(\lambda )L_w(\lambda )$ (6)where _{$L_R(\lambda )$}is the radiance from Rayleigh scattering, _{$L_a(\lambda )$} is the radiance from aerosol scattering, _{$T(\lambda )$} is the direct transmittance and _{$t(\lambda )$} is the total (direct plus diffuse) transmittance. The atmospheric correction procedure [2, 25] includes isolation of the _{$L_w(\lambda )$}term from Eq. (6) through the pre-calculation of the Rayleigh radiance _{$L_R(\lambda )$} [1,3] and creating lookup tables (LUT) for various aerosol models to determine aerosol radiance _{$L_a(\lambda )$}. The surface effects are included in both the Rayleigh component (calculated independently of aerosol parameters) and in the radiances for aerosol models [2, 25, 26]. For above water and satellite observations, the determination of _$\rho$for specific illumination and viewing conditions the Cox-Munk model [6] is used, which defines wave slopes as a function of the wind speed. Following [27] and the manual of the RayXP code [11] used in this work, the reflection matrix by the wind-roughened surface when the interface is illuminated from the atmosphere is

$R_{surf}({\theta }_{i1},{\theta }_{v1},{\phi }_v)=P_r({\mu }_{i1},{\mu }_{v1},{\phi }_v)Y_{sur}({\theta }_{i1},{\theta }_{v1},{\phi }_v)$ (7)where

$P_r({\mu }_i,{\mu }_v,{\phi }_v)=\frac{1}{{\mu }_{i1}{\mu }_{v1}{({\mu }_{i1}+{\mu }_{v1})}^4}\frac{1}{2{\sigma }^2}\exp \left[\frac{\frac{2(1-x)}{{({\mu }_{i1}+{\mu }_{v1})}^2}-1}{2{\sigma }_{}^2}\right]$ (8)characterizes effects of the surface, _{$Y_{sur}({\theta }_{i1},{\theta }_{v1},{\phi }_v)$}is the reflection matrix from the flat surface with rotational terms applied, _{$x=-{\mu }_{v1}{\mu }_{i1}+\sqrt{1-\mu {}_{v1}{}^2}\sqrt{1-{\mu }_{i1}^2}\cos ({\phi }_v)$}. The term _{${\sigma }_{}^2$}is the variance [6] of the slope distribution related to the wind speed W at 12.5m above the surface level (to the first order independent of the wind direction):

${\sigma }^2=0.003+0.00512W$ (9)Incident and viewing angles in Eq. (7) and (8) and corresponding cosines are shown with an additional subscript ‘1’, since the angles in RayXP are measured from nadir and thus they are supplementary to the _{${\theta }_i$} and _{${\theta }_v$} defined above. No shadowing effects were considered.

For above surface measurements, the impact of the wavy surface is handled in a similar way; other terms in Eq. (6) are modified in accordance with the sensor position.

It was noticed in [1, 10, 13] that there are significant differences in the simulations using scalar and vector RT models. The RayXP model [11] (version 2.04) used in this work is a 1-D, VRT model and polarization effects are fully taken into account.

The reflectance coefficients were calculated from the VRT model by simulating a fully absorbing ocean. Under this condition, _{$L_w({\theta }_v,{\phi }_v,\lambda )=0$}, and _{$\rho ({\theta }_v,{\phi }_v,\lambda )$} may be computed from Eq. (4). Simulations were conducted for solar zenith angles (SZA) of _{${\theta }_i$} = 30° and 60°, an azimuth viewing angle of _{${\phi }_v$} = 90°, sensor viewing angles of _{${\theta }_v$} = 10° to 60°, wind speeds of W = 0 to 15 m/s (white caps effects were not considered) and aerosol optical thicknesses of_{${\tau }_a(440)$} = 0, 0.1, 0.4, and 1.00. Simulations were performed in the 400-900 nm wavelength range with steps of 20 nm. The wavelength dependence of _{${\tau }_a$} was determined using an Angstrom coefficient (_$\gamma$) of 1.4 for all simulations with a few exceptions described below. Three atmospheric layers were included into the model: top layer making up 64.95% of the total Rayleigh optical thickness_{${\tau }_R$}, a middle layer containing 35% of _{${\tau }_R$} and the full optical thickness of aerosols _{${\tau }_a$}, and third layer with 0.05% of _{${\tau }_R$} between the sensor and the surface. The aerosol single scattering albedo (SSA) was 0.99. Several scattering matrices of aerosols were considered but results are presented only for the oceanic and partially continental aerosol models from the RayXP library [28]. RayXP processing is extremely fast, due to the novel approaches to the solution of the radiative transfer equation [29, 30]; one configuration with the given W and _{${\tau }_a$}is processed in a few seconds for one wavelength on a regular PC. This fast evaluation allows many radiance distributions to be computed in a short amount of time, therefore simulations were performed in the plane of interest for _{${\theta }_v$} = 0-180° with 1° step. An example of such a radiance distribution in _{${\phi }_v$} = 90° plane for the windy surface is shown in Fig. 1 where higher radiance values at _{${\theta }_v$} = 0-25° are due to the Sun glint.

 

Fig. 1 Radiance distribution for four wavelengths in the plane with _{${\phi }_v$} = 90°, _{${\theta }_i$} = 30° and W = 5 m/s with the sensor near the surface.

Download Full Size | PDF

The RayXP code was benchmarked against well-established Monte Carlo [11] and VRT programs [31] for the ocean-atmosphere system and demonstrated excellent agreement with polarimetric measurements in the open ocean and coastal environments [32–35] where the sensors were positioned below water, above water on the ship, on an airplane and on a satellite. However, in this study special attention is required for the accuracy of the modeling of the ocean surface and Fig. 2 shows a very good match between _{$\rho ({\theta }_v)$}values simulated by RayXP and Mobley’s VRT code [10,21] for the windy surface without aerosols at 550nm. A discrepancy exists for _{${\theta }_i$} = 30° and _{${\theta }_v$} = 10 and 20° where there is a strong increase of the reflectance coefficient due to additional Sun glint effects which were not masked anyhow in our simulations.

 

Fig. 2 Comparison of reflectance coefficient simulated by RayXP code for windy surface with Fresnel reflectance and with results of other VRT simulations [10]: a) _{${\phi }_v$} = 90°, _{${\theta }_i$} = 30° and W = 5m/s, b) _{${\phi }_v$} = 90°, _{${\theta }_i$} = 60° and W = 5m/s.

Download Full Size | PDF

In addition, for validation of the RayXP computations, the reflectance coefficients were compared with those obtained using the OSOAA code [36] which numerically solves the radiative transfer equation for the coupled system of atmosphere and water body with a rough air-water interface considering the linear polarization state of light. The results are shown in Fig. 3 for a rough sea following the Cox-Munk distribution with a wind speed of 5 m/s and a pure Rayleigh atmosphere at three different wavelengths. The comparison demonstrates the consistency between the two codes and results are comparable with those obtained at only one given wavelength by [10].

 

Fig. 3 Values of the reflectance coefficient calculated from RayXP (top row) and OSOAA (bottom row) codes using isotropic Cox-Munk slopes distribution for wind speed of 5 m/s. Values are given for a pure Rayleigh atmosphere at 410, 550 and 870 nm. The viewing zenith angle is set to 40°. The concentric circles represent the solar zenith angles by step of 10° (from 0° to 70°) and angles in the polar diagrams represent the relative azimuth with the Sun. (The azimuth is equal to 0° when the Sun and the sensor are in opposition.)

Download Full Size | PDF

3. Results of simulations

For W = 0 and _{${\tau }_a$} = 0, _$\rho$ is almost spectrally flat for all viewing angles, but becomes spectrally dependent in the presence of aerosols as shown in Fig. 4. Changes in aerosol amount and type induce changes in the sky radiance distribution, which in turn changes the _$\rho$value. These changes are shown to be partially associated with polarization; repeating the simulations without considering polarization (an optional feature of the RayXP code) shows significant differences in _$\rho$ if polarization is ignored (Fig. 4(a) vs. Figure 4(c)). For all viewing angles except _{${\theta }_v$} = 40° in Fig. 4(a), an increase of _{${\tau }_a$}leads to lower_$\rho$ values for all wavelengths. For _{${\theta }_v$} = 40°,_$\rho$increases about 0.005 at 400nm and 0.01 at 870nm. Spectral dependence is well visible for all aerosol _{${\tau }_a$}values and _$\rho$ is quite sensitive even to small aerosol additions. The dependence is different for the continental aerosol model (Fig. 4(b)). Scalar mode computations for the oceanic aerosol model and viewing angles of _{${\theta }_v$} = 50° and 60° show that_$\rho$becomes almost spectrally flat. For viewing angles closer to nadir (_{${\theta }_v\le$} 40°), scalar and vector computations exhibit only small differences (Fig. 4(c)). A similar effect is observed for the continental aerosol model (not shown).

 

Fig. 4 Reflectance coefficient spectra for Rayleigh scattering only and with additions of aerosols: _{${\theta }_i$} = 30°, W = 0. Aerosol optical thicknesses _{${\tau }_a$} listed in all figures are specified at 440nm: a) oceanic aerosol model, vector mode; in insert M11 element of the phase matrix (PM) for this model [11]; b) continental aerosol model, vector mode, in insert M11 element of PM for this model; c) oceanic model, scalar mode.

Download Full Size | PDF

Simulated spectra of _{$\rho ({\theta }_v,{\tau }_a)$}for SZA _{${\theta }_i$} = 30° and 60° are shown in Fig. 5(a) and (b) respectively for W = 5m/s. There are very significant changes, especially for _{${\theta }_i$} = 30° in comparison with no wind conditions (W = 0).

 

Fig. 5 Reflectance coefficient spectra for Rayleigh scattering only and with additions of aerosols, W = 5m/s: a) _{${\theta }_i$} = 30°, b) _{${\theta }_i$} = 60°.

Download Full Size | PDF

As can be seen in Fig. 5(a) for _{${\theta }_i$} = 30° and _{${\tau }_a$} = 0, there is a strong increase of the reflectance coefficient with the wavelength due to the additional Sun glint effects for the viewing angles of _{${\theta }_v$} = 10, 20 and 30° (also seen in Fig. 1). The slopes of these curves decrease rapidly with increasing aerosol load leading to the strong spectral dependence for all_{${\tau }_a$}. (Such regimes were used for the retrieval of additional properties of aerosols [37]). The effect is much more moderate for _{${\theta }_v$} = 40° but there is still a spectral change of _$\rho$from 0.03 to 0.04 in the absence of aerosols with smaller changes for larger aerosol loads. In all cases, the reflectance coefficient is strongly biased in the presence of aerosols, in agreement with previous studies [16].

For satellite observations, some areas exhibiting Sun glint will be masked. The normalized threshold for Sun glint (_{$l_{gn}$}) currently used in satellite data processing is _{$l_{gn}$} = 0.005 sr⁻¹ with additional adjustment [38] based on the model of [1]. The maximal expected reflectance coefficients, which correspond to this threshold, were estimated as:

${\rho }_{\max }(\lambda ,{\theta }_v)=(l_{gn}/T_v(\lambda ))*T_0(\lambda )*F_0(\lambda )\cos ({\theta }_i)/L_s(\lambda ,{\theta }_v)$ (10)where _{$T_v(\lambda )$} is the direct transmittance of the sunlight from the surface to TOA which adjusts _{$l_{gn}$}to the surface conditions as a function of the viewing angle, _{$T_0(\lambda )$} is the direct transmittance of the sunlight through the atmosphere. It was found that the thresholds lines for_{$l_{gn}=0.005s{r}^{-1}$} and all viewing angles are in the very upper left corner in Fig. 5(a) (not shown) and do not mask any data in this figure.

Additional glint corrections [38] are more complex and were not applied, but the cases of _{${\theta }_v$} = 30° through 60° are usually not masked and can be further analyzed in terms of satellite observations. In Fig. 5(b) there are no visible Sun glint effects but there are significant changes of _$\rho$as a function of _{${\tau }_a$}for most of the viewing angles. It should be also noted that for _{${\theta }_v$} = 30° and 40° without aerosols,_$\rho$ becomes smaller than the flat surface Fresnel coefficients, as was already observed in Fig. 2(b), and then becomes much higher with increasing _{${\tau }_a$}. Smaller than Fresnel values occur probably due to the redistribution of the skylight on the rough surface which results in smaller reflected radiances at some viewing angles than in the case of the flat surface.

The differences in radiances corresponding to Fig. 4(a) and 5(a) simulated through the full VRT are shown in Fig. 6 where a significant impact of _{${\tau }_a$}is also seen for _{${\theta }_v$} = 30° - 60°, which are primarily not affected by Sun glint. While above water radiometric measurements are conducted, as discussed above, primarily at _{${\theta }_v$} = 40°, in above water imaging a full range of viewing angles need to be considered. Figure 6 provides the differences in radiances for the rough surface (W = 5m/s) compared with W = 0 with Rayleigh atmosphere (_{${\tau }_a$}) and several aerosol loads, which characterize possible errors due to variability of _$\rho$. Further comparing the values in Fig. 6 for no aerosol case (_{${\tau }_a$} = 0, lines with circles) and different aerosol loads (other lines) for example for _{${\theta }_v$} = 40°, we see that possible errors are of the order of 0.01-0.02 mW/cm²/µm/sr depending on_{${\tau }_a$}. In coastal areas, this can result in errors in retrieved water leaving radiances of about 10% in the blue, 2-5% in the green and 20% or more in the near infra-red (NIR). The numbers are much smaller for the open ocean atmosphere with_{${\tau }_a\approx 0$}.

 

Fig. 6 Radiance difference above water surface corresponding to the _$\rho$in Fig. 4a and 5a.

Download Full Size | PDF

In addition to simulations, the ocean surface was studied experimentally using a novel snapshot hyperspectral imager (Cubert UHD285, Germany) which allows the recording of a panchromatic image in the wavelength range of 450-1000nm with 1000x1000 pixels and simultaneously a hyperspectral data cube of 4nm channel spacing for 50x50 pixels. Examples of the panchromatic images are shown in Fig. 7; there are obvious areas of strong glint for viewing angles of 20-30°, smaller for 40° in Fig. 7(a) and no glint for larger _{${\theta }_v$}. It should be emphasized that even in the area of no glint in Fig. 5(a) there is a strong change of _$\rho$with changes of the aerosol load. All glint issues discussed above become small or non-existent for higher solar zenith angles when the wind speed is W = 5m/s, as shown in Fig. 5(b) where there are practically no Sun glint effects with small exceptions. A corresponding image is shown in Fig. 7(b) for _{${\theta }_i\approx$} 50° with very few visible points of Sun glint. However, the changes in the reflectance coefficient in Fig. 5(b) reach 50% at 900 nm compared to the values at 400 nm for some aerosol loads. This is probably due to the “oceanic” aerosol model chosen for the simulations, which is a largely coarse mode aerosol and therefore has a greater effect on NIR wavelengths.

 

Fig. 7 a) Brooklyn pier 05/03/2017, 2:01pm, W = 5 m/s, _{${\theta }_i$} = 29.7° _{${\theta }_v$} = 20°-60°, exposure time _{$t_{\exp }$} = 20ms, _{${\tau }_a\left(440\right)$} = 0.146 b) Duck, NC pier 10/15/17, 2:14pm, W = 5 m/s, _{${\theta }_i$} = 49.3° _{${\theta }_v$} = 20°-60°, exposure time _{$t_{\exp }$} = 22ms, _{${\tau }_a\left(440\right)$} = 0.243

Download Full Size | PDF

With increasing wind speed, as expected, all effects become stronger and more pronounced even for large Sun zenith angles: see Fig. 8 for W = 15m/s, with all curves exhibiting spectral dependence with different spectral slopes depending on aerosol_{${\tau }_a$}.

 

Fig. 8 Reflectance coefficient spectra for Rayleigh scattering only and with additions of aerosols, W = 15m/s, _{${\theta }_i$} = 60°.

Download Full Size | PDF

Analyzing results of the simulations one can see that presence of aerosols changes values of the reflectance coefficient and its spectral dependence and this impact can even have different signs depending on the illumination and observation geometry and_{${\tau }_a$}. In addition, there are many ‘grey’ conditions where small amount of Sun glint is present; it is not eliminated by the threshold even in the satellite processing and should be considered together with ‘pure’_$\rho$to determine _{$L_w$} from_{$L_t$}since it is impossible to separate both. This can also affect black-pixel assumptions in NIR and UV part of the spectra [39, 40] used in the atmospheric correction.

4. Comparison of simulations and measurements, implications to ocean color

4.1 Shipborne observations

Separation of sky and Sun glints is difficult to do experimentally, so results from the simulations above should be probably taken into account to improve the quality of above water measurements. Thus, Lee et al. [14], Fig. 3 and 4, calculated the reflectance coefficient based on measurements in very clear waters in Hawaii with _{${\theta }_i$} = 30°, _{${\theta }_v$} = 30° and W = 8 m/s and claimed that reflectance coefficient can change strongly with wavelength, with_$\rho$values at 700 nm reaching 0.1-0.15. The atmosphere in Hawaii is usually very clear with a typical _{${\tau }_a(440)$}of about 0.05. Comparing these results with Fig. 5(a) (even for a slightly lower wind speed and Rayleigh atmosphere), we see a very good match of experimental [14] and our simulated data. While measurements with _{${\theta }_v$} = 30° look attractive as being more stable to pointing errors than for larger viewing angles, it is clear that this configuration is also quite sensitive to Sun glint and probably should not be recommended for the field measurements. Although there are some complications for each viewing angle, radiometric measurements should always be accompanied by simultaneous measurements of_{${\tau }_a$}. If _{${\theta }_v$} = 40° is used, as recommended by [4], in accordance with Fig. 5(a), for open ocean shipborne measurements with low aerosol load,_{${\theta }_i$} = 30° and W = 5 m/s, spectrally increasing _$\rho$ should be used in Eq. (5), which will lead to the lower water leaving radiance and _{$R_{rs}$}in the NIR part of the spectra than in the case of the spectrally constant _$\rho$ = 0.028. This may eliminate (or minimize) the often used biasing of the _{$R_{rs}$}spectrum through the subtraction of _{$R_{rs}(750)$} [4]. Figure 9 illustrates the difference in _{$R_{rs}$}with each method for two different spectra measured in the field (Stations LS6 and LS9, SABOR cruise, 2014 [35]), where the spectrum with _{${\tau }_a\approx 0.1$} (Fig. 9(a)) almost reaches 0 at 900nm. The impact of variations of the reflectance coefficient is minimal for blue water but becomes more noticeable with increasing chlorophyll-a and colored dissolved organic matter (CDOM) concentrations and where _{$R_{rs}$}increases in the 400-450nm range. In coastal waters (Fig. 9(b)) with typically higher aerosol loads (_{${\tau }_a(440)$} = 0.4), following Fig. 4(a), a spectrally increasing _$\rho$ coefficient should be used leading to_{$R_{rs}(750)\approx 0$}. As already discussed in Section 3, the main impact is on radiance and _{$R_{rs}$} in the NIR. This produces different values in the visible part of the spectrum than use of the spectrally constant _$\rho$ = 0.028 and subtraction of _{$R_{rs}(750)$}. The results can depend on aerosol scattering matrix as shown in Fig. 4, and since the matrix cannot be easily measured in the field, additional studies should characterize _{$R_{rs}$}sensitivity to the variations of aerosol models typical for atmospheres over the ocean before recommending this approach.

 

Fig. 9 Determination of remote sensing reflectance using three different approaches: a) open ocean spectra, _{${\theta }_i$} = 36°, W = 1.5m/s. b) coastal water spectra _{${\theta }_i$} = 31°, W = 3.6 m/s.

Download Full Size | PDF

A viable alternative may also be the use of a polarizer with vertical orientation [5, 15, 41], which can block sky glint, thus potentially avoiding sky measurements at all; this methodology should be revisited since water-leaving light exhibits also specific polarization patterns.

4.2 AERONET-OC observations

SeaPRISM instruments on ocean platforms also make measurements of the water and sky at 40° from nadir and zenith following the recommendations of [4] with a very small field of view 1.5° [42]. To minimize Sun glint the moving radiometer is always positioned at _$\pm$90° from the Sun. As can be seen from Fig. 4 and 5, and as expected due to the chosen observational geometry, Sun glint conditions are not typical for _{${\theta }_v$} = 40°, especially in coastal areas with higher _{${\tau }_a$}, except for small SZA, very clear atmosphere and high wind speeds. Additional filtering is further applied in protocols [42], which practically eliminates Sun glint effects. However, a sky glint correction methodology, which is valid for a broad range of wind speeds, should take into account the variability of the reflectance coefficient. In accordance with established AERONET-OC protocols, only two of eleven $L_t$ measurements are considered for further processing which minimizes sky glint effects followed by the application of a spectrally constant$\rho$[4]. The lowest $L_t$values, due to the instantaneous slope of the sea surface and small integration times (75ms) of SeaPRISM instrument, are often associated with effective viewing angles smaller than 40°. In accordance with Fig. 4, for low wind speeds $\rho$can be smaller than the constant value currently applied and spectrally dependent resulting in a sky overcorrection. Improved processing schemes may be considered based on the spectral dependence of_$\rho$and aerosol characteristics including_{${\tau }_a$}, particle size distribution and phase functions, as these characteristics are available from the same SeaPRISM instrument and AERONET retrieval.

4.3 Implications for satellite observations

To evaluate the impact of _$\rho$variability on the radiances expected on the satellite sensors, several VRT simulations were conducted with the sensor at the top of the atmosphere (TOA) level. The thin layer near the surface was removed, water was still considered fully absorbing. As discussed above, in the current atmospheric correction algorithm [24] Rayleigh component of TOA radiance _{$L_R(\lambda )$} is pre-calculated and corrected for the surface effects without considering aerosol parameters, as they are not known at that stage of the processing.

In the presence of aerosols, the correction for the surface effects should be different depending on the wind speed and characteristics of aerosols. The differences in radiances at TOA for the wind- roughened surface with W = 10 m/s and flat surface (W = 0) for different _{${\tau }_a$}and _{${\theta }_i$} = 30° are shown in Fig. 10(a) and for _{${\theta }_i$} = 60° in Fig. 11(a). Further, assuming that for _{${\tau }_a=0$} effects of windy surface are fully accounted for (lines with circles), in Fig. 10(b) and 11(b) the differences of radiances with _{${\tau }_a>0$} and _{${\tau }_a=0$} from Fig. 10(a) and 11(a) respectively are presented, which estimates possible errors in radiances due to the variation of the surface effects in the presence of aerosols. They were characterized before through the variation of _$\rho$in Fig. 4 and 5. The effects increase with increasing wind speeds and are smaller for lower wind speeds (shown in Fig. 12). In Fig. 10 results for the viewing angles of 20° and 30° are removed since they are strongly affected by the Sun glint. For _{${\tau }_a\le 0.1$}, the errors are in the range of 0.01-0.02 mW/cm²/μm/sr which is up to 5-10% of the TOA radiance in the NIR and thus can affect atmospheric correction in the open ocean and through it satellite vicarious calibration [43,44].

 

Fig. 10 a) Differences in radiances at TOA between wind roughened and flat surface for Rayleigh scattering only and with different aerosol loads, _{${\theta }_i$} = 30°, W = 10m/s, b) Expected errors in radiances due to the non-including of aerosols in the processing of the surface effects.

Download Full Size | PDF

 

Fig. 11 a) Differences in radiances at TOA between wind roughened and flat surface for Rayleigh scattering only and with different aerosol loads, _{${\theta }_i$} = 60°, W = 10m/s, b) Expected errors in radiances due to the non-including of aerosols in the processing of the surface effects.

Download Full Size | PDF

 

Fig. 12 Expected errors in radiances at TOA due to the non-including of aerosols in the processing of the surface effects, W = 5m/s: a)_{${\theta }_i$} = 30°, b) _{${\theta }_i$} = 60°.

Download Full Size | PDF

In coastal waters with similar wind speeds but higher aerosol load, the impact on the water leaving radiance and atmospheric correction will be strong with the errors of up to 20% in the NIR at TOA and up to 25% and more for the water leaving radiance in the blue (typically 0.2-0.3 mW/cm²/μm/sr). Estimated errors for W = 5m/s are presented in Fig. 12 for the same two incident angles showing smaller effects for _{${\tau }_a\le 0.1$} and thus smaller impact on the open ocean processing but with still strong impact for coastal waters (_{${\tau }_a>0.1$}). The radiance errors in Fig. 10(b) and 11(b) will be partially mitigated by the non-optimal choice of aerosol model but the complex spectra like shown in Fig. 11(b) and 12(b) cannot be mitigated fully and thus can at least partially explain the vulnerability of the blue bands to the atmospheric correction. This is in line with recent findings [18] where the comparison of aerosol’s phase functions measured at AERONET-OC sites and used in the satellite atmospheric correction match only in about 10% of cases and the difference between satellite and AERONET-OC _{${\tau }_a$}depends on the wind speed.

As mentioned before, all simulations (with few exceptions described above) were carried out with the “oceanic” aerosol model from RayXP library. Obviously, aerosol models with different fractions of coarse and fine mode aerosols and thus different phase matrices [26, 44] will also have different impacts on the reflectance coefficient and corresponding radiances. These effects should be further studied.

In coastal areas aerosols are often absorbing. To analyze this effect, additional simulations of the TOA radiances were conducted with the single scattering albedo (SSA) reduced from 0.99 in all previous simulations to 0.8. Results for _{${\theta }_i$} = 30° and W = 5m/s are shown in Fig. 13 demonstrating similar impact of aerosols as in Fig. 10a with SSA = 0.99.

 

Fig. 13 Differences in radiances at TOA between wind roughened and flat surface for Rayleigh scattering only and with different aerosol loads for the case with absorbing aerosols (SSA = 0.8)

Download Full Size | PDF

All previous simulations were carried out with an Angstrom coefficient of _$\gamma$ = 1.4. Comparison of these results with radiances simulated with _$\gamma$ = 1.0 and 2.0 revealed, as expected, only small differences with the previous cases in the red and NIR part of the spectrum which nevertheless can affect the atmospheric correction procedures.

5. Conclusions

Comprehensive simulations of the spectral characteristics of the reflectance coefficient of the sky from the ocean surface in the presence of various aerosol loads and the corresponding radiances at the near surface level and at TOA were carried out by a well benchmarked VRT code. The results reveal a profound impact of the aerosol _{${\tau }_a$}and other parameters, when included in the surface processing, on the radiances measured by the sensors at all levels, which greatly affects the quality of retrieval of the water leaving radiance.

It is shown that due to the very significant effect of aerosols on the reflectance coefficient, in any above water radiometric observations, radiance measurements should be always complemented by simultaneous measurements or simulations of _{${\tau }_a$}.

For shipborne observations, the viewing angle of _{${\theta }_v$} = 40° for water and sky measurements is still probably the best choice but _$\rho$spectra for the typical Sun zenith angles, wind speeds and _{${\tau }_a$} can be tabulated and used accordingly to determine the proper water leaving radiances. The subtraction of _{$R_{rs}(750)$} from the whole spectrum is probably not necessary but this requires further experimental tests due to the variability of aerosol parameters. Above water measurements with a polarizer with vertical orientation for the blocking of Sun and sky glints looks attractive and should be further tested through comparison with in-water measurements.

In the AERONET-OC data processing where measurements are also carried out with a viewing angle of _{${\theta }_v$} = 40°, similar to shipborne observations, tabulated _$\rho$ spectra based on measured _{${\tau }_a$} should improve data quality and in some conditions sky overcorrection with modified algorithms for sky glint removal.

In the processing of current data from satellite sensors it is expected that at high wind speeds (W_$\approx$10m/s) even with low aerosol loads _{${\tau }_a\le 0.1$} inaccuracies of processing of surface effects without including aerosols can affect atmospheric correction and vicarious calibration; Rayleigh lookup tables (LUT) used in the atmospheric correction include corrections for non-flat surfaces [43, 44] but without additions of aerosols. The effects are smaller for lower wind speeds, so the best satellite retrievals occur in areas with low aerosol loadings and moderate wind speeds similar to the ones at the Marine Optical Buoy (MOBY) site. In coastal areas with generally larger aerosol loads, the atmospheric correction procedure can be more strongly affected by the errors in the NIR and resulting in large impact on the water leaving radiance in the blue part of the spectra. It may be improved with correction for the variability of _$\rho$ spectra for various Sun and viewing angles, wind speeds and aerosol loads. Technically that will require multiple LUT with different _{${\tau }_a$} correcting main LUTs for Rayleigh scattering and/or an additional iteration(s) in the atmospheric correction when/if _{${\tau }_a$}was accurately estimated.

Obviously the presented results provide only mean values of the reflectance coefficients and radiances are calculated only for fully absorbing waters. The whole picture becomes much more complex when the variances of wind speeds, observational geometries and water conditions are taken into account. First results in such possible variations were provided in [16] from Monte Carlo simulations. Continuation of this work as well as field studies with the snapshot imager are underway which are expected to shed more light on the impact of the wind- roughened surface on the measured radiances.