华南沿海赤潮灾害应急监测与预警方法

技术领域

本发明涉及一种赤潮水体分析方法，具体涉及一种赤潮灾害应急监测与预警的方法。

背景技术

基于卫星遥感影像MODIS反演的水质参数分布数据产品可以用于分析赤潮爆发几率判别，水体富营养化情况以及初级生产力水平等。赤潮爆发的准确判别对研究海水水质、渔业资源、旅游业等发展有着重要的意义。

近年来，随着人类无节制的经济活动对自然界产生的耦合效应，我国沿岸赤潮呈频发势态，给当地渔业造成了重大损失。2010年中国沿海共发生赤潮69次，累计面积10892平方公里。南海14次，累计面积223平方公里。2011年我国沿海共发生赤潮55次，累计面积6076平方千米。南海11次，累计面积190平方千米。2012年，我国沿海共发现赤潮73次，12次造成灾害，直接经济损失20.15亿元，其中南海区10次。2013年，我国沿海共发现赤潮46次，其中有毒赤潮7次。2013年最大面积大于等于100平方千米。2014年南海区共监测到赤潮16次，赤潮累计面积约684平方千米，赤潮生物种类为夜光藻、球形棕囊藻、红色赤潮藻、多纹膝沟藻、中肋骨条藻、赤潮异弯藻和条纹环沟藻。2014年两次较大面积的赤潮均发生在广东湛江港附近海域，赤潮种类为中肋骨条藻和夜光藻。从2001年以来南海区赤潮调查报告来看，2014年赤潮发生的面积较大，次数较多。赤潮直接威胁着人类的生存环境，已经成为我国目前主要的海洋灾害之一。因此，南海海域的赤潮监测研究工作已经迫在眉睫。赤潮监测预警体系亟需完善，赤潮防治、应急处理能力有待提高。

研究表明，利用遥感影像进行赤潮水体判别，可以即时地对赤潮灾害进行预警以及采取及时有效地灾害决策。因此，基于遥感反演数据产品对赤潮水体进行有效判别，是进一步提高海域资源管理的主要发展方向。

发明内容

本发明的目的是提供一种赤潮水体判别方法，用于华南沿海赤潮灾害应急监测与预警，该方法基于多因素综合分析，提高了赤潮判别的时效性，能够进一步提高赤潮应急灾害管理效率。

为实现上述目的，本发明提供以下技术方案：

一种华南沿海赤潮灾害应急监测与预警方法，其步骤包括：

A：按照所需时间、区域测量得到相应水体光谱数据后计算得到相应水体的反射率数据；

B：利用步骤A所得反射率数据计算荧光高度(FLH)，其中FLH的计算采用667nm、678nm、748nm三个波段进行比值计算；计算FLH的具体方法如下式(1)所示：

$>FLH={\mathrm{nL}}_w\left(678\right)-{\mathrm{nL}}_w\left(667\right)-[{\mathrm{nL}}_w\left(748\right)-{\mathrm{nL}}_w\left(667\right)]\frac{678-667}{748-667}---\left(1\right)$>

C：利用步骤A所得反射率数据采取6波段算法提取海水颗粒物后向散射系数，具体步骤如下：

步骤A得到的已经过大气校正的反射率数据r_rs表达为：

$>r_{rs}=g_1\frac{b_b}{a+b_b}+g_2{\left(\frac{b_b}{a+b_b}\right)}^2---\left(2\right)$>

其中g₁、g₂分别为0.0949和0.0797；a是海水总吸收系数，b_b是海水总后向散射系数；

其中a包含纯海水吸收、CDM吸收和浮游植物吸收的总和；其中浮游植物吸收可以表达为：

a_ph(λ)＝a^*[chl](3)

其中，a^*是叶绿素比吸收系数(m²/mg)；[chl]是叶绿素浓度(mg/m³)；(这里，a^*在412,443,488,531,555,667nm)，参数选取参见表1

表1：模型参数

CDOM吸收与波长存在一个指数关系，表达为：

a_CDM(λ)＝a_CDM(443)exp[-S(λ-443)](4)

S是a_CDM光谱斜率，可以通过反射率比值计算得出：

