地类图斑与多源在线地图坐标粗匹配研究

黄亮平，李德平,2*，周 亮,2，刘师师

（1.湖南师范大学 资源与环境科学学院，湖南 长沙 410081；2.湖南师范大学地理空间大数据挖掘与应用湖南省重点实验室，湖南 长沙 410081）

地类图斑与多源在线地图集成不仅可以促进地类图斑和多源在线地图资源的社会化应用，还可为政府相关部门和社会大众提供一种获取地类图斑与多源在线地图集成信息的手段。本文采用的地类图斑数据来源于某地区第三次土地调查数据库，多源在线地图以天地图、百度地图、Google Earth为例，分别将地类图斑与这3种在线地图进行集成。在集成过程中，由于地类图斑的坐标系为CGCS2000高斯投影坐标系，而天地图为CGCS2000坐标系、百度地图为百度坐标系（BD09）、Google Earth为1984世界大地坐标系（WGS84），因此需进行坐标匹配。本文研究的是坐标粗匹配，误差控制在10 m以内即可。

目前针对CGCS2000高斯投影坐标与CGCS2000坐标、BD09坐标、WGS84坐标的匹配均有研究。对于投影坐标与地理坐标的转换问题，无论是CGCS2000、WGS84、西安1980还是北京1954，均可利用高斯反算公式进行坐标转换，得到的坐标精度不低于0.0001″，虽然存在一定的误差，但能满足日常生产需求[1-5]。CGCS2000坐标与BD09坐标之间的转换，可利用多项式变换模型进行分区坐标转换，区域大小基本一致的坐标转换精度也基本一致[6]。对于百度地图坐标（BD09II）与WGS84坐标存在较大偏移，且BD09II坐标无法转换为WGS84坐标的现象，利用坐标转换算法可实现WGS84、GCJ02和BD09三种坐标之间的相互转换[7-8]；利用百度地图API或依托控制点库进行离线转换可实现WGS84坐标到BD09坐标的转换[9-10]。BD09坐标向WGS84坐标的转换方式主要包括二阶差分公式法、等量偏移法、格网法以及BP神经网络法[7,11]。对于CGCS2000坐标与WGS84坐标之间的转换，若需精确转换，目前推荐从政府部门和CORS服务中实现[12]。CGCS2000坐标与WGS84坐标的转换需明确所对应的历元和框架，可根据具体精度要求，分别利用历元转换、框架转换、布尔莎七参数模型或结合使用进行坐标转换[12-14]。历元转换影响远大于框架转换，历元引起的坐标差别可达dm级甚至m级，而框架引起的坐标差别仅为cm级[12-13]，因此进行框架转换时，必须先统一历元[13]，当精度要求为cm级时，只需进行历元转换[13]。CGCS2000坐标系与WGS84坐标系因扁率的差别不会影响某点的经度，对纬度和高度的影响则与该点的纬度有关，且最大差值远小于1 mm，在目前坐标系的实现精度范围内是可以忽略的[12,15-16]，因此当精度要求为cm级及以上时，认为CGCS2000坐标系和WGS84坐标系是基本相容的[12-13,16]。

1 地类图斑与3种在线地图的坐标系

地类图斑为CGCS2000高斯投影坐标系，按3°分带投影，位于37带内，中央子午线为111°；天地图为CGCS2000坐标系，百度地图为BD09坐标系，Google Earth为WGS84坐标系。CGCS2000坐标系是我国目前最新的国家大地坐标系，是一种地心坐标系，即坐标系原点为地球质心。我国自2008年7月1日起，全面启用CGCS2000坐标系，国家测绘局授权组织实施[15]。CGCS2000高斯投影坐标系是基于CGCS2000坐标系利用高斯投影方法进行投影得到的投影坐标系。WGS84坐标系是美国国防部制图局建立的坐标系，是国际上采用的一种地心坐标系，WGS84坐标系与CGCS2000坐标系椭球参数仅在扁率上略有差别[16-17]。BD09坐标系是在WGS84坐标系的基础上经过两次加密而得到的，第一次是从WGS84坐标系加密为火星坐标系（CGJ02），第二次是从GCJ02坐标系加密为BD09坐标系，且不同区域加密算法不一样[6]。

2 地类图斑与3种在线地图坐标粗匹配方法

2.1 地类图斑与天地图坐标粗匹配方法

2.1.1 CGCS2000坐标系的椭球参数

CGCS2000坐标系的椭球参数主要包括长半轴、短半轴、扁率、第一偏心率平方和第二偏心率平方，如表1所示。

表1 CGCS2000坐标系的椭球参数

2.1.2 高斯反算

高斯平面坐标(x,y)向大地坐标(L,B)转换的过程称为高斯反算，计算公式为[1]：

2.2 地类图斑与百度地图坐标粗匹配方法

