Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Estimation of near-infrared water-leaving reflectance for satellite ocean color data processing

Open Access Open Access

Abstract

The atmospheric correction algorithm employed by the NASA Ocean Biology Processing Group requires an assumption of negligible water-leaving reflectance in the near-infrared region of the spectrum. For waters where this assumption is not valid, an optical model is used to estimate near-infrared water-leaving reflectance. We describe this optical model as implemented for the sixth reprocessing of the SeaWiFS mission-long time-series (September 2009). Application of the optical model resulted in significant reductions in the number of negative water-leaving reflectance retrievals in turbid and optically complex waters, and improved agreement with in situ chlorophyll-a observations. The incidence of negative water-leaving reflectance retrievals at 412 nm was reduced by 40%, while negative reflectance at 490 nm was nearly eliminated.

©2010 Optical Society of America

1. Introduction

The retrieval of accurate geophysical data products (e.g., spectral remote sensing reflectance Rrs(λ), and concentrations of the phytoplankton pigment chlorophyll-a, Ca) from satellite-borne radiometers requires the use of an ‘atmospheric correction’ algorithm [1]. A critical component of this atmospheric correction algorithm is estimation of spectral aerosol reflectance, ρa(λ), which are subtracted from the total measured signal as part of the process of determining Rrs(λ). The algorithm implemented by the NASA Ocean Biology Processing Group (OBPG) [1] requires an assumption of negligible water-leaving reflectance in the near infrared (NIR) region of the spectrum (i.e., Rrs(NIR) = 0 sr−1). The NIR bands of the NASA Sea-viewing Wide Field-of-view Sensor (SeaWiFS) are centered on 765 and 865 nm. With this ‘black-pixel’ assumption, the measured top-of-atmosphere (TOA) reflectance in two NIR bands can be used to estimate both the magnitude and spectral dependence of ρa(NIR). Unfortunately, this assumption is rarely valid for waters with significant particulate (e.g., algal and mineral) backscattering [2]. If Rrs(NIR) can be modeled, however, it can be removed from the TOA signal, allowing the atmospheric correction to proceed in estimating ρa(NIR) as if Rrs(NIR) were negligible.

The OBPG implemented a bio-optical model to estimate Rrs(NIR) as part of the third reprocessing of the SeaWiFS mission in May 2000 (R2000) [2,3]. This approach required an initial run through the atmospheric correction process (ignoring Rrs(NIR) > 0 sr−1) to estimate Ca. This preliminary Ca was used to estimate spectral particulate backscattering, bbp(NIR), which was used in turn to reconstruct Rrs(NIR). With this modeled Rrs(NIR) removed from the TOA signal, the atmospheric correction process was repeated, Ca was recalculated, the process was iterated upon until convergence in Ca was achieved [2,4]. With the fourth SeaWiFS reprocessing in July 2002 (R2002), this Ca-driven algorithm was replaced with a reflectance-based algorithm [5]. The OBPG switched from [2] to [5] after observing that [2] depressed ratios of ρa(NIR) in optically complex waters. This often resulted in the selection of spectrally flat aerosol models within the atmospheric correction process, which artificially depressed the final Ca retrievals. The initial estimates of Ca were unrealistically high in these complex water masses, which led to over-estimates of bbp(NIR) and, subsequently, Rrs(NIR). The reflectance-based approach [5] was less prone to such over-correction. In this revised model, Rrs(λ) from an initial run through the atmospheric correction process and empirical estimates of the absorption coefficients for particles and dissolved materials at 670-nm were used to estimate bbp(670) via an ocean surface reflectance model [6]. A spectral scattering function [7] was then used to derive bbp(NIR) from bbp(670). Finally, Rrs(NIR) was reconstructed from bbp(NIR), the atmospheric correction process was repeated, and the process was iterated upon until convergence in Rrs(NIR) was achieved.

Despite such efforts to account for non-negligible Rrs(NIR), the atmospheric correction algorithm continued to demonstrate sub-optimal performance in highly scattering waters, as manifest by negative Rrs(λ), particularly in the blue region of the spectrum, and erroneous spectral dependence in visible (VIS) part of the spectrum, with the latter leading to inaccurate Ca retrievals. To address these issues, a re-evaluation of the implementation of [5] was initiated prior to the sixth SeaWiFS reprocessing in September 2009 (R2009) and several modifications were made to the Rrs(NIR) model, which are presented below. First, an inconsistency in the relationship between Rrs(λ) and the inherent optical properties (IOP) was eliminated. Second, alternative estimates of the spectral dependence of bbp(λ) and of the absorption coefficients for particles and dissolved materials at 670 nm were adopted. Finally, a revised iteration scheme was implemented.

