一种基于高光谱遥感的水体富营养化状态测算方法

1.本发明属于水体测算方法技术领域，具体涉及一种基于高光谱遥感的水体富营养化状态测算方法。


背景技术：

2.水体富营养化是指水体中n、p等营养盐过多引起的水质污染现象，随着营养盐的增加，将会导致水体中的悬浮植物爆发式生长，从而引起水体中溶解氧的含量下降，进而导致水体恶化，鱼类和其他水生物死亡，因此对水体富营养化的监测和测算内陆水体的营养状态等级具有极其重要的生态意义。
3.传统的水体营养状态监测方法通常为在目标水体实地采集水体样本后再进行水质的参数测量，最后再根据计算出的营养状态指数来对目标水体的营养状态进行分类，但这样的方式会耗费大量的人力、物力，不能够大范围的、动态的对目标水体进行监测，因此现在逐渐开始采用遥感技术并配合水体颜色的营养状态分类方法，来对水体的营养状态进行分类，但该方法仅仅只能够将水体的营养状态粗略的分为贫、中、富三个营养状态，但对富营养化水体的治理是一个长期的过程，在面对不同程度的富营养化水体时，往往需要优先考虑治理富营养化较为严重的区域，在这样的情况下，富营养化的精细化分类就显得尤为重要。


技术实现要素：