$>s=0.015+\frac{0.002}{0.6+\frac{r_{rs}\left(443\right)}{r_{rs}\left(555\right)}}---\left(5\right)$>

所述的海水总后向散射系数由纯水后向散射与颗粒物后向散射系数计算得出，其中颗粒物后向散射系数通过下式计算得出：

b_bp(λ)＝b_bp(λ₀)(λ/λ₀)^-Y(6)

Y是后向散射光谱指数，在这里取常数1(λ₀＝443nm)

固定参数S与Y，将式(3)-(6)代入式(2)，解非线性方程，估算出叶绿素浓度[chl]、b_bp(443)和a_CDM(443)；

D：将步骤C中提取得到的b_bp(443)作为输入参数，同时固定S，将所述的式(3)-(6)再次代入式(2)，再次解非线性方程，估算出叶绿素浓度[chl]’、a_CDM(443)’和Y’；

E：根据步骤B所得的荧光高度值和步骤D所得的叶绿素浓度值判读赤潮爆发几率，具体包括：

E1.根据步骤B得到的结果做第一级判读：FLH大于等于0.05mWcm-²μm^-1sr^-1的区域判别为高度疑似，反之判别为低度疑似；

E2.根据步骤D得到的叶绿素浓度[chl]’值对步骤F1判别为高度疑似的区域进一步做第二级判读：[chl]’值高于5mg/L的区域判别为赤潮爆发几率高，反之判别为赤潮爆发几率低。

本发明所述的方法中，步骤A所述的测量得到相应水体光谱数据的测量方法可以是实地测量方法，也可以是遥感测量方法；本发明优选遥感测量方法。

本发明一种优选的方案中，步骤A的具体方法是按照所需时间、区域获取相应光学遥感数据，对其进行辐射定标和几何校正，以及大气校正处理；然后提取已经过大气校正的反射率数据。

进一步优选的方案中，所述的光学遥感数据是MODISL1遥感影像数据，经过辐射定标、几何校正、以及大气校正处理后得到经过大气校正的MODISL2反射率数据。

更进一步优选的方案中，步骤B计算荧光高度前先利用所述的经过大气校正的MODISL2反射率数据采用ENVI软件进行强化真彩影像合成(ERGB)，其中红色波段选取555nm(R)波段，绿色波段选取488nm(G)波段，蓝色波段选取443nm(B)波段；且步骤E1在做所述第一级判读前先根据强化真彩影像合成结果做第零级判读：合成的影像中呈现黑褐色的区域判别为赤潮疑似，未呈现黑褐色的区域判别为非赤潮疑似。

本发明优选的另一种方案中，进一步将步骤C得到的443nm颗粒物后向散射系数b_bp(443)和步骤B得到的荧光高度(FLH)进行比较，计算b_bp/FLH的值，记作“赤潮爆发几率判别指数”；并根据b_bp/FLH的值对步骤F3判别为赤潮爆发几率高的情况做进一步判读：b_bp/FLH小于0.2mW^-1cm²μmm^-1sr的判别为赤潮爆发几率很高，反之判别为赤潮爆发几率略高。

本发明的方法运用多因素的综合分析手段对于海水赤潮爆发几率进行判读，具体通过对荧光高度(FLH)、叶绿素浓度反演(chl)、后向散射系数反演(b_bp)、赤潮水体判别指数b_bp/FLH这些不同角度、不同表达方式的因素进行综合分析，最终实现了对水体赤潮爆发几率的综合地、准确地判别。该方法中的指标选择、数据处理和判读标准都是建立在对大量的实践经验和自然规律的总结的基础上，保证了判断水体赤潮爆发几率的高效性和准确性，特别是在运用遥感技术获取遥感影像数据的方案中还有效地利用了遥感数据，进一步利用强化真彩数据合成(ERGB)技术为提高水色遥感反演算法提供了技术支持。

附图说明

图1为本发明实施例1中步骤C得到的荧光高度计算结果图。

图2为本发明实施例1中步骤E得到的叶绿素浓度[chl]’反演结果图。

图3为本发明实施例1中步骤F4得到的b_bp/FLH的遥感影像数据分析结果图。

图4体现了本发明实施例2的6波段法对叶绿素浓度的反演精度。

图5是本发明实施例3步骤A得到的观测样本水体反射率光谱曲线图。

