Suppression of local haze variations in MERIS images over turbid coastal waters for retrieval of suspended sediment concentration

Atmospheric correction over turbid waters can be problematic if atmospheric haze is spatially variable. In this case the retrieval of water quality is hampered by the fact that haze variations could be partly mistaken for variations in suspended sediment concentration (SSC). In this study we propose the suppression of local haze variations while leaving sediment variations intact. This is accomplished by a multispectral data projection (MDP) method based on a linear spectral mixing model, and applied prior to the actual standard atmospheric correction. In this linear model, the hazesediment spectral mixing was simulated by a coupled water-atmosphere radiative transfer (RT) model. As a result, local haze variations were largely suppressed and transformed into an approximately homogenous atmosphere over the MERIS top-of-atmosphere (TOA) radiance scene. The suppression of local haze variations increases the number of satellite images that are still suitable for standard atmospheric correction processing and subsequent water quality analysis. © 2010 Optical Society of America OCIS codes: (010.1285) Atmospheric correction; (010.1690) Atmospheric and oceanic opticsColor; (100.0100) Image processing References and links: 1. Y. J. Kaufman, “The atmospheric effect on remote sensing and its correction,” In Theory and applications of optical remote sensing, G. Asrar, ed., (Wiley, New York, 1998), pp.734. 2. J. Lavreau, “De-hazing Landsat Thematic Mapper images,” Photogramm. Eng. Remote Sensing 57, 1297–1302 (1991). 3. Y. Zhang, B. Guindon, and J. Cihlar, “An image transform to characterize and compensate for spatial variations in thin cloud contamination of Landsat images,” Remote Sens. Environ. 82(2-3), 173–187 (2002). 4. Y. Zhang, and B. Guindon, “Quantitative Assessment of a Haze Suppression Methodology for Satellite Imagery: Effect on Land Cover Classification Performance,” IEEE Trans. Geosci. Rem. Sens. 41(5), 1082–1089 (2003). 5. Z. Zhang, “Studies on basic characteristics of fine sediment in Changjiang Estuary,” J. Sediment. Res. 1, 67–73 (1996) (in Chinese with English abstract). 6. P. M. Teillet, N. T. O'Neill, A. Kalinauskas, D. Stugeon, and G. Fedosejevs, “A dynamic regression algorithm for incorporating atmospheric models into image correction procedures,” in Proceedings of the IGARSS'87 Symposium, Ann Arbor, MI, pp. 913−918 (1987). 7. R. J. Kauth, and G. S. Thomas, “The Tasseled Cap—a graphic description of the spectral-temporal development of the agricultural crops as seen by Landsat,” in Proceedings of the Symposium on Machine processing of Remotely Sensed Data, Purdue University, pp. 41−51 (1976). 8. R. Richter, “Atmospheric correction of satellite data with haze removal including a haze/clear transition region,” Comput. Geosci. 22(6), 675–681 (1996). 9. 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). 10. C. Hu, K. L. Carder, and F. E. Muller-Karger, “Atmospheric correction of SeaWiFS imagery over turbid coastal waters: A practical method,” Remote Sens. Environ. 74(2), 195–206 (2000). #123815 $15.00 USD Received 4 Feb 2010; revised 22 Mar 2010; accepted 13 May 2010; published 28 May 2010 (C) 2010 OSA 7 June 2010 / Vol. 18, No. 12 / OPTICS EXPRESS 12653 11. Th. Schroeder, I. Behnert, M. Schaale, J. Fischer, and R. Doerffer, “Atmospheric correction algorithm for MERIS above case-2 waters,” Int. J. Remote Sens. 28(7), 1469–1486 (2007). 12. C. Y. Ji, “Haze reduction from the visible bands of LANDSAT TM and ETM+ images over a shallow water reef environment,” Remote Sens. Environ. 112(4), 1773–1783 (2008). 13. European Space Agency, “MERIS user guide,” http://envisat.esa.int/handbooks/meris. 14. F. Shen, W. Verhoef, Y-X. Zhou, Mhd. S. Salama, and X–L. Liu, “Satellite estimates of wide-range suspended sediment concentrations in estuarine and coastal waters using MERIS data,” (submitted). 15. A. Berk, G. P. Anderson, P. K. Acharya, J. H. Chetwynd, L. S. Bernstein, E. P. Shettle, et al. MODTRAN4 USERS MANUAL, Air Force Research Laboratory, Space Vehicles Directorate, Air Force Materiel Command, Hanscom AFB, MA 01731–3010, pp.97 (2000).


