GPS衛(wèi)星位置解算-創(chuàng)新互聯(lián)

本文介紹了基于C語(yǔ)言的GPS衛(wèi)星位置解算原理與程序設(shè)計(jì)。針對(duì)每個(gè)原理、公式、代碼設(shè)計(jì)進(jìn)行了詳細(xì)講解，希望能夠給測(cè)繪學(xué)子們帶來(lái)幫助。
創(chuàng)新互聯(lián)公司是一家集網(wǎng)站建設(shè),弋江企業(yè)網(wǎng)站建設(shè),弋江品牌網(wǎng)站建設(shè),網(wǎng)站定制,弋江網(wǎng)站建設(shè)報(bào)價(jià),網(wǎng)絡(luò)營(yíng)銷,網(wǎng)絡(luò)優(yōu)化,弋江網(wǎng)站推廣為一體的創(chuàng)新建站企業(yè)，幫助傳統(tǒng)企業(yè)提升企業(yè)形象加強(qiáng)企業(yè)競(jìng)爭(zhēng)力。可充分滿足這一群體相比中小企業(yè)更為豐富、高端、多元的互聯(lián)網(wǎng)需求。同時(shí)我們時(shí)刻保持專業(yè)、時(shí)尚、前沿，時(shí)刻以成就客戶成長(zhǎng)自我，堅(jiān)持不斷學(xué)習(xí)、思考、沉淀、凈化自己，讓我們?yōu)楦嗟钠髽I(yè)打造出實(shí)用型網(wǎng)站。
參考書籍：
李征航黃勁松：GPS測(cè)量與數(shù)據(jù)處理（第三版）

基礎(chǔ)原理

1）衛(wèi)星位置解算流程

2）衛(wèi)星有關(guān)參數(shù)

3）衛(wèi)星位置計(jì)算過(guò)程?

程序設(shè)計(jì)

數(shù)據(jù)預(yù)處理

1）GPS時(shí)轉(zhuǎn)換(TimetoGPSsec)

2）選擇最佳歷元(select_epoch)

代碼與解析

BDS衛(wèi)星位置

基礎(chǔ)原理 1）衛(wèi)星位置解算流程

如流程圖所示，在程序設(shè)計(jì)中對(duì)于一顆衛(wèi)星的位置解算，是要經(jīng)過(guò)多次迭代，直到信號(hào)傳播時(shí)間收斂。

2）衛(wèi)星有關(guān)參數(shù)

在正式開始程序設(shè)計(jì)之前，我們還需要了解一下衛(wèi)星的相關(guān)參數(shù)及其含義。

平近點(diǎn)角M?

衛(wèi)星平均運(yùn)動(dòng)角速度為 $n=\frac{2\Pi \ }{T}\$ ，根據(jù)開普勒第三定律得：

$\frac{T^{2}}{a^{3}}=\frac{4\Pi ^{2}}{GM}$

所以衛(wèi)星平均運(yùn)動(dòng)角速度n：

$n= \sqrt{\frac{GM}{a^{3}}}$

因此，平近點(diǎn)角M定義為：

$M= n\left ( t-t_{0} \right )$

t為信號(hào)發(fā)射時(shí)刻，t0為衛(wèi)星過(guò)近地點(diǎn)時(shí)刻。

偏近點(diǎn)角E

如下圖所示：以衛(wèi)星橢圓軌道的焦點(diǎn)為原點(diǎn)，以指向近地點(diǎn)的方向?yàn)閄軸，Y軸平行于橢圓短半軸，建立軌道平面坐標(biāo)系。以衛(wèi)星橢圓中心為圓心，以衛(wèi)星軌道的長(zhǎng)半軸a為半徑作一個(gè)輔助圓。

過(guò)衛(wèi)星S作X軸的垂線，其延長(zhǎng)線與輔助圓交于S‘，圓心到S’的向徑與X軸的夾角定義為偏近點(diǎn)角E。?

平近點(diǎn)角與偏近點(diǎn)角的關(guān)系可用開普勒方程表示：

$M= E-e\sin E$

