CN101571923B - 基于智能蚁群算法的遥感图像水系网半自动的提取方法 - Google Patents

Abstract

本发明公开提出了一种基于智能蚁群算法的遥感图像水系网半自动的提取方法，该方法包括如下步骤：A、利用多光谱遥感数据初次提取水系网，并得到水体掩膜图初步结果；B、分别建立5*1矩阵、南北矩阵、东西矩阵，判定大型湖泊并同时去除大型湖泊；C、在对水系网全局搜索的基础上，自动运行蚁群算法，优化小断流；D、建立一个3*3搜索矩阵，采用人工交互方式，并利用智能蚁群算法优化大断流；E、利用遥感图像处理系统软件ENVI中矢量转换算法完成水系网矢量图的绘制。该方法生成了水系网图像，提高了水系网提取的准确度；避免湖泊等大型水体对河流等线状地物提取的干扰；并可准确完成部分小断流的识别和优化，比人工判定更加快速、准确；最后智能蚁群算法转为人工交互方式计算，结合人工判定的主观能动性，完成大断流的识别和优化，并可降低此算法的计算复杂度。

Description

基于智能蚁群算法的遥感图像水系网半自动的提取方法

技术领域

本发明涉及一种利用智能蚁群算法实现在中低分辨率遥感图像中水系网半自动的提取方法，属于计算机图像处理和模式识别技术领域。

背景技术

一般在世界上，独流入海的河流极少，大多数河流是由大大小小的支流汇合到一条主流中，从而形成一个比较复杂水系统，称为水系网。因为水系网是一个非常复杂的研究系统，所以水系网的特征提取工作一般采用比较常规手段---即依靠人工目视判读。随着遥感卫星数据的海量获取，目视手段已经不能满足大量解译工作的要求，因此传统的目视解译遥感图像的方法正慢慢被计算机自动解译所取代。随着计算机技术和图形学的发展，利用计算机从遥感影像中自动提取目标地物已具有一定的理论基础，使遥感信息的自动化提取成为可能。因此，利用计算机自动提取水系网将成为遥感影像特征提取的必然趋势。目前，对于水系网提取方法主要以基于光谱信息水系网提取法、基于空间特征的提取方法和基于数学形态学的特征提取方法为主。然而，上述三种方法分别存在不足，对于常用的水系网提取方法，基于光谱信息水系网提取方法提取水系网的结果具有明显的断点和干扰点，尤其在中低分辨率遥感图像上更为明显；对于基于空间特征的水系网提取方法，其数学建模比较复杂，不利于计算机参与智能解译；对于基于数学形态学的特征提取方法，其在断点、大型湖泊的判读具有明显的优势，但其在中低分辨率遥感图像应用中具有计算复杂度高，缺少智能性等缺点，不利于在比较大的地理区域内独立进行水系网的判读。

发明内容

为了解决以上诸多技术问题，本发明的目的是在充分利用多光谱判读和数学形态学提取方法的优点并结合智能蚁群算法的自组织性和分布式计算的特点的基础上，提供一种基于智能蚁群算法的遥感图像水系网半自动的提取方法。进而，克服了随着遥感卫星数据的海量获取，人工目视判读无法完成准确、高效判读工作的缺点。

本发明的技术解决方案如下：

为了达到上述目的，本发明的技术方案如下：

上述基于智能蚁群算法的遥感图像水系网半自动的提取方法包括如下几个步骤：

A、利用多光谱遥感数据初次提取水系网，得到水体掩膜图的初步结果：

首先利用近红外波段和红光波段数据计算归一化植被指数NDVI，再对其直方图进行分析，确定水体判断阈值K₁；其次，根据水体的物理特性，利用中分辨率成像光谱仪Moderate-resolution Imaging Spectroradiometer地面反射率产品MODIS09的第三波段BAND3和第二波段BAND2波段数据，计算其直方图，并分别取得第三波段BAND3直方图的陡降拐点坐标K₂及第二波段BAND2直方图的陡升拐点坐标K₃，确定水体提取的阈值，进而利用布尔公式，RIVER_MASK＝(NDVI≤K₁)AND(BAND3≥K₂)AND(BAND2≤K₃)，得到水体掩膜图的初步结果，RIVER_MASK：提取水体掩膜图的初步结果图像。