2. NIR reflectance model development

2.1 Inversion of remote-sensing reflectance

Remote sensing reflectance is proportional to the IOPs, specifically the marine absorption and backscattering coefficients, through:

Rrs(λ)=G(λ)X(λ)
where

X(λ)=bb(λ)a(λ)+bb(λ),

G(λ) is a function of the illumination conditions, water constituents, and sea-state, bb(λ) is the total backscattering coefficient [ bbw(λ) + bbp(λ)], and a(λ) is the total absorption coefficient [ aw(λ) + ap(λ) + ag(λ)]. The additional subscripts w, p, and g indicate specific contributions by water, particles, and dissolved material, respectively.

The NIR reflectance model [5] implemented for R2002 and used through SeaWiFS reprocessing 5.2 in 2007 (R2007) suffered from an inconsistent application of Eq. (1). First, an approximation of this relationship, where G(λ) was fixed at 0.51 and X(λ) was simplified to bb(λ) / a(λ), was used to obtain near-infrared IOPs from Rrs(VIS). The quadratic form of [6], where G(λ) = [g1 + g2 X(λ)], g1 = 0.0949, and g2 = 0.0794, was later used to obtain Rrs(NIR) from the modeled near-infrared IOPs. While a minor inconsistency, using dissimilar definitions for G(λ) and X(λ) was found to have significant effects in highly scattering waters.

Furthermore, the approach was inconsistent with the methods of [9] that are used elsewhere (and independently) in the atmospheric correction process to normalize Rrs(λ) for bi-directional effects. In [9], G(λ) = f(λ) / Q(λ), where f(λ) is a function of the illumination conditions, water constituents, and sea-state, Q(λ) is the ratio of spectral upwelling irradiance to spectral upwelling radiance, and ratios of f(λ) to Q(λ) are pre-tabulated from extensive radiative transfer simulations. For the Rrs(NIR) model implemented in R2009, G(λ) = f(λ) / Q(λ), using the look up tables (LUTs) of [9], for both the conversion of retrieved Rrs(VIS) to IOPs and the conversion of modeled IOPs to Rrs(NIR). While these LUTs were originally generated for waters whose light field is predominantly influenced by phytoplankton, they proved suitable for use in this model, as will be demonstrated in Section 4. A full discussion of the applicability of [9] to all water types is beyond the scope of this paper.

2.2 Spectral backscattering

Given the relationship in Eq. (1), a measured Rrs0), and knowledge of a(λ0), an estimate of bb0) can be obtained, where λ0 is any reference wavelength (e.g., 670 nm). This retrieved bb0) can then be extrapolated to the NIR through a spectral model. The spectral dependence of bb(λ) can be expressed as:

bb(λ)=bb(λ0)Y,
bb(λ)=bbw(λ)+bbp(λ0)Y,
where Y can take a number of forms, but is commonly parameterized as a power law function:

Y=[λ0λ]η.

In Eq. (2)b, the η term is typically in the range of [0, 2], with values around 2 observed in open ocean waters and turbid waters approaching zero [12]. In contrast, the model employed in R2002-R2007 [3] used a spectral scattering relationship for Eq. (2)a developed by [7]:

Y=(αλ+β)/(αλ0+β),
where the empirically derived coefficients α and β are 0.00113 and 1.62517, respectively. For use in the Rrs(NIR) model, the latter relationship required an assumption of equivalent spectral dependency for backscattering and scattering. In the Rrs(NIR) model employed in R2009, the more common power-law form of Eq. (2)b was adopted. Rather than relying on a fixed η value, however, the empirical relationship developed in [12] was adopted:

η=2 .0 * [1 . - 1.2 * e(-0 .9*Rrs(443)/Rrs(555))].

The NASA bio-Optical Marine Algorithm Data set (NOMAD) [13] includes measurements of bb(VIS) at a number of wavelengths spanning a range of values from 0.0007 to 0.0133 m−1. Unfortunately, this data set does not include measured bb(NIR). Nevertheless, this data set was used to evaluate the various forms of Eq. (2) and 3 using a reference wavelength of 670 nm and a retrieval wavelength of 555 nm. The formulation of [12] had the minimal bias and median percent difference relative to in situ measurements (Table 1 , Fig. 1 ). and performed better than the linear, spectral scattering function [7].