Introduction
For ocean colour and water quality mapping one needs atmospherically clear satellite scenes, but these are seldom available.Haze, a kind of partly transparent clouds of varying optical thickness, arising from a variety of atmospheric constituents including water droplets, ice crystals or fog/smog particles [1], presents a major problem in remote sensing, not only because it masks the objects being imaged but also because they alter their spectral signatures [2] in the visible through near infrared spectral ranges [3,4].Contamination by spatially varying haze on satellite ocean scenes results in the wrong identification or overestimation of suspended particulate matters over turbid coastal waters, so as to primarily impair the accuracy of retrieval of water constituents and concentrations (e.g.concentrations of chlorophyll and suspended sediment which are principal products of ocean-colour observations).Terrestrial suspended particulate matters or resuspension from the bottom in shallow waters carried by river plume to the ocean, along with amounts of nutrients and pollutant materials, can impact oceanic environments, for instance, in the Changjiang (Yangtze River) Estuary and Coast.Such matters consist mainly of suspended finest-grain sediment particles with an average diameter of 8.6 µm [5], hereafter referred as "sediment".
Regarding correction of spatially varying haze, an extension of the image-based dark object subtraction method proposed by [6] is to partition a scene into sub-regions which are treated individually, the Tasseled-Cap Transforms [7], where the haze component (TC4) is applied pixel-by-pixel to estimate the radiometric contribution from haze, and correction of haze is made by subtracting the amount of additive shift for each pixel (e.g.see [2]), or by histogram matching (e.g.see [8]).The Haze Optimized Transform (HOT) reported by [3] is to detect haze and its effect in conjunction with a "dark target" subtraction methodology to radiometrically adjust visible-band imagery.These methods were applied to land scenes.
For ocean colour applications, SeaWiFS and MODIS standard algorithms for haze correction applied to oceanic waters (e.g.see [9]) and a few studies for turbid coastal waters (e.g.see [10,11]) are based on the assumption of spatially constant haze over the whole scene.For spatially varying haze, Ji (2008) reported that neither the Tasseled Cap haze component approach nor the HOT transform method is suited for haze removal in the case of shallow water marine environment [12].Therefore, the author proposed an image-based linear regression model for the haze correction.This method requires determining a proper threshold for the near-infrared (NIR) band in order to separate haze contaminated areas from haze-free areas by empirical operations [12].However, this may be not applicable to turbid coastal waters.The reason is that suspended sediment concentration (SSC) in the turbid waters, similar to haze levels, spatially varies and surface reflectance leaving the waters in the NIR spectral range is non-zero, so that it is difficult to find the proper threshold for the separation.
This study results from a pressing need for a methodology to effectively suppress the haze variations and to improve the accuracy of SSC retrieval in the presence of locally distributed haze over turbid waters.Such a capability would give a significant improvement on the accuracy of satellite-derived products.However, few publications on suppression of local varying haze mixed with sediment in turbid coastal waters have appeared so far.Since the spectral effects of haze variation in the atmosphere and sediment variation in the water on topof-atmosphere (TOA) radiances are different, in principle both effects can be separated, provided sufficient spectral information is available.The MERIS sensor has fifteen spectral bands in the VNIR (visible/near infrared) range and a good signal to noise ratio (see [13]), and therefore linear unmixing of haze and sediment components over turbid waters is considered feasible.A method called multispectral data projection (MDP), based on a linear spectral mixing model, is proposed to suppress the haze variations while leaving sediment variations intact, and which is applied prior to the actual atmospheric correction.The haze-sediment spectral mixing was simulated by a coupled water-atmosphere radiative transfer (RT) model.For turbid waters, we used a semi-empirical radiative transfer (SERT) model [14], and for the atmosphere we used MODTRAN [15] simulation results.The advantage of employing RT models is that one has full control over the selection of the haze and sediment end members, and these spectral basis vectors can be determined with high precision.Thereafter, the actual atmospheric correction can be conducted under the situation of spatially constant atmosphere and then the obtained water-leaving reflectance used as input for the SSC retrieval.

