具体实施方式

下面对本发明的具体实施方式进行描述，以便于本技术领域的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。

参考图1，图1示出了地形勘测中的地形勘测点布测优化方法；如1所示，该方法包括步骤S1至步骤S5。

在步骤S1中，获取待研究区域的多光谱影像和全色影像，并对多光谱影像和全色影像进行预处理，之后将预处理后的影像进行镶嵌拼接得到研究区影像。

实施时，本方案优选步骤S1进一步包括：

S11、获取待研究区域的多光谱影像和全色影像，并对多光谱影像依次进行辐射定标、大气校正和正射校正及对全色影像进行辐射定标和正射校正；

S12、对正射校正后的多光谱和全色影像进行图像融合，以得到更高精度的待研究区域影像，之后采用分辨率高于多光谱影像和全色影像分辨率的遥感影像对融合后的影像进行几何校正；采用多光谱影像和全色波段影像进行波段融合，进而提升遥感影像的分辨率。

S13、对几何校正后的影像进行镶嵌拼接，得到最终的研究区影像。

在步骤S2中，采用Sobel算子对研究区影像进行边缘检测，得到待研究区域内的地物轮廓线。

在本发明的一个实施例中，所述步骤S2进一步包括：

S21、将研究区影像转换为由R、G、B三波段组成的jp2格式的图片，该种格式的图片可以使转换后图片分辨率不变。

S22、将图片转化为灰度图后，进行归一化处理；

S23、调用sobel算子分别对归一化后的图像进行0°、45°、90°、135°四个方向上的卷积，分别得到图像中每个像素点在四个方向的梯度值：

其中，A为归一化后的像素点的像素值；

S23、根据像素点在四个方向的梯度值，计算每个像素点卷积后的梯度值：

其中，G_0°、G_45°、G_90°、G_135°分别为像素点在0°、45°、90°、135°方向上的梯度值；G_x和G_y分别为像素点在x和y方向上的梯度值；G为像素点卷积后的梯度值；

S24、判断像素点的梯度值是否大于预设梯度阈值；若是，将像素点标记为边缘点，否则将像素点标记为非边缘点；

本方案的预设梯度阈值优选为待研究区域所有像素点卷积后的梯度值的平均值。

为了进一步提高标记的准确性，本方案的预设梯度阈值还可以做如下优化：

S241、选取研究区域的遥感影像中地物类型大于预设数量的区域作为试验区域，并采用人工标记的方式描绘地物轮廓线作为验证样本；

S242、获取待研究区域中与试验区域对应区块提取的地物轮廓线，并采用提取的地物轮廓线构成的试验样本；

S243、分别计算验证样本和试验样本构成的区域的覆盖面积F1和F2，并判断覆盖面积F2与覆盖面积F1的比值是否大于等于面积阈值；

S244、当比值大于等于面积阈值时，输出预设梯度阈值；当比值小于面积阈值时，令初始预设梯度阈值＝2*初始预设梯度阈值，之后返回步骤S243。

S25、采用标记的所有边缘点作为待研究区域内的地物轮廓线，为了便于后续查看，可以把边缘点置为白点，非边缘点置为黑点。

在步骤S3中，获取由地物轮廓线围成的面矢量多边形和地物轮廓线构成的面矢量多边形，并结合分辨率为30m的DEM数据(30m的DEM数据的数据来源可以参考网址http：//www.gscloud.cn/)统计每个多边形的特征值。

在本发明的一个实施例中，所述步骤S3进一步包括：

S31、采用S13获得的遥感影像进行地理配准，并定义地理坐标系为GCS_WGS_1984，得到具有空间参考的边缘影像；

S32、将边缘影像中地物轮廓线的像素值赋值为1，将非边缘点的像素赋值为0，并将赋值的边缘影像由栅格转换为面矢量，得到由地物轮廓线围成的面矢量多边形和地物轮廓线构成的面矢量多边形；

S33、根据30m分辨率的DEM数据，统计每个多边形内部高程的特征值均值、方差、标准差、中位数、众数和最值：

X_min＝min(X₁，X₂...X_i)，X_max＝max(X₁，X₂...X_i)

其中，X_min为多边形内高程值的最小值；X_max为多边形内高程值的最大值；M_O为多边形内高程值的众数；为多边形内高程值的平均值；X_i为多边形内的第i个最小高程值，M_d为多边形内高程值的中位数，σ为多边形内高程的标准差，n为多边形内高程的数量；

当n为奇数时，当N为偶数时，

在步骤S4中，计算多边形之间特征值的相似度，并将相似度大于预设阈值时对应的多边形合并为同一分区块；

在本发明的一个实施例中，所述步骤S4进一步包括：

S41、将所有多边形构成一个集合，在集合中随机选取一个多边形；

S42、计算选取的多边形与集合中余下多边形间相应特征值的比值：

Q_Ak2＝σ_A/σ_k，Q_Ak3＝M_{O_A}/M_{O_k}

Q_Ak4＝X_{min_A}/X_{min_k}、Q_Ak5＝X_{max_A}/X_{max_k}、Q_Ak6＝M_{d_A}/M_{d_k}

其中，σ_A、M_{O_A}、X_{min_A}、X_{max_A}和M_{d_A}分别为选取的多边形A高程值的平均值、标准差、众数、最大值、最小值和中位数；σ_k、M_{O_k}、X_{min_k}、X_{max_k}和M_{d_k}分别为余下多边形中第k个多边形高程值的平均值、标准差、众数、最大值、最小值和中位数；Q_Ak1、Q_Ak2、Q_Ak3、Q_Ak4、Q_Ak5和Q_Ak6分别为多边形A与余下多边形中第k个多边形的平均值比值、标准差比值、众数比值、最大值比值、最小值比值；

之后将两个多边形的6个比值(Q_Ak1、Q_Ak2、Q_Ak3、Q_Ak4、Q_Ak5和Q_Ak6)作为一组数据。

S43、判断每组数据中6个比值是否存全部大于预设阈值，若是，则标记余下多边形中的当前多边形，否则，不标记；

S44、当每组数据均与预设阈值比较后，判断是否存在被标记的多边形，若是，进入步骤S45，否则进入步骤S46；

S45、将被标记的多边形与选取的多边形合并成一个多边形，采用未进行合并的多边形更新集合，并随机选取一个多边形进入步骤S42；

S46、采用集合中的余下多边形更新集合，并随机选取一个多边形进入步骤S42。

本方案通过统计特征值比较，对待研究区域的分区进行二次验证；由于同一分区范围内的地貌较为平缓，且其在地形高程上相差较小。如芦苇的起伏程度相似，划分至一个区域后，在进行芦苇地形勘测时，由于地貌相同，选取少量地形勘测即可，从而减轻整体间隔数据采集所需的经济成本。

在步骤S5中、采用碎部测量原则将每个分区块划分成若干图块，当所述图块中的特征值大于等于预设特征值时，表明地形相对比较复杂，比如洼地、山地等，此时根据数字地形测量规范进行地形勘测点布测，其中地形勘测点尽量布在图块的边界轮廓线上。

当图块中的特征值小于预设特征值时，表明图块的地形相对平坦，此时可以采用米字形进行地形勘测点布测，且米字形的边点位于图块的边界轮廓线上。

综上所述，本方案的地形勘测点的布测优化方法将地形相似的区域划分至同一区块，并通过对区块内地形勘测点的合理布局以减少地形勘测点的布局数量，本方案通过在减少地形勘测点的数量同时，还能保证后续通过这些点采集回的高程数据插值得到的DEM的精度，同时还能大幅降低采集高程数据的成本。