Tables Icon

Table 1. Statistics of Modeled versus Measured Backscattering at 443 nm

 figure: Fig. 1

Fig. 1 Modeled vs. measured bb(555 nm) from the NOMAD data set. Gray symbols are the power law formulation of [12]. Red symbols are the linear function of [7].

Download Full Size | PDF

2.3 Estimating absorption

In the NIR, a(λ) is dominated by the absorption due to water, aw(λ) [10,11], such that absorption coefficients for particles, ap(NIR), and dissolved material, ag(NIR), can be ignored. However, a reasonable estimation of the absorption coefficient is crucial for the inversion of Rrs(λ) to retrieve bb(670), where ap(670) and ag(670) can be a significant portion of the total. In [5], the estimate of a(670) involved two functions, one based on Ca [14] and the other on a green-to-red reflectance ratio [5]. In the modified Rrs(NIR) model, this was updated to a single, Ca-based relationship derived from NOMAD [13] (Type II linear regression using spectroscopy only, limited to Ca > 0.2 mg m−3; Fig. 2 ):

a(670)=e(ln(Ca)0.93893.7589)+aw(670)
The ap(670) and ag(670) values contained within the NOMAD data set span a range of 0.00001–0.03 m−1 for ap(670) and 0.00057 0.8413 m−1 for ag(670).

 figure: Fig. 2

Fig. 2 Particulate and dissolved absorption at 670 nm vs. Ca derived from the NOMAD data set.

Download Full Size | PDF

3. Iteration scheme

Since the retrieval of Rrs(λ) in highly scattering waters requires knowledge of Rrs(NIR), and estimation of Rrs(NIR) requires knowledge of Rrs(VIS), the Rrs(NIR) algorithm must be implemented iteratively. The iteration begins by assuming Rrs(NIR) = 0 sr−1, so that initial Rrs(λ) can be derived and a first estimate of Rrs(NIR) can be achieved. The iteration converges when Rrs(λ) changes by less than 2%. In the majority of cases, the method converges after 3-4 iterations, but the algorithm allows up to 10 iterations.

If, however, the initial atmospheric correction results in negative or an otherwise non-physical spectral distribution for Rrs(λ) (criteria based on in situ observations), the iteration is re-initialized assuming ρa(NIR) = 0 (i.e., all reflectance except that from Rayleigh scattering and Sun glint originates from the water mass). This is effectively the opposite extreme from the initial condition of Rrs(NIR) = 0 sr−1. The iteration process is then allowed to proceed to convergence as previously described. For the rare cases where the iteration fails to converge after this reset, an additional iteration is forced, again with an assumption of ρa(NIR) = 0. This pixel is then flagged with an ‘atmospheric correction warning’ to alert users to the questionable nature of the pixel, although a qualitatively useful retrieval may still result.

The controls on this iteration scheme differ slightly from previous implementations [8], especially with regards to the how the iteration is reset, which is currently more robust than that of past methods. It should also be noted that the Rrs(NIR) is forced to zero if the initial iteration results in a Ca < 0.3 mg m−3. This is done to avoid the addition of noise into clear-water pixels through the estimation of otherwise negligible Rrs(NIR) contributions. Furthermore, to avoid the introduction of transitional artifacts, the Rrs(NIR) estimation is linearly weighted from 0 to 1 for 0.3 ≤ Ca ≤ 0.7 mg m−3 (Fig. 3 ). This Ca-based weighting of Rrs(NIR) was first introduced in R2002 and is effectively unchanged for R2009.

 figure: Fig. 3

Fig. 3 Global map of likely application of the NIR correction. This was created from the SeaWiFS mission cumulative climatology. Black indicates land; grey Ca < = 0.3 mg m3; white Ca > 0.3 mg m3. The white area indicates where the NIR correction is likely to be applied.

Download Full Size | PDF

4. Impact

A long-term satellite and in situ data set was recently compiled to study Ca variability in Chesapeake Bay [15]. These data were used to evaluate the modifications to the NIR correction described above, as implemented for the most recent SeaWiFS reprocessing. Specifically, the NIR correction implemented with R2002 (and used through R2007) was compared to the method described in this work and employed for R2009. The satellite and in situ data acquisition and processing methods are described in [15].

