ITRF与CGCS2000坐标系统转换

1.引言

北京54坐标系和西安1980坐标系的建立极大的促进了新中国测绘的发展,然而随着空间大地测量技术的兴起,这两种经典的局部大地坐标系已不能满足中国经济建设和国防建设的要求,建立以地球质量中心为原点的地心大地坐标系势在必行.建在ITRF基础上的2000中国大地坐标系(Chian Geodetic Coordinate System 2000,CGCS200)在21世纪初正式开始实施.中国自主研发的卫星导航和授时系统也采用了CGCS2000坐标系.然而国际上统一采用的是ITRF的坐标系,在处理测量数据是难免会遇到坐标系统的转换问题.本文在介绍ITRF和CGCS2000的基础上,分析如何使用七参数法来实现ITRF到CGCS2000的坐标转换.

2.ITRF与CGCS2000简介

2.1 ITRF

由于地球的大小和形状不断的发生变化,故建立在地球上坐标系统也应当随之而变动.1988看国际地球自转和参考系统服务(International Earth Rotation and Reference Systems Service,IERS)成立,联合所有可用的监测手段和数据(包括GPS,VLBI,SLR和DORIS)计算出实时的地球参考系统的参数以及参数变化速度,从而实现精准的地球参考框架.ITRF为国际公认的应用最广泛、精度最高的地心坐标框架。ITRF是ITRS(国际地球参考系统)的具体实现.ITRF每过一段时间就会更新,目前最新的参考框架为ITRF2014.

2.2 CGCS2000

CGCS2000是ITRF1997都为参考框架,而以2000.0作为参考历元.从某种意义上说,CGCS2000就是ITRS的一个具体实现.故CGCS2000的笛卡尔坐标系原点和轴的定义为:

原点: 地球的质量中心
Z 轴: 指向IERS参考极方向
X 轴: 指向IERS参考子午面与通过原点且同Z轴正交的赤道面的交线
Y 思: 与Z轴和X轴共同构成右手坐标系

CGCS2000的参考椭球为一等位旋转椭球。等位椭球(或水准椭球)定义为其椭球面是一等位面的椭球。CGCS2000的参考椭球的几何中心与坐标系的原点重合,旋转轴与坐标系的z轴一致。参考椭球既是几何应用的参考面,又是地球
表面上及空间正常重力场的参考面。
CGCS2000参考椭球四个常数为:

长半轴       a=6378137.0m
扁率         f=1/298.257222101
地球引力常数 GM=3986004.418×E8 m3s-2
地球角速度   w=7292115.0E-11 rad/S

3.转换原理

因为CGCS2000就是ITRS的一个具体实现,所以ITRF与CGCS2000的转换就可以理解为ITRF框架的相互转换.一般来说,坐标转换可以采用下面三种方式来完成

  1. 通过IGS的联测,平差与计算获得CGCS2000的坐标
  2. 通过公共点解算,获得七参数法中的相应参数,应用这些参数进行转换.
  3. 通过使用ITRF提供的参数结合历元间变换实现.

方案一很明显计算量较大,且需要大量IGS站点的数据,这很多时候要求难以达到.而方案二很难保证公共点的精度问题,也就是说,使用这种解算方式无法保证转换精度.本文主要来探讨使用方案三来实现坐标框架之间的转换.
使用方案三有三大核心要求必须要满足:

  1. 控制点的坐标
  2. 坐标框架间的转换参数
  3. 当前测区的速度场

其中控制点的坐标是已知的要求,而ITRF之间的转换参数IERS已提供,而当前测区的速度场在实际计算过程难以获得,一般采用加权平均法和NNR-NUVEL-1A运动板块模型获得.

在ITRF网站上能够获得取到14个参数,其中七个参数(Helmert parameters,赫尔默特参数)为在指定初始历元下的具体值,而另外七个参数则为参数的变化速率.如果想要获得指定历元下的参数值,计算方式为:
\quad P=P_0 + V_P(t - t_0) \qquad  (1)
其中P为指定历元的参数值,而P_0则为初始参数值, V_P为该参数的变化速度,而t_0则为指定的初始历元.通过上面计算,即可获得指定历元和坐标框架下的七个转换参数.
然而再使用七参数转换法中布尔沙模型(亦称为赫尔默特变换)来完成坐标框架之前的转换.

| XS |    | X |   | Tx |   |  D   -Rz   Ry | | X |
|    |    |   |   |    |   |               | |   |
| YS |  = | Y | + | Ty | + |  Rz   D   -Rx | | Y |  (2)
|    |    |   |   |    |   |               | |   |
| ZS |    | Z |   | Tz |   | -Ry   Rx   D  | | Z |

式(2)中(X,Y,Z)为被转换的坐标,而(XS,YS,ZS)则为转换后的坐标,(Tx,Ty,Tz,D,Rx,Ry,Rz)为之前计算出来的指定历元和框架下的坐标转换参数.地面上的观测点会随着时间进行移和升降,在参考历元的时间跨度较大情况下,需要考虑板块漂移和地壳形变特多方面因素的影响,故不同历元基准的坐标需要进行历元基准的归算.

为t1时刻下的速度场年变化率.低精度的速度场场对于转换精度影响较大,然而该速度场平时难以获得,由NNR-NUVEL-IA的运动板块模型计算出来的点位精度误差一般较大,也不宜采用.所以在实际应用当中将ITRF与CGCS2000高精度转换还是比较困难的.

4. 实现范例