在計(jì)算偏近點(diǎn)角E時(shí)，E0可用平近點(diǎn)角M做初始值，然后進(jìn)行迭代計(jì)算，因?yàn)镋非常小，因此迭代幾次即可收斂，可將收斂標(biāo)準(zhǔn)設(shè)為1.0e-5。?

真近點(diǎn)角V

由上圖所示，真近點(diǎn)角V與偏近點(diǎn)角E存在著一定關(guān)系，真近點(diǎn)角V可又下述式子得出：

$V=\arctan \frac{\sqrt{1-e^{2}}\sin E}{\cos E-e}$

其中，e為第一偏心率。

如果已知真近點(diǎn)角，那么衛(wèi)星在軌道坐標(biāo)系中的位置向量可以表示為：

$r= \frac{a\left ( 1-e^{2} \right )}{1+e\cos V}\binom{cosV}{sinV}$

衛(wèi)星的運(yùn)動(dòng)速度? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?

由上式得得衛(wèi)星的運(yùn)動(dòng)速度為：

$r= \frac{na}{1-e\cos E}\binom{-sinE}{\sqrt{1-e^{2}\cos E}}$

由開普勒方程可得：

$\widehat{M}=\widehat{E}\left ( 1-e\sin E \right )= n$

因此，可以推導(dǎo)出衛(wèi)星的運(yùn)動(dòng)方程：

$\widehat{r}=\frac{GM}{r^{3}}r$

衛(wèi)星鐘差

衛(wèi)星鐘差方程為：

$\Delta t_{s}=a_{0}+a_{1}(t-t_{oc})+a_{2}(t-t_{oc})^{2}+\Delta t_{r}$

a0，a1，a2為N文件中的鐘頻，鐘偏和鐘漂。（廣播星歷中的衛(wèi)星鐘差誤差為10ns，等效距離為3m）

$\Delta t_{r}=-\frac{2\sqrt{A\cdot GM}}{c^{2}}e^{2}sinE$

其中 $\sqrt{A}$ 為長(zhǎng)半軸a的平方根，e為軌道偏心率，二者都由廣播星歷播發(fā)；GM為引力常數(shù)，大小為：398600500000000；E為衛(wèi)星的偏近點(diǎn)角。(因?yàn)槠c(diǎn)角E很小，該誤差在偽距中可忽略不計(jì))

3）衛(wèi)星位置計(jì)算過(guò)程?

（1）求平均角速度n：

$n_{0}=\sqrt{\frac{GM}{a^{3}}}$

$n=n_{0}+\Delta n$

$\Delta n$ 為廣播星歷(N文件)中衛(wèi)星平均運(yùn)動(dòng)速率與計(jì)算值之差（ $\Delta n$ 代碼中用deltan表示）。?

（2）求規(guī)劃時(shí)刻 $t_{k}$

$t_{k}=t-t_{oe}$

t為信號(hào)發(fā)射時(shí)刻， $t_{oe}$ 為N文件中的參考?xì)v元（ $t_{oe}$ 代碼中用TOE表示）。?

信號(hào)發(fā)射時(shí)刻為：O文件觀測(cè)時(shí)刻-信號(hào)傳播時(shí)間。

近似信號(hào)傳播時(shí)間=偽距觀測(cè)值/光速-接收機(jī)鐘差+衛(wèi)星鐘差+電離層改正+對(duì)流層改正

上式中，用參考時(shí)刻toe時(shí)刻代替了衛(wèi)星過(guò)近地點(diǎn)t0時(shí)刻。
這是因?yàn)椋簭V播星歷2h更新一次，間給參考時(shí)刻設(shè)在中央時(shí)刻時(shí)，外推間隔小于1h。而衛(wèi)星的運(yùn)行周期為12h，采取衛(wèi)星過(guò)近地點(diǎn)時(shí)刻t0，外推間隔大可能達(dá)到6h。用toe代替t0，模型簡(jiǎn)單，并且可以保證精度。
計(jì)算t和toe時(shí)，需要保證二者之差在兩小時(shí)內(nèi)，計(jì)算才能生效。因此，在計(jì)算衛(wèi)星位置前，需要根據(jù)O文件中的觀測(cè)選擇N文件中的最佳觀測(cè)歷元（代碼見后文）。