Multispectral data projection (MDP)
The effect of atmospheric haze on satellite imagery of ocean scenes is that, compared to a clear atmosphere, the radiance signal is increased by haze due to scattering by the aerosol particles.This increase is wavelength-dependent and, in turn, this wavelength-dependence (as expressed by the Ångström coefficient) can be related to the aerosol type (rural, urban, maritime, etc.).However, in all cases one may assume in a first approximation that the increase is linear with the amount of haze.In the following it is assumed that the aerosol type is constant over a scene, so that for the wavelength-dependence this applies as well.However, moderate spatial variations in aerosol optical thickness (AOT) may occur.
Similar considerations apply to sediment in ocean waters.An increase in sediment concentration results into an approximately linear increase of the radiance all over the spectrum, which is again wavelength-dependent.However, the wavelength-dependencies for increases of sediment and haze are not the same, and therefore these are in principle separable by a spectral unmixing technique.
The linear approximations hold only to a certain extent.That is why it will not be possible to remove clouds by this technique.But one may expect that moderate variations of AOT can be suppressed by correcting for the haze variations, in such a way that the effect of sediment variations over the scene stays intact.This is called the multispectral data projection (MDP) method.By means of model simulations with the atmospheric RT code MODTRAN [15] it can be shown that decreasing the atmospheric visibility (corresponding to an increasing AOT) in a number of steps gives an increase of the TOA radiance which is wavelength-dependent, but in a consistent spectral direction.
In the following it is assumed that TOA radiance spectra of ocean scenes over coastal waters can often be decomposed into two major constituents, atmospheric haze and sediment in the water.The linear mixing model applied for the combination haze-sediment can be formulated as follows: Suppose that the TOA spectral radiance vector of a pixel is approximated by The real TOA radiance spectrum of the pixel is supposed to be given by = + p p ε where ε is the residual error radiance spectrum.Suppression of haze variations is possible by transforming the TOA radiance spectrum to one which corresponds to the standard haze level.This spectrum is called ' p and it is created by subtracting the increased haze contribution, which is equal to 1 a h , so we have The operation to go from p to ' p is called a projection, since by removing the increased haze contribution; the data are projected onto the s vector along the direction of the h vector.For a two-dimensional case this is illustrated in Fig. 1.In order to enable this transformation, the actual haze amount in arbitrary spectra p , given by the coefficient 1 a , has to be estimated.For this, in Eq. ( 1) p is replaced by p , and, after subtracting the reference vector r , multiplied by the transposed vectors of h and s , which yields a a a a

− = + − = + h p r h h h s s p r s h s s
From this, one finds both coefficients from

h h h s h p r h p r s h s s s p r s p r
One can write where I is the identity matrix.The matrix T − I hu in Eq. ( 2) contains the coefficients applied to the radiance data in all bands, and T u rh is an offset vector.For the image processing it is important to establish which quantities are fixed and which ones are spatially variable.Of the right-hand side of Eq. ( 2), only p is spatially variable, and the other quantities are only dependent on the end members selected.However, Eq. ( 2) involves a matrix multiplication, which might be timeconsuming.Therefore, it is probably more efficient to perform the two-step operation In this case, no matrices are involved, only vectors.The two-step algorithm can be summarized as follows: 1.The reference spectrum r is subtracted from the pixel's radiance spectrum p .A vector inner product is formed by means of the transposed of vector u , which results into the image of scalar 1 a , expressing the spatial haze variations.
2. The image of the haze variations magnitude expressed in 1 a is multiplied by the haze variation spectrum h and the resulting spectrum is subtracted from the original pixel's spectrum p .In this second step the haze variations are removed while sediment variations remain intact.Regarding the error of the linear unmixing approximation, it should be mentioned that in the two-dimensional case the decomposition into haze and sediment components is exact, i.e. there will be no difference between the actual and the modelled TOA radiance spectrum.This is because in a two-dimensional space a given vector can always be decomposed into two basis-vectors.In a multidimensional space the approximation according Eq. (1) will not be perfect, but only the minimum sum of squared errors solution.In section 3.1, Fig. 5(a) gives a visual impression of how accurate the approximation is for MERIS data in the spectral range from 400 to 800 nm.This can also be regarded as a first step in the validation of the algorithm.Note that the error can only be estimated from model-generated results, and that it also depends on the choice of end member spectra.One may expect that the approximation will only be fair in the spectral space spanned by the reference spectral vector and the chosen vectors expressing the haze and sediment variations.However, investigating these aspects was considered to fall outside the scope of the present paper.

