本博客是转载,原文请看下面的连接
感谢:RTKLIB单点定位处理流程之一(postpos/后处理)_wuwuku123的博客-CSDN博客_rtklib伪距单点定位
版本:
#define VER_RTKLIB "2.4.3" /* library version */
#define PATCH_LEVEL "b29" /* patch level */
处理参数设置:主要是三个结构体
- prcopt_t opt = prcopt_default; // processing options
- solopt_t sopt = solopt_default; // solution options
- filopt_t fopt = { 0 }; // file options
#include "rtklib.h"
int main() {
char* infile[] = { {"移动站观测数据"},
{"移动站星历文件"} };
char* ofile = { "输出结果文件/myposbyVS.pos" };
gtime_t ts = { 0 }, te = { 0 };
double st[] = { 2019,12,5,0,0,0 }; //开始处理时间
double ed[] = { 2019,12,5,1,0,0 }; //结束处理时间
ts = epoch2time(st);
te = epoch2time(ed);
double ti = 0.0; // processing interval(s) (0:all) 时间间隔
double tu = 0.0; // processing unit time(s) (0:all)
prcopt_t opt = prcopt_default; //processing options
solopt_t sopt = solopt_default; //solution options
filopt_t fopt = { 0 }; //file options
int n = 2; //number of input files
//GPS单点定位参数设置
opt.navsys = SYS_GPS;
opt.mode = PMODE_SINGLE; // PMODE_SINGLE; // PMODE_MOVEB; // ; PMODE_FIXED
opt.nf = 1; /* number of frequencies (1:L1,2:L1+L2,3:L1+L2+L5) */
opt.elmin = 10 * D2R;
opt.ionoopt = IONOOPT_BRDC; /* ionosphere option (IONOOPT_???) */
opt.tropopt = TROPOPT_SAAS; /* troposphere option (TROPOPT_???) */
opt.dynamics = 0; /* dynamics model (0:none,1:velociy,2:accel) */
opt.tidecorr = 0; /* earth tide correction (0:off,1:solid,2:solid+otl+pole) */
opt.niter = 1; /* number of filter iteration */
opt.codesmooth = 0; /* code smoothing window size (0:none) */
opt.rovpos = 0; /* rover position for fixed mode */
opt.refpos = 0; /* base position for relative mode */
/* (0:pos in prcopt, 1:average of single pos, */
/* 2:read from file, 3:rinex header, 4:rtcm pos) */
//opt.eratio[] = {100.0,100.0}; /* code/phase error ratio */
//opt.err[] = { 100.0,0.003,0.003,0.0,1.0 }; /* err[] */
/* measurement error factor */
/* [0]:reserved */
/* [1-3]:error factor a/b/c of phase (m) */
/* [4]:doppler frequency (hz) */
//opt.std[] = { 30.0,0.03,0.3 }; /* std[] */
/* initial-state std [0]bias,[1]iono [2]trop */
//opt.prn[] = { 1E-4,1E-3,1E-4,1E-1,1E-2 ,0.0 };
/* process-noise std [0]bias,[1]iono [2]trop [3]acch [4]accv [5] pos */
//求解选项定义
sopt.posf = SOLF_LLH; // SOLF_XYZ; // SOLF_ENU;
sopt.times = TIMES_GPST; /* time system: gps time */
sopt.timef = 1; /* time format (0:sssss.s,1:yyyy/mm/dd hh:mm:ss.s) */
sopt.timeu = 3; /* time digits under decimal point */
sopt.degf = 0; /* latitude/longitude format (0:ddd.ddd,1:ddd mm ss) */
sopt.outhead = 1; /* output header (0:no,1:yes) */
sopt.outopt = 1; /* output processing options (0:no,1:yes) */
sopt.outvel = 1; /* output velocity options (0:no,1:yes) */
sopt.datum = 0; /* datum (0:WGS84,1:Tokyo) */
sopt.height = 0; /* height (0:ellipsoidal,1:geodetic) */
sopt.geoid = 0; /* geoid model (0:EGM96,1:JGD2000) */
sopt.solstatic = 0; /* solution of static mode (0:all,1:single) */
sopt.sstat = 1; /* solution statistics level (0:off,1:states,2:residuals) */
sopt.trace = 3; /* debug trace level (0:off,1-5:debug) */
//sopt.nmeaintv[] = { 0.0,0.0 }; /* nmea output interval (s) (<0:no,0:all) */
/* nmeaintv[0]:gprmc,gpgga,nmeaintv[1]:gpgsv */
//sopt.sep[] = " "; /* field separator */
//sopt.prog[] = ""; /* program name */
//sopt.maxsolstd = 0; /* max std-dev for solution output (m) (0:all) */
//文件设置选项 /* file options type */
// #define MAXSTRPATH 1024 /* max length of stream path */
//fopt.satantp[MAXSTRPATH] = { 0 }; /* satellite antenna parameters file */
//fopt.rcvantp[MAXSTRPATH] = { 0 }; /* receiver antenna parameters file */
//fopt.stapos[MAXSTRPATH] = { 0 }; /* station positions file */
//fopt.geoid[MAXSTRPATH] = { 0 }; /* external geoid data file */
//fopt.iono[MAXSTRPATH] = { 0 }; /* ionosphere data file */
//fopt.dcb[MAXSTRPATH] = { 0 }; /* dcb data file */
//fopt.eop[MAXSTRPATH] = { 0 }; /* eop data file */
//fopt.blq[MAXSTRPATH] = { 0 }; /* ocean tide loading blq file */
//fopt.tempdir[MAXSTRPATH] = { 0 }; /* ftp/http temporaly directory */
//fopt.geexe[MAXSTRPATH] = { 0 }; /* google earth exec file */
//fopt.solstat[MAXSTRPATH] = { 0 }; /* solution statistics file */
//fopt.trace[MAXSTRPATH] = { 0 }; /* debug trace file */
postpos(ts,te,ti,tu,&opt,&sopt,&fopt,infile,n,ofile,"","");
printf("*********postpos over*********\n");
getchar();
return 0;
}
函数调用关系:
首先通过函数:postpos(ts,te,ti,tu,&opt,&sopt,&fopt,infile,n,ofile,"","");
进入子函数,对其内部的函数解释如下:
1、openses函数:预处理
目的:通过读取 filopt_t fopt 中的参数,其中filopt_t结构体定义如下,实现对‘文件的读取’,读取子函数为 openses!
typedef struct { /* file options type */
char satantp[MAXSTRPATH]; /* satellite antenna parameters file */
char rcvantp[MAXSTRPATH]; /* receiver antenna parameters file */
char stapos [MAXSTRPATH]; /* station positions file */
char geoid [MAXSTRPATH]; /* external geoid data file */
char iono [MAXSTRPATH]; /* ionosphere data file */
char dcb [MAXSTRPATH]; /* dcb data file */
char eop [MAXSTRPATH]; /* eop data file */
char blq [MAXSTRPATH]; /* ocean tide loading blq file */
char tempdir[MAXSTRPATH]; /* ftp/http temporaly directory */
char geexe [MAXSTRPATH]; /* google earth exec file */
char solstat[MAXSTRPATH]; /* solution statistics file */
char trace [MAXSTRPATH]; /* debug trace file */
} filopt_t;
/* open procssing session-------------------------------*/
static int openses(const prcopt_t *popt, const solopt_t *sopt, const filopt_t *fopt, nav_t *nav, pcvs_t *pcvs, pcvs_t *pcvr)
- 然后为:判断处理时间是否正确,否则程序终止;
- 然后为:通过malloc函数为ifile[],动态分配内存,否则程序终止;
- 然后为:进行文件的读取,对于sp3/SP3、.eph/EPH文件的时间,起使和终止相差1小时;brdc文件的开始和终止时间差2小时。
- 然后为:路径替换,实现语句:
nf+=reppaths(infile[j],ifile+nf,MAXINFILE-nf,tts,ttte,"",""); 子函数内实现语句: reppath(path,rpath[n],time,rov,base);
经过两次读取循环,已经将观测、导航文件读取到ifie数组中:
将输出结果文件路径保存到ofile数组中:
if (!reppath(outfile,ofile,tts,"","")&&i>0) flag=0; 子函数实现代码: strcpy(rpath,path);
2、execses_b
在main中,输入文件中:第一个为观测文件,第二个为导航文件;
2.1 readpreceph
/* 读取文件为:read prec ephemeris, sbas data, lex data, tec grid and open rtcm -*/
static void readpreceph(char **infile, int n, const prcopt_t *prcopt,nav_t *nav, sbs_t *sbs, lex_t *lex)读取文件的函数依次为:
/* read precise ephemeris files */
/* read precise clock files */
/* read satellite fcb files */
/* read solution status files for ppp correction */
/* read sbas message files */
/* read lex message files */
/* allocate sbas ephemeris */
/* set rtcm file and initialize rtcm struct */
以函数:readsp3为例
/* read precise ephemeris files */ for (i=0;i<n;i++) { if (strstr(infile[i],"%r")||strstr(infile[i],"%b")) continue; readsp3(infile[i],nav,0); }
第一次循环读入的是观测文件,因此进入函数之后,动态分配内存、释放内存,结束!
第二次循环读入的是导航文件,因此进入函数之后,动态分配内存、释放内存,结束!
在函数readrnxc中:读取钟差文件时,读取导航文件时,将导航文件,文件头的相关参数保存到nav中:
/* read precise clock files */ for (i=0;i<n;i++) { if (strstr(infile[i],"%r")||strstr(infile[i],"%b")) continue; readrnxc(infile[i],nav); }
导航文件,文件头为:
nav:
2.2、execses_r
stat=execses_r(ts,te,ti,popt,sopt,fopt,flag,infile,index,n,outfile,rov);
2.2.1 execsesstat=execses(ts,te,ti,popt,sopt,fopt,flag,infile,index,n,outfile);
如果main函数中,设置了调试等级:
sopt.sstat = 1; /* solution statistics level (0:off,1:states,2:residuals) */ sopt.trace = 3; /* debug trace level (0:off,1-5:debug) */
则在此处新建生成的调试数据,保存的路径、文件名
/* open debug trace */ if (flag&&sopt->trace>0) { if (*outfile) { strcpy(tracefile,outfile); strcat(tracefile,".trace"); //输出文件名,后加.trace } else { strcpy(tracefile,fopt->trace); } traceclose(); //关闭记录文件 traceopen(tracefile); //开始对 trace文件进行写处理 tracelevel(sopt->trace); //读取调试等级 }
顺序执行,子函数如下:
1、/* read ionosphere data file */ 读取数据放到如下结构体中:
typedef struct { /* TEC grid type */ gtime_t time; /* epoch time (GPST) */ int ndata[3]; /* TEC grid data size {nlat,nlon,nhgt} */ double rb; /* earth radius (km) */ double lats[3]; /* latitude start/interval (deg) */ double lons[3]; /* longitude start/interval (deg) */ double hgts[3]; /* heights start/interval (km) */ double *data; /* TEC grid data (tecu) */ float *rms; /* RMS values (tecu) */ } tec_t;
2、/* read erp data */ 读取数据放到如下结构体中:
typedef struct { /* earth rotation parameter type */ int n,nmax; /* number and max number of data */ erpd_t *data; /* earth rotation parameter data */ } erp_t;
3、/* read obs and nav data */
其中观测数据结构体的使用:
typedef struct { /* observation data record */ gtime_t time; /* receiver sampling time (GPST) */ unsigned char sat,rcv; /* satellite/receiver number */ unsigned char SNR [NFREQ+NEXOBS]; /* signal strength (0.25 dBHz) */ unsigned char LLI [NFREQ+NEXOBS]; /* loss of lock indicator */ unsigned char code[NFREQ+NEXOBS]; /* code indicator (CODE_???) */ double L[NFREQ+NEXOBS]; /* observation data carrier-phase (cycle) */ double P[NFREQ+NEXOBS]; /* observation data pseudorange (m) */ float D[NFREQ+NEXOBS]; /* observation data doppler frequency (Hz) */ } obsd_t; typedef struct { /* observation data */ int n,nmax; /* number of obervation data/allocated */ obsd_t *data; /* observation data records */ } obs_t;
- readobsnav
函数: static int readobsnav(gtime_t ts, gtime_t te, double ti, char **infile,const int *index, int n, const prcopt_t *prcopt,obs_t *obs, nav_t *nav, sta_t *sta)
程序中,每次都是将相关的变量置0、指针变量赋值为null
obs->data=NULL; obs->n =obs->nmax =0; nav->eph =NULL; nav->n =nav->nmax =0; nav->geph=NULL; nav->ng=nav->ngmax=0; nav->seph=NULL; nav->ns=nav->nsmax=0; nepoch=0;
- readrnxt
/* read rinex obs and nav files ------------------------------------------------
extern int readrnxt(const char *file, int rcv, gtime_t ts, gtime_t te,double tint, const char *opt, obs_t *obs, nav_t *nav,sta_t *sta)
- readrnxfile
调用子函数:
stat=readrnxfile(files[i],ts,te,tint,opt,0,rcv,&type,obs,nav,sta);
调用子函数:
if (sta) init_sta(sta);
- readrnxfp
stat=readrnxfp(fp,ts,te,tint,opt,flag,index,type,obs,nav,sta);
调用子函数:先读头文件,判断.o或.n;然后进行body读取
- readrnxh
static int readrnxh(FILE *fp, double *ver, char *type, int *sys, int *tsys,
char tobs[][MAXOBSTYPE][4], nav_t *nav, sta_t *sta)
先读取文件头1部分,再通过调用子函数读取2部分(见下图)
读取Rinex文件头的代码实现:while (fgets(buff,MAXRNXLEN,fp)) { if (strlen(buff)<=60) continue; else if (strstr(label,"RINEX VERSION / TYPE")) { *ver=str2num(buff,0,9); *type=*(buff+20); /* satellite system */ switch (*(buff+40)) { case ' ': case 'G': *sys=SYS_GPS; *tsys=TSYS_GPS; break; case 'R': *sys=SYS_GLO; *tsys=TSYS_UTC; break; case 'E': *sys=SYS_GAL; *tsys=TSYS_GAL; break; /* v.2.12 */ case 'S': *sys=SYS_SBS; *tsys=TSYS_GPS; break; case 'J': *sys=SYS_QZS; *tsys=TSYS_QZS; break; /* v.3.02 */ case 'C': *sys=SYS_CMP; *tsys=TSYS_CMP; break; /* v.2.12 */ case 'I': *sys=SYS_IRN; *tsys=TSYS_IRN; break; /* v.3.03 */ case 'M': *sys=SYS_NONE; *tsys=TSYS_GPS; break; /* mixed */ default : trace(2,"not supported satellite system: %c\n",*(buff+40)); break; } continue; } else if (strstr(label,"PGM / RUN BY / DATE")) continue; else if (strstr(label,"COMMENT")) { /* opt */ /* read cnes wl satellite fractional bias */ if (strstr(buff,"WIDELANE SATELLITE FRACTIONAL BIASES")|| strstr(buff,"WIDELANE SATELLITE FRACTIONNAL BIASES")) { block=1; } else if (block) { /* cnes/cls grg clock */ if (!strncmp(buff,"WL",2)&&(sat=satid2no(buff+3))&& sscanf(buff+40,"%lf",&bias)==1) { nav->wlbias[sat-1]=bias; } /* cnes ppp-wizard clock */ else if ((sat=satid2no(buff+1))&&sscanf(buff+6,"%lf",&bias)==1) { nav->wlbias[sat-1]=bias; } } continue; } /* file type */ switch (*type) { case 'O': decode_obsh(fp,buff,*ver,tsys,tobs,nav,sta); break; case 'N': decode_navh (buff,nav); break; case 'G': decode_gnavh(buff,nav); break; case 'H': decode_hnavh(buff,nav); break; case 'J': decode_navh (buff,nav); break; /* extension */ case 'L': decode_navh (buff,nav); break; /* extension */ } if (strstr(label,"END OF HEADER")) return 1; if (++i>=MAXPOSHEAD&&*type==' ') break; /* no rinex file */ }
如下示例一个观测文件的头文件:
其中,数据类型为:C1 L1 P1 C2 L2 P2 C5 L5;rinex版本为:2.11;在程序中数据类型保存为:
调用子函数:
/* decode obs header ---------------------------------------------------------*/
static void decode_obsh(FILE *fp, char *buff, double ver, int *tsys,
char tobs[][MAXOBSTYPE][4], nav_t *nav, sta_t *sta)
- readrnxobs
读取Rinex Body /* read rinex body */
具体实现代码如下:
switch (*type) { case 'O': return readrnxobs(fp,ts,te,tint,opt,index,ver,&tsys,tobs,obs, sta); case 'N': return readrnxnav(fp,opt,ver,sys ,nav); case 'G': return readrnxnav(fp,opt,ver,SYS_GLO,nav); case 'H': return readrnxnav(fp,opt,ver,SYS_SBS,nav); case 'J': return readrnxnav(fp,opt,ver,SYS_QZS,nav); /* extension */ case 'L': return readrnxnav(fp,opt,ver,SYS_GAL,nav); /* extension */ case 'C': return readrnxclk(fp,opt,index,nav); }
读取观测文件,首先为: obsd_t *data; 分配内存:
if (!(data=(obsd_t *)malloc(sizeof(obsd_t)*MAXOBS))) return 0;
其中; #define MAXOBS 64 /* max number of obs in an epoch */
- readrnxobsb // 循环读取观测数据
static int readrnxobsb(FILE *fp, const char *opt, double ver, int *tsys,char tobs[][MAXOBSTYPE][4], int *flag, obsd_t *data,
sta_t *sta)/* read rinex obs data body */ while ((n=readrnxobsb(fp,opt,ver,tsys,tobs,&flag,data,sta))>=0&&stat>=0) { for (i=0;i<n;i++) { /* utc -> gpst */ if (*tsys==TSYS_UTC) data[i].time=utc2gpst(data[i].time); /* save cycle-slip */ saveslips(slips,data+i); } /* screen data by time */ if (n>0&&!screent(data[0].time,ts,te,tint)) continue; for (i=0;i<n;i++) { /* restore cycle-slip */ restslips(slips,data+i); data[i].rcv=(unsigned char)rcv; /* save obs data */ if ((stat=addobsdata(obs,data+i))<0) break; } }
读入Rinex观测文件每一个历元的时间后,进入子函数(decode_obsepoch):读取每一个历元的卫星个数,即下图中的第一个红框
/* decode obs epoch ----------------------------------------------------------*/
static int decode_obsepoch(FILE *fp, char *buff, double ver, gtime_t *time,int *flag, int *sats)
读取上图中的第二个红框:
子函数: /* decode obs data */
if (decode_obsdata(fp,buff,ver,mask,index,data+n)&&n<MAXOBS) n++;每一历元的观测数据读取完成后,对LLI标识的周跳进行处理:
输入文件中,第一个为移动站观测文件;当读取完第一个观测文件后,又开始如下图的循环:
本例程中,存在两个输入文件分别为:移动站观测文件、移动站导航文件;
注意点:当读取到导航文件头文件中的GPS电离层α、β参数时,
程序跳转到此处:
将观测、导航文件的数据读取出来,并保存到相应的结构体中:obs_t *obs, nav_t *nav
下一步就是:
/* sort observation data */
nepoch=sortobs(obs);
4、/* read dcb parameters */
保存在:nav_t *nav 里面的cbias、rbias中:
double cbias[MAXSAT][3]; /* satellite dcb (0:p1-p2,1:p1-c1,2:p2-c2) (m) */ double rbias[MAXRCV][2][3]; /* receiver dcb (0:p1-p2,1:p1-c1,2:p2-c2) (m) */
5、/* set antenna paramters 设置天线参数*/
6、/* read ocean tide loading parameters 设置潮汐参数*/
7、/* rover/reference fixed position 获取基准站坐标*/
#define PMODE_SINGLE 0 /* positioning mode: single */ #define PMODE_DGPS 1 /* positioning mode: DGPS/DGNSS */ #define PMODE_KINEMA 2 /* positioning mode: kinematic */ #define PMODE_STATIC 3 /* positioning mode: static */ #define PMODE_MOVEB 4 /* positioning mode: moving-base */ #define PMODE_FIXED 5 /* positioning mode: fixed */ #define PMODE_PPP_KINEMA 6 /* positioning mode: PPP-kinemaric */ #define PMODE_PPP_STATIC 7 /* positioning mode: PPP-static */ #define PMODE_PPP_FIXED 8 /* positioning mode: PPP-fixed */
#define POSOPT_POS 0 /* pos option: LLH/XYZ */ #define POSOPT_SINGLE 1 /* pos option: average of single pos */ #define POSOPT_FILE 2 /* pos option: read from pos file */ #define POSOPT_RINEX 3 /* pos option: rinex header pos */ #define POSOPT_RTCM 4 /* pos option: rtcm station pos */ #define POSOPT_RAW 5 /* pos option: raw station pos */
8、/* open solution statistics */输出统计信息:solution statistics level (0:off,1:states,2:residuals)
建立路径和保存文件名!
if (flag&&sopt->sstat>0) { strcpy(statfile,outfile); strcat(statfile,".stat"); rtkclosestat(); rtkopenstat(statfile,sopt->sstat); }
* $POS,week,tow,stat,posx,posy,posz,posxf,posyf,poszf * week/tow : gps week no/time of week (s) * stat : solution status * posx/posy/posz : position x/y/z ecef (m) float * posxf/posyf/poszf : position x/y/z ecef (m) fixed * * $VELACC,week,tow,stat,vele,veln,velu,acce,accn,accu,velef,velnf,veluf,accef,accnf,accuf * week/tow : gps week no/time of week (s) * stat : solution status * vele/veln/velu : velocity e/n/u (m/s) float * acce/accn/accu : acceleration e/n/u (m/s^2) float * velef/velnf/veluf : velocity e/n/u (m/s) fixed * accef/accnf/accuf : acceleration e/n/u (m/s^2) fixed * * $CLK,week,tow,stat,clk1,clk2,clk3,clk4 * week/tow : gps week no/time of week (s) * stat : solution status * clk1 : receiver clock bias GPS (ns) * clk2 : receiver clock bias GLO-GPS (ns) * clk3 : receiver clock bias GAL-GPS (ns) * clk4 : receiver clock bias BDS-GPS (ns) * * $ION,week,tow,stat,sat,az,el,ion,ion-fixed * week/tow : gps week no/time of week (s) * stat : solution status * sat : satellite id * az/el : azimuth/elevation angle(deg) * ion : vertical ionospheric delay L1 (m) float * ion-fixed: vertical ionospheric delay L1 (m) fixed * * $TROP,week,tow,stat,rcv,ztd,ztdf * week/tow : gps week no/time of week (s) * stat : solution status * rcv : receiver (1:rover,2:base station) * ztd : zenith total delay (m) float * ztdf : zenith total delay (m) fixed * * $HWBIAS,week,tow,stat,frq,bias,biasf * week/tow : gps week no/time of week (s) * stat : solution status * frq : frequency (1:L1,2:L2,...) * bias : h/w bias coefficient (m/MHz) float * biasf : h/w bias coefficient (m/MHz) fixed * * $SAT,week,tow,sat,frq,az,el,resp,resc,vsat,snr,fix,slip,lock,outc,slipc,rejc * week/tow : gps week no/time of week (s) * sat/frq : satellite id/frequency (1:L1,2:L2,...) * az/el : azimuth/elevation angle (deg) * resp : pseudorange residual (m) * resc : carrier-phase residual (m) * vsat : valid data flag (0:invalid,1:valid) * snr : signal strength (dbHz) * fix : ambiguity flag (0:no data,1:float,2:fixed,3:hold,4:ppp) * slip : cycle-slip flag (bit1:slip,bit2:parity unknown) * lock : carrier-lock count * outc : data outage count * slipc : cycle-slip count * rejc : data reject (outlier) count * *-----------------------------------------------------------------------------*/
9、 /* write header to output file */
输出结果数据格式如下:
#define SOLF_LLH 0 /* solution format: lat/lon/height */ #define SOLF_XYZ 1 /* solution format: x/y/z-ecef */ #define SOLF_ENU 2 /* solution format: e/n/u-baseline */ #define SOLF_NMEA 3 /* solution format: NMEA-183 */
输出 结果文件头:
输出结果类型结构体:sol_t /* solution type */
typedef struct { /* solution type */ gtime_t time; /* time (GPST) */ double rr[6]; /* position/velocity (m|m/s) */ /* {x,y,z,vx,vy,vz} or {e,n,u,ve,vn,vu} */ float qr[6]; /* position variance/covariance (m^2) */ /* {c_xx,c_yy,c_zz,c_xy,c_yz,c_zx} or */ /* {c_ee,c_nn,c_uu,c_en,c_nu,c_ue} */ float qv[6]; /* velocity variance/covariance (m^2/s^2) */ double dtr[6]; /* receiver clock bias to time systems (s) */ unsigned char type; /* type (0:xyz-ecef,1:enu-baseline) */ unsigned char stat; /* solution status (SOLQ_???) */ unsigned char ns; /* number of valid satellites */ float age; /* age of differential (s) */ float ratio; /* AR ratio factor for valiation */ float thres; /* AR ratio threshold for valiation */ } sol_t;
RTK配置信息结构体:rtk_t; /* RTK control/result type */
typedef struct { /* RTK control/result type */ sol_t sol; /* RTK solution */ double rb[6]; /* base position/velocity (ecef) (m|m/s) */ int nx,na; /* number of float states/fixed states */ double tt; /* time difference between current and previous (s) */ double *x, *P; /* float states and their covariance */ double *xa,*Pa; /* fixed states and their covariance */ int nfix; /* number of continuous fixes of ambiguity */ ambc_t ambc[MAXSAT]; /* ambibuity control */ ssat_t ssat[MAXSAT]; /* satellite status */ int neb; /* bytes in error message buffer */ char errbuf[MAXERRMSG]; /* error message buffer */ prcopt_t opt; /* processing options */ } rtk_t;
ambc_t /* ambiguity control type */
typedef struct { /* ambiguity control type */ gtime_t epoch[4]; /* last epoch */ int n[4]; /* number of epochs */ double LC [4]; /* linear combination average */ double LCv[4]; /* linear combination variance */ int fixcnt; /* fix count */ char flags[MAXSAT]; /* fix flags */ } ambc_t;
卫星状态信息: ssat_t; /* satellite status type */
typedef struct { /* satellite status type */ unsigned char sys; /* navigation system */ unsigned char vs; /* valid satellite flag single */ double azel[2]; /* azimuth/elevation angles {az,el} (rad) */ double resp[NFREQ]; /* residuals of pseudorange (m) */ double resc[NFREQ]; /* residuals of carrier-phase (m) */ unsigned char vsat[NFREQ]; /* valid satellite flag */ unsigned char snr [NFREQ]; /* signal strength (0.25 dBHz) */ unsigned char fix [NFREQ]; /* ambiguity fix flag (1:fix,2:float,3:hold) */ unsigned char slip[NFREQ]; /* cycle-slip flag */ unsigned char half[NFREQ]; /* half-cycle valid flag */ int lock [NFREQ]; /* lock counter of phase */ unsigned int outc [NFREQ]; /* obs outage counter of phase */ unsigned int slipc[NFREQ]; /* cycle-slip counter */ unsigned int rejc [NFREQ]; /* reject counter */ double gf; /* geometry-free phase L1-L2 (m) */ double gf2; /* geometry-free phase L1-L5 (m) */ double mw; /* MW-LC (m) */ double phw; /* phase windup (cycle) */ gtime_t pt[2][NFREQ]; /* previous carrier-phase time */ double ph[2][NFREQ]; /* previous carrier-phase observable (cycle) */ } ssat_t;
11、单点定位&前向处理
if (popt_.mode==PMODE_SINGLE||popt_.soltype==0)
执行子函数:procpos
/* process positioning -------------------------------------------------------*/
static void procpos(FILE *fp, const prcopt_t *popt, const solopt_t *sopt,int mode)
3、procpos
rtkinit(&rtk,popt);
1、RTKINIT 初始化
/* initialize rtk control ------------------------------------------------------
* initialize rtk control struct
* args : rtk_t *rtk IO rtk control/result struct
* prcopt_t *opt I positioning options (see rtklib.h)
* return : none
*-----------------------------------------------------------------------------*/
extern void rtkinit(rtk_t *rtk, const prcopt_t *opt)
2、当前历元观测卫星多少、观测值;/* input obs data, navigation messages and sbas correction -------------------*/
static int inputobs(obsd_t *obs, int solq, const prcopt_t *popt)循环执行:每一个历元、每一个历元进行:
while ((nobs=inputobs(obs,rtk.sol.stat,popt))>=0) { /* exclude satellites */ for (i=n=0;i<nobs;i++) { if ((satsys(obs[i].sat,NULL)&popt->navsys)&& popt->exsats[obs[i].sat-1]!=1) obs[n++]=obs[i]; //去除配置中不参与的卫星的观测数据 } if (n<=0) continue; /* carrier-phase bias correction */ if (navs.nf>0) { //nf 为 number of satellite fcb data corr_phase_bias_fcb(obs,n,&navs); } else if (!strstr(popt->pppopt,"-DIS_FCB")) { corr_phase_bias_ssr(obs,n,&navs); } /* disable L2 */ #if 0 if (popt->freqopt==1) { for (i=0;i<n;i++) obs[i].L[1]=obs[i].P[1]=0.0; } #endif 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++; } }
定位、测速程序:
3.1 rtkpos
extern int rtkpos(rtk_t *rtk, const obsd_t *obs, int n, const nav_t *nav)
3.1.1 pntpos 单点定位
/* single-point positioning ----------------------------------------------------
* compute receiver position, velocity, clock bias by single-point positioning
* with pseudorange and doppler observablesextern int pntpos(const obsd_t *obs, int n, const nav_t *nav,const prcopt_t *opt, sol_t *sol, double *azel, ssat_t *ssat,char *msg)
- satposs 计算卫星位置和钟差
/* satellite positions and clocks --------compute satellite positions, velocities and clocks
extern void satposs(gtime_t teph, const obsd_t *obs, int n, const nav_t *nav,int ephopt, double *rs, double *dts, double *var, int *svh)
- ephclk // 利用广播星历计算卫星钟差
/* satellite clock with broadcast ephemeris ----------------------------------*/
static int ephclk(gtime_t time, gtime_t teph, int sat, const nav_t *nav,double *dts)
调用子函数:
1、首先选择与当前历元观测量对应的卫星星历:
/* select ephememeris --------------------------------------------------------*/
static eph_t *seleph(gtime_t time, int sat, int iode, const nav_t *nav)2、计算观测量对应卫星的钟差
/* broadcast ephemeris to satellite clock bias ---------------------------------
extern double eph2clk(gtime_t time, const eph_t *eph)
- satpos /* satellite position and clock at transmission time */
/* satellite position and clock -----------------* compute satellite position, velocity and clock
extern int satpos(gtime_t time, gtime_t teph, int sat, int ephopt,const nav_t *nav, double *rs, double *dts, double *var,
int *svh)注意:satpos中,根据系统对星历的配置,进行不同的卫星位置、速度求解!
#define EPHOPT_BRDC 0 /* ephemeris option: broadcast ephemeris */ #define EPHOPT_PREC 1 /* ephemeris option: precise ephemeris */ #define EPHOPT_SBAS 2 /* ephemeris option: broadcast + SBAS */ #define EPHOPT_SSRAPC 3 /* ephemeris option: broadcast + SSR_APC */ #define EPHOPT_SSRCOM 4 /* ephemeris option: broadcast + SSR_COM */ #define EPHOPT_LEX 5 /* ephemeris option: QZSS LEX ephemeris */
switch (ephopt) { case EPHOPT_BRDC : return ephpos (time,teph,sat,nav,-1,rs,dts,var,svh); case EPHOPT_SBAS : return satpos_sbas(time,teph,sat,nav, rs,dts,var,svh); case EPHOPT_SSRAPC: return satpos_ssr (time,teph,sat,nav, 0,rs,dts,var,svh); case EPHOPT_SSRCOM: return satpos_ssr (time,teph,sat,nav, 1,rs,dts,var,svh); case EPHOPT_PREC : if (!peph2pos(time,sat,nav,1,rs,dts,var)) break; else return 1; case EPHOPT_LEX : if (!lexeph2pos(time,sat,nav,rs,dts,var)) break; else return 1; }
调用子函数:
1、根据设置的星历类型,进行对应的卫星位置、钟差计算:
/* satellite position and clock by broadcast ephemeris -----------------------*/
static int ephpos(gtime_t time, gtime_t teph, int sat, const nav_t *nav,int iode, double *rs, double *dts, double *var, int *svh)1.1 选择观测量对应的卫星星历:
/* select ephememeris --------------------------------------------------------*/
static eph_t *seleph(gtime_t time, int sat, int iode, const nav_t *nav)2、计算卫星位置、速度、钟差
/* broadcast ephemeris to satellite position and clock bias --------------------
extern void eph2pos(gtime_t time, const eph_t *eph, double *rs, double *dts,double *var)
节对应的函数调用顺序如下:
- estpos 利用伪距计算接收机位置
/* estimate receiver position ------------------------------------------------*/
static int estpos(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)估计参数7个:3个位置、1个钟差、3个系统间钟差参数;
x[NX]={0},dx[NX],Q[NX*NX]
v=mat(n+4,1); H=mat(NX,n+4); var=mat(n+4,1);
下面进行循环迭代:求解接收机位置、速度;迭代终止条件为:(norm(dx,NX)<1E-4)
伪距残差/* pseudorange residuals -----------------------------------------------------*/
static int rescode(int iter, const obsd_t *obs, int n, const double *rs,
const double *dts, const double *vare, const int *svh,
const nav_t *nav, const double *x, const prcopt_t *opt,
double *v, double *H, double *var, double *azel, int *vsat,
double *resp, int *ns)
调用如下子函数:
1、计算卫星-接收机之间几何距离:
主要去除地球自转影响!
/* geometric distance ----------------------------------------------------------
* compute geometric distance and receiver-to-satellite unit vectorextern double geodist(const double *rs, const double *rr, double *e)
2、对伪距进行码偏差改正: 单频点定位考虑TGD
/* psendorange with code bias correction -------------------------------------*/
static double prange(const obsd_t *obs, const nav_t *nav, const double *azel,
int iter, const prcopt_t *opt, double *var)3、剔除预先设置中不参与定位的卫星
/* test excluded satellite -----------------------------------------------------
extern int satexclude(int sat, int svh, const prcopt_t *opt)
4、电离层改正:
/* ionospheric correction ------------------------------------------------------
* compute ionospheric correctionextern int ionocorr(gtime_t time, const nav_t *nav, int sat, const double *pos,const double *azel, int ionoopt, double *ion, double *var)
调用子函数:广播星历中电离层α、β模型:
5、对流层改正:
/* tropospheric correction -----------------------------------------------------
* compute tropospheric correctionextern int tropcorr(gtime_t time, const nav_t *nav, const double *pos, const double *azel, int tropopt, double *trp, double *var)
调用子函数:compute tropospheric delay by standard atmosphere and saastamoinen model
6、伪距残差:
/* pseudorange residual */ //nv表示方程的个数
v[nv]=P-(r+dtr-CLIGHT*dts[i*2]+dion+dtrp);7、设计/构造矩阵 H
/* design matrix */ //
for (j=0;j<NX;j++) H[j+nv*NX]=j<3?-e[j]:(j==3?1.0:0.0);
8、如果是多系统参与定位,需要估计系统间时间偏差
/* time system and receiver bias offset correction */
if (sys==SYS_GLO) {v[nv]-=x[4]; H[4+nv*NX]=1.0; mask[1]=1;}
else if (sys==SYS_GAL) {v[nv]-=x[5]; H[5+nv*NX]=1.0; mask[2]=1;}
else if (sys==SYS_CMP) {v[nv]-=x[6]; H[6+nv*NX]=1.0; mask[3]=1;}
else mask[0]=1;其具体参数为:
sol->dtr[0]=x[3]/CLIGHT; /* receiver clock bias (s) */ sol->dtr[1]=x[4]/CLIGHT; /* glo-gps time offset (s) */ sol->dtr[2]=x[5]/CLIGHT; /* gal-gps time offset (s) */ sol->dtr[3]=x[6]/CLIGHT; /* bds-gps time offset (s) */
9、伪距测量中各部分误差项中方差和
/* pseudorange measurement error variance ------------------------------------*/
static double varerr(const prcopt_t *opt, double el, int sys)函数代码:
/* pseudorange measurement error variance ------------------------------------*/ static double varerr(const prcopt_t *opt, double el, int sys) { double fact,varr; fact=sys==SYS_GLO?EFACT_GLO:(sys==SYS_SBS?EFACT_SBS:EFACT_GPS); varr=SQR(opt->err[0])*(SQR(opt->err[1])+SQR(opt->err[2])/sin(el)); if (opt->ionoopt==IONOOPT_IFLC) varr*=SQR(3.0); /* iono-free */ return SQR(fact)*varr; }
对应的项为:
结构体:prcopt_t
double err[5]; /* measurement error factor */
/* [0]:reserved */
/* [1-3]:error factor a/b/c of phase (m) */
/* [4]:doppler frequency (hz) */默认参数配置为:{100.0,0.003,0.003,0.0,1.0}, /* err[] */
为了保持“卫星个数”大于参数个数(7个),保持为超定方程:设置H矩阵及对应的方差,代码如下所示:
/* constraint to avoid rank-deficient */ for (i=0;i<4;i++) { if (mask[i]) continue; v[nv]=0.0; for (j = 0; j < NX; j++) { H[j + nv * NX] = j == i + 3 ? 1.0 : 0.0; printf("%f ", H[j + nv * NX]); if (j == NX - 1) printf("\n"); } var[nv++]=0.01; //printf("%f ",H[j + nv * NX]); }
某一个历元观测到卫星个数为5,故需要增加3个“伪卫星”,构成超定方程;代码运行如下:
对伪距残差进行定权,高度角定权:
/* weight by variance */ printf("伪距误差标准差:\n"); for (j=0;j<nv;j++) { sig=sqrt(var[j]); v[j]/=sig; printf("sig[%d]=%f%t ",j,sig); for (k=0;k<NX;k++) H[k+j*NX]/=sig; }
加权最小二乘:
/* least square estimation 用最小二乘法,解Ax=b-----------------------------------------------------
* least square estimation by solving normal equation (x=(A*A')^-1*A*y)* least square estimation by solving normal equation (x=(A*A')^-1*A*y)
* args : double *A I transpose of (weighted) design matrix (n x m)
* double *y I (weighted) measurements (m x 1)
* int n,m I number of parameters and measurements (n<=m) m为方程个数,测量个数
* double *x O estmated parameters (n x 1)
* double *Q O esimated parameters covariance matrix (n x n)
* return : status (0:ok,0>:error)
* notes : for weighted least square, replace A and y by A*w and w*y (w=W^(1/2))
* matirix stored by column-major order (fortran convention)
*-----------------------------------------------------------------------------*/保存结果:接收机位置、钟差、方差、协方差、参与定位卫星个数:
if (norm(dx,NX)<1E-4) { sol->type=0; sol->time=timeadd(obs[0].time,-x[3]/CLIGHT); //求解位置的时刻为:扣除接收机钟差影响 sol->dtr[0]=x[3]/CLIGHT; /* receiver clock bias (s) */ sol->dtr[1]=x[4]/CLIGHT; /* glo-gps time offset (s) */ sol->dtr[2]=x[5]/CLIGHT; /* gal-gps time offset (s) */ sol->dtr[3]=x[6]/CLIGHT; /* bds-gps time offset (s) */ for (j=0;j<6;j++) sol->rr[j]=j<3?x[j]:0.0; for (j=0;j<3;j++) sol->qr[j]=(float)Q[j+j*NX]; sol->qr[3]=(float)Q[1]; /* cov xy */ sol->qr[4]=(float)Q[2+NX]; /* cov yz */ sol->qr[5]=(float)Q[2]; /* cov zx */ sol->ns=(unsigned char)ns; sol->age=sol->ratio=0.0; /* validate solution */ if ((stat=valsol(azel,vsat,n,opt,v,nv,NX,msg))) { sol->stat=opt->sateph==EPHOPT_SBAS?SOLQ_SBAS:SOLQ_SINGLE; }
同时还包括:有效性判断:
1、chisqr
2、最大DOP值
/* validate solution ---------------------------------------------------------*/
static int valsol(const double *azel, const int *vsat, int n,
const prcopt_t *opt, const double *v, int nv, int nx,
char *msg)estpos 对应的函数调用层级如下:
- raim_fde //有效性判断,完好性 !
进入raim_fde前进行的判断:
/* 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);}1、参与定位卫星个数6个及以上;
2、opt->posopt[4] 为 1 ?表示啥意思?
/* 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)
- estvel 计算接收机速度:
/* estimate receiver velocity ------------------------------------------------*/
static void estvel(const obsd_t *obs, int n, const double *rs, const double *dts,
const nav_t *nav, const prcopt_t *opt, sol_t *sol,
const double *azel, const int *vsat)调用子函数:
1、多普勒残差:
/* doppler residuals ------------------------*/
static int resdop(const obsd_t *obs, int n, const double *rs, const double *dts,const nav_t *nav, const double *rr, const double *x,
const double *azel, const int *vsat, double *v, double *H)2、最小二乘计算速度:
ssat 保存相关信息:ssat_t ssat;
typedef struct { /* satellite status type */ unsigned char sys; /* navigation system */ unsigned char vs; /* valid satellite flag single */ double azel[2]; /* azimuth/elevation angles {az,el} (rad) */ double resp[NFREQ]; /* residuals of pseudorange (m) */ double resc[NFREQ]; /* residuals of carrier-phase (m) */ unsigned char vsat[NFREQ]; /* valid satellite flag */ unsigned char snr [NFREQ]; /* signal strength (0.25 dBHz) */ unsigned char fix [NFREQ]; /* ambiguity fix flag (1:fix,2:float,3:hold) */ unsigned char slip[NFREQ]; /* cycle-slip flag */ unsigned char half[NFREQ]; /* half-cycle valid flag */ int lock [NFREQ]; /* lock counter of phase */ unsigned int outc [NFREQ]; /* obs outage counter of phase */ unsigned int slipc[NFREQ]; /* cycle-slip counter */ unsigned int rejc [NFREQ]; /* reject counter */ double gf; /* geometry-free phase L1-L2 (m) */ double gf2; /* geometry-free phase L1-L5 (m) */ double mw; /* MW-LC (m) */ double phw; /* phase windup (cycle) */ gtime_t pt[2][NFREQ]; /* previous carrier-phase time */ double ph[2][NFREQ]; /* previous carrier-phase observable (cycle) */ } ssat_t;
outsolstat
/* output solution status ----------------------------------------------------*/
static void outsolstat(rtk_t *rtk)调用子函数:
rtkoutstat
/* write solution status to buffer -------------------------------------------*/
extern int rtkoutstat(rtk_t *rtk, char *buff)保存的为如下信息到文件:myposbyVS.pos.stat
/* write ppp solution status to buffer */ if (rtk->opt.mode>=PMODE_PPP_KINEMA) { return pppoutstat(rtk,buff); } est=rtk->opt.mode>=PMODE_DGPS; nfreq=est?nf:1; tow=time2gpst(rtk->sol.time,&week); /* receiver position */ if (est) { for (i=0;i<3;i++) xa[i]=i<rtk->na?rtk->xa[i]:0.0; p+=sprintf(p,"$POS,%d,%.3f,%d,%.4f,%.4f,%.4f,%.4f,%.4f,%.4f\n",week,tow, rtk->sol.stat,rtk->x[0],rtk->x[1],rtk->x[2],xa[0],xa[1], xa[2]); } else { p+=sprintf(p,"$POS,%d,%.3f,%d,%.4f,%.4f,%.4f,%.4f,%.4f,%.4f\n",week,tow, rtk->sol.stat,rtk->sol.rr[0],rtk->sol.rr[1],rtk->sol.rr[2], 0.0,0.0,0.0); } /* receiver velocity and acceleration */ if (est&&rtk->opt.dynamics) { ecef2pos(rtk->sol.rr,pos); ecef2enu(pos,rtk->x+3,vel); ecef2enu(pos,rtk->x+6,acc); if (rtk->na>=6) ecef2enu(pos,rtk->xa+3,vela); if (rtk->na>=9) ecef2enu(pos,rtk->xa+6,acca); p+=sprintf(p,"$VELACC,%d,%.3f,%d,%.4f,%.4f,%.4f,%.5f,%.5f,%.5f,%.4f,%.4f,%.4f,%.5f,%.5f,%.5f\n", week,tow,rtk->sol.stat,vel[0],vel[1],vel[2],acc[0],acc[1], acc[2],vela[0],vela[1],vela[2],acca[0],acca[1],acca[2]); } else { ecef2pos(rtk->sol.rr,pos); ecef2enu(pos,rtk->sol.rr+3,vel); p+=sprintf(p,"$VELACC,%d,%.3f,%d,%.4f,%.4f,%.4f,%.5f,%.5f,%.5f,%.4f,%.4f,%.4f,%.5f,%.5f,%.5f\n", week,tow,rtk->sol.stat,vel[0],vel[1],vel[2], 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0); } /* receiver clocks */ p+=sprintf(p,"$CLK,%d,%.3f,%d,%d,%.3f,%.3f,%.3f,%.3f\n", week,tow,rtk->sol.stat,1,rtk->sol.dtr[0]*1E9,rtk->sol.dtr[1]*1E9, rtk->sol.dtr[2]*1E9,rtk->sol.dtr[3]*1E9); /* ionospheric parameters */ if (est&&rtk->opt.ionoopt==IONOOPT_EST) { for (i=0;i<MAXSAT;i++) { ssat=rtk->ssat+i; if (!ssat->vs) continue; satno2id(i+1,id); j=II(i+1,&rtk->opt); xa[0]=j<rtk->na?rtk->xa[j]:0.0; p+=sprintf(p,"$ION,%d,%.3f,%d,%s,%.1f,%.1f,%.4f,%.4f\n",week,tow, rtk->sol.stat,id,ssat->azel[0]*R2D,ssat->azel[1]*R2D, rtk->x[j],xa[0]); } } /* tropospheric parameters */ if (est&&(rtk->opt.tropopt==TROPOPT_EST||rtk->opt.tropopt==TROPOPT_ESTG)) { for (i=0;i<2;i++) { j=IT(i,&rtk->opt); xa[0]=j<rtk->na?rtk->xa[j]:0.0; p+=sprintf(p,"$TROP,%d,%.3f,%d,%d,%.4f,%.4f\n",week,tow, rtk->sol.stat,i+1,rtk->x[j],xa[0]); } } /* receiver h/w bias */ if (est&&rtk->opt.glomodear==2) { for (i=0;i<nfreq;i++) { j=IL(i,&rtk->opt); xa[0]=j<rtk->na?rtk->xa[j]:0.0; p+=sprintf(p,"$HWBIAS,%d,%.3f,%d,%d,%.4f,%.4f\n",week,tow, rtk->sol.stat,i+1,rtk->x[j],xa[0]); } }
- outsol
/* output solution body --------------------------------------------------------
* output solution body to file
* args : FILE *fp I output file pointer
* sol_t *sol I solution
* double *rb I base station position {x,y,z} (ecef) (m)
* solopt_t *opt I solution options
* return : none
*-----------------------------------------------------------------------------*/
extern void outsol(FILE *fp, const sol_t *sol, const double *rb, const solopt_t *opt)
outsols
/* output solution body --------------------------------------------------------
* output solution body to buffer
* args : unsigned char *buff IO output buffer
* sol_t *sol I solution
* double *rb I base station position {x,y,z} (ecef) (m)
* solopt_t *opt I solution options
* return : number of output bytes
*-----------------------------------------------------------------------------*/
extern int outsols(unsigned char *buff, const sol_t *sol, const double *rb,const solopt_t *opt)
判断语句:根据设置的输出位置形式,进行选择:switch (opt->posf) { case SOLF_LLH: p+=outpos (p,s,sol,opt); break; case SOLF_XYZ: p+=outecef(p,s,sol,opt); break; case SOLF_ENU: p+=outenu(p,s,sol,rb,opt); break; case SOLF_NMEA: p+=outnmea_rmc(p,sol); p+=outnmea_gga(p,sol); break; }
/* output solution as the form of lat/lon/height -----------------------------*/
static int outpos(unsigned char *buff, const char *s, const sol_t *sol,const solopt_t *opt)2.2 节函数调用层级,如下:
GPS伪距单频/L1 单点定位整个过程包括:
1、读取移动站观测数据、导航数据;将整个文件的相关数据放入结构体中;进行排序和重复数据删除;
2、进行单点定位,测速;保存数据.pos、.pos.trace、.pos.stat文件;
3、最终freeobsnav(&obss,&navs);