（3）求平近角

$M=M_{0}+nt_{k}$

$M_{0}$ 為N文件中的參考時(shí)間的平近點(diǎn)角。

（4）求偏近點(diǎn)角E

$E_{0}=M$

$E_{n}=M+e\sin E_{n-1}$

e為N文件中軌道偏心率；偏近點(diǎn)角需要迭代計(jì)算，將平近角M作為E的初始值帶入，因?yàn)镋很小，一般迭代幾次即可收斂，收斂條件為

$\left | E_{n}-E_{n-1} \right |< 1e-12$

1.0e-12為科學(xué)計(jì)數(shù)法，表示10的負(fù)12次方。

（5）求真近點(diǎn)角

$V=\arctan \frac{\sqrt{1-e^{2}}\sin E_{n}}{\cos E_{n}-e}$

（6）求升交角距u

$U=V_{k}+\omega$

$\omega$ 為N文件中的近地點(diǎn)角距（ $\omega$ 代碼中用omega表示）。

（7）攝動(dòng)改正

升交角距改正項(xiàng)： $\delta _{u}=C_{uC}\cos 2U+C_{uS}\sin 2U$

軌道向徑改正項(xiàng)： $\delta _{r}=C_{rC}\cos 2U+C_{rS}\sin 2U$

軌道向徑改正項(xiàng)： $\delta _{i}=C_{iC}\cos 2U+C_{iS}\sin 2U$

$C_{uC}C_{uS}C_{rC}C_{rS}C_{iC}C_{iS}$ 分別為N文件中提供的6個(gè)攝動(dòng)改正參數(shù)。（ $\delta _{u}\delta _{r}\delta _{i}$ 在代碼中分別用Epsilon_u、Epsilon_r、Epsilon_i表示）

（8）計(jì)算改正后的升交角距、軌道向徑、軌道向徑

$u_{k}=t_{k}+\delta _{u}$

$r_{k}=a\left ( 1-e\cos E_{k} \right )\delta _{r}$

$i_{k}=i_{0}+t_{k}\left ( IDOT \right )+\delta _{i}$

（9）衛(wèi)星在升交點(diǎn)軌道直角坐標(biāo)系的坐標(biāo)

$x_{k}=r_{k}\cos u_{k}$

$y_{k}=r_{k}\sin u_{k}$

（10）計(jì)算升交點(diǎn)經(jīng)度L

$L=\Omega _{0}+\left ( \Omega -\omega _{e} \right )t_{k}-\omega_{e} t_{oe}$

$\Omega _{0}$ 為N文件中GPS周開始時(shí)刻的升交點(diǎn)經(jīng)度， $\Omega$ 為N文件中升交點(diǎn)赤經(jīng)變化率， $t_{oe}$ 為N文件中參考?xì)v元時(shí)刻。

（11）計(jì)算衛(wèi)星在地固坐標(biāo)系中的空間直角坐標(biāo)系

$X=x_{k}\cos L-y_{k}\cos i_{k}sinL$

$Y=x_{k}\sin L-y_{k}\cos i_{k}cosL$

$Z_{k}=y_{k}sini_{k}$

（12）計(jì)算Z軸旋轉(zhuǎn)變換，得到的新坐標(biāo)

$x=\cos \left ( e_{E}\Delta t \right )X+\sin \left ( e_{E} \Delta t\right )Y$

$y=-\sin\left ( e_{E}\Delta t \right )X+\cos \left ( e_{E} \Delta t\right )Y$

z=Z

$e_{E}$ 為地球自轉(zhuǎn)角速度， $\Delta t$ 為信號(hào)傳播時(shí)間。

（13）重新計(jì)算衛(wèi)地距

$R=\sqrt{\left ( x-x_{r} \right )^{2}+\left ( y-y_{r} \right )^{2}+\left ( z-z_{r} \right )^{2}}$

$x_{r}y_{r}z_{r}$ 為接收機(jī)坐標(biāo)，可將O文件中接收機(jī)概略坐標(biāo)當(dāng)作接收機(jī)坐標(biāo)初始值，減少迭代次數(shù)。

