在光学区，雷达目标的散射回波可视为有限个强散射点合成响应的结果，这些强散射点又可称为目标的散射中心^[1]。作为雷达目标的主要电磁散射特性之一，散射中心在雷达目标识别^[2]、逆合成孔径雷达(Inverse Synthetic Aperture Radar，ISAR)超分辨率成像^[3-4]、雷达散射截面(Radar Cross Section，RCS)外推^[5-6]等军事领域有着广泛的应用。因此，如何快速、有效提取散射中心成为****们致力于解决的一类重要课题。
随着隐身技术的高速发展，边缘绕射等散射类型逐渐替代镜面散射，成为隐身目标的主要散射源。基于几何绕射理论(GTD)的模型^[7-8]可有效地描述包括边缘绕射、镜面散射等多种散射类。
因此，GTD散射中心模型可准确刻画隐身目标的电磁散射特性，并成为目前描述隐身目标的主要散射中心模型之一。
准确估计提取GTD散射中心模型参数对描述隐身目标的电磁散射特性具有重要意义。目前，估计GTD散射中心模型参数的方法主要可分为两大类：一类为非参量法，如基于逆傅里叶变换(IFFT)的CLEAN算法^[9]、基于射线跟踪的CLEAN算法^[10]等；另一类为参量法，如MUSIC类算法^[11-12]、ESPRIT类算法^[13-14]、MEMP类算法^[15]等。由于非参量法分辨率较低，参数估计精度不高，因而大多数研究工作均是针对参量法开展的。
本文基于经典MUSIC算法，对一维GTD散射中心模型参数进行估计提取。针对经典MUSIC算法噪声鲁棒性较差、参数估计精度不高这一问题，本文提出了一类改进算法，有效提高了原算法的噪声鲁棒性与参数估计性能。此类改进算法通过构建一个置换矩阵J，可得到原始回波数据的共轭信息；进而将原始回波数据的自相关协方差矩阵R_EeEe，共轭数据的自相关协方差矩阵R_YY叠加取平均处理，得到新的总协方差矩阵R；最后对矩阵R作偶次方处理，得到最终的总协方差矩阵R₁。仿真实验验证了改进算法的有效性与先进性，并综合考虑参数估计精度与算法运算量，给出最合理的偶次方次数。
1 一维GTD散射中心模型及其简化 基于GTD，用来表征雷达目标散射回波的一维GTD散射中心模型的表达式为

(1)

式中: {A_i, α_i, r_i}分别为第i个散射中心的散射强度、散射类型、位置参数；I为散射中心的数目；f_m=f₀+mΔf，f₀为初始频率，Δf为步进频率，m为频率序号, 总频率步进数为M；c=3×10⁸ m/s为电磁波传播速度；w(f_m)为复高斯白噪声。α_i为0.5的整数倍^[7]，共有5种，不同散射类型对应不同的α值，如表 1所示。
表 1 典型散射结构的α值 Table 1 α values of typical scattering structures

典型散射结构	α值
二面角、三面角、平面法向反射	1
单曲面反射、圆柱面反射	0.5
双曲面反射、球面反射	0
边缘绕射	-0.5
尖顶绕射	-1

表选项






为了简化参数估计的复杂度，当雷达步进频率相对于初始频率较小时，即Δf/f₀ < < 1时，有如下近似：

(2)

将式(2)代入式(1)，可得

(3)

式中：

(4)

(5)

因此，求得式(5)中的P_i，即可求解得到散射中心的类型参数α_i、位置参数r_i，具体求解如下：

(6)

(7)

2 一维GTD散射中心模型参数估计 2.1 经典MUSIC算法 利用经典MUSIC算法估计一维GTD散射中心模型各类参数的步骤如下：
步骤1 ?为使GTD散射中心模型符合经典MUSIC算法估计参数时的模型，对式(3)作如下向量化重构：

(8)

式中：

(9)

(10)

(11)

(12)

(13)

步骤2 ?对式(8)进行数据重构，得到具有汉克尔形式的散射回波数据矩阵X₀：

(14)

式中：L=M－P+1，I < P < M/2。
步骤3 ?求X₀的自相关矩阵：

(15)

步骤4 ?对R_X0X0作特征值分解：

(16)

式中：U_S和U_N分别为信号子空间矩阵和噪声子空间矩阵；Σ_S和Σ_N两矩阵对角线上的元素分别为信号的特征值和噪声的特征值。
步骤5 ?根据信号子空间与噪声子空间的正交特性，定义一个参数估计算子P_w(α, r)：

(17)

