林草湿数据与三调对接融合图斑边界重合方法的研究

梁志强，丛喜东

(1.黑龙江省自然资源权益调查监测院，黑龙江 哈尔滨 150080；2.黑龙江省生态研究所，黑龙江 哈尔滨 150080)

1)2021年5月国家林业和草原局发布《林草湿数据与第三次全国国土调查数据对接融合技术指南》(以下简称“指南”)中指出，充分发挥第三次全国国土调查(以下简称“三调”)成果在国土空间管理中的“统一底版”作用，依据“三调”成果，理清林地、草地、湿地的现状范围界线，解决地类交叉重叠问题，融合林地、草地、湿地等资源信息，优化国家级公益林范围，将国家级公益林落到山头地块，形成与“三调”无缝衔接的林草资源图(含湿地)，支撑森林生态效益补偿制度，为林草生态网络感知系统建设，开展林草资源及生态状况监测，加强林草资源科学保护、系统修复和合理利用，推进林草治理体系和治理能力现代化提供支撑[1]。

2)指南中明确了目的、任务、原则要求、技术流程、工作成果、质量管理等内容，并对落界经营界线和小班细化做了具体说明。指南中指出经营界线存在误差性质偏移，无实质改变的情况下，应保持与三调界线重合；小班进行分割及碎班合并，对于区划尺度不同或区划误差，引起小班界线相互间里出外进等情况，可根据影像特征和现地情况，对小班边界进行修正或合并小班，原则上应以较高区划精度结果为准。

3)在具体的对接融合过程中由于林草湿数据与三调数据在区划精度上存在一定的差异，常规GIS软件，对容差范围内的现象不视为错误，故存在边线或折点无法重叠情况，现将具体问题进行分析，并针对问题提出解决思路与方法，最终形成满足对接融合要求的修正数据。

1 对接融合方法概述

1)采用GIS软件，将林草湿图层分别以县为单位，进行关键字赋值，用于图斑融合后挂接属性数据[2]，为加快图形融合数据处理速度，可保留适量相应图层的影响图斑边界处理的重要属性因子，其他属性因子则在对接融合后进行挂接，如三调图斑保留标识码、地类、权属单位名称、坐落单位名称、图斑面积等，林地保留单位名称、林班、小班、地类等，湿地保留斑块号、斑块名称、湿地类型等，草地保留斑块编号、草地类型等；

2)将统一空间参考的三调数据图斑与林地、湿地、草地图斑分别进行联合操作，形成相应的对接融合图层；对不满足最小区划面积条件的细碎图斑进行合并，合并后图斑不可跨三调图斑，这部分可通过技术手段自动合并；对大于最小区划面积的图斑，需要采用人工方式根据实际情况进行合并、切割、保留等处理；图斑融合完毕后进行空间拓扑、质检、修正，挂接重要属性因子，进行数据统计分析，形成空间数据、专题图、统计分析表、分析报告等验收成果。

2 对接融合数据边线及折点不重叠分析

在空间数据对接融合过程中，由于第三次全国国土调查数据矢量图在GIS软件中XY容差为0.0001m，通过分析三调数据图斑折点坐标发现，所有折点坐标均进行了处理，采用了4舍6入的方式，保留5位小数。而以往林草湿数据图斑的生产过程中采用GIS软件默认XY容差为0.001m，折点坐标采用默认保留小数位，二者在容差及坐标精度方面存在一定差异。在图斑对接融合过程中，在容差范围内二者边界虽然是一致的，但在放大比例尺到一定范围时，二者将出现不完全重合的现象，即边线或折点存在位置偏移，且该偏移值在GIS软件的容差范围内，图斑进行编辑操作时无法实现国土图斑的折点捕捉，从而修改操作难度大，且修改量大。如图1对接融合图斑与三调数据图斑边界及折点不重合现象示例所示。

图1 对接融合图斑与三调图斑边界及折点不重合示例Fig.1 Example of the Mismatch between the Boundaries and Vertices of the Butt Fusion Patch and Tritone Patch

3 解决问题的思路和方法

针对上述对接融合数据所产生的情况，可以通过遍历对接融合数据图斑折点获取距离三调图斑的最近折点或最近边线的垂点的方式实现对接融合数据图斑折点坐标的修改，使对接融合数据图斑边线及折点与国土数据图斑边线及折点保持完全一致。具体方法步骤如下：

1)获取对接融合数据图斑边线的所有折点。将对接融合数据增加一个唯一值字段并赋值，用于图斑修正后数据验证及属性数据挂接。

2)获取该对接融合数据图斑对应的国土图斑边线折点。对接融合数据图斑的标识码(BSM)与三调数据图斑标识码应一致，可通过BSM筛选三调数据图斑。如标识码不一致会导致修正折点坐标值计算有误。

3)获取对接融合数据图斑折点与对应三调数据图斑距离最近的折点。逐点遍历对接融合数据图斑折点，查找与该折点距离最近的对应的三调数据图斑折点，且该距离小于0.0005m，即视为该对接融合数据图斑折点应与三调数据图斑折点为同一点，对接融合数据图斑折点的坐标用该查找到的三调数据图斑折点坐标进行替换。

