这个学期一直在做卫星导航算法重构,重写一套算法过程中借鉴学习了各类开源程序,对目前用过的程序进行总结。
1、开源程序总结
卫星导航算法整体架构推荐阅读:上交GICI、与武大中心的POSGO ,架构合理优秀,他们的数据结构会更适用于多频多系统。强烈建议卫星导航的一些数据结构一定要按照C++ map容器设计,卫星号为健,对应相应的卫星信息。摆脱传统的数组或者结构体,不然在设计多频多系统时会吃亏。 (后面会单独详细总结)。
多频多系统程序设计推荐阅读:武大中心PRIDE-PPP-AR,适应了所有频点,模糊度固定,功能强大。我们一般可以通过卫星位置与sunpos信息,大致确定卫星的姿态,并且对于GPS与GLONASS 的eclip 现象可进一步精确计算卫星的姿态。但是在此程序中还添加了att.obx文件,我们能获取高精度卫星姿态信息,在进行windup修正时,会有进一步的提升,这样类似的提升在此开源程序还有很多很多
2、卫星导航程序算法数据结构总结
学会使用C++ IO流
如果要重构一份卫星导航程序,个人建议摆脱Rtklib的C语言文件读取,随着文件版本的迭代、例如CLK文件3,04等文件版本更新,出现了许多空格等因素造成的文件读取失败,建议用C++ IO流方式进行文件读取,方便高效快速。
观测文件数据设计
利用C++ STL库设计观测数据,以下是一个简单的例子
struct Carrierphase {
double value;//载波相位值
uint8_t LLi; //周跳标识符
int snr; //载噪比
double freq;
string code;
};
//伪距
struct Pseudorange{
double freq;
double value;
string code;
};
struct MeasData
{
map<int,Carrierphase> L;
map<int,Pseudorange> P;
map<int,double> D;
map<int,double> S;
};
//每一个历元的存储的观测数据(根据不同的系统分开存储,为了后续可以快速修改)
struct MeasDatas {
map<int,MeasData> GpsObservationData;
map<int,MeasData> GlonassObservationData;
map<int,MeasData> GalieoObservationData;
map<int,MeasData> Bds2ObservationData;
map<int,MeasData> Bds3ObservationData;
//对应的prn,每颗卫星存储的观测值,因为不存在同一个历元出现相同的卫星好的数据,所以用map存储
};
后面发现在GICI、POSGO 开源代码都有类似的数据结构,(英雄所见略同) 并且由于map是二叉树,查找真的很方便。
使用多线程多进程读取多个文件、设计程序
1、以下是我文件读取的程序,简单利用了C++17的特性 future。2、并且我在设计程序中,打算对每一种卫星导航单独开进程,四大卫星导航在计算卫星位置、对流层改正量、电离层改正量等....其他信息都是分开的。3、我们在处理同一个历元相同卫星也可以采取多进程的方式,因为每一颗卫星计算其改正值方式都应该是相同的,例如算GPS 1号卫星的方式与算GPS 3号卫星是可以同时进行的,整体上北斗导航系统与GPS导航系统也是可以同时进行的,利用C++ 多线程多进程的特性,提高代码的效率。
// 定义互斥锁
std::mutex obs_mtx;
std::mutex dcb_mtx;
std::mutex osb_mtx;
std::mutex atx_mtx;
std::mutex ephemeris_mtx;
std::map<std::string, std::optional<std::future<bool>>> futures;
for (const auto &entry : file_map) {
const auto &file_type = entry.first;
const auto &file_name = entry.second;
if (file_type == OBSFILE) {
futures[file_type] = std::async(std::launch::async, [&]() {
std::lock_guard<std::mutex> lock(obs_mtx);
return obs_data_.loadobservation(file_name);
});
} else if (file_type == EPHEMERISFILE) {
futures[file_type] = std::async(std::launch::async, [&]() {
std::lock_guard<std::mutex> lock(ephemeris_mtx);
return ephemeris_data_.loadEphemeris(file_name);
});
} else if (file_type == OSBFILE) {
futures[file_type] = std::async(std::launch::async, [&]() {
std::lock_guard<std::mutex> lock(osb_mtx);
return osb_data_.loadOSB(file_name);
});
} else if (file_type == DCBFILE) {
futures[file_type] = std::async(std::launch::async, [&]() {
std::lock_guard<std::mutex> lock(dcb_mtx);
return dcb_data_.loadDCB(file_name);
});
} else {
throw std::runtime_error("SPPSever Unsupported file type: " + file_type);
}
}
for (auto &future_pair : futures) {
// if (future_pair.second.has_value() && !future_pair.second.value().get()) {
// all_files_loaded = false;
// break;
// }
if ((future_pair.first == OBSFILE || future_pair.first == EPHEMERISFILE) &&
future_pair.second.has_value() && !future_pair.second.value().get()) {
all_files_loaded = false;
break;
}
}
设计一套多频多系统的程序的需求
我认为满足多频多系统的数据结构一定是:对于任意code优先级、任意频率、任意卫星系统都应该是可选择的、并且当写好其中的一个频点的程序后,整个数据结构不需要大改就满符合上述所有要求。因此我们在读完程序后一定要设计一个数据切割模块,此某块需求如下:输入为起始终止时间、卫星系统、卫星频率、以及卫星频率下的code(如果没有输出code,则按照一定的优先级排列),这样的切割数据模块与设计的配置文件对应。
bool Observation::findgpscodeByPriority(const OBSdata& sysobsdata, const string& stringfreq, const int freqnum, const vector<string> &setcode, MeasData& sysmeasdata)
{
std::map<string, std::vector<string>> gpsCodePriority_P = {
{"L1", {"C1C", "C1P", "C1Y", "C1W", "C1M", "C1S", "C1L", "C1X"}},
{"L2", {"C2P", "C2Y", "C2W", "C2C", "C2M", "C2D", "C2L", "C2S", "C2X"}},
{"L5", {"C5I", "C5Q", "C5X"}}
};
std::map<string, std::vector<string>> gpsCodePriority_L = {
{"L1", {"L1C", "L1P", "L1Y", "L1W", "L1N", "L1M", "L1S", "L1L", "L1X"}},
{"L2", {"L2P", "L2Y", "L2W", "L2C", "L2M", "L2D", "L2L", "L2S", "L2X"}},
{"L5", {"L5I", "L5Q", "L5X"}}
};
std::map<string, std::vector<string>> gpsCodePriority_S = {
{"L1", {"S1C", "S1P", "S1Y", "S1W", "S1N", "S1M", "S1S", "S1L", "S1X"}},
{"L2", {"S2P", "S2Y", "S2W", "S2C", "S2M", "S2D", "S2L", "S2S", "S2X"}},
{"L5", {"S5I", "S5Q", "S5X"}}
};
std::map<string, std::vector<string>> gpsCodePriority_D = {
{"L1", {"D1C", "D1P", "D1Y", "D1W", "D1N", "D1M", "D1S", "D1L", "D1X"}},
{"L2", {"D2P", "D2Y", "D2W", "D2C", "D2M", "D2D", "D2L", "D2S", "D2X"}},
{"L5", {"D5I", "D5Q", "D5X"}}
};
if((!setcode.empty())&&(setcode.size() == 4))
{
sysmeasdata.P[freqnum].value = sysobsdata.P.at(setcode[0]);
sysmeasdata.P[freqnum].freq = getgpsfreq(setcode[0]);
sysmeasdata.P[freqnum].code = setcode[0];
sysmeasdata.L[freqnum].value = sysobsdata.L.at(setcode[1]).value;
sysmeasdata.L[freqnum].LLi = sysobsdata.L.at(setcode[1]).LLi;
sysmeasdata.L[freqnum].snr = sysobsdata.L.at(setcode[1]).snr;
sysmeasdata.L[freqnum].freq = getgpsfreq(setcode[1]);
sysmeasdata.L[freqnum].code = setcode[1];
sysmeasdata.D[freqnum] = sysobsdata.D.at(setcode[2]);
sysmeasdata.S[freqnum] = sysobsdata.D.at(setcode[3]);
return true;
}
if ((gpsCodePriority_P.find(stringfreq) != gpsCodePriority_P.end())||(gpsCodePriority_L.find(stringfreq) != gpsCodePriority_L.end())) {
const auto& PCodePriority = gpsCodePriority_P[stringfreq];
const auto& lCodePriority = gpsCodePriority_L[stringfreq];
const auto& SCodePriority = gpsCodePriority_S[stringfreq];
const auto& DCodePriority = gpsCodePriority_D[stringfreq];
for (const string& p : PCodePriority) {
auto P_code = sysobsdata.P.find(p);
if (P_code != sysobsdata.P.end()) {
sysmeasdata.P[freqnum].value = P_code->second;
sysmeasdata.P[freqnum].freq = getgpsfreq(p);
sysmeasdata.P[freqnum].code = p;
break;
}
}
for (const string& L : lCodePriority) {
auto L_code = sysobsdata.L.find(L);
if (L_code != sysobsdata.L.end()) {
sysmeasdata.L[freqnum].value = L_code->second.value;
sysmeasdata.L[freqnum].LLi = L_code->second.LLi;
sysmeasdata.L[freqnum].snr = L_code->second.snr;
sysmeasdata.L[freqnum].freq = getgpsfreq(L);
sysmeasdata.L[freqnum].code = L;
break;
}
}
for (const string& D : DCodePriority) {
auto D_code = sysobsdata.D.find(D);
if (D_code != sysobsdata.D.end()) {
sysmeasdata.D[freqnum] = D_code->second;
break;
}
}
for (const string& S : SCodePriority) {
auto S_code = sysobsdata.D.find(S);
if (S_code != sysobsdata.S.end()) {
sysmeasdata.S[freqnum] = S_code->second;
break;
}
}
return true;
}
return false;
}
以上是我目前的重构中的一些经验,后续会持续更新,如果有空也可以更新下PPP-b2b踩过的坑。同时欢迎大佬们批评指正,如有其他想法可留言讨论。