Diffuse Attenuation Coefficient Retrieval in CDOM Dominated Inland Water with High Chlorophyll-a Concentrations

The kd(490) is a vertical light attenuation coefficient and an important parameter for water quality. The kd estimates are often based on empirical and semi-analytical algorithms, designed for oceanic and coastal waters. However, there is a lack of information about the performances of these models to inland waters dominated by chromophoric dissolved organic matter (CDOM). Therefore, to contribute to this investigation, nine empirical models based on the blue-to-green and blue-to-red ratios and chlorophyll-a (Chl-a) concentration were evaluated, as well as three semi-analytical models using bands from the Operational Land Imager (OLI)/Landsat-8. The errors (mean absolute percentage error, MAPE > 80%) presented by the empirical models confirmed that the blue-to-green ratio failed in retrieving kd(490) in an environment dominated by CDOM. Similar failures occurred with the models using the Chl-a concentration (MAPE ~60%) as input. A semi-analytical approach showed the lowest error (MAPE = 41.04%) in the estimate of the inherent optical properties for complex waters in order to reduce the errors above. After retrieval of kd(490) using the semi-analytical model, seasonal patterns were observed, and high values of kd(490) were detected in the dry season possibly due to the increase of the concentration of the optically-significant substances (OSS).


Introduction
The diffuse attenuation coefficient (k d ) is an apparent optical property (AOP) and is defined as the exponential decrease of the ambient downwelling irradiance (E d ) with depth and is, therefore, related to light penetration and availability and can be used to predict the euphotic depth. AOPs are those properties that depend both on the medium and on the geometric structure of the radiance distribution [1]. k d is largely determined by the inherent optical properties (IOPs), absorption (a), and backscattering (b b ) coefficients of the first-order, and are less dependent on the incident radiation field as the sun angle. IOPs are properties of the medium and do not depend on the ambient light field [1].
Non-algal particles (NAP), phytoplankton, chromophoric dissolved organic matter (CDOM), and water itself are considered the four optically-significant substances (OSS) that control the k d. The k d at 490 nm, k d (490), is generally considered and classified as a parameter of water quality and, therefore, is essential for monitoring the eutrophication process due to light attenuation by phytoplankton growth or suspended matter [2]. The capacity of the inverse of k d (490) to provide the depth from which 90% of the upwelling signal at 490 nm originates allows for the remote estimation. For primary production models, the k d (490) was used as a surrogate for the attenuation coefficient of photosynthetically-active radiation, k d (PAR); the significant positive correlation between k d (490) and k d (PAR) allowed the k d (PAR) inference from k d (490) [3,4]. For monitoring eutrophication and the associated decrease in light attenuation, the accurate estimation of k d is critical to better understand and model the primary productivity, heat, and gas transfer in aquatic systems.
In order to estimate k d (490), empirical algorithms were developed using a direct conversion from normalized water-leaving radiance (nL w ) or remote sensing reflectance (R rs ) ratios or indirect conversion through products based on chlorophyll-a (Chl-a) concentrations [3][4][5][6][7][8][9]. However, these models are considered site-specific, which is a limiting factor for broader applications. Wang et al. [9] and Lee et al. [10] came up with a semi-analytical approach covering a wide range of waters. Wang et al. [9] used a combination of two models suitable for ocean and turbid coastal waters, aiming to improve the application range of k d (490) estimations. The model from Lee et al. [9] considered the IOPs derived from R rs through a quasi-analytical algorithm (QAA) [11] and the solar zenith angle using the radiative transfer theory. The algorithm makes it possible the estimation of spectral k d , as well as k d (490), for water bodies ranging from the clearest ocean to turbid coastal waters.
Even so, studies showing the use of remote sensing for k d (490) estimation in tropical productive inland waters are hampered by the lack of in situ data and theoretical framework to predict and interpret water color data [12]. The optical complexity of turbid waters (mineral suspended solids, algae, and associated organic particles) produce multiple scattering, hampering the mathematical analysis of radiative transfer [12]. From our understanding there are very few investigations of k d (490) estimation through the oceanic-based algorithms in inland waters dominated by CDOM, which is an optically-active indicator of dissolved organic matter and plays an important role in cycling carbon.
CDOM is considered an important water quality indicator due to its impact on drinking water, carbon balance, and aquatic ecosystems; mainly because CDOM affects the penetration of light into the water column, which affects the primary productivity. Most models use R rs or nL w from the blue-green spectral region, in which the CDOM has high absorption. Due to this, the following questions arise: (1) Will the k d (490) models that use the blue/green ratio be impacted in such a way that the errors make it impossible to use? (2) How will these models be impacted in the presence of high Chl-a concentrations? (3) Will the semi-analytical models, which use IOPs as the input data, perform better than the empirical ones? (4) What are the potential for k d monitoring from space in inland waters dominated by CDOM? Therefore, the main goal of this paper was to assess the performance of algorithms to estimate k d (490) in inland waters dominated by CDOM, which are widely available to end users. We tested nine empirical and three semi-analytical algorithms. The selected models will help us to answer the above questions.

