RTKLIB专题学习(五)—单点定位实现进阶(二)
今天我们继续来了解一下,RTKLIB中的单点定位是如何实现的:
上一篇RTKLIB专题学习(五)—单点定位实现进阶(一)讲到了得到了最小二乘计算结果之后,需要进行结果有效性检验。
那么今天我们一起来看看,在调用完estpos之后,pntpos又做了哪些工作吧!
1.当最小二乘结果无效的时候,进行错误的探测与剔除
/* estimate receiver position with pseudorange */
stat=estpos(obs,n,rs,dts,var,svh,nav,&opt_,sol,azel_,vsat,resp,msg);
/* RAIM FDE */
if (!stat&&n>=6&&opt->posopt[4]) {
stat=raim_fde(obs,n,rs,dts,var,svh,nav,&opt_,sol,azel_,vsat,resp,msg);
}
使用到的是raim_fde
函数:
/* RAIM FDE (failure detection and exclution) -------------------------------*/
static int raim_fde(const obsd_t *obs, int n, const double *rs,
const double *dts, const double *vare, const int *svh,
const nav_t *nav, const prcopt_t *opt, sol_t *sol,
double *azel, int *vsat, double *resp, char *msg)
{
obsd_t *obs_e;
sol_t sol_e={{0}};
char tstr[32],name[16],msg_e[128];
double *rs_e,*dts_e,*vare_e,*azel_e,*resp_e,rms_e,rms=100.0;
int i,j,k,nvsat,stat=0,*svh_e,*vsat_e,sat=0;
trace(3,"raim_fde: %s n=%2d\n",time_str(obs[0].time,0),n);
if (!(obs_e=(obsd_t *)malloc(sizeof(obsd_t)*n))) return 0;
rs_e = mat(6,n); dts_e = mat(2,n); vare_e=mat(1,n); azel_e=zeros(2,n);
svh_e=imat(1,n); vsat_e=imat(1,n); resp_e=mat(1,n);
for (i=0;i<n;i++) {
/* satellite exclution */
for (j=k=0;j<n;j++) {
if (j==i) continue;
obs_e[k]=obs[j];
matcpy(rs_e +6*k,rs +6*j,6,1);
matcpy(dts_e+2*k,dts+2*j,2,1);
vare_e[k]=vare[j];
svh_e[k++]=svh[j];
}
/* estimate receiver position without a satellite */
if (!estpos(obs_e,n-1,rs_e,dts_e,vare_e,svh_e,nav,opt,&sol_e,azel_e,
vsat_e,resp_e,msg_e,i)) {
trace(3,"raim_fde: exsat=%2d (%s)\n",obs[i].sat,msg);
continue;
}
for (j=nvsat=0,rms_e=0.0;j<n-1;j++) {
if (!vsat_e[j]) continue;
rms_e+=SQR(resp_e[j]);
nvsat++;
}
if (nvsat<5) {
trace(3,"raim_fde: exsat=%2d lack of satellites nvsat=%2d\n",
obs[i].sat,nvsat);
continue;
}
rms_e=sqrt(rms_e/nvsat);
trace(3,"raim_fde: exsat=%2d rms=%8.3f\n",obs[i].sat,rms_e);
if (rms_e>rms) continue;
/* save result */
for (j=k=0;j<n;j++) {
if (j==i) continue;
matcpy(azel+2*j,azel_e+2*k,2,1);
vsat[j]=vsat_e[k];
resp[j]=resp_e[k++];
}
stat=1;
*sol=sol_e;
sat=obs[i].sat;
rms=rms_e;
vsat[i]=0;
strcpy(msg,msg_e);
}
if (stat) {
time2str(obs[0].time,tstr,2); satno2id(sat,name);
trace(2,"%s: %s excluded by raim\n",tstr+11,name);
}
free(obs_e);
free(rs_e ); free(dts_e ); free(vare_e); free(azel_e);
free(svh_e); free(vsat_e); free(resp_e);
return stat;
}
该函数的意义就是,每次排除一颗卫星,然后重新进行单点定位,若是计算成功,则剔除的这颗卫星就是有故障的,否则继续寻找故障卫星。成功寻找卫星并结算成功,函数值返回1;否则返回0,即该历元结算失败。
2.在进行完错误检测与排除后,pntpos
进行如下工作,并返回1(计算成功)或0(计算失败)
/* estimate receiver velocity with Doppler */
if (stat) {
estvel(obs,n,rs,dts,nav,&opt_,sol,azel_,vsat);
}
if (azel) {
for (i=0;i<n*2;i++) azel[i]=azel_[i];
}
if (ssat) {
for (i=0;i<MAXSAT;i++) {
ssat[i].vs=0;
ssat[i].azel[0]=ssat[i].azel[1]=0.0;
ssat[i].resp[0]=ssat[i].resc[0]=0.0;
ssat[i].snr[0]=0;
}
for (i=0;i<n;i++) {
ssat[obs[i].sat-1].azel[0]=azel_[ i*2];
ssat[obs[i].sat-1].azel[1]=azel_[1+i*2];
ssat[obs[i].sat-1].snr[0]=obs[i].SNR[0];
if (!vsat[i]) continue;
ssat[obs[i].sat-1].vs=1;
ssat[obs[i].sat-1].resp[0]=resp[i];
}
}
free(rs); free(dts); free(var); free(azel_); free(resp);
return stat;
}
3.接着在调用pntpos
之后,rtkpos
会进行下列操作,完成rtkpos
中的单点定位操作
/* rover position by single point positioning */
if (!pntpos(obs,nu,nav,&rtk->opt,&rtk->sol,NULL,rtk->ssat,msg)) {
errmsg(rtk,"point pos error (%s)\n",msg);
if (!rtk->opt.dynamics) {
outsolstat(rtk);
return 0;
}
}
if (time.time != 0) rtk->tt = timediff(rtk->sol.time, time); opt->tt = rtk->tt;
/* single point positioning */
if (opt->mode==PMODE_SINGLE) {
outsolstat(rtk);
return 1;
}
4.继续进行,在rtkpos
处理完之后,若前面解算成功,则函数返回值为1,在procpos
里面继续进行:
if (!rtkpos(&rtk,obs,n,&navs)) continue;
if (mode==0) { /* forward/backward */
if (!solstatic) {
outsol(fp,&rtk.sol,rtk.rb,sopt);
}
else if (time.time==0||pri[rtk.sol.stat]<=pri[sol.stat]) {
sol=rtk.sol;
for (i=0;i<3;i++) rb[i]=rtk.rb[i];
if (time.time==0||timediff(rtk.sol.time,time)<0.0) {
time=rtk.sol.time;
}
}
}
else if (!revs) { /* combined-forward */
if (isolf>=nepoch) return;
solf[isolf]=rtk.sol;
for (i=0;i<3;i++) rbf[i+isolf*3]=rtk.rb[i];
isolf++;
}
else { /* combined-backward */
if (isolb>=nepoch) return;
solb[isolb]=rtk.sol;
for (i=0;i<3;i++) rbb[i+isolb*3]=rtk.rb[i];
isolb++;
}
}
if (mode==0&&solstatic&&time.time!=0.0) {
sol.time=time;
outsol(fp,&sol,rb,sopt);
}
rtkfree(&rtk);
可以从这句看出,若前面解算失败,也就是rtkpos
返回值为0,则不进行后面的操作
5.这样,单点定位就算是结算完毕啦,后面的操作其实就是输出结果、关闭一些文件流之类的操作了!(还会有调用procpos
的函数,前面几篇讲过啦,这里不再重复)