[0001] 本发明涉及遥感技术领域，尤其涉及一种富营养化湖泊水体叶绿素a的MODIS卫星高精度监测方法。

[0002] 遥感技术广泛地应用于监测水体富营养化和藻华，大大改善了水质的实时监测能力(Matthews,2010；Odermatt et al.,2012)。叶绿素a(Chla,μg/L)是藻华强度和水质营养状况的重要标志，成为水质遥感最普遍的参数，(Chaeffer et al.,2012；Matthews,2010；Qi et al.,2014)。近年来，由于频繁出现的藻华，对饮用水和灌溉用水水质产生了严重的威胁(Duan et al.,2009；Guo,2007)，利用遥感监测水体叶绿素a迅速成为研究热点(Paerl & Huisman,2008；Paerl & Huisman,2009；Matthews,2010；Odermatt et al.,2012)。然而，往往出现类似问题的水体通常都是光学特性复杂的二类水(Morel & Prieur,1977)，因此，对于光学特性复杂的内陆和沿海水域进行精准的水色遥感反演(特别是Chla)一直是一个严峻的挑战(IOCCG,2000)。

[0003] 通过全面回顾关于二类水体的文献，四个特性通常被用来评估Chla浓度：色素在442nm波段处吸收最大值(Bricaud et al.,1995；Gitelson et al.,1992)；色素在665nm波段处第二吸收最大值(Bricaud et al.,1995)；572nm波段附近由于色素吸收最小值和SPM的散射而产生的反射峰(Gitelson et al.,1992；Schalles et al.,1998)；700nm波段附近由于浮游植物、SPM、CDOM和纯水吸收的总值最小而形成的反射峰(Gitelson et al.,1992；Vasilkov & Kopelevich,1982；Vos et al.,1986)；以及685nm波段附近的荧光峰(Gitelson et al.,1994；Gower,1980；Gower et al.,1999)。基于第一特性，OC2,OC3和OC4算法得到发展并且可以应用于一类水以及低Chla浓度的二类水(D'Sa & Miller,2003；Horion et al.,2010；Mélin et al.,2007；Witter et al.,2009)，因为所有的中等分辨率海洋水色光谱分析仪都能获取Chla的第一特征的吸收波段(Matthews,2010)。但是OC2-OC4算法只适用于反演Chla浓度小于10μg/L的光学特性复杂水体，无法了应用于富营养化水体。近年来有些算法发展尝试使用了第二和第三特性(绿光和红光波段)，并在改善Chla数据结果精确性方面取得了重大的进步(Dall'Olmo & Gitelson,2005；Gitelson et al.,2008；Le et al.,2013；Le et al.,2009；Shen et al.,2010；Tassan & Ferrari,2003；Thiemann & Kaufmann,2000)。然而，除了shen和le的研究，大多数这方面研究都是基于实地测量的光谱数据(Le et al.,2013；Shen et al.,2010)。目前，基于700nm左右和670nm左右的波段比值算法被广泛应用于评估富营养化水体的Chla浓度(Duan et al.,2007；Gitelson et al.,1993；Gons,1999；Moses & Gitlson,2009)。Chla与700/670nm比值之间的相关性主要来源于叶绿素a光谱的第四个特性。后来，这个比值算法又演变和发展为“三波段算法”(Dall'Olmo & Gitelson,2005；
Duan et al.,2010；Gitelson et al.,2008；Moses et al.,2009；Zimba & Gitelson,2006)和“四波段算法”，这些算法大大提高了浑浊水体中Chla浓度的反演精度(Le et al.,2013；
Le et al.,2009)。FLH(Dierberg & Carriker,1994；Giardino et al.,2005)、MCI(Binding et al.,2011；Gower et al.,2005)以及MPH(Matthews et al.,2012)则是基于685nm波段附近的荧光性最大值的线性基线算法，FLH适用于浓度小于30μg/L的水体，其他两个适用于浓度小于100μg/L的水体(Gower et al.,2005；Matthews et al.,2012)。除了上面提到的Chla算法外，人工神经网络算法(Keiner,1999)，多元回归分析(Tyler et al.,2006)，EOF(Craig et al.,2005；Qi et al.,2014)，SDA(Oyama et al.,2010；Oyama et al.,2009)，以及LUT方法(Yang et al.,2011)也被应用于Chla的估算。尽管这些算法比其他算法改善了部分相关性，但是这些改善与其复杂的数学过程而言则显得微不足道(Matthews,2010)。

