摘要
本发明提供一种视觉导航完好性监测的计算方法。该方法利用合适的视觉定位模型、数学算法和丰富的导航测量数据,提高了定位精度和定位结果的可用性。解决了GNSS在复杂环境下无法保证可用性而导致卫星完好性算法性能不足的问题,有助于实现飞机紧密进近和自主着陆,对确保航空飞行安全具有重要意义。
图1 视觉完好性监测流程图
?
图2 人工目标的定位模型
?
图3 机场跑道的定位模型
?
图4 保护水平计算流程图
交叉参考相关专利申请
略
技术领域
本发明涉及卫星导航技术领域,具体地说,涉及一种用于视觉导航完好性监测的计算方法。
1 背景
随着卫星导航技术的发展和广泛应用,人们越来越认识到全球卫星导航系统完好性的重要性。一份调查报告显示,在航空应用中,致命事故大多发生在最后的进近和降落阶段。尽管飞机的进近和降落阶段虽然最短,但也是保障生命安全最关键的环节。
全球定位系统(GPS)接收机的自主完好性监测通过检查卫星定位观测的一致性来提供导航系统的完好性。然而,在飞行路径的进近和着陆阶段,由于卫星定位观测的不足,现有的接收机自主完好性监测方法的性能一般不能满足严格航空应用对完好性和可用性的要求。因此,开发辅助GNSS的视觉导航设备已成为航空卫星导航领域的一个热点。
目前,北京航空航天大学和西安导航技术研究所等高校和研究机构对卫星导航的完好性进行了理论和技术研究。从GPS到差分GPS,再到其它导航系统,如惯性导航系统和组合导航系统,对卫星导航的完好性进行了越来越多的研究。与此同时,卫星定位完好性监测方法已经相对成熟,但卫星定位信号在所有环境或操作条件下的永久可用性无法得到保证,导致卫星完好性算法的性能不足。因此,本发明提出了一种用于监测视觉辅助导航系统完好性的研究方法。
2 总结
为解决现有技术中由于无法保证复杂环境或操作条件下卫星定位的性能可用性而导致卫星完好性算法性能不足的问题,本发明提供了一种视觉导航完整性监测计算方法。其中包括:
- 根据具体的视觉定位场景建立合适的视觉定位模型,计算出视觉定位解;
- 修正视觉系统参数的误差,并计算误差状态向量;
- 建立视觉观测模型,计算视觉系统观测矩阵
H
H
H;
- 引入视觉故障偏差,分析故障偏差向量
b
b
b随视觉观测值的变化;
- 结合视觉系统观测方程和卫星定位方程计算组合导航定位方程;
- 计算不同故障模式下的保护水平。
步骤(1)中的视觉定位模型在具体的视觉定位场景中可分为两类,第一类为手动设置地标作为着陆参考,第二类为机场跑道作为着陆参考。手动设置的地面标志为“H”或“I”标志;而使用机场跑道作为着陆参考就是借用机场环境附近的机场跑道等高线,这样相机就可以可靠地识别机场跑道目标。进一步,在导航坐标系下构建视觉定位方程如下:
t
m
n
=
s
m
n
+
p
n
+
q
n
t_m^n=s_m^n+p^n+q^n
tmn?=smn?+pn+qn 当路标在相机系下的视线向量
s
m
c
s_m^c
smc?、相机在机体系的位置
q
b
=
[
x
1
?
y
2
?
z
3
]
T
q^b=[x_1\ y_2\ z_3]^T
qb=[x1??y2??z3?]T和着陆路标的位置向量
t
m
n
=
[
x
2
?
y
2
?
z
2
]
T
t_m^n=[x_2\ y_2\ z_2]^T
tmn?=[x2??y2??z2?]T已知时,飞机在三维世界的位置
p
n
=
[
x
3
?
y
3
?
z
3
]
T
p^n=[x_3\ y_3\ z_3]^T
pn=[x3??y3??z3?]T可以被计算。其中
t
m
n
t_m^n
tmn?表示第
m
m
m个路标点在导航系中的位置向量;
s
m
n
s_m^n
smn?表示第
m
m
m个目标在导航系下的视线向量;
p
n
p^n
pn表示在导航系下的飞机位置向量;
q
n
q^n
qn表示相机中心到飞机质心的距离,它也被称为相机和飞机之间的最小角度偏移向量;
q
b
q^b
qb表示机体系中相机的位置,它可以通过相机进行校准。
进一步,对视觉系统参数的误差进行修正,并计算误差状态向量的步骤具体如下:
- 将基于上述数学参数的误差状态向量表示为:
δ
x
=
[
β
p
n
?
δ
t
m
n
]
\delta x=[\beta p^n \ \delta t_m^n]
δx=[βpn?δtmn?] 其中
δ
p
n
\delta p^n
δpn指飞机三维位置的误差向量:
δ
p
n
=
p
~
n
?
p
n
\delta p^n=\tilde{p}^n-p^n
δpn=p~?n?pn;
δ
t
m
n
\delta t_m^n
δtmn?是路标误差向量:
δ
t
m
n
=
t
~
m
n
?
t
m
n
\delta t_m^n=\tilde{t}_m^n-t_m^n
δtmn?=t~mn??tmn?,其中标称值用波浪线表示。 - 当一个着陆路标进入和离开相机视野时,路标的数量实际上会随着时间而变化,但必须有一个合理的最大路标跟踪数量,以保持计算效率。
进一步,计算视觉系统观测矩阵
H
i
s
H_{is}
His?的步骤具体如下:
-
m
g
mg
mg路标视觉观测方程表示为:
z
m
(
t
i
)
=
T
c
p
i
x
s
 ̄
m
c
(
t
i
)
+
v
(
t
i
)
z_m(t_i)=T_c^{pix}\underline{s}_m^c(t_i)+v(t_i)
zm?(ti?)=Tcpix?s?mc?(ti?)+v(ti?) 其中
z
m
(
t
i
)
z_m(t_i)
zm?(ti?)是着陆路标点
m
m
m的像素位置;
T
c
p
i
x
T_c^{pix}
Tcpix?是像素坐标系和相机坐标系之间的变换矩阵;
s
m
c
s_m^c
smc?指像素位置对应于从相机焦点向外延伸的视线向量;
s
m
z
c
s_{m_z}^c
smz?c?是相机在
s
z
s_z
sz?方向上到光学中心的距离;
s
 ̄
