RTKLIB中的卫星天线与接收机天线修正

文章目录

接收机天线相位中心模型
卫星天线相位中心模型
RTKLIB中相关源码
- 一、读取天线文件、存储相关参数
- 二、卫星端PCO修正
- 三、卫星端PCV修正
- 四、接收机端PCO和PCV修正
- 五、PPP中的修正流程
- 六、rtkpost中的修正流程
参考文献

天线相位中心即天线接收信号的电气中心，其空间位置在出厂时往往不在天线的几何中心上。天线所辐射出的电磁波在离开天线一定的距离后，其等相位面会近似为一个球面，该球面的球心即为该天线的等效相位中心，即天线相位中心（Antenna Phase Center ）。

NGS提供的ANTEX格式天线模型，包含了卫星天线模型以及部分接收机天线模型。使用天线模型的目的包括：

修正天线参考点和天线相位中心的之间的偏差；
修正和仰角有关的误差；
修正L1和L2之间的相位中心偏差（这个误差可能对整周模糊度固定造成影响）

接收机天线相位中心模型

一般选取接收机天线底部与天线中轴的交点作为参考点（称天线参考点，ARP）：

ARP与实际相位中心的几何偏差值称为天线相位中心偏差（PCO）
由不同高度角、方位角测得的距离产生系统性的测量偏差为天线相位中心变化（PCV）

目前RTKLIB支持NGS PCV以及ANTEX格式的天线模型，其中包括了PCO和PCV修正参数。通过手册E.8章节可知，接收机天线修正如下：

PCO通常是当地坐标系NED参数，因此需要利用转换矩阵转到ECEF坐标系
PCV修正则通过对高度角进行插值得到

卫星天线相位中心模型

卫星天线偏移（Satellite Antenna Offsets）是指卫星的天线相位中心与卫星质心的偏移。在PPP中我们需要用到精密星历。IGS 提供的精密星历(卫星位置)是相对卫星质心的，而载波相位、伪距测量值是相对天线相位中心的，因此我们需要计算PCV和PCO值，将所有值统一转换成相对天线相位中心的。

具体计算公式可以参考手册的E.8章节。

RTKLIB中相关源码

一、读取天线文件、存储相关参数

openses函数：位于postpos.c中，调用readpcv函数读取卫星天线参数和接收机天线参数

    /* read satellite antenna parameters */if (*fopt->satantp&&!(readpcv(fopt->satantp,pcvs))) {showmsg("error : no sat ant pcv in %s",fopt->satantp);trace(1,"sat antenna pcv read error: %s\n",fopt->satantp);return 0;}/* read receiver antenna parameters */if (*fopt->rcvantp&&!(readpcv(fopt->rcvantp,pcvr))) {showmsg("error : no rec ant pcv in %s",fopt->rcvantp);trace(1,"rec antenna pcv read error: %s\n",fopt->rcvantp);return 0;}

readpcv函数代码如下：

/* read antenna parameters ------------------------------------------------------*/
extern int readpcv(const char *file, pcvs_t *pcvs)
{pcv_t *pcv;char *ext;int i,stat;trace(3,"readpcv: file=%s\n",file);if (!(ext=strrchr(file,'.'))) ext="";if (!strcmp(ext,".atx")||!strcmp(ext,".ATX")) {stat=readantex(file,pcvs);}else {stat=readngspcv(file,pcvs);}for (i=0;i<pcvs->n;i++) {pcv=pcvs->pcv+i;trace(4,"sat=%2d type=%20s code=%s off=%8.4f %8.4f %8.4f  %8.4f %8.4f %8.4f\n",pcv->sat,pcv->type,pcv->code,pcv->off[0][0],pcv->off[0][1],pcv->off[0][2],pcv->off[1][0],pcv->off[1][1],pcv->off[1][2]);}return stat;
}

setpcv函数（设置天线参数）：位于postpos.c中，查找各颗卫星的天线参数，并放在nav->pcvs中；查找接收机的天线参数，并放在popt->pcvr中，需要注意的是，只有当popt->anttype[i]是*(auto)的时候，通过RINEX头识别到的接收机天线类型，以及天线的H/E/N偏差才会起作用。