Study Area and Sampling Planning
The Bariri hydroelectric reservoir (22 • 9 49.260" S 48 • 44 21.420" W) is in the middle of São Paulo State and is part of the Cascading Reservoir System (CRS) of the Tietê River ( Figure 1). Bariri is the second of a total of six reservoirs and presents the smallest flooded area of 63 km 2 at an average altitude of 450 m. The Bariri reservoir is situated in a tropical climate with a dry period (April-September) and a wet period (October-March) according to the Köppen classification. Water retention time varies between 7 and 24 days [13,14].
The amount of wastewater coming from the metropolitan region of the Tietê River characterizes the Bariri Reservoir as highly productive water with high average concentrations of total nitrogen (2750 µg·L −1 ), phosphorus (87 µg·L −1 ), and Chl-a (55.8 µg·L −1 ) with the phytoplankton community dominated by cyanobacteria.
Two surveys were carried out where, in the first case (BAR1), 30 samples were collected from 15-18 August 2016 and, in the second case (BAR2), 18 samples were taken from 23-24 June 2017 (see Figure 1 for sample locations). The days of the field campaigns were planned with respect to the temporal resolution of OLI/Landsat-8, aiming to match the data collection with the satellite overpass of the study site.

Water Quality and Optical Data
Wind speed (m s −1 ), pH, depth (m), Secchi disk depth (Z SD ; m), and turbidity (NTU) were measured in the field with an anemometer, pH meter, depth gauge, Secchi disk, and turbidity meter, respectively. The water samples were collected in the field to derive the OSS, such as the Chl-a and the suspended particulate matter (SPM), as well as the inorganic (IPM) and organic particulate matter (OPM) concentrations according to Golterman et al. [15] and the American Public Health Association Protocol [16], respectively.
The IOP data analysis followed the methodology proposed in Bricaud et al. [17] for the absorption coefficient of the colored dissolved organic matter (a CDOM ). A baseline correction was performed by subtracting the average value between 700 and 750 nm from all the spectrum values; while for the absorption coefficients of non-algal particle (a NAP ) and the total particulate (a p ), that consists in the sum of a NAP and phytoplankton absorption coefficient (a ϕ ), the acquisition and analyses were made in accordance with the transmittance-reflectance method, according to Tassan and Ferrari [18][19][20].
The radiometric data were obtained from the two hyperspectral RAMSES sensors (TriOS, Rastede, Germany). The sky and total radiance data (L sky and L t , respectively, both in W m −2 sr −1 ) and the downwelling and sky irradiance data (E d and E sky , respectively, both in W m −2 ) were taken with the sensors fixed by a steel frame in the boat, in a configuration of 40 • from nadir (zenith) and to azimuthal angle of 90 • in order to minimize the specular reflection [21]. The radiometric quantities were used to calculate the remote sensing reflectance [22] according to Equation (1): where L SR (λ) is the surface-reflected light and consists of a fraction (ρ sky ≈0.02-0.05) of L S (λ) and the sun glint, residual glint, whitecaps and foam-reflected light (∆). The ∆ value can be assumed to be zero for oceanic water in the near-infrared (NIR); however, in turbid water, this element must be retrieved.
In agreement with that, Lee et al. [22] proposed a calculation approach of R rs as a function of the total remote-sensing reflectance (T rs , ratio of L t to E d ) and sky remote-sensing reflectance (S rs , ratio of L sky to E d ) for each L t and L sky scan from the Equation (2): where ∆ is spectrally constant, settled before R rs can be derived [22]. F is the surface Fresnel reflectance based on the viewing geometry (~0.021), and S rs (λ) is the sky remote sensing reflectance. The value of ∆ is then numerically derived by spectral optimization, which minimizes the error between the modeled R rs and optimized R rs [22]. In sequence, the R rs used in estimate models of k d must be simulated by the satellite signal at each spectral channel centered at a wavelength (λ = 443; 482; 561; 655 nm). The band simulation consists in the convolution of the radiation signal of the hyperspectral sensor and the spectral response function [F r (λ)] of the OLI sensor in a wavelength interval of the spectral resolution [23]: where R c rs (λ) is the remote-sensing reflectance simulated at center wavelength; and x max and x min are, respectively, the maximum and minimum values of the sensor spectral channel.
According to E d (490, 0-) data, the k d (490, 0-) in both field campaigns were determined as the slope of -ln[E d (490, z)] at subsurface depth (z) [1]: The choice of depths used in derivation of k d (490) has important impact in its accuracy. The wave-introduced fluctuations in the subsurface light field generates about 20-40% of the variability in the signal [24,25]. Therefore, in the subsurface, the E d profile was visually analyzed in order to choose the upper bound with the minimal effect of focusing and defocusing sunlight on the water's surface [26]. For determination of the lower bound, the depth of 1 optical depth (OD), where the E d (z) signal achieves 1% of the subsurface value, was considered [27].
The E d measurements were taken in an integration time of 100 s generating 10 spectral curves which the most representative one was chosen according the median method. In BAR1, the lower bound used for k d (490) derivation was of 4 m and in BAR2, the depth considered was of 6 m. Therefore, five and seven E d measurements, respectively, were used in BAR1 and BAR2 for the k d spectrum at each sampling point.
In order to eliminate the noise in E sky (λ) due to changes in the sun illumination condition caused by cloud cover during the E − d (z) measurements, a normalization factor was required in all scans.
The E d normalization consists of a division of E sky at first scan t(z 1 ) to E sky at subsequent scans t(z i ) as the factor normalization of E − d t(z i ) according to Mishra et al. and Mueller [5,26]:

