基于多元非线性回归的CORS 实时平面坐标转换方法

张永峰，李 琪，程晓晖，邹 璇，崔健慧

（1.武汉攀达时空科技有限公司，武汉 430079；2.广州市城市规划勘测设计研究院，广州 510060；3.武汉大学 卫星导航定位技术研究中心，武汉 430079； 4.山东科技大学 测绘科学与工程学院，山东 青岛 266510）

0 引言

卫星定位连续运行参考站（continuously operating reference stations, CORS）是1 种将全球卫星导航系统（global navigation satellite system, GNSS）的导航定位技术、现代通讯技术、计算机技术、测绘技术等多种技术集成的实用性网络系统[1]，系统全自动无需人工干预，代替了分散的、短期的、静态的定位测量，极大地提高了测绘工作的效率[2]。但是CORS 测量所得到的测绘成果一般为1984 世界大地坐标系（world geodetic coordinate system 1984,WGS84）或2000 国家大地坐标系（China geodetic coordinate system 2000, CGCS2000）坐标系下的大地坐标，而测绘生产工作中使用的坐标系通常为国家坐标系统或城市独立坐标系统，因此就出现了CORS 的坐标转换问题[3]。

根据国家相关法律法规的规定，坐标转换参数属于机密级以上保密范围测绘成果，不能直接交给用户使用。因此在目前的CORS 系统中，通常采用以下3 种坐标转换方式[4]进行坐标转换：

1）提供人工或在线坐标转换服务，用户自行将WGS84 或CGCS2000 测量成果发送至相关部门或网络服务器，待转换完成后再从相关部门或服务器上接收坐标转换结果。该方法保密性强，但增加了用户及相关部门的工作量，转换效率及实时性低。

2）在用户终端上进行参数加密，在终端上安装坐标转换软件并进行加密设置。该方法实时性较强，但由于转换工作均在用户端完成，因此保密性一般，而且用户终端种类繁杂，难以进行操作。

3）通过参数分解向用户提供伪参数，将4 参数或7 参数分解为2 部分，在CORS 数据播发服务器与用户之间架设中间服务器，使用第1 部分参数对虚拟参考站（virtual reference station, VRS）进行坐标转换，再通过用户终端的第2 部分参数进行2 次转换后，最终得到真实坐标。该方法保密性好，实时性强，广泛应用于7 参数坐标转换。

第3 种方式是目前解决CORS 实时坐标转换常用的方式，可以满足7 参数坐标转换参数的需求，在广东CORS（Guangdong CORS, GDCORS）等系统中已得到实际应用；但在其他方式的坐标转换模型中不一定适用，如广州等大型城市地方坐标系。由于目前国内城市的地域范围逐渐扩大；导致东西跨度太大，因此在使用地方坐标系时，为了减少投影变形带来的误差，往往不会按标准3°带或6°带选取中央子午线，而是将中央子午线取至广州市中央，而在某些高原地区也会将投影面高程取至测区高程均值[5-6]。这些修改后的投影参数同样作为坐标转换参数的一部分。当4 参数与投影参数同时作为保密参数时，由于存在2 种转换方式，无法像独立的7 参数一样分解为2 套参数进行处理。本文针对CORS 坐标转换4 参数与投影参数模型，介绍了1 种基于多元回归方程的CORS 实时平面坐标转换服务方法。

1 基于多元回归方程的实时平面坐标转换

首先，对4 参数与投影参数进行修改生成伪参数。然后，对于使用伪参数产生的坐标偏差，本文提出利用多元回归方程计算补偿值，在VRS 虚拟参考站上进行坐标补偿。最后，将改正后的坐标编码成国际海运事业无线电技术委员会（Radio Technical Commission for Maritime Services,RTCM）电文播发给流动站用户，从而对流动站坐标进行偏差改正，使用户得到准确的地方坐标。本文提出的方法流程如图1 所示。

图1 方法流程

1.1 伪参数生成

在计算多元回归方程之前，首先需要进行伪参数的生成。由于VRS 的坐标偏差会导致基线解算出现尺度和方向上的偏差，因此在实际工作中，需要尽量提高基准站坐标与观测值间的一致性，以避免这种情况发生[7]。因此，出于对基线解算精度与保密的双重要求，可以在原始4 参数的平移参数上各加上数米的偏移量，在实际操作时，可以加入偏移量使平移参数变为整数，方便用户使用；将投影参数的中央子午线偏移数10″，投影高程可修改为0 m。

1.2 多元回归方程的建立

目前很多采用格网法来对坐标偏差值进行改正[8-11]；但对本文伪参数坐标偏差值进行补偿时，由于补偿值间变化趋势较大，为了满足CORS 的厘米级精度要求，需要将格网划分得很密，必然会导致格网数据量变得很大，严重影响了CORS 的处理速度，降低了实时性。