Scheme for suppression of haze variations
Ocean colour remote sensing usually requires two major steps.First, atmospheric correction has to be conducted for ocean scenes to obtain surface reflectance spectra leaving the water, because the atmosphere, in most cases, contributes more than 90% of TOA radiance obtained by sensors.Second, water constituents and concentrations have to be derived from the atmospherically corrected reflectance spectra.The first step is extremely critical.Usual schemes of atmospheric correction basically take the assumption of a spatially constant haze over the whole scene.To reduce the problem of local haze variations over the ocean scene, suppression of the haze variations by the MDP transform of the TOA radiance prior to the actual atmospheric correction is proposed here.The scheme for the processing procedure is shown in Fig. 2. The reference spectrum was taken to be the one for 40-km visibility and a minimum SSC level of 10 mg l −1 .The increased haze spectrum was represented by a 10-km visibility and the increased sediment spectrum by a SSC level of 50 mg l −1 .These three basis spectra (end members) were selected to represent the spectral effects of haze and sediment increases, relative to a reference situation of good visibility (i.e.40 km) and low SSC (i.e. 10 mg l −1 ) and shown in Fig. 4. Transforming all 36 spectra (shown in Fig. 3) by the proposed MDP method gives a result as displayed in Fig. 5(a).This may be compared to the subset of 9 spectra corresponding to a standard level of haze (i.e., in the case of 40-km visibility), which is given in Fig. 5(b).Comparison reveals that the radiance spectra (shown in Fig. 5(a)) resulting from the MDP transform resemble quite well those for the real radiance case (shown in Fig. 5(b)) with the standard level of haze.This implies that the proposed MDP method can effectively suppress the haze variations on the spectra of TOA radiance.All 36 modelled spectra (shown in Fig. 3) are transformed by the MDP method and shown in (a).A subset of 9 real spectra at 40 km visibility is shown in (b) for comparison.