地类图斑采用CGCS2000参考椭球，百度地图采用WGS84参考椭球，且在WGS84的基础上经过了两次加密。虽然可利用严密的坐标转换公式将CGCS2000参考椭球转换为WGS84参考椭球，再利用百度地图API、坐标转换算法或依托控制点库进行离线转换，将WGS84坐标转换为BD09坐标；但该方式比较复杂、需要控制点、百度地图API，一次只能转换100个点且不可逆。由于百度地图在不同区域采用不同的加密算法，且多项式变换模型阶数越高、校正变形越复杂、非控制点转换误差越显著，因此不宜采用三阶以上的多项式变换模型进行坐标变换[18]。本文试图利用不同大小的格网，结合平均位移法（零阶多项式）、一阶多项式、二阶多项式、三阶多项式以及平面四参数法对两套坐标系之间的匹配关系进行对比研究，如图1所示。

图1 地类图斑与百度地图坐标粗匹配技术流程图

2.2.1 坐标匹配方法

本文分别采用平均位移法、一阶多项式、二阶多项式、三阶多项式和平面四参数法进行对比研究，从而挑选出最适合本文的方法。5种方法的适用情况、需求参数个数、需要的最少非相关公共点数如表2所示，可以看出，除平均位移法外，其他方法若公共点个数大于最少非相关公共点数，则可利用最小二乘法求解。

表2 不同坐标匹配方法的适用特征

设任意一点A在地类图斑中的坐标为(x,y)，利用坐标匹配方法转换后的坐标为(x′,y′)，在百度地图中对应的真实坐标为(X,Y)。

1）平均位移法的数学模型为：

式中，Δx、Δy分别为x和y方向的平均位移量。

2）一阶多项式的数学模型为：

3）二阶多项式的数学模型为：

4）三阶多项式的数学模型为：

5）平面四参数法的数学模型为：

式中，Δx和Δy分别为x和y方向的平移参数；m为尺度参数；α为旋转参数。

2.2.2 精度评定标准

通过比较检验点预测值与真实值之间的均方根误差（RMSE）得到坐标匹配方法的精度。

式中，n为检验点个数。

2.3 地类图斑与Google Earth坐标粗匹配方法

为了验证在本文精度条件下地类图斑坐标与Google Earth坐标基本相容，获取了大量地类图斑与Google Earth相同位置特征点的坐标，再分别计算X方向和Y方向的差值进行验证，如图2所示。

图2 地类图斑与Google Earth坐标粗匹配技术流程图

3 地类图斑与3种在线地图坐标粗匹配过程

3.1 地类图斑与天地图坐标粗匹配

地类图斑为CGCS2000高斯投影坐标系，天地图为CGCS2000坐标系，因此基于CGCS2000坐标系椭球参数，利用高斯反算方法，将地类图斑的高斯投影坐标系转换为CGCS2000坐标系，即可进行地类图斑与天地图的坐标粗匹配。

3.2 地类图斑与百度地图坐标粗匹配

为了比较不同格网划分和不同坐标匹配方法对坐标匹配精度的影响，以便挑选出最适合的格网大小和坐标匹配方法，本文以某地区第三次土地调查数据库中的地类图斑为实验区，面积约为1 549 km2（未分区），首先在ArcGIS中依次添加20×20 km（2×2分区）、10×10 km（4×4 分区）、5×5 km（8×8分区）、2.5×2.5 km（16×16分区）的格网；再按“分布均匀，边界、四角有点”的原则[19]，分别利用ArcGIS和百度地图拾取坐标系统获取地类图斑和百度地图相同位置特征点的坐标（16对公共点、10对检验点）；然后对于每一对特征点，设置地类图斑平面坐标为(X地,Y地)，百度地图地理坐标为(L百,B百)，利用ArcGIS将百度地图地理坐标进行3°带高斯37带投影，得到平面坐标(X百,Y百)，并将地类图斑平面坐标与百度地图平面坐标整理成一个表格文件（表3）；最后分别采用平均位移法、一阶多项式、二阶多项式、三阶多项式、平面四参数法对16对公共点坐标进行参数计算，得到转换参数，基于转换参数利用10对检验点坐标进行坐标匹配，并进行精度评定。

表3 地类图斑与百度地图相同位置特征点的平面坐标/m

3.3 地类图斑与Google Earth坐标粗匹配

地类图斑为CGCS2000高斯投影坐标系，Google Earth为WGS84坐标系，分别利用ArcGIS和Google Earth客户端获取地类图斑和Google Earth相同位置的425对特征点坐标。对于每一对特征点，设置地类图斑平面坐标为(X地,Y地)，Google Earth地理坐标为（L谷,B谷)，利用ArcGIS将Google Earth地理坐标进行3°带高斯37带投影，得到平面坐标(X谷,Y谷)，并将地类图斑平面坐标与Google Earth平面坐标整理成一个表格文件（表4）。基于此，分别计算X和Y方向差值，得到其散点图（图3、4）。