4.本发明的目的是：旨在提供一种基于高光谱遥感的水体富营养化状态测算方法，通过高光谱成像仪构建出初步辨别水体营养化状态的检测参数和进一步辨别水体富营养化状态的检测参数，从而进一步达到精细化评价内陆水体富营养化等级的方法。
5.为实现上述技术目的，本发明采用的技术方案如下：
6.一种基于高光谱遥感的水体富营养化状态测算方法，包括设于目标水体的地基高光谱成像仪，该方法包含以下步骤：
7.s1：首先利用遥感图像处理软件，对地基高光谱成像仪采集得到的影像数据进行处理，分别提取出计算遥感反射率的三个参量；
8.s2：利用步骤s1中已提取出的三个参量，构建敏感波长的反射率曲线r
rs
；
9.s3：利用步骤s2中构建出的敏感波长反射率曲线，提取特定波长的遥感反射率，反演出cie颜色空间中的三刺激值，并计算目标水体颜色在cie色度图上对应的坐标；其次分别反演出叶绿素的浓度和目标水体的浊度；
10.s4：构建初步辨别目标水体营养化等级的检测参数fuci，利用步骤s3中的坐标点计算出对应的色度角，并运用查表法找出色度角对应的检测参数值；
11.s5：利用步骤s4中得到的检测参数对目标水体的营养化等级进行辨别，若目标水体为重度富营养化状态，则转至步骤s6，否则则将其辨别为检测参数所对应的营养化等级；
12.s6：构建用单个水质参数表征水体富营养化程度的检测参数itlc，利用步骤s3得
到的水质参数分别计算出基于叶绿素浓度或者目标水体浊度的itlc；
13.s7：将目标水体的浊度与透明度做相关性分析，判别该目标水体的浊度和透明度是否高度相关，若是存在高度相关性，则转至步骤s8，否则令该目标水体中浊度的营养状态指数的相关权重为0，并根据步骤s3中反演得到的叶绿素浓度和浊度，重新计算目标水体中浊度和叶绿素浓度的营养状态指数的相关权重；
14.s8：利用基于不同水质参数计算出的itlc以及其营养状态指数的相关权重，构建基于叶绿素浓度和目标水体浊度的内陆水体富营养化综合评价指数itlc(∑)；
15.s9：利用内陆水体富营养化综合评价参数itlc(∑)对目标水体进行判定，最终辨别水体的富营养化等级。
16.其中，所述步骤s1中提取的三个参量分别为标准板测量值l
p
、天空漫散射光l
sky
、离水辐亮度l
w
，具体提取方法如下：
17.l
p
：使用遥感图像处理软件打开对地基高光谱成像仪采集得到的影像数据源文件，选取任意一端波段形成监测区域影像，将该影像中的整个白板选取为感兴趣区域，利用roi工具读取出该感兴趣区域的最大值，形成各个波段的标准板测量最大值集合，并将其以text文档形式导出保存。
18.l
sky
：使用遥感图像处理软件打开对地基高光谱成像仪采集得到的影像数据源文件，选取任意一波段形成监测区域影像，选取该影像中最远处天空做为感兴趣区域，并利用roi工具读取出该感兴趣区域的平均值，形成包含各个波段的光谱曲线，并将其同样以text文档形式导出保存。
19.l
w
：使用遥感图像处理软件打开对地基高光谱成像仪采集得到的影像数据源文件，选取任意一波段形成监测区域影像，从目标水域中选择无白帽、无波纹、无倒影、无船只的水体做为感兴趣区域，同样利用roi工具读取出该感兴趣区域的平均值，形成包含各个波段的光谱曲线，并将其同样以text文档形式导出保存。
20.其中，所述步骤s2的具体步骤为，将步骤s1得到的三种参数数据导入excel表格中，并进行通过如下两步计算形成敏感波长的反射率曲线：
21.(1)
22.(2)
23.ρ
p
为标准板的反射率，e
d
(0
+
)表示水表面入射总辐照度，ρ表示气—水界面对天空光的反射率。
24.其中，所述步骤s3提取的特定波长遥感反射率r
rs
(r)、r
rs
(g)、r
rs
(b)、r
rs
(λ
a
)、r
rs
(λ
b
)、r
rs
(λ
c
)和r
rs
(λ
d
)，其对应波长为λ
r
、λ
g
、λ
b
、λ
a
、λ
b
、λ
c
和λ
d
，然后利用地基高光谱红绿蓝三个波段反射率值计算cie颜色空间中的三刺激值x、y、z，并进一步计算cie色度图中颜色对应的坐标(m，n)，三刺激值和坐标的计算公式如下：
25.cie颜色空间三刺激值：
26.cie色度图中颜色所对应的坐标：
27.分别利用叶绿素浓度和水体浊度的反演公式，反演得到叶绿素的浓度以及目标水体浊度。叶绿素浓度和水体浊度的反演公式如下：
28.叶绿素浓度反演公式：
29.浊度反演公式：
30.其中，所述步骤s4中所构建的检测参数为水体色度fuci，cie色度图中坐标对应的色度角为θ，其计算公式如下：
31.色度角θ：
32.将计算得到的色度角θ与fuci查找表中的数值相对应，找到查找表中与色度角θ最接近的色度值，其值对应的fuci即为该水体色度。
33.其中，所述步骤s5中对目标水体的营养化等级判别标准如下：
34.贫营养状态：fuci＜7
35.中营养状态：7≤fuci＜10
36.富营养状态：10≤fuci且r
rs
(r)≥h，否则仍为中营养状态
37.将步骤s4中得到的fuci与标准进行对比，若目标水体属于富营养状态的水体，则转至s6，否则该目标水体则为fuci所对应的营养状态。
38.其中，所述步骤s6中所构建的检测参数为营养等级分类指数itlc，其计算公式如下：
39.营养等级分类指数itlc：itlc(x)＝a+bln
x
40.其中x为叶绿素浓度或水体浊度中的一种；
41.所述步骤s6具体为：利用步骤s3所得到的叶绿素浓度和目标水体浊度分别代入营养等级分类指数公式，经计算得到基于叶绿素浓度的营养等级分类指数和基于目标水体浊度的营养等级分类指数。
42.其中，所述步骤s7将目标水体的实测浊度和透明度进行相关性分析，辨别该目标水体浊度和透明度的相关系数是否大于0.8，如果r2>0.8,则说明在该水体中，浊度和透明度有很强的相关性，可以用浊度代替透明度构建富营养化综合评价指数，转至步骤s8，否则则表示在该水体中，不能用浊度代替透明度来构建富营养化综合评价指数，于是令其浊度的营养状态指数的相关权重为0；
43.对浊度或叶绿素的营养状态指数的相关权重重新计算，其计算公式如下：
44.相关权重计算公式：
45.r
ij
表示第j种参数与基准参数i(叶绿素浓度)的相关系数，在这里只有两种参数，第一种为叶绿素浓度，第二种为水体浊度。
46.其中，所述步骤s8中，利用基于不同水质参数计算出的itlc以及其营养状态指数的相关权重，构建基于叶绿素浓度和目标水体浊度的内陆水体富营养化综合评价指数itlc(∑)，其计算公式如下：
47.内陆水体富营养化综合评价指数itlc(∑)：
[0048][0049]
步骤s9具体为：将在步骤s6中所得到的叶绿素浓度和目标水体浊度带上其相关权重求和，最终得到基于叶绿素浓度和目标水体浊度的内陆水体富营养化综合评价指数itlc(∑)。
[0050]
其中，所述步骤s9对目标水体的富营养化程度进行了判别，其判别标准如下：
[0051]
轻度富营养状态：0＜itlc(∑)≤60
[0052]
中度富营养状态：60＜itlc(∑)≤70
[0053]
重度富营养状态：70＜itlc(∑)
[0054]
将步骤s8中计算得到的itlc(∑)进行判别，最终得到目标水体的富营养化状态。
[0055]
通过地基高光谱成像仪采集水体数据并处理，进而得到目标水体的敏感波长反射率曲线，进一步演算得到cie颜色空间中的三刺激值和坐标，配合fuci查找表和反演公式初步得到水体的色度和富营养化等级，而为了水体治理的时效性和准确性，通过对水体的浊度与透明度关系的判断，进一步区分富营养化程度下水体的等级，从而更加利于对水体的分类，及时地治理富营养化更加严重的水体。
附图说明
[0056]
本发明可以通过附图给出的非限定性实施例进一步说明。
[0057]
图1为本发明方法的初步判断流程图示意图。
[0058]
图2为本发明方法的进一步判断流程图示意图。
[0059]
图3为本发明方法的fuci查找表。
具体实施方式
[0060]
为了使本领域的技术人员可以更好地理解本发明，下面结合附图和实施例对本发明技术方案进一步说明。
[0061]
如图1
‑
3所示，一种基于高光谱遥感的水体富营养化状态测算方法，包括设于目标水体的地基高光谱成像仪，该方法包含以下步骤：
[0062]
s1：首先利用遥感图像处理软件，对地基高光谱成像仪采集得到的影像数据进行处理，分别提取出计算遥感反射率的三个参量；
[0063]
s2：利用步骤s1中已提取出的三个参量，构建敏感波长的反射率曲线r
rs
；
[0064]
s3：利用步骤s2中构建出的敏感波长反射率曲线，提取特定波长的遥感反射率，反演出cie颜色空间中的三刺激值，并计算目标水体颜色在cie色度图上对应的坐标；其次分别反演出叶绿素的浓度和目标水体的浊度；
[0065]
s4：构建初步辨别目标水体营养化等级的检测参数fuci，利用步骤s3中的坐标点计算出对应的色度角，并运用查表法找出色度角对应的检测参数值；
[0066]
s5：利用步骤s4中得到的检测参数对目标水体的营养化等级进行辨别，若目标水
体为重度富营养化状态，则转至步骤s6，否则则将其辨别为检测参数所对应的营养化等级；
[0067]
s6：构建用单个水质参数表征水体富营养化程度的检测参数itlc，利用步骤s3得到的水质参数分别计算出基于叶绿素浓度或者目标水体浊度的itlc；
[0068]
s7：将目标水体的浊度与透明度做相关性分析，判别该目标水体的浊度和透明度是否高度相关，若是存在高度相关性，则转至步骤s8，否则令该目标水体中浊度的营养状态指数的相关权重为0，并根据步骤s3中反演得到的叶绿素浓度和浊度，重新计算目标水体中浊度和叶绿素浓度的营养状态指数的相关权重；
[0069]
s8：利用基于不同水质参数计算出的itlc以及其营养状态指数的相关权重，构建基于叶绿素浓度和目标水体浊度的内陆水体富营养化综合评价指数itlc(∑)；
[0070]
s9：利用内陆水体富营养化综合评价参数itlc(∑)对目标水体进行判定，最终辨别水体的富营养化等级。
[0071]
在本实施例中，地基高光谱成像仪的传感器为cmos成像仪，分辨率为200万像素，光谱范围包含420
‑
1000nm，采集速度为8
‑
20秒一次，其辅助成像为500w的高清相机，该成像仪的拍摄方向与水面法线的夹角约40
°
，且其观测平面与太阳光入射平面的夹角为135
°
，且背对太阳方向安装；
[0072]
在对水体营养化程度进行测算时，先通过地基高光谱成像仪采集目标水体的影像数据，并通过遥感图像处理软件进行处理，得到计算遥感反射率的三个参量，然后通过三个参量构建敏感波长的反射率曲线，再利用敏感波长反射率曲线提取cie颜色空间中三刺激值的特定波长反射率，并根据这三刺激值计算对应的坐标，同时通过三波段算法和波段比值法反演出目标水体中叶绿素的浓度和浊度，此时根据目标水体色度角θ与fuci查找表进行对比，得到最匹配的水体色度；
[0073]
当初步得到水体色度时则对水体进行初步判断其水体营养化程度，当水体为富营养状态时则再通过水体浊度与透明度的关系判断，对其状态进行进一步细化分类；
[0074]
若水体浊度与透明度的相关性较差，则无法用浊度代替透明度的方法来构建富营养化综合评价指数，若水体浊度与透明度的相关性较高，则可以用浊度代替透明度的方法来构建富营养化程度的辨别，即重新计算水体浊度的营养状态指数的相关权重和水体叶绿素浓度的营养状态指数的相关权重，最后进行综合评价后对富营养化状态下的水体进一步细分，从而更加利于对水体的治理分类，优先治理富营养化程度较高的水体；
[0075]
这样通过地基高光谱成像仪采集水体数据并处理，进而得到目标水体的敏感波长反射率曲线，进一步演算得到cie颜色空间中的三刺激值和坐标，配合fuci查找表和反演公式初步得到水体的色度和富营养化等级，而为了水体治理的时效性和准确性，通过对水体的浊度与透明度关系的判断，进一步区分富营养化程度下水体的等级，从而更加利于对水体的分类，及时地治理富营养化更加严重的水体。
[0076]
在本实施例中具体测算步骤如下：
[0077]
步骤s1的具体步骤为，提取的三个参量分别为标准板测量值、天空漫散射光、离水辐亮度，利用高光谱成像仪，对目标水域水体进行连续拍摄，形成包含水体信息的影像数据文件，在本实施例中遥感图像处理软件采用envi软件，再利用envi软件提取该影像数据raw源文件任意选取一波段后形成监测区域影像，并可从该影像上读取出水面上的三个参量，具体提取方法如下：
[0078]
l
p
：首先在该影像上可以观测到固定在成像仪视野范围内的标准板，打开roi面板将这整个标准板选为感兴趣区域，利用roi工具自动读取出该感兴趣区域的一系列数值，包括最大值、最小值、平均值、中位数值等(标准板测量值一般选取最大值)，同时由于通过roi开启影像数据的raw源文件，将使得读取的一系列数据包含整个波段，因此便形成了各个波段标准板测量值的最大值集合，并将其以text文档形式导出保存。
[0079]
l
sky
：首先在该影像上找到最远处的天空区域，然后打开roi面板将整个影像区域中靠近影像最顶端的天空区域选为感兴趣区域，利用roi工具自动读取出该感兴趣区域的一系列数值，包括最大值、最小值、平均值、中位数值等(天空漫散射光一般选取平均值)，同时由于通过roi开启影像数据的raw源文件，将使得读取的一系列数据包含整个波段，因此便形成了包含各个波段的天空漫散射光光谱曲线，同样将其以text文档形式导出保存。
[0080]
l
w
：首先在该影像上找到无倒影、无白帽、无波纹、无船只的水体区域，然后打开roi面板将其选择为感兴趣区域，利用roi工具自动读取出该感兴趣区域的一系列数值，包括最大值、最小值、平均值、中位数值等(离水辐亮度一般选取平均值)，同理将形成包含各个波段的离水辐亮度光谱曲线，将其以text文档形式导出保存。
[0081]
步骤s2的具体步骤为，将步骤s1得到的三种参数数据导入excel表格中，并进行通过如下两步计算形成敏感波长的反射率曲线：
[0082]
(1)
[0083]
(2)
[0084]
其中，ρ
p
为标准板的反射率；e
d
(0
+
)表示水表面入射总辐照度；ρ表示气—水界面对天空光的反射率。
[0085]
步骤s3的具体步骤为，利用步骤s2中的敏感波长反射率曲线提取计算cie颜色空间中三刺激值的特定波长反射率，其中红、绿、蓝三个波段的中心波长分别为λ
r
＝645nm、λ
g
＝555nm、λ
b
＝469nm，其对应波长反射率为r
rs
(r)、r
rs
(g)、r
rs
(b)，利用红绿蓝三波段的中心波长反射率计算出cie颜色空间中的三刺激值，并进一步计算出cie色度图中颜色所对应的坐标(m，n)，其计算公式如下，其中该公式中的所有系数值均由国际照明委员会建立的cie色空间标准得到：
[0086]
cie颜色空间三刺激值：
[0087]
cie色度图中颜色所对应的坐标：
[0088]
利用敏感波长反射率曲线，提取叶目标水体绿素浓度和浊度的特定波长反射率，其中提取叶绿素浓度和水体浊度的波长分别为λ
a
＝850nm、λ
b
＝660nm、λ
c
＝830nm、λ
d
＝660nm，其对应波长反射率为r
rs
(λ
a
)、r
rs
(λ
b
)、r
rs
(λ
c
)、r
rs
(λ
d
)，分别利用三波段算法和波段
比值法反演出目标水体叶绿素的浓度和浊度，其计算公式如下：
[0089]
三波段算法：
[0090]
波段比值法：
[0091]
步骤s4的具体步骤为，利用经验模型构建初步辨别水体营养化程度的检测参数为水体色度fuci对照表，然后利用步骤s3中得到的颜色坐标，计算出对应的色度角θ，最后将色度角与fuci对照表进行对比，查找出表中与色度角θ最接近的值，即目标水体的fuci，其中色角度θ计算公式如下：
[0092]
色度角θ：
[0093]
步骤s5的具体步骤为，首先设置二次限定条件的阈值h为0.00625，然后对目标水体的营养状态进行初步分级，分级标准如下：
[0094]
贫营养状态：fuci＜7
[0095]
中营养状态：7≤fuci＜10
[0096]
富营养状态：10≤fuci且r
rs
(r)≥0.00625，否则水体仍为中营养状态；
[0097]
将步骤s4中所得到的fuci值与步骤s5中的标准进行对比，如果10≤fuci并且目标水体在红波段的反射率值大于0.00625，则目标水体属于富营养状态，转至步骤s6，否则目标水体则为该辨别标准所判定的营养状态(贫营养状态或中营养状态)。
[0098]
步骤s6的具体步骤为，当目标水体已判别为富营养状态时，此时便对富营养状态的水体进行进一步的判别，首先构建基于单个参数评价水体富营养化等级的检测参数itlc，分别计算出基于叶绿素浓度和基于水体浊度的itlc，基于叶绿素浓度和基于水体浊度的计算公式分别如下，其中以下两个基于不同水质参数的itlc计算公式来源于《湖泊富营养化评价方法及分级标准》：
[0099]
基于叶绿素浓度的营养等级分类指数itlc：itlc(c
chl
)＝25+10.861ln
x
[0100]
基于水体浊度的营养等级分类指数itlc：itlc(turbidity)＝51.18
‑
19.4ln
x
[0101]
步骤s7的具体步骤为，将目标水体中浮漂式在线监测设备监测到的实测目标水体浊度和透明度进行相关性分析，将目标水体样本集中的水体浊度和透明度同时导入excel表格中，利用excel表格算法工具，判别水体浊度和透明度的相关性是否强烈(r2>0.8)，如果r2>0.8，则表示在目标水体中，浊度和透明度的相关性较好，可以用浊度代替透明度进行水体富营养化程度的辨别，转至步骤八，否则则表示不能用浊度代替透明度来构建富营养化综合评价指数，于是令其浊度的营养状态指数的相关权重为0；
[0102]
若可以用目标水体浑浊度代替透明度，则需要对水体浊度的营养状态指数和水体叶绿素浓度的营养状态指数的相关权重进行重新计算，其计算公式如下：
[0103]
浊度营养状态指数的相关权重计算公式：
[0104]
叶绿素浓度营养状态指数的相关权重计算公式：
[0105]
将步骤三计算得到的叶绿素浓度和浊度同时导入excel表格中，利用excel算法工具，分别计算出r
ii
和r
ij
，其中，r
ii
表示目标水体的叶绿素浓度的自相关系数，r
ij
表示目标水体叶绿素浓度与浊度的相关系数。
[0106]
步骤s8的具体步骤为，利用基于单个水质参数的营养状态指数和其对应的相关权重，构建基于叶绿素浓度和目标水体浊度的内陆水体富营养化综合评价指数itlc(∑)，其计算公式如下：
[0107]
内陆水体富营养化综合评价指数itlc(∑)：
[0108][0109]
然后将步骤s6中得到的基于单个水质参数的营养等级分类指数(itlc(c
chl
)、itlc(turbidity))带上其对应的相关权重(w
turbidity
)求和，最终得到内陆水体的富营养化综合评价指数。
[0110]
步骤s9的具体步骤为，设定内陆水体富营养化综合评价指数的判别标准：
[0111]
轻度富营养状态：0＜itlc(∑)≤60
[0112]
中度富营养状态：60＜itlc(∑)≤70
[0113]
重度富营养状态：70＜itlc(∑)
[0114]
最后将步骤s8得到的内陆水体的富营养化综合评价指数与判别标准进行对比，则得到目标水体在富营养化状态下的进一步细化分级状态。
[0115]
上述实施例仅示例性说明本发明的原理及其功效，而非用于限制本发明。任何熟悉此技术的人士皆可在不违背本发明的精神及范畴下，对上述实施例进行修饰或改变。因此，凡所属技术领域中具有通常知识者在未脱离本发明所揭示的精神与技术思想下所完成的一切等效修饰或改变，仍应由本发明的权利要求所涵盖。