图6是本发明实施例3步骤B得到的观测样本水体FLH分布图。

图7是本发明实施例3步骤D得到的观测样本水体叶绿素浓度分布图。

图8是本发明实施例3步骤E3得到的观测样本水体b_bp/FLH分布图。

具体实施方式

以下结合实施例及附图对本发明进行详细的描述。

实施例1

本实施以2014年11月24日发生在粤西近海发生的赤潮事件为例，以MODISL1b数据作为遥感影像，详细说明利用所述方法进行华南沿海赤潮水体判别的流程。所述方法包括以下步骤：

A：接收2011年11月24日的MODISL1数据，利用定标参数对其进行辐射定标和几何校正，以及大气校正处理；利用Matlab软件通过hdf头文件参数对MODISL1b数据进行转换处理，估算大气分子瑞利散射。

F₀(λ)是大气层外对应不同波长的太阳辐照度，其大小随着日地距离而变化，μ和μ₀分别是太阳天顶角以及卫星天顶角的余弦，α_↓和α_↑，分别表示入射光和反射光的散射相位角，P_r(α_↓)是入射光的散射相函数，P_r(α_↑)反射光的散射相函数，T_OZ(λ)为经过臭氧校正的大气透过率。τ_r(λ)是瑞利光学厚度。

以氧气与水汽吸收波段组合进一步估算大气气溶胶散射。

ε(λ_n，λ_NIR)＝exp[C(λ_NIR-λ_n)](8)

其中λ_n为待估波长，λ_NIR为参考波长。

利用上述所得结果，提取已经过大气校正的MODISL2反射率数据；

根据文件参数，进行反射率计算。

$>r_{rs}=\frac{L_w}{E_d}---\left(9\right)$>

L_W是经过大气校正后的离水辐亮度，E_d是入射辐照度。

B：利用步骤A所得MODISL2反射率数据，进一步采用ENVI软件进行强化真彩影像合成(ERGB)，其中红色波段选取555nm(R)波段，绿色波段选取488nm(G)波段，蓝色波段选取443nm(B)波段，得到ERGB结果图。

C：利用步骤A所得MODISL2反射率数据，计算荧光高度(FLH)，其中FLH的计算采用667nm，678nm，748nm三个波段进行比值计算。计算FLH的具体方法如下：

$>FLH={\mathrm{nL}}_w\left(678\right)-{\mathrm{nL}}_w\left(667\right)-[{\mathrm{nL}}_w\left(748\right)-{\mathrm{nL}}_w\left(667\right)]\frac{678-667}{748-667}---\left(10\right)$>

计算结果显示FLH分布在0.005到0.2mWcm-²μm^-1sr^-1之间，相应的遥感影像数据分析结果见图1所示。

D：利用步骤A所得MODISL2反射率数据，采取6波段算法计算海水颗粒物后向散射系数，具体步骤如下：

412-665nm遥感反射率表达r_rs(已经过校正的MODISL2数据)为：

$>r_{rs}=g_1\frac{b_b}{a+b_b}+g_2{\left(\frac{b_b}{a+b_b}\right)}^2---\left(11\right)$>

其中g₁和g₂分别为0.0949和0.0797。a是海水总吸收系数，b_b是海水总后向散射系数。

其中a包含纯海水吸收、CDM吸收和浮游植物吸收的总和。

其中叶绿素浓度与浮游植物的吸收关系为：

a_ph(λ)＝a^*[chl](12)

其中，a^*是叶绿素比吸收系数(m²/mg)；[chl]是叶绿素浓度(mg/m³)；(这里，a^*在412,443,488,531,555,667nm)，参数选取参见表1

表1：模型参数

CDOM吸收与波长存在一个指数关系，表达为：

a_CDM(λ)＝a_CDM(443)exp[-S(λ-443)](13)

S是a_CDM光谱斜率，可以通过反射率比值计算得出：

$>s=0.015+\frac{0.002}{0.6+\frac{r_{rs}\left(443\right)}{r_{rs}\left(555\right)}}---\left(14\right)$>

所述的海水总后向散射系数由纯水后向散射与颗粒物后向散射系数计算得出，其中颗粒物后向散射系数通过下式计算得出：