B、建立5*1矩阵、101*1南北矩阵、1*101东西矩阵，判定大型湖泊和完成去除大型湖泊，其具体步骤如下：

B1、建立一个小型矩阵，即5*1的矩阵，然后对全幅水体掩膜图进行全局搜索边界点；

B2、建立一个101*1南北矩阵和一个1*101东西矩阵，分别承担从南北方向和东西方向搜索大型湖泊的任务；

B3、遍历整幅图像，分别设立南北、东西阈值。利用101*1南北矩阵和1*101东西矩阵，分别在两个方向上检测所有像元等于1的个数，记录数据，用其数据除以101，若这个结果值大于0.9，且小于1.0，则分别记录101*1南北矩阵和1*101东西矩阵中心点的位置，并将其像素值赋值为2；再次检验由101*1南北矩阵和1*101东西矩阵赋值为2的像素位置，选择位置公共的像素点，并重新赋值为2，其他值按布尔公式得出，即若只有一个为2的值，则恢复像素点值为1，而值都为0的像素点的值不变；

B4、以像素值为2的像素点(湖面中心点)为圆心，以100为半径，其圆形内最北边，最南边，最东边，最西边的边界点按下式计算：E_col ：最东边像素点的列数；W_col ：最西边像素点的列数；N_row：最北边像素点的行数；S_row；最南边像素点的行数。

如果Rate大于高宽比的阈值，则算法继续，否则，跳出算法，重设阈值，Rate为圆内高宽比值。

B5、最后以像素值为2的像素点为圆心，以100为半径，其圆形内所有值清零，即完成大型湖泊的去除。

C、大型湖泊去除以后，重新对水系网进行全局搜索，自动运行智能蚁群算法，优化小断流。

D、然后，在建立一个3*3搜索矩阵基础上，采用人工交互方式标定未被智能蚁群算法所处理大断流，并在小范围内运行智能蚁群算法，并在蚁群信息素指引下，优化大断流。

E、最后，利用遥感图像处理系统软件ENVI中的矢量转换算法，完成水系网矢量图的绘制。

本发明的基于智能蚁群算法的遥感图像水系网半自动的提取方法的优点在于：(1)充分利用了光谱信息和图像代数的优点，取其长处，提高了水系网的提取准确度；(2)利用数学图形学原理，引入了两个简单的单向矩阵，实现了大型内陆湖泊准确判定即对大面积片状物体的标定，从而避免了在后续步骤中大型湖泊等对河流等线状地物提取的干扰；(3)利用智能蚁群算法的自组织性和自启发性，可准确代替图像代数中的膨胀和腐蚀算法并且可自动完成部分小断流的识别和优化，比传统算法具有智能性且比人工判定更加快速、准确；(4)与传统算法相比，保留其人工判断的主观能动性，即在最后几步中算法转成半自动型计算，即人为引入标定点，同时再次利用智能蚁群算法，使算法充分融合了计算机的准确性和人工的主观能动性，避免了中低分辨率遥感图像中由于云层等因素干扰而产生的大间隔的断点因智能蚁群算法自主寻优性而忽略的缺点，且引入部分人工判定可降低一定的算法的计算复杂度。

附图说明

图1为基于智能蚁群算法的遥感图像水系网半自动的提取方法的总流程框图；

图2为利用归一化植被指数NDVI和多光谱图像计算初步提取水体掩膜图的流程框图；

图3为利用南北矩阵、东西矩阵判定大型湖泊并消除大型湖泊的流程框图；

图4为利用智能蚁群算法优化水系网的算法流程框图。

具体实施方式

