孙延鹏^1,2, 张石¹, 屈乐乐², 白文静²
1. 东北大学 计算机科学与工程学院, 辽宁 沈阳 110169;
2. 沈阳航空航天大学 电子信息工程学院, 辽宁 沈阳 110136
收稿日期：2016-01-15
基金项目：国家自然科学基金资助项目(61302172, 61671310)；辽宁省自然科学基金资助项目(2014024002, 201602565)；航空科学基金资助项目(2016ZC54013)。
作者简介：孙延鹏(1973-)，男，山东聊城人，东北大学博士研究生，沈阳航空航天大学教授;
张石(1963-)，男，辽宁抚顺人，东北大学教授，博士生导师。

摘要：针对传统压缩感知SFGPR成像重建算法在强杂波测量环境中往往会失效的问题, 提出一种基于子空间投影杂波抑制技术的SFGPR压缩感知成像重建算法.该算法首先在每个天线测量位置通过压缩感知测量模型重建所有的频域原始均匀采样数据, 然后采用子空间投影杂波抑制技术滤除较强的地面回波, 最后结合稀疏重建算法对地下目标图像进行压缩感知重建.实验数据处理结果验证了所提方法的有效性和准确性.
关键词：频率步进探地雷达压缩感知子空间投影成像算法
Subspace Projection Based Compressive Sensing SFGPR Imaging Algorithm
SUN Yan-peng^1,2, ZHANG Shi¹, QU Le-le², BAI Wen-jing²
1. School of Computer Science & Engineering, Northeastern University, Shenyang 110169, China;
2. School of Electronic and Information Engineering, Shenyang Aerospace University, Shenyang 110136, China
Corresponding author: SUN Yan-peng, E-mail: 781074263@qq.com
Abstract: The traditional compressive sensing (CS) stepped frequency ground penetrating radar (SFGPR) imaging algorithm usually loses effect in strong clutter environment. To alleviate this problem, a CS SFGPR imaging algorithm based on the subspace projection clutter suppression technique was proposed. The original uniform frequency sampling data at each measurement position were reconstructed from the reduced set of randomly measured data using CS measurement model. Then the subspace projection clutter suppression technique was employed to suppress the strong ground reflection signal. Finally the sparse reconstruction algorithm was used to reconstruct the image of underground targets. The experimental data has verified the validity and effectiveness of the proposed imaging method.
Key Words: stepped frequency ground penetrating radar (SFGPR)compressive sensingsubspace projectionimaging algorithm
由于频率步进探地雷达(SFGPR, stepped frequency ground penetrating radar)具有较高的动态范围和较强的射频抗干扰能力等优点, 近年来在探地雷达领域受到越来越多的关注^[1-3].目前, 探地雷达不断朝着多通道、多极化、多波段和高分辨率方向发展, 使得SFGPR存在数据采集量大和数据采集时间长等缺点.近年来, 美国学者Donoho等提出了压缩感知(CS, compressive sensing)理论^[4].目前, 国内外学者将CS理论应用到SFGPR的系统构建和成像重建等领域^[5-7].
但在SFGPR的实际测量过程中, 由于强地面杂波的存在, CS成像算法的重建性能会急剧下降^[8].所以需采用合适的杂波抑制方法对原始频域采样数据进行预处理以去除较强的地面反射波.传统的SFGPR杂波抑制算法大都是基于均匀采样数据^[9], 而在CS成像过程中, 考虑到成像数据的大幅度减少且在空频域表现为非均匀采样数据, 经典杂波抑制算法往往都不适用.针对上述问题, 本文提出了一种基于子空间投影杂波抑制的压缩感知SFGPR成像重建算法, 该算法首先在每个天线测量位置基于接收信号的稀疏性重建所有频点的信号, 然后通过经典的子空间投影杂波抑制算法去除地面强反射波, 最后基于CS成像模型对地下目标进行高分辨成像重建, 从而实现在强杂波环境对地下目标进行准确成像定位.
1 SFGPR信号模型与压缩采样对于传统SFGPR系统而言, 设系统的初始频率为f₀, 频率步长为Δf, 带宽B=(N-1)Δf.在发射带宽B内, 第n个频点的发射频率为

(1)

设测量过程中共有M个天线测量位置, 如果忽略地下媒质的衰减效应, 则系统在第m(m=0, 1, …, M-1) 个天线位置第n个工作频点的接收信号可表示为

(2)

其中:σ_s是对应空气与地面分界面的反射率; τ_s为系统到地面的双程传输延时; P为地下点目标的个数; σ_p(σ_p?σ_s)是第p个点目标对应的反射率; τ_{p, m}为在第m个测量位置系统到第p个点目标的双程传输延时.
为得到地下待成像区域的2维反射率分布r(k, l), 其中k=0, 1, …, N_x-1；l=0, 1, …, N_y-1, N_x和N_y分别表示水平和深度方向划分的网格数, 首先可将r(k, l)通过列堆叠转换成N_xN_y×1维反射率向量r, 然后根据r的估计结果通过重排即可得到r(k, l).这样, 式(2) 可表示成矩阵形式:

(3)