b_bp(λ)＝b_bp(λ₀)(λ/λ₀)^-Y(15)

Y是后向散射光谱指数，在这里取常数1(λ₀＝443nm)。

固定参数S与Y，将式(12)-(15)代入式(11)，解非线性方程，估算出叶绿素浓度[chl]、b_bp(443)和a_CDM(443)。

E：将步骤D中提取得到的b_bp(443)作为输入参数，同时固定S，将所述的式(12)-(15)再次代入式(11)，再次解非线性方程，估算出叶绿素浓度[chl]’、a_CDM(443)’和Y’；计算结果显示[chl]’分布在0到7mg/L之间；相应的叶绿素浓度[chl]’反演结果见图2所示。

F：根据步骤B所得的强化真彩影像合成结果、步骤C所得的荧光高度值、步骤D得到的443nm颗粒物后向散射系数b_bp(443)和步骤E所得的叶绿素浓度值判读赤潮爆发几率，具体包括：

F1.根据结果初步判别为疑似或非疑似：由于赤潮水体在ERGB影像中主要呈现黑褐色，而步骤B的得到的ERGB影像中部分区域出现明显的黑褐色，因此出现黑褐色的区域判别为疑似，而未出现黑褐色的区域判别为非疑似；

F2.根据步骤C得到的结果对步骤F1中判读为疑似的情况做进一步判读：判别标准是FLH值大于等于0.05mWcm-²μm^-1sr^-1判别为高度疑似，反之判别为低度疑似；步骤C计算结果如图1所示，FLH分布在0.005到0.2mWcm-²μm^-1sr^-1之间，因此，在被判别为疑似的区域内，FLH值大于等于0.005且小于0.05mWcm-²μm^-1sr^-1的区域(图1中方框标出的区域)进一步判别为低度疑似，而FLH值在0.05～0.2mWcm-²μm^-1sr^-1的区域进一步判别为高度疑似；

F3.根据步骤E得到的叶绿素浓度[chl]’值对步骤F2判别为高度疑似的情况做进一步判读：判别标准是[chl]’值高于5mg/L的判别为赤潮爆发几率高，反之判别为赤潮爆发几率低；步骤E计算结果如图2所示[chl]’分布在0到7mg/L之间，因此，在被判别为高度疑似的区域内，[chl]’分布在0到5mg/L之间的区域(图2中方框标出的区域)更进一步判别为赤潮爆发几率低，而[chl]’高于5mg/L直至7mg/L的区域更进一步判别为赤潮爆发几率高；

F4.进一步将步骤D得到的443nm颗粒物后向散射系数b_bp(443)和步骤C得到的荧光高度(FLH)进行比较，计算b_bp/FLH的值，记作“赤潮爆发几率判别指数”，相应的遥感影像数据分析结果见图3；并根据图3中b_bp/FLH的值对步骤F3判别为赤潮爆发几率高的情况做进一步判读：判别标准是b_bp/FLH小于0.2mW^-1cm²μmm^-1sr的判别为赤潮爆发几率很高，反之判别为赤潮爆发几率略高；因此，被判别为赤潮爆发几率高的区域内，b_bp/FLH值小于0.2mW^-1cm²μmm^-1sr的广东阳江市南部海湾地区(图3中方框标出的区域)被最终判读为赤潮爆发几率很高。另外基于该海域优势藻种为赤潮藻种夜光藻，并且该时间段为该海域赤潮多发时间段，因此判断赤潮发生可能性非常大。

2014年11月25日，南海航空支队执法组3808飞行快报报告在执行102号线巡航任务途中，在阳江闸坡以及茂名博贺港西南侧放鸡岛附近海域(21°23′N，111°16′E)以及沿线至下川岛以南海域发现疑似赤潮，呈砖红色、品红色，条带状、块状分布，累计分布面积约300平方公里。在南海航空支队执法组3808飞行快报报告中，实测的疑似赤潮区域的FLH分布在0.05到0.1mWcm-²μm^-1sr^-1之间；实测的疑似赤潮区域的[chl]’分布在5到7mg/L之间。上述预警判别方法中的预测结果与实测结果高度吻合，由此表明本发明所述的赤潮预警判别方法对赤潮爆发几率的判别能力较好，并可以为应对赤潮事件提供决策支持。

