1 引言

干涉合成孔径雷达(InSAR)利用两幅SAR图像中的相位信息可以获取大范围、高精度的地表3维信息和变化信息，在军事侦察和民用测绘等方面有广阔的应用前景^[1]。但是InSAR只能获取各个分辨单元中单个散射点的高程信息，无法分辨出同一分辨单元内的不同散射点，即无法获取叠掩处的多个散射点的高程及散射系数等信息^[2]。同时随着高分辨率SAR的出现，叠掩问题显得尤为突出，对高精度城市遥感、民用测绘及军事侦察产生了很大影响，也成为了近些年国内外专家研究的热点之一^[3]。

在已有的解决叠掩问题的方法中，TomoSAR技术是解决叠掩问题的有效手段，其利用多基线技术，在高度维形成合成孔径从而实现SAR图像高度向的分辨能力，可有效获取SAR图像中同一分辨单元内多个散射点的高程及散射特性等信息。已有的TomoSAR方法主要有谱估计方法^[4–9]、最大似然估计方法^[10,11]、CS^[12–17]算法等。在传统的谱估计方法中，Capon, MUSIC, ESPRIT等方法具有高度向超分辨率的能力，而实际的机载多基线InSAR系统中基线大多为非均匀基线，使得谱估计的高度-谱函数中谱峰周围形成较大的旁瓣，从而影响叠掩处散射点的高程估计精度及分辨性能^[18]，同时机载多基线InSAR系统基线数目较少(3～4个左右)，使得估计叠掩处散射点的数目较少，高度估计精度也大幅降低。此外，传统的谱估计算法均基于信号的2阶统计特性，即阵列接收数据的协方差矩阵，但实际中机载SAR回波信号并不能保证被2阶统计特性完全描述，此时谱估计方法会产生较大误差。

基于以上问题，本文将高阶统计量中的4阶累积量用于信号的统计特性分析中，扩展信号的统计特性维度，利用4阶累积量对高斯噪声的不敏感性，提高算法对噪声的抗干扰性和散射点高度估计精度，同时4阶累积量在非均匀阵列条件下有良好的虚拟阵列扩展性能，从而提高散射点高度向估计精度及散射点的分辨能力^[19,20]。该算法将4阶累积量与已有的谱估计方法(Capon, MUSIC)相结合来解决叠掩问题，通过仿真实验表明，基于4阶累积量的谱估计方法在叠掩处散射点高度估计精度及分辨能力方面较传统的谱估计方法均有所提高，也通过实测数据验证了所提方法的有效性。

本文结构如下，第2节给出相应的信号模型及基于4阶累积量的相关介绍；第3节介绍基于4阶累积量的谱估计方法，引出基于4阶累积量的Capon方法和MUSIC方法；第4节利用4组仿真实验和实测数据验证本文方法与传统的谱估计方法相比的优越性；最后，第5节对本文方法进行总结。

2 信号模型 2.1 模型

机载非均匀多基线SAR系统叠掩成像原理模型如图1所示。在一发多收模式下，机载SAR各天线所获取的SAR图像经过图像配准和去斜处理后，第 $k$ 个天线在高程向的成像模型可表示为^[21]：

其中， $\Delta s$ 表示高程向的范围， ${\xi _k}$ 表示沿着高程向的空间频率， $\gamma (s)$ 表示各散射点散射系数在高程向的分布。考虑到噪声，高程向的离散模型简化为：

其中， ${g}$ 表示天线数为 $K$ 的观测向量， ${A}$ 中的矩阵形式为 ${{A}_{kl}} = \exp ( - {\rm j}2{{π}}{\xi _k}{s_l})$ ，为 $K \times L$ 的矩阵， $L$ 为散射点数目， ${\xi _k} = -{{{B_k}\cos ({\theta _l} - \alpha )}}/\left({{\lambda R}}\right)$ , ${\theta _l}$ 为散射点 $l$ 对应的下视角， $R$ 为参考斜距， ${B_k}$ 为第 $k$ 个天线对应的基线长度， $\alpha $ 为基线倾角， ${λ} $ 为信号波长。高程向 $s$ 和高度向 $h$ 的关系可表示为 $s = {h}/{{\sin ({\theta _l})}}$ ，此时 ${{A}_{kl}} = \exp \left( {{\rm j}\displaystyle\frac{{2{{π}}}}{\lambda } \cdot \frac{{{B_k}\cos ({\theta _l} - \alpha )}}{{R\sin ({\theta _l})}} \cdot {h_l}} \right)$ ，同时令 ${A} = \left[ {{a}({\varphi _1}),{a}({\varphi _2}), ·\!·\!· ,{a}({\varphi _L})} \right]$ ，其中 ${a}({\varphi _l}) = \left[ {{A}_{1l}},\right.$ $ \left. {{A}_{2l}}, ·\!·\!· ,{{A}_{Kl}} \right]^{\rm T}$ 称为阵列方向矢量， ${\rm T}$ 表示转置操作。