/* set antenna parameters ----------------------------------------------------*/
static void setpcv(gtime_t time, prcopt_t *popt, nav_t *nav, const pcvs_t *pcvs,const pcvs_t *pcvr, const sta_t *sta)
{pcv_t *pcv,pcv0={0};double pos[3],del[3];int i,j,mode=PMODE_DGPS<=popt->mode&&popt->mode<=PMODE_FIXED;char id[64];/* set satellite antenna parameters */for (i=0;i<MAXSAT;i++) {nav->pcvs[i]=pcv0;if (!(satsys(i+1,NULL)&popt->navsys)) continue;if (!(pcv=searchpcv(i+1,"",time,pcvs))) {satno2id(i+1,id);trace(3,"no satellite antenna pcv: %s\n",id);continue;}nav->pcvs[i]=*pcv;}for (i=0;i<(mode?2:1);i++) {popt->pcvr[i]=pcv0;if (!strcmp(popt->anttype[i],"*")) { /* set by station parameters */strcpy(popt->anttype[i],sta[i].antdes);if (sta[i].deltype==1) { /* xyz */if (norm(sta[i].pos,3)>0.0) {ecef2pos(sta[i].pos,pos);ecef2enu(pos,sta[i].del,del);for (j=0;j<3;j++) popt->antdel[i][j]=del[j];}}else { /* enu */for (j=0;j<3;j++) popt->antdel[i][j]=stas[i].del[j];}}if (!(pcv=searchpcv(0,popt->anttype[i],time,pcvr))) {trace(2,"no receiver antenna pcv: %s\n",popt->anttype[i]);*popt->anttype[i]='\0';continue;}strcpy(popt->anttype[i],pcv->type);popt->pcvr[i]=*pcv;}
}

二、卫星端PCO修正

satantoff函数：位于preceph.c中，分别计算L1和L2的卫星端PCO。由于PPP中使用的是无电离层组合，因此会计算无电离层组合后的PCO修正值。之后在peph2pos函数中，会在由精密星历计算的卫星位置上进行PCO修正。

extern void satantoff(gtime_t time, const double *rs, int sat, const nav_t *nav,double *dant)
{const double *lam=nav->lam[sat-1];const pcv_t *pcv=nav->pcvs+sat-1;double ex[3],ey[3],ez[3],es[3],r[3],rsun[3],gmst,erpv[5]={0};double gamma,C1,C2,dant1,dant2;int i,j=0,k=1;trace(4,"satantoff: time=%s sat=%2d\n",time_str(time,3),sat);/* sun position in ecef */sunmoonpos(gpst2utc(time),erpv,rsun,NULL,&gmst);/* unit vectors of satellite fixed coordinates */for (i=0;i<3;i++) r[i]=-rs[i];if (!normv3(r,ez)) return;for (i=0;i<3;i++) r[i]=rsun[i]-rs[i];if (!normv3(r,es)) return;cross3(ez,es,r);if (!normv3(r,ey)) return;cross3(ey,ez,ex);if ((NFREQ>=3)&&(satsys(sat,NULL)&SYS_SBS)) k=2;  /* always use L5 for SBAS */if ((NFREQ>=3)&&(satsys(sat,NULL)&SYS_GAL)&&(lam[1]==0)) k=2; /* use E5a if no E5b for GAL */if (NFREQ<2||lam[j]==0.0||lam[k]==0.0) return;gamma=SQR(lam[k])/SQR(lam[j]);C1=gamma/(gamma-1.0);C2=-1.0 /(gamma-1.0);/* iono-free LC */for (i=0;i<3;i++) {dant1=pcv->off[j][0]*ex[i]+pcv->off[j][1]*ey[i]+pcv->off[j][2]*ez[i];dant2=pcv->off[k][0]*ex[i]+pcv->off[k][1]*ey[i]+pcv->off[k][2]*ez[i];dant[i]=C1*dant1+C2*dant2;}
}

