为接收阵中心在地心直角坐标系X轴方向的位置，Ζ Λ为接收 阵中心在地心直角坐标系Z轴方向的位置；当k = 1时，获得采样时刻的位置参数X(l)= [X⑴，y⑴，z⑴]τ，当k = 2时，获得采样时刻的位置参数X (2) = [X⑵，y⑵，z⑵]T;k 取自然数；
[0057] b)对前两个采样时刻的位置参数进行差分计算目标在第2个采样时刻的近似速 度
[0059] 其中，T表不米样周期；
[0060] c)根据上述计算结果构造滤波初值，
[0061] X(0) .. [x<2}> ^ 2), (1 /(2), K2), 0;^(2}; z{ 2), OJ'
[0062] P(0) = 106Ig
[0063] 其中，19表示9维的单位矩阵。
[0064] 进一步地，在步骤9中，所述完成对空间机动目标跟踪的具体步骤为：
[0065] a)状态向量采样
[0066] 对状态估计向量、、、!?.进行sigma采样，k取自然数，X p(k_l)表示第p个状态 采样点，P = 〇，…，2nx，nx=9;采样点为：
[0068] 其中，κ为比例参数：：表示矩阵
的第i列，i = 1，2，···，Πχ
[0069] b)时间更新
[0070] 利用CS模型对目标状态采样点进行一步预测得到x p (k I k-i)，一步预测状态转移 矩阵为φa，根据目标测元与目标运动状态之间的映射函数h (·)计算观测一步预测的采样 点Tp(k| k-i)，再对其加权求和，计算观测量的一步预测fU…丨）；
[0071] xp(k|k-l)表示第ρ个状态采样点从k-i时刻到k时刻的预测值；T p(k|k-1)表 示第p个观测一步预测的采样点；表示第p个采样点的均值权重为第k个 采样时刻的观测量的一步预测；
[0078] G(k)表示第k个采样时刻目标状态更新增益矩阵；p(k|k-l)表示目标状态协方差 矩阵的一步预测；Pyy(k|k-1)表示观测协方差矩阵的一步预测；
分 别表示第k-ι到第k个采样时刻的目标状态一步预测和测元一步预测的均值；GT (k)表示第 k个采样时刻的目标状态更新增益矩阵的转置。
[0079] 有益效果：
[0080] (1)本发明采用发射机、接收机均为均匀圆阵的配置，不但能够估计目标的四个角 度和一个多普勒频率信息，还能够对空间目标进行三维定位。
[0081] (2)在进行目标运动状态参数解算时，传统的几何解算方法只能利用角度测元，无 法利用多普勒频率测元，造成了测元信息浪费，而本发明采用了基于CS模型的UKF算法估 计目标的运动状态参数，一方面提高了对测量噪声的抑制能力，另一方面能够融合利用所 有测元信息，这两方面共同提高了目标运动状态参数的解算精度。
[0082] (3)本发明采用了适用于机动目标跟踪的CS模型，基于此可以对空间机动目标的 运动状态进行一定时间的预测，实现对空间机动目标的有效稳定跟踪。
【附图说明】
[0083] 图1为本发明的一种双基地ΜΜ0雷达空间机动目标跟踪方法实现流程图
[0084] 图2为本发明的双基地Μ頂0雷达配置示意图。
[0085] 图3为用本发明方法仿真信噪比20dB情况下某临近空间运动目标相对于发射机 的俯仰角和方位角的估计误差图。
[0086] 图4为用本发明方法仿真信噪比20dB情况下某临近空间运动目标相对于发射机 的俯仰角和方位角的估计误差图。
[0087] 图5为用本发明方法仿真信噪比20dB情况下某临近空间运动目标多普勒频率的 估计误差图。
[0088] 图6为用本发明方法仿真信噪比为20dB时的目标位置滤波误差图。
[0089] 图7为用本发明方法仿真信噪比为20dB时的目标速度滤波误差图。
【具体实施方式】
[0090] 下面结合附图并举实施例，对本发明进行详细描述。
[0091] 图1为本发明的一种双基地ΜΙΜΟ雷达空间机动目标跟踪方法实现流程图，实现本 发明的具体步骤如下：
[0092] 步骤1，分别将双基地ΜΜ0雷达的发射机配置为Μ个阵元的均匀圆阵、接收机配置 为Ν个阵元的均匀圆阵，并使发射机中Μ个阵元发射相互正交的波形信号；其中，Μ表示发 射机阵元个数，Ν表示接收机阵元个数，且Μ、Ν均为自然数；
[0093] 图2为本发明的双基地ΜΜ0雷达配置示意图。坐标系采用地心直角坐标系 ，原点％为地球中心，0 轴的正向为沿起始天文子午面与地球赤道平面的交线指 向外方向，0A轴沿地球自转轴指向北极，0CYC轴与0CXC、0 CZC轴构成右手直角坐标系。其 中，图2中(\为发射圆阵的中心点，01^为接收圆阵的中心点。
[0094] 发射圆阵与接收圆阵均布设在当地水平面内。设发射圆阵的圆点0,的大地坐 标（经度Lt，炜度Bt，高程ht)为（Lt，B t，ht)，接收圆阵的圆点4的大地坐标（经度Lp炜 度民，高程1〇为（1^，81^，1〇。发射圆阵半径为心其阵元个数1=2打〇证(2 3^/\)+1， fl〇〇r( ·)表示向下取整运算，λ表示发射波的波长；接收圆阵半径为&，其阵元个数N = 2f loor (2 n rr/ λ ) +1 〇
[0095] 假设空间有一目标P，在地心系下的位置坐标为（x，y，z)、速度坐标为1????发 射圆阵的圆心与目标的连线〇tP与发射圆阵平面的夹角为俯仰角，记为〇t，〇tP在发射圆阵 平面的投影与当地正北方向的夹角为方位角，记为0t;接收圆阵的圆心至目标连线〇^与 接收圆阵平面的夹角为俯仰角，记为Φρ 〇^在接收圆阵平面的投影与当地正北方向的夹 角为方位角，记为Θ」由目标运动引起的接收信号多普勒频率为fd。
[0096] 步骤2,利用所述发射机中Μ个阵元发射相互正交的波形信号，接收机中的N个阵 元分别接收所述发射机中Μ个阵元发射的相互正交的波形信号，并进行匹配滤波，获得匹 配滤波后的匪X 1维雷达回波信号s和L次快拍积累得到的匪XL维雷达回波信号矩阵S， 进而获得矩阵S中第η个接收阵元的MX L维切片矩阵形式Sn，其中，L表示快拍次数，η为 自然数，η = 1，2, · · ·，Ν ;
[0097] 具体地，利用发射机中Μ个阵元发射相互正交的波形信号，接收机中的Ν个阵元分 别接收该发射机中Μ个阵元的发射信号，并进行匹配滤波，依次获得匹配滤波后的ΝΜΧ 1维 雷达回波信号s，考察第1个脉冲匹配滤波后的回波信号s (1)，其表达式为
[0099] 其中，Ρ为目标雷达波反射系数，假定目标对不同正交发射信号的反射系数在相 同脉冲内是完全相关的，在不同脉冲间是相互独立的，Φ)为目标接收方向向量， at(0t，Φ,)为目标发射方向向量，#表示Kronecker积，爾i}为第1个脉冲匹配滤波后的 接收噪声，1 = 1，2，…乂；
[0101] 其中，β m为第η个接收阵元的方位角，β m= 2 31 (η-1)/Ν，η = 1，2，···，N，符号 ?急表示"简记为";(^)表示第n个接收阵元的接收方向因子。
[0103] 其中，β tm为第m个发射阵元的方位角，β tm= 2 Jr (m_l)/M，m = 1，2，…，Μ ; atni(0t，（i>t)表示第m个发射阵元的的发射方向因子。
[0104] 经过L次快拍积累得到匪X L维矩阵S，在L次快拍内认为目标运动状态参数不 变，S的表达式为
[0106] 其中：
、上标T表示转置。兔 表示由.3(/..)形成的矩阵，Β表示归一化多普勒频率方向向量。
[0107] 进而得到L次快拍积累中第η个接收阵元的MX L维切片矩阵Sn，其表达式为：
[0108] Sn= A R (n) ATBT+Wn，η = 1，2, · · ·，N
[0109] 其中，AR(η)表示向量^的第η个元素，Wn表示第η个接收阵元L次快拍积累中的 噪声矩阵;Ατ表示发射方向向量，A R表示接收方向向量。
[0110] 步骤3,根据N个接收阵元的MXL维切片矩阵形SSi~SN，利用平行因子算法分别 获得发射方向向量的估计接收方向向量的估计,?和多普勒频率方向向量的估计
[0111] 步骤3的具体子步骤为：
[0112] 3a)第η个接收阵元的MX L维切片矩阵Sn，根据N个接收阵元的MX L维切片矩 阵Si~S N，获得MX NX L的三维数据集，进而得到第m个发射阵元的第η个接收阵元的第1 次快拍的平行因子三线性模型形式，其具体表达式为
[0113] Sm>njl= A T(m)AR(n)B(l)+wmjnjl, m e {1, e {1, {1, ···,〇,
[0114] 其中，AT(m)表示发射方向向量AT的第m个元素，A R(n)表示接收方向向量AR的第 η个元素，B(l)表示多普勒频率方向向量B的第1个元素 ，表示三维噪声数据集；
[0115] 3b)根据Sni, n, i，分别得到第m个发射阵元的LXN维切片矩阵形式F"、第1次快拍 的NXM维切片矩阵形式Zi，进而分别获得Μ个发射阵元的LXN维切片矩阵形式匕~F M，L 次快拍的NXM维切片矩阵形式Zi~Z p
[0116] 具体地，根据第m个发射阵元的第η接收阵元的第1次快拍的平行因子三线性模 型形式Sm, n, i的对称性，对其第二维和第三维分别进行切片，分别获得第m个发射阵元的 LXN维切片矩阵形式Fm、第1次快拍的NXM维切片矩阵形式Ζ:，其表达式分别为：
[0119] 3c)将N个接收阵元的MXL维切片矩阵形KSi~SN按列形成匪XL维矩阵S，将 Μ个发射阵元的L X N维切片矩阵形式Fi~F M按列形成ML X N维矩阵F，将L次快拍的NX Μ 维切片矩阵形式Zi~Zj安列形成NLXM维矩阵Ζ，具体地，
[0123] 其中，。表示Khatri-Rao积，Ws表示NMXL维矩阵S的噪声矩阵，W F表示MLXN 维矩阵F的噪声矩阵