（14）重新計(jì)算信號(hào)傳播時(shí)間

$t=\frac{R}{Cv}$

R為新的衛(wèi)地距，Cv表示光速。?

（15）重復(fù)步驟1~14，直到信號(hào)傳播時(shí)間收斂，收斂條件為

$\left | t_{n}-t_{n-1} \right |<0.0001$

程序設(shè)計(jì)

單顆衛(wèi)星位置解算流程：

將年月日轉(zhuǎn)換成GPS周內(nèi)秒
為衛(wèi)星選擇最佳歷元
計(jì)算信號(hào)近似傳播時(shí)間
衛(wèi)星位置解算

數(shù)據(jù)預(yù)處理 1）GPS時(shí)轉(zhuǎn)換(TimetoGPSsec)

GPS時(shí)，也稱GPS time，GPST，由GPS星載原子鐘和地面監(jiān)控站原子鐘組成的一種原子時(shí)基準(zhǔn)，其起點(diǎn)為1980年1月6日0時(shí)00分00秒。

觀測(cè)值文件中的參考時(shí)刻是以年月日進(jìn)行記錄的，不利于計(jì)算，需轉(zhuǎn)換成周內(nèi)秒的形式。?

程序設(shè)計(jì)主要思想：?

先計(jì)算到1980年1月1號(hào)0時(shí)有多少天
考慮閏年（4年一潤(rùn)，百年不潤(rùn)，四百年一潤(rùn)），閏月（潤(rùn)年潤(rùn)月為29天），以及大小月。
1980年1月6日為星期日，一周的開始，解算出來(lái)的天數(shù)需要減6（1980.1.6為GPS時(shí)起算時(shí)刻）
將得到天數(shù)除以7，得到周內(nèi)的天數(shù)，將周內(nèi)天數(shù)轉(zhuǎn)換成周內(nèi)秒數(shù)

TimetoGPSsec代碼如下：

double TimetoGPSsec(int year, int month, int day, int hour, int min, double sec)
{
	int DayofYear = 0, DayofMonth = 0, SecofWeek = 0;
	int i = 0;
	for (i = 1980; i< year; i++)
	{
        //判斷閏年
		if ((i % 4 == 0 && i % 100 != 0) || i % 400 == 0)
			DayofYear += 366;
		else
			DayofYear += 365;
	}
	for (i = 1; i< month; i++)
	{
         //判斷大小月
		if (i == 1 || i == 3 || i == 5 || i == 7 || i == 8 || i == 10 || i == 12)
			DayofMonth += 31;
		else if (i == 4 || i == 6 || i == 9 || i == 11)
			DayofMonth += 30;
		else
		{
			if ((year % 4 == 0 && year % 100 != 0) || year % 400 == 0)
				DayofMonth += 29;
			else
				DayofMonth += 28;
		}
	}
	int Day = DayofYear + DayofMonth + day - 6;//計(jì)算到1980年1月6號(hào)（星期日）有多少天
	SecofWeek = (double)(Day % 7) * 24 * 60 * 60 + (double)hour * 3600 + (double)min * 60 + sec;
	return SecofWeek;
}

其中：

傳入的參數(shù)為O、N文件中得到的年、月、日、時(shí)、分、秒。
返回值為計(jì)算得到的周內(nèi)秒。

2）選擇最佳歷元(select_epoch)

根據(jù)衛(wèi)星的PRN號(hào)以及觀測(cè)時(shí)間，選擇與N文件中時(shí)間間隔最小的數(shù)據(jù)塊進(jìn)行解算。

主要思想：?

遍歷整個(gè)N文件數(shù)據(jù)塊部分
先判斷衛(wèi)星的PRN號(hào)是否相等
設(shè)置最小時(shí)間差初始值為10000秒，每次更新最小值，并記錄最小值所對(duì)應(yīng)的歷元數(shù)
返回時(shí)間差最小的歷元數(shù)

select_epoch代碼如下：