[0004] 通常，遥感仪器所测量的可见光和近红外光谱(400—1000nm)部分，通常被用于监测水体。因此，典型的海洋水色传感器SeaWiFS、MODIS和MERIS更适用于监测水色参数。从Chla遥感反演的文献可知，SeaWiFS和MODIS可用于OC算法和低Chla浓度，MERIS可通过红光和近红外波段比值、MCI算法和NN算法等可用于高浓度Chla反演(Odermatt et al.,2012)。MERIS在测量浑浊富营养糊水体的Chla方面比其他传感器有明显的优势，但是其16天的测量周期却差强人意。相比之下，MODIS有1-2天的高时间分辨率以及比SeaWiFS更好的空间分辨率(MODIS为250/500m，SeaWiFS为1000m)，更加适合实时的水色参数监测。遗憾的是，用MODIS用于实时Chla监测还存在一些尚未需要解决的问题。首先，现有的MODIS的Chla算法(如OC2,OC3,OC4和OCI等等)不适用于高浑浊的富营养化水体；第二，MODIS用于水体研究的海洋水色波段在沿海水域和内陆水域都是饱和的，几乎没有数据可以使用，因此，只能利用不饱和陆地宽波段的算法来反演Chla；最后，目前没有较为可靠的大气校正办法帮助MODIS获取到内陆浑浊湖泊的Rrs数据(经过完全的准确的大气校正的数据)。因此，利用MODIS经过局部大气校正的不饱和光谱波段，发展一个新的方法，来实现高混浊富营养化水体中Chla的反演，成为下一步研究的方向。

[0005] 本发明的目的在于提供一种富营养化湖泊水体叶绿素a的MODIS卫星高精度监测方法，可精确获取浅水湖泊叶绿素a浓度时空分布，准确分析蓝藻水华发生、发展状况及趋势，科学评估湖泊污染治理与生态修复效果，为水利、环保等部门的水资源管理、水环境保护的科学决策提供科技支撑。

[0006] 本发明的上述目的通过独立权利要求的技术特征实现，从属权利要求以另选或有利的方式发展独立权利要求的技术特征。

[0007] 为达成上述目的，本发明所采用的技术方案如下：

[0008] 一种富营养化湖泊水体叶绿素a的MODIS卫星高精度监测方法，所述方法包括以下步骤：

[0009] 1)筛选对叶绿素a浓度变化敏感而又不受高悬浮物、薄云、太阳耀斑影响的叶绿素a评价指数BNDBI；

[0010] 其中，所述对叶绿素a浓度变化敏感而又不受高悬浮物、薄云、太阳耀斑影响的叶绿素a评价指数BNDBI是指基于叶绿素a、无机悬浮物和黄色物质光谱响应特征，以蓝光与近红外波段为基线，选择红、绿波段间的差值表达形式，可避免高悬浮物、薄云、太阳耀斑对叶绿素a浓度估算的不利影响，并以该藻类指数作为叶绿素a浓度遥感监测指数；

[0011] 2)基于生物光学模型模拟，明确BNDBI与叶绿素a浓度间的定量关系；

[0012] 在生物光学模型的基础上，结合湖泊的实测数据，进行不同情景下的数值模拟，确定BNDBI与叶绿素a浓度的定量关系，同时从理论上确定该指数对高浑浊水体不敏感性；

[0013] 3)获取地面监测遥感反射比Rrs与模拟的瑞利散射矫正后Rrc之间的定量关系；