对式(17)进行谱峰搜索，则峰值对应的信号参数r和α即为散射中心位置参数和类型参数。
步骤6?基于估计得到的r和α，利用最小二乘法^[16]可得到散射中心强度参数A_i，求解公式如下：

(18)

2.2 改进的MUSIC算法 利用改进MUSIC算法估计一维GTD散射中心模型各类参数的步骤如下：
步骤1 ?定义一个维度为P×L的置换矩阵J，其反对角线上的元素为1，其他位置上的元素为0，如下：

(19)

步骤2?基于含有目标原始回波信息的矩阵X₀及置换矩阵J，构造包含目标原始回波共轭信息的矩阵Y：

(20)

式中：X₀^*为X₀的共轭转置矩阵。
步骤3 ?分别计算矩阵束X₀及Y的自相关协方差矩阵：

(21)

步骤4 ?对求得的矩阵R_X0X0、R_YY相加并取平均值，得到一新的总协方差矩阵R：

(22)

步骤5?分析可知，新的总协方差矩阵R为P×P的Hermite矩阵，因此可得到如下关系式：

(23)

由式(23)可知，对(RR^H)ⁿ进行特征值分解，所得的特征值为原特征值的2n倍，且其对应的特征向量与原特征向量相同，即

(24)

式中：λ_(RR^H)ⁿ和υ_(RR^H)ⁿ分别为新构建矩阵(RR^H)ⁿ的特征值和特征向量；λ_R和υ_R分别为原矩阵R的特征值和特征向量。
因此，通过构建如式(25)中的矩阵R₁替代矩阵R，可增大信号特征值与噪声特征值之间的差距，从而可有效区分低信噪比情况下的信号与噪声，等效为增大了信噪比，改善了算法的噪声鲁棒性与参数估计性能。

(25)

从数学表达式来定量分析，各散射中心参数的方差为

(26)

式中：和z分别为第i次蒙特卡罗估计的参数和设定参数；σ²和γ_i分别为噪声特征值和信号特征值；(U_S, U_N)为信号子空间和噪声子空间的合并矩阵。
则由式(26)可知，通过增大σ²与γ_i之间的差距，会相应减小。因此，构建式(25)的矩阵R₁用以代替矩阵R，可增大等效信噪比，从而可以有效提高参数的估计精度。
步骤6 ?对得到的矩阵R₁作2.1节中的步骤3~步骤6处理，即可估计得到散射中心各类参数。
3 仿真实验 3.1 均方根误差比较 为验证本文提出的此类改进MUSIC算法的有效性，并与经典MUSIC算法的参数估计性能作比较，设置如下仿真实验。
本文设定雷达目标的后向散射回波由4个散射中心构成，具体散射参数值如表 2所示。步进雷达的初始频率f₀为10 GHz，步进频率Δf为16 MHz，总频率步进数M为100。
表 2 4个散射中心参数 Table 2 Parameters of four scattering centers

位置参数ri/m	散射类型αi	散射强度Ai
1.21	1.0	4.20
1.45	0.5	3.50
1.63	-0.5	2.20
1.84	-1	1.36

表选项






在目标电磁散射回波中加入信噪比SNR为-10~20 dB的一维复高斯白噪声，且每个信噪比下进行200次蒙特卡罗实验。信噪比SNR具体定义式^[17]为

(27)

式中：(f_m)和w(f_m)分别为信号对应的表达式及噪声对应的表达式，如下：

(28)

为定量比较算法的参数估计性能，定义均方根误差(Root Mean Square Error, RMSE):

(29)

式中：z_i、z和K分别为第i次仿真实验估计的散射中心参数、原散射中心参数和每个信噪比下的蒙特卡罗实验次数。
基于经典MUSIC算法与不同偶次方下的改进MUSIC算法进行一维GTD散射中心模型参数估计，各类参数均方差的仿真结果如图 1~图 3所示。

图 1 不同算法下，r1~r4的RMSE比较 Fig. 1 comparison of r1-r4 among different algorithms

图选项

图 2 不同算法下，α1~α4的RMSE比较 Fig. 2 RMSE comparison of α1-α4 among different algorithms

图选项

图 3 不同算法下，A1~A4的均方差比较 Fig. 3 RMSE comparison of A1-A4 among different algorithms

图选项