//挑選合適的星歷
extern int select_epoch(double SecofWeek, int sPRN, pnav_body nav_b, int n_n)
{
	int n_epoch = 0;
	double min=10000;//假設(shè)最小值是1萬(wàn)秒
	double Min;
	int i;
	for (i = 0; i< n_n / 8; i++)
	{
		if (sPRN == nav_b[i].sPRN)
		{
			 Min = fabs(SecofWeek - nav_b[i].TOE);
			 if (Min<= min)
			 {
				 n_epoch = i;
				 min = Min;
			 }
		}
	}
	return n_epoch;
}

其中：

SecofWeek為觀測(cè)時(shí)刻的GPS周內(nèi)秒，sPRN為衛(wèi)星的PRN號(hào)，nav_b為導(dǎo)航文件數(shù)據(jù)塊結(jié)構(gòu)體，n_n為導(dǎo)航文件數(shù)據(jù)塊總行數(shù)（不包含頭部）
需要傳入的參數(shù)有：觀測(cè)的周內(nèi)秒時(shí)間、衛(wèi)星的PRN號(hào)、N文件中數(shù)據(jù)塊的結(jié)構(gòu)體以及N文件數(shù)據(jù)塊部分的行數(shù)。
假設(shè)最小值min是一萬(wàn)秒，當(dāng)比min小時(shí)，記錄新的min，并記下該數(shù)據(jù)塊的位置i到n_epoch中
循環(huán)所有N文件中所有的數(shù)據(jù)塊，返回最佳值

代碼與解析

在程序設(shè)計(jì)時(shí)，將單次解算衛(wèi)星位置的過(guò)程封裝成一個(gè)函數(shù)，然后放入while循環(huán)中迭代計(jì)算，當(dāng)信號(hào)傳播時(shí)間收斂后，輸出衛(wèi)星位置結(jié)果。

定義一個(gè)結(jié)構(gòu)體，用來(lái)當(dāng)作函數(shù)的傳參

//衛(wèi)星位置傳參
typedef struct Pos_t
{
	double X;
	double Y;
	double Z;
	double deltat;//改正前的接收機(jī)鐘差
	double delta_t;//改正后的接收機(jī)鐘差
}Pos_t;

#define C_V 299792458//光速（m）
#define GM 398600500000000//地心引力常數(shù)
#define math_e 2.718281828459 //e值
#define PI 3.141592653589793
#define Earth_e 7.2921151467e-5 //地球自轉(zhuǎn)角速度
#define C1 4 

//時(shí)間轉(zhuǎn)換
double GPSsec = TimetoGPSsec(obs_e[i].y, obs_e[i].m, obs_e[i].d, obs_e[i].h, obs_e[i].min, obs_e[i].sec);
//根據(jù)O文件歷元參考時(shí)間選擇合適的N文件數(shù)據(jù)
int n_epoch = select_epoch(GPSsec, obs_e[i].sPRN[j], nav_b, n_n);
//觀測(cè)時(shí)刻 - 參考時(shí)刻
double detat_toc = GPSsec - nav_b[n_epoch].TOE;
//計(jì)算近似的信號(hào)傳播時(shí)間,接收機(jī)鐘差已初始化為0(偽距/光速-接收機(jī)鐘差+衛(wèi)星鐘差)
double delta_t = (obs_b[i].obs[j][C1] / C_V) - sta[i].delta_TR + nav_b[n_epoch].sa0 + nav_b[n_epoch].sa1 * detat_toc + nav_b[n_epoch].sa2 * pow(detat_toc, 2);
//計(jì)算衛(wèi)星位置
double deltat = 0.0;//判斷收斂
while (fabs(delta_t - deltat) >0.0001)
{
	//臨時(shí)結(jié)構(gòu)體變量tmp，用來(lái)傳參
	Pos_t tmp = { 0 };
    //計(jì)算衛(wèi)星位置函數(shù)
	tmp = sat_pos(nav_b[n_epoch].deltan, nav_b[n_epoch].sqrtA, nav_b[n_epoch].TOE,
		nav_b[n_epoch].M0, nav_b[n_epoch].e, nav_b[n_epoch].omega, nav_b[n_epoch].OMEGA,
		nav_b[n_epoch].deltaomega, nav_b[n_epoch].Cuc, nav_b[n_epoch].Crc, nav_b[n_epoch].Cic,
		nav_b[n_epoch].Cus, nav_b[n_epoch].Crs, nav_b[n_epoch].Cis, nav_b[n_epoch].i0, nav_b[n_epoch].IDOT,
		delta_t, deltat, sta[i].X, sta[i].Y, sta[i].Z, GPSsec);
    //傳參，更新信號(hào)傳播時(shí)間
	deltat = tmp.deltat;
	delta_t = tmp.delta_t;
	sat[i][j].X = tmp.X;
	sat[i][j].Y = tmp.Y;
	sat[i][j].Z = tmp.Z;
}