4)获取距离对接融合数据图斑折点最近的国土数据图斑边线上的线段，取点到线段的垂点作为折点。当未找到满足具体误差范围(<0.0005m)的最近折点，说明该对接融合数据图斑折点误差限范围内无三调图斑折点，此时需查找与该折点最近的三调图斑边线上的线段，且该点与线段两个端点的距离与线段长度的误差也在0.0005m范围内，即该点位于线段的两个端点之间。获取该对接融合数据图斑折点到三调图斑线段的垂足点，且垂足点与折点距离在误差限范围内，将垂点坐标替换对接融合数据图斑的折点坐标，也可以获取对接融合数据折点处线段与三调图斑最近线段的交点作为替换坐标值，通过分析垂点与交点误差很小，故取垂点坐标即可。

5)创建新的矢量图层，将对接融合数据的每个图斑折点按照重新计算的折点坐标值重新生成图斑要素。

4 软件工具实现及成果验证

通过图斑的折点修正方法实现对接融合数据与三调图斑边界的重合操作，图斑的折点遍历、点间距、点线距、垂点等计算采用Visual Studio开发环境、C#语言、GDAL开源技术进行实现。

1)GDAL开源技术。GDAL(Geospatial Data Abstraction Library，地理空间数据抽象库)是一个在X/MIT许可协议下的开源空间数据(栅格数据和矢量数据)转换库，包括读取、写入、转换、处理各种栅格和矢量数据格式[3]，可通过GDAL读取矢量图斑要素属性及几何图形坐标节点等空间及属性信息。在Visual Studio开发环境中，可加载GDAL的DLL类库文件实现对矢量图层及图斑的具体操作。

2)通过C#编码的方式实现软件功能的开发。软件工具的主要功能是加载对接融合和三调成果矢量数据，自动读取对接融合数据图斑的边线折点及对应三调图斑的边线或折点，通过计算获取对接融合数据图斑的每个折点的修正值，并自动创建新的图斑要素，要素包含对接融合数据图斑的关键字及三调图斑标识码字段。对接融合数据修正后结果可通过关键字及BSM进行相应林草湿数据与三调数据的属性连接。

软件实现的部分关键代码如下：

//计算点与点距离

private double GetMinDistanceP2P(double p1x, double p1y, double p2x, double p2y)

{

return Math.Pow((Math.Pow((p1x - p2x), 2) + Math.Pow((p1y - p2y), 2)), 0.5);

}

//计算点到线段的距离

private double GetMinDistanceP2L(double pX, double pY, double segAX, double segAY, double segBX, double segBY)

{

double product_ap_ab = (pX - segAX) * (segBX - segAX) + (pY - segAY) * (segBY - segAY);

double product_ab_ab = (segBX - segAX) * (segBX - segAX) + (segBY - segAY) * (segBY - segAY);

double r = product_ap_ab / product_ab_ab;

if (r >=1)

{

return Math.Sqrt((pX-segBX)*(pX-segBX) + (pY - segBY) * (pY - segBY));

}

else if (r <=0)

{

return Math.Sqrt((pX-segAX)*(pX-segAX) + (pY - segAY) * (pY - segAY));

}

else

{

double vector_ac_X = r * (segBX - segAX);

double vector_ac_Y = r * (segBY - segAY);

double vector_cp_X = pX - segAX - vector_ac_X;

double vector_cp_Y = pY - segAY - vector_ac_Y;

return Math.Sqrt(vector_cp_X * vector_cp_X + vector_cp_Y * vector_cp_Y);

}

}

//获取点到线段的垂足坐标，x1、y1、x2、y2为线段的端点坐标

private double[] GetVerticalAixP2L (double x0, double y0, double x1, double y1, double x2, double y2)

{

double[] ret = new double[] { 0, 0 };

if (x2 - x1 == 0)

{

ret[0] = x2;

ret[1] = y0;

}

else

{

//计算斜率

double k1 = (y2 - y1) / (x2 - x1);

// 垂足x坐标

double x = (Math.Pow(k1, 2) * x1 + (k1 * (y0 - y1)) + (x0)) / (Math.Pow(k1, 2) + 1);

ouble y = (k1) * (x - (x1)) + (y1);

ret[0] = x;

ret[1] = y;

}

return ret;

}

(3)通过修正前后数据比对分析。如图3对接融合图斑边界及折点修正前后对比图所示，红色边线为对接融合数据的修正前图斑边线，紫色边线为修正后边线，黑色边线为三调图斑边线，通过对比图可以看出修正后边线与三调数据图斑边线及折点完全一致，并通过对某乡镇16000多个图斑进行验证，修正前后对接融合数据的面积、周长等均与对接融合数据图斑未修正前保持基本一致，其修正后结果满足数据生产和成果质量要求。

图3 对接融合图斑边界及折点修正前后对比Fig.3 Contrast before and after the Patch Boundary and Vertex Correction of the Butt Fusion Pattern