K d (490) Models
The Table 1 summarizes the k d (490) algorithms tested in this work.
Global Ocean Waters * k w (490) = 0.016 m −1 is a constant of the diffuse attenuation coefficient for pure sea water.
We tested three types of models: (1) empirical relationships between the k d and Chl-a concentration; (2) empirical relationships between the water-leaving radiance (nL w ) and k d ; and (3) semi-analytical models which are based on radiative transfer models.
The nL w is calculated through the conversion of R rs as nL w = F 0 (λ)R rs (λ) where F 0 (λ) is the extraterrestrial solar irradiance [9,28]. The R(λ) is the irradiance reflectance beneath the water-surface and is calculated as a function of extraterrestrial solar irradiance [F 0 (λ)] and nL w [9].
The equations from Morel et al. [3] are based just on Chl-a concentration, whereas the other empirical algorithms use a simple ratio of nL w (λ). Zhang and Fell [29] and Wang et al. [9] are empirical models based on a R rs (λ) ratio. The semi-analytical model of Wang et al. [9] was developed using the R(λ) ratios, and the Lee et al. [10] model was preliminarily used for the estimatiion of a(λ) and b b (λ), which are obtained from the QAA_v5 (the version 5 was used in this work) [30] with R rs as the input at the reference wavelength (λ 0 ) of 561 nm, the nearest OLI center band of λ 0 = 55× of QAA_v5.