其中: 为N×1维频域测量数据向量; Ψ_m为N×N_xN_y维字典矩阵, 其第j列可表示为

(4)

根据M个测量孔径可得到MN×1维复合测量数据向量和MN×N_xN_y维复合字典矩阵.复合测量数据向量e和反射率向量r之间的关系可表示为

(5)

由于P往往远远小于网格数N_xN_y, 所以r具有较强的稀疏性.根据CS理论, 可以用Q₁Q₂个测量数据重建稀疏向量r.CS降采样过程可视为一个测量矩阵Φ作用到原始测量数据向量e上.则压缩采样数据向量可表示为

(6)

其中Φ为Q₁Q₂×MN维测量矩阵, 考虑到在空频域实现联合压缩采样, 则Φ可表示为

(7)

其中:?表示Kronecker积；I_Q1为Q₁×Q₁维单位矩阵；Φ_s表示Q₂×M维空域测量矩阵, 其构造可从M×M维单位矩阵中随机选取Q₂行得到.Φ_m(m=0, 1, …, M-1) 表示Q₁×N维频域测量矩阵, 其构造可从N×N维单位矩阵中随机选取Q₁行得到.
图 1给出了传统SFGPR和压缩感知SFGPR的空频域采样示意图, 可以看出, 与传统SFGPR相比, 压缩感知SFGPR的空频域采样数据将大幅减少.对于给定的压缩采样测量数据向量, 可通过求解l₁范数最优化问题即可重建反射率向量r：
图 1(Fig. 1)

图 1 SFGPR空频域采样示意图Fig.1 Views of SFGPR space frequency sampling pattern (a)—传统SFGPR空频域采样示意图；(b)—压缩感知SFGPR空频域采样示意图.

(8)

式(8) 表示的稀疏最优化问题可利用正交匹配追踪算法得到求解^[10].
2 压缩采样子空间投影杂波抑制2.1 单个天线测量位置的频域数据重建在CS测量模型下, 由于测量数据大幅度减少且表现为空频域随机采样数据, 使得传统杂波抑制算法不再适用.为了能够使用经典杂波抑制算法, 本文提出在每个天线测量位置基于CS测量数据重建原始发射频点的均匀采样数据.首先将双程探测时窗τ_max均匀划分为L个时延网格, 则可以得到L×1维双程传输时延向量τ=[τ₁, τ₂, …, τ_L]^T.则在第m个测量位置, 压缩采样频域数据向量可表示为

(9)

其中A是N×L维字典矩阵, 其第l列可表示为

(10)

b_m是L×1维回波幅度向量, 其第l个元素对应第l个时延网格的幅度, b_m表示为

(11)

在SFGPR实际探测过程中, 由于目标回波个数通常远小于划分的时延网格数, 所以b_m是稀疏向量. b_m可通过求解下面最优化问题得到：

(12)

式(12) 可通过正交匹配追踪算法求解得到b_m, 从而得到：

(13)

当所有的天线位置处频域均匀采样数据得到重建后, 即可使用合适的杂波抑制方法去除地面强反射波.
2.2 子空间投影杂波抑制当在每个测量位置恢复N个频点测量信号后, 假设共有Q₂个天线测量位置, 那么对应Q₂个天线位置N个发射频点的所有测量数据可构成N×Q₂维矩阵E.由于测量数据可分解为地面回波子空间和地下目标回波子空间^[11], 所以可对矩阵E进行奇异值分解:

(14)

其中H代表共轭转置, 正交矩阵U和V分别包含左右奇异向量, Λ为对角矩阵, 可表示为

(15)

其中：λ₁, λ₂, …, λ_Q2为奇异值, 奇异值按λ₁≥λ₂≥…≥λ_Q2降序排列.
由于地面反射波幅度比目标回波幅度强很多, 通常可以用E的前K个较大的奇异值对应的奇异向量来构造地面回波子空间A_ground.

(16)

其中u_i和v_i分别代表前K个较大奇异值对应的左、右奇异向量.
地面回波子空间的正交子空间A_ground^⊥为

(17)

其中I是单位矩阵.将E投影到A_ground^⊥上即可得到杂波抑制后的空频域采样数据：

(18)

至此就能得到地面杂波抑制后的目标回波空频域测量数据.然后再采用CS成像算法对地下目标进行成像重建.
3 实测数据处理利用本文所提算法对佐治亚理工大学公布的SFGPR系统实测数据^[1]进行处理, 系统工作频段为60 MHz~8.06 GHz, 频率步进间隔为20 MHz, 在每个天线测量位置共测得401个频点采样数据.实验在沙箱里进行, 收发天线相位中心距离地面的高度为27.8 cm, 收发天线间隔为12 cm, 移动步长为2 cm, 两个直径为12.7 cm的金属球埋入沙坑中, 沙的介电常数为4.两个金属球在地下的真实位置如图 2所示.选取天线位置从-60 cm到60 cm之间的61道数据, 对应61×401=24 461个空频域均匀采样数据.
图 2(Fig. 2)

图 2 金属球在地下的真实位置示意图Fig.2 True location of the buried metal balls