其中：

obs_e：觀測(cè)文件歷元結(jié)構(gòu)體，記錄了每個(gè)觀測(cè)值歷元第一行時(shí)間及衛(wèi)星數(shù)等信息
obs_b：觀測(cè)文件數(shù)據(jù)塊結(jié)構(gòu)體，記錄了每個(gè)觀測(cè)歷元觀測(cè)值的數(shù)據(jù)
nav_b：導(dǎo)航文件數(shù)據(jù)塊結(jié)構(gòu)體，記錄了每個(gè)歷元的衛(wèi)星參數(shù)
觀測(cè)值文件和導(dǎo)航文件的讀取可參考：RINEX O文件讀取和RINEX N文件讀取
detat_toc：計(jì)算衛(wèi)星鐘差所用的時(shí)間
delta_t：近似信號(hào)傳播時(shí)間
obs_b[i].obs[j][C1] / C_V：表示第i個(gè)歷元第j顆衛(wèi)星的偽距值除以光速，C1和C_V為宏定義
sta[i].delta_TR：接收機(jī)鐘差，i表示第i個(gè)歷元，接收機(jī)鐘差初始化時(shí)需為0
nav_b[n_epoch].sa0 + nav_b[n_epoch].sa1 * detat_toc + nav_b[n_epoch].sa2 * pow(detat_toc, 2)：該部分表示衛(wèi)星鐘差的計(jì)算，對(duì)照上文公式，省略了末項(xiàng)
deltat = 0.0：判斷信號(hào)收斂，第一次計(jì)算時(shí)，需要初始化為0
Pos_t tmp：上文中用于傳參的結(jié)構(gòu)體
sat_pos：計(jì)算衛(wèi)星位置函數(shù)，下文會(huì)詳細(xì)展示

計(jì)算衛(wèi)星位置函數(shù)