三、卫星端PCV修正

satantpcv函数：计算卫星天线端PCV修正值，具体计算方法可以参考手册E.8：

/* satellite antenna phase center variation ----------------------------------*/
static void satantpcv(const double *rs, const double *rr, const pcv_t *pcv,double *dant)
{double ru[3],rz[3],eu[3],ez[3],nadir,cosa;int i;for (i=0;i<3;i++) {ru[i]=rr[i]-rs[i];rz[i]=-rs[i];}if (!normv3(ru,eu)||!normv3(rz,ez)) return;cosa=dot(eu,ez,3);cosa=cosa<-1.0?-1.0:(cosa>1.0?1.0:cosa);nadir=acos(cosa);antmodel_s(pcv,nadir,dant);
}

四、接收机端PCO和PCV修正

antmodel函数计算接收机端修正值，需要注意的是，代码中实际计算方法和手册不太一致。代码先计算line of sight在NED中的单位矢量e，然后计算PCO在line of sight上的投影长度，并没有将PCO从NED转换到ECEF这一过程。

extern void antmodel(const pcv_t *pcv, const double *del, const double *azel,int opt, double *dant)
{double e[3],off[3],cosel=cos(azel[1]);int i,j;trace(4,"antmodel: azel=%6.1f %4.1f opt=%d\n",azel[0]*R2D,azel[1]*R2D,opt);e[0]=sin(azel[0])*cosel;e[1]=cos(azel[0])*cosel;e[2]=sin(azel[1]);for (i=0;i<NFREQ;i++) {for (j=0;j<3;j++) off[j]=pcv->off[i][j]+del[j];dant[i]=-dot(off,e,3)+(opt?interpvar(90.0-azel[1]*R2D,pcv->var[i]):0.0);}trace(4,"antmodel: dant=%6.3f %6.3f\n",dant[0],dant[1]);
}

五、PPP中的修正流程

ppp_res函数：位于ppp.c中，通过satantpcv函数计算卫星天线端PCV修正，通过antmodel函数计算接收机天线端修正值，最后在corr_meas函数中进行修正。

        /* satellite and receiver antenna model */if (opt->posopt[0]) satantpcv(rs+i*6,rr,nav->pcvs+sat-1,dants);antmodel(opt->pcvr,opt->antdel[0],azel+i*2,opt->posopt[1],dantr);/* phase windup model */if (!model_phw(rtk->sol.time,sat,nav->pcvs[sat-1].type,opt->posopt[2]?2:0,rs+i*6,rr,&rtk->ssat[sat-1].phw)) {continue;}/* corrected phase and code measurements */corr_meas(obs+i,nav,azel+i*2,&rtk->opt,dantr,dants,rtk->ssat[sat-1].phw,L,P,&Lc,&Pc);

在corr_meas函数中进行修正，如下所示，在伪距、载波相位上减去了dants和dantr。

        /* antenna phase center and phase windup correction */L[i]=obs->L[i]*CLIGHT/freq[i]-dants[i]-dantr[i]-phw*CLIGHT/freq[i];P[i]=obs->P[i]       -dants[i]-dantr[i];

六、rtkpost中的修正流程

zdres函数:位于rtkpos.c中，调用antmodel函数计算接收机端修正值，在zdres_sat函数中进行修正。需要注意的是，相对定位只需要计算接收机端的天线相位中心修正值。这是由于，相对定位进行单差时，已经将卫星端的天线误差消除了。但是PPP则需要同时处理接收机和卫星端误差。

        /* calc receiver antenna phase center correction */antmodel(opt->pcvr+index,opt->antdel[index],azel+i*2,opt->posopt[1],dant);/* calc undifferenced phase/code residual for satellite */zdres_sat(base,r,obs+i,nav,azel+i*2,dant,opt,y+i*nf*2,freq+i*nf);

参考文献

[1] GPS从入门到放弃（二十） — 天线偏移
[2] 天线高改正和天线相位中心改正
[3] 101-RTKLIB中的天线相位中心改正