水体提取方法、装置、设备、存储介质及程序产品

1.本发明涉及图像处理技术领域，尤其涉及一种水体提取方法、装置、设备、存储介质及程序产品。


背景技术：

2.地表水是指河流、湖泊、池塘与冰川等存在于地球表面的水体总称，在农业灌溉、人类生存和发展以及维持生态系统稳定等方面发挥着重要作用，所以对地表水精细地、高频率地动态监测显得尤为重要。及时、准确地获取地表水的动态信息对于政府部门制定政策、调整措施方面十分重要，在灾害预防等方面也有重要应用。随着卫星遥感技术的发展，利用遥感影像进行地表水监测成为一个重要的方式，包括光学遥感和雷达遥感。
3.现有技术中，光学遥感受到云层的影响，大多数影像由于云量较大而无法使用，所以难以实现水体连续的动态监测，更重要的是在水体识别时还会受到雪、建筑物、土壤和山体阴影等的影响，难以对水体进行精细地刻画。雷达遥感具有穿云透雾的能力，不受气候条件和光照水平的影响，可以全天时、全天候地获取数据，但由于其图像处理较复杂，且由于其侧视成像的特点会因地形的起伏造成几何畸变，如叠掩、透视收缩、阴影等，阴影的后向散射系数与水体相近，同样导致进行水体提取时阴影与水体混淆。水体提取时，如何减小阴影的影响成为一个亟待解决的问题。


技术实现要素：