[0014] 模拟湖泊地区在不同气溶胶类型及厚度、不同太阳高度角、卫星观测角以及方位角情况下，基于地面监测的遥感反射比Rrs的BNDBI与基于模拟的瑞利散射矫正的Rrc的BNDBI之间的定量关系；

[0015] 4)获取基于MODIS卫星影像的叶绿素a浓度反演算法；

[0016] 基于前述步骤和方法，将基于地面实测光谱数据的叶绿素a反演算法应用至经过瑞利散射矫正的卫星影像数据，基于前述的反演算法方法，在对多幅时间序列的卫星影像处理后获取富营养化湖泊叶绿素a浓度的年际、月际变化规律及其空间分布。

[0017] 作为本发明的进一步改进，所述步骤1)中，叶绿素a、无机悬浮物和黄色物质光谱响应特征来自于湖泊野外实测的光谱数据Rrs，采用的监测仪器为美国ASD公司的双通道地面光谱监测仪。

[0018] 作为本发明的进一步改进，所述步骤2)中，基于地面实测光谱数据的叶绿素a评价指数BNDBI表达形式为：

[0019] BNDBI＝(Rrs555'-Rrs645')/(Rrs555'-Rrs645')

[0020]

[0021]

[0022] 作为本发明的进一步改进，所述步骤2)中“进行不同情景的数值模拟”具体包括：

[0023] 首先在无机悬浮物浓度和黄色物质保持不变的情况下，获取BNDBI与叶绿素a浓度间的定量关系；

[0024] 其次，模拟叶绿素a和黄色物质浓度不变时，BNDBI对无机悬浮物浓度的响应；

[0025] 最后模拟叶绿素a和无机悬浮物浓度保持不变时，黄色物质浓度变化对BNDBI的影响。

[0026] 作为本发明的进一步改进，所述步骤3)中，气溶胶类型参照SeaDas的LUT的结果，气溶胶厚度参照湖泊地区常年监测结果范围，卫星观测角则依据太阳、卫星和湖泊的相对位置确定。

[0027] 作为本发明的进一步改进，所述步骤4)中，应用于MODIS影像的BNDBI指数表达形式为：

[0028]

[0029]

[0030]

[0031] 并且，MODIS影像的BNDBI指数建立在MODIS卫星影像的辐射定标、几何纠正和大气瑞利散射校正的基础上。

[0032] 由以上本发明的技术方案可知，本发明的富营养化湖泊水体叶绿素a的MODIS卫星高精度监测方法，基于湖泊叶绿素a遥感监测基本原理，确定对水体中无机悬浮物以及黄色物质不敏感的叶绿素a监测指数BNDBI，通过叶绿素a浓度与地面实测光谱的BNDBI及经过瑞利散射矫正后的BNDBI数据的定量关系这一核心，来实现对MODIS影像中叶绿素a的高精度估算，更加客观真实地反映湖泊叶绿素a的时空分布。叶绿素a的高精度监测可以有效的实现对湖泊藻华风险以及对水源区进行有效评估；湖泊叶绿素a浓度的长期高精度监测，有助于科学评估年际间叶绿素a的变化及其发展趋势，有效评估湖泊污染治理和生态修复的绩效，为水利、环保等部门的水资源管理、水环境保护的科学决策提供科技支撑。

[0033] 应当理解，前述构思以及在下面更加详细地描述的额外构思的所有组合只要在这样的构思不相互矛盾的情况下都可以被视为本公开的发明主题的一部分。另外，所要求保护的主题的所有组合都被视为本公开的发明主题的一部分。

[0034] 结合附图从下面的描述中可以更加全面地理解本发明教导的前述和其他方面、实施例和特征。本发明的其他附加方面例如示例性实施方式的特征和/或有益效果将在下面的描述中显见，或通过根据本发明教导的具体实施方式的实践中得知。

