一、调试运行RTKLIB
话不多说,开启rtklib学习第一步,放上链接:
测绘工程本科生如何入门GNSS算法 (一)- 在Visual Studio 2017中跑通RTKLIB开源代码 - 知乎 (zhihu.com)
RTKlib代码调试教程-Visual Studio 2017_rtklib调试无法打开源文件-CSDN博客
RTKLIB源码——如何在VS2017中配置、调试_rtklib在vs中运行运行-CSDN博客
这些博主的方法最后得到一个.exe命令行程序,对于学习rtklib来说远远不够,为此需要改写main文件代码,使其能够在代码中读取文件,调试运行。
二、对代码的先行知识
在main.c文件的main函数中,上图的postpos函数是通往各类定位功能的函数入口,因此,改写程序,必须要了解此函数的参数的定义与作用。
1、rtklib.h
rtklib.h存储了各种常量、结构体和函数原型的定义
先在main.c文件找到上图代码,确定&prcopt,&solopt,&filopt结构体的类型,在到rtklib.h文件中找到定义,如下:
typedef struct { /* processing options type */
int mode; /* positioning mode (PMODE_???) */
int soltype; /* solution type (0:forward,1:backward,2:combined) */
int nf; /* number of frequencies (1:L1,2:L1+L2,3:L1+L2+L5) */
int navsys; /* navigation system */
double elmin; /* elevation mask angle (rad) */
snrmask_t snrmask; /* SNR mask */
int sateph; /* satellite ephemeris/clock (EPHOPT_???) */
int modear; /* AR mode (0:off,1:continuous,2:instantaneous,3:fix and hold,4:ppp-ar) */
int glomodear; /* GLONASS AR mode (0:off,1:on,2:auto cal,3:ext cal) */
int bdsmodear; /* BeiDou AR mode (0:off,1:on) */
int maxout; /* obs outage count to reset bias */
int minlock; /* min lock count to fix ambiguity */
int minfix; /* min fix count to hold ambiguity */
int ionoopt; /* ionosphere option (IONOOPT_???) */
int tropopt; /* troposphere option (TROPOPT_???) */
int dynamics; /* dynamics model (0:none,1:velociy,2:accel) */
int tidecorr; /* earth tide correction (0:off,1:solid,2:solid+otl+pole) */
int niter; /* number of filter iteration */
int codesmooth; /* code smoothing window size (0:none) */
int intpref; /* interpolate reference obs (for post mission) */
int sbascorr; /* SBAS correction options */
int sbassatsel; /* SBAS satellite selection (0:all) */
int rovpos; /* rover position for fixed mode */
int refpos; /* base position for relative mode */
/* (0:pos in prcopt, 1:average of single pos, */
/* 2:read from file, 3:rinex header, 4:rtcm pos) */
double eratio[NFREQ]; /* code/phase error ratio */
double err[5]; /* measurement error factor */
/* [0]:reserved */
/* [1-3]:error factor a/b/c of phase (m) */
/* [4]:doppler frequency (hz) */
double std[3]; /* initial-state std [0]bias,[1]iono [2]trop */
double prn[5]; /* process-noise std [0]bias,[1]iono [2]trop [3]acch [4]accv */
double sclkstab; /* satellite clock stability (sec/sec) */
double thresar[4]; /* AR validation threshold */
double elmaskar; /* elevation mask of AR for rising satellite (deg) */
double elmaskhold; /* elevation mask to hold ambiguity (deg) */
double thresslip; /* slip threshold of geometry-free phase (m) */
double maxtdiff; /* max difference of time (sec) */
double maxinno; /* reject threshold of innovation (m) */
double maxgdop; /* reject threshold of gdop */
double baseline[2]; /* baseline length constraint {const,sigma} (m) */
double ru[3]; /* rover position for fixed mode {x,y,z} (ecef) (m) */
double rb[3]; /* base position for relative mode {x,y,z} (ecef) (m) */
char anttype[2][MAXANT]; /* antenna types {rover,base} */
double antdel[2][3]; /* antenna delta {{rov_e,rov_n,rov_u},{ref_e,ref_n,ref_u}} */
pcv_t pcvr[2]; /* receiver antenna parameters {rov,base} */
unsigned char exsats[MAXSAT]; /* excluded satellites (1:excluded,2:included) */
char rnxopt[2][256]; /* rinex options {rover,base} */
int posopt[6]; /* positioning options */
int syncsol; /* solution sync mode (0:off,1:on) */
double odisp[2][6*11]; /* ocean tide loading parameters {rov,base} */
exterr_t exterr; /* extended receiver error model */
} prcopt_t;
typedef struct { /* solution options type */
int posf; /* solution format (SOLF_???) */
int times; /* time system (TIMES_???) */
int timef; /* time format (0:sssss.s,1:yyyy/mm/dd hh:mm:ss.s) */
int timeu; /* time digits under decimal point */
int degf; /* latitude/longitude format (0:ddd.ddd,1:ddd mm ss) */
int outhead; /* output header (0:no,1:yes) */
int outopt; /* output processing options (0:no,1:yes) */
int datum; /* datum (0:WGS84,1:Tokyo) */
int height; /* height (0:ellipsoidal,1:geodetic) */
int geoid; /* geoid model (0:EGM96,1:JGD2000) */
int solstatic; /* solution of static mode (0:all,1:single) */
int sstat; /* solution statistics level (0:off,1:states,2:residuals) */
int trace; /* debug trace level (0:off,1-5:debug) */
double nmeaintv[2]; /* nmea output interval (s) (<0:no,0:all) */
/* nmeaintv[0]:gprmc,gpgga,nmeaintv[1]:gpgsv */
char sep[64]; /* field separator */
char prog[64]; /* program name */
} solopt_t;
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;
2、inputfile,outfile
inputfile包含了.o,.n,.clk,.sp3文件
outfile文件则为程序输出文件.pos
在postpos.c文件的postpos函数中也可以获取这些参数的信息
了解完这些参数,就可以动手去修改main.c文件了。
三、改写代码示例
1、spp运行
参考链接:RTKLIB专题学习(四)---单点定位实现初识(二)_rtklib processing unit time_十八与她的博客-CSDN博客
此时只需要在
①filopt filopt结构体这设置对应参数文件(文件是有顺序的)
②char infile_[MAXFILE][MAXSTRPATH]处设置观测值文件和星历文件
③double eps[],double epe[]处设置起始历元和结束历元时间
输出文件如下:
2、PPP代码简单但输入文件麻烦版本
main.c函数代码给出(同SPP的代码区别不大,就是多了参数文件的路径):
/*------------------------------------------------------------------------------
* rnx2rtkp.c : read rinex obs/nav files and compute receiver positions
*
* Copyright (C) 2007-2009 by T.TAKASU, All rights reserved.
*
* version : $Revision: 1.1 $ $Date: 2008/07/17 21:55:16 $
* history : 2007/01/16 1.0 new
* 2007/03/15 1.1 add library mode
* 2007/05/08 1.2 separate from postpos.c
* 2009/01/20 1.3 support rtklib 2.2.0 api
* 2009/12/12 1.4 support glonass
* add option -h, -a, -l, -x
* 2010/01/28 1.5 add option -k
* 2010/08/12 1.6 add option -y implementation (2.4.0_p1)
* 2014/01/27 1.7 fix bug on default output time format
*-----------------------------------------------------------------------------*/
#include <stdarg.h>
#include "rtklib.h"
#pragma once
static const char rcsid[] = "$Id: rnx2rtkp.c,v 1.1 2008/07/17 21:55:16 ttaka Exp $";
#define PROGNAME "rnx2rtkp" /* program name */
#define MAXFILE 8 /* max number of input files */
/* show message --------------------------------------------------------------*/
static int showmsg(char *format, ...)
{
va_list arg;
va_start(arg,format); vfprintf(stderr,format,arg); va_end(arg);
fprintf(stderr,"\r");
return 0;
}
static void settspan(gtime_t ts, gtime_t te) {}
static void settime(gtime_t time) {}
/* print help ----------------------------------------------------------------*/
static void printhelp(void)
{
int i;
for (i=0;i<(int)(sizeof(help)/sizeof(*help));i++) fprintf(stderr,"%s\n",help[i]);
exit(0);
}
/* rnx2rtkp main -------------------------------------------------------------*/
int main(int argc, char **argv)
{
int i, n, ret, x;
double tint = 0; /* 求解时间间隔(0:默认) */
gtime_t ts = { 0 }, te = { 0 }; /* 历元时段始末控制变量 */
char* infile[MAXFILE], outfile[MAXSTRPATH] = { '\0' };
/* 结果输出路径 */
char resultpath[MAXSTRPATH] = "D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\result";
char sep = (char)FILEPATHSEP;
solopt_t solopt = solopt_default; /* 默认求解格式设置 */
prcopt_t prcopt = prcopt_default;
sol_t solt;
filopt_t filopt = { /* 参数文件路径设置 */
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\igs14.atx", /* 卫星天线参数文件 */
"", /* 接收机天线参数文件 */
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\site.crd", /* 测站位置文件 */
"", /* 扩展大地水准面数据文件 */
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\CODG2440.17I", /* 电离层数据文件 */
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\CAS0MGXRAP_20172440000_01D_01D_DCB.BSX",/* DCB数据文件 */
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\igs19647.erp"/* 地球自转参数文件 */
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\ocnload.blq", /* 海洋潮汐负荷文件 */
};
//观测数据文件
char infile_[MAXFILE][MAXSTRPATH] = {
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\brdm2440.17p",
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\jfng2440.17o",
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\wum19644.clk",
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\wum19644.sp3",
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\wum19645.clk",
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\wum19645.sp3",
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\wum19646.clk",
"D:\\visual_studio2017\\project\\rtklib_ppp\\rtklib\\testdata\\test1\\wum19646.sp3",
};
long t1, t2;
double eps[] = { 2017,9,1,0,0,0 }, epe[] = { 2017,9,2,0,0,0 }; /* 设置计算的历元时段 */
ts = epoch2time(eps); te = epoch2time(epe);
for (i = 0, n = 0; i < MAXFILE; i++)
if (strcmp(infile_[i], "")) infile[n++] = &infile_[i][0];
sprintf(outfile, "%s%c", resultpath, sep);//设置输出路径
/* 自定义求解格式 --------------------------------------------------------*/
solopt.posf = SOLF_XYZ; /* 选择输出的坐标格式,经纬度或是XYZ坐标等 */
solopt.times = TIMES_GPST; /* 控制输出解的时间系统类型 */
solopt.degf = 0; /* 输出经纬度格式(0:°, 1:°′″) */
solopt.outhead = 1; /* 是否输出头文件(0:否,1:是) */
solopt.outopt = 1; /* 是否输出prcopt变量(0:否,1:是) */
solopt.height = 1; /* 高程(0:椭球高,1:大地高) */
solopt.sstat = 2;
/* 自定义处理选项设置-PPP ----------------------------------------------------*/
prcopt.navsys = SYS_ALL;// GPS+GAL //SYS_ALL;
prcopt.mode = PMODE_PPP_STATIC;
prcopt.soltype = 0; // 求解类型(0:向前滤波,1:向后滤波,2:混合滤波)
prcopt.nf=2; /* 参与计算的载波频率个数 (1:L1,2:L1+L2,3:L1+L2+L5) */
prcopt.navsys = SYS_GPS | SYS_CMP; // 处理的导航系统
prcopt.elmin = 10.0 * D2R;//卫星截止高度角
prcopt.sateph = EPHOPT_PREC;
prcopt.modear=1; /* AR mode (0:off,1:continuous,2:instantaneous,3:fix and hold,4:ppp-ar) */
prcopt.elmin = 7 * D2R;
prcopt.ionoopt = IONOOPT_TEC;
prcopt.tropopt = TROPOPT_ESTG;
prcopt.tidecorr = 2; //地球潮汐改正选项(0:关闭,1:固体潮,2:固体潮+?+极移)
sprintf(outfile, "%s%cAJAC-ppp-static.pos", resultpath, sep); // 输出结果名称
/*rtk*/
solt.stat = SOLQ_SINGLE;
t1 = clock();
ret = postpos(ts, te, tint, 0.0, &prcopt, &solopt, &filopt, infile, n, outfile, "", "");
t2 = clock();
if (!ret) fprintf(stderr, "%40s\r", "");
printf("\n * The total time for running the program: %6.3f seconds\n%c", (double)(t2 - t1) / CLOCKS_PER_SEC, '\0');
printf("Press any key to exit!\n");
getchar();
return ret;
}
与spp一样的思路,但是当文件过多时 ,一个个输入很繁琐。那就看下面的!
3、PPP代码升级但输入文件简单版本
在一位研究生学长指导下完成,对他和CSDN博主十八与她表达感谢!
思路是手动输入一个目录,其中有ppp定位所需文件和一个子目录result存放输出文件。
因为在读取目录内文件时涉及到entry入口,因此将main.c改为main.cc,并添加头文件。
此时只需在代码中更改文件夹路径和历元时间及MAXFILE宏定义即可
//rnx2rtkp.c : read rinex obs/nav files and compute receiver positions
#include <stdarg.h>
#include "rtklib.h"
#include<stdio.h>
#include<assert.h>
#include <iostream>
#include <thread>
#include "dirent.h"
#include<string.h>
#include<stdlib.h>
#pragma once
static const char rcsid[] = "$Id: rnx2rtkp.c,v 1.1 2008/07/17 21:55:16 ttaka Exp $";
#define PROGNAME "rnx2rtkp" /* program name */
#define MAXFILE 14 /* max number of input files */
/* show message --------------------------------------------------------------*/
static int showmsg(char *format, ...);
static void settspan(gtime_t ts, gtime_t te) {}
static void settime(gtime_t time) {}
// read_filedir -----------------------------------------------------
//读取filepath目录中的文件
static void read_filedir(char* filePath, char file[MAXFILE][MAXSTRPATH],char * resultpath);
//substrend -----------------------------------------------------
//获取字符串的后n位
static char* substrend(const char* str, int n);
/* rnx2rtkp main -------------------------------------------------------------*/
int main(int argc, char **argv)
{
int i, n, ret, x;
double tint = 0; /* */
gtime_t ts = { 0 }, te = { 0 }; /* */
char* infile[MAXFILE], outfile[MAXSTRPATH] = { '\0' };
char resultpath[MAXSTRPATH] = {'\0'};/*the output file's path */
//input dir path
char filepath[MAXSTRPATH] = "D:\\visual_studio2017\\project\\rtklib_ppp\\auto_rtklib\\testdata\\test1";
double eps[] = { 2017,9,1,0,0,0 }, epe[] = { 2017,9,2,0,0,0 }; /*the start time and end time */
ts = epoch2time(eps); te = epoch2time(epe);
char sep = (char)FILEPATHSEP;
solopt_t solopt = solopt_default;
prcopt_t prcopt = prcopt_default;
sol_t solt;
char *infile_[MAXFILE];
filopt_t filopt;
//auto acquire infile_ and filopt
n = 0;
char file[MAXFILE][MAXSTRPATH] = { "" };//all files from input dir
read_filedir(filepath, file,resultpath);
if (resultpath[0] == '\0') strcpy(resultpath,filepath);//无输出文件夹时在当前目录输出
for (int i = 0; i < MAXFILE; i++)
{
char* substr = substrend(*(file + i), 4);
char* substr2 = substrend(*(file + i), 1);
if (strcmp(substr, ".atx") == 0 || strcmp(substr, ".ATX") == 0)
{
sprintf(filopt.satantp, *(file + i)); /* satellite antenna parameters file */
//sprintf(filopt.rcvantp, *(file + i)); /* receiver antenna parameters file */
continue;
}
else if (strcmp(substr, ".snx") == 0 || strcmp(substr, ".crd") == 0 || strcmp(substr, ".SNX") == 0) { //
sprintf(filopt.stapos, *(file + i)); /* station positions file */
continue;
}
else if (strcmp(substr, ".BSX") == 0 || strcmp(substr, ".bsx") == 0 || strcmp(substr, ".DCB") == 0) {//----------------
sprintf(filopt.dcb, *(file + i)); /* dcb data file */
continue;
}
else if (strcmp(substr2, "I") == 0 || strcmp(substr2, "i") == 0) {//----------------
sprintf(filopt.iono, *(file + i)); /* ION data file */
continue;
}
else if (strcmp(substr, ".erp") == 0 || strcmp(substr, ".ERP") == 0) {//--------------------------
sprintf(filopt.eop, *(file + i)); /* eop data file */
continue;
}
else if (strcmp(substr, ".blq") == 0) { //----------------------------
sprintf(filopt.blq, *(file + i)); /* ocean tide loading blq file */
continue;
}
else {
//infile obs,nav,clk and sp3
infile_[n] = *(file + i);
n++;
}
}
long t1, t2;
for (i = 0; i < n; i++) {
if (strcmp(infile_[i], "")) infile[i] = infile_[i];
printf("infile[%d]:%s\n", i, infile[i]);
}
sprintf(outfile, "%s%c", resultpath, sep);//
/* --------------------------------------------------------*/
solopt.posf = SOLF_XYZ; /* */
solopt.times = TIMES_GPST; /* */
solopt.degf = 0; /* */
solopt.outhead = 1; /* */
solopt.outopt = 1; /* */
solopt.height = 1; /* */
solopt.sstat = 2;
/* the prcopt of PPP ----------------------------------------------------*/
prcopt.navsys = SYS_ALL;// GPS+GAL //SYS_ALL;
prcopt.mode = PMODE_PPP_STATIC;//rtklib.h/319line
prcopt.soltype = 0; //
prcopt.nf=2;
prcopt.navsys = SYS_GPS | SYS_CMP; //
prcopt.elmin = 10.0 * D2R;//
prcopt.sateph = EPHOPT_PREC;
prcopt.modear=1; /* AR mode (0:off,1:continuous,2:instantaneous,3:fix and hold,4:ppp-ar) */
prcopt.elmin = 7 * D2R;
prcopt.ionoopt = IONOOPT_TEC;
prcopt.tropopt = TROPOPT_ESTG;
prcopt.tidecorr = 2; //
sprintf(outfile, "%s%cAJAC-ppp-static.pos", resultpath, sep); // the outfile name
/*rtk*/
solt.stat = SOLQ_SINGLE;
t1 = clock();
ret = postpos(ts, te, tint, 0.0, &prcopt, &solopt, &filopt, infile, n, outfile, "", "");
t2 = clock();
if (!ret) fprintf(stderr, "%40s\r", "");
printf("\n * The total time for running the program: %6.3f seconds\n%c", (double)(t2 - t1) / CLOCKS_PER_SEC, '\0');
printf("Press any key to exit!\n");
getchar();
return ret;
}
static char* substrend(const char* str, int n)
{
int k = 0;
char* substr = (char*)malloc(n + 1);
//char a = *str;
int length = strlen(str);
if (n >= length)
{
//*substr = a;
strcpy(substr, str);
//printf("substr:%s", substr);
return substr;
}
for (int i = length - n; i < length; i++)
{
substr[k] = str[i];
k++;
}
substr[k] = '\0';
//printf("substr:%s\n", substr);
return substr;
}
static void read_filedir(char* filePath, char file[MAXFILE][MAXSTRPATH],char *resultpath)
{
int i = 0, j = 0;
DIR* dir = NULL;
struct dirent* entry;
char temp[MAXSTRPATH];
if ((dir = opendir(filePath)) == NULL)
{
printf("opendir failed!");
}
else
{
while (entry = readdir(dir))
{
if (strcmp(entry->d_name, ".") == 0 || strcmp(entry->d_name, "..") == 0 || !strchr(entry->d_name, '.')) {
if (!strchr(entry->d_name, '.')) {//找到result路径
temp[0] = { '\0' };
strcat(temp, filePath); strcat(temp, "\\");
strcat(temp, entry->d_name); strcpy(resultpath, temp);
}
continue;
}
strcpy(temp, filePath);
strcat(temp, "\\");
strcat(temp, entry->d_name);
strcpy(file[i], temp);
i++;
printf("filename%d = %s\n", i, entry->d_name);
}
}
for (i = 0; i < MAXFILE; i++) printf("filename%d = %s\n", i, file[i]);
}
static const char *help[] = {""};
static int showmsg(char *format, ...)
{
va_list arg;
va_start(arg, format); vfprintf(stderr, format, arg); va_end(arg);
fprintf(stderr, "\r");
return 0;
}
static void printhelp(void)
{
int i;
for (i = 0; i < (int)(sizeof(help) / sizeof(*help)); i++) fprintf(stderr, "%s\n", help[i]);
exit(0);
}
至此开启rtklib学习新篇章
四、 致读者与自己
作为一名测绘工程本科生,尤其在前不久的大二暑期实习结束后,我深刻体验了测绘工程传统作业的艰辛。因此我坚定了要快速掌握测绘3S技术的想法,真的不想在接触外业了。在暑期我开始入手rtklib,从看代码和看论文两手齐抓,起初进度很慢,无论是晦涩难懂的代码还是论文都犹如天书很难进行下去。但一天天积累一点点,终于在国庆假期,我给自己安排一段时间去攻克代码调试难关,尤其这几天,彻底解决了代码调试问题,附加的对rtklib的结构,尤其是main函数里的各种参数都可以算的是很了解的了。
再次郑重感谢以上博主和我的研究生师兄,他们在我舍不得开CSDN会员的情况下,在一个个难关面前给我提供供了宝贵的经验与灵感。学习就是在巨人的肩膀上到下一个巨人的肩膀上,为此我坚定的免费开源CSDN上所有文章。并希望读者在评论区一同分享自己觉得给自己帮助很大的链接。
愿大家一起努力,尤其在GNSS领域用知识武装自己,成为高级人才,建设伟大中国!