//計(jì)算衛(wèi)星位置
 Pos_t sat_pos( double deltan,double sqrtA,double TOE,double M0,double e,
	double omega,double OMEGA,double deltaomega,double Cuc,double Crc, double Cic, 
	double Cus,double Crs,double Cis,double i0,double IDOT,double delta_t, double deltat, double X_R,
	double Y_R, double Z_R, double GPSsec)
{
	Pos_t pos_t={0};
	//計(jì)算信號(hào)發(fā)射時(shí)刻=O文件觀測(cè)時(shí)刻-信號(hào)傳播時(shí)間
	double T_S = GPSsec - delta_t;
	//計(jì)算衛(wèi)星的平均角速度n
	double n0 = sqrt(GM) / pow(sqrtA, 3);
	double n = n0 + deltan;
	//計(jì)算歸劃時(shí)間tk
	double tk = T_S - TOE;
	if (tk >32400)tk -= 604800;//規(guī)劃時(shí)間改正
	else if (tk< -32400)tk += 604800;
	//計(jì)算衛(wèi)星發(fā)射時(shí)衛(wèi)星的平近角M
	double M = M0 + n * tk;
	//迭代計(jì)算偏近點(diǎn)角
	double E = M, E0;
	do
	{
		E0 = E;
		E = M + e * sin(E);
	} while (fabs(E - E0) >1.0e-5);
	//計(jì)算真近點(diǎn)角V
	double cosv = (cos(E) - e) / (1 - e * cos(E));//cosV
	double sinv = sqrt(1 - pow(e, 2)) * sin(E) / (1 - e * cos(E));//sinV
	double Vk = atan2(sinv, cosv);//注意atan與atan2的不同
	//計(jì)算升交角距u
	double u = Vk + omega;
	//計(jì)算攝動(dòng)改正項(xiàng)
	double Epsilon_u = Cuc * cos(2 * u) + Cus * sin(2 * u);
	double Epsilon_r = Crc * cos(2 * u) + Crs * sin(2 * u);
	double Epsilon_i = Cic * cos(2 * u) + Cis * sin(2 * u);
	//計(jì)算改正后的升交角距、軌道向徑、軌道傾角
	double uk = u + Epsilon_u;
	double rk = pow(sqrtA, 2) * (1 - e * cos(E)) + Epsilon_r;
	double ik = i0 + Epsilon_i + IDOT * tk;
	//計(jì)算衛(wèi)星在升交點(diǎn)軌道直角坐標(biāo)系的坐標(biāo)
	double xk = rk * cos(uk);
	double yk = rk * sin(uk);
	//計(jì)算升交點(diǎn)經(jīng)度
	double L = OMEGA + (deltaomega - Earth_e) * tk - Earth_e * TOE; 
	//計(jì)算衛(wèi)星在地固系下的直角坐標(biāo)
	double X = xk * cos(L) - yk * cos(ik) * sin(L);
	double Y = xk * sin(L) + yk * cos(ik) * cos(L);
	double Z = yk * sin(ik);
	//通過(guò) Z 軸的旋轉(zhuǎn)變換，對(duì)該坐標(biāo)進(jìn)行地球自轉(zhuǎn)改正，得到新坐標(biāo)
	pos_t.X = cos(Earth_e * delta_t) * X + sin(Earth_e * delta_t) * Y;
	pos_t.Y = -sin(Earth_e * delta_t) * X + cos(Earth_e * delta_t) * Y;
	pos_t.Z = Z;
	//重新計(jì)算信號(hào)傳播的時(shí)間
	double R = sqrt(pow(X - X_R, 2) + pow(Y - Y_R, 2) + pow(Z - Z_R, 2));
	pos_t.deltat = delta_t;
	pos_t.delta_t = R / C_V;
	return pos_t;
}

（上述代碼為了更加清晰的體現(xiàn)哪個(gè)是新變量、哪個(gè)是舊變量，將開辟變量的工作放在函數(shù)中間；讀者自己寫代碼時(shí)，建議將開辟變量的工作統(tǒng)一放到函數(shù)的起始位置，增加運(yùn)行效率；有些語(yǔ)言不支持程序運(yùn)行后再開辟變量，例如：Fortran）

傳入?yún)?shù)：從deltann到IDOT，均為N文件中讀取的參數(shù)，delta_t和deltat為信號(hào)傳播時(shí)間，用于判斷信號(hào)傳播時(shí)間是否收斂，X_R，Y_R，Z_R為測(cè)站坐標(biāo)，首次迭代時(shí)，可初始化為0或概略坐標(biāo)（O文件頭），GPSsec為轉(zhuǎn)化成GPS時(shí)的觀測(cè)時(shí)刻

上述代碼中有關(guān)時(shí)間的參數(shù)比較多，這里再次梳理一遍，防止混淆：

觀測(cè)時(shí)刻（GPSsec）：接收機(jī)收到信號(hào)的時(shí)刻，O文件每個(gè)歷元的第一行記錄，需轉(zhuǎn)換成GPS周內(nèi)秒

衛(wèi)星鐘差時(shí)間（detat_toc）：由觀測(cè)時(shí)刻減去衛(wèi)星信號(hào)發(fā)射時(shí)刻TOE，TOE以GPS周內(nèi)秒的形式存儲(chǔ)，這個(gè)時(shí)間只用于衛(wèi)星鐘差改正

信號(hào)近似傳播時(shí)間（delta_t）：由偽距進(jìn)行計(jì)算，需要進(jìn)行衛(wèi)星鐘差、接收機(jī)鐘差、電離層和對(duì)流層改正（上述代碼沒(méi)有加入電離層和流層改正，推薦兩種簡(jiǎn)單的改正模型：電離層經(jīng)驗(yàn)?zāi)Ｐ蚄lobuchar，對(duì)流層經(jīng)驗(yàn)?zāi)Ｐ虶CAT）

