一种对遥感影像统计回归降尺度并提高水体清晰度的方法

1.本发明涉及一种对遥感影像统计回归降尺度并提高水体清晰度的方法，涉及遥感图像处理领域，尤其涉及水体识别指数、遥感计算河道流量、统计回归降尺度等相结合的技术领域。


背景技术：

2.利用遥感影像来获得河道流量数据是目前科研和生产中十分热门的手段，其拥有快速、便捷、成本低廉的特征。从1996年有学者提出归一化水体指数(ndwi)以来，许多研究和工作围绕精确测量地表的水体分布展开。过去大量的工作重点在于通过测量间接参数来推算水的分布，现在更多的研究工作集中在直接测量例如河道中的水量、流速等，通过对于卫星得到的遥感影像进行光谱信息的运算能够精确地分离出代表水体的信息，帮助确定水信息数据。目前拥有多种不同属性的卫星，其得到的数据也不同，每种卫星最小能分辨的单位用空间分辨率来表示，空间分辨率越高的卫星数据越清晰，一年中能够对同一个地方扫过获得的数据量更多，但是对同一个地方重新扫过的时间周期也更长。因此高分辨率数据重访周期长，低分辨率数据重访周期短，不同的数据优劣势不同，将这些数据结合使用能够提升数据使用效率。
3.同时，大部分对地面的卫星主要依靠光的反射也就是电磁波，非常容易受到大气层中的云、雾等干扰，出现数据缺失或看不清图像，不同的地面物体之间数值差别不大的情况，同样影响数据的使用效率。可以采用计算机算法程序来降低云的影响，这种算法被称为去云算法，但是它对厚云层效果好，对于薄云薄雾效果不佳。
4.一般云雾主要出现在热带、亚热带地区，比如云贵高原等，由于云雾分布不均匀，有时候有遮挡，有时候没有，因此光使用一种遥感影像会导致数据不足，结合多种影像可以互相补充，利用各种数据的优势。在统计学领域有一种回归方法，叫做统计回归降尺度，可以将多种数据结合起来，保留各自的优势，例如高空间分辨率数据的高清晰度，和地空间数据的高数据量。借助这种方法能够将数量更多的低空间分辨影像提高清晰度，从而提高数据利用效率。
5.目前能够降低云雾等对图像清晰度影响的技术方法按图像特征可分为：(1)单幅影像：通过去云算法，获取受到云层影响的图像上的光谱信息，筛选出云层的光谱特征，或直接采用滤波方法，对光谱中特定的波段去除，来进行对遥感影像或者ndwi图像计算处理，降低云层的影响；(2)多幅影像：利用多幅影像中对同一地点的不同遥感影像，将没有云雾的部分提取出来，合称为一幅没有云雾的图像。【其中，ndwi即归一化水体指数：1996年，mcfeeters提出了归一化水体指数(ndwi)，能够精确有效地获取遥感图像中水体分布信息，如水体面积、平均河道宽度、河流长度等，进而用于遥感演算河道流量。其原理是通过将卫星影像的绿色波段与近红外波段进行线性运算，能够凸显水体与其他地物光谱特性的差异，从而更容易在图像中识别出水体。】
6.然而上述方法仍存以下问题：(1)去云算法能够解决部分有云层分布的图像，对于
整个图像都存在云雾或者云雾对整个图像影响较为接近的情况不能很好的提升清晰度；(2)多幅影像只能得到一个影像结果，对于识别河流水体等目的来说结果数据减少了，不利于进一步工作的展开；(3)去除的影像虽然云雾影响降低了，但是仍然是不同影像的单独使用，没有合理利用不同数据各自的优势。
7.本发明中涉及的概念和相关技术，简要概括如下：
8.统计回归：通过对两组不同的数据之间建立统计学关系，得到两组数据之间的数学公式，可以输入其中一组数据中的某一个数值，来得到另一组数据中对应的数值。
9.水利几何关系：根据河流在实际中的形状，用数学函数对其进行模拟，并模拟得到其中物理量如河宽、水深、流量等之间的关系。


技术实现要素：