下面根据图1至4给出本发明的基于智能蚁群算法的遥感图像水系网半自动的提取方法的一个实施例，要指出的是，所给出的实施例是为了说明本发明方法的技术特点和功能特点，使能更易于理解本发明，而不是用来限制本发明的范围。

先请参照图1，它是基于智能蚁群算法的遥感图像水系网半自动的提取方法的总流程框图，运行步骤如下：

如图2，利用多光谱遥感数据初次提取水系网，得到水体掩膜图的初步结果，其具体步骤如下：

A1、输入中分辨率成像光谱仪Moderate-resolution Imaging Spectroradiometer地面反射率产品MODIS09的第一波段BAND1、第二波段BAND2和第三波段BAND3遥感数据，同时计算归一化植被指数NDVI，即

得到归一化植被指数NDVI遥感图像；

A2、分别对归一化植被指数NDVI、第二波段BAND2和第三波段BAND3数据进行直方图分析，确定水体提取的阈值，利用归一化植被指数NDVI，选定其直方图第一个陡升拐点坐标值K₁，第二波段BAND2同样选定直方图第一个陡升拐点坐标值K₃，而第三波段BAND3选定直方图第一个陡降拐点坐标值K₂，保存三个初始阈值K₁，K₂，K₃；

A3、根据水体的物理特性，即水体在近红外波段具有低反射率及在蓝光波段具有高反射率的特性，利用数学布尔公式得出了初始的水体掩膜图，即RIVER_MASK＝(NDVI≤K₁)AND(BAND3≥K₂)AND(BAND2≤K₃)。若满足条件，则将单像元值赋值为1，否则赋值为0，从而得到一个初步的水体掩膜图；

分别建立5*1矩阵、101*1南北矩阵和1*101东西矩阵，判定大型湖泊和去除大型湖泊，如图3所示，其具体步骤如下：

B1、首先，输入初步提取水体掩膜图，同时建立一个小型矩阵，即5*1矩阵，然后对全幅水体掩膜图进行全局搜索边界点，若满足以下条件：对于5*1矩阵，如果前四个像元值为0，而最后一个像元为1，则标记最后一个像元位置并将其赋值为2；同时，如果前四个像元值都为1，而最后一个值为0，则也标记最后一个像元位置并将其赋值为2，即利用5*1矩阵粗略地完成水体地物的边界检测并得到相应的位置信息；

B2、建立一个101*1南北矩阵和一个1*101东西矩阵，分别承担从南北方向和东西方向搜索湖泊的任务：

大型湖泊绝大部分在遥感图像上呈现出片状结构，相对于狭长的内陆河流，其南北宽度与东西宽度之比更接近于1，基于上述原理，分别建立一个南北方向上搜索矩阵，即101*1的纵向矩阵，与一个东西方向上的搜索矩阵，即1*101横向矩阵。

B3、遍历整幅图像，分别设立南北方向、东西方向阈值：已经建立两个方向上的检测矩阵，在101*1和1*101矩阵中检测所有像元等于1的个数，记录数据，并用其数据除以101，若这个结果值大于0.9，且小于1.0，则分别记录101*1南北矩阵和1*101东西矩阵中心点的位置，并将其像素值赋值为2，再次检验由101*1南北矩阵和1*101东西矩阵赋值为2的像素位置，选择位置公共的像素点，并重新赋值为2，其他值利用下面布尔计算得出，即若只有一个为2的值，则恢复像素点值为1，而值都为0的像素点值不变。即：

(1)ifcount1≥0.9 and count1≤1.0 So mask1＝2；if count2≥0.9 and count2≤1.0 So mask2＝2

(2)ifmask1＝2 and mask2＝2 So mask＝2；if mask1＝2 and mask2＝1 So mask＝1；if mask1＝1 andmask2＝2 So mask＝1；if mask 1＝0 and mask2＝0 So mask＝0

count1：南北矩阵等于1的总个数除以101的值；count2：东西矩阵等于1的总个数除以101的值；