表4 地类图斑与Google Earth相同位置特征点平面坐标/m

4 匹配结果与分析

4.1 地类图斑与天地图坐标粗匹配结果

对于地类图斑与天地图坐标的粗匹配，地类图斑与天地图均采用CGCS2000参考椭球，只是地类图斑为投影坐标系，而天地图为大地坐标系。因此，只需利用高斯反算就能实现地类图斑与天地图的坐标粗匹配，且坐标精度不低于0.000 1″，符合本文坐标粗匹配精度要求。

4.2 地类图斑与百度地图坐标粗匹配结果

对于地类图斑与百度地图坐标的粗匹配，针对不同格网分区，本文分别利用5种方法对地类图斑与百度地图的坐标匹配关系进行了对比研究，得到不同分区不同方法坐标匹配的RMSE（表5）。

表5 不同分区不同方法坐标匹配的RMSE/m

由表5可知，无论利用哪种方法，区域越小，则RMSE越小；在所有分区中，利用二阶多项式得到的RMSE最小；当区域较大时（未分区、2×2分区），利用5种方法得到的RMSE基本上一致，考虑到方法的简便性，宜选用平均位移法；当区域较小时（4×4分区、8×8分区、16×16分区），二阶多项式优于其他方法，因此宜选用二阶多项式方法。本文研究地类图斑与百度地图坐标的粗匹配，精度要求控制在10 m以内，同时考虑到区域越小带来的计算量越大，因此宜选用8×8分区格网结合二阶多项式方法进行坐标匹配。

4.3 地类图斑与Google Earth坐标粗匹配结果

对于地类图斑与Google Earth坐标的粗匹配，由图3、4可知，地类图斑与Google Earth在X、Y方向的差值范围约为(-6,10)和(-10,10)，符合本文坐标粗匹配精度要求，因此可将CGCS2000坐标系与WGS84坐标系当成同一坐标系看待。同时，在精度要求为cm级及以上时，可认为CGCS2000坐标与WGS84坐标基本一致[13]。

图3 地类图斑与Google Earth在X方向的差值图

图4 地类图斑与Google Earth在Y方向的差值图

5 实例验证

以VS2013为开发平台，以C#为编程语言，本文首先利用FME API获取地类图斑的空间信息，将Web Browser控件嵌入HTML，并在HTML中分别利用各自地图API调用3种在线地图；再将地类图斑以覆盖物的形式分别叠加在3种在线地图上，实现地类图斑与3种在线地图的集成；然后利用CGCS2000参考椭球参数进行高斯反算实现地类图斑与天地图的坐标粗匹配；最后利用8×8分区格网结合二阶多项式方法将地类图斑转换为百度地图平面坐标，由于百度地图为地理坐标，因此还需利用WGS84参考椭球参数进行高斯反算才能实现地类图斑与百度地图的坐标粗匹配。由于本文的精度要求在m级，CGCS2000坐标系与WGS84坐标系可通用，因此只需利用CGCS2000参考椭球参数将地类图斑进行高斯反算即可与Google Earth相匹配。

地类图斑与3种在线地图能很好地进行坐标粗匹配（图5），规避了传统地理信息保密性、参数保密性以及坐标精匹配。由图5可知，本文选用的8×8分区格网结合二阶多项式方法可较好地解决地类图斑与百度地图集成时出现的坐标匹配问题，是一种更简单、快速、直观的适用于普通社会大众的坐标粗匹配途径。

图5 地类图斑与天地图、百度地图、Google Earth坐标粗匹配示例

6 结 语

针对地类图斑与多源在线地图集成时出现的坐标不匹配问题，以天地图、百度地图、Google Earth为例，对地类图斑与3种在线地图的坐标匹配进行了研究，并通过实例对坐标匹配结果进行了验证。

1）对于采用同一参考椭球的地类图斑与天地图，利用高斯反算方法即可实现二者的粗匹配，坐标精度不低于0.000 1″。

2）对于地类图斑与百度地图的粗匹配，通过对比不同格网大小和不同方法对地类图斑与百度地图坐标匹配精度的影响发现，格网越小，坐标匹配精度越高。在5种方法中，二阶多项式方法的效果最好。区域较大时，宜选用平均位移法；区域较小时，宜选用二阶多项式方法；但实际应用时，可根据具体情况选择，应同时考虑精度与计算量，本文中宜选用8×8分区格网结合二阶多项式方法实现地类图斑与百度地图的粗匹配。

3）对于地类图斑与Google Earth坐标的粗匹配，本文验证了在精度要求为m级时，CGCS2000坐标系和WGS84坐标系可视为同一坐标系。