本次主要实现2016年北斗系统由精密星历内插出来的卫星位置坐标(IGb08)转换到CGCS2000坐标系下.IGb08与ITRF08相比,虽然仅仅只采用了GPS的数据(ITRF采用多种手段),不过在精度上两者相差无几,即使在高精度的定位中,也可以认为两种是同一的.而由于卫星在空中不受地面速度场的影响,所以本次从标转换实质是将当前历元下的
ITRF08坐标转化为ITRF97的坐标(注:此处存在疑问,可能理解有误)首先从ITRF获取转换参数,目前ITRF网站已提供ITRF各个框架之前转换参数(而之前只提供相邻框架的转换参数,往往需要转换数次才能达到想想要坐标框架下的坐标).如下所示:

--------------------------------------------------------------------
SOLUTION      Tx    Ty     Tz     D     Rx      Ry      Rz    EPOCH
UNITS-------> mm    mm     mm    ppb    .001"   .001"   .001"
              .     .      .      .     .       .       .
       RATES  Tx    Ty     Tz     D     Rx      Ry      Rz
UNITS-------> mm/y  mm/y   mm/y  ppb/y .001"/y .001"/y .001"/y
---------------------------------------------------------------------
  ITRF2005    -2.0  -0.9   -4.7   0.94   0.00    0.00    0.00  2000.0
       rates   0.3   0.0    0.0   0.00   0.00    0.00    0.00
  ITRF2000    -1.9  -1.7  -10.5   1.34   0.00    0.00    0.00  2000.0
       rates   0.1   0.1   -1.8   0.08   0.00    0.00    0.00
  ITRF97       4.8   2.6  -33.2   2.92   0.00    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
  ITRF96       4.8   2.6  -33.2   2.92   0.00    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
  ITRF94       4.8   2.6  -33.2   2.92   0.00    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
  ITRF93     -24.0   2.4  -38.6   3.41  -1.71   -1.48   -0.30  2000.0
       rates  -2.8  -0.1   -2.4   0.09  -0.11   -0.19    0.07
  ITRF92      12.8   4.6  -41.2   2.21   0.00    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
  ITRF91      24.8  18.6  -47.2   3.61   0.00    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
  ITRF90      22.8  14.6  -63.2   3.91   0.00    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
  ITRF89      27.8  38.6 -101.2   7.31   0.00    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
  ITRF88      22.8   2.6 -125.2  10.41   0.10    0.00    0.06  2000.0
       rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02
_____________________________________________________________________

因为CGCS2000采用的是ITRF1997的框架,故采用的十四个参数如下,其中参考初始历元为2010.0

ITRF97       4.8   2.6  -33.2   2.92   0.00    0.00    0.06  2000.0
     rates   0.1  -0.5   -3.2   0.09   0.00    0.00    0.02

然后再利用公式(1)和(2)计算即可得到该坐标在ITRF97框架下的坐标.(参考代码见附件)

4.总结

通过之前的摸索和研究,主要收获有可总结下面几点:

  1. ITRF是ITRS的具体实现,而相同的ITRF框架在不同的参考历元下的坐标系是不相同的;
  2. CGCS2000本质就是ITRS一个具体实现;CGCS2000与ITRF之间的转化实际上就是同框架间的转化;
  3. ITRF的坐标框架框架是随时间变化的,也就意味着两个坐标系在某历元下即便精度完全相同(坐标差异为0),那么随着时间的变化,同一参考点在两坐标系下的坐标也会有偏差;
  4. WGS84经过多改修正后,与IGS和ITRF之前差异已经到了毫米级,很多情况可以认为三者坐标是同一的.

虽然收获颇多,然而本次研究并没有完全解决ITRF与CGCS2000的转化问题,存在的问题主要有:

  1. 对于同一ITRF框架不同历元基准的转换,地面参考点选用的是速度场的变化,而对于空中坐标如何使用速度场问题没有弄清;
  2. 以2000.0为基准参考历元的ITRF97框架与ITRF00框架的差异没有弄清;

附件

利用七参数转换法将ITRF2008坐标转换成ITRF97的C/C++函数代码

/*Transformat ITRF2008 to ITRF97
*I: (x,y,z) the coordinates in ITRF2014
*I: epoch   the epoch of the coordiantes
*O: (xs,ys,zs) the coordinates in CRCS2000
* This function just transformat  ITRF08 to  ITRF97
*/
void ITRF08toCRCS2000(const double& x, const double&y,const double
                 const double& epoch,double& xs,double& ys, double
{
    double P[] = {4.8,2.6,-33.2,2.92,0.00,0.00,0.06};
    double vP[] = {0.1, -0.5, -3.2, 0.08, 0.00, 0.00, 0.02};
    double ref_ehoch = 2000.0;
    double tP[7] = {0};
    for(int i = 0; i < 7; ++i)
    {
        tP[i] = P[i] + vP[i]*(epoch - ref_ehoch);
    }
    xs = x + tP[0] + tP[3]*x - tP[6]*y + tP[5]*z;
    ys = y + tP[1] + tP[6]*x + tP[3]*y - tP[4]*z;
    zs = z + tP[2] - tP[5]*x + tP[4]*y - tP[3]*z;
    return;
}

参考文献

1.苗龙.ITRF2008与CGCS2000坐标系的转换.地理信息系统.2011年12月(http://www.gissky.net/paper/UploadFiles_4495/201208/2012081721174571.pdf) 

此条目发表在GNSS分类目录,贴了标签。将固定链接加入收藏夹。

发表评论

电子邮件地址不会被公开。