2.2 4阶累积量

与2阶累积量相比，4阶累积量具有更好的信号表征特性，其作为高阶统计特性，广泛应用于阵列信号处理领域，其优点在于：高斯色噪声的高阶累积量为0，具有自动抑制加性高斯白噪声和高斯色噪声的能力；同时能够增加虚拟阵元，扩展阵列孔径，从而使得阵列处理方法的分辨率和测向性能得到提高^[19,22]。

由式(2)中观测向量g得到的4阶累积量表达式为：

其中， $ \otimes $ 表示Kronecker积， $ * $ 表示共轭操作， ${\rm H}$ 表示共轭转置操作。

由式(3)可知，4阶累积量的计算需要多视来进行处理。如果各信源间相互独立，则有^[20]

其中，

${{C}\!_s}$ 为信号源的4阶累积量矩，可表示为：

${{C}\!_{ε}}$ 为噪声的4阶累积量矩，可表示为：

当噪声为高斯白噪声或高斯色噪声时，有 ${{C}\!_{ε}} \!=\! 0$ ，此时

从式(4)，式(5)可以看出4阶累积量矩阵经过阵列扩展后的阵列导引矢量为：

而相应的2阶统计量的阵列方向矢量为 ${a}(\varphi )$ ，其2阶统计特性为：

其中， ${P}\,$ 表示信号的协方差矩阵， $\sigma \!_R^{\ 2}$ 表示加性噪声功率。

比较式(3)和式(10)可以发现，基于4阶累积量的空间维度进行了扩展，由原来的2阶统计量的 $K \!\times\! K$ 维变为 ${K^2} \times {K^2}$ 维。同时由式(9)定义的导向矢量可以看出，该矢量除了包含 ${a}(\varphi )$ 矢量外，也增加了许多虚拟阵元产生的矢量，与均匀阵列相比，非均匀阵列模式下产生的虚拟阵元更多，虚拟基线长度也更长，而阵元数的增加和基线长度的扩展能够增强对叠掩处散射点高度处的识别能力，从而提高散射点的高度测量精度及散射点的分辨能力^[23]。

3 基于4阶累积量的谱估计方法 3.1 基于4阶累积量的Capon

由2.2节的分析可得，基于4阶累积量的Capon方法可表示为：

式(11)中的 ${\hat P_C}$ 表示随高度变化的功率谱函数，因此基于4阶累积量的Capon方法可通过搜索谱函数 ${\hat P_C}$ 中谱峰对应的高度获取叠掩处散射点的高度信息。

3.2 基于4阶累积量的MUSIC

为保证方法的稳定，提高信源的检测能力，有效保护4阶累积量的信号子空间的信息，现将 ${{C}\!_{g}}$ 进行重构^[24]，即

其中， ${{I}\!_{v}} = {\left[ {\begin{array}{*{20}{c}} 0&0& ·\!·\!· &1 \\ 0& ·\!·\!· &1&0 \\ \vdots & \vdots & \ddots & \vdots \\ 1& ·\!·\!· &0&0 \end{array}} \right]_{{L^2} \times {L^2}}}$ 。

由式(4)可知，可以通过对4阶累积量矩阵进行特征值分解，来获取4阶累积量的信号子空间和噪声子空间，此时有

其中， ${{Σ} _s}$ 表示大特征值组成的对角矩阵， ${{U}\!_s}$ 表示大特征值(4阶累积量的信号空间)对应的特征向量组成的矩阵， ${{Σ}_{ε} }$ 表示小特征值组成的对角矩阵， ${{U}\!_{ε}}$ 表示小特征值(4阶累积量的噪声空间)对应的特征向量组成的矩阵。

在获取叠掩处散射点的数目后，可确定相应的4阶累积量的噪声子空间，然后根据噪声子空间和导向矢量的正交性，基于4阶累积量的MUSIC可表示为：

式(15)中的 ${\hat P\!_M}$ 表示随高度变化的伪谱函数，因此基于4阶累积量的MUSIC方法可通过搜索谱函数 ${\hat P_M}$ 中谱峰对应的高度获取叠掩处散射点的高度信息。此方法需要预先确定叠掩处散射点的数目。

