GPS从入门到放弃（八） — GPS卫星速度解算
在阅读这一篇强烈建议先阅读
GPS卫星位置解算
。
为了计算卫星速度，需要对卫星的位置求导。
计算各参变量对时间的导数如下：
M
˙
k
=
n
\dot{M}_k = n
M
˙
k
=
n
E
˙
k
=
M
˙
k
1
−
e
cos
E
k
\dot{E}_k = \frac{\dot{M}_k}{1-e\cos E_k}
E
˙
k
=
1
−
e
cos
E
k
M
˙
k
ν
˙
k
=
1
−
e
2
E
˙
k
1
−
e
cos
E
k
\dot{\nu}_k = \frac{\sqrt{1-e^2}\dot{E}_k}{1-e\cos E_k}
ν
˙
k
=
1
−
e
cos
E
k
1
−
e
2
E
˙
k
Φ
˙
k
=
ν
˙
k
\dot{\Phi}_k = \dot{\nu}_k
Φ
˙
k
=
ν
˙
k
δ
u
˙
k
=
2
Φ
˙
k
(
C
u
s
cos
2
Φ
k
−
C
u
c
sin
2
Φ
k
)
\delta\dot{u}_k = 2\dot{\Phi}_k(C_{us}\cos2\Phi_k – C_{uc}\sin2\Phi_k)
δ
u
˙
k
=
2
Φ
˙
k
(
C
u
s
cos
2
Φ
k
−
C
u
c
sin
2
Φ
k
)
δ
r
˙
k
=
2
Φ
˙
k
(
C
r
s
cos
2
Φ
k
−
C
r
c
sin
2
Φ
k
)
\delta\dot{r}_k = 2\dot{\Phi}_k(C_{rs}\cos2\Phi_k – C_{rc}\sin2\Phi_k)
δ
r
˙
k
=
2
Φ
˙
k
(
C
r
s
cos
2
Φ
k
−
C
r
c
sin
2
Φ
k
)
δ
i
˙
k
=
2
Φ
˙
k
(
C
i
s
cos
2
Φ
k
−
C
i
c
sin
2
Φ
k
)
\delta\dot{i}_k = 2\dot{\Phi}_k(C_{is}\cos2\Phi_k – C_{ic}\sin2\Phi_k)
δ
i
˙
k
=
2
Φ
˙
k
(
C
i
s
cos
2
Φ
k
−
C
i
c
sin
2
Φ
k
)
Ω
˙
k
=
Ω
˙
−
Ω
˙
e
\dot{\Omega}_k = \dot{\Omega} – \dot{\Omega}_e
Ω
˙
k
=
Ω
˙
−
Ω
˙
e
i
˙
k
=
i
˙
+
δ
i
˙
k
\dot{i}_k = \dot{i} + \delta\dot{i}_k
i
˙
k
=
i
˙
+
δ
i
˙
k
r
˙
k
=
A
e
E
˙
k
sin
E
k
+
δ
r
˙
k
\dot{r}_k = Ae\dot{E}_k\sin E_k + \delta\dot{r}_k
r
˙
k
=
A
e
E
˙
k
sin
E
k
+
δ
r
˙
k
u
˙
k
=
Φ
˙
k
+
δ
u
˙
k
\dot{u}_k = \dot{\Phi}_k + \delta\dot{u}_k
u
˙
k
=
Φ
˙
k
+
δ
u
˙
k
计算卫星在轨道平面内的速度和WGS-84坐标系中的速度如下：
x
˙
k
′
=
r
˙
k
cos
u
k
−
r
k
u
˙
k
sin
u
k
\dot{x}’_k = \dot{r}_k\cos u_k – r_k \dot{u}_k \sin u_k
x
˙
k
′
=
r
˙
k
cos
u
k
−
r
k
u
˙
k
sin
u
k
y
˙
k
′
=
r
˙
k
sin
u
k
+
r
k
u
˙
k
cos
u
k
\dot{y}’_k = \dot{r}_k\sin u_k + r_k \dot{u}_k \cos u_k
y
˙
k
′
=
r
˙
k
sin
u
k
+
r
k
u
˙
k
cos
u
k
x
˙
k
=
−
y
k
Ω
˙
k
−
(
y
˙
k
′
cos
i
k
−
z
k
i
˙
k
)
sin
Ω
k
+
x
˙
k
′
cos
Ω
k
\dot{x}_k = -y_k\dot{\Omega}_k – (\dot{y}’_k\cos i_k – z_k\dot{i}_k)\sin\Omega_k + \dot{x}’_k\cos\Omega_k
x
˙
k
=
−
y
k
Ω
˙
k
−
(
y
˙
k
′
cos
i
k
−
z
k
i
˙
k
)
sin
Ω
k
+
x
˙
k
′
cos
Ω
k
y
˙
k
=
x
k
Ω
˙
k
+
(
y
˙
k
′
cos
i
k
−
z
k
i
˙
k
)
cos
Ω
k
+
x
˙
k
′
sin
Ω
k
\dot{y}_k = x_k\dot{\Omega}_k + (\dot{y}’_k\cos i_k – z_k\dot{i}_k)\cos\Omega_k + \dot{x}’_k\sin\Omega_k
y
˙
k
=
x
k
Ω
˙
k
+
(
y
˙
k
′
cos
i
k
−
z
k
i
˙
k
)
cos
Ω
k
+
x
˙
k
′
sin
Ω
k
z
˙
k
=
y
˙
k
′
sin
i
k
+
y
k
′
i
˙
k
cos
i
k
\dot{z}_k = \dot{y}’_k\sin i_k + y’_k\dot{i}_k\cos i_k
z
˙
k
=
y
˙
k
′
sin
i
k
+
y
k
′
i
˙
k
cos
i
k
在之前的位置解算的过程中，我们已经求得了很多中间变量，利用以上公式，即可求得卫星速度。
代码如下：
M_k_Dot = n
E_k_Dot = M_k_Dot/(1-e*math.cos(E_k))
nu_k_Dot = math.sqrt(1-e**2)*E_k_Dot/(1-e*math.cos(E_k))
Phi_k_Dot = nu_k_Dot
delta_u_k_Dot = 2*Phi_k_Dot*(Cus*math.cos(2*Phi_k) - Cuc*math.sin(2*Phi_k))
delta_r_k_Dot = 2*Phi_k_Dot*(Crs*math.cos(2*Phi_k) - Crc*math.sin(2*Phi_k))
delta_i_k_Dot = 2*Phi_k_Dot*(Cis*math.cos(2*Phi_k) - Cic*math.sin(2*Phi_k))
Omega_k_Dot = Omega_Dot - Omega_e_Dot
i_k_Dot = i_Dot + delta_i_k_Dot
r_k_Dot = A*e*E_k_Dot*math.sin(E_k) + delta_r_k_Dot
u_k_Dot = Phi_k_Dot + delta_u_k_Dot
x_p_k_Dot = r_k_Dot*math.cos(u_k) - r_k*u_k_Dot*math.sin(u_k)
y_p_k_Dot = r_k_Dot*math.sin(u_k) + r_k*u_k_Dot*math.cos(u_k)
x_k_Dot = -y_k*Omega_k_Dot - (y_p_k_Dot*math.cos(i_k) - z_k*i_k_Dot)*math.sin(Omega_k) + x_p_k_Dot*math.cos(Omega_k)
y_k_Dot = x_k*Omega_k_Dot + (y_p_k_Dot*math.cos(i_k) - z_k*i_k_Dot)*math.cos(Omega_k) + x_p_k_Dot*math.sin(Omega_k)
z_k_Dot = y_p_k_Dot*math.sin(i_k) + y_p_k*i_k_Dot*math.cos(i_k)
print("x_k_Dot={}".format(x_k_Dot))
print("y_k_Dot={}".format(y_k_Dot))
print("z_k_Dot={}".format(z_k_Dot))
用
GPS卫星位置解算
中的数据进行解算，可以得到卫星在WGS-84坐标系中的速度为：
[
x
˙
k
y
˙
k
z
˙
k
]
=
[
1088.4093298528
2189.0005213284817
−
1527.1796186285617
]
(
米
/
秒
)
\left[ \begin{array}{r} \dot{x}_k \\ \dot{y}_k \\ \dot{z}_k \end{array} \right]= \left[ \begin{array}{c} 1088.4093298528 \\ 2189.0005213284817 \\ -1527.1796186285617 \end{array} \right](米/秒)
⎣
⎡
x
˙
k
y
˙
k
z
˙
k
⎦
⎤
=
⎣
⎡
1
0
8
8
.
4
0
9
3
2
9
8
5
2
8
2
1
8
9
.
0
0
0
5
2
1
3
2
8
4
8
1
7
−
1
5
2
7
.
1
7
9
6
1
8
6
2
8
5
6
1
7
⎦
⎤
(
米
/
秒
)