[0035] 附图不意在按比例绘制。在附图中，在各个图中示出的每个相同或近似相同的组成部分可以用相同的标号表示。为了清晰起见，在每个图中，并非每个组成部分均被标记。现在，将通过实施例并参考附图来描述本发明的各个方面的实施例，其中：

[0036] 图1是BNDBI指数监测叶绿素a的基本原理示意图。

[0037] 图2是理论模拟下BNDBI与叶绿素a浓度间定量关系。

[0038] 图3是不同气溶胶类型及其厚度，不同太阳高度角、卫星观测角以及方位角时，基于Rrs的BNDBI与基于Rrc的BNDBI间的定量关系。

[0039] 图4是巢湖叶绿素a的MODIS卫星高精度监测空间分布结果(2012年9月25日)。

[0040] 图5是BNDBI算法对于高悬浮物影响的应用示意图。

[0041] 图6是BNDBI算法对于薄云影响的应用示意图。

[0042] 图7是BNDBI算法对于太阳耀斑影响的应用示意图。

[0043] 前述图示1-7中，作为英文形式表达的各坐标、标识或其他表示，均为本领域所公知的，并不在本例中再做赘述。

[0044] 为了更了解本发明的技术内容，特举具体实施例并配合所附图式说明如下。

[0045] 在本公开中参照附图来描述本发明的各方面，附图中示出了许多说明的实施例。本公开的实施例不必定意在包括本发明的所有方面。应当理解，上面介绍的多种构思和实施例，以及下面更加详细地描述的那些构思和实施方式可以以很多方式中任意一种来实施，这是应为本发明所公开的构思和实施例并不限于任何实施方式。另外，本发明公开的一些方面可以单独使用，或者与本发明公开的其他方面的任何适当组合来使用。

[0046] 本发明给予MODIS卫星数据对富营养化湖泊叶绿素a浓度进行高精度监测，上述目的是这样实现的：

[0047] 筛选对叶绿素a浓度变化敏感而又不受高悬浮物、薄云、太阳耀斑影响的叶绿素a评价指数(BNDBI)；

[0048] 在生物光学模型模拟的基础上，从理论上论证BNDBI与叶绿素a浓度间的定量关系；

[0049] 结合2013-2014年巢湖野外监测的水体光谱信息及相应的水体叶绿素a浓度，获取基于地面实测光谱和BNDBI的叶绿素a反演算法；

[0050] 通过模拟不同气溶胶类型及厚度、不同太阳高度角、卫星观测角及方位角，获取基于Rrs的BNDBI与基于Rrc的BNDBI间的定量关系；进而将基于地面实测光谱数据的叶绿素a反演算法推广至经过瑞利散射矫正的卫星影像数据，进而估算湖泊全水域叶绿素a浓度及其空间分布。

[0051] 作为示例性的描述，下面结合附图所示，对前述方法的实施进行具体说明[0052] 步骤1、确定叶绿素a的监测评价指数BNDBI

[0053] 对叶绿素a浓度变化敏感而又不受高悬浮物、薄云、太阳耀斑影响的叶绿素a评价指数BNDBI是指基于叶绿素a以及悬浮物光谱响应特征，以蓝光波段和近红外波段为基线，选择红、绿波段标准差值形式，可避免高悬浮物、薄云、太阳耀斑对叶绿素a浓度估算的不利影响，并以该藻类指数作为叶绿素a浓度遥感监测指数。

[0054] 具体地，基于水体光学活性物质(叶绿素a、无机悬浮物、黄色物质)的基本监测原理，研究水体中三种光学活性物质的光谱特征，并结合现有叶绿素a在二类水体中的反演算法，比较分析各自优缺点，选择准确估算叶绿素a浓度的同时，不受水体中无机悬浮物以及黄色物质的影响的监测指数作为蓝藻水华MODIS卫星监测的基本指数，以克服水体中其他光学活性物质对叶绿素a浓度监测的不利情况。