4 实验验证

为了验证基于4阶累积量的算法在叠掩处散射点高度估计的有效性和优越性，本文进行4组仿真实验，分别从信噪比、视数、散射点高度差及高度-归一化谱函数对比等4个方面进行仿真实验，与传统的谱估计解叠掩方法进行对比，并结合相关实测数据进行处理验证。本节首先介绍仿真的系统参数，然后进行4组仿真实验验证算法的有效性，最后利用实测数据进行分析。

4.1 仿真系统参数

仿真实验以德国MEMPHIS的三基线系统为例，详细参数如表1所示。由文献[21]可知，单航过机载多基线InSAR系统高度向的理论分辨率为 ${\rho _h} = $ $ {{\lambda R\sin (\theta )}}/{{\Delta B}}$ ，对应的高度向的克拉美罗界^[25]为 ${\sigma _h} = \displaystyle\frac{{\lambda R\sin (\theta )}}{{2{{π}}\sqrt K \cdot \sqrt {2 \cdot{\rm SNR}} \cdot {\sigma\!_B}}}$ ，高度向最大不模糊距离为 ${h_{\rm amb}} = \displaystyle\frac{{\lambda R\sin (\theta )}}{{\Delta {b_{\min }}}}$ 。其中， $\theta $ 为雷达下视角， $R$ 为斜距， $ {λ} $ 代表波长， $\Delta B$ 为基线总长度， $\Delta {b_{\min }}$ 为最短基线长度， ${\sigma _B}$ 为基线的均方差， $K$ 代表天线数目。由表1可知，仿真系统的高度向理论分辨率约为42 m，高度向理论最小均方差约为1.2 m，高度向最大不模糊距离为210 m。

4.2 仿真结果分析

由多基线SAR成像模型可知，信噪比、视数、高度差对叠掩处散射点高度估计产生重要影响，因此从这3方面分别进行仿真实验，并利用高度-归一化谱分布与传统的谱估计方法(Capon, MUSIC方法)进行直观对比仿真实验，验证算法的有效性。在仿真实验中，计算叠掩处散射点高度估计的均值、均方根误差(RMSE)，每次试验均通过1000 次Monte Carlo仿真来确保实验可靠性，最后通过高度-归一化谱函数分布图对本文方法进行直观说明。

仿真实验1 首先考察信噪比对基于4阶累积量的谱估计方法的高度估计性能。仿真实验1主要为了验证基于4阶累积量的谱估计方法随信噪比(SNR)的变化情况，并与传统的谱估计解叠掩方法进行比较。由于仿真系统参数的高度向理论分辨率为42 m，为了确保其超分辨性能，现选取叠掩处高度为 ${h_1} = 0 \ {\rm m}$ , ${h_2}{\rm{ = }}30 \ {\rm m}$ ，幅度一致的两个散射点来进行实验，并与传统的谱估计方法进行对比验证，仿真实验视数设为8视，SNR变化范围设置为 ${\rm SNR} = \left[ {0:2:30} \right]$ dB。由于 ${h_1}$ , ${h_2}$ 的结果类似，现将 ${h_2}{\rm{ = }}30 \ {\rm m}$ 的仿真结果(均值、RMSE)显示如图2，其中图2(a)中绿色实线表示理想高度 $h = 30 \ {\rm m}$ 。

由图2(a)可以看出，传统的Capon, MUSIC方法随SNR变化不大，且均值较真实值差别较大，而基于4阶累积量的Capon, MUSIC方法测量的高度向均值随SNR的变化有减小的趋势，与传统的Capon, MUSIC方法相比其均值更为趋近真实值。图2(b)基于4阶累积量的Capon, MUSIC方法中散射点高度测量均方差小于传统的谱估计方法，其随SNR变化不断下降，并在高信噪比(>20 dB)下趋于稳定。

从仿真实验1可以看出，基于4阶累积量的谱估计方法(Capon, MUSIC方法)在叠掩处散射点高度向测量均值、均方差等方面均优于传统的谱估计方法(Capon, MUSIC方法)，充分验证了4阶累积量对高斯噪声的抑制性能。