Suppression of haze variations in MERIS TOA radiance images
As an example, the MERIS TOA radiance scene of 23 June 2005, displaying local haze variations over the Yangtze delta and the adjacent coasts, shown in Fig. 6(a), was taken.
Varying concentration levels of haze and sediment components occurred in the left part of the scene (see Fig. 6(c) and the ellipse area in Fig. 6(a)).After the performance of the MDP transform, the local haze variations are pronouncedly suppressed on the MERIS TOA radiance We took a region of interest on the scene (i.e. the ellipse area shown in Fig. 6(a)) and generated a scatter plot of radiance at MERIS band 442 nm verses 620 nm within the area (see Fig. 7 shown in the (b).

Retrieval of suspended sediment concentrations
Since local haze variations on the ocean TOA radiance scene have been suppressed, the usual processing approach of atmospheric correction can simply be performed under the condition of a homogenous atmosphere.This part is not studied here.A more important purpose was to obtain mappings of water constituents and concentrations.The study on SSC retrieval from MERIS ocean scenes over turbid waters of the Changjiang Estuary and adjacent coasts (see [14]) indicated that the atmospherically corrected water-leaving reflectance at shorter wavelength (e.g.560 or 620 nm) could be appropriate for the retrieval of low SSC, and longer wavelength (e.g.708 or 778 nm) should be used for the retrieval of high SSC.
Figure 8 shows SSC retrieved from the water-leaving reflectance of the MERIS ocean scene, by using the SERT model coupled with a multi-conditional algorithm scheme that has been established in the study of [14].Figure 8(b) shows the MERIS-derived SSC that resulted from applying the proposed scheme of suppression of local haze variations (termed as hazesuppressed SSC).By contrast, Fig. 8(a) shows the MERIS-derived SSC resulting from usual processing procedures, i.e. without applying suppression of local haze variations (termed as haze-unsuppressed SSC).It is evident that wrong identifications of sediment occur in the two dashed-circled areas indicated in Fig. 8(a), caused by the increased reflectance over the clear water due to the contamination by local haze.Moreover, the concentrations of sediment shown in Fig. 8  which was largely overestimated.The haze-suppressed SSC at the Pin 1 point (see Fig. 8(b)) was 0.1 g l −1 , which approximates the real SSC.Furthermore, the Gehu Lake (e.g.Pin 2 point in Fig. 8(b)) was investigated and the water stays clear there with SSC below 0.02 g l −1 .Obviously, the haze-unsuppressed SSC at the Pin 2 point (above 0.3 g l −1 in Fig. 8(a)) was greatly overestimated too.The validation of the results from actual atmospheric correction and the algorithm of SSC retrieval in the study area was discussed in detail in Ref [14].

Conclusion
It can be concluded that the spectral projection method transforms given TOA radiance spectra into spectra which correspond to a standard atmospheric condition with a visibility of 40 km.These transformed spectra resemble those for the real 40-km visibility case quite well, and multispectral images of similar data can next be used as input for a simple and standard atmospheric correction.
The present results are partly dependent on the choice of basis vectors.It should also be possible to apply the method with varying basis vectors by using the data from the look-up table in a more flexible way.

1 a
where r = reference spectrum for standard haze level and zero sediment h = spectrum for increased haze level and zero sediment, relative to r s = spectrum for standard haze level and increased sediment, relative to r p = modelled TOA radiance spectrum #123815 -$15.00USD Received 4 Feb 2010; revised 22 Mar 2010; accepted 13 May 2010; published 28 May 2010 (C) 2010 OSA 7 June 2010 / Vol. 18, No. 12 / OPTICS EXPRESS 12655 = amount of haze above the standard level in arbitrary units 2 a = amount of sediment in arbitrary units

Fig. 1 .
Fig. 1.Two-dimensional vector diagram on the modelling of haze and sediment variations by linear spectral mixing and the suppression of haze variations by subtracting the increased haze contribution.

Fig. 3 .
Fig. 3. Model-simulated TOA radiance spectra at MERIS bands for 36 combinations of haze and sediment

Fig. 4 .
Fig. 4. Three basis spectra (end members) at MERIS bands for the MDP method

Fig. 5 .
Fig. 5. TOA radiance spectra at MERIS bands at standard level of haze with a 40-km visibility.

Fig. 6 .
Fig. 6.MERIS scene of 23 June 2005 over the Changjiang Estuary and adjacent coast.Original TOA radiance scene at MERIS 620 nm is shown in (a).The scene where spatial variations of haze were suppressed by the MDP transform is shown in (b).The spatially varying haze on the scene is shown in (c).scene(shown in Fig.6(b)) which has been transformed to the one at standard level of haze with 40-km visibility.Visually, the scene is clear and haze appears uniform (Fig.6(b)).We took a region of interest on the scene (i.e. the ellipse area shown in Fig.6(a)) and generated a scatter plot of radiance at MERIS band 442 nm verses 620 nm within the area (see Fig.7(a)).It is shown in Fig.7(a) that the original haze and sediment components in various

Fig. 7 .
Fig. 7. Scatter plot of original TOA radiance at MERIS 442 nm vs. the one at 620 nm in the region of interest (red dashed-circled area in Fig. 6(a)) shown in the (a) and of transformed by the MDP method in the same region of interest (green solid-circled area in Fig. 6(b)) (a) were overestimated, due to the same reason.This can be proved by comparing retrieved SSC to real SSC values.The real SSC measured at the South Branch (e.g.Pin 1 point in Fig. 8(b)) on 23 June 2005 was 0.09 g l −1 , coincident with the day of MERIS sampling.The haze-unsuppressed SSC (see Fig. 8(a)) at the Pin 1 point was above 0.4 g l −1 , #123815 -$15.00USD Received 4 Feb 2010; revised 22 Mar 2010; accepted 13 May 2010; published 28 May 2010

Fig. 8 .
Fig. 8. Mappings of SSC retrieved from MERIS TOA radiance scene of 23 June 2005.The MERIS haze-unsuppressed SSC is shown in (a) and the haze-suppressed SSC in (b).