由图 1~图 3可知，当信噪比为-10~0 dB时基于不同次方处理的改进算法估计得到的各类散射中心参数的均方差均要低于经典MUSIC算法；而当信噪比为0~20 dB时，改进算法的参数均方差与经典算法相差不大或持平。换言之，在信噪比较低，即-10~0 dB时，此类改进算法的参数估计性能要优于经典MUSIC算法，改进算法的抗噪性能更优越。由图 1可知，随着偶次方次数的增加，改进算法的位置参数估计性能稍有提高；而比较图 2、图 3可得，基于不同次方下的改进算法估计得到的类型参数与强度参数的估计精度彼此相差不大。
分析原因可知，本文提出的改进算法首先利用原始回波数据共轭信息，有效提高了原始回波数据的利用率，其次通过将协方差矩阵R作偶次方处理，增大了噪声特征值与信号特征值之间的差距，避免了由于信号特征值与噪声特征值接近而出现的谱峰混叠问题，从而改善了算法的噪声鲁棒性，提高了算法的参数估计性能。随着偶次方次数的增大，信号特征值与噪声特征值差距越来越大；但较低的偶次方即可分辨出信号特征值与噪声特征值，因此随着偶次方次数的增大，算法估计性能改进得并不明显。
3.2 RCS重构精度比较 基于估计得到的模型各参数与一维GTD散射中心模型，可反演得到目标的后向散射电场E^s，进而由远场条件下目标电场与RCS之间的关系式，如式(30)所示，对目标的远场RCS进行拟合重构。

(30)

式中：E^s和Eⁱ分别为目标的散射电场和电磁波的入射电场；R_f为目标距探测雷达的远场距离。
在信噪比为20 dB的噪声条件下，利用经典MUSIC算法与二次方下的改进MUSIC算法，分别估计一维GTD散射中心模型参数并进行RCS拟合重构，与理论RCS作比较，分别得到目标的重构RCS曲线及RCS差值曲线如图 4、图 5所示。

图 4 重构RCS比较 Fig. 4 Comparison of reconstructed RCS

图选项

图 5 RCS差值比较 Fig. 5 Comparison of errors of RCS

图选项

由图 4可得，基于改进MUSIC算法重构得到的RCS与理论RCS吻合程度要高于经典MUSIC算法，因此可以更加准确地反映目标电磁散射特性。
从图 5可更为直观地看出，改进MUSIC算法重构的RCS与理论RCS的差值大多均在0 dBsm上下浮动，而经典MUSIC算法重构的RCS与理论RCS的差值多数要远大于0 dBsm。
由此可以验证，基于改进MUSIC算法估计得到的模型参数精度要高于经典MUSIC算法，即本文改进MUSIC算法的参数估计性能更为优越，噪声鲁棒性更好。
3.3 主要改进步骤运算量与运算时间比较 由于此类改进MUSIC算法是针对总协方差矩阵构造的步骤进行改进，因此，本文仅比较2次方、4次方、6次方对应的改进算法与经典MUSIC算法在总协方差矩阵构造这一步的计算量。并在MATLAB 2017A、计算机内存为8 GB的仿真环境下，比较不同算法进行200次蒙特卡罗实验估计散射中心3类参数所需的运行时间，如表 3所示。
表 3 不同算法运算量与运算时间比较 Table 3 Comparison of computational complexity and running time among different algorithms

算法	计算协方差矩阵	200次运行时间/s
MUSIC	P2L	6.177181
改进MUSIC(2次方)	3P2L+PL	6.919602
改进MUSIC(4次方)	6P2L+PL	7.269313
改进MUSIC(6次方)	8P2L+PL	7.617376

表选项






分析表 3可知，由于改进算法MUSIC需要构建置换矩阵J及新的协方差矩阵R，并增加了对协方差矩阵作偶次方处理这一步骤，因而运算量与运算时间相比于经典MUSIC算法要稍有增加，但可以看出，算法运算量增加不大。
4 结论 本文在经典MUSIC算法的基础上，提出了一类改进的MUSIC算法，主要创新点如下：
1) 通过对经典MUSIC算法提取参数过程中的协方差矩阵作偶次方处理，得到一类改进的MUSIC算法，仿真实验验证了改进的MUSIC算法可以有效地提高散射中心参数估计精度。
2) 通过设置仿真实验，探究了不同偶次方次数对参数估计精度的影响，并权衡考虑运算复杂度和参数估计精度，得到最佳的偶次方次数。
3) 根据散射中心模型及目标远场条件下电场值与RCS之间的关系，对目标的RCS进行了频率外推，仿真结果验证了外推结果的可靠性。

参考文献