4.本发明提供一种水体提取方法、装置、设备、存储介质及程序产品，用以解决现有技术中水体图像中的阴影造成水体提取结果精度低的缺陷，实现高精度精准地提取地表水体。
5.本发明提供一种水体提取方法，包括：获取目标区域的第一水体图像，以及获取所述目标区域的第二水体图像，其中，所述第一水体图像包括至少一个第一水体单元区域，所述第二水体图像包括至少一个第二水体单元区域，非水体造成的第一阴影区域分布于所述至少一个第一水体单元区域，所述非水体造成的第二阴影区域分布于所述至少一个第二水体单元区域，所述第一阴影区域和所述第二阴影区域的分布不同；将所述第一水体图像和所述第二水体图像进行叠加，获得每一个所述第一水体单元区域分别与对应的所述第二水体单元区域的重叠区域；根据至少一个所述重叠区域，确定所述目标区域的最终水体图。
6.根据本发明提供的一种水体提取方法，所述根据至少一个所述重叠区域，确定所述目标区域的最终水体图，包括：针对每一个所述重叠区域分别作如下处理：计算所述重叠区域占对应的所述第一水体单元区域的第一比例，以及计算所述重叠区域占对应的所述第二水体单元区域的第二比例；当所述第一比例大于第一阈值，且所述第二比例大于第二阈值时，提取所述第一水体图像中，所述重叠区域对应的所述第一水体单元区域；基于提取出的至少一个所述第一水体单元区域，确定所述目标区域的最终水体图。
7.根据本发明提供的一种水体提取方法，所述将所述第一水体图像和第二水体图像
进行叠加，获得每一个所述第一水体单元区域分别与对应的所述第二水体单元区域的重叠区域之后，还包括：基于每一个所述重叠区域，形成第一中间结果；叠加所述第一中间结果与预设的河湖矢量数据，将所述第一中间结果中属于河流水体的所述重叠区域剔除，其中，所述河湖矢量数据包括所述河流水体的区域分布数据；所述针对每一个所述重叠区域分别作如下处理：计算所述重叠区域占对应的所述第一水体单元区域的第一比例，以及计算所述重叠区域占对应的所述第二水体单元区域的第二比例；当所述第一比例大于第一阈值，且所述第二比例大于第二阈值时，提取所述第一水体图像中，所述重叠区域对应的所述第一水体单元区域，包括：针对剔除所述河流水体的所述第一中间结果中的、每一个所述重叠区域分别作如下处理：计算所述重叠区域占对应的所述第一水体单元区域的第一比例，以及计算所述重叠区域占对应的所述第二水体单元区域的第二比例；当所述第一比例大于第一阈值，且所述第二比例大于第二阈值时，提取所述第一水体图像中，所述重叠区域对应的所述第一水体单元区域。
8.根据本发明提供的一种水体提取方法，所述基于提取出的至少一个所述第一水体单元区域，确定所述目标区域的最终水体图，包括：基于提取出的至少一个所述第一水体单元区域，形成第二中间结果；叠加所述第二中间结果与预设的土地覆盖数据，将所述第二中间结果中属于预设土地类型的所述第一水体单元区域剔除，形成第三中间结果，其中，所述土地覆盖数据包括所述预设土地类型的区域分布数据；基于所述第三中间结果，确定所述目标区域的最终水体图。
9.根据本发明提供的一种水体提取方法，所述基于所述第三中间结果，确定所述目标区域的最终水体图，包括：叠加所述第三中间结果和所述河流水体的区域分布数据，将叠加结果作为所述目标区域的所述最终水体图。
10.根据本发明提供的一种水体提取方法，所述第一水体图像为雷达图像，所述雷达图像通过雷达设备采集所述目标区域的水体形成；所述第二水体图像为光学图像，所述光学图像通过光学设备采集所述目标区域的水体形成。
11.本发明还提供一种水体提取装置，包括：获取模块，用于获取目标区域的第一水体图像，以及获取所述目标区域的第二水体图像，其中，所述第一水体图像包括至少一个第一水体单元区域，所述第二水体图像包括至少一个第二水体单元区域，非水体造成的第一阴影区域分布于所述至少一个第一水体单元区域，所述非水体造成的第二阴影区域分布于所述至少一个第二水体单元区域，所述第一阴影区域和所述第二阴影区域的分布不同；叠加模块，用于将所述第一水体图像和所述第二水体图像进行叠加，获得每一个所述第一水体单元区域分别与对应的所述第二水体单元区域的重叠区域；确定模块，用于根据至少一个所述重叠区域，确定所述目标区域的最终水体图。
12.本发明还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现如上述任一种所述水体提取方法的步骤。
13.本发明还提供一种计算机程序产品，包括计算机程序，所述计算机程序被处理器执行时实现如上述任一种所述水体提取方法的步骤。
14.本发明提供的水体提取方法、装置、设备、存储介质及程序产品，获取目标区域的第一水体图像，以及获取目标区域的第二水体图像之后，将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域，然
后根据至少一个重叠区域，确定目标区域的最终水体图。
15.上述过程中，第一水体图像中存在由非水体造成的第一阴影区域。第二水体图像中存在由非水体造成的第二阴影区域。第一阴影区域分布在至少一个第一水体单元区域中，第二阴影区域分布在至少一个第二水体单元区域中。而由于第一阴影区域和第二阴影区域分布不同，当同一个目标区域的第一水体图像和第二水体图像叠加时，同一个水体的实际分布区域会重叠，即会存在重叠区域，但是两个水体图像中不同分布的阴影区域则不会重叠。而基于重叠区域确定最终水体图，则降低剔除第一水体图和/第二水体图中阴影的影响，提高水体提取的精度。
附图说明
16.为了更清楚地说明本发明或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作一简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
17.图1是本发明提供的水体提取方法的流程示意图；
18.图2是本发明提供的重叠区域示例图；
19.图3是本发明提供的水体提取方法的具体实施示例图；
20.图4是本发明提供的雷达影像水体初步提取结果示例图；
21.图5是本发明提供的光学影像水体初步提取结果示例图；
22.图6是本发明提供的sentinel-1与sentinel-2初步水体提取结果的交集示例图；
23.图7是本发明提供的交集剔除河流的结果示例图；
24.图8是本发明提供的sentinel-1初步水体剔除河流的结果示例图；
25.图9是本发明提供的sentinel-2初步水体剔除河流的结果示例图；
26.图10是本发明提供的percent_s1》37％、percent_s2》40％时的水体结果示例图；
27.图11是本发明提供的与林地叠加的结果示例图；
28.图12是本发明提供的与林地、草地叠加的结果示例图；
29.图13是本发明提供的与林地、草地、建筑叠加的结果示例图；
30.图14是本发明提供的sentinel-1与去除小斑块后的结果有交集的水体示例图；
31.图15是本发明提供的添加河流水体，得到最终水体图的示例图；
32.图16是本发明提供的水体提取装置的结构示意图；
33.图17是本发明提供的电子设备的结构示意图。
具体实施方式
34.为使本发明的目的、技术方案和优点更加清楚，下面将结合本发明中的附图，对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
35.随着卫星遥感技术的发展，利用遥感影像进行地表水监测成为一个重要的方式，即利用水体在不同波段上与其他地物之间的光谱差异来进行水体提取，主要有阈值法(单
波段阈值和多波段阈值)和分类器法(面向对象法、决策树法和支持向量机法等)，阈值法操作简单但精度较低，分类器法精度较高但操作复杂。光学遥感受到云层的影响，大多数影像由于云量较大而无法使用，所以难以实现水体连续的动态监测；同时在水体识别时还会受到雪、建筑物、土壤和山体阴影等的影响，难以对水体进行精细地刻画。雷达遥感具有穿云透雾的能力，不受气候条件和光照水平的影响，可以全天时、全天候地获取数据，但由于其图像处理较复杂，且由于其侧视成像的特点会因地形的起伏造成几何畸变，如叠掩、透视收缩、阴影等，阴影的后向散射系数与水体相近，导致进行水体提取时阴影与水体混淆。
36.例如，区域z地理位置偏南，其位于北纬20
°
54
′‑
26
°
24
′
，东经104
°
26
′‑
112
°
04
′
之间，总面积约23.67万km2，位于云贵高原往两广丘陵的过渡地带，属于典型的喀斯特地貌，总体地势呈现西北高东南低的走势，四周多被山地、高原环绕，呈盆地状，且山地多，平原少，具有独特的喀斯特地貌。区域z以山地丘陵为主，山地丘陵主要分布在区域z的西部、北部和东部，山地以中低山为主，山地地形约占区域z总面积的40％。区域z的山地，由于受到印度洋板块和亚欧板块的碰撞挤压，山脉多呈弧形脉状分布，走向上以东北西南走向和西北东南走向为主。区域z最小高程为-134m，最大高程为2111m，地形起伏较大。光学影像在水体识别时受到山体阴影的影响较大，简单阈值法无法分离水体和阴影；合成孔径雷达(synthetic aperture radar，sar)采用侧视成像方式，高的山峰背面由于雷达照射不到而没有回波信号，在sar图像中表现为非常黑的阴影，即雷达阴影，而阴影的后向散射系数与水体的相近，导致提取水体时水体与阴影混淆，增加了水体提取的难度。而光学影像容易受到云的影响，区域z所在区域云雨天气较多，大部分影像云量较大而无法使用。
37.因此，在同类型的多云雨山区，对地表水体的提取受山区的影响大。如何减小阴影的影响成为一个亟待解决的问题。下面结合图1-图15描述本发明的水体提取方法。
38.一个实施例中，如图1所示，水体提取方法的实现过程如下：
39.步骤101，获取目标区域的第一水体图像，以及获取目标区域的第二水体图像，其中，第一水体图像包括至少一个第一水体单元区域，第二水体图像包括至少一个第二水体单元区域，非水体造成的第一阴影区域分布于至少一个第一水体单元区域，非水体造成的第二阴影区域分布于至少一个第二水体单元区域，第一阴影区域和第二阴影区域的分布不同。
40.本实施例中，目标区域指的是需要提取水体信息的区域，例如，目标区域可以指上述中的区域z。提取目标区域的水体信息时，需要通过目标区域的水体图像来实现。但是水体图像由于收到非水体的影响，会出现阴影区域，该阴影区域与水体的实际区域混淆，则会导致水体图像呈现的各个水体单元区域不准确，例如，一个水体单元区域在实际地形上其实为一个山体的阴影，而不是一个实际水体。又例如，一个水体单元区域中的一部分在实际地形上为实际水体，而剩下的部分其实为草地阴影。
41.一个目标区域可以通过多种手段获得不同的水体图像，例如，雷达采集、光学设备采集、红外线设备采集以及其他任意一种能够采集水体分布信息的手段，第一水体图像和第二水体图像可以为上述手段中的任意两种采集得到的图像。而不同采集手段受外界因素的影响不同，使获得的水体图像中阴影区域的分布不同。
42.本实施例中，通过阴影区域分布不同的第一水体图像和第二水体图像，采用本发明提供的方法，更好地避免阴影的影响，实现目标区域水体的精细提取。
43.一个实施例中，第一水体图像为雷达图像，雷达图像通过雷达设备采集目标区域的水体形成；第二水体图像为光学图像，光学图像通过光学设备采集目标区域的水体形成。
44.更具体的，第一水体图像为哨兵1号(sentinel-1)卫星观测得到的雷达图像。sentinel-1为主动微波遥感卫星，搭载c波段合成孔径雷达传感器，能够全天时、全天候地对地面进行观测并进行雷达成像。哨兵1号(sentinel-1)卫星是欧空局(european space agency，esa)哥白尼计划的地球观测卫星，由a、b两颗卫星组成，单颗卫星的重访周期为12天，两颗互补，重访周期为6天。sentinel-1共有4种成像模式：条带成像模式(strip map mode，sm)、干涉宽幅模式(interferometric wide swath，iw)、超宽幅模式(extra wide swath，ew)、波浪模式(wave mode，wv)。
45.本实施例中，采用sentinel-1a iw成像模式的地距影像(ground range detected，grd)产品，包含vv和vh两种极化方式，其中vv指的是一种同向极化，vh指的是一种交叉极化。sentinel-1覆盖区域z一共需要12景影像，一个月能实现2～3次的重复观测。
46.第二水体图像为哨兵2号(sentinel-2)卫星观测得到的光学图像。sentinel-2属于高分辨率多光谱成像卫星，主要由2a和2b两颗卫星组成。一颗卫星的重访周期为10天，双星的重访周期为5天。sentinel-2携带有多光谱仪，其13个光谱带跨越可见光、近红外和短波红外，蓝、绿、红和近红外波段为10m分辨率。sentinel-2具有高时间分辨率和光谱分辨率，在监视地球表面方面具有集中且连续的优势。
47.本实施例中，光学影像受到云的影响，区域z云雨天气较多，选择sentinel-2过境拍摄的大部分影像云量较大而无法使用，且sentinel-2影像完整覆盖区域z需要44景影像，所以选取每三个月一张晴空的数据进行水体提取。
48.本实施例中，sentinel-2采集图像时容易收到云层的影响，尤其是类似于区域z的多云雨天气，需要经过三个月甚至六个月的时间，才能确定一张能够使用的影像来实现水体提取。由于可用影像采集消耗的时间过长，单独采用sentinel-2很难实现地表水体的动态监测。而雷达数据在多云雨的山地区域，具有全天候全天时的监测能力，可以实现在空间分辨率条件下，实现多颗卫星雷达数据达到高时间频次观测的能力，实现山区水体的动态监测能力，满足最终用户的需求。如sentinel 1a在12天内可以实现重访能力。
49.将sentinel 1采集得到的雷达图像和sentinel 2采集得到的光学图像结合处理，例如，sentinel 1每12天得到一张雷达影像，均与该时段内sentinel 2采集的光学图像，分别作为第一水体图像和第二水体图像进行处理，在保证高精度提取地表水体的同时，还能够实现动态监测。
50.本实施例中，一段时间内，sentinel-1每次采集得到第一水体图像后，即可以与该时间段sentinel-2采集得到的第二水体图像作为基础，采用本发明提供的方法进行水体提取，实现对目标区域中水体的动态观测。
51.一个实施例中，在进行水体提取之前，可以对sentinel-1采集的原始雷达数据，和sentinel-2采集的原始光学数据进行预处理，滤除图像噪音等因素的影响，提高水体提取结果的精确度。
52.一个实施例中，对原始雷达数据的预处理包括但不限于轨道校正、辐射定标、辐射地形校正、滤波、几何地形校正和/或分贝化。
53.轨道校正。卫星在运行过程中由于外界干扰或者自身原因会使飞行姿态发生轻微
抖动、偏移，从而为数据获取引入一定的系统误差。而sar元数据中提供的轨道状态参数通常不准确，所以为了更好地去除轨道引起的系统性误差，需要利用精密轨道信息数据进行轨道精确校正，以获取精确的卫星几何位置和速度信息。
54.辐射定标。辐射定标(radiometric calibration)是将图像原始记录的是遥感影像像元亮度值(digital number，dn)，转换为辐射亮度值的过程。sentinel-1的辐射定标公式如下所示：
[0055][0056]
式(1)中，dn为像素点(又称像元)的亮度值，i表示像素点位置，a为定标参数，为sigma-naught(标准后向散射系数)、gamma-naught(归一化的后向散射系数)或beta-naught(雷达反射率系数)的其中之一。本实施例中将图像校准为beta-naught产品。
[0057]
辐射地形校正。辐射地形校正是为了消除地形对遥感反射率的影响，在复杂的山区地形中，不同的地形坡度使得影像辐射失真，增加了地物光谱信息的提取难度。辐射地形校正可以消除与地形相关的辐射变化，同时保留与土地覆盖相关的辐射变化，但是入射角的影响仍然存在。这一步骤在一定程度上消除了透视收缩、叠掩和阴影的影响。
[0058]
滤波。由于sar系统本身是一个相干系统，相干波的加强和减弱叠加会引起雷达图像上出现周期性的斑点噪声，而且在影响sar图像质量的各种噪声因素中，相干噪声的影响远比其他各类噪声大，所以为了保持尽量多的图像细节信息，抑制噪声对影像解译的影响，需要进行滤波处理，本方法利用refined lee滤波器(一种改进的图像去噪滤波器)对斑点噪声进行处理。
[0059]
几何地形校正。几何地形校正是消除和改正遥感影像几何误差的过程。由于场景的地形变化和卫星传感器的倾斜，sar图像中的距离可能会发生畸变，因此会造成sar影像的一些失真，如南北方向颠倒等。本方法利用距离多普勒算法进行地形校正，减少遥感影像地形畸变的误差影响。
[0060]
线性转换为分贝
[0061]
为了避免后向散射系数值过小的问题，将预处理结果转换为对数值，使可视化和分析变得容易。转换公式如下：
[0062][0063]
式(2)中，γ0表示针对归一化的后向散射系数进行分贝化处理。
[0064]
至此，sentinel-1原始雷达数据的预处理完成。
[0065]
一个实施例中，sentinel-1原始雷达数据的预处理完成后，对预处理后的sentinel-1雷达图像进行水体初步提取。
[0066]
具体的，利用双极化水体指数(简称sdwi)对预处理后的sentinel-1雷达图像进行水体初步提取，得到sentinel-1每12天的水体图。其基本原理是基于vv和vh后向散射系数，进行波段运算，以此扩大水体与其他地物之间的差异。划分水体的分界点根据实际情况来确定，例如，一种情况下，以0为分界点，把sdwi大于0的部分划分为水体，sdwi小于0则划分为非水体。又例如，根据实验数据，选取-0.2为阈值分割点(即分界点)，将sdwi影像二值化，提取水体。sdwi表达式如下：
[0067]
sdwi＝ln(10
×
vv
×
vh)-8
ꢀꢀ
(3)；
[0068]
式(3)中，sdwi为水体提取指数，vv和vh分别表示vv极化模式和vh极化模式。
[0069]
至此，将sentinel-1雷达图像进行水体初步提取的结果作为第一水体图像。
[0070]
一个实施例中，对原始光学数据的预处理包括但不限于辐射定标、大气校正、重采样和/或波段合成。
[0071]
辐射定标和大气校正。辐射定标是消除传感器自身造成的误差，将dn值转换为辐射亮度值的过程。大气校正是将大气上层表观反射率转化为地表反射率以消除大气散射、大气吸收和大气反射等产生的误差的过程。
[0072]
重采样。大气校正后的sentinel-2数据存在三种分辨率(10m、20m和60m)，采用双线性内插法将所有波段的分辨率统一重采样为10m。
[0073]
波段合成。重采样的处理结果将sentinel-2数据存储为一个波段对应一个图像，采用波段合成工具将其合成为多波段图像，方便后续分析。
[0074]
至此，sentinel-2原始光学数据的预处理完成。
[0075]
一个实施例中，sentinel-2原始光学数据的预处理完成后，对预处理后的sentinel-2光学图像进行水体初步提取。
[0076]
具体的，本实施例中采用单波段阈值法对预处理后的sentinel-2光学图像进行水体初步提取。单波段阈值法是利用遥感影像在近红外波段上的水体光谱特性与其他地物明显不同的特征，在这个波段上水体表现出高吸收率而其他地物表现出高反射率，通过对比不同地物在该波段上的光谱曲线差异，利用试错法找到能将水体与背景地物分开的最佳阈值，其表达式一般为：
[0077]
ρ
nir
《t
ꢀꢀ
(4)；
[0078]
式(4)中，ρ
nir
代表水体在近红外波段的反射率；t代表水体提取中选取的阈值。根据实验数据，选取阈值为800。
[0079]
至此，将sentinel-2光学图像进行水体初步提取的结果作为第二水体图像。
[0080]
步骤102，将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域。
[0081]
本实施例中，当第一水体图像为雷达图像，第二水体图像为光学图像时，由于太阳照射角度的不同，以及光学和雷达成像方式不同，所以在遥感图像中，光学的阴影(对应第二阴影区域)和雷达的阴影(对应第一阴影区域)不在同一个方位，两者重叠较少，甚至不重叠，而实际水体的位置大体相同。将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域，可以实现剔除阴影区域的目的。
[0082]
本实施例中，目标区域内的单个水体，例如河、湖等，呈现分散式分布的特点。所以，第一水体图像包括至少一个第一水体单元区域，各个第一水体单元区域呈分散式分布。同理，第二水体图像包括至少一个第二水体单元区域，各个第二水体单元区域呈分散式分布。
[0083]
而对于同一个目标区域而言，理论上来说，同一个单个水体在第一水体图像和第二水体图像中，呈现的图像分布应该一致。但由于阴影的存在，导致第一水体图像和第二水体图像中，同一个单个水体的分布存在差异。
[0084]
将第一水体图像和第二水体图像叠加，同一个单个水体对应的第一水体单元区域和第二水体单元区域会出现重叠的现象。如图2所示，对于一个单个水体而言，实线圈内的区域a表示该单个水体在第一水体图像中的第一水体单元区域，虚线圈内的区域b表示该单个水体在第二水体图像中的第二水体单元区域。将第一水体图像和第二水体图像叠加后，虚线和实现相交得到的区域c则为一个重叠区域。
[0085]
将第一水体图像和第二水体图像叠加后，则可以得到至少一个重叠区域。
[0086]
步骤103，根据至少一个重叠区域，确定目标区域的最终水体图。
[0087]
本实施例中，根据至少一个重叠区域，实现目标区域水体的精细提取，获得目标区域的最终水体图。
[0088]
一个实施例中，根据至少一个重叠区域，确定目标区域的最终水体图，具体实现过程如下：针对每一个重叠区域分别作如下处理：计算重叠区域占对应的第一水体单元区域的第一比例，以及计算重叠区域占对应的第二水体单元区域的第二比例；当第一比例大于第一阈值，且第二比例大于第二阈值时，提取第一水体图像中，重叠区域对应的第一水体单元区域；基于提取出的至少一个第一水体单元区域，确定目标区域的最终水体图。
[0089]
本实施例中，以基于图2所示的一个重叠区域为例，计算区域c占区域a的第一比例，以及计算区域c占区域b的第二比例。当第一比例大于第一阈值，例如，第一阈值取值为37％，且第二比例大于第二阈值，例如，第二阈值取值为40％，则认为区域a为水体，提取出区域a。若第一比例小于或等于第一阈值，和/或，第二比例小于或等于第二阈值，则认为区域a为阴影，剔除区域a。
[0090]
基于提取出来的每一个第一水体单元区域，确定目标区域的最终水体图。
[0091]
本实施例中，第一阈值和第二阈值的比例，可以根据实验数据和/或实际情况进行设定。例如，以样本雷达图像和样本光学图像为基础，第一阈值和第二阈值分别以2％为间隔，在0至100％之间进行取值，直至预测的最终样本水体图达到预设的误差范围，则确定第一阈值和第二阈值的最终取值。当然，还可以采用其他方式确定第一阈值和第二阈值的取值，本发明的保护范围不以第一阈值和第二阈值的具体数值和取值方式为限制。
[0092]
一个实施例中，地面上的水体中，河流水体相对于湖泊来说为流动水体，会存在季节性断流、河道改变等情况。而在对目标区域的水体进行动态观测的过程中，为避免河流水体对观测的影响，采用河湖矢量数据对河流水体进行处理。具体的，将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域之后，基于每一个重叠区域，形成第一中间结果；叠加第一中间结果与预设的河湖矢量数据，将第一中间结果中属于河流水体的重叠区域剔除，其中，河湖矢量数据包括河流水体的区域分布数据。针对每一个重叠区域分别进行处理，具体的，针对剔除河流水体的第一中间结果中的、每一个重叠区域分别进行处理。
[0093]
本实施例中，河湖矢量数据以全球多方河流湖泊数据集作为源数据，结合国产高分影像数据，产生的一套全球河流、湖泊数据集，是具有可编辑性的较高精度的湖、河矢量数据集。
[0094]
本实施例中，预先将河湖矢量数据中的河流水体的区域分布数据提取出来，即预设的河湖矢量数据为河流水体的区域分布数据，叠加重叠区域形成的第一中间结果和预设的河湖矢量数据，剔除属于河流水体的重叠区域，则可以避免河流的影响。
[0095]
一个实施例中，基于第一比例和第二比例，提取出至少一个第一水体单元区域后，则可以去除大部分阴影。但是一些面积小于一定值的小斑块阴影无法去除。为了避免剔除小斑块阴影，基于提取出的至少一个第一水体单元区域，确定目标区域的最终水体图，具体实现过程如下：基于提取出的至少一个第一水体单元区域，形成第二中间结果；叠加第二中间结果与预设的土地覆盖数据，将第二中间结果中属于预设土地类型的第一水体单元区域剔除，形成第三中间结果，其中，土地覆盖数据包括预设土地类型的区域分布数据；基于第三中间结果，确定目标区域的最终水体图。
[0096]
本实施例中，土地覆盖数据可以为esa土地覆盖数据，是欧空局联合全球多家科研机构，基于sentinel-1和sentinel-2数据制作的2020年全球10米土地覆盖产品，地物类别一共分为林地、灌木、草地、耕地、建筑、荒漠、冰雪、水体、湿地、红树林、苔藓和地衣11类，该产品的总体精度为74％。
[0097]
本实施例中，小斑块阴影主要分布在预设土地类型的区域中，将第一水体单元区域形成的第二中间结果与预设的土地覆盖数据叠加，剔除属于预设土地类型的第一水体单元区域，则去除了小斑块阴影。
[0098]
本实施例中，预设土地类型根据实验数据和/或实际情况来确定。例如，通过叠加sentinel-1样本雷达图像进行水体初步提取的样本结果与esa土地覆盖数据，发现sentinel-1样本结果中误提的阴影大多分布在esa土地覆盖数据的林地、草地、建筑和荒漠类型上，所以根据esa土地覆盖数据分别提取出林地、草地和建筑类型作为掩膜，依次与sentinel-1的水体结果进行叠加分析，去除水体误提的部分。esa土地覆盖数据的荒漠定义为一年中任何时候植被覆盖面积都没有超过10％的裸露土壤、沙子和岩石，而其水体定义为一年中超过9个月被水体覆盖的地理区域，包括湖泊、水库和河流。某些细小的河流或水塘可能在汛期才有水，为了保留水体的完整性，所以不利用荒漠类型作为掩膜去除阴影。又例如，只确定一种预设土地类型，例如林地，来去除小斑块阴影。本技术的保护范围不以预设土地类型的具体内容和确定方式为限制。本实施例中，通过土地覆盖数据去除了小斑块阴影的影响，进一步提高水体提取的精确度。
[0099]
一个实施例中，当上述实施例中通过河湖矢量数据剔除了河流水体，则生成最终水体图时，将河流水体再加入最终结果中。具体的，基于第三中间结果，确定目标区域的最终水体图，实现过程如下：叠加第三中间结果和河流水体的区域分布数据，将叠加结果作为目标区域的最终水体图。
[0100]
一个具体的实施例中，如图3所示，水体提取的具体实施方式如下：
[0101]
对采用sentinel-1a iw成像模式的地距影像(ground range detected，grd)数据(采用vv和vh两种极化方式)，进行水体初步提取。具体包括预处理的轨道校正、辐射定标、辐射地形校正、滤波、几何地形校正和分贝化处理，以及利用sdwi进行阈值分割完成水体初步提取，得到第一水体图像。该第一水体图像是由sentinel-1得到的每12天的水体图。
[0102]
对采用sentinel-2得到的影响数据，进行水体初步提取。具体包括预处理的辐射定标、大气校正、重采样和波段合成，以及利用阈值法完成水体初步提取，得到第二水体图像。该第二水体图像是由sentinel-2得到的每季度(三个月)的水体图。
[0103]
将第一水体图像和第二水体图像进行叠加，得到sentinel-1和sentinel-2相交水体(即上述实施例中的第一中间结果)。通过预先从河湖矢量数据中提取的河流水体的分布
区域数据，剔除第一中间结果中的河流水体。然后计算各个静止单个水体的重叠度，例如湖泊、水库、坑塘。该重叠度通过第一比例和第二比例表示。当通过第一比例和第二比例确定一个重叠区域为阴影时，剔除该重叠区域对应的第一水体单元区域；当确定一个重叠区域为水体时，保留该重叠区域对应的第一水体单元区域。
[0104]
然后通过esa土地覆盖数据，剔除提取出的各个第一水体单元区域中，与林地、草地和建筑重叠的部分，最终得到需要的最终水体图。需要的话，可以采用预设的方法对最终水体图的精度进行验证。
[0105]
一个具体的实施例中，如图4所示，图4为表示雷达影像(sentinel-1)水体(sentinel-1_water)初步提取结果；如图5所示，图5为表示光学影像(sentinel-2)水体(sentinel-2_water)初步提取结果。叠加图4和图5，得到如图6所示的结果，图6表示sentinel-1与sentinel-2初步水体提取结果的交集(water_intersect)。
[0106]
接下来剔除河流水体，该过程可以采用河湖矢量数据来实现。图7表示交集剔除河流的结果(water_lake)，图8表示sentinel-1初步水体剔除河流的结果(s1_water_lake)，图9表示sentinel-2初步水体剔除河流的结果(s2_water_lake)。
[0107]
如图10所示，剔除河流水体后，计算单个水体分别占sentinel-1、sentinel-2初步水体的第一比例和第二比例。percent_s1表示第一比例，percent_s2表示第二比例，以percent_s1》37％、percent_s2》40％为分割点进行提取。图10表示percent_s1》37％、percent_s2》40％的水体提取结果(result)。
[0108]
如图11至图13所示，接下来去除小斑块阴影。具体的，利用esa土地覆盖数据的林地、草地、建筑类型进行叠加分析，去除小斑块阴影。以上一步处理结果为基础，图11表示与林地叠加的结果(result_tree)，图12表示与林地、草地叠加的结果(result_tree_grass)，图13表示与林地、草地、建筑叠加的结果(result_tree_grass_built)。
[0109]
如图14所示，将去除小斑块后的结果与sentinel-1水体结果做相交分析，提取出sentinel-1水体结果与去除小斑块后的结果有交集的水体，图14表示提取sentinel-1与去除小斑块后的结果有交集的水体(s1_lake_result)。
[0110]
如图15所示，将得到的sentinel-1湖泊水体与sentinel-1的河流水体合并，得到最后结果。图15表示添加河流水体，得到最终水体图(last_result)。
[0111]
本发明提供的水体提取方法，获取目标区域的第一水体图像，以及获取目标区域的第二水体图像之后，将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域，然后根据至少一个重叠区域，确定目标区域的最终水体图。上述过程中，第一水体图像中存在由非水体造成的第一阴影区域。第二水体图像中存在由非水体造成的第二阴影区域。第一阴影区域分布在至少一个第一水体单元区域中，第二阴影区域分布在至少一个第二水体单元区域中。而由于第一阴影区域和第二阴影区域分布不同，当同一个目标区域的第一水体图像和第二水体图像叠加时，同一个水体的实际分布区域会重叠，即会存在重叠区域，但是两个水体图像中不同分布的阴影区域则不会重叠。而基于重叠区域确定最终水体图，则降低剔除第一水体图和/第二水体图中阴影的影响，提高水体提取的精度。
[0112]
尤其对于以山地丘陵地形为主的山区而言，例如，上述实施例中提到的区域z，山体会导致雷达图像和光学图像中存在大量阴影，而本发明提供的方法，能够有效提取水体
区域，避免阴影的影响，大大降低了山地地形对水体提取的影响，实现针对山区高精度精准地提取地表水体。
[0113]
当然，对于一个地域由于其他因素，例如，云层、林地、草地等，导致水体提取困难时，若能够获得阴影区域分布不同的第一水体图像和第二水体图像，同样可以采用本发明提供的方法实现地表水体的精细提取。
[0114]
下面对本发明提供的水体提取装置进行描述，下文描述的水体提取装置与上文描述的水体提取方法可相互对应参照。如图16所示，水体提取装置包括：
[0115]
获取模块1601，用于获取目标区域的第一水体图像，以及获取目标区域的第二水体图像，其中，第一水体图像包括至少一个第一水体单元区域，第二水体图像包括至少一个第二水体单元区域，非水体造成的第一阴影区域分布于至少一个第一水体单元区域，非水体造成的第二阴影区域分布于至少一个第二水体单元区域，第一阴影区域和第二阴影区域的分布不同；
[0116]
叠加模块1602，用于将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域；
[0117]
确定模块1603，用于根据至少一个重叠区域，确定目标区域的最终水体图。
[0118]
一个实施例中，确定模块1603，具体用于针对每一个重叠区域分别作如下处理：计算重叠区域占对应的第一水体单元区域的第一比例，以及计算重叠区域占对应的第二水体单元区域的第二比例；当第一比例大于第一阈值，且第二比例大于第二阈值时，提取第一水体图像中，重叠区域对应的第一水体单元区域；基于提取出的至少一个第一水体单元区域，确定目标区域的最终水体图。
[0119]
一个实施例中，水体提取装置还包括精处理模块1604。精处理模块1604，具体用于将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域之后，还包括：基于每一个重叠区域，形成第一中间结果；叠加第一中间结果与预设的河湖矢量数据，将第一中间结果中属于河流水体的重叠区域剔除，其中，河湖矢量数据包括河流水体的区域分布数据；针对每一个重叠区域分别作如下处理，包括：针对剔除河流水体的第一中间结果中的、每一个重叠区域分别作如下处理。
[0120]
一个实施例中，确定模块1603，具体用于基于提取出的至少一个第一水体单元区域，形成第二中间结果；叠加第二中间结果与预设的土地覆盖数据，将第二中间结果中属于预设土地类型的第一水体单元区域剔除，形成第三中间结果，其中，土地覆盖数据包括预设土地类型的区域分布数据；基于第三中间结果，确定目标区域的最终水体图。
[0121]
一个实施例中，确定模块1603，具体用于叠加第三中间结果和河流水体的区域分布数据，将叠加结果作为目标区域的最终水体图。
[0122]
一个实施例中，第一水体图像为雷达图像，雷达图像通过雷达设备采集目标区域的水体形成；第二水体图像为光学图像，光学图像通过光学设备采集目标区域的水体形成。
[0123]
图17示例了一种电子设备的实体结构示意图，如图17所示，该电子设备可以包括：处理器(processor)1701、通信接口(communications interface)1702、存储器(memory)1703和通信总线1704，其中，处理器1701，通信接口1702，存储器1703通过通信总线1704完成相互间的通信。处理器1701可以调用存储器1703中的逻辑指令，以执行水体提取方法，该方法包括：获取目标区域的第一水体图像，以及获取目标区域的第二水体图像，其中，第一
水体图像包括至少一个第一水体单元区域，第二水体图像包括至少一个第二水体单元区域，非水体造成的第一阴影区域分布于至少一个第一水体单元区域，非水体造成的第二阴影区域分布于至少一个第二水体单元区域，第一阴影区域和第二阴影区域的分布不同；将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域；根据至少一个重叠区域，确定目标区域的最终水体图。
[0124]
此外，上述的存储器1703中的逻辑指令可以通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：u盘、移动硬盘、只读存储器(rom，read-only memory)、随机存取存储器(ram，random access memory)、磁碟或者光盘等各种可以存储程序代码的介质。
[0125]
另一方面，本发明还提供一种计算机程序产品，所述计算机程序产品包括计算机程序，计算机程序可存储在非暂态计算机可读存储介质上，所述计算机程序被处理器执行时，计算机能够执行上述各方法所提供的水体提取方法，该方法包括：获取目标区域的第一水体图像，以及获取目标区域的第二水体图像，其中，第一水体图像包括至少一个第一水体单元区域，第二水体图像包括至少一个第二水体单元区域，非水体造成的第一阴影区域分布于至少一个第一水体单元区域，非水体造成的第二阴影区域分布于至少一个第二水体单元区域，第一阴影区域和第二阴影区域的分布不同；将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域；根据至少一个重叠区域，确定目标区域的最终水体图。
[0126]
又一方面，本发明还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现以执行上述各方法提供的水体提取方法，该方法包括：获取目标区域的第一水体图像，以及获取目标区域的第二水体图像，其中，第一水体图像包括至少一个第一水体单元区域，第二水体图像包括至少一个第二水体单元区域，非水体造成的第一阴影区域分布于至少一个第一水体单元区域，非水体造成的第二阴影区域分布于至少一个第二水体单元区域，第一阴影区域和第二阴影区域的分布不同；将第一水体图像和第二水体图像进行叠加，获得每一个第一水体单元区域分别与对应的第二水体单元区域的重叠区域；根据至少一个重叠区域，确定目标区域的最终水体图。
[0127]
以上所描述的装置实施例仅仅是示意性的，其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性的劳动的情况下，即可以理解并实施。
[0128]
通过以上的实施方式的描述，本领域的技术人员可以清楚地了解到各实施方式可借助软件加必需的通用硬件平台的方式来实现，当然也可以通过硬件。基于这样的理解，上述技术方案本质上或者说对现有技术做出贡献的部分可以以软件产品的形式体现出来，该计算机软件产品可以存储在计算机可读存储介质中，如rom/ram、磁碟、光盘等，包括若干指
令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行各个实施例或者实施例的某些部分所述的方法。
[0129]
最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。