把待成像区域划分成61×31个空间网格, 每个网格尺寸为1 cm×1 cm.对应CS算法, 随机抽取30个天线位置, 每个天线位置抽取80个频点采样数据, 共有2 400个空频域随机采样数据, 只占传统SFGPR系统数据采集量的10.1%.图 3给出了未采用杂波抑制技术预处理的后向投影(back projection, BP)成像重建结果和CS成像重建结果.为了便于比较, 成像结果均采用归一化dB值显示, 显示动态范围为40 dB.可以看出在有地面强反射波情况下, CS成像算法已无法重建准确的目标像.
图 3(Fig. 3)

图 3 未经杂波抑制技术预处理的成像结果Fig.3 Imaging result without clutter suppression technology (a)—BP成像结果；(b)—CS成像结果.

图 4给出了经子空间投影杂波抑制预处理后的BP成像结果和CS成像结果.BP成像算法中直接采用子空间投影杂波抑制算法去除地面强反射波, 而CS成像算法采用本文方法去除地面强反射波.采用的计算环境为Intel(R) Core(TM) i7-4510U CPU @ 2.0 GHz, 4 GB内存, BP成像结果的重建时间为203.8 s, CS成像结果的重建时间为23.8 s.从图 4可以看出, 尽管两种成像算法都能得到准确的目标像, 但这时CS成像算法不仅所需数据量少, 而且成像结果的目标旁瓣较小, 分辨率高, 更利于对地下目标的辨识.
图 4(Fig. 4)

图 4 经子空间投影杂波抑制预处理的成像结果Fig.4 Imaging result with clutter suppression by the subspace projection (a)—BP成像结果；(b)—CS成像结果.

4 结论由于CS成像重建算法在强杂波测量环境中会失效, 而在CS测量模式下测量数据表现为空频域随机采样数据, 使得传统的杂波抑制算法不再适用.本文提出一种基于子空间投影杂波抑制的压缩感知SFGPR成像重建算法, 该算法首先在每个天线测量位置重建所有频域均匀采样数据, 然后基于正交子空间投影技术得到目标回波子空间, 最后利用CS成像模型重建高质量的地下目标像.实测数据处理结果表明所提成像重建算法能够实现在强杂波环境对地下目标进行准确成像定位, 而且成像所需数据量少, 成像结果具有旁瓣低和分辨率高等优点.
参考文献

[1]	Counts T, Gurbuz A C, Scott W R. Multistatic ground-penetrating radar experiments[J].IEEE Transactions on Geoscience and Remote Sensing, 2007, 45(8): 2544–2553.DOI:10.1109/TGRS.2007.900677
[2]	屈乐乐, 方广有, 杨天虹. 频率步进探地雷达系统小型化的设计与实现[J].仪器仪表学报, 2010, 31(3): 699–703. ( Qu Le-le, Fang Guang-you, Yang Tian-hong. Design and implementation of miniature stepped-frequency ground penetrating radar system[J].Chinese Journal of Scientific Instrument, 2010, 31(3): 699–703.)
[3]	Nicolaescu I. Improvement of stepped-frequency continuous wave ground-penetrating radar cross-range resolution[J].IEEE Transactions on Geoscience and Remote Sensing, 2013, 51(1): 85–92.DOI:10.1109/TGRS.2012.2198069
[4]	Donoho D. Compressed sensing[J].IEEE Transactions on Information Theory, 2006, 52(4): 1289–1306.DOI:10.1109/TIT.2006.871582
[5]	Yang J, Jin T, Huang X, et al. Sparse MIMO array forward-looking GPR imaging based on compressed sensing in clutter environment[J].IEEE Transactions on Geoscience and Remote Sensing, 2014, 52(7): 4480–4494.DOI:10.1109/TGRS.2013.2282308
[6]	Krueger K R, McClellan J H, Scott W R. Efficient algorithm design for GPR imaging of landmines[J].IEEE Transactions on Geoscience and Remote Sensing, 2015, 53(7): 4010–4021.DOI:10.1109/TGRS.2015.2388786
[7]	Ambrosanio M, Pascazio V. A compressive-sensing-based approach for the detection and characterization of buried objects[J].IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2015, 8(7): 3386–3395.DOI:10.1109/JSTARS.2015.2421812
[8]	Qu L L, Yang T H. Investigation of air/ground reflection and antenna beam width for compressive sensing SFCW GPR migration imaging[J].IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(8): 3143–3149.DOI:10.1109/TGRS.2011.2179049
[9]	Wei X, Zhang Y. Interference removal for autofocusing of GPR data from RC bridge decks[J].IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2015, 8(3): 1145–1151.DOI:10.1109/JSTARS.2015.2402211
[10]	Tropp J A, Gilbert A C. Signal recovery from random measurements via orthogonal matching pursuit[J].IEEE Transactions on Information Theory, 2007, 53(12): 4655–4666.DOI:10.1109/TIT.2007.909108
[11]	Solimene R, Cuccaro A, Aversano A, et al. Ground clutter removal in GPR surveys[J].IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2014, 7(3): 792–1151.DOI:10.1109/JSTARS.2013.2287016