为了降低CORS 数据处理中心的数据处理量，提高CORS 处理速度，本文针对坐标补偿值建立了多元回归方程，通过方程及用户位置计算坐标偏差改正值。

对于方程样本点，使用伪参数进行坐标转换带来的坐标偏差值为

式中：( Δx , Δy) 为坐标偏差值；( x,y) 为真实地方坐标，( x ′,y)′为使用伪参数得到的地方坐标。

根据测区实际情况选择合适的回归方程形式，本文对比了1 阶、2 阶和3 阶多项式的效果为：

式中：( Δxi、 Δyi)为坐标偏差值；( B0、 L0)为测区中央方程原点的经纬度；( Bi、 Li)为各样本点的经纬度；a、 b、 c、 d、 e、 f 、 g为方程参数。

根据最小二乘估计方法求解回归参数[12]，得出多元回归方程及残差值，计算出中误差为

式中：V 为残差值阵；n 为样本总数；m 为方程参数个数。

2 多元回归方程建立实例

本文以广州地区的实测数据为例建立多元回归方程。在选取样本点时，在全广州市范围内以经度50″、纬度40″为间距进行采样，共计得到4 267 个样本点，采用相同的伪坐标转换参数，以这些坐标点的经度、纬度、北方向差值、东方向差值作为样本，分别使用式（2）～式（4）来建立多元回归方程。由于涉及基础测绘坐标数据和转换参数的保密性，因此本文不给出具体的计算数据及多元回归方程，直接给出最后的精度评估结果。

样本点的精度统计结果如表1 所示，其中：2 阶、3 阶结果的北方向、东方向精度均高于1 阶，尤其在东方向上，提升效果明显；3 阶与2 阶相比，残差值相近，同时在3 阶项上的参数值均接近于0，而且采用2 阶多项式可以减少运算量。因此本文将采用式（3）得到的2 阶多元回归方程结果进行CORS 坐标转换。

表1 样本点残差值统计结果

2 阶多元回归方程得到的样本点残差值统计结果中：东方向中误差为2 mm，最大偏差为7 mm，考虑到CORS 平面精度为3 cm，坐标转换带来的东方向误差可以忽略不计；北方向中误差接近1 cm，最大偏差4.5 cm，已经超过CORS 平面精度要求。为进一步探讨北方向误差的影响区域，得到北方向误差空间分布情况如图2 所示。

图2 样本点残差值分布

由精度统计结果和误差分布可知：北方向多元回归精度在广州市中央大部分区域优于1 cm，带来的误差可以忽略不计；在东南角和西北角的误差值会大于2 cm，但考虑到东南角为海域，西北角为山区，基本不存在测绘作业任务，因此北方向多元回归误差也基本不会对广州CORS 服务性能带来影响。而西南区域测量作业量较大，因此必须对这片区域的残差值进行修正。

由于该区域残差值均在2 ～3.5 cm 之间，且往同一方向单向递增，因此可以将其作为系统差进行消除，根据CORS 服务精度要求，将该部分残差值统一改正2 cm。在进行修正时，可根据用户上传位置经纬度进行判断，若用户处于该区域范围内，则在多元回归方程的基础上再进行1 次改正。经过以上方式处理后，得到新的样本点残差值分布情况如图3 所示。

图3 改正后样本点残差值分布

3 坐标转换服务精度分析

为了验证本文提出的CORS 实时平面坐标转换服务精度，在广州市范围内选取了13 个广州市二等GNSS 控制点和57 个广州市二等平面控制网点作为测试点，测试过程中使用了天宝、攀达时空等多种品牌的接收机，这些流动站终端设备能够接收北斗卫星导航系统（BeiDou navigation satellite system, BDS）及其他GNSS 的卫星信号。

测试时，首先在流动站终端上设置用户端参数，然后连接广州CORS 进行实时动态差分（realtime kinematic, RTK）测量，在采集数据时，以控制点方式采集20 个历元进行平滑处理得到测试结果。然后将结果与控制点广州市地方坐标系成果进行对比。控制点精度统计结果如表2 及图4 所示，北方向中误差为1.4 cm，东方向中误差为1.3 cm，最大值均小于3 cm，满足CORS 服务的平面精度要求[13-14]。

表2 控制点统计结果 单位：m

图4 控制点误差散点图

4 结束语

在CORS 系统日益发展的今天，只有将工程应用与CORS 更加紧密地结合起来，才能更好地利用CORS 提高测绘工作的效率，而坐标系统是所有工程应用都需要考虑到的问题。本文提出的基于多元回归方程的CORS 实时平面坐标转换服务系统，能够在保证测绘成果资料保密性的前提下，让用户直接通过CORS 获取地方平面坐标。经过精度分析表明，通过该方式得到的地方平面坐标精度优于3 cm，能够满足CORS 服务在工程应用中的精度要求。