The revised correction reduced the percentage of negative Rrs(412) for the time series of SeaWiFS retrievals over the Bay by nearly half (18.5% vs. 31.9%). The incidence of negative Rrs(490) was reduced to negligible levels (Table 2 ). The satellite estimates of Ca showed a substantial improvement in the agreement with ground truth measurements in the turbid and highly productive regions of the Middle and Upper Bay (Fig. 4 ).The revised correction also improved reflectance retrievals relative to coincident in situ measurements. Following the method of [16], and using a global validation data set compiled by the OBPG, we compared Rrs(λ) retrievals for SeaWiFS after limiting the data set to where the NIR correction was used in the processing of the satellite data (e.g., Ca ≥ 0.3 mg m−3). Improvement in the satellite agreement with in situ Rrs was evident for all bands, with the greatest improvement seen for the 412 and 443 nm bands (Table 3 ).

Tables Icon

Table 2. Percent Negative Remote Sensing Reflectance

 figure: Fig. 4

Fig. 4 Histograms of in situ and satellite derived Ca for the Lower, Middle and Upper Chesapeake Bay

Download Full Size | PDF

Tables Icon

Table 3. Validation Statistics for SeaWiFS Rrs(λ)

5. Conclusion

The revised Rrs(NIR) model that is presented in this paper has been implemented in the operational processing of satellite ocean color sensor data by the NASA Ocean Biology Processing Group. This revised correction eliminated a number of inconsistencies present in the previous version. A variable estimate for the spectral slope of the backscattering coefficient was implemented in place of a fixed linear relationship. A simplified estimate of absorption in the red region of the visible spectrum was also implemented. This resulted in a significant reduction in the incidence of negative Rrs for satellite retrievals in the 412 to 490 nm range, while improving the agreement of satellite-derived Rrs(λ) and Ca with in situ measurements.

The results and discussion presented here have been focused on ocean color retrievals from the SeaWiFS sensor. It should be noted, however, that the same algorithm has been applied in the processing of MODIS-Aqua and other ocean color sensor data sets archived by NASA’s Ocean Biology Processing Group.

The approach presented, while a significant improvement over the previous implementation, cannot correct all cases where the black pixel assumption is invalid. The use of empirical optical models (e.g. spectral backscattering, particulate and dissolved absorption) limits the applicability of this algorithm to waters that are similar to those over which these empirical models were developed. In addition, the practical implementation of the iteration scheme can result in conditions of non-convergence. Fortunately, the optical domains covered by the data sets used in the generation of these algorithms are quite diverse and the cases of convergence failure are small. Therefore, these limitations do not detract from the benefits realized by the use of this algorithm.

References and Links

1. H. R. Gordon and M. Wang, “Retrieval of water-leaving radiance and aerosol optical thickness over the oceans with SeaWiFS: a preliminary algorithm,” Appl. Opt. 33(3), 443–452 (1994). [CrossRef]   [PubMed]  

2. D. A. Siegel, M. Wang, S. Maritorena, and W. Robinson, “Atmospheric correction of satellite ocean color imagery: the black pixel assumption,” Appl. Opt. 39(21), 3582–3591 (2000). [CrossRef]  

3. C. R. McClain, E. Ainsworth, R. Barnes, R. Eplee, Jr., F. Patt, W. Robinson, M. Wang, and S. W. Bailey, “SeaWiFS Postlaunch Calibration and Validation Analyses, Part 1,” NASA Tech. Memo. 206892, National Aeronautics and Space Administration, Goddard Space Flight Center, Greenbelt, MD (2000).

4. A. Morel, “Optical modeling of the upper ocean in relation to its biogenous matter content (case 1 waters),” J. Geophys. Res. 93(C9), 10749–10768 (1988). [CrossRef]  

5. R. P. Stumpf, R. A. Arnone, J. R. W. Gould, P. M. Martinolich, and V. Ransibrahmanakul, “A partially coupled ocean-atmosphere model for retrieval of water-leaving radiance from SeaWiFS in coastal waters,” in Patt, F.S., et al., 2003: Algorithm Updates for the Fourth SeaWiFS Data Reprocessing. NASA Tech. Memo. 206892, National Aeronautics and Space Administration, Goddard Space Flight Center, Greenbelt, MD (2003).

6. H. R. Gordon, O. B. Brown, R. H. Evans, J. W. Brown, R. C. Smith, K. S. Baker, and D. K. Clark, “A semianalytic radiance model of ocean color,” J. Geophys. Res. 93(D9), 10909–10924 (1988). [CrossRef]  