仿真实验2 仿真实验2主要为了验证基于4阶累积量的谱估计方法随视数的变化情况，并与传统的谱估计解叠掩方法进行比较。与仿真实验1类似，仍选取叠掩处高度为 ${h_1} = 0 \ {\rm m}$ , ${h_2}{\rm{ = }}30 \ {\rm m}$ ，幅度一致的两个散射点来进行实验。同时由仿真实验1可以看出两类方法在20 dB后表现趋于稳定，因此仿真实验的信噪比设置为SNR=20 dB。此外，视数变化范围设置为 $N = \left[ {2:2:32} \right]$ 视。由于 ${h_1}$ , ${h_2}$ 的结果类似，将 ${h_2}{\rm{ = }}30 \ {\rm m}$ 的仿真结果(均值、RMSE)显示如图3，其中图3(a)中绿色实线表示理想高度 $h = 30 \ {\rm m}$ 。

由图3(a)可以看出，4种方法得到的散射点高度均值随视数变化不大，但基于4阶累积量的Capon, MUSIC方法估计的高度均值小于传统的Capon, MUSIC方法，更为接近真实值。图3(b)中随着视数的增加，高度均方差(RMSE)也相应地变小，但基于4阶累积量的Capon, MUSIC方法得到的RMSE小于传统的Capon, MUSIC方法，且在较小视数(6视)下RMSE趋于稳定。

从仿真实验2可以看出，基于4阶累积量的谱估计方法(Capon, MUSIC方法)能在较小视数下(6视)高度测量趋于稳定，且其在叠掩处散射点高度估计均值、RMSE等方面均优于传统的Capon, MUSIC方法。

仿真实验3 接着从叠掩处散射点高度差方面来验证本文所提方法的有效性和稳定性，并与传统的谱估计解叠掩方法进行比较。由仿真实验1和仿真实验2可以看出，两类方法在20 dB, 8视后高度测量较为稳定，因此本次仿真实验的信噪比和视数分别设置为20 dB和8视。同时叠掩处两个散射点中的一个高度为 ${h_1} = 0 \ {\rm m}$ ，另一个散射点高度变化范围设置为 ${h_2} = \left[ {2:2:50} \right]{\rm m}$ ，两个散射点幅度一致。由于 ${h_1}$ 的高度固定，现将高度变化的第2个散射点 ${h_2}$ 的仿真结果(均值、RMSE)显示如图4，其中图4(a)中绿色斜线表示随高度差变化的散射点的理想高度范围 $h = [2,50]\;{\rm m}$ 。

从图4(a)可以看出4种方法均具有一定的超分辨能力，其中基于4阶累积量的Capon, MUSIC方法散射点高度向分辨能力分别优于传统的Capon, MUSIC方法。图4(b)中4种方法各自在其高度差范围内的散射点高度估计RMSE趋势为先增大后减小，未能成功地分离两个散射点，也能较好地验证图4(a) 4种方法的高度分辨能力，而随着高度差的增大能够成功区分两个散射点，且基于4阶累积量的Capon, MUSIC方法所得到的RMSE分别小于传统的Capon, MUSIC方法。

从仿真实验3可以看出，基于4阶累积量的谱估计方法(Capon, MUSIC方法)的高度向分辨率较传统的谱估计方法(Capon, MUSIC方法)有所提高，且在相同高度差下的散射点的高度估计均值误差、均方差更小，因而有更高的稳定性。

仿真实验4 本次仿真实验从高度-谱函数图方面直观说明本文算法与传统谱估计算法的差别，验证本文算法的有效性。从前面仿真实验可以看出系统在信噪比为20 dB、视数为8视、高度差为30 m参数下，有较强的稳定性，本次仿真采用以上3个参数进行仿真实验，并对谱函数进行归一化，结果如图5所示。图5中的绿色虚线代表 ${h_1} = 0 \ {\rm m}$ , ${h_2} = 30 \ {\rm m}$ 的理想高度。

从图5可以看出基于4阶累积量的谱估计方法比传统的谱估计方法(Capon, MUSIC方法)峰值更为尖锐，由于4阶累积量对高斯噪声的抑制作用，使其归一化谱函数变化更为迅速，即在峰值附近下降更为明显，提高了叠掩处散射点高度估计的稳定性(均值、均方差)，与仿真实验1—仿真实验3的相关结论相吻合，同时由于其非均匀阵列下的虚拟阵元拓展作用，扩展虚拟阵列天线数目、统计空间维度、基线长度，使得其高度-归一化谱函数分布图谱峰主瓣更窄，提高了主瓣的分辨率，也能在一定程度上提高叠掩处散射点高度向的分辨率和估计精度，验证了仿真实验3的结果。

4.3 实测数据处理