10.本发明的目的在于提供一种对遥感影像统计回归降尺度并提高水体清晰度的方法，以解决目前大量的遥感影像都受到云雾等的影响导致清晰度不够，同时这些影像各自存在数量不足、空间分辨率低等问题，以及数据之间没有互相补充，充分发挥各自的优势，导致用来识别计算地表水体时的数据来源不足的问题。
11.本发明的目的在于提供一种对遥感影像统计回归降尺度并提高水体清晰度的方法，具体通过在两种不同空间分辨率的遥感影像之间构建统计回归关系，并借助构建得到的关系将粗空间分辨率的遥感影像代入计算，计算得到新的细空间分辨率的影像，来补充原有的细分辨率影像，并根据两种数据受到云雾影响的不同程度，提高受云雾影响更严重的数据的清晰度，提高水体识别程度。
12.本发明的另一个目的在于提供一种河道流量计算的方法，是基于上述得到的细空间分辨率的影像，对遥感图像中河道流量计算的方法。
13.所述一种对遥感影像统计回归降尺度并提高水体清晰度的方法，具体包括：
14.s101、图像统计回归关系计算
15.该步骤的目的在于：在两种不同分辨率的数据之间建立关系，使得粗分辨率数据能够通过计算公式得到细分辨率的数据值。具体过程如下：
16.读取两种不同分辨率的ndwi数据：1)粗分辨率数据集a与细分辨率数据集b，将粗分辨率数据集a中的每一份数据a’与细分辨率数据集b中每一份具有同一时间属性的数据b’一一对应组成一对数据；2)把粗分辨率数据集a中每一份数据的每一个大格点按照细分辨率数据集b的分辨率分解成和细分辨率数据集b一样大小的小格点，并把原粗分辨率数据集a的大格点的值赋给分解后的小格点；3)将新的每一对数据a与数据b格点一一对应。借助统计线性回归函数(具体计算公式见公式(1)至(3))模拟两个函数之间的线性相关系数coe，把数据a与coe相乘后再与数据b相减得到残差数据res。
17.obj
i,j
＝coe
i,j
×
ori
i,j
res
i,j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(1)
18.res
i,j
＝obj
i,j-coe
i,j
×
ori
i,j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2)
19.reg
i,j
＝coe
i,j
×
ori
i,j
res
i,j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
20.式中obj为数据b，ori为数据a，coe为数据b关于数据a的相关系数，res为残差影像，reg为回归后影像，i,j分别对应影像中格点的行、列序号。
21.s102、图像降尺度处理：该步骤的目的在于：借助步骤s101确定的粗、细分辨率数
据之间的关系，将同样的粗分辨率数据代入计算公式中，即可获得其对应的细分辨率数据，即将粗分辨率的遥感影像数据处理为精细分辨率数据，以提高同一面积内数据的丰富程度，便于进一步计算。具体过程如下：
22.根据步骤s101的统计回归关系计算程序，将需要降尺度的粗分辨率数据集c每一份数据处理成与细分辨率数据集b统一分辨率后代入公式(1)中，计算得到需要降尺度的粗分辨率数据集c对应的新的细分辨率数据集d，并对新的细分辨率数据集d内的数据重新统计所有格点上的数值，按照每个格点数值与所有格点的平均值的差除以统计值中最大值与最小值的差，得到的商作为每个格点的新数值，来对结果归一化处理，便于计算与查看结果，具体归一化公式如下：
[0023][0024]
式中p为格点新数值，v为格点原始数值，v
ave
为每份数据中的所有格点影像值的平均值，v
max
为统计值中最大值，v
min
为统计值中最小值。
[0025]
本发明方法进一步包括：
[0026]
s103、河道流量计算：该步骤的目的在于：通过现有的遥感计算流量原理，构建出河流的模拟形态和流量计算公式，完成从水面面积到河道流量的计算过程。借助现有的遥感测流计算方法，输入无人机测量得到的数字表面模型数据、实测水深、坡降数据，构建水力几何关系，并输入糙率与流速数据，能够拟合出得到河流断面形状。通过谷歌地球引擎计算平台(gee，https://earthengine.google.com/)获取站点的历史遥感影像并计算成ndwi图像值，对图像中选取水体、陆地训练区，批量获取水体阈值，在训练区中得到河道内的表示水体的格点以及分布区域，并统计计算得到在河道内的水体宽度，带入所述的水利几何关系，计算得到河道流量值。
[0027]
s104、水体清晰度评估
[0028]
计算以下两个数值来评估图像中水体清晰度：
①
通过无人机在检测站点拍摄生成的高精度低空影像来框选出河流范围，提取训练区学习的水体ndwi阈值，对范围内的水体像元数进行统计，计算得到河道范围内的水体面积，并从河道范围内提取平均河道宽度，与水体面积相除得到河道水体长度，计算其在河道整体长度中的相对河长占比l。
②
通过统计在河道内的水体数量，除以总水体面积，与1做差值，计算得到河道中水体连续性值c，具体计算公式如下：
[0029][0030]
li＝ai÷
(rwii×
li)
ꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0031]
式中i为对应河流序号，ni表示i河的连续水体斑块数，ai表示i河的水体总格点数，rwi表示河宽，lhi表示河道整体长度。
[0032]
通过依次上述四个步骤后，可以将粗分辨率的遥感影像数据处理为精细分辨率数据，并计算得到图像中河道内的流量数据。
[0033]
本发明的有益效果：
[0034]
本发明方法通过建立两种遥感影像数据之间的统计回归关系，可以借助受云雾影响较小的数据来计算出对应区域内另一种数据的值，提升图像中受云雾影响的部分区域的清晰度；图像之间一一对应得到统计回归关系，每一幅需要处理的遥感影像在应用本发明
的方法后都能得到一个对应的结果，不会减少数据量；应用本发明方法后原本分辨率较粗的遥感影像数据分辨率会提升，而同时其结果也能与细分辨率遥感影像数据互相补充使用，保留了两种数据各自的优势。
附图说明
[0035]
图1所示为本发明方法流程示意图
[0036]
图2a、b河流一使用本发明方法前后水体清晰度情况，a为原始影像，b计算处理后影像；
[0037]
图3a、b河流二使用本发明方法前后水体清晰度情况，a为原始影像，b计算处理后影像；
[0038]
图4a、b河流三使用本发明方法前后水体清晰度情况，a为原始影像，b计算处理后影像；
[0039]
图5a、b河流四使用本发明方法前后水体清晰度情况，a为原始影像，b计算处理后影像；
[0040]
图6a、b河流五使用本发明方法前后水体清晰度情况，a为原始影像，b计算处理后影像；
[0041]
图7a、b河流六使用本发明方法前后水体清晰度情况，a为原始影像，b计算处理后影像；
[0042]
图8a、b河流七使用本发明方法前后水体清晰度情况，a为原始影像，b计算处理后影像；
[0043]
图9a为上述7个河流的水体连续性；
[0044]
图9b为上述7个河流的相对河长占比。
具体实施方式
[0045]
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合附图对本发明的各实施方式进行详细的阐述。然而，本领域的普通技术人员可以理解，在本发明各实施方式中，为了使读者更好地理解本发明而提出了许多技术细节。但是，即使没有这些技术细节和基于以下各实施方式的种种变化和修改，也可以实现本发明所要求保护的技术方案。
[0046]
一种对遥感影像统计回归降尺度并提高水体清晰度的方法，具体包括：
[0047]
s101、图像统计回归关系计算：本发明对两种具有不同空间分辨率的ndwi等遥感数据进行统计回归降尺度：其中按照空间分辨率高低将其分为粗分辨率数据集a与细分辨率数据集b，每种数据集中按照时间顺序排列有若干遥感影像图片，两种数据集中遥感影像的空间范围均一致。首先按照时间顺序，依次选择a中的图像a1’、a2’、a3’至an’，并从b中挑选与其时间属性相同的数据b1’、b2’、b3’至bn’，将a1’、b1’组成第1组数据，a2’、b2’组成第2组数据，直至an’、bn’组成第n组数据。其次将每组数据中来自a的数据a1’～an’在空间上重新采样，采样格点大小细分辨率数据集b的空间最小单元，得到新数据a1～an。随后对每一组数据都进行如下操作：按空间分布顺序，从第一行第一列格点开始直至最后一行最后一列格点，依次选择a与b两张遥感影像的对应格点中的影像值，分别记来自影像a的格点值为ori
i,j
，记来自影像b的格点值为obj
i,j
，其中i与j分别代表行序号与列序号。对于每个格
点，每一组数据组中均可得到一个ori值与一个obj值，n组数据一共能获得n个ori值与n个obj值，将这n个对应的ori值与n个obj值带入线性回归公式中，组成n维方程组，以计算得到方程组的回归参数coe
i,j
与残差参数res
i,j
，每个格点均能得到各自的回归参数与残差参数，具体回归方程见下公式(1)～(2)。根据每个格点各自的回归参数与残差参数，即可带入降尺度方程中，构成统计回归降尺度计算公式，按格点一共能得到i
×
j个计算公式，具体降尺度方程见下公式(3)。
[0048]
obj
i,j
＝coe
i,j
×
ori
i,j
res
i,j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(1)
[0049]
res
i,j
＝obj
i,j-coe
i,j
×
ori
i,j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2)
[0050]
reg
i,j
＝coe
i,j
×
ori
i,j
res
i,j
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
[0051]
式中obj为数据b，ori为数据a，coe为数据b关于数据a的相关系数，res为残差影像，reg为回归后影像，i,j分别对应影像中格点的行、列序号。
[0052]
s102、图像降尺度处理：在建立统计回归降尺度公式后，图像中的每一个格点都有一个对应的公式，通过将需要计算的影像值代入公式即可得到新的影像值。选取存在图像不清晰，水体模糊的遥感影像数据集c，按时间顺序依次选取其中的影像c1，c2，至cn。每一幅影像同s101步骤中一样，按照细分辨率影像的空间最小单元重新采样。随后将重新采样完成后的影像按照各点顺序，依次将其影像值代入对应格点的公式中，得到新的影像值，按时间顺序计算得到新的影像d1，d2至dn，组成降尺度后的数据集d。对新数据集d中的每一张影像中所有格点的影像值进行统计，计算出其平均值v
ave
、最大值v
max
、最小值v
min
，并通过归一化公式将每个格点数值代入公式计算得到归一化结果，将数值分布更为平均，具体归一化公式见下公式(4)。
[0053][0054]
式中p为格点新数值，v为格点原始数值，v
ave
为每份数据中的所有格点影像值的平均值，v
max
为每份数据中的最大值，v
min
为每份数据中的最小值。
[0055]
s103、河道流量计算：借助归一化后的细分辨率新遥感影像数据集e，采用遥感图像测流方法来获取河流相关参数。利用事先测量的、与数据集中的遥感图像空间范围一致的高精度无人机影像，来确定无人机影像中河流与陆地所在范围，定为water与land。将河流范围water与陆地范围land同样叠加到数据集e上，框选出河流与陆地应当存在的区域，来读取该影像中陆地和河流对应的平均影像值v
land
与v
water
，并通过v
land
与v
water
这两个河流与陆地的影像阈值来提取出图像中河流的分布面积area与平均宽度rivwidth。将事先实际测量得到的对应区域的水深、流速、糙率、坡降参数带入水力学流量计算公式曼宁公式，见下公式(5)，通过当天测量已知的流量结果qi，与清晰的无人机图像测量得到的riverwidthi，构建出河道流量与河流宽度的关系。最终将对数据集e中的每张影像e1，e2至en，依次提取出对应的水体宽度riverwidth1至riverwidthn，带入构建好的河宽-流量关系式，计算得到每张影像对应的流量值q。
[0056][0057]
式中k是转换系数，这里取1，n是糙率，a是过水断面面积，p是湿周长度，s是水力坡度，这里用河道坡降代替。
[0058]
s104、水体清晰度评估
[0059]
计算以下两个数值来评估图像中水体清晰度：
①
通过无人机在检测站点拍摄生成的高精度低空影像来框选出河流范围，提取训练区学习的水体ndwi阈值，对范围内的水体像元数进行统计，计算得到河道范围内的水体面积，并从河道范围内提取平均河道宽度，与水体面积相除得到河道水体长度，计算其在河道整体长度中的相对河长占比l。
②
通过统计在河道内的水体数量，除以总水体面积，与1作差值，计算得到河道中水体连续性值c，具体计算公式如下：
[0060][0061]
li＝ai÷
(rwi×
lhi)
ꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0062]
式中i为对应河流序号，ni表示i河的连续水体斑块数，ai表示i河的水体总格点数，rwi表示河宽，lhi表示河道整体长度。
[0063]
图2a、b至图8a、b为名为0501、0502、0401、0402、0404、0405、0406的7个河流使用本发明方法前后水体清晰度情况；每组中a图代表原始影像，b图代表计算处理后影像，河道中浅色为所识别出的水体范围，深色为陆地范围，图中l与c分别代表相对河长占比与水体连续性，可以代表水体的分布情况；图9a代表7个河流的水体连续性，图9b代表7个河流的相对河长占比。其中深色柱代表原始影像中得到的水体连续性指数及相对河长占比，浅色柱代表应用本方法后的影像得到的相应结果。