m
c
(
t
i
)
\underline{s}_m^c(t_i)
s?mc?(ti?)是相机到目标
m
m
m的视线向量的标准化形式;
v
(
t
i
)
v(t_i)
v(ti?)是一个均值为0协方差为
R
R
R的独立的加性高斯白噪声,它满足以下关系:
E
[
v
(
t
i
)
v
t
(
t
j
)
]
=
{
R
,
??
t
i
=
t
j
0
,
???
t
i
≠
t
j
E[v(t_i)v^t(t_j)]=\begin{cases} R, \ \ t_i=t_j \\ 0, \ \ \ t_i \neq t_j \end{cases}
E[v(ti?)vt(tj?)]={R,??ti?=tj?0,???ti??=tj?? - 设在一个高
H
H
H宽
W
W
W的矩形
(
M
×
N
)
(M\times N)
(M×N)像素网格中,像素坐标系与相机坐标系的变换关系为
T
c
p
i
x
T_c^{pix}
Tcpix?,它是已知的,
T
c
p
i
x
T_c^{pix}
Tcpix?和
T
p
i
x
c
T_{pix}^c
Tpixc?之间的关系如下:
T
p
i
x
c
=
[
?
H
f
M
0
H
(
M
+
1
)
2
f
M
0
W
f
N
?
W
(
N
+
1
)
2
f
N
0
0
1
]
=
(
T
c
p
i
x
)
?
1
T_{pix}^c = \begin{bmatrix} -\frac{H}{fM} & 0 & \frac{H(M+1)}{2fM} \\ 0 & \frac{W}{fN} & -\frac{W(N+1)}{2fN} \\ 0 & 0 & 1 \end{bmatrix} = (T_c^{pix})^{-1}
Tpixc?=?????fMH?00?0fNW?0?2fMH(M+1)??2fNW(N+1)?1?????=(Tcpix?)?1 一个非线性方程
z
m
(
t
i
)
=
T
c
p
i
x
s
 ̄
m
c
(
t
i
)
+
v
(
t
i
)
z_m(t_i)=T_c^{pix}\underline{s}_m^c(t_i)+v(t_i)
zm?(ti?)=Tcpix?s?mc?(ti?)+v(ti?)被线性化,一个非线性函数被设为
h
[
x
(
t
i
)
]
=
T
c
p
i
x
s
 ̄
m
c
(
t
i
)
h[x(t_i)]=T_c^{pix}\underline{s}_m^c(t_i)
h[x(ti?)]=Tcpix?s?mc?(ti?),计算了一阶泰勒展开式,得到视觉系统观测矩阵
H
H
H,由下式得到的视觉系统观测矩阵
H
H
H中元素的偏导数为:
H
=
[
?
h
?
p
n
?
0
]
H=[\frac{\partial h}{\partial p^n}\ 0]
H=[?pn?h??0] 因为,
q
n
=
C
b
n
C
c
b
s
m
c
+
C
b
n
q
d
=
C
b
n
(
C
c
b
s
m
c
+
q
d
)
q^n=C_b^nC_c^bs_m^c+C_b^nq^d=C_b^n(C_c^bs_m^c+q^d)
qn=Cbn?Ccb?smc?+Cbn?qd=Cbn?(Ccb?smc?+qd)
h
h
h相对于
p
n
p^n
pn的偏导数为:
?
h
?
p
n
=
1
s
m
z
c
T
c
p
i
x
(
s
 ̄
m
c
C
n
b
C
b
c
?
C
n
b
C
b
c
)
\frac{\partial h}{\partial p^n}=\frac{1}{s_{m_z}^c}T_c^{pix}(\underline{s}_m^cC_n^bC_b^c-C_n^bC_b^c)
?pn?h?=smz?c?1?Tcpix?(s?mc?Cnb?Cbc??Cnb?Cbc?) 因为视觉系统观测矩阵
H
H
H与
t
m
n
t_m^n
tmn?无关,
H
H
H对
δ
t
m
n
\delta t_m^n
δtmn?的偏导数为零。
进一步,分析故障偏差向量
b
b
b对视线观测向量
z
z
z的变化过程如下:
- 状态向量
x
x
x的线性化视觉观测模型表示为:
z
=
H
x
+
v
+
b
z=Hx+v+b
z=Hx+v+b 其中
z
z
z是视觉观测向量;状态向量
x
x
x是
n
n
n维的列向量;
v
v
v是一个
m
m
m维的观测噪声向量,服从均值为零的高斯分布,它的协方差矩阵是对角阵
R
R
R,表示为
R
=
σ
2
I
m
?
n
R=\sigma^2I_{m*n}
R=σ2Im?n?;
b
b
b是一个
m
m
m维的故障偏差向量,在无故障条件下,
b
b
b是一个零向量,当第
i
i
i个测量中出现单故障时,
b
b
b中的第
i
i
i个元素仍表示为
b
b
b。 - 当视觉观测向量变化时,视觉观测误差
δ
z
\delta z
δz表示为:
δ
z
=
H
δ
x
\delta z=H\delta x
δz=Hδx;因此,得到
z
=
H
x
+
v
+
b
z=Hx+v+b
z=Hx+v+b的最小二乘解为:
δ
x
=
(
H
T
H
)
?
1
H
T
δ
z
\delta x=(H^TH)^{-1}H^T\delta z
δx=(HTH)?1HTδz 如果不考虑噪声向量,则有状态向量
x
x
x的线性观测模型可得:
δ
x
=
x
^
?
x
=
H
ˉ
z
=
H
ˉ
b
\delta x = \hat{x}-x=\bar{H}z=\bar{H}b
δx=x^?x=Hˉz=Hˉb 其中
H
ˉ
=
(
H
T
H
)
?
1
H
T
\bar{H}=(H^TH)^{-1}H^T
Hˉ=(HTH)?1HT - 由于故障偏差向量
b
b
b对像素点
(
x
,
y
)
(x,y)
(x,y)的每个分量
x
x
x和
y
y
y都有影响,因此故障偏差向量在极坐标下表示为:
b
i
=
∣
∣
b
∣
∣
c
o
s
θ
b_i=||b||cos\theta
bi?=∣∣b∣∣cosθ和
b
j
=
∣
∣
b
∣
∣
s
i
n
θ
b_j=||b||sin\theta
bj?=∣∣b∣∣sinθ;则
b
i
=
∣
∣
b
∣
∣
c
o
s
θ
b_i=||b||cos\theta
bi?=∣∣b∣∣cosθ为故障偏差向量
b
b
b对分量
x
x
x的影响,
b
j
=
∣
∣
b
∣
∣
s
i
n
θ
b_j=||b||sin\theta
bj?=∣∣b∣∣sinθ为故障偏差向量
b
b
b对分量
y
y
y的影响;因此,位置误差的范数为:
∣
∣
δ
x
∣
∣
2
=
(
H
ˉ
b
)
T
H
ˉ
b
||\delta x||^2=(\bar{H}b)^T\bar{H}b
∣∣δx∣∣2=(Hˉb)THˉb;因为,
(
H
ˉ
b
)
T
H
ˉ
b
=
b
T
H
ˉ
T
H
ˉ
b
=
b
i
b
i
(
H
ˉ
T
H
ˉ
)
i
i
+
b
i
b
j
(
H
ˉ
T
H
ˉ
)
i
j
+
b
i
b
j
(
H
T
H
)
j
i
+
b
j
b
j
(
H
ˉ
T
H
ˉ
)
j
j
(\bar{H}b)^T\bar{H}b=b^T\bar{H}^T\bar{H}b=b_ib_i(\bar{H}^T\bar{H})_{ii}+b_ib_j(\bar{H}^T\bar{H})_{ij}+b_ib_j(H^TH)_{ji}+b_jb_j(\bar{H}^T\bar{H})_{jj}
(Hˉb)THˉb=bTHˉTHˉb=bi?bi?(HˉTHˉ)ii?+bi?bj?(HˉTHˉ)ij?+bi?bj?(HTH)ji?+bj?bj?(HˉTHˉ)jj? 将其代入极坐标公式,可得故障偏差向量
b
b
b对水平位置误差的影响公式为:
δ
x
=
(
∣
∣
b
∣
∣
2
(
s
i
n
2
θ
(
H
ˉ
T
H
ˉ
)
i
i
+
s
i
n
2
θ
(
H
ˉ
T
H
ˉ
)
i
j
+
c
o
s
2
θ
(
H
ˉ
T
H
ˉ
)
j
j
)
)
1
2
\delta x=\big( ||b||^2(sin2\theta(\bar{H}^T\bar{H})_{ii}+sin2\theta(\bar{H}^T\bar{H})_{ij}+cos^2\theta(\bar{H}^T\bar{H})_{jj}) \big)^{\frac{1}{2}}
δx=(∣∣b∣∣2(sin2θ(HˉTHˉ)ii?+sin2θ(HˉTHˉ)ij?+cos2θ(HˉTHˉ)jj?))21? 因此,因此,当角度值
θ
\theta
θ发生变化时,水平位置误差也会发生相应的变化,特别是当对
θ
\theta
θ取一些特殊的角度值时,水平位置误差
δ
x
\delta x
δx的公式可以简化。
进一步介绍了组合导航定位方程的求解过程:
- 假设卫星定位方程是
z
=
G
x
z=Gx
z=Gx,结合视觉系统的观测矩阵
H
H
H和卫星的观测矩阵
G
G
G来获得扩展的观测矩阵
H
′
H'
H′,它可以表示为,
H
′
=
[
H
G
]
H'=\begin{bmatrix} H \\ G \end{bmatrix}
H′=[HG?]
H
′
H'
H′是一个
m
×
n
m\times n
m×n的线性观测矩阵。 - 当一个视觉系统是由两个相机或三个相机组成时,单个相机的变换矩阵
M
1
M_1
M1?、
M
2
M_2
M2?和
M
i
M_i
Mi?可以合并成
M
M
M,
M
M
M可以表示为:
M
=
[
M
1
0
0
0
M
2
0
0
0
M
i
]
M=\begin{bmatrix} M_1 & 0 & 0\\ 0 & M_2 & 0 \\ 0 & 0 & M_i \end{bmatrix}
M=???M1?00?0M2?0?00Mi?????
进而,将保护水平的计算方法分为三个假设:假设1、假设2和假设3。假设1是在观测中没有故障。假设2是观测中有单个故障。假设3是在观测中存在两个故障。
进而,在假设1的情况下,可以根据无故障放大因子
K
f
f
K_{ff}
Kff?和水平位置标准差计算保护水平,其中无故障放大因子
K
f
f
K_{ff}
Kff?由无故障概率
P
f
f
P_{ff}
Pff?确定。因此,无故障保护水平的计算公式为:
H
P
L
1
=
K
f
f
(
k
)
?
σ
x
HPL_1=K_{ff}(k)*\sigma_x
HPL1?=Kff?(k)?σx? 其中
K
f
f
K_{ff}
Kff?是无故障放大因子,
σ
x
\sigma_x
σx?是飞机水平方向的标准差。
在假设2的情况下,单个视觉观测故障或单个卫星观测故障是可能的。当只有一个观测故障时,故障偏差向量
b
b
b中只包含一个偏差分量,现在水平位置误差
δ
x
\delta x
δx可以表示为:
δ
x
=
x
^
?
x
=
H
ˉ
z
=
H
ˉ
b
\delta x = \hat{x} - x=\bar{H}z=\bar{H}b
δx=x^?x=Hˉz=Hˉb 其中
x
^
\hat{x}
x^表示最小二乘估计值,
H
=
(
H
T
H
)
?
1
H
T
H=(H^TH)^{-1}H^T
H=(HTH)?1HT。每个观测的斜率被设置为
H
s
l
o
p
e
Hslope
Hslope,表示如下:
H
s
o
l
p
e
i
=
δ
x
2
D
=
∣
∣
δ
x
∣
∣
2
∣
∣
p
∣
∣
2
=
δ
x
T
δ
x
p
T
p
=
b
T
H
ˉ
T
H
ˉ
b
b
T
P
T
P
b
Hsolpe_i=\frac{\delta x^2}{D}=\frac{||\delta x||^2}{||p||^2}=\frac{\delta x^T\delta x}{p^Tp}=\frac{b^T\bar{H}^T\bar{H}b}{b^TP^TPb}
Hsolpei?=Dδx2?=∣∣p∣∣2∣∣δx∣∣2?=pTpδxTδx?=bTPTPbbTHˉTHˉb? 其中
b
T
H
ˉ
T
H
ˉ
b
=
b
i
b
i
(
H
ˉ
T
H
ˉ
)
i
i
+
b
i
b
j
(
H
ˉ
T
H
ˉ
)
i
j
+
b
i
b
j
(
H
ˉ
T
H
ˉ
)
j
i
+
b
j
b
j
(
H
ˉ
T
H
ˉ
)
j
j
b^T\bar{H}^T\bar{H}b=b_ib_i(\bar{H}^T\bar{H})_{ii}+b_ib_j(\bar{H}^T\bar{H})_{ij}+b_ib_j(\bar{H}^T\bar{H})_{ji}+b_jb_j(\bar{H}^T\bar{H})_{jj}
bTHˉTHˉb=bi?bi?(HˉTHˉ)ii?+bi?bj?(HˉTHˉ)ij?+bi?bj?(HˉTHˉ)ji?+bj?bj?(HˉTHˉ)jj?;
D
D
D是统计检验量;
p
p
p表示奇偶向量;
P
P
P表示奇偶变换矩阵,
p
=
P
b
p=Pb
p=Pb。因此,水平定位误差的保护水平
H
P
L
2
HPL_2
HPL2?为:
H
P
L
2
=
H
s
o
l
p
e
m
a
x
×
s
b
i
a
s
HPL_2=Hsolpe_{max}\times s_{bias}
HPL2?=Hsolpemax?×sbias? 其中
s
b
i
a
s
=
σ
λ
n
?
4
s_{bias}=\sigma \frac{\sqrt{\lambda}}{\sqrt{n-4}}
sbias?=σn?4
?λ
?? 其中
σ
\sigma
σ是噪声协方差;
λ
\lambda
λ是分散参数。当无故障时,
H
P
L
2
HPL_2
HPL2?的值最小;当观测中存在单故障时,
H
P
L
2
HPL_2
HPL2?的值会增加。
在假设3的情况下,有可能两个视觉观测发生故障,或者两个卫星定位观测发生故障,又或者一个视觉观测发生故障和一个卫星定位观测发生故障。当两个观测发生故障时,第
i
i
i次观测和第
j
j
j次观测同时故障的故障偏差向量为
b
(
i
)
b^{(i)}
b(i)和
b
(
j
)
b^{(j)}
b(j),可以表示为:
b
=
[
0
?
?
?
0
?
b
(
i
)
?
b
(
j
)
?
0
?
?
?
0
]
T
b=[0\ \cdots \ 0 \ b^{(i)} \ b^{(j)} \ 0 \ \cdots \ 0]^T
b=[0???0?b(i)?b(j)?0???0]T 因此,定位误差可以表示如下:
Δ
=
H
ˉ
b
=
A
b
=
[
b
(
i
)
A
1
i
+
b
(
j
)
A
1
j
b
(
i
)
A
2
i
+
b
(
j
)
A
2
i
b
(
i
)
A
3
i
+
b
(
j
)
A
3
i
b
(
i
)
A
4
i
+
b
(
j
)
A
4
j
]
T
\Delta = \bar{H}b=Ab=[b^{(i)}A_{1i}+b^{(j)}A_{1j}b^{(i)}A_{2i}+b^{(j)}A_{2i}b^{(i)}A_{3i}+b^{(j)}A_{3i}b^{(i)}A_{4i}+b^{(j)}A_{4j}]^T
Δ=Hˉb=Ab=[b(i)A1i?+b(j)A1j?b(i)A2i?+b(j)A2i?b(i)A3i?+b(j)A3i?b(i)A4i?+b(j)A4j?]T 水平定位误差表示如下:
δ
x
2
=
(
b
(
i
)
A
1
i
+
b
(
j
)
A
1
j
)
2
+
(
b
(
i
)
A
2
i
+
b
(
j
)
A
2
j
)
2
\delta x^2=(b^{(i)}A_{1i}+b^{(j)}A_{1j})^2+(b^{(i)}A_{2i}+b^{(j)}A_{2j})^2
δx2=(b(i)A1i?+b(j)A1j?)2+(b(i)A2i?+b(j)A2j?)2 因此,斜率计算公式
H
s
l
o
p
e
Hslope
Hslope为:
H
s
l
o
p
e
i
=
δ
x
2
D
=
∣
∣
δ
x
∣
∣
2
∣
∣
p
∣
∣
2
=
δ
x
T
δ
x
p
T
p
=
(
b
(
i
)
A
1
i
+
b
(
j
)
A
1
j
)
2
+
(
b
(
i
)
A
2
i
+
b
(
j
)
A
2
j
)
2
b
T
P
T
P
b
Hslope_i=\frac{\delta x^2}{D}=\frac{||\delta x||^2}{||p||^2}=\frac{\delta x^T\delta x}{p^Tp}=\frac{(b^{(i)}A_{1i}+b^{(j)}A_{1j})^2+(b^{(i)}A_{2i}+b^{(j)}A_{2j})^2}{b^TP^TPb}
Hslopei?=Dδx2?=∣∣p∣∣2∣∣δx∣∣2?=pTpδxTδx?=bTPTPb(b(i)A1i?+b(j)A1j?)2+(b(i)A2i?+b(j)A2j?)2? 因此,计算
H
P
L
3
HPL_3
HPL3?的公式如下:
H
P
L
3
=
H
s
l
o
p
e
m
a
x
×
s
b
i
a
s
HPL_3=Hslope_{max}\times s_{bias}
HPL3?=Hslopemax?×sbias? 将假设1、假设2和假设3结合起来,将所有情况下的误差都包含在保护水平中,得到总保护水平的计算公式如下:
H
P
L
=
m
a
x
(
H
P
L
1
,
H
P
L
2
)
2
+
H
P
L
3
2
HPL=\sqrt{max(HPL_1,HPL_2)^2+HPL_3^2}
HPL=max(HPL1?,HPL2?)2+HPL32?
? 进而,第
m
m
m个路标点的视觉测量方程
z
m
(
t
i
)
=
T
c
p
i
x
s
 ̄
m
c
(
t
i
)
+
v
(
t
i
)
z_m(t_i)=T_c^{pix}\underline{s}_m^c(t_i)+v(t_i)
zm?(ti?)=Tcpix?s?mc?(ti?)+v(ti?)也称为伪距离定位方程。
与现有技术相比,本发明具有以下有益效果。本发明提供一种视觉导航完好性监测的计算方法。针对不同的实际情况,提出了合适的视觉定位测量模型和数学算法,丰富了导航观测数据,提高了定位精度和完好性监测性能。解决了现有技术中由于无法保证复杂环境或操作下卫星定位的可用性而导致卫星完好性算法性能不足的缺点。
3 图的简要说明
为了使本发明的内容更容易理解,下面参照本发明的具体实施例及其附图进一步详细地描述本发明。图1是本实例中视觉完好性监测的流程图。图2显示了本实例中人工目标的定位模型。图3显示了本实例中机场跑道的定位模型。图4是本实例中保护水平计算的流程图。
4 详细描述
本发明的目标和功能以及用于实现这些目标和功能的方法将参照实例加以阐明。然而,本发明不限于以下所公开的实例。它可以以不同的形式实施。本说明书的本质仅仅是帮助该领域的技术人员全面地理解本发明的具体细节。
下面将参照附图描述本发明的实例。在图纸中,相同的参考数字代表相同或相似的部分,或相同或相似的步骤。下面通过具体实例描述了一种根据本发明的视觉导航完好性监测的计算方法。
该方法通过引入视觉导航着陆系统的视觉观测来提高导航完好性监测的性能。如图1所示,是本实例中视觉完好性监测的流程图。首先,根据相机捕捉到的图像进行图像预处理,检测并保留带有着陆路标点的图像,用于视觉系统的定位。其次,为了保证定位精度,减少系统参数引起的误差,需要对误差进行修正和估计。然后,通过构建视觉定位方程,得到视觉观测模型,并计算观测矩阵,进一步计算飞机的位置和姿态信息。最后进行故障检测和排除,当检验统计量小于阈值时,系统没有发生故障,可以直接输出定位结果。否则,当检验统计量大于或等于该阈值时,表示系统发生故障,可以排除故障,并分别计算不同故障情况下的保护水平,进一步判断系统的可用性。
该方法包括以下具体步骤:
(1)根据具体的视觉定位场景建立合适的视觉定位模型,计算出视觉定位方案。
通过对飞机视觉着陆场景的研究,可以将视觉定位模型归纳为两种场景。第一种场景中,使用具有明显特征的人工目标作为着陆标志,如“H”或“I”标志,如图2所示。在第二种场景中,借用机场环境附近的机场跑道等高线,以便相机能够可靠地识别机场跑道目标,如图3所示。
如图2和图3所示,在一个导航坐标系中,根据向量相加,得到一个飞行器位置向量
p
n
p^n
pn,相机和机身之间的小角度偏移向量
q
n
q^n
qn,视线向量
s
m
n
s_m^n
smn?,它们相加以获得着陆路标向量
t
m
n
t_m^n
tmn?,如下公式被建立:
t
m
n
=
s
m
n
+
p
n
+
q
n
t_m^n=s_m^n+p^n+q^n
tmn?=smn?+pn+qn 因此,当相机系下路标点的视线向量
s
m
c
s_m^c
smc?,机体系下相机的位置
q
b
=
[
x
1
?
y
2
?
z
3
]
q^b=[x_1\ y_2\ z_3]
qb=[x1??y2??z3?]和着陆路标点的位置向量
t
m
n
=
[
x
2
?
y
2
?
z
2
]
T
t_m^n=[x_2\ y_2\ z_2]^T
tmn?=[x2??y2??z2?]T已知时,可以通过计算获得飞机在三维世界中的位置
p
n
=
[
x
3
?
y
3
?
z
3
]
T
p^n=[x_3\ y_3\ z_3]^T
pn=[x3??y3??z3?]T。其中
p
n
p^n
pn表示飞机在导航系下的位置向量;
t
m
n
t_m^n
tmn?表示第
m
m
m个路标点在导航系下的位置向量;
C
b
n
C_b^n
Cbn?表示机体系到导航系的方向余弦矩阵;
s
m
n
s_m^n
smn?表示第
m
m
m个目标在导航系中的视线向量;
s
m
b
s_m^b
smb?表示第
m
m
m个目标在机体系下的视线向量,
s
m
n
=
C
b
n
s
m
b
s_m^n=C_b^ns_m^b
smn?=Cbn?smb?;
q
n
q^n
qn表示相机光心到飞机质心之间的距离;
q
b
q^b
qb表示相机在机体系下的位置,这可以通过标定获得。
(2)修正视觉系统参数误差,并计算误差状态向量。
由于系统参数容易产生误差,因此需要对视觉参数的误差进行估计和校正,以减少视觉观测计算中线性化误差的影响。基于上述数学参数的误差状态向量可以表示为:
δ
x
=
[
β
p
n
?
δ
t
m
n
]
\delta x = [\beta p^n \ \delta t_m^n]
δx=[βpn?δtmn?] 其中
δ
p
n
\delta p^n
δpn是飞机位置的三维误差向量,可以表示为:
δ
p
n
=
p
~
n
?
p
n
\delta p^n=\tilde{p}^n-p^n
δpn=p~?n?pn
δ
t
m
n
\delta t_m^n
δtmn?是路标点位置误差,可以表示为:
δ
t
m
n
=
t
^
m
n
?
t
m
n
\delta t_m^n=\hat{t}_m^n-t_m^n
δtmn?=t^mn??tmn? 其中标称值用波浪线表示。
需要注意的是,当一个着陆路标进入和离开相机视野时,路标的数量实际上会随着时间而变化,但必须有一个合理的最大路标跟踪数量,以保持计算效率。
(3)建立视觉观测模型,计算视觉系统观测矩阵H。
如图2所示,第
m
m
m个路标点的视觉观测方程可以表示为:
z
m
(
t
i
)
=
T
c
p
i
x
s
 ̄
m
c
(
t
i
)
+
v
(
t
i
)
z_m(t_i)=T_c^{pix}\underline{s}_m^c(t_i)+v(t_i)
zm?(ti?)=Tcpix?s?mc?(ti?)+v(ti?) 该公式也称为伪距离定位方程,是基于着陆路标点和对应图像中的一个像素点得到的。其中
z
m
(
t
i
)
z_m(t_i)
zm?(ti?)是着陆路标点
m
m
m的像素位置;
T
c
p
i
x
T_c^{pix}
Tcpix?是像素坐标系和相机坐标系之间的变换矩阵;
s
m
c
s_m^c
smc?表示像素位置对应于从相机焦点向外延伸的视线向量;
s
m
z
c
s_{m_z}^c
smz?c?是相机在
s
z
s_z
sz?方向上到光心的距离;
s
 ̄
m
c
(
t
i
)
\underline{s}_m^c(t_i)
s?mc?(ti?)是相机到目标
m
m
m的视线向量的归一化形式;
v
(
t
i
)
v(t_i)
v(ti?)是独立加性高斯白噪声,其均值为0,方差为
R
R
R,它满足如下方程:
E
[
v
(
t
i
)
v
t
(
t
j
)
]
=
{
R
,
??
t
i
=
t
j
0
,
???
t
i
≠
t
j
E[v(t_i)v^t(t_j)]=\begin{cases} R, \ \ t_i=t_j \\ 0, \ \ \ t_i \neq t_j \end{cases}
E[v(ti?)vt(tj?)]={R,??ti?=tj?0,???ti??=tj??
T
p
i
x
c
T_{pix}^c
Tpixc?是已知的,设在一个高
H
H
H宽
W
W
W的矩形像素网格
(
M
×
N
)
(M\times N)
(M×N)中,像素坐标系与相机坐标的变换关系为
T
c
p
i
x
T_c^{pix}
Tcpix?,它满足如下公式:
T
p
i
x
c
=
[
?
H
f
M
0
H
(
M
+
1
)
2
f
M
0
W
f
N
?
W
(
N
+
1
)
2
f
N
0
0
1
]
=
(
T
c
p
i
x
)
?
1
T_{pix}^c=\begin{bmatrix} -\frac{H}{fM} & 0 & \frac{H(M+1)}{2fM} \\ 0 & \frac{W}{fN} & -\frac{W(N+1)}{2fN} \\ 0 & 0 & 1 \end{bmatrix}=(T_c^{pix})^{-1}
Tpixc?=?????fMH?00?0fNW?0?2fMH(M+1)??2fNW(N+1)?1?????=(Tcpix?)?1
为了便于计算,对非线性方程
z
m
(
t
i
)
=
T
c
p
i
x
s
 ̄
m
c
(
t
i
)
+
v
(
t
i
)
z_m(t_i)=T_c^{pix}\underline{s}_m^c(t_i)+v(t_i)
zm?(ti?)=Tcpix?s?mc?(ti?)+v(ti?)进行了线性化,一个非线性函数被设为
h
[
x
(
t
i
)
]
=
T
c
p
i
x
s
 ̄
m
c
(
t
i
)
h[x(t_i)]=T_c^{pix}\underline{s}_m^c(t_i)
h[x(ti?)]=Tcpix?s?mc?(ti?),计算了一阶泰勒级数展开式,得到视觉系统观测矩阵
H
H
H,其中
H
H
H中的元素由下式计算的偏导数给出:
H
=
[
?
h
?
p
n
?
0
]
H=[\frac{\partial h}{\partial p^n} \ 0]
H=[?pn?h??0] 因为,
q
n
=
C
b
n
C
c
b
s
m
c
+
C
b
n
q
d
=
C
b
n
(
C
c
b
s
m
c
+
q
d
)
q^n=C_b^nC_c^bs_m^c+C_b^nq^d=C_b^n(C_c^bs_m^c+q^d)
qn=Cbn?Ccb?smc?+Cbn?qd=Cbn?(Ccb?smc?+qd)
h
h
h相对于
p
n
p^n
pn的偏导数为:
?
h
?
p
n
=
1
s
m
z
c
T
c
p
i
x
(
s
m
c
C
n
b
C
b
c
?
C
n
b
C
b
c
)
\frac{\partial h}{\partial p^n}=\frac{1}{s_{m_z}^c}T_c^{pix}(s_m^cC_n^bC_b^c-C_n^bC_b^c)
?pn?h?=smz?c?1?Tcpix?(smc?Cnb?Cbc??Cnb?Cbc?) 因为
H
H
H与
t
m
n
t_m^n
tmn?无关,
H
H
H对
δ
t
m
n
\delta t_m^n
δtmn?的偏导数为零。
(4)引入状态向量
x
x
x的视觉观测模型,表示为:
z
=
H
x
+
v
+
b
z = Hx+v+b
z=Hx+v+b 其中
z
z
z是视觉观测向量;状态向量
x
x
x是
n
n
n维列向量;
v
v
v是
m
m
m维观测噪声向量,它服从零均值协方差为
R
R
R的高斯分布,表示为
R
=
σ
2
I
m
?
m
R=\sigma^2I_{m*m}
R=σ2Im?m?;
b
b
b是
m
m
m维故障偏差向量,在无故障情况下,
b
b
b是一个零向量,当第
i
i
i个测量发生故障时,
b
b
b中第
i
i
i个元素仍然表示为
b
b
b。
当视觉观测向量变化时,视觉观测误差表示为:
δ
z
=
H
δ
x
\delta z=H \delta x
δz=Hδx 因此,关于
z
=
H
x
+
v
+
b
z=Hx+v+b
z=Hx+v+b的最小二乘解为,
δ
x
=
(
H
T
H
)
?
1
H
T
δ
z
\delta x =(H^TH)^{-1}H^T\delta z
δx=(HTH)?1HTδz
如果不考虑噪声向量,则由状态向量
x
x
x的线性观测模型可得,
δ
x
=
x
^
?
x
=
H
ˉ
z
=
H
ˉ
b
\delta x=\hat{x}-x=\bar{H}z=\bar{H}b
δx=x^?x=Hˉz=Hˉb 其中
H
ˉ
=
(
H
T
H
)
?
1
H
T
\bar{H}=(H^TH)^{-1}H^T
Hˉ=(HTH)?1HT。
由于故障偏差向量
b
b
b对像素点
(
x
,
y
)
(x,y)
(x,y)的每个分量
x
x
x和
y
y
y都有影响,因此故障偏差向量
b
b
b在极坐标下表示为:
b
i
=
∣
∣
b
∣
∣
c
o
s
θ
b_i=||b||cos\theta
bi?=∣∣b∣∣cosθ
b
j
=
∣
∣
b
∣
∣
s
i
n
θ
b_j=||b||sin\theta
bj?=∣∣b∣∣sinθ
其中
b
i
=
∣
∣
b
∣
∣
c
o
s
θ
b_i=||b||cos\theta
bi?=∣∣b∣∣cosθ是故障偏差向量
b
b
b对分量
x
x
x的影响;
b
j
=
∣
∣
b
∣
∣
s
i
n
θ
b_j=||b||sin\theta
bj?=∣∣b∣∣sinθ是故障偏差向量
b
b
b对分量
y
y
y的影响。
因此,位置误差的范数表示为:
∣
∣
δ
x
∣
∣
2
=
(
H
ˉ
b
)
T
H
ˉ
b
||\delta x||^2=(\bar{H}b)^T\bar{H}b
∣∣δx∣∣2=(Hˉb)THˉb 依据如下公式:
(
H
ˉ
b
)
T
H
ˉ
b
=
b
T
H
ˉ
T
H
ˉ
b
=
b
i
b
i
(
H
ˉ
T
H
ˉ
)
i
i
+
b
i
b
j
(
H
ˉ
T
H
ˉ
)
i
j
+
b
i
b
j
(
H
T
H
)
j
i
+
b
j
b
j
(
H
ˉ
T
H
ˉ
)
j
j
(\bar{H}b)^T\bar{H}b=b^T\bar{H}^T\bar{H}b=b_ib_i(\bar{H}^T\bar{H})_{ii}+b_ib_j(\bar{H}^T\bar{H})_{ij}+b_ib_j(H^TH)_{ji}+b_jb_j(\bar{H}^T\bar{H})_{jj}
(Hˉb)THˉb=bTHˉTHˉb=bi?bi?(HˉTHˉ)ii?+bi?bj?(HˉTHˉ)ij?+bi?bj?(HTH)ji?+bj?bj?(HˉTHˉ)jj?
将其代入极坐标公式,可得到故障偏差向量
b
b
b对水平位置误差的影响公式,具体公式如下:
δ
x
=
(
∣
∣
b
∣
∣
2
(
s
i
n
2
θ
(
H
ˉ
T
H
ˉ
)
i
i
+
s
i
n
2
θ
(
H
ˉ
T
H
ˉ
)
i
j
+
c
o
s
2
θ
(
H
ˉ
T
H
ˉ
)
j
j
)
)
1
2
\delta x=(||b||^2(sin^2\theta(\bar{H}^T\bar{H})_{ii}+sin2\theta(\bar{H}^T\bar{H})_{ij}+cos^2\theta(\bar{H}^T\bar{H})_{jj}))^{\frac{1}{2}}
δx=(∣∣b∣∣2(sin2θ(HˉTHˉ)ii?+sin2θ(HˉTHˉ)ij?+cos2θ(HˉTHˉ)jj?))21? 因此,当角度值
θ
\theta
θ发生变化时,水平位置误差也会发生相应的变化,特别是当取一些特殊的角度值时,水平位置误差
δ
x
\delta x
δx公式可以简化。
(5)结合视觉系统观测方程和卫星定位方程,计算组合导航定位方程。
假设卫星定位方程为
z
=
G
x
z=Gx
z=Gx,合并视觉系统观测矩阵
H
H
H和卫星观测矩阵
G
G
G可以获得扩展观测矩阵
H
′
H'
H′,表示如下:
H
′
=
[
H
G
]
H'=\begin{bmatrix} H \\ G \end{bmatrix}
H′=[HG?] 其中
H
′
H'
H′是
m
×
n
m\times n
m×n的线性观测矩阵。
当仍然有有限数量的卫星可用时,可以利用有限数量的卫星观测来加强视觉导航设备。如果视觉系统是由两个或三个相机组成的系统,则视觉系统中单个相头传感器的变换矩阵
M
1
M_1
M1?、
M
2
M_2
M2?和
M
i
M_i
Mi?可以组合为
M
M
M,
M
M
M可以具体表示为:
M
=
[
M
1
0
0
0
M
2
0
0
0
M
i
]
M=\begin{bmatrix} M_1 & 0 & 0\\ 0 & M_2 & 0 \\ 0 & 0 & M_i \end{bmatrix}
M=???M1?00?0M2?0?00Mi?????
(6)在不同的故障模式下计算保护水平。
保护水平的计算方法分为三个假设:假设1、假设2和假设3。
假设1为观测无故障的情况,此时保护水平可根据无故障放大因子
K
f
f
K_{ff}
Kff?和水平位置标准差
σ
x
\sigma_x
σx?计算,其中无故障放大因子
K
f
f
K_{ff}
Kff?由无故障概率
P
f
f
P_{ff}
Pff?确定。因此,无故障保护水平的计算公式为:
H
P
L
1
=
K
f
f
(
k
)
?
σ
x
HPL_1=K_{ff}(k)*\sigma_x
HPL1?=Kff?(k)?σx? 其中
H
P
L
1
HPL_1
HPL1?表示无故障情况下的保护水平;
K
f
f
K_{ff}
Kff?是无故障放大因子;
σ
x
\sigma_x
σx?是飞机在水平方向上的标准差。
假设2是观测中出现单故障的情况,在这种情况下,单个视觉观测或GNSS观测都有可能发生故障。当只有一个观测故障时,故障偏差向量
b
b
b中只包含一个偏差分量,此时水平位置误差
δ
x
\delta x
δx可表示为:
δ
x
=
x
^
?
x
=
H
ˉ
z
=
H
ˉ
b
\delta x=\hat{x}-x=\bar{H}z=\bar{H}b
δx=x^?x=Hˉz=Hˉb 其中
x
^
\hat{x}
x^表示最小二乘估计值,
H
ˉ
=
(
H
T
H
)
?
1
H
T
\bar{H}=(H^TH)^{-1}H^T
Hˉ=(HTH)?1HT。每个观测的斜率设为
H
s
l
o
p
e
Hslope
Hslope,表示为:
H
s
l
o
p
e
i
=
δ
x
2
D
=
∣
∣
δ
x
∣
∣
2
∣
∣
p
∣
∣
2
=
δ
x
T
δ
x
p
T
p
=
b
T
H
ˉ
T
H
ˉ
b
b
T
P
T
P
b
Hslope_i=\frac{\delta x^2}{D}=\frac{||\delta x||^2}{||p||^2}=\frac{\delta x^T \delta x}{p^Tp}=\frac{b^T\bar{H}^T\bar{H}b}{b^TP^TPb}
Hslopei?=Dδx2?=∣∣p∣∣2∣∣δx∣∣2?=pTpδxTδx?=bTPTPbbTHˉTHˉb? 其中
b
T
H
ˉ
T
H
ˉ
b
=
b
i
b
i
(
H
ˉ
T
H
ˉ
)
i
i
+
b
i
b
j
(
H
ˉ
T
H
ˉ
)
i
j
+
b
i
b
j
(
H
ˉ
T
H
ˉ
)
j
i
+
b
j
b
j
(
H
ˉ
T
H
ˉ
)
j
j
b^T\bar{H}^T\bar{H}b=b_ib_i(\bar{H}^T\bar{H})_{ii}+b_ib_j(\bar{H}^T\bar{H})_{ij}+b_ib_j(\bar{H}^T\bar{H})_{ji}+b_jb_j(\bar{H}^T\bar{H})_{jj}
bTHˉTHˉb=bi?bi?(HˉTHˉ)ii?+bi?bj?(HˉTHˉ)ij?+bi?bj?(HˉTHˉ)ji?+bj?bj?(HˉTHˉ)jj?;
D
D
D是统计检验量;
p
p
p是奇偶向量;
P
P
P是奇偶变换矩阵,
p
=
P
b
p=Pb
p=Pb。因此,水平定位误差保护水平
H
P
L
2
HPL_2
HPL2?为:
H
P
L
2
=
H
s
l
o
p
e
m
a
x
×
s
b
i
a
s
HPL_2=Hslope_{max}\times s_{bias}
HPL2?=Hslopemax?×sbias? 其中
s
b
i
a
s
=
σ
λ
n
?
4
s_{bias}=\sigma \frac{\sqrt{\lambda}}{\sqrt{n-4}}
sbias?=σn?4
?λ
??
σ
\sigma
σ是噪声方差,它是一个分布参数。当无故障时,
H
P
L
2
HPL_2
HPL2?的值最小;当观测中出现单故障时,
H
P
L
2
HPL_2
HPL2?的值会增加。
假设3为观测出现两个故障的情况,出现两个视觉观测故障、两个GNSS观测故障或一个视觉观测和一个GNSS观测故障的情况。当两个观测故障时,第
i
i
i个观测与第
j
j
j个观测同时故障的故障偏差向量为
b
(
i
)
b(i)
b(i)和
b
(
j
)
b(j)
b(j),可表示为:
b
=
[
0
?
?
?
0
?
b
(
i
)
?
b
(
j
)
?
?
?
0
?
0
?
?
0
]
T
b=[0\ \cdots \ 0\ b^{(i)} \ b^{(j)} \ \cdots \ 0\ 0 \cdots \ 0] ^T
b=[0???0?b(i)?b(j)???0?0??0]T 因此,定位误差表示如下:
Δ
=
H
ˉ
b
=
A
b
=
[
b
(
i
)
A
1
i
+
b
(
j
)
A
1
j
b
(
i
)
A
2
i
+
b
(
j
)
A
2
i
b
(
i
)
A
3
i
+
b
(
j
)
A
3
i
b
(
i
)
A
4
i
+
b
(
j
)
A
4
j
]
T
\Delta=\bar{H}b=Ab=[b^{(i)}A_{1i}+b^{(j)}A_{1j}b^{(i)}A_{2i}+b^{(j)}A_{2i}b^{(i)}A_{3i}+b^{(j)}A_{3i}b^{(i)}A_{4i}+b^{(j)}A_{4j}]^T
Δ=Hˉb=Ab=[b(i)A1i?+b(j)A1j?b(i)A2i?+b(j)A2i?b(i)A3i?+b(j)A3i?b(i)A4i?+b(j)A4j?]T 水平定位误差表示如下:
δ
x
2
=
(
b
(
i
)
A
1
i
+
b
(
j
)
A
1
j
)
2
+
(
b
(
i
)
A
2
i
+
b
(
j
)
A
2
j
)
2
\delta x^2=(b^{(i)}A_{1i}+b^{(j)}A_{1j})^2+(b^{(i)}A_{2i}+b^{(j)}A_{2j})^2
δx2=(b(i)A1i?+b(j)A1j?)2+(b(i)A2i?+b(j)A2j?)2 因此,计算
H
s
l
o
p
e
Hslope
Hslope如下:
H
s
l
o
p
e
i
=
δ
x
2
D
=
∣
∣
δ
x
∣
∣
2
∣
∣
p
∣
∣
2
=
δ
x
T
δ
x
p
T
p
=
(
b
(
i
)
A
1
i
+
b
(
j
)
A
1
j
)
2
+
(
b
(
i
)
A
2
i
+
b
(
j
)
A
2
j
)
2
b
T
P
T
P
b
Hslope_i=\frac{\delta x^2}{D}=\frac{||\delta x||^2}{||p||^2}=\frac{\delta x^T \delta x}{p^Tp}=\frac{(b^{(i)}A_{1i}+b^{(j)}A_{1j})^2+(b^{(i)}A_{2i}+b^{(j)}A_{2j})^2}{b^TP^TPb}
Hslopei?=Dδx2?=∣∣p∣∣2∣∣δx∣∣2?=pTpδxTδx?=bTPTPb(b(i)A1i?+b(j)A1j?)2+(b(i)A2i?+b(j)A2j?)2? 因此,计算
H
P
L
3
HPL_3
HPL3?的公式如下:
H
P
L
3
=
H
s
l
o
p
e
m
a
x
×
s
b
i
a
s
HPL_3=Hslope_{max}\times s_{bias}
HPL3?=Hslopemax?×sbias?
结合假设1、假设2和假设3,所有情况下的误差在保护水平中都被考虑,计算保护水平的公式如下,
H
P
L
=
m
a
x
(
H
P
L
1
,
H
P
L
2
)
2
+
H
P
L
3
2
HPL=\sqrt{max(HPL_1,HPL_2)^2+HPL_3^2}
HPL=max(HPL1?,HPL2?)2+HPL32?
?
在所有情况下,保护水平的值越准确,在给定风险值的情况下,用户越能避免所有偏差,系统可用性也就越高。根据本实例,通过所提出的视觉导航完好性监测算法,可以有效地排除故障,从而提高导航定位的准确性和可用性。上述实例仅仅是用于清楚说明的示例,而不是对实现方式的限制。那些具有一般技能的人也可以在上述描述的基础上做出其它不同形式的改变。没有必要也没有办法用尽所有的实现方式。其发生的变化或变更,仍在本发明的保护范围内。
本发明的其它实例可以很容易地结合本发明的描述和实践而被本领域中的技术人员理解。本描述和实例仅被认为是示例性的,本发明的真正范围和主题由权利要求书定义。
|