mask1：南北矩阵检测的结果；mask2：东西矩阵检测的结果；mask：最终检测的水体掩膜图。

B4、以像素值为2的像素点为圆心，以100为半径，其圆形内最北边，最南边，最东边，最西边的边界点按下式计算：

$\mathrm{Rate}=\frac{E\_\mathrm{col}-W\_\mathrm{col}}{S\_\mathrm{row}-N\_\mathrm{row}}\mathrm{or}\frac{S\_\mathrm{row}-N\_\mathrm{row}}{E\_\mathrm{col}-W\_\mathrm{col}}$

E_col：最东边像素点的列数；W_col：最西边像素点的列数；N_row：最北边像素点的行数；S_row：最南边像素点的行数；

如果Rate大于高宽比的阈值，则算法继续，否则，跳出算法，重设阈值，Rate为圆内高宽比值；

B5、以像素值为2的像素点为圆心，以100为半径，其圆形内所有值清零，即完成大型湖泊的去除；

C、大型湖泊去除以后，重新对水系网进行全局搜索，自动运行智能蚁群算法，优化水系网中的小断流，如图4所示，其具体步骤如下：

C1、首先，输入去除大型湖泊的水体掩膜图以及相对应的归一化植被指数NDVI图像，同时建立一个3*3搜索矩阵，并搜索独立在水系网外的各个孤立点，并标记它们，若3*3搜索矩阵中心像素为1，且周围8个点为0，标定其中心点的位置，并赋值为2；

C2、接着，以所述孤立点为智能蚁群算法的起点，并将其作为3*3搜索矩阵的中心点，利用智能蚁群算法的启发性寻找下一个最有可能为水系网的像素点，并将其标定成1，并将3*3搜索矩阵中至少存在3个标定为1的值作为终止条件，其具体步骤如下：

C21、以标定为2的像素点为3*3搜索矩阵中心，建立一个3*3搜索矩阵，同时建立一个与之相对应的信息素矩阵并将其初始化。另外，利用中分辨率成像光谱仪Moderate-resolution Imaging Spectroradiometer地面反射率产品MODIS09的第一波段BAND1和第二波段BAND2数据，计算得出与之对应的归一化植被指数NDVI图像，并输入算法；

C22、启动智能蚁群算法，以3*3搜索矩阵为单位，每个起点放置10只蚂蚁，按下面公式完成计算：

左式为一个像素点的蚂蚁转移概率的计算公式，

其中：α，β为随机参数，作用为调节τi，ηi对蚂蚁的影响度；

τ_i为第i个像素点的归一化植被指数NDVI的倒数；

η_i为第i个像素点上蚂蚁的信息素浓度，

由于信息素随着时间的流逝，其浓度会减弱，而对于蚂蚁爬过的地方，由于蚂蚁本身会释放信息素，所以其爬过的地方浓度也会相应的增加，从而出现了正负反馈现象，所以可有效避免蚂蚁在寻找水体像素点的过程中陷入局部最优化，信息素更新公式如下：

P_i＝P_i-1+Δτ-0.1*T

其中，P_i，P_i-1为蚂蚁爬过前后的蚂蚁转移概率，P_i为当前状态，P_i-1为前一状态；

Δτ为蚂蚁爬过后信息素的增值，且由前一状态的信息素决定；

T为前一状态转换成当前状态所消耗的时间；

C23、完成上一步以后，计算3*3搜索矩阵中转移概率最大的值，并得到它在3*3搜索矩阵中具体位置，蚂蚁总是向着信息素最大的地方爬行，它最可能为遥感图像遗漏的水体像素点，并标定它的值为1，若发现最大值为起点位置，则随机赋值给其他8个点中的任意一个点为1；再者，当10只蚂蚁都爬行完毕，更新信息素矩阵，同时将新发现的水体像素点作为新的起点，继续寻找，直到满足终止条件为止；