实施例2

本实施例以光学辐射传输模拟软件-Hydrolight为基础,应用了二类水体模型生成了350组模拟反射率数据，应用反射率数据提取了叶绿素浓度数据。并进行了算法精度验证。所述方法包括以下步骤：

A：应用Hydrolight辐射传输模拟软件模拟不同叶绿素浓度的海水反射率数据。

B：利用步骤A所得反射率数据，进一步采用采取6波段算法计算海水颗粒物后向散射系数，具体步骤如下：

412-665nm遥感反射率表达r_rs为：

$>r_{rs}=g_1\frac{b_b}{a+b_b}+g_2{\left(\frac{b_b}{a+b_b}\right)}^2---\left(11\right)$>

其中g₁和g₂分别为0.0949和0.0797。a是海水总吸收系数，b_b是海水总后向散射系数。

其中a包含纯海水吸收、CDM吸收和浮游植物吸收的总和。

其中叶绿素浓度与浮游植物的吸收关系为：

a_ph(λ)＝a^*[chl](12)

其中，a^*是叶绿素比吸收系数(m²/mg)；[chl]是叶绿素浓度(mg/m³)；(这里，a^*在412,443,488,531,555,667nm)，参数选取参见表1

表1：模型参数

CDOM吸收与波长存在一个指数关系，表达为：

a_CDM(λ)＝a_CDM(443)exp[-S(λ-443)](13)

S是a_CDM光谱斜率，可以通过反射率比值计算得出：

$>s=0.015+\frac{0.002}{0.6+\frac{r_{rs}\left(443\right)}{r_{rs}\left(555\right)}}---\left(14\right)$>

所述的海水总后向散射系数由纯水后向散射与颗粒物后向散射系数计算得出，其中颗粒物后向散射系数通过下式计算得出：

b_bp(λ)＝b_bp(λ₀)(λ/λ₀)^-Y(15)

Y是后向散射光谱指数，在这里取常数1(λ₀＝443nm)。

固定参数S与Y，将式(12)-(15)代入式(11)，解非线性方程，估算出叶绿素浓度[chl]、b_bp(443)和a_CDM(443)；

C：将步骤B中提取得到的b_bp(443)作为输入参数，同时固定S，将所述的式(12)-(15)再次代入式(11)，再次解非线性方程，估算出叶绿素浓度[chl]’、a_CDM(443)’和Y’；该模拟数据计算结果显示[chl]’分布在0到50mg/L之间；而实际输入Hydrolight数据也在此范围之内(如图4所示)，说明所述6波段法反演叶绿素浓度的精度较高。

本实施例利用Hydrolight模拟数据，用matlab实现优化算法，反演水质分布产品。将6波段叶绿素提取方法应用于该模拟数据。结果证明该算法可以有效地反演叶绿素浓度，并为赤潮遥感监测判断做出可靠的技术支持。

实施例3

本实施以粤西近海发生的赤潮事件为例，以现场观测数据作为输入参数，详细说明利用所述方法进行华南沿海赤潮水体判别的流程。所述方法包括以下步骤：

A：在粤西近海某观测现场随机选择10个观测点，现场采集光谱数据，进行反射率计算。

$>r_{rs}=\frac{L_w}{E_d}---\left(9\right)$>

L_W是离水辐亮度，E_d是入射辐照度。

反射率计算结果如图5所示。

B：利用步骤A所得反射率数据计算荧光高度(FLH)，其中FLH的计算采用667nm、678nm、748nm三个波段进行比值计算；计算FLH的具体方法如下式(10)所示：

$>FLH={\mathrm{nL}}_w\left(678\right)-{\mathrm{nL}}_w\left(667\right)-[{\mathrm{nL}}_w\left(748\right)-{\mathrm{nL}}_w\left(667\right)]\frac{678-667}{748-667}---\left(10\right)$>

计算得到的结果见图6所示，图中小圆圈代表观测点样本，横坐标对应样本的编号，纵坐标对应样本测量计算得到的FLH值。

C：利用步骤A所得反射率数据，进一步采取6波段算法计算海水颗粒物后向散射系数，具体步骤如下：

412-665nm遥感反射率表达r_rs为：

