武汉大学卫星导航原理课程作业——广播星历与卫星位置计算
发布时间
阅读量：
阅读量
武汉大学卫星导航原理课程作业
采用C++语言编译，实现rinex3.04格式的广播星历文件读取和卫星的位置解算。
- 作业内容
- 根据所给定的广播星历文件，计算当天多颗GPS卫星（3颗及以上）在24h内每隔30s在瞬时地心地固系下的坐标。
- 将上述计算得到的GPS卫星位置与精密星历对比（15min间隔），评估广播星历计算出的卫星位置的误差，绘制出卫星的轨道误差时序图。
- 给定文件说明
- 广播星历文件：BRDC00IGS_R_20213480000_01D_MN.rnx
- 精密星历文件：WUM0MGXFIN_20213480000_01D_15M_ORB.SP3
- rinex 04格式说明文档。给出的广播星历文件格式为rinex 04，参数含义参照说明文档附录A5、A6。
一、广播星历文件的读取
1.首先读取所有内容，存储到lines容器里面，并记录下头文件结束的位置，便于后续获取数据
vector<string> lines;
string end_of_header = " END OF HEADER";
int end_of_line = 0;//文件头结束的位置
void readfile(string filename)
{
int i = 0;
ifstream infile(filename);
if (!infile.is_open())
{
cerr << "open error!" << endl;
}
string line="0";
while (getline(infile, line))
{
if (line == end_of_header)
{
end_of_line = i;
}
lines.push_back(line);
i = i + 1;
}
infile.close();
}//读取广播星历文件
AI写代码2.定义结构体并进行数据分配
class Data
{
private:
struct data_body
{
string sPRN;//卫星号
int TOC_Y;
int TOC_M;
int TOC_D;
int TOC_H;
int TOC_Min;
int TOC_Sec;//卫星钟的参考时刻，年月日...
double sa0;//钟差
double sa1;//钟速
double sa2;//钟漂
//广播轨道-1
double IODE;//星历发布时间
double Crs;
double deltan;
double M0;
//广播轨道-2
double Cuc;
double e;//轨道偏心角
double Cus;
double sqrtA;//长半轴平方根
//广播轨道-3
double TOE;//星历的参考时刻（GPS周内秒）
double Cic;
double OMEGA;//升交点赤经
double Cis;
//广播轨道-4
double i0;//参考时刻的轨道倾角
double Crc;
double omega;//近地点角距
double deltaomega;//
//广播轨道-5
double IDOT;
double L2code;
double GPSweek;//GPS周数
double L2flag;
//广播轨道-6
double aACC;//卫星精度
double sHEA;//卫星健康状态
double TGD;//TGD
double IODC;//钟的数据龄期
//广播轨道-7
double TTN;//电文发送时刻
double fit;//拟合区间
double none1;
double none2;//备用
};
struct position
{
int GPSweek;
double GPSsec;
double x;
double y;
double z;//地心地固坐标系下的坐标
};
public:
vector <data_body> data;
vector <position> data_P;
void getdata()
{
for (int j = end_of_line + 1; j <= size(lines) - 1; j = j + 8)
{
data_body data1;
istringstream iss(lines[j]);
iss >> data1.sPRN;
iss >> data1.TOC_Y >> data1.TOC_M >> data1.TOC_D >> data1.TOC_H >> data1.TOC_Min >> data1.TOC_Sec;
iss>> data1.sa0 >> data1.sa1 >> data1.sa2;
istringstream iss1(lines[j + 1]);
iss1 >> data1.IODE >> data1.Crs >> data1.deltan >> data1.M0;
istringstream iss2(lines[j + 2]);
iss2 >> data1.Cuc >> data1.e >> data1.Cus >> data1.sqrtA;
istringstream iss3(lines[j + 3]);
iss3 >> data1.TOE >> data1.Cic >> data1.OMEGA >> data1.Cis;
istringstream iss4(lines[j + 4]);
iss4 >> data1.i0 >> data1.Crc >> data1.omega >> data1.deltaomega;
istringstream iss5(lines[j + 5]);
iss5 >> data1.IDOT >> data1.L2code >> data1.GPSweek >> data1.L2flag;
istringstream iss6(lines[j + 6]);
iss6 >> data1.aACC >> data1.sHEA >> data1.TGD >> data1.IODC;
istringstream iss7(lines[j + 7]);//逐行解析数据
vector<double> row;//存储数据的临时容器
double value=0;
while (iss7 >> value)
{
row.push_back(value);
}
while (row.size() < 4)
{
row.push_back(0);
}//数据不足4个时，补齐为0
data1.TTN = row[0];
data1.fit = row[1];
data1.none1 = row[2];
data1.none2 = row[3];
data.push_back(data1);//将解析完的数据添加到容器里面
}
}
AI写代码二、卫星位置解算
1.计算步骤：
2.时间转换和最佳历元选取
需要将世界时转换成GPS时，同时需要查找时间差最近的历元，以此历元的数据进行卫星位置解算，代码如下：
double get_GPST(int Y, int M, int D, int H, int min, double S)
{
double y = 0;
double m = 0;
if (M <= 2)
{
y = Y - 1;
m = M + 12;
}
else
{
y = Y;
m = M;
}
double UT = 0;//世界时
double MJD = 0;//简化儒略日
UT = H + min / 60.0 + S / 3600.0;
MJD = int(365.25*y) + int(30.6001*(m + 1)) + D + UT / 24 - 679019;
int GPSweek = 0;
double secofweek = 0;
GPSweek = int((MJD - 44244) / 7);//GPS周
secofweek = (MJD - 44244 - GPSweek * 7) * 86400;//GPS周秒
return secofweek;
}//计算GPS周内秒的函数
int select_epoch(double secofweek,string PRN)
{
int n = 0;
double min_s = 10000;
double Min;
for (int i = 0; i <= size(data) - 1; i++)
{
if (data[i].sPRN == PRN)
{
Min = fabs(secofweek - data[i].TOE);
if (Min <= min_s)
{
n = i;
min_s = Min;
}
}
}
return n;
}//选择时间差最小的历元计算
void caculate(int k,double secofweek)
{
double A = pow(data[k].sqrtA,2);
double n0 = sqrt(u / pow(A, 3));
double n = n0 + data[k].deltan;//计算卫星的平均角速度n
double tk = secofweek - data[k].TOE;//计算规划时间
if (tk > 302400)
{
tk = tk - 604800;
}
if (tk < -302400)
{
tk = tk + 604800;
}//考虑跨周的影响
double Mk = data[k].M0 + n * tk;//计算平近点角
double E0 = 0;
double E1 = Mk;
while (abs(E1-E0)>1e-6)
{
E0 = E1;
E1 = Mk + data[k].e*sin(E0);
}//计算偏近点角
double cosv = (cos(E1) - data[k].e) / (1 - data[k].e * cos(E1));
double sinv = sqrt(1 - pow(data[k].e, 2))*sin(E1) / (1 - data[k].e * cos(E1));
double Vk = atan2(sinv, cosv);
double w = Vk + data[k].omega;
double Ou = data[k].Cuc*cos(2 * w) + data[k].Cus*sin(2 * w);
double Or = data[k].Crc*cos(2 * w) + data[k].Crs*sin(2 * w);
double Oi = data[k].Cic*cos(2 * w) + data[k].Cis*sin(2 * w);
double uk = w + Ou;
double rk = A * (1 - data[k].e*cos(E1)) + Or;
double ik = data[k].i0 + Oi + data[k].IDOT*tk;
double xk1 = rk * cos(uk);
double yk1 = rk * sin(uk);
double Ok = data[k].OMEGA + (data[k].deltaomega - Oe)*tk - Oe * data[k].TOE;
double xk = xk1 * cos(Ok) - yk1 * cos(ik)*sin(Ok);
double yk = xk1 * sin(Ok) + yk1 * cos(ik)*cos(Ok);
double zk = yk1 * sin(ik);
position P;
P.GPSweek = data[k].GPSweek;
P.GPSsec = secofweek;
P.x = xk;
P.y = yk;
P.z = zk;
data_P.push_back(P);//把坐标数据添加到容器里面
}
AI写代码三、结果文件保存
将计算所得结果按照指定要求保存在文本文件中：
void savefile()
{
ofstream outfile("广播星历.txt");
outfile <<fixed<< setprecision(3);
for (int i = 0; i <= size(data_P) - 1; i++)
{
outfile << data_P[i].GPSweek << " " << data_P[i].GPSsec << " " << data_P[i].x << " " << data_P[i].y << " " << data_P[i].z << endl;
}
}//将广播星历的数据写入文件
AI写代码武汉大学测绘学院在读本科生一枚，本次的内容就分享到这里了，希望大家多多关注，后面会持续给大家分享更多的优质内容。本次作业的所有测试用例和源代码都可以免费下载，希望大家点点关注支持一下。
全部评论 (0)
还没有任何评论哟~