OLI/Landsat-8 Data Processing and Acquisition
The Landsat 8 Surface Reflectance on-demand data generated by the Landsat Surface Reflectance Code (LaSRC) were obtained in the U.S. Geological Survey platform (http://earthexplorer.usgs.gov/). The LaSRC algorithm for atmospheric correction was developed using the Second Simulation of the Satellite Signal in the Solar Spectrum Vectorial (6SV) model. The algorithm uses the OLI Coastal Aerosol Band (0.433-0.450 µm) which works as cover for shorter wavelength as the blue band in previous Landsat and helpful to retrieve aerosol properties [31]. Pahlevan et al. [32] and Bernardo et al. [33] showed that LaSRC is a consistent product to derive aquatic estimates. A total of 10 atmospherically corrected images (Path/Row 221/75) were acquired during the year of 2016, covering the months from February to November. The criteria used for the image choice was attributed to cloud-free data over the reservoir.
To obtain the nL w (λ) and R(λ), and to test the k d (490) models, the surface reflectance images were divided by π to convert them into R rs [34]. The k d (490) model with the best performance was applied to the time-series of LaSRC images selected to obtain the spatial distribution of the vertical attenuation coefficient. The k d (490) model best fitted for Bariri reservoir was validated through the use of the satellite image from 15 August 2016, correspondent to the first day of BAR1. In the same way, the LaSRC accuracy in k d (490) estimations was evaluated from the analysis between the R rs converted images with the R rs calculated from TriOS data in BAR1.

Statistical Analysis and Accuracy Assessment
Statistical analyses, including calculations of the maximum, minimum, and average values, and linear and non-linear regressions, were performed. The k d (490) analyzed models were applied for BAR1 and BAR2 datasets. The root mean square error (RMSE), mean absolute percentage error (MAPE), and the bias were used to assess the accuracy of the k d (490) models: where x est is the k d (490) estimated value and x mea is the in situ measurement of k d (490).

Remote Sensing Reflectance Spectra
In BAR1 spectral curves ( Figure 2a) include a high absorption feature at about 680 nm, and a reflectance peak at approximately 710 nm that can be associated with the Chl-a pigment. The sampling sites with the highest Chl-a concentrations presented a significant peak at 550 nm, expected for water with a great amount of algae, as seen in the field. The increase of reflectance in longer wavelengths is indicative of the increase of the total suspended matter concentration. The Chl-a concentration in BAR1 was approximately 37 times higher than in BAR2 dataset (see Table 2), which explains why pigment absorption was more evident during the first survey. Still, in BAR1, the green and dark red curves (Figure 2a) represent two sampling points with expressive reflectance peaks at 550 nm and 700 nm and the highest values of Chl-a concentration of 709.89 µg L −1 and 623.70 µg L −1 , respectively, among the other sampling points. Due to this high concentration of Chl-a some spectra presented the red edge feature as found in healthy vegetation, which explain the rise in the reflectance near 689 nm. In both datasets, the organic particles were predominant (Table 2), therefore, the spectra of BAR1 and BAR2 were quite similar with less reflectance in BAR2 mainly because the lowest OSS concentrations were verified in relation to BAR1. In the BAR2 spectra, the features characteristic of phytoplankton absorption were not evident such as that for BAR1 due to the decrease of Chl-a concentration, with the exception of the feature of the sampling point that presented the maximum value of 19.11 µg L −1 of the Chl-a concentration in BAR2 (navy blue curve in Figure 2b).
Light absorption by CDOM was highest in both field campaigns (1.57 ± 0.37 m −1 in BAR1 and 1.98 ± 0.76 m −1 in BAR2) in relation to light absorption by other OSSs. In CDOM-rich lakes, the spectra shapes vary according to CDOM levels and the influence of other substances in absorption processes [35]. In BAR1, the light absorption was dominated by CDOM, however, the reflectance spectra assumed a shape of Chl-a due to the high concentrations verified, and presented low reflectance values at 400-500 nm. In BAR2, with the predominance of CDOM in the light absorption process and the reduction of OSS concentrations, reflectance values assumed the nearly flat feature at~600 nm.  Table 2 presents the descriptive statistics for optical and water parameters for the two field campaigns, BAR1 and BAR2.

Optical Properties and Water Constituents
In BAR1, it was sunny with some overcast moments. In BAR2, the sky was more favorable with sunny days during all field campaigns. In all datasets, the average wind speed was 3.3 ± 1.8 m s −1 producing some small waves on the water surface. The minimum values for both field campaigns were 0 m s −1 and the maximum value was 8 m s −1 in BAR1 and 6.5 m s −1 in BAR2.
The water collected in BAR1 was green in most of the samples, which is explained by the presence of phytoplankton pigments due to high Chl-a concentration averaging 119.8 ± 96.4 µg L −1 , and ranging between 25.7 and 709.9 µg L −1 . In BAR2, the Chl-a variability reduced significantly, with an average of 8.0 ± 3.3 µg L −1 , varying from 3.8 to 19.1 µg L −1 . Therefore, the mixed data resulted in a range from 3.8 to 709.9 µg L −1 of Chl-a.
The SPM concentration (average of 8.4 ± 4.6 mg L −1 ) was predominately of organic particles with an average of 6.06 ± 4.57 mg L −1 in BAR1, as well as in BAR2, corresponding to a major part of SPM (1.59 ± 0.44 mg L −1 ) with an average of 1.11 ± 0.32 mg L −1 . The SPM presented a low concentration in BAR2. With respect to absorption, a NAP (482) was higher in BAR1 (average of 0.4 ± 0.1 m −1 ) than BAR2 (average of 0.1 ± 0.04 m −1 ). The mean of a CDOM (482) was 1.00 ± 0.3 m −1 in BAR1 and 1.1 ± 0.6 m −1 in BAR2. The mean of a ϕ (482) was 0.7 ± 1.0 m −1 in BAR1 and 0.07 ± 0.02 m −1 in BAR2. Figure 3 shows the average spectra of IOPs for BAR1 and BAR2 datasets. The a ϕ spectra was higher in BAR1 than BAR2; it is possible to see two peaks of absorption at 450 and 670 nm in BAR1 (Figure 3a) whereas in BAR2 (Figure 3b) the characteristic a ϕ spectra was not evident. In BAR2, the OSS concentration was less than in BAR1, therefore, the absorption processes were reduced in relation to BAR1, with the exception of a CDOM . In BAR2, the a ϕ and a NAP spectra were similar with light absorption by NAP slightly higher than the absorption by phytoplankton, with the absorption curve decay at 550 nm (Figure 3b). The absorption curve of NAP in the blue portion of the visible spectrum represents the typical spectra of mineral and/or detrital particles, verified in BAR1 and BAR2 absorption spectra. The atypical a ϕ spectra in BAR2, with the linear increase at wavelengths shorter than 500 nm, indicate an environment with elevated CDOM absorption [36]. Among the three absorption coefficients previously mentioned, a CDOM (482) showed the highest values in both field campaigns. Zhang et al. [37] proved that high Chl-a concentration and high a CDOM indicate that the accumulation and degradation of phytoplankton are a source of CDOM in eutrophic waters. Organic matter likely originates from phytoplankton degradation in BAR1. In BAR2, the a CDOM is also predominant with significant reduction of Chl-a concentration, which can suggest an allochthonous source of CDOM. The Bariri reservoir has been experiencing impacts due to sugarcane production and the input of urban and industrial wastewaters [38,39]. The proportional contribution of the IOPs in BAR1 and BAR2 among the total absorption (a t ) budget in the Bariri aquatic environment was computed considering the OLI bands (443, 482, 561, and 655 nm) (Figure 4a-d).
The average of k d (490) in BAR1 and BAR2 was 2.8 ± 0.3 m −1 (ranging from 1.9 to 4.0 m −1 ) and 1.8 ± 0.1 m −1 (ranging from 1.5 to 2.3 m −1 ), respectively. The water column was more attenuated in August 2016 than June 2017. Both months are in the dry season; however, 2016 presented the highest mean rate of precipitation from January to June when compared with 2017 over the Bariri reservoir area (Figure 1d). The rain carries organic matter, mainly from industrial effluents, wastewater, and particulate matter to the Bariri reservoir, as well as facilitates the increase of organic matter fluxes along the Tietê River from upstream to downstream. The Tietê River is in one of the most industrialized basins of São Paulo State and the reservoir's construction promoted the rapid transformation in the land use and land cover, facilitating the pollution problems and accelerating the sedimentation process [38]. The concentration of water constituents was higher in August, increasing the process of absorption and the backscattering of light and, consequently, showing the highest k d (490) values. Figure 5 presents the spectra of k d derivation of BAR1 and BAR2.  After normalization of Ed measurements, it was possible to derive k d in an accurate way and reduce the noise introduced by the wave-focusing. In BAR1, the k d (490) showed higher magnitude in relation to BAR2.
In both field surveys, the SPM concentration ([SPM]) was dominated by organic matter. The relationship between [SPM] and Chl-a concentration ([Chl-a]) presented a significant positive correlation (R = 0.97; p < 0.001) for the mixed data (Figure 6a). High proportions of [Chl-a] in water indicate excessive input of nutrients which increase the attenuation of light in water [39], verified by the positive correlation (R = 0.74; p < 0.001) between k d (490) and [Chl-a] for the Bariri mixed data (Figure 6b). The nutrient enrichment in water accelerated the primary productivity with a great amount of phytoplankton production, confirmed by the strong relationship between [Chl-a] and a φ (482) (R = 0.89; p < 0.001) in Figure 6c.

Assessment of the Vertical Attenuation Coefficient Models
Due to a lack of data from turbid inland waters used in the calibration of k d (490) algorithms, the application for the Bariri reservoir can be compromised. The algorithm coefficients developed for a wide range of waters tried to minimize this limitation although it is not a guarantee of success [25]. For all k d (490) algorithms, bands were defined according to the center bands from OLI/Landsat-8, including the blue band at 490 nm that was replaced by 482 nm. Table 3 summarized the statistical results yielded after the application of k d (490) models. The empirical models of k d (490) based on Chl-a concentration underestimated the results for almost the entire dataset, except for two samples with high values of Chl-a concentration which were overestimated. The coefficients of the model from Morel et al. [3] were developed for Case I waters, with the optical properties dominated by phytoplankton, therefore, the accuracy of the model in aquatic environments with high rates of CDOM and/or suspended solids is expected to fail [28]. The models from Morel et al. [3] were calibrated with open ocean waters with values of Chl-a < 2.4 µg L −1 , much lower of the values found in Bariri reservoir, resulting in k d (490) errors of MAPE~61% and RMSE~1.5 m −1 .
The empirical models using the nL w spectral ratio developed for ocean open and global waters presented the worst results with MAPE around 80% and RMSE~2 m −1 . Mueller [5], Mueller and Tress [6], Chauhan et al. [7], and Werdell [40]  The blue-green ratio nL w (490)/nL w (555) or R rs (490)/R rs (555) have large uncertainties when k d (490) is greater than 0.2 m −1 . The spectral ratio presents an asymptotic value with increasing OSS concentration [5]. Thus, this spectral ratio is not sensitive for the variations in OSS that occurred in turbid inland waters, resulting in significant underestimation of k d (490). In addition, this spectral ratio does not consider the effects of solar angle changes that decrease the accuracy of empirical algorithms for estimating k d (490) [25].
Aiming to investigate the estimation of k d (490) from MERIS, Kratzer et al. [8] developed the empirical algorithm using the ratio of 490 nm and 620 nm for Case 2 water optically dominated by CDOM, which showed an inverse relationship to salinity in the Baltic Sea. The Kratzer et al. [8] model applied to all datasets from Bariri presented a MAPE of 60.6% and an RMSE of 1.6 m −1 .
The empirical algorithm developed by Wang et al. [9] combined open ocean and turbid coastal waters (from 0.3 to 0.6 m −1 ) through a linear regression equation of the R rs (670)/R rs (490) ratio that showed better matching of k d (490) data when the k d (490) values of calibrated data were higher than 0.3 m −1 . The application of this model in the Bariri reservoir data resulted in a MAPE of 54.0% and an RMSE of 1.4 m −1 , yielding underestimated values (bias of −1.1) of k d (490).
The semi-analytical algorithms of Wang et al. [9] were developed considering the MODIS satellite data and the absorption and backscattering coefficients derived from Lee et al. [11]. The algorithm based on the ratio R(667)/R(488) presented a MAPE of 46.3% and an RSME of 1.3 m −1 , while the ratio R(645)/R(488) showed a MAPE around 63.7% and an RMSE of 1.7 m −1 . The OLI sensor does not have a spectral band centered at 488, neither at 667 nm, nor at 645 nm, but at 482 and 655 nm, respectively. Therefore, in both semi-analytical equations, R(482) and R(655) were used and the differences found in the results above are explained by the coefficients, which, in the first model, presented better adjustment than the second model in the Bariri reservoir.
The Zhang and Fell [29] empirical algorithm and the semi-analytical algorithm developed by Lee et al. [10] presented similar performances corresponding to a MAPE of 42.90%, an RMSE of 1.22 m −1 , and a MAPE of 41.0%, and an RMSE of 1.1 m −1 , respectively. Calibrated with a wide range of environments (0.01-4.6 m −1 ), the Zhang and Fell [29] model have showed greater correlation for clear waters, whereas the Lee et al. [10] algorithm showed the greatest application for turbid waters, such as what occurs in this study site. The difference in errors among these empirical and semi-analytical algorithms could be explained by the use of the 655 nm wavelength in the Zhang and Fell [29] model, which introduces noise due to the high absorption rate of water in the red wavelength, even after corrections in data processing.
The Lee et al. [10] semi-analytical algorithm showed the lowest errors (MAPE = 41.0%; RMSE = 1.1 m −1 ) for k d (490) estimate in the Bariri reservoir. In inland waters, the regional variability of the OSSs generate significant changes in the attenuation of the light process, therefore, the Lee et al. [10] algorithm in considering a(λ) and b b (λ), estimated via QAA_v5 [31], and the sun angle resulted in fewer uncertainties regarding data matching. Additionally, the algorithms' coefficients derived from numerical simulations reduce the dependence of a specific data range.
Empirical algorithms were developed (results not shown here) from Bariri using all datasets, however, the results were not satisfactory. Therefore, the Lee et al. [10] model was chosen for mapping the variability of k d (490) in the study site.

Assessment of the OLI Atmospheric Correction Product
The evaluation of LaSRC accuracy consisted in the comparison between the OLI-derived R rs product from August 2016 and the R rs calculated from TriOS data in BAR1. The relative percentage errors were 50.1% at 443 nm, 20.7% at 482 nm, 17.2% at 561 nm, and 29.6% at 655 nm. The highest errors in the coastal and blue regions are related to the increase of scattering in those regions. The lowest error is verified in the green region with 17.2% and, in sequence, in the red region, the error increased again to 29.6%. The LaSRC accuracy results encountered here are compatible with the atmospheric correction analysis presented in Rodrigues, T. et al. [41] and Bernardo, N. et al. [33].

Application of k d (490) on OLI/Landsat-8 Images
The variability of k d (490) in the Bariri reservoir expressed in OLI/Landsat-8 is presented in Figure 6. The spatial distribution of Lee et al. [10] in the Bariri reservoir was validated from the k d (490) estimates obtained from the image of 15 August and presented a MAPE of 32.5% and an RMSE of 0.9 m −1 . The errors of k d (490) estimates from the Lee et al. [10] model, applied for all datasets, using the in situ R rs , were 41.0% of MAPE and 1.1 m −1 of RMSE. The differences between the k d (490) estimations from in situ R rs and R rs from LaSRC are related to the decrease of the variability among the sampling points through the sample size reduction that, in the former case, was 48 (all datasets) and, in the latter case, was 26 (four sampling points with missing information in the image from August in the region highlighted as zone B in Figure 7) that presented the lowest errors.
The increase of the percentage error between the OLI-derived R rs at 482 nm used as the input of k d (490) estimation and the k d (490) estimations for BAR1 was 11.9% and could be related to IOP estimations via QAA_v5. The application with no modifications in some empirical steps can affect the QAA_v5 performance in optically-complex waters [42,43]. The maximum value of k d (490) was 5.60 m −1 in April and the minimum value was 0.89 m −1 in September, both in the dry season. In the image from April, specifically in zone B highlighted by the arrow, an expressive increase of the k d (490) was observed, which is probably related to the reduction of water velocity favoring the increase of Chl-a concentration and, in a dramatic way, an algal bloom, intensifying the attenuation of light in the water. The decomposition of algae within the water can characterize a source of CDOM, which also intensifies the attenuation of light in the water [35]. Along the reservoir, except for zone B, the k d (490) values remained below the mean (Figure 7), which was expected for a period of low precipitation, confirming an atypical event in April. In September, the k d (490) was uniformly distributed around the mean, with values closer to the maximum (4.40 m −1 ) in the northeast region of the reservoir which presents agricultural activities and bare soil areas that facilitates the input of external material into the reservoir even in periods of low precipitation. Other regions with high values of k d (490) were in zone B and another bend where water flux is slower than in zone A.
The boxplot in Figure 8 summarizes the variability of k d (490) for each month. In the dry season, the k d (490) values were closer to the maximum in May and June, with most values above the mean and, concomitantly, values above 2.5 m −1 in zone A. In July and August, the k d (490) started to decrease with most of the values closer to the minimum and below the mean, with a decrease of k d (490) in zone A. The precipitation rates were quite smaller in July (~4.0 mm) and August (~53.8 mm) than in May (~125.2 mm) and June (~125.4 mm), in conjunction with a reduction of k d (490) in zone A, which could be revealed as a point source of wastewater. Unfortunately, the image from August, coincident with the BAR1 field campaign, has no information in the region highlighted as zone B that prevents showing the highest Chl-a concentrations in these sampling points, which generated high values of k d (490) as mentioned in Section 3.1. The mean was 1.5 m −1 , with values closer to the minimum, ranging from 1.3 to 1.7 m −1 , not corresponding with in situ k d (490). The precipitation rate in August was very low, resulting in a low k d (490) in zone A. The months of October, November, February, and March are representative of the wet season, however, in October, the mean k d (490) almost coincided with the minimum resulting in low values with some peaks northeast of the reservoir and in a winding region.
In November, the increase of k d (490) was observed with the mean closer to the maximum value, presenting some points above 2.5 m −1 . In February, with higher precipitation rates than the months previously mentioned, the k d (490) was around 3.0 m −1 even with approximately 30% of the reservoir covered by clouds. In March, the k d (490) was well distributed around the mean, ranging from 1.7 to 4.9 m −1 , with the highest values in the northeast and extremely north of the reservoir and values above 2.5 m −1 in some parts of zone B and in zone A.

Conclusions
After the application of several tests to estimate k d (490) through the available models developed to a wide range of waters, it was verified that empirical models do not estimate the k d (490) accurately in the Bariri reservoir. The uncertainty presented by the empirical models (~80% of MAPE) confirms that the blue/green wavelength ratio is not able to express the k d (490) variability in an environment dominated by CDOM, which was verified in the two surveys. The empirical models with Chl-a concentration as the input presented intermediary errors (~60% of MAPE) when applied to the aggregated data from the two surveys realized. In BAR1, the Chl-a concentration was significantly higher than in BAR2; the light absorption by CDOM was maintained predominantly in both surveys. Among the tested empirical models, the model developed by Zhang and Fell [29] was an exception in achieving lower errors and statistical results close to those obtained from the semi-analytical models.
The semi-analytical model developed by Lee et al. [9] presented the lowest error (MAPE of 41.0%) among all twelve models analyzed. The coefficients present in the Lee et al. [9] model were calibrated from a wide range of waters, filling the lack of the empirical models that are based on a specific dataset. However, the IOPs' estimates through the original QAA_v5 can introduce some errors in inland waters due to the optical complexity.
The Lee et al. [9] model was applied in atmospherically-corrected OLI/Landsat-8 images (LaSRC) with a relative error percentage difference of 8.5% in relation to the application to all datasets from the R rs obtained in situ. Therefore, the atmospheric correction was appropriate to retrieve the k d (490) and the Lee et al. [9] model was able to highlight the main variations of k d (490) in an environment dominated by CDOM in the Bariri reservoir, allowing the identification of the regions with more or less light attenuation and, consequently, the biota modifications that are dependent of the photic conditions in the water column.
The application of the k d (490) model on OLI/Landsat-8 images exhibited important temporal and spatial variability. The spatial variability was verified in points where the sinuosity of the reservoir promoted the reduction of water velocity and facilitated the OSS concentrations, increasing the attenuation of light. The temporal variability was linked to significant activities around the reservoir that increase the runoff and the input of external material into the reservoir mainly in periods of high precipitation rates. Additionally, the precipitation in the wet season promotes the resuspension and carriage of sediments, increasing the attenuation of light.
Author Contributions: A.C.C.G., N.B., A.C.d.C., and T.R. carried out the field campaigns and the laboratory analysis; E.A. supervised the experiment and managed the projects that funded the dataset acquisition and processing, as well as supervised the first author in their Master's degree; A.C.C.G. conducted the experiments and wrote the paper; and the other authors corrected the manuscript.