$>r_{rs}=g_1\frac{b_b}{a+b_b}+g_2{\left(\frac{b_b}{a+b_b}\right)}^2---\left(11\right)$>

其中g₁和g₂分别为0.0949和0.0797。a是海水总吸收系数，b_b是海水总后向散射系数。

其中a包含纯海水吸收、CDM吸收和浮游植物吸收的总和。

其中叶绿素浓度与浮游植物的吸收关系为：

a_ph(λ)＝a^*[chl](12)

其中，a^*是叶绿素比吸收系数(m²/mg)；[chl]是叶绿素浓度(mg/m³)；(这里，a^*在412,443,488,531,555,667nm)，参数选取参见表1

表1：模型参数

CDOM吸收与波长存在一个指数关系，表达为：

a_CDM(λ)＝a_CDM(443)exp[-S(λ-443)](13)

S是a_CDM光谱斜率，可以通过反射率比值计算得出：

$>s=0.015+\frac{0.002}{0.6+\frac{r_{rs}\left(443\right)}{r_{rs}\left(555\right)}}---\left(14\right)$>

所述的海水总后向散射系数由纯水后向散射与颗粒物后向散射系数计算得出，其中颗粒物后向散射系数通过下式计算得出：

b_bp(λ)＝b_bp(λ₀)(λ/λ₀)^-Y(15)

Y是后向散射光谱指数，在这里取常数1(λ₀＝443nm)。

固定参数S与Y，将式(12)-(15)代入式(11)，解非线性方程，估算出叶绿素浓度[chl]、b_bp(443)和a_CDM(443)；

D：将步骤C中提取得到的b_bp(443)作为输入参数，同时固定S，将所述的式(12)-(15)再次代入式(11)，再次解非线性方程，估算出叶绿素浓度[chl]’、a_CDM(443)’和Y’；该数据计算结果显示[chl]’分布在0到9mg/L之间；计算得到的[chl]’结果见图7所示，图7中小圆圈代表观测点样本，横坐标对应样本的编号，纵坐标对应样本测量计算得到的[chl]’值。

E：根据步骤B所得的荧光高度值、步骤C所得的b_bp(443)和步骤D所得的叶绿素浓度值判读赤潮爆发几率，具体包括：

E1.根据步骤B得到的结果赤潮水体做初步判读：判别标准是FLH大于等于0.05mWcm-²μm^-1sr^-1判别为高度疑似，反之判别为低度疑似；因此，采集数据对应的区域内，如图6所示，FLH值分布在大于等于0.05mWcm-²μm^-1sr^-1范围内的水体被划分为“赤潮水体”，其对应的观测地区被判别为赤潮高度疑似，其余水样被划分为“非赤潮水体”，其对应的观测区域为赤潮低度疑似；

E2.根据步骤D得到的叶绿素浓度[chl]’值对步骤E1判别为高度疑似的情况做进一步判读：判别标准是[chl]’值高于5mg/L的判别为赤潮爆发几率高，反之判别为赤潮爆发几率低；步骤D计算结果显示所述地区[chl]’分布在0到9mg/L之间，因此，被判别为高度疑似的区域内，如图7所示，[chl]’高于5mg/L直至9mg/L的水样被划分为“赤潮水体”，其对应的观测区域被进一步判别为赤潮爆发几率高，而其余水样被划分为“非赤潮水体”，其对应的观测区域被进一步判别为赤潮爆发几率低。

E3.进一步将步骤C得到的b_bp(443)和步骤B得到的FLH进行比较，计算b_bp/FLH的值，记作“赤潮爆发几率判别指数”；并根据b_bp/FLH的值对步骤E2判别为赤潮爆发几率高的情况做进一步判读：如图8所示，b_bp/FLH小于0.2mW^-1cm²μmm^-1sr的水样被划分为“赤潮水体”，其对应的观测地区判别为赤潮爆发几率很高，反之判别为赤潮爆发几率略高。

本实施例利用现场赤潮监测数据，用matlab实现优化算法，反演水质分布产品。将6波段叶绿素提取方法应用于该模拟数据。结果证明该赤潮判别方法可以有效地识别赤潮水体，并为赤潮遥感监测、预判做出可靠的技术支持。