捷联惯导系统学习4.3(静基座误差)

捷联惯导的动态误差传播本身是一个复杂的时变系统，与特定航行轨迹密切相关，很难去求解。
但在静基座条件下，误差传播退化成一个线性定常系统，通过获得其解析，解获得惯导误差传播特性。

静基座误差方程

静基座条件：vn=0;p=[Lλh]T(固定已知);(roll,yaw,pitch,不变)v^n=0;p=\left[\begin{matrix} L&\lambda&h\end{matrix}\right]^T(固定已知);(roll,yaw,pitch,不变)vn=0;p=[Lλh]T(固定已知);(roll,yaw,pitch,不变)
已知：加速度计测得比力fsfn=[00g]Tf_{sf}^n=\left[\begin{matrix} 0&0&g\end{matrix}\right]^Tfsfn=[00g]TRMh=RNh=RR_{Mh}=R_{Nh}=RRMh=RNh=R
可以简化并解耦为高度通道和水平通道，如下
{δv˙U=2wnδvE−g0sin2L(β−4β1cos2L)δL+β2δh+▽U北向速度误差分量δh˙=δvU高度误差ϕ˙E=wUϕN+wNϕU−δvN/R−ξE东向姿态误差ϕ˙N=−wUϕE+δvE/R−wUδL−ξN北向姿态误差ϕ˙U=wNϕE+δvEtanL/R+wNδL−ξU天向姿态误差δv˙E=−gϕN+2wUδvN+▽E东向速度误差δv˙N=gϕE−2wUδvE+▽N北向速度误差δL˙=δvN/R纬度误差分量δλ˙=δvEsecL/R经度误差分量\begin{cases} \delta \dot v_U=2w_n\delta v_E-g_0sin2L(\beta-4\beta_1cos2L)\delta L+\beta_2\delta h+\bigtriangledown_U&北向速度误差分量 \\ \delta \dot h=\delta v_U&高度误差\\ \dot \phi_E=w_U\phi_N+w_N\phi_U-\delta v_N/R-\xi_E&东向姿态误差 \\ \dot \phi_N=-w_U\phi_E+\delta v_E/R-w_U\delta L-\xi_N&北向姿态误差\\ \dot \phi_U=w_N\phi_E+\delta v_E tanL/R+w_N\delta L-\xi_U&天向姿态误差\\ \delta \dot v_E=-g\phi_N+2w_U\delta v_N+\bigtriangledown_E&东向速度误差\\ \delta \dot v_N=g\phi_E-2w_U\delta v_E+\bigtriangledown_N&北向速度误差\\ \delta \dot L=\delta v_N/R&纬度误差分量\\ \delta \dot \lambda=\delta v_E secL/R&经度误差分量 \end{cases}⎩⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎧δv˙U=2wnδvE−g0sin2L(β−4β1cos2L)δL+β2δh+▽Uδh˙=δvUϕ˙E=wUϕN+wNϕU−δvN/R−ξEϕ˙N=−wUϕE+δvE/R−wUδL−ξNϕ˙U=wNϕE+δvEtanL/R+wNδL−ξUδv˙E=−gϕN+2wUδvN+▽Eδv˙N=gϕE−2wUδvE+▽NδL˙=δvN/Rδλ˙=δvEsecL/R北向速度误差分量高度误差东向姿态误差北向姿态误差天向姿态误差东向速度误差北向速度误差纬度误差分量经度误差分量

高度通道

再高度通道，水平速度不大(vE≈vN≈0v_E\approx v_N\approx0vE≈vN≈0),且平稳运动(fE≈fN≈0f_E \approx f_N\approx0fE≈fN≈0)
{δv˙U=2wnδvE−g0sin2L(β−4β1cos2L)δL+β2δh+▽U北向速度误差分量δh˙=δvU高度误差\begin{cases} \delta \dot v_U=2w_n\delta v_E-g_0sin2L(\beta-4\beta_1cos2L)\delta L+\beta_2\delta h+\bigtriangledown_U&北向速度误差分量 \\ \delta \dot h=\delta v_U&高度误差\\ \end{cases}{δv˙U=2wnδvE−g0sin2L(β−4β1cos2L)δL+β2δh+▽Uδh˙=δvU北向速度误差分量高度误差
等效天向加速度计零偏为：▽U′=2wnδvE−g0sin2L(β−4β1cos2L)δL+▽U\bigtriangledown_U'=2w_n\delta v_E-g_0sin2L(\beta-4\beta_1cos2L)\delta L+\bigtriangledown_U▽U′=2wnδvE−g0sin2L(β−4β1cos2L)δL+▽U

▽U′\bigtriangledown_U'▽U′:在时域中为常值
δh(s)=1s2−β2∗▽U′s=▽U′2β2(1s+β2+1s−β2−2s)\delta h(s)=\frac{1}{s^2-\beta_2}*\frac{\bigtriangledown_U'}{s}=\frac{\bigtriangledown_U'}{2\beta_2}(\frac{1}{s+\sqrt{\beta_2}}+\frac{1}{s-\sqrt{\beta_2}}-\frac{2}{s})δh(s)=s2−β21∗s▽U′=2β2▽U′(s+β21+s−β21−s2)
进行拉普拉时变换得到频域模型：
δh(t)=▽U′β2((β2t)22!+(β2t)44!+(β2t)66!)+...(β2=3.08×10−6)s−2\delta h(t)=\frac{\bigtriangledown_U'}{\beta_2}(\frac{(\sqrt{\beta_2}t)^2}{2!}+\frac{(\sqrt{\beta_2}t)^4}{4!}+\frac{(\sqrt{\beta_2}t)^6}{6!})+...(\beta_2=3.08×10^{-6})s^{-2}δh(t)=β2▽U′(2!(β2t)2+4!(β2t)4+6!(β2t)6)+...(β2=3.08×10−6)s−2
由频域模型可知：β2t>1,即t>1β2≈1755s\sqrt{\beta_2t}>1,即t>\frac{1}{\sqrt{\beta_2}}\approx1755sβ2t>1,即t>β21≈1755s 以后2次项将会对结果产生显著明显，产生发散。
纯惯导系统高度通道不能长时间使用单独使用，往往需要高度计进行高度阻尼。
可以在短时间内（半小时内），或者固定高度的场景，此时忽略2此方以上影响使用，即高度误差为$δh(t)≈▽U′t22（物理意义上为等效天向零偏的2次积分）\delta h(t)\approx \frac{\bigtriangledown_U' t^2}{2}（物理意义上为等效天向零偏的2次积分）δh(t)≈2▽U′t2（物理意义上为等效天向零偏的2次积分）

水平通道

{ϕ˙E=wUϕN+wNϕU−δvN/R−ξE东向姿态误差ϕ˙N=−wUϕE+δvE/R−wUδL−ξN北向姿态误差ϕ˙U=wNϕE+δvEtanL/R+wNδL−ξU天向姿态误差δv˙E=−gϕN+2wUδvN+▽E东向速度误差δv˙N=gϕE−2wUδvE+▽N北向速度误差δL˙=δvN/R纬度误差分量δλ˙=δvEsecL/R经度误差分量\begin{cases} \dot \phi_E=w_U\phi_N+w_N\phi_U-\delta v_N/R-\xi_E&东向姿态误差 \\ \dot \phi_N=-w_U\phi_E+\delta v_E/R-w_U\delta L-\xi_N&北向姿态误差\\ \dot \phi_U=w_N\phi_E+\delta v_E tanL/R+w_N\delta L-\xi_U&天向姿态误差\\ \delta \dot v_E=-g\phi_N+2w_U\delta v_N+\bigtriangledown_E&东向速度误差\\ \delta \dot v_N=g\phi_E-2w_U\delta v_E+\bigtriangledown_N&北向速度误差\\ \delta \dot L=\delta v_N/R&纬度误差分量\\ \delta \dot \lambda=\delta v_E secL/R&经度误差分量 \end{cases}⎩⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎧ϕ˙E=wUϕN+wNϕU−δvN/R−ξEϕ˙N=−wUϕE+δvE/R−wUδL−ξNϕ˙U=wNϕE+δvEtanL/R+wNδL−ξUδv˙E=−gϕN+2wUδvN+▽Eδv˙N=gϕE−2wUδvE+▽NδL˙=δvN/Rδλ˙=δvEsecL/R东向姿态误差北向姿态误差天向姿态误差东向速度误差北向速度误差纬度误差分量经度误差分量
将上述公式转化为矩阵：
X=[ϕ˙Eϕ˙Nϕ˙Uδv˙Eδv˙NδL˙]TX=\left[\begin{matrix} \dot \phi_E&\dot \phi_N&\dot \phi_U&\delta \dot v_E&\delta \dot v_N&\delta \dot L \end{matrix}\right]^TX=[ϕ˙Eϕ˙Nϕ˙Uδv˙Eδv˙NδL˙]T
U=[−ξE−ξN−ξU▽E▽N0]TU=\left[\begin{matrix}-\xi_E&-\xi_N&-\xi_U&\bigtriangledown_E&\bigtriangledown_N&0 \end{matrix}\right]^TU=[−ξE−ξN−ξU▽E▽N0]T
F=[0wUwN0−1/R0−wU001/R0−wUwN00tanL/R0wN00g002wU0g00−2wU0000001/R0]F=\left[\begin{matrix}0&w_U&w_N&0&-1/R&0\\ -w_U&0&0&1/R&0&-w_U\\ w_N&0&0&tanL/R&0&w_N\\ 0&0g&0&0&2w_U&0\\ g&0&0&-2w_U&0&0\\ 0&0&0&0&1/R&0 \end{matrix}\right]F=⎣⎢⎢⎢⎢⎢⎢⎡0−wUwN0g0wU000g00wN0000001/RtanL/R0−2wU0−1/R002wU01/R0−wUwN000⎦⎥⎥⎥⎥⎥⎥⎤

得到:{X˙=FX+U公式一δλ˙=δvEsecL/R公式二\begin{cases} \dot X=FX+U&公式一\\ \delta \dot \lambda=\delta v_E secL/R&公式二\\ \end{cases}{X˙=FX+Uδλ˙=δvEsecL/R公式一公式二
假设公式一、公式二均为定常系统，进行拉式变换，可以得到：
{X(s)=(sI−F)−1[X(0)+U(s)]公式一δλ(s)=1s[(δvE(s)secL)/R+δλ(0)]公式二\begin{cases}X(s)=(sI-F)^{-1}[X(0)+U(s)]&公式一\\ \delta \lambda(s)=\frac{1}{s}[(\delta v_E(s) secL)/R+\delta\lambda(0)]&公式二\\ \end{cases}{X(s)=(sI−F)−1[X(0)+U(s)]δλ(s)=s1[(δvE(s)secL)/R+δλ(0)]公式一公式二
公式一通过其伴随矩阵，求解可以得到：
Δ(s)=∣sI−F∣=(s2+wie2)[(s2+ws2)2+4s2wf2]\Delta (s)=|sI-F|=(s^2+w_{ie}^2)[(s^2+w^2_s)^2+4s^2w_f^2]Δ(s)=∣sI−F∣=(s2+wie2)[(s2+ws2)2+4s2wf2]
ws=g/R≈84.4min(g=9.8m/s2,R=6371km)w_s=\sqrt{g/R}\approx 84.4 min(g=9.8m/s^2,R=6371km)ws=g/R≈84.4min(g=9.8m/s2,R=6371km):休拉(Schuler)角频率
wf=wiesinL:w_f=w_{ie}sinL:wf=wiesinL:傅科(Foucault)角频率
运行在3-4小时的惯导系统无需考虑傅科周期影响
其他略：详见《捷联惯导算法与组合导航原理》98页

失准角误差通道

令δvE=δvN=0,δL=0\delta v_E=\delta v_N=0,\delta L=0δvE=δvN=0,δL=0

{ϕ˙E=wUϕN+wNϕU−ξE东向姿态误差ϕ˙N=−wUϕE−ξN北向姿态误差ϕ˙U=wNϕE−ξU天向姿态误差\begin{cases} \dot \phi_E=w_U\phi_N+w_N\phi_U-\xi_E&东向姿态误差 \\ \dot \phi_N=-w_U\phi_E-\xi_N&北向姿态误差\\ \dot \phi_U=w_N\phi_E-\xi_U&天向姿态误差\\ \end{cases}⎩⎪⎨⎪⎧ϕ˙E=wUϕN+wNϕU−ξEϕ˙N=−wUϕE−ξNϕ˙U=wNϕE−ξU东向姿态误差北向姿态误差天向姿态误差
可以得到误差方程：
Δ=∣sI−[0wU−wN−wU00wN00]∣=s(s2+wie2)\Delta=|sI-\left[\begin{matrix} 0&w_U&-w_N\\ -w_U&0&0\\ w_N&0&0\\ \end{matrix}\right]|=s(s^2+w_{ie}^2)Δ=∣sI−⎣⎡0−wUwNwU00−wN00⎦⎤∣=s(s2+wie2)
失准角的震荡频率为地球自转频率，周期为24h

水平北向/水平东向通道

令:ϕE=ϕU=0,δvN=0,δL=0\phi_E=\phi_U=0,\delta v_N=0,\delta L=0ϕE=ϕU=0,δvN=0,δL=0
{ϕ˙N=δvE/R−ξNϕ˙E=−δvN/R−ξEδλ˙=δvEsecL/Rδv˙E=−gϕN+▽Eδv˙N=gϕE+▽NδL˙=δvN/R\begin{cases} \dot \phi_N=\delta v_E/R-\xi_N \\ \dot \phi_E=-\delta v_N/R-\xi_ E\\ \delta \dot \lambda=\delta v_E secL/R\\ \delta \dot v_E=-g\phi_N+\bigtriangledown_E\\ \delta \dot v_N=g\phi_E+\bigtriangledown_N\\ \delta \dot L=\delta v_N/R\\ \end{cases}⎩⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎧ϕ˙N=δvE/R−ξNϕ˙E=−δvN/R−ξEδλ˙=δvEsecL/Rδv˙E=−gϕN+▽Eδv˙N=gϕE+▽NδL˙=δvN/R
得到水平东向与与天向双通道的特征根均为：
Δ=s2+g/R=s2+ws2=0\Delta=s^2+g/R=s^2+w_s^2=0Δ=s2+g/R=s2+ws2=0
即水平通道的无阻尼震荡周期均为休拉频率ws=g/Rw_s=\sqrt{g/R}ws=g/R
约束条件物理意义上表示：在东向和天向通道上无运动。
在惯导低速、短时、小加速度条件下休拉震荡明显。

水平北向与水平东向双通道

令：ϕU=0\phi_U=0ϕU=0
{ϕ˙E=wUϕN+wNϕU−δvN/R−ξE东向姿态误差ϕ˙N=−wUϕE+δvE/R−wUδL−ξN北向姿态误差δv˙E=−gϕN+2wUδvN+▽E东向速度误差δv˙N=gϕE−2wUδvE+▽N北向速度误差δL˙=δvN/R纬度误差分量δλ˙=δvEsecL/R经度误差分量\begin{cases} \dot \phi_E=w_U\phi_N+w_N\phi_U-\delta v_N/R-\xi_E&东向姿态误差 \\ \dot \phi_N=-w_U\phi_E+\delta v_E/R-w_U\delta L-\xi_N&北向姿态误差\\ \delta \dot v_E=-g\phi_N+2w_U\delta v_N+\bigtriangledown_E&东向速度误差\\ \delta \dot v_N=g\phi_E-2w_U\delta v_E+\bigtriangledown_N&北向速度误差\\ \delta \dot L=\delta v_N/R&纬度误差分量\\ \delta \dot \lambda=\delta v_E secL/R&经度误差分量 \end{cases}⎩⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎧ϕ˙E=wUϕN+wNϕU−δvN/R−ξEϕ˙N=−wUϕE+δvE/R−wUδL−ξNδv˙E=−gϕN+2wUδvN+▽Eδv˙N=gϕE−2wUδvE+▽NδL˙=δvN/Rδλ˙=δvEsecL/R东向姿态误差北向姿态误差东向速度误差北向速度误差纬度误差分量经度误差分量
特征方程为：
Δ=s[s2+(ws+62wU)2][s2−(ws−62wU)2]\Delta=s[s^2+(w_s+\frac{\sqrt{6}}{2}w_U)^2][s^2-(w_s-\frac{\sqrt{6}}{2}w_U)^2]Δ=s[s2+(ws+26wU)2][s2−(ws−26wU)2]
无阻尼震荡周期均为休拉频率ws=g/Rw_s=\sqrt{g/R}ws=g/R和62wU\frac{\sqrt{6}}{2}w_U26wU

水平北向及方位双通道

令ϕN=0,δvE=0\phi_N=0,\delta v_E=0ϕN=0,δvE=0
可以得到：
{ϕ˙E=wNϕU−δvN/R−ξE东向姿态误差ϕ˙U=wNϕE+wNδL−ξU天向姿态误差δv˙N=gϕE+▽N北向速度误差δL˙=δvN/R纬度误差分量\begin{cases} \dot \phi_E=w_N\phi_U-\delta v_N/R-\xi_E&东向姿态误差 \\ \dot \phi_U=w_N\phi_E+w_N\delta L-\xi_U&天向姿态误差\\ \delta \dot v_N=g\phi_E+\bigtriangledown_N&北向速度误差\\ \delta \dot L=\delta v_N/R&纬度误差分量\\ \end{cases}⎩⎪⎪⎪⎨⎪⎪⎪⎧ϕ˙E=wNϕU−δvN/R−ξEϕ˙U=wNϕE+wNδL−ξUδv˙N=gϕE+▽NδL˙=δvN/R东向姿态误差天向姿态误差北向速度误差纬度误差分量
可以得到特征根为：
Δ=s2+(wN2+g/R)s2+wN2g/R=(s2+ws2)(s2+wN2)=0\Delta=s^2+(w_N^2+g/R)s^2+w_N^2g/R=(s^2+w_s^2)(s^2+w_N^2)=0Δ=s2+(wN2+g/R)s2+wN2g/R=(s2+ws2)(s2+wN2)=0
罗经效应：当存在方位误差ϕU\phi_UϕU,通过地球自转的北分量wNw_NwN耦合引起东向失准角ϕE\phi_EϕE,失准角ϕE\phi_EϕE在耦合重力引起北向速度误差，δvN\delta v_NδvN.
δvN(s)=−gwns(s2+w2)(s2+w2)ϕU0\delta v_N(s)=\frac{-gw_ns}{(s^2+w^2)(s^2+w^2)}\phi_{U0}δvN(s)=(s2+w2)(s2+w2)−gwnsϕU0
δvN(t)=−ϕU0gwNcoswNt−coswst(ws+wN)(ws−wN)\delta v_N(t)=-\phi_{U0}g_{w_N}\frac{cosw_Nt-cosw_st}{(w_s+w_N)(w_s-w_N)}δvN(t)=−ϕU0gwN(ws+wN)(ws−wN)coswNt−coswst