C24、进一步在算法外围设置时间步阈值(COUNT)，其中规定每个时间步内蚂蚁只能移动一个像素长度，从而可控制蚂蚁的转移范围；

C3、完成步骤C2后，输出步骤C2的最终优化完小断流的水系网图像；

D、然后，建立一个3*3搜索矩阵，采用人工交互方式，并利用智能蚁群算法处理大断流，优化大断流，参照图4，利用人工交互方式标定未被智能蚁群算法所处理大断流，在小范围内运行智能蚁群算法，并在蚂蚁信息素指引下，优化大断流。即当完成智能蚁群算法自动优化小断流以后，算法状态转为人工交互方式(即半自动方式)。因为在中低分辨率遥感图像中云层等因素影响而产生比较大的断流会被自动运行的智能蚁群算法所忽略，所以加入部分人工判读可以使算法计算复杂度降低并使水系网的提取更加准确，如图4所示，其具体步骤如下：

D1、首先，输入去除大型湖泊的水体掩膜图以及相对应的归一化植被指数NDVI图像(因计算机程序源码相同，图4中C1与D1为共用，说明书为叙述方便，采用分别介绍方式)，同时建立一个3*3搜索矩阵，通过人工判定，设置小范围内智能蚁群算法的起点和终点，其中起点和终点可并称为“种子点”或“标志点”；

D2、随后，启动智能蚁群算法，在小范围内进行大断流的优化，输入归一化植被指数NDVI遥感图像，且同时初始化智能蚁群算法中的信息素矩阵，同样以3*3搜索矩阵为一单位，其中在每个起点上放置10只蚂蚁，按步骤C22中的蚂蚁转移概率公式和信息素更新公式完成计算；

D3、在人工交互方式中，智能蚁群算法具有算法的终止点，且与全自动智能蚁群算法有所不同，为防止智能蚁群算法陷入局部最优解而失去寻找终止点的能力，设置一个信息素阈值K，即在智能蚁群算法运行中，如果在信息素矩阵中有元素大于K，则该像素的信息素归零，从而避免陷入局部最优解。在完成步骤D2后，计算3*3搜索矩阵中转移概率最大的值，并得到它在3*3搜索矩阵中具体位置，标定它的值为1，若发现最大值为起点位置，则随机赋值给其他8个点中的任意一个点为1，进而更新信息素矩阵。同时将新发现的水体点作为新的起点，继续寻找，并设置一个终点判定值，即若此蚂蚁在规定的时间步阈值(COUNT)内无法到达终点，则将此只蚂蚁爬过的地方重新归零，并在起点重新启动一只新的蚂蚁，转步骤D2；若蚂蚁在规定的时间步阈值内到达终点，则保留此蚂蚁爬过地方的所有值，同时终止起点再次发送蚂蚁，记录蚂蚁爬过地方的所有值，从而完成水系网大断流的优化；

D4、最后，输出步骤D3的最终优化完大断流的水系网图像(因计算机程序源码相同，图4中C3与D4为共用，说明书为叙述方便，采用分别介绍方式)；

E、最后，利用遥感图像处理系统软件ENVI中的矢量转换算法，完成水系网矢量图的绘制。

利用矢量转换算法，完成水系网的绘制：利用遥感图像处理系统软件ENVI中矢量转换算法完成计算，其中包括细化河流、栅格图像转换成矢量图等步骤，最终输出优化后的水系网矢量图。

Claims (3)

1.基于智能蚁群算法的遥感图像水系网半自动的提取方法，其特征在于它包括如下步骤：

A、利用多光谱遥感数据和归一化植被指数NDVI数据初次提取水系网，得到水体掩膜图的初步结果；

B、分别建立5*1矩阵、101*1南北矩阵、1*101东西矩阵，判定大型湖泊并同时去除大型湖泊；

C、大型湖泊去除以后，重新对水系网进行全局搜索，自动运行智能蚁群算法，优化小断流，其具体步骤为：