[1]	聂在平, 方大纲. 目标与环境电磁散射特性建模[M]. 北京: 国防工业出版社, 2009: 278-286. NIE Z P, FANG D G. Modeling of target and environment electromagnetic scattering characteristics[M]. Beijing: National Defence Industry Press, 2009: 278-286. (in Chinese)
[2]	DING B Y, WEN G J, HUANG X D, et al. Target recognition in synthetic aperture radar images via matching of attributed scattering centers[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2017, 10(7): 3334-3347. DOI:10.1109/JSTARS.2017.2671919
[3]	ZHU W, JIANG M, HE X. Comprehensive sensing based multilevel fast multipole acceleration for fast scattering center extraction and ISAR imaging[J]. Sensors, 2018, 18(7): 2024. DOI:10.3390/s18072024
[4]	CUI S, LI S, YAN H, et al.A method of 3D scattering center extraction based on ISAR images[C]//2016 IEEE International Conference on Electronic Information and Communication Technology.Piscataway: IEEE Press, 2016: 439-442.
[5]	QUAN X Y, ZHANG B C, WANG Z D, et al. An efficient data compression technique based on BPDN for scattered fields from complex targets[J]. Science China Information Sciences, 2017, 60(10): 109302. DOI:10.1007/s11432-016-0529-6
[6]	邱志强.基于空间谱估计的雷达目标散射中心提取研究[D].成都: 电子科技大学, 2016: 116-130. QIU Z Q.Research on radar target scattering center extraction on the spatial estimation algorithm[D].Chengdu: University of Electronic Science and Technology, 2016: 116-130(in Chinese).
[7]	POTTER L C, CHANG D M, CARRIER R, et al. A GTD-based parametric model for radar scattering[J]. IEEE Transactions on Antennas and Propagation, 1995, 43(10): 1058-1066. DOI:10.1109/8.467641
[8]	CUI S M, JIN D G, FANG D G. Scattering center analysis of radar targets with the use of the matrix pencil combined with the GTD-based method[J]. Microwave Optical Technology Letters, 1997, 14(4): 244-247. DOI:10.1002/(SICI)1098-2760(199703)14:4<244::AID-MOP15>3.0.CO;2-8
[9]	TSAO J, STEIBERG B D. Reduction of sidelobe and speckleartifacts in microwave imaging:The CLEAN technique[J]. IEEE Transactions on Antennas and Propagation, 1988, 36(4): 543-556.
[10]	RAJAN B, LING H. Three dimensional scattering center extraction using the shooting and bouncing ray technique[J]. IEEE Transactions on Antennas and Propagation, 1996, 44(11): 1445-1453. DOI:10.1109/8.542068
[11]	郑舒予, 张小宽, 宗彬锋. 基于改进MUSIC算法的散射中心参数提取及RCS重构[J]. 系统工程与电子技术, 2020, 42(1): 76-82. ZHENG S Y, ZHANG X K, ZONG B F. Extraction of scattering center parameters and reconstruction of RCS based on improved MUSIC algorithm[J]. Systems Engineering and Electronics, 2020, 42(1): 76-82. (in Chinese)
[12]	ZHANG H M, ZHANG H Y. Research on DOA estimation method of sonar radar target based on MUSIC algorithm[J]. Journal of Physics:Conference Series, 2019, 1176(3): 1-5.
[13]	LEYMEN S, KONSTANTIN U, PIERRE C. Multidimensional ESPRIT for damped and undamped signals:Algorithm, computations, and perturbation analysis[J]. IEEE Transactions on Signal Processing, 2017, 65(22): 5897-5910. DOI:10.1109/TSP.2017.2736512
[14]	郑舒予, 张小宽, 宗彬锋, 等. 基于PQ-FB-2D-ESPRIT算法的二维GTD模型参数估计[J]. 系统工程与电子技术, 2020, 42(6): 1283-1289. ZHENG S Y, ZHANG X K, ZONG B F, et al. Parameter estimation of 2D-GTD model based on the PQ-FB-2D-ESPRIT algorithm[J]. Systems Engineering and Electronics, 2020, 42(6): 1283-1289. (in Chinese)
[15]	HUA Y B. Estimating two-dimensional frequencies by matrix enhancement and matrix pencil[J]. IEEE Transactions on Signal Processing, 1992, 40(9): 2267-2280. DOI:10.1109/78.157226
[16]	张贤达. 现代信号处理[M]. 北京: 清华大学出版社, 2002: 56-60. ZHANG X D. Modern signal processing[M]. Beijing: Tsinghua University Press, 2002: 56-60. (in Chinese)
[17]	王菁, 周建江. 一种基于GTD模型的目标散射中心提取方法[J]. 系统工程与电子技术, 2008, 30(11): 2147-2150. WANG J, ZHOU J J. Method of scattering center extraction based on GTD model[J]. Systems Engineering and Electronics, 2008, 30(11): 2147-2150. (in Chinese)