由于实际情况的限制，本次实验所用的实测数据为机载毫米波三天线的SAR数据，由文献[21]求得的高度向理论分辨率为 ${\rho _h} \approx 72.7\ {\rm{m}}$ 。原始图像幅度图如图6(a)所示，其中纵向为距离向，横向为方位向，距离方位像素值为1024 $* $ 1536。长基线干涉条纹图和相干系数图分别如图6(b)和图6(c)所示，从图6(b)中可以看出中间部分存在“相位反偏”现象，同时结合图6(a)和图6(c)可以看出相干系数较低的区域包含阴影区域和叠掩区域。依据信息论方法^[10]获取的图像中散射点个数分布如图7(a)所示，可以看出图7(a)能分别对阴影、叠掩、单散射点区域进行较好地区分，验证了图6的分析，同时地形的复杂性导致叠掩处有可能出现多于两个散射点的现象存在，而本文实测数据所用的机载三天线SAR数据在叠掩处至多检测出两个散射点，因此多于两个散射点的情况不在本文讨论范围之内。

本节对原始三天线数据进行SAR图像配准、去斜等操作后，截取数据图像中方位向第700 像素、距离像素为600～800的像素点，并分别用传统的Capon, MUSIC方法和基于4阶累积量的Capon方法、MUSIC方法对该剖面图进行处理，实测数据的处理视数为5 $ *$ 5，得到相应的高度-归一化谱函数，如图8所示，图8中右侧色柱代表dB。该剖面图对应的散射点个数分布如图7(b)所示，可以看出该剖面图后半段包含两个散射点，利用该散射点剖面图，求得图8中对应的散射点高度剖面图如图9所示。

从图8中高度-归一化谱函数可以看出，基于4阶累积量的Capon, MUSIC方法分别较传统的Capon, MUSIC方法其谱峰更为尖锐明显，旁瓣更低。图9中距离向前半段只有一个散射点的高度存在，且高度连续性较好，而随着距离的增大，两个散射点随之出现，且随着距离向的增大，叠掩处两个散射点的高度在逐渐减小，能较好地印证叠掩处散射点高度随斜距的变化情况，且叠掩处散射点的高度差小于高度向理论分辨率，说明4种方法具有较好的叠掩处散射点超分辨能力。从图9中可以看出基于4阶累积量的Capon, MUSIC方法分别较传统的Capon, MUSIC方法求得的叠掩处散射点的高度连续性更好，更为符合实际情况，说明其高度测量更为稳定。

下面对所提算法的散射点高度向分辨能力进行进一步分析。图7(a)中叠掩的像素点个数为103182，占总数据的6.56%。现统计数据中每个叠掩处的两个散射点的高度，并对其高度差进行统计分析，如图10所示。图10(a)和图10(b)中的黑色竖线代表理论高度分辨率(约72.7 m)，统计4种方法高度差小于理论分辨率的像素点与数据中所有叠掩处像素点的比例，分别为 ${P_{\rm Capon}} = 82.88\% $ , ${P_{{4阶}\rm Capon}} = $ 84.87%, ${P_{\rm MUSIC}} = 84.06\% $ , ${P_{{4阶}\rm MUSIC}} = 84.18\% $ ，可以看出4种方法求得的高度差大部分小于高度理论分辨率，说明4种方法均具有一定的超分辨能力，并且基于4阶累积量的Capon, MUSIC方法所占的比例分别略大于传统的Capon, MUSIC方法。同时从图10(a)和图10(b)可以看出基于4阶累积量的Capon, MUSIC方法与传统的Capon, MUSIC方法在较小的高度差(<40 m)下统计像素点数更多，说明本文所提方法的叠掩处散射点高度向分辨能力强于传统的Capon, MUSIC方法。在基于4阶累积量的MUSIC方法中，需要确定散射点数目、信号及噪声子空间，而4阶累积量的方法使得信号子空间和噪声子空间造成了一定程度的耦合，因此对噪声子空间的阶数选取造成了一定的困难。本文噪声子空间的阶数可以选取为天线数与信源数差值的平方和。

5 结束语

针对机载多基线InSAR系统中非均匀基线和基线数目少导致传统的谱估计解叠掩方法中散射点高度估计不足的问题，本文将4阶累积量的统计特性与谱估计方法相结合，在有效去除高斯噪声的同时，利用非均匀阵列虚拟阵元扩展特性提高叠掩处散射点的分辨性能及高度估计精度，同时利用实际机载多基线InSAR仿真系统参数进行仿真试验，仿真实验验证了本文方法的有效性，也对实际机载多基线InSAR系统解决叠掩问题提供了良好的依据，实测数据实验初步说明本文方法在单航过机载多基线SAR系统下具有较好的叠掩处散射点分离能力，散射点高度向测量也更为稳定。