C1、首先，输入去除大型湖泊的水体掩膜图以及相对应的归一化植被指数NDVI图像，同时建立一个3*3搜索矩阵，并搜索独立在水系网外的各个孤立点，并标记它们，若3*3搜索矩阵中心像素为1，且周围8个点为0，标定其中心点的位置，并赋值为2；

C2、接着，以所述孤立点为智能蚁群算法的起点，并将其作为3*3搜索矩阵的中心点，利用智能蚁群算法的启发性寻找下一个最有可能为水系网的像素点，并将其标定成1，并将3*3搜索矩阵中至少存在3个标定为1的值作为终止条件，其具体步骤如下：

C21、以标定为2的像素点为3*3搜索矩阵中心，建立一个3*3搜索矩阵，同时建立一个与之相对应的信息素矩阵并将其初始化；另外，利用中分辨率成像光谱仪Moderate-resolutionImaging Spectroradiometer地面反射率产品MODIS09的第一波段BAND1和第二波段BAND2数据，计算得出与之对应的归一化植被指数NDVI图像，并输入算法；

C22、启动智能蚁群算法，以3*3搜索矩阵为单位，每个起点放置10只蚂蚁，按下面公式完成计算：

左式为一个像素点的蚂蚁转移概率的计算公式，

其中，α，β为随机参数，作用为调节τ_i，η_i对蚂蚁的影响度；

τ_i为第i个像素点的归一化植被指数NDVI的倒数；

η_i为第i个像素点上蚂蚁的信息素浓度，

信息素具有正负反馈现象，可有效避免蚂蚁在寻找水体像素点的过程中陷入局部最优化，信息素更新公式如下：

P_i＝P_i-1+Δτ-0.1*T

其中，P_i，P_i-1为蚂蚁爬过前后的蚂蚁转移概率，P_i为当前状态，P_i-1为前一状态；

Δτ为蚂蚁爬过后信息素的增值，且由前一状态的信息素决定；

T为前一状态转换成当前状态所消耗的时间；

C23、完成上一步骤C22后，计算3*3搜索矩阵中转移概率最大的值，并得到它在3*3搜索矩阵中具体位置，蚂蚁总是向着信息素最大的地方爬行，它最可能为遥感图像遗漏的水体像素点，并标定它的值为1，若发现最大值为起点位置，则随机赋值给其他8个点中的任意一个点为1；再者，当10只蚂蚁都爬行完毕，更新信息素矩阵，同时将新发现的水体像素点作为新的起点，继续寻找，直到满足终止条件为止；

C24、进一步在算法外围设置时间步阈值，其中规定每个时间步内蚂蚁只能移动一个像素长度，从而可控制蚂蚁的转移范围；

C3、完成步骤C2后，输出步骤C2的最终优化完小断流的水系网图像；

D、然后，建立一个3*3搜索矩阵，采用人工交互方式，并利用智能蚁群算法优化大断流，续接大断流，其具体步骤如下：

D1、首先，输入去除大型湖泊的水体掩膜图以及相对应的归一化植被指数NDVI图像，同时建立一个3*3搜索矩阵，通过人工判定，设置小范围内智能蚁群算法的起点和终点，其中起点和终点可并称为“种子点”或“标志点”；

D2、随后，启动智能蚁群算法，在小范围内进行大断流的续接，输入归一化植被指数NDVI遥感图像，且同时初始化智能蚁群算法中的信息素矩阵，同样以3*3搜索矩阵为一单位，其中在每个起点上放置10只蚂蚁，按步骤C22中的蚂蚁转移概率公式和信息素更新公式完成计算；