7. R. W. Gould Jr, R. A. Arnone, and P. M. Martinolich, “Spectral dependence of the scattering coefficient in case 1 and case 2 waters,” Appl. Opt. 38(12), 2377–2383 (1999). [CrossRef]  

8. NIR, html, NIR Correction, http://oceancolor.gsfc.nasa.gov/REPROCESSING/SeaWiFS/R4/NIR.html

9. A. Morel, D. Antoine, and B. Gentili, “Bidirectional reflectance of oceanic waters: accounting for Raman emission and varying particle scattering phase function,” Appl. Opt. 41(30), 6289–6306 (2002). [CrossRef]   [PubMed]  

10. L. Kou, D. Labrie, and P. Chylek, “Refractive indices of water and ice in the 0.65-2.5 m spectral range,” Appl. Opt. 32(19), 3531–3540 (1993). [CrossRef]   [PubMed]  

11. R. M. Pope and E. S. Fry, “Absorption spectrum (380-700 nm) of pure water. II. Integrating cavity measurements,” Appl. Opt. 36(33), 8710–8723 (1997). [CrossRef]  

12. Z. Lee, R. A. Arnone, C. Hu, P. J. Werdell, and B. Lubac, “Uncertainties of optical parameters and their propagations in an analytical ocean color inversion algorithm,” Appl. Opt. 49(3), 369–381 (2010). [CrossRef]   [PubMed]  

13. P. J. Werdell and S. W. Bailey, “An improved in-situ bio-optical data set for ocean color algorithm development and satellite data product validation,” Remote Sens. Environ. 98(1), 122–140 (2005). [CrossRef]  

14. A. Bricaud, A. Morel, M. Babin, K. Allali, and H. Claustre, “Variations of light absorption by suspended particles with the chlorophyll a concentration in oceanic (case 1) waters: analysis and implications for bio-optical models,” J. Geophys. Res. 103(C13), 31033–31044 (1998). [CrossRef]  

15. P. J. Werdell, S. W. Bailey, B. A. Franz, L. W. Harding Jr, G. C. Feldman, and C. R. McClain, “Regional and seasonal variability of chlorophyll-a in Chesapeake Bay as observed by SeaWiFS and MODIS- Aqua,” Remote Sens. Environ. 113(6), 1319–1330 (2009). [CrossRef]  

16. S. W. Bailey and P. J. Werdell, “A multi-sensor approach for the on- orbit validation of ocean color satellite data products,” Remote Sens. Environ. 102(1-2), 12–23 (2006). [CrossRef]  

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (4)

Fig. 1
Fig. 1 Modeled vs. measured bb(555 nm) from the NOMAD data set. Gray symbols are the power law formulation of [12]. Red symbols are the linear function of [7].
Fig. 2
Fig. 2 Particulate and dissolved absorption at 670 nm vs. Ca derived from the NOMAD data set.
Fig. 3
Fig. 3 Global map of likely application of the NIR correction. This was created from the SeaWiFS mission cumulative climatology. Black indicates land; grey Ca < = 0.3 mg m3; white Ca > 0.3 mg m3. The white area indicates where the NIR correction is likely to be applied.
Fig. 4
Fig. 4 Histograms of in situ and satellite derived Ca for the Lower, Middle and Upper Chesapeake Bay

Tables (3)

Tables Icon

Table 1 Statistics of Modeled versus Measured Backscattering at 443 nm

Tables Icon

Table 2 Percent Negative Remote Sensing Reflectance

Tables Icon

Table 3 Validation Statistics for SeaWiFS Rrs(λ)

Equations (8)

Equations on this page are rendered with MathJax. Learn more.

R r s ( λ ) = G ( λ ) X ( λ )
X ( λ ) = b b ( λ ) a ( λ ) + b b ( λ )
b b ( λ ) = b b ( λ 0 ) Y ,
b b ( λ ) = b b w ( λ ) + b b p ( λ 0 ) Y ,
Y = [ λ 0 λ ] η .
Y = ( α λ + β ) / ( α λ 0 + β ) ,
η = 2 .0 *  [ 1 . - 1 .2 * e ( -0 .9* R rs ( 443 ) / R rs ( 555 ) ) ] .
a ( 670 ) = e ( ln ( C a ) 0.9389 3.7589 ) + a w ( 670 )
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.