本发明提出一种雷达组网系统的空间配准方法，主要解决现有雷达组网系统信号融合时空间配准的能量损失较大、配准效果不佳问题。其技术方案是：1)设定监视区域，选定监视区域中心点作为初始栅格单元的中心点；2)根据雷达组网系统内各雷达波束指向选出照射到监视区域的雷达；3)计算所述各雷达在所述中心点处形成的初始空间分辨单元信息；4)计算所述各雷达初始空间分辨单元的融合区域的协方差矩阵；5)计算初始栅格单元的大小和对应的旋转矩阵；6)根据初始栅格单元对监视区域进行栅格划分；7)将划分出的网格单元分别对应到各雷达的接收通道和距离单元，完成空间配准。本发明能量损失小，性能优异，可用于目标检测中的信号融合。
1.一种雷达组网系统的空间配准方法，所述雷达组网系统内包含多部雷达，其实现步骤包括如下：
1)设定雷达组网系统的监视区域，选定监视区域中心点作为初始栅格单元的中心点；
2)根据雷达组网系统内各雷达波束指向选出照射到监视区域的雷达；
3)计算雷达组网系统内各雷达在监视区域中心点处形成的初始空间分辨单元信息以及初始栅格单元中心点相对于指向监视区域的各雷达的方位角θi和俯仰角 其中1≤i≤N,N表示所述指向监视区域的雷达总数；
4)将每个雷达的初始空间分辨单元近似为椭球形状，计算该椭球区域的协方差矩阵Ci，其中1≤i≤N,N表示所述指向监视区域的雷达总数；
5a)根据各椭球区域的协方差矩阵Ci，计算各椭球相交区域的协方差矩阵C：其中，(·)-1表示矩阵求逆运算，ki表示所述指向监视区域的第i部雷达对应的椭球区域的协方差矩阵的权重系数，ki的取值依据行列式最小准则确定；
5b)对各椭球相交区域的协方差矩阵C进行特征值分解，C＝VDV-1，得到特征值矩阵D和特征向量矩阵V；
其中，的对角线元素表示初始栅格单元三维尺寸，m表示松弛因子，m的取值由雷达参数和布站情形确定；
5c)根据特征向量矩阵V的正交特性，将C＝VDV-1变换为：C＝(V-1)TDV-1，得到初始栅格单元对应的旋转矩阵为：P＝V-1其中(·)Τ表示矩阵的转置运算；
6)根据初始栅格单元的大小 和对应的旋转矩阵P，将监视区域划分为E个格网单元；
根据指向监视区域的每个雷达与E个格网单元的角度关系，将每个格网单元对应到该雷达的接收通道；
根据指向监视区域的每个雷达与E个格网单元的距离关系，将每个划分出的格网单元对应到该雷达的距离单元。
2.根据权利要求1所述的方法，其中步骤2)中根据雷达组网系统内各雷达波束指向选出照射到监视区域的雷达，是对雷达组网系统内每个雷达分别作判别，根据监视区域中心点坐标和雷达位置坐标计算出监视区域中心相对于该雷达的方位角ɑ和俯仰角β，若ɑ、β分别在该雷达当前时刻发射波束的方位、俯仰角度范围内，则将该雷达归为指向监视区域的雷达。
3.根据权利要求1所述的方法，其中步骤3)中计算雷达组网系统内各雷达在监视区域中心点处形成的初始空间分辨单元信息以及初始栅格单元中心点相对于指向监视区域的各雷达的方位角θi和俯仰角 按如下步骤进行：
3a)为了简化计算，将雷达在监视区域中心点处形成的初始空间分辨单元近似为长方体；
3b)设初始空间分辨单元的长方体中心与雷达站的距离为r，计算得出该长方体在距离维上的长度Δr，在方位维上的长度rΔθ，在俯仰维上长度 其中，r通过距离计算公式得到，Δr、Δθ和 分别为已知的雷达距离分辨率、方位角分辨率和俯仰角分辨率；
3c)根据初始栅格单元中心点坐标和雷达位置坐标即可计算出初始栅格单元相对于指向监视区域的各雷达的方位角θi和俯仰角 其中1≤i≤N,N表示所述指向监视区域的雷达总数。
4.根据权利要求1所述的方法，其中步骤4)中计算椭球区域的协方差矩阵Ci，按如下步骤进行：
4a)构建关于俯仰角的旋转矩阵R1和关于方位角的旋转矩阵R2：其中 为初始栅格单元相对于指向监视区域的第i部雷达的俯仰角，θi为初始栅格单元相对于指向监视区域的第i部雷达的方位角，1≤i≤N,N表示所述指向监视区域的雷达总数；
4c)根据综合旋转矩阵R，得到初始空间分辨单元的协方差矩阵Ci为：其中(·)Τ表示矩阵的转置运算，ai表示第i部雷达在监视区域中心点处形成的初始空间分辨单元所代表的长方体在距离维上的长度，bi表示该长方体在方位维上的长度,ci表示该长方体在俯仰维上的长度,1≤i≤N,N表示所述指向监视区域的雷达总数。
5.根据权利要求1所述的方法，其中步骤6)中根据初始栅格单元的大小 和对应的旋转矩阵P，将监视区域划分为E个格网单元，按如下步骤划分：
6a)在笛卡尔坐标系中绘制一立方体，其中，立方体的中心位于坐标系原点，各边平行于坐标轴，且边长等于所述监视区域中心点到所述监视区域边缘的最大距离；
6b)设定基准长方体的中心在坐标系原点，各边平行于坐标轴，且各边在x,y,z轴上长度分别为 的对角线元素值；
以基准长方体为中心对6a)的立方体用多个长方体进行无缝填充，其中所有长方体的形状、尺寸均与基准长方体相同；
记各长方体中心点的坐标集合为O＝{ri,i＝1,2,3,...,M}，其中，ri代表第i个长方体的中心点坐标行向量，M表示长方体的个数；
6c)对各长方体中心点的坐标集合O中的点作如下变换：
记坐标平移旋转变化后格网单元中心点的坐标集合 其中，qi表示变换后的第i个元素，s表示所述初始栅格单元的中心点坐标行向量；
6d)将坐标平移旋转变化后格网单元中心点的坐标集合 中，不在所述监视区域中的多余点剔除掉，并记剔除多余点后格网单元中心点的坐标集合为 其中ti为最终划分出的各格网单元中心点的坐标行向量，E为最终划分出的格网单元个数。
雷达组网系统的空间配准方法
技术领域
[0001] 本发明属于雷达技术领域，具体的说是一种雷达组网系统的空间配准方法，可用于目标检测中的信号融合。
背景技术
[0002] 雷达组网系统的信息融合是雷达领域的热点问题。雷达组网系统的信息融合分为三个层次：航迹融合、判决融合和信号融合，信号融合的难度最大，但其能够提高雷达对目标的检测性能。雷达组网系统进行信号融合的前提是雷达组网系统内各雷达的信号能够在空间中配准，但由于各雷达观察警戒模式不同、位置不同、雷达参数不同以及所处环境不同等因素，导致对雷达组网系统内各雷达接收信号成功进行空间配准有较大的难度。实际上，空间配准问题几乎成了信号融合检测领域普遍回避的问题，大部分信号融合的研究都是在假设各雷达站信号已经完成空间配准的条件下进行的。
[0003] 已有的空间配准方法主要分为两类：离线处理和在线处理。离线处理方法以解决传感器系统偏差为主要目的。离线处理方法主要包括：最小二乘法、广义最小二乘法、最大似然法和精确极大似然估计等。在线处理方法主要用于实时估计系统偏差，能较好的表现复杂环境、噪声影响及特殊运动状态下的误差动态变化，具有更好的灵活性和适用性。在线处理方法主要包括：基于卡尔曼滤波算法、改进的卡尔曼滤波算法以及序贯卡尔曼滤波算法。
[0004] 上述已有的空间配准方法中，离线处理方法仅适用于固定误差因素，对非固定误差因素无法给出有效评测；在线处理方法尽管灵活性强，但操作不够简便，不利于工程应用。
发明内容
[0005] 本发明的目的在于提出一种雷达组网系统的空间配准方法，以解决上述现有技术操作复杂和不能对非固定误差因素进行有效评测的不足，实现雷达组网系统内信号融合时的空间配准。
[0006] 本发明对于一个包含多部雷达的雷达组网系统，其进行空间配准的技术方案包括如下：
[0007] 1)设定雷达组网系统的监视区域，选定监视区域中心点作为初始栅格单元的中心点；
[0008] 2)根据雷达组网系统内各雷达波束指向选出照射到监视区域的雷达；
[0009] 3)计算雷达组网系统内各雷达在监视区域中心点处形成的初始空间分辨单元信息以及初始栅格单元中心点相对于指向监视区域的各雷达的方位角θi和俯仰角 其中1≤i≤N,N表示所述指向监视区域的雷达总数；
[0010] 4)将每个雷达的初始空间分辨单元近似为椭球形状，计算该椭球区域的协方差矩阵Ci，其中1≤i≤N,N表示所述指向监视区域的雷达总数；
[0011] 5)计算初始栅格单元的大小 和对应的旋转矩阵P；
[0012] 5a)根据各椭球区域的协方差矩阵Ci，计算各椭球相交区域的协方差矩阵C：
[0013]
[0014] 其中，(·)-1表示矩阵求逆运算，ki表示所述指向监视区域的第i部雷达对应的椭球区域的协方差矩阵的权重系数，ki的取值依据行列式最小准则确定；
[0015] 5b)对各椭球相交区域的协方差矩阵C进行特征值分解，C＝VDV-1，得到特征值矩阵D和特征向量矩阵V；
[0016] 根据特征值矩阵D得到初始栅格单元的大小 为：
[0017]
[0018] 其中， 的对角线元素表示初始栅格单元三维尺寸，m表示松弛因子，m的取值由雷达参数和布站情形确定；
[0019] 5c)根据特征向量矩阵V的正交特性，将C＝VDV-1变换为：C＝(V-1)TDV-1，得到初始栅格单元对应的旋转矩阵为：P＝V-1；
[0020] 6)根据初始栅格单元的大小 和对应的旋转矩阵P，将监视区域划分为E个格网单元；
[0021] 7)完成空间配准：
[0022] 根据指向监视区域的每个雷达与E个格网单元的角度关系，将每个格网单元对应到该雷达的接收通道；
[0023] 根据指向监视区域的每个雷达与E个格网单元的距离关系，将每个划分出的格网单元对应到该雷达的距离单元；
[0024] 将雷达的回波信号放入对应的格网单元完成空间配准。
[0025] 本发明方法根据多源信息融合理论，对多部雷达在监视区域中心的初始分辨单元的相交区域拟合得到初始栅格单元，并利用初始栅格单元对监视区域划分成格网，能够完成组网雷达系统空间配准的要求，与传统航迹起始方法相比能够保证在监视区域中能量损失较小。
[0026] 为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
附图说明
[0027] 图1是本发明的实现流程示意图；
[0028] 图2是本发明仿真使用的雷达布站示意图；
[0029] 图3是在松弛因子m＝2/3时，用本发明仿真4部雷达在栅格划分后的监视区域能量损失示意图；
[0030] 图4是在松弛因子m＝1时，用本发明仿真4部雷达在栅格划分后的监视区域能量损失示意图；
[0031] 图5是在松弛因子m＝3/4时，用本发明仿真4部雷达在栅格划分后的监视区域能量损失示意图；
[0032] 图6是在松弛因子m＝1/2，时用本发明仿真4部雷达在栅格划的分后监视区域能量损失示意图。
具体实施方式
[0033] 下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述。
[0034] 参照图1，本发明的实现步骤如下：
[0035] 步骤1，设定监视区域，并选取初始栅格单元中心点。
[0036] 设定雷达组网系统的监视区域，并选定所述监视区域的中心点作为初始栅格单元的中心点，所述初始栅格单元用于对所述监视区域进行栅格化分，本实例中，取初始栅格单元中心点坐标为(50,50)km。
[0037] 步骤2，初始化雷达参数，判断出指向监视区域的雷达。
[0038] 初始化雷达组网系统内的雷达参数，并根据雷达参数中当前时刻雷达发射波束的方位角和俯仰角确定指向监视区域的雷达，其中，所述雷达参数至少包括雷达站的坐标，雷达的距离分辨率，雷达的方位角分辨率，雷达的俯仰角分辨率，雷达接收通道的方位角和俯仰角以及当前时刻雷达发射波束的方位角和俯仰角。
[0039] 步骤3，计算每部指向监视区域的雷达的初始空间分辨单元信息，及初始栅格单元中心点相对于指向监视区域的各雷达的方位角和俯仰角。
[0040] 计算初始空间分辨单元信息具体表现为计算初始空间分辨单元的大小，为了简化计算，将雷达在监视区域中心点处形成的初始空间分辨单元近似为长方体。
[0041] 3.1)设初始空间分辨单元的长方体中心与雷达站的距离为r，该长方体在距离维上的长度为a＝Δr，在方位维上的长度为b＝rΔθ，在俯仰维上长度为 其中，r通过距离计算公式得到，Δr、Δθ和 分别为已知的雷达距离分辨率、方位角分辨率和俯仰角分辨率，本实例中，采用二维雷达进行仿真，取距离分辨率Δr＝150m，方位角分辨率Δθ＝π/180；
[0042] 3.2)根据初始栅格单元中心点坐标和雷达位置坐标即可计算出初始栅格单元相对于指向监视区域的各雷达的方位角和俯仰角。
[0043] 步骤4，计算各雷达在监视单元中心形成的初始空间分辨单元相交区域的协方差矩阵。
[0044] 为了提取初始空间分辨单元的相交区域信息，本发明对初始空间分辨单元重新构建计算模型：以初始空间分辨单元的中心点为球心，对指向监视区域的每部雷达构建一个用以拟合初始空间分辨单元的椭球区域，即取初始空间分辨单元的内切椭球作为用于拟合的椭球区域：
[0045] 4.1)记第i部指向监视区域的雷达对应的椭球区域的协方差矩阵为Ci，其计算方法为：
[0046]
[0047] 其中(·)Τ表示矩阵的转置运算，R表示由步骤3.2)中的初始栅格单元相对于指向监视区域的第i部雷达的方位角和俯仰角构建的综合旋转矩阵，ai表示第i部雷达在监视区域中心点处形成的初始空间分辨单元所代表的长方体在距离维上的长度，bi表示该长方体在方位维上的长度,ci表示该长方体在俯仰维上的长度，1≤i≤N,N表示所述指向监视区域的雷达总数。
[0048] 4.2)根据行列式最小原则和步骤4.1)得到的椭球区域的协方差矩阵Ci，计算包含步骤4.1)所述椭球区域的相交区域的协方差矩阵C：
[0049]
[0050] 其中，(·)-1表示矩阵求逆运算，N表示所述指向监视区域的雷达总数，ki表示所述指向监视区域的第i部雷达对应的椭球区域的协方差矩阵的权重系数；ki的取值依据行列式最小准则确定。
[0051] 步骤5，计算初始栅格单元的大小 和对应的旋转矩阵P。
[0052] 步骤4中所述的初始空间分辨单元的相交区域是一个椭球形的区域，直接在监视区域无空隙生长比较困难；本发明以步骤4中得到的初始空间分辨单元的相交区域的信息构建初始栅格单元，将所述相交区域重构为长方体，即取该椭球区域的外切长方体，以简化栅格划分计算模型，具体方法如下：
[0053] 5.1)对各椭球相交区域的协方差矩阵C进行特征值分解，C＝VDV-1，得到特征值矩阵D和特征向量矩阵V；
[0054] 根据特征值矩阵D得到初始栅格单元的大小 为：
[0055]
[0056] 其中，的对角线元素表示初始栅格单元三维尺寸，m表示松弛因子，m的取值由雷达参数和布站情形确定，本实例中，取松弛因子m＝2/3；
[0057] 5.2)根据特征向量矩阵V的正交特性，将C＝VDV-1变换为：C＝(V-1)TDV-1，得到初始栅格单元对应的旋转矩阵为：P＝V-1。
[0058] 步骤6，对监视区域进行栅格划分。
[0059] 由步骤5得到的长方体是对相交椭球区域的近似，该长方体的边与该椭球的三个长半轴方向相同，在大多数情况下与笛卡尔坐标系的轴不平行，故本发明中所述的栅格划分即为从监视区域的中心开始对任意分布的长方体状初始栅格单元进行生长直至完成对监视区域的无缝填充；长方体的生长存在两种方案：1、直接在以长方体的中心为原点，以与长方体的边平行的三个方向为轴的坐标系上进行生长，2、先在在笛卡尔坐标系中进行生长，而后根据旋转矩阵对坐标系进行旋转；方案2在实际操作中可以简化计算量，在本发明中采取方案2的思想，具体步骤为：
[0060] 6.1)在笛卡尔坐标系中绘制一立方体，其中，立方体的中心位于坐标系原点，立方体各边平行于坐标轴，且边长等于所述监视区域中心点到所述监视区域边缘的最大距离；
[0061] 6.2)设定基准长方体的中心在坐标系原点，各边平行于坐标轴，且各边在x,y,z轴上长度分别为 的对角线元素值；以基准长方体为中心对6.1)的立方体用多个长方体进行无缝填充，其中所有长方体的形状、尺寸均与基准长方体相同；记各长方体中心点的坐标集合为：
[0062] O＝{ri,i＝1,2,3,...,M}，
[0063] 其中，ri代表第i个长方体的中心点坐标行向量，M表示长方体的个数；
[0064] 6.3)对各长方体中心点的坐标集合O中的点作如下变换：
[0065] qi＝riP+s；
[0066] 记坐标平移旋转变化后格网单元中心点的坐标集合为：
[0067]
[0068] 其中，qi表示上述变换后的第i个元素，s表示所述初始栅格单元的中心点坐标行向量；
[0069] 6.4)在坐标平移旋转变化后格网单元中心点的坐标集合 中，剔除掉不在所述监视区域中的多余点，并记剔除多余点后格网单元中心点的坐标集合为：
[0070]
[0071] 其中ti为最终划分出的各格网单元中心点的坐标行向量，E为最终划分出的格网单元个数。
[0072] 步骤7，空间配准。
[0073] 本发明所述的空间配准是指：针对划分出的格网，对各雷达站的回波信号在格网单元处进行采样，其实现途径为：
[0074] 根据指向监视区域的每个雷达与步骤6得到的E个格网单元的角度关系，将每个格网单元与该雷达的接收通道形成对应关系；
[0075] 根据指向监视区域的每个雷达及步骤6得到的E个格网单元的距离关系，将每个划分出的格网单元与该雷达的距离单元形成对应关系；
[0076] 将每个雷达的回波信号放入与该雷达的接收通道和距离单元相对应的格网单元内，完成空间配准。
[0077] 本发明的效果可通过以下仿真对比试验进一步说明：
[0078] 1.仿真场景：实验采用二维雷达，认为雷达1～雷达4分别坐落于(0,0)km，(50,0)km、(100,0)km和(0,50)km，波束指向分别为π/4、π/2、3π/4和0，设光速c＝3×108m/s，发射信号带宽均为B＝1MHz，方位角分辨率均为θBeam＝π/180，监视区域为点(49.7,50.3)km、(50.3,50.3)km、(49.7,49.7)km和(50.3,49.7)km组成的矩形区域，监视区域中心点为(50,50)km，松弛因子m＝2/3，各雷达在监视区域处的能量均为1，雷达布站示意图如图2所示。
[0079] 2.仿真内容：
[0080] 仿真1：采用以上实验场景，给出在松弛因子m＝2/3时，4部雷达在栅格划分后的监视区域能量损失图，结果如图3；
[0081] 仿真2：采用以上实验场景，给出在松弛因子m＝1时，4部雷达在栅格划分后的监视区域能量损失图，结果如图4；
[0082] 仿真3：采用以上实验场景，给出在松弛因子m＝3/4时，4部雷达在栅格划分后的监视区域能量损失图，结果如图5；
[0083] 仿真4：采用以上实验场景，给出在松弛因子m＝1/2时，4部雷达在栅格划分后的监视区域能量损失图，结果如图6；
[0084] 每幅图的纵横坐标单位均为米，用x(m)和y(m)表示的，每幅图中黑色*表示栅格单元中心，颜色的对应值表示能量值，监视区域能量的最大值等于照射到监视区域中的雷达站的能量值之和。
[0085] 3.实验结果分析：
[0086] 在图3中，栅格单元个数为19个，在栅格单元中心点处，能量损失为0，在各栅格单元交汇处，存在少量能量损失，且损失的能量为最大值的1/4＜1/2，符合工程要求，可见，该方法性能较好，可用于信号融合时的空间配准。
[0087] 比较图3～图6，栅格单元的个数分别为19个、9个、13个、和32个，可见在松弛因子不同时，监视区域内栅格点数量也不同。且各图中，在栅格单元中心点处，融合后能量都不存在损失，但在栅格单元的交汇处，随着松弛因子的减小，最大的能量损失分别为3/4、1/2、1/4和0，当能量损失超过1/2时，就会严重影响信号融合的性能。可见松弛因子影响着空间配准的性能，针对不同的雷达布站方式，松弛因子需重新确定。
[0088] 以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到变化或替换，都应涵盖在本发明的保护范围之内。