D3、为防止智能蚁群算法陷入局部最优解，设置一个信息素阈值K，在智能蚁群算法运行中，如果在信息素矩阵中有元素大于K，则该像素的信息素归零，从而避免陷入局部最优解，在完成步骤D2后，计算3*3搜索矩阵中转移概率最大的值，并得到它在3*3搜索矩阵中具体位置，标定它的值为1，若发现最大值为起点位置，则随机赋值给其他8个点中的任意一个点为1，进而更新信息素矩阵，同时将新发现的水体点作为新的起点，继续寻找，并设置一个终点判定值，即若此蚂蚁在规定的时间步阈值内无法到达终点，则将此只蚂蚁爬过的地方重新归零，并在起点重新启动一只新的蚂蚁，转步骤D2；若蚂蚁在规定的时间步阈值内到达终点，则保留此蚂蚁爬过地方的所有值，同时终止起点再次发送蚂蚁，记录蚂蚁爬过地方的所有值，从而完成水系网大断流的续接；

D4、输出步骤D3的最终优化完大断流的水系网图像；

E、将上述输出优化大断流的水系网图像，利用遥感图像处理系统软件ENVI中的矢量转换算法，完成水系网矢量图的绘制。

2.根据权利要求1所述的基于智能蚁群算法的遥感图像水系网半自动的提取方法，其特征在于所述的步骤A中利用多光谱遥感数据和归一化植被指数NDVI数据初次提取水系网，得到水体掩膜图的初步结果的具体步骤如下：

首先利用近红外波段和红光波段数据计算归一化植被指数NDVI，再对其直方图进行分析，确定水体判断阈值K₁；其次，根据水体的物理特性，利用中分辨率成像光谱仪Moderate-resolution Imaging Spectroradiometer地面反射率产品MODIS09的第三波段BAND3和第二波段BAND2波段数据，计算其直方图，并分别取得BAND3直方图的陡降拐点坐标K₂及BAND2直方图的陡升拐点坐标K₃，确定水体提取的阈值，进而利用布尔公式，RIVER_MASK＝(NDVI≤K₁)AND(BAND3≥K₂)AND(BAND2≤K₃)，得到水体掩膜图的初步结果。

3.根据权利要求2所述的基于智能蚁群算法的遥感图像水系网半自动的提取方法，其特征在于所述的步骤B中分别建立5*1矩阵、101*1南北矩阵、1*101东西矩阵，判定大型湖泊并同时完成去除大型湖泊的具体步骤如下：

B1、首先，输入初步提取水体掩膜图，同时建立一个小型矩阵，即5*1矩阵，然后对全幅水体掩膜图进行全局搜索边界点；

B2、建立一个101*1南北矩阵和一个1*101东西矩阵，分别承担从南北方向和东西方向搜索大型湖泊的任务；

B3、遍历整幅图像，分别设立南北、东西阈值，利用101*1南北矩阵和1*101东西矩阵，分别在两个方向上检测所有像元等于1的个数，并记录数据，用其数据除以101，若这个结果值大于0.9，且小于1.0，则分别记录101*1南北矩阵和1*101东西矩阵中心点的位置，并将其像素值赋值为2；再次检验由101*1南北矩阵和1*101东西矩阵赋值为2的像素位置，选择位置公共的像素点，并重新赋值为2，其他值按布尔公式得出，即若只有一个为2的值，则恢复像素点值为1，而值都为0的像素点的值不变；

B4、以像素值为2的像素点为圆心，以100为半径，其圆形内最北边，最南边，最东边和最西边的边界点按下式计算： $\mathrm{Rate}=\frac{E\_\mathrm{col}-W\_\mathrm{col}}{S\_\mathrm{row}-N\_\mathrm{row}}\mathrm{or}\frac{S\_\mathrm{row}-N\_\mathrm{row}}{E\_\mathrm{col}-W\_\mathrm{col}}$

E_col：最东边像素点的列数；W_col：最西边像素点的列数；N_row：最北边像素点的行数；S_row：最南边像素点的行数，

如果Rate大于高宽比的阈值，则算法继续，否则，跳出算法，重设阈值，Rate为圆内高宽比值；

B5、最后以像素值为2的像素点为圆心，以100为半径，其圆形内所有值清零，即完成大型湖泊的去除。

Images