信號(hào)發(fā)射時(shí)刻（T_S）：由觀測(cè)時(shí)刻減去信號(hào)近似傳播時(shí)間

規(guī)劃時(shí)間（tk）：將信號(hào)發(fā)射時(shí)刻規(guī)劃到以TOE為基準(zhǔn)的時(shí)刻

上述代碼與公式一一對(duì)應(yīng)，需要注意以下幾點(diǎn)：

1）if (tk >32400)tk -= 604800;else if (tk< -32400)tk += 604800;規(guī)劃時(shí)間改正，604800為一周的秒數(shù)，32400為一周半的秒數(shù)。該部分內(nèi)容可參考其他博主：衛(wèi)星發(fā)射時(shí)刻歸化

2）計(jì)算偏近點(diǎn)角E時(shí)，因?yàn)镋極小，迭代幾次即可收斂，因此收斂條件只要保證精度即可,在計(jì)算衛(wèi)星鐘差時(shí)，末項(xiàng)的改正數(shù)用到了偏近點(diǎn)角E，如需改正，可將E保留進(jìn)行傳參

3）C語(yǔ)言中角度計(jì)算均是以弧度為單位；atan和atan2均是求反正切函數(shù)，后者的使用范圍是是四象限角

4）上述中用到的sqrt，pow，fabs分別是：開方函數(shù)，冪函數(shù)，取絕對(duì)值函數(shù)，均需包含頭文件

BDS衛(wèi)星位置

BDS的衛(wèi)星位置解算與GPS衛(wèi)星有很大相似之處，其中MEO與IGSO衛(wèi)星位置解算與GPS一致，GEO衛(wèi)星從計(jì)算升交點(diǎn)赤經(jīng)開始有所變化，在解算位置前，需要根據(jù)衛(wèi)星的PRN號(hào)來(lái)判斷衛(wèi)星的種類，GEO計(jì)算公式如下（公式選自本科時(shí)期的論文，如有錯(cuò)誤多多包涵）：

注意：（4-27）中-5°表示角度，為GEO衛(wèi)星的傾角

你是否還在尋找穩(wěn)定的海外服務(wù)器提供商？創(chuàng)新互聯(lián)www.cdcxhl.cn海外機(jī)房具備T級(jí)流量清洗系統(tǒng)配攻擊溯源，準(zhǔn)確流量調(diào)度確保服務(wù)器高可用性，企業(yè)級(jí)服務(wù)器適合批量采購(gòu)，新人活動(dòng)首月15元起，快前往官網(wǎng)查看詳情吧

當(dāng)前標(biāo)題：GPS衛(wèi)星位置解算-創(chuàng)新互聯(lián)
本文路徑：http://www.rwnh.cn/article2/cejhoc.html

成都網(wǎng)站建設(shè)公司_創(chuàng)新互聯(lián)，為您提供關(guān)鍵詞優(yōu)化、響應(yīng)式網(wǎng)站、網(wǎng)站收錄、網(wǎng)站策劃、動(dòng)態(tài)網(wǎng)站、外貿(mào)建站

聲明：本網(wǎng)站發(fā)布的內(nèi)容（圖片、視頻和文字）以用戶投稿、用戶轉(zhuǎn)載內(nèi)容為主，如果涉及侵權(quán)請(qǐng)盡快告知，我們將會(huì)在第一時(shí)間刪除。文章觀點(diǎn)不代表本網(wǎng)站立場(chǎng)，如需處理請(qǐng)聯(lián)系客服。電話：028-86922220；郵箱：631063699@qq.com。內(nèi)容未經(jīng)允許不得轉(zhuǎn)載，或轉(zhuǎn)載時(shí)需注明來(lái)源：創(chuàng)新互聯(lián)

猜你還喜歡下面的內(nèi)容

中文字幕日韩精品一区二区免费_精品一区二区三区国产精品无卡在_国精品无码专区一区二区三区_国产αv三级中文在线

GPS衛(wèi)星位置解算-創(chuàng)新互聯(lián)