基于UKF和IUFIR的机动目标跟踪方法 |
|||||||
申请号 | CN201610154301.9 | 申请日 | 2016-03-17 | 公开(公告)号 | CN105891820A | 公开(公告)日 | 2016-08-24 |
申请人 | 四川长虹电器股份有限公司; | 发明人 | 何娟; 冯霞; 王必娇; | ||||
摘要 | 本 发明 公开了一种基于UKF和IUFIR的机动目标 跟踪 方法,它包括以下步骤:建立系统模型,设定目标的 状态方程 和测量方程;对目标初始状态和误差协方差进行初始化;使用UKF滤波方法进行线性转换;使用IUFIR滤波方法对线性化的数据进行滤波初始化,得到s时刻的广义噪声功率增益分别为和Gs;更新IUFIR滤波方法中的广义噪声功率增益G;对IUFIR滤波增益K更新;使用IUFIR滤波方法进行滤波估计,本发明解决在非线性雷达目标跟踪系统中,状态估计方面因为过程噪声的不确定性而增加预测估计的误差甚至发散的问题。 | ||||||
权利要求 | 1.一种基于UKF和IUFIR的机动目标跟踪方法,其特征在于,它包括以下步骤: |
||||||
说明书全文 | 基于UKF和IUFIR的机动目标跟踪方法技术领域[0001] 本发明涉及空间机动目标跟踪技术领域,具体涉及基于UKF和IUFIR的机动目标跟踪方法。 背景技术[0002] 空间机动目标跟踪技术很早就被提出来了,但是它一直是雷达目标跟踪中的难点,随着数据关联等技术在多目标跟踪理论方面取得了开创性的突破,由数据关联技术和滤波技术相结合提出了很多新的目标跟踪算法。但是由于机动目标的机动情况时未知的,这就增加了状态方程的过程噪声难以确定,一般模型的过程噪声均为高斯白噪声,后来也有一些方法来解决这个问题,例如可调白噪声方法、变维滤波方法等,这些方法虽然通过各种方法来降低过程噪声的不确定性来提高跟踪的准确度,但是没有从根本上解决这个问题。 发明内容[0003] 本发明克服了现有技术的不足,提供一种基于UKF和IUFIR的机动目标跟踪方法。 [0004] 为实现上述目的,本发明采用以下技术方案: [0005] 一种基于UKF和IUFIR的机动目标跟踪方法,它包括以下步骤: [0006] 步骤1、建立系统模型,设定目标的状态方程和测量方程; [0007] 步骤2、对目标初始状态和误差协方差进行初始化; [0008] 步骤3、使用UKF滤波方法进行线性转换; [0009] 步骤4、使用IUFIR滤波方法对线性化的数据进行滤波初始化,得到s时刻的广义噪声功率增益分别为 和Gs; [0010] 步骤5、更新IUFIR滤波方法中的广义噪声功率增益G; [0011] 步骤6、对IUFIR滤波增益K更新; [0012] 步骤7、使用IUFIR滤波方法进行滤波估计。 [0013] 更进一步的技术方案是,所述步骤1中建立的非线性跟踪系统模型中状态方程和测量方程分别为: [0014] xn=f(xn-1,θn-1)+wn-1 [0015] zn=hn(xn,θn)+vn [0016] 其中,xn为目标状态向量,f为状态转移非线性函数,wn-1为过程噪声,其协方差为Qn,h为量测函数,vn为测量噪声,zn为系统测量向量。 [0017] 更进一步的技术方案是,所述步骤3使用UKF滤波中的无迹转换进行线性转换得到离散的系统模型为: [0018] xn=Fnxn-1+Bnwn [0019] zn=Hnxn+vn [0020] xn∈RK×1为n时刻目标状态向量,Fn为n时刻非线性状态函数,wn∈Rp×1为n时刻过程噪声,vn∈Rp×1为n时刻量测噪声,zn∈Rm×1为n时刻系统的量测向量,Hn为n时刻非线性量测函数。 [0021] 更进一步的技术方案是,所述步骤4使用IUFIR滤波方法对线性化的数据进行滤波初始化,得到s时刻的广义噪声功率增益分别为 和Gs,得到: [0022] [0023] [0024] 其中, [0025] [0026] [0027] [0028] 更进一步的技术方案是,所述滤波增益和状态估计,具体公式表示为: [0029] [0030] [0031] 本技术方案中,无迹卡尔曼滤波方法中的无迹转换过程具体如下,无迹变换是对状态概率密度函数(Probability Density Function,PDF)做近似。它首先选择有限个近似高斯分布离散点(称为σ点),对每个σ点施以非线性变换,得到一簇UT变换后的线性点,再使用UFIR滤波方法来对转换的线性离散系统进行线性滤波,实现目标的定位和跟踪。 [0032] (1)UT变换步骤为: [0033] 设 为一个随机变量,g代表任意一个非线性函数,实现 令y=g(x), Cov(x,x)=Px。UT变换的具体过程如下: [0034] 先按照以下方法得到(2nx+1)个带权值的粒子,mi为粒子权值,取(i=1…2nx),选取原则如下: [0035] [0036] [0037] 其中,λ为尺度参数,用于调节δ点和先验均值的距离远近,一般取0或3-nx;是指矩阵平方根的第i行或第i列。 [0038] 这些δ点通过非线性函数得到,公式如下: [0039] ξi=g(χi) i=0,1,…2nx (3) [0040] y的均值和方差可通过这些的线性化的δ点得出: [0041] [0042] [0043] (2)IUFIR滤波方法 [0044] IUFIR滤波方法主要针对线性滤波模型,假设离散线性模型系统的状态方程和观测方程如下: [0045] xn=Fnxn-1+Bnwn (6) [0046] zn=Hnxn+vn (7) [0047] xn∈RK×1为n时刻目标状态向量,Fn为n时刻非线性状态函数,wn∈Rp×1为n时刻过程噪声,vn∈Rp×1为n时刻量测噪声,zn∈Rm×1为n时刻系统的量测向量,Hn为n时刻非线性量测函数。 [0048] 设在系统中的vn、wn为零均值噪声序列,它们互不相关,其协方差分别为:Qn和Rn。对所有的i、j时刻,噪声的基本统计性 质为 :E{vn }=0,E{wn }=0,[0049] 当利用IUFIR滤波方法进行滤波时,对多个时刻抽样值取平均值进行滤波(即批处理IUFIR滤波),此特性能够在忽视噪声统计特性的情况下对信号进行滤波。具体步骤如下: [0050] 初始化:通过批处理IUFIR滤波得到初始的s时刻的目标状态初值和广义噪声功率增益分别为 和Gs,此时的批处理滤波器窗长为N。 [0051] [0052] [0053] 其中: [0054] [0055] [0056] [0057] IUFIR滤波的迭代更新思想为:此时批处理的窗长发生变化为M,且M□N,其中s=m+M-1,当迭代变量从m+M到n,当l=n时,就得到n时刻的状态预测估计 [0058] [0059] 广义噪声功率增益的更新:在机动目标跟踪中,机动会导致量测数据和滤波的数据产生很大的偏差,此偏差能够目标反映机动的情况。但是广义噪声功率增益G只与状态转移矩阵和观测矩阵有关系,并且不同时刻的G之间是相互独立的,所以本发明IUFIR滤波方法通过调整系数γ来自适应调整G的值,提高滤波精度。 [0060] 在i时刻,选取量测与滤波结果的偏差的平方根ηi来描述机动的基本信息,用公式表示为: [0061] [0062] 其中,zi为第i时刻的量测数据,为第i时刻的状态预测结果,Hi为第i时刻的状态观测矩阵,κ为目标运动的维数。 [0063] 由于目标的机动情况只能通过不同时刻的偏差的比值才能观察出来,所以本发明将此时刻的偏差平方根与上一时刻的偏差的平方根的比值来反映机动情况,第i时刻与第i-1时刻的偏差的平方根比值用数学公式表示如下: [0064] [0065] 同时考虑到野值可能会造成滤波发散的原因,以及一般机动会持续几秒到几十秒,所以G的调整系数γ取比值λ的均值的方法来消除野值影响,同时也不影响对目标机动的描述,在第i时刻得到G得调整系数γ为: [0066] [0067] 其中, 由于必须满足l0≥m+M+1,所以得到: 说明对广义噪声功率增益的调整只能从 时刻开始,第l时刻的广义噪声功率增益的计算表达式为: [0068] [0069] IUFIR滤波方法的滤波增益为: [0070] [0071] IUFIR滤波方法的l时刻的状态预测估计为: [0072] [0073] 当l=n时,得到第n时刻的滤波状态估计为 [0074] 与现有技术相比,本发明的有益效果是: [0077] 图2为本发明一种实施例中跟踪路径和目标实际运行的路径的对比图。 具体实施方式[0078] 下面结合附图对本发明作进一步阐述。 [0079] 如图1所示,一种基于UKF和IUFIR的机动目标跟踪方法,首选引入涉及本实施例的相关参数,并详细描述如下: [0080] g代表任意一个非线性函数,y=g(x); [0081] E(x)为x的期望值,Cov(x,x)=Px,Px为协方差; [0082] mi为粒子权值,λ为尺度参数,和Py为y的均值和方差: [0083] 在k-1时刻的线性化状态和均方差为 Pk-1|k-1; [0084] wm、wc为权值; [0085] xn∈RK×1为n时刻目标状态向量; [0086] Fn为n时刻非线性状态函数; [0087] wn∈Rp×1为n时刻过程噪声,vn∈Rp×1为n时刻量测噪声; [0088] zn∈Rm×1为n时刻系统的量测向量; [0089] Hn为n时刻非线性量测函数; [0090] Qn和Rn分别为噪声序列vn、wn的协方差分别; [0091] IUFIR滤波方法中s时刻的目标状态初值和广义噪声功率增益初始值分别为 和Gs; [0092] 批处理的窗长M,且M□N,其中s=m+M-1; [0093] n时刻的状态预测估计 [0094] 在i时刻的量测与滤波结果的偏差的平方根ηi; [0095] zi为第i时刻的量测数据, 为第i时刻的状态预测结果,Hi为第i时刻的状态观测矩阵;κ为目标运动的维数; [0096] λi为第i时刻与第i-1时刻的偏差的平方根比值; [0097] γ为广义噪声功率增益G的调整系数; [0098] Gl为第l时刻的广义噪声功率增益; [0099] Kl为IUFIR滤波方法的滤波增益; [0100] 为IUFIR滤波方法的l时刻的状态预测估计。 [0101] 本发明模拟的环境参数为实验仿真参数:目标在雷达直角坐标系中初始位置为[10km,10km,10km],初始速度为[100m/s,100m/s,100m/s],加入随机噪声,仿真帧数为100,采样周期为T=2s,蒙特卡罗采样次数位100。 [0102] 步骤1、建立系统模型,设定本方法的目标的状态方程和测量方程,雷达测量系统为非线性的,本发明中假设非线性模型为: [0103] xn=f(xn-1,θn-1)+wn-1 (19) [0104] zn=hn(xn,θn)+vn (20) [0105] 步骤2、对目标初始状态和误差协方差进行初始化,如图1的步骤101对无迹卡尔曼滤波初始化:k=0,得到如下公式: [0106] [0107] cov(x0)=P0 (22) [0108] 步骤3、使用UKF滤波方法进行线性转换,如图1的步骤102、103,假设在k-1时刻的状态为 pk-1|k-1,首先计算在k-1时刻采样的样本点集σ点(i=1,…2n): [0109] [0110] [0111] [0112] 通过计算样本点集σ点中的 i=1,2…,2n来计算 pk|k-1,计算如下: [0113] [0114] [0115] [0116] 其中wm、wc表达式如下: [0117] wm=λ/(n+λ) (29) [0118] wc=λ(n+λ)+(1-α2+β) (30) [0119] 等式(29)中的β表示状态的分布,对于高斯分布β=2是最优的,如果状态变量是单变量,则最佳的选择是β=0。 [0120] 步骤4,使用改进的IUFIR滤波方法对线性化的数据进行滤波初始化,如图1的步骤104,得到s时刻的广义噪声功率增益分别为 和Gs,经UKF的UT转换线性化的系统模型如式(6)和式(7),此时的批处理滤波器窗长为N=55。 [0121] [0122] [0123] [0124] [0125] 步骤5、更新IUFIR滤波方法中的广义噪声功率增益G,如图1的步骤105,此时批处理的窗长发生变化为M=5,且M□N,其中s=m-M-1,当迭代变量从m+M到n,当l=n时,得到Gl为 [0126] 步骤6、如图1的步骤106对IUFIR滤波增益K更新,得到: [0127] 步骤7、如图1的步骤107,使用IUFIR滤波方法进行滤波估计,得到l时刻的滤波估计为: [0128] 通过上述给定的参数,利用本发明去跟踪目标,虽然计算量增大了,但是滤波精度大大提高,大大降低了漏跟目标的概率,跟踪路径和目标实际运行的路径对比结果如图2所示,表明本发明能稳定精确的跟踪目标。 [0129] 以上具体实施方式对本发明的实质进行详细说明,但并不能对本发明的保护范围进行限制,显而易见地,在本发明的启示下,本技术领域普通技术人员还可以进行许多改进和修饰,需要注意的是,这些改进和修饰都落在本发明的权利要求保护范围之内。 |