[0055] 本实施例中，由于高叶绿素a的水体在绿光波段(570nm)有个反射峰，而因叶绿素a在665nm的强吸收而造成红光波段的反射谷(图1)，因此，可以从MODIS对应的叶绿素a特征波段来估算叶绿素a的含量。图1是MODIS波段设置下高叶绿素a、高浑浊以及一般水体的光谱以及三者的差别，可以看出如果以555nm、645nm波段为两端基点，高叶绿素a水体与高浑浊水体具有最大的差别，此外，基于蓝光波段(469nm)和近红外波段(859nm)为基线，可以滤除一定的薄云和太阳耀斑造成的影响。根据此特点，提出了BNDBI(Baseline Normalized difference bloom index)指数：

[0056] BNDBI＝(Rrs555'-Rrs645')/(Rrs555'-Rrs645')

[0057]

[0058]

[0059] 其中，Rrs(λ)是地面测量获得的λ波长处水体遥感反射比。

[0060] 步骤2、生物光学模型模拟证明BNDBI与叶绿素a浓度间的定量关系[0061] 在生物光学模型的基础上，结合巢湖的实测数据，开展不同情景下的数值模拟，从理论上证明BNDBI与叶绿素a浓度的定量关系，以及水体中其他光学活性物质对该算法影响。

[0062] 本实施例中，对于一般水体，水体的遥感反射比与水体的固有光学属性成正比关系，

[0063]

[0064] a(λ)＝aw(λ)+aph(λ)+ad(λ)+ag(λ)

[0065] bb(λ)＝bbw(λ)+bbp(λ) (2)

[0066] 其中aw和bbw对应着纯水的吸收系数和后向散射系数；而aph、ad和ag则是浮游植物色素、无机悬浮物和黄色物质的吸收系数，它们都与水体中相应物质的量存在着密切关系，bbp是水体中颗粒物的后向散射系数，在藻类含量不高的水体中，该系数与无机悬浮物有着密切关系。其中，

[0067]

[0068] 根据式(1)，BNDBI与叶绿素a浓度间存在以下关系，

[0069]

[0070] 根据式(4)，BNDBI和叶绿素a浓度间存在着单调关系，也就是BNDBI随叶绿素a浓度的增加而增加。因此，假定水体中无机悬浮物浓度为50mg/L，在忽略黄色物质影响的情况下，图2为基于生物光学模型模拟的BNDBI和叶绿素a浓度间定量关系。

[0071] 根据我们2013-2014年在巢湖的野外实测的光谱数据以及相应的叶绿素a浓度数据，我们构建了基于实测光谱数据Rrs的巢湖叶绿素a的反演算法。4 3 2

[0072] Chla＝982.3*BNDBI+71.86*BNDBI+562.4*BNDBI+79.05*BNDBI+6.6(5)[0073] 步骤3、获取基于地面实测光谱数据Rrs的BNDBI与基于模拟的瑞利散射校正后的Rrc的BNDBI之间的定量关系

[0074] 考察巢湖地区在不同气溶胶类型及厚度、不同太阳高度角、卫星观测角以及方位角对基于地面实测光谱数据Rrs的BNDBI与基于模拟的瑞利散射校正后Rrc的BNDBI之间的定量关系的影响，并通过模拟数据确定两者之间的定量模型。

[0075] 本实施例中，要将基于实测光谱数据得到叶绿素a的反演算法推广至卫星影像数据，大气矫正不可忽略。但是目前仍然缺乏针对高浑浊水体有效的精确大气矫正算法，本次采用MODIS影像的瑞利散射矫正，也就是通过这种矫正，大气层顶的光学信息去除了瑞利散射的影响，依然包含着气溶胶信息以及地面信息。基于瑞利散射校正后的数据，BNDBI表达为：

[0076]

[0077]

[0078]

[0079] 其中，Rrc(λ)是经过瑞利校正的λ波长处的反射率。Rrc是MODIS数据进行瑞利散射校正，然后基于Hu等(2004)的研究将其转换为瑞利散射校正后的反射率：

[0080]

[0081] 式中， 是校正臭氧和其他气体吸收效应后的传感器辐射率，F0是获取数据时的大气圈外太阳辐照度，θ0是太阳天顶角，Rr是采用6S(Vermote等，1997)预测的瑞丽反射率。

[0082] 基于辐射传输理论以及假定一个非耦合的海洋—大气系统，Rrc可以表达为：

[0083] Rrc＝Ra+t0tRtarget (8)

[0084] 式中，Ra是气溶胶反射率(包括来自于气溶胶分子的相互作用)，Rtarget是野外实测目标(藻类或者水体)的表面反射率，t0是从太阳到目标物的大气透射率，t是从目标物到卫星传感器的大气透射率。由于受风和水流的影响，浮游藻类通常呈现一种水面浮油的形态，因此t可以看作浮游藻类的光透射率。

[0085] 为了考察不同气溶胶类型及其厚度，以及卫星观测造成的影响，我们根据巢湖地区在不同气溶胶类型及厚度、不同太阳高度角、卫星观测角以及方位角对基于地面实测光谱数据Rrs的BNDBI与基于模拟的瑞利散射校正后Rrc的BNDBI之间的定量关系的影响(图3)，并通过模拟数据确定两者之间的定量模型，

[0086] BNDBI(Rrc)＝1.051BNDBI(Rrs)-0.007

[0087] (9)

[0088] 步骤4、获取基于MODIS卫星影像的叶绿素a浓度反演算法

[0089] 基于前述步骤和方法，将基于地面实测光谱数据的叶绿素a反演算法应用至经过瑞利散射矫正的卫星影像数据，基于前述的反演算法方法，在对多幅时间序列的卫星影像处理后获取富营养化湖泊叶绿素a浓度的年际、月际变化规律及其空间分布。

[0090] 基于式(5)和式(9)，即可获取巢湖叶绿素a的MODIS卫星高精度反演模型。

[0091] 依据基于MODIS影像的瑞利散射矫正，结合式(5)和式(9)即可实现全影像中叶绿素a浓度的高精度估算。具体流程主要如下：①对获取的MODIS影像进行了几何纠正和辐射定标计算。几何纠正采用Geographic Lat/Lon投影，结合1B数据中的经纬度信息进行校正,校正后的位置精度达到0.5个像素。在ERDAS中利用湖泊矢量边界，通过掩膜技术提取湖泊水域，除去岛屿植被的影响，利用最近邻法，将MODIS 500m影像数据重采样为250m；②MODIS影像中逐一像元计算其在band 1(645nm)、band 2(859nm)、band 3(469nm)和band 4(555nm)的Rrc值；③根据公式(6)逐一像元计算BNDBI值；④继而根据式(5)和式(9)，即可得到计算后的叶绿素a空间分布结果(见图4)。

[0092] 此外，我们还将该基线算法进行了高浑浊水体、薄云以及太阳耀斑条件下叶绿素a浓度反演对比，结果分别如图5-图7所示，相比较发现该模型依然具有较高的精度，具有一定的普适性。

[0093] 通过上述方法即可实现对某一MODIS影像中藻华混合像元的藻华面积估算，更加客观真实地反映湖泊藻华面积的高精度估算及其时空分布。蓝藻面积的高精度监测可以有效的实现对湖泊藻华风险以及对水源区进行有效评估；此外，将MODIS历史影像通过上述方法逐一计算后，即可实现湖泊蓝藻面积的长期高精度监测，有助于科学评估年际间藻华实际强度的变化及其发展趋势，有效评估湖泊污染治理和生态修复的绩效，为水利、环保等部门的水资源管理、水环境保护的科学决策提供科技支撑。

[0094] 虽然本发明已以较佳实施例揭露如上，然其并非用以限定本发明。本发明所属技术领域中具有通常知识者，在不脱离本发明的精神和范围内，当可作各种的更动与润饰。因此，本发明的保护范围当视权利要求书所界定者为准。