CN102879782A - 基于分数阶傅里叶变换的压缩感知sar成像方法 - Google Patents

Abstract

本发明提出一种基于分数阶傅里叶变换的压缩感知SAR成像方法，属于微波成像技术领域，对稀疏目标场景，***发射线性调频脉冲信号，采用将回波信号向低维测量矩阵投影的方式，获取比Nyquist采样定理所需测量数据量少得多的测量数据。结合回波信号在简化分数阶傅里叶变换域稀疏形式和数据获取方式构建重构矩阵之后，通过压缩感知信号重构方法优化求解距离脉压结果，并对距离脉压结果进行距离徙动校正和方位脉压获得目标场景成像结果。本发明所提出的成像方法针对具有稀疏特性的目标场景可实现大幅降低测量数据量，有效缓解数据存储和传输的压力。

Description

基于分数阶傅里叶变换的压缩感知SAR成像方法

技术领域

本发明属于微波成像技术领域，具体涉及一种基于分数阶傅里叶变换的压缩感知合成孔径雷达(SAR)成像的信号处理方法。 

背景技术

合成孔径雷达(SAR)通过发射天线对被观测区域发射电磁脉冲，并利用接收天线接收雷达回波获得目标场景电磁散射特性的空间分布。作为一种有源***，SAR主动发射微波段电磁波而不依赖太阳光照，而且云、雾和降水对微波的影响都不大，确保了SAR全天时全天候成像。在环境保护、灾害监测、海洋观测、资源勘查、地质测绘等方面有着广泛的应用，目前已成为高分辨率对地观测和全球资源管理的重要手段。 

合成孔径雷达是二维高分辨率成像雷达，一维是电磁波波束指向的距离向，一维是与平台运动方向一致的方位向。高分辨率合成孔径雷达提供了丰富的目标场景信息，不仅适用于地形测绘和目标成像，而且也适用于在难度不断增加的监视环境中进行目标检测。因此高分辨率是合成孔径雷达发展的重要方向之一。合成孔径雷达采用综合孔径原理提高方位向分辨率，而距离分辨率的提高则需借助脉冲压缩（简称脉压）技术。由于雷达的理论距离分辨率与信号的带宽成反比，信号带宽越大，分辨率越高。要想获得距离向的高分辨率，就要求雷达的发射信号具有大带宽。对于距离分辨率达到厘米级的雷达来说，信号带宽要达到几GHz。以Nyquist采样定理为基础的现代雷达成像体制面临采样率过高、数据量过大，海量数据存储以及传输难以实现等问题的挑战。处理过程中又存在测量数据和处理结果间数据量上不平衡的问题，即对于比较稀疏的场景仍需按Nyquist采样定理的要求获取庞大的数据，而处理结果却只保留少量位置和散射系数等有用信息。 

压缩感知(Compressed Sensing，简称CS)是近几年发展起来的一种充分利用信号稀疏性或可压缩性的全新信号采集、编解码技术，为解决SAR中采样率过高、数据量过大以及所获取的数据和观测场景及处理结果间不平衡等问题提供了解决方案。压缩感知理论表明，只要信号是稀疏的或可压缩的，或者在某个变换域是稀疏的或可压缩的，就可用一个与变换基不相关的测量矩阵将变换所得高维信号投影到一个低维空间上，测量值并非信号本身，而是从高维到低维的投影值。解码过程不是编码的简单逆过程，而是利用信号稀疏分解中已有的重构方法，通过优化算法在概率意义上从少量的投影中有效恢复和逼近原始信号，所需测量值 的数目远小于传统理论下的样本数。压缩感知理论的详细内容可参阅以下文献： 

[1]D.L.Donoho,“Compressed sensing,”IEEE Trans.Inf.Theory,vol.52,no.4,pp.1289-1306,2006. 

[2]D.L.Donoho,M.Elad,and V.N.Temlyakov,“Stable recovery of sparse overcomplete representations in the presence ofnoise,”IEEE Trans.Inf．Theory,vol.52,no.1,pp.6-18,2006. 

[3]E.J.Candes,and T.Tao,“Decoding by linear programming,”IEEE Trans.Inf.Theory,vol.51,no.12,pp.4203-4215,2005. 

[4]E.J.Candès,“Compressive sampling,”Proc.Int.Congr.Math.,vol.3,pp.1433–1452,2006. 

[5]E.J.Candès,and M.B.Wakin,“An introduction to compressive sampling,”IEEE,Signal Process.Mag.,vol.25,no.2,pp.21-30,2008. 

[6]R.Baraniuk,“Compressive sensing,”IEEE Signal Process.Mag.,vol.24,no.4,pp.118–121,2007 

压缩感知依赖于两个限定条件： 

1.被测量信号具有稀疏性（或可压缩性）或者在某个变换域具有稀疏性（或可压缩性）； 

2.感知矩阵具有不相关性（或约束等距性RIP）。 

当N维目标信号x＝(x₁,x₂,…x_N)^T所依赖的自由度的维数远远小于N时被称为是稀疏的。如果在x中仅有K个非零元(K＜＜N)，则称x为K-稀疏的，K为x的稀疏度。在实际应用中，向量x可能是一个信号f在正交基矩阵Ψ下的系数，即： 

$f\left(n\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{x}_i{\mathrm{\ψ}}_i\left(n\right),n=1,\·\·\·,N,$ 或者f=Ψ^Hx，      (1) 

其中f=[f(1),f(2),…,f(N)]^T,T代表转置变换，Ψ^H是Ψ的共轭转置矩阵。如果f的系数序列x是K-稀疏的，则称f在Ψ域稀疏，Ψ称为稀疏基矩阵；如果x仅有K项较大项，其它项的值都很小，则称f是可压缩的。 

在压缩感知测量中，并不直接测量信号f，而是将f向M×N(M＜＜N)维测量矩阵Φ投影（即f与Φ的行向量做内积），得到一组测量值y_m=<f,ф_m>,m=1,…,M，其中ф_m是Φ的第m个行向量，写成矩阵形式为y=фf。将(1)代入y=фf并考虑到实际中噪声的存在，得 

y=ΦΨ^Hx+e=Θx+e，        (2) 

其中Θ=ΦΨ^H，e为随机性或确定性噪声，噪声能量满足||e||₂≤ε，ε＞0为噪声水平。 

Φ的行向量不能够稀疏表示Ψ的列向量，反之Ψ的行向量不能够稀疏表示Φ的列向量时，Φ和Ψ满足不相关性原则，通常测量矩阵Φ选择随机矩阵能与正交矩阵Ψ在很大程度上满足不相关性原则。当系数向量x是K-稀疏或可压缩的、测量模式服从不相关性原则并且线性测量值y的测量数目满足M=O(Klog(N/K))时，通过求解下面这个优化问题，系数向量x能被高概率地重构。 

$\arg \underset{x}{\min }{\left|\right|x\left|\right|}_{0}s.t.{\left|\right|y-\mathrm{\&Theta;x}\left|\right|}_2\&le;\mathrm{\&epsiv;},---\left(3\right)$

其中||·||₀为所谓的0-范数，即x中非零元素的个数，s.t.表示使得满足，Θ=ΦΨ^H为M×N维的重构矩阵，ε是噪声存在时，优化收敛的门限，与噪声水平相同。 

求解(3)式的算法有凸优化算法/线性规划方法中的去噪基追踪算法(BPDN)、梯度追踪算法(GPSR)以及贪婪算法中的匹配追踪(MP)、正交匹配追踪(OMP)和分段正交匹配追踪(StOMP)等，另外还有一些非凸优化算法。在国防科学技术大学刘吉英的博士论文“压缩感知理论及在成像中的应用”一文中对这些重构算法有较详细论述。 

压缩感知为降低接收机A／D转换速率提供了可能，2007年莱斯大学的R.Baraniuk等人首次将压缩感知引入高分辨雷达，所提雷达***的发射机与传统雷达相同，接收端为一低速率A／D转换器，将昂贵的接收机硬件设计转移到信号恢复算法上。自此以来压缩感知理论在雷达成像中的应用受到越来越多的关注，相应的基础研究逐步展开。 

发明内容

本发明的目的是提出一种基于分数阶傅里叶变换的压缩感知合成孔径雷达(SAR)成像的信号处理方法，对于高分辨率和宽测绘带SAR***，能有效降低测量数据量，缓解信号存储和传输的压力。 

为达到上述目的，本发明所述基于分数阶傅里叶变换的压缩感知SAR成像方法的包括如下步骤： 

一种基于分数阶傅里叶变换的压缩感知SAR成像方法，包括以下几个步骤： 

步骤一：雷达天线向具有稀疏特性的目标场景发射以下线性调频脉冲信号： 

$s_t\left(t\right)=\mathrm{rect}\left(\frac{t}{T_p}\right)\exp \left({j2\mathrm{\&pi;f}}_{c}t\right)\exp \left({-\mathrm{j\&pi;K}}_r{t}^2\right),---\left(4\right)$

其中t为时间变量，T_p、f_c、K_r分别为发射的线性调频脉冲信号的脉冲宽度、载频和调频率，rect(·)为矩形窗函数，定义为

步骤二：构建测量矩阵，获取测量数据； 

发射信号被场景中的目标反射形成回波，将场景中的目标假设为点目标，获得每个信号接收点处的回波信号s⁽ⁿ⁾，其中上角标(n)表示第n个信号接收点，n=1,…,N_a，N_a为信号接收点数目，即方位向采样点数；采用测量矩阵对每个信号接收点处的回波信号进行测量，第n个 信号接收点处的测量结果为M×1维复向量： 

$s_{\mathrm{cs}}^{\left(n\right)}={\mathrm{\&Phi;s}}^{\left(n\right)};$

对回波进行测量所用M×N_r维测量矩阵Φ是M×N_r维随机矩阵Θ与N_r×N_r维简化分数阶傅里叶阵的乘积，即

所述N_r×N_r维简化分数阶傅里叶阵

的第k行第l列元素为： 

${\left[T_F^{\mathrm{\&alpha;}}\right]}_{\mathrm{kl}}=\exp (-j\frac{2\mathrm{\&pi;}(k-(N_r/2)-1)(l-(N_r/2)-1)}{N_r})\&CenterDot;\exp \left(\frac{j}{2}{(l-\frac{N_r}{2}-1)}^2\left(\cot \mathrm{\&alpha;}\right){\left({\mathrm{\&Delta;}}_t\right)}^2\right),---\left(5\right)$

其中α=arccot(-2πK_r)，k,l＝1,…,N_r，Δ_t为Nyquist采样率下采样间隔，即Nyquist采样率的倒数，N_r为以Nyquist采样率对回波进行测量所需测量数量点数；测量数量M由场景的稀疏度K决定，K表征场景中目标数目的稀疏度；M满足K＜M＜＜N_r； 

作为优选，所述的随机矩阵Θ选取M×N_r维的随机高斯矩阵、二值矩阵、部分傅里叶阵或部分哈达玛阵中的一种随机矩阵。 

作为优选，步骤二所述用测量矩阵对回波进行测量获得测量结果前，先以Nyquist采样率对第n个信号接收点处的连续回波信号进行采样，所获得的第n个信号接收点处的离散回波信号s⁽ⁿ⁾的点数为N_r，s⁽ⁿ⁾是N_r×1的行向量，再与测量矩阵Φ做内积获得测量结果

步骤三：通过压缩感知信号重构方法，优化求解距离脉压结果： 

对步骤二获得的每一个接收点处的测量信号

n=1,…,N_a，求解优化问题： 

$\arg \min {\left|\right|s_{\mathrm{ofr}}^{\left(n\right)}\left|\right|}_{0}s.t.{\left|\right|s_{\mathrm{cs}}^{\left(n\right)}-{\mathrm{\&Theta;s}}_{\mathrm{ofr}}^{\left(n\right)}\left|\right|}_2\&le;\mathrm{\&epsiv;}---\left(6\right)$ 即将满足条件

的

作为距离向在简化分数阶傅里叶域脉压的重构结果，将重构结果写成列向量形式并依次排列形成二维信号

N_a为信号接收点数目；其中||·||₀为0-范数，即x中非零元素的个数；||·||₂表示2-范数；s.t.表示使得满足条件；ε是预设的噪声存在时优化收敛的门限； 

作为优选，所述噪声存在时优化收敛的门限ε≥||e⁽ⁿ⁾||₂，e⁽ⁿ⁾＝Φu⁽ⁿ⁾；u⁽ⁿ⁾为第n个信号接收点处噪声信号。 

步骤四：距离徙动校正； 

根据雷达参数包括：场景中心到雷达天线中心的最短距离R₀、雷达发射脉冲重复频率PRF、SAR平台飞行速度v、载频f_c，构建距离徙动校正RMC矩阵T_rmc，T_rmc是N_r×N_a维的矩阵，它的第m列为： 

$T_{\mathrm{rmc}}\left(m\right)={(\exp \left(j\frac{4\mathrm{\&pi;}{R}_0{K}_r}{\mathrm{c\&beta;}\left(m\right)}\mathrm{\&tau;}\left(1\right)\right),\&CenterDot;\&CenterDot;\&CenterDot;,\exp \left(j\frac{4\mathrm{\&pi;}{R}_0{K}_r}{\mathrm{c\&beta;}\left(m\right)}\mathrm{\&tau;}\left(N_r\right)\right))}^T,m=1,\&CenterDot;\&CenterDot;\&CenterDot;,N_a---\left(7\right)$

其中

为方位向离散化频率，τ(1),…,τ(N_r)是在Nyquist采样率下的距离向离散化时间； 

将步骤三获得的二维信号

做二维傅里叶变换，结果为S_2f；然后通过将距离徙动校正RMC矩阵T_rmc与S_2f做点乘来完成距离徙动校正，即距离徙动校正的结果为： 

Srmc=T_rmc.*S_2f， 

其中.*表示矩阵点乘； 

步骤五：根据所选用的方位向脉压方式构建相应的方位向脉压矩阵，对步骤四中经过距离脉压和距离徙动校正的数据S_rmc实施方位向脉压，得到场景最终成像结果。 

作为优选，步骤五采用频域匹配滤波的方法实现方位向脉压，根据雷达参数包括：场景中心距离R₀、雷达发射脉冲重复频率PRF、SAR平台飞行速度v、载频f_c，构建方位向上脉压函数形成方位向脉压矩阵H_a；将步骤五中经过距离脉压和距离徙动校正的数据S_rmc与H_a相乘，得到频域匹配滤波结果为S=S_rmcH_a，最后将S进行二维逆傅里叶变换得到场景最终成像结果。 

进一步地，若步骤五采用频域匹配滤波的方法实现方位向脉压时，方位向脉压矩阵  $H_a=\mathrm{diag}(\exp \left(j\frac{4\mathrm{\&pi;}}{c}{R}_0{f}_c\mathrm{\&beta;}\left(1\right)\right),\&CenterDot;\&CenterDot;\&CenterDot;,\exp \left(j\frac{4\mathrm{\&pi;}}{c}{R}_0{f}_c\mathrm{\&beta;}\left(N_a\right)\right)).$

进一步地，所述步骤二所述第n个信号接收点处的连续回波s⁽ⁿ⁾(t)具有如下形式： 

$s^{\left(n\right)}\left(t\right)=s_0^{\left(n\right)}\left(t\right)+u^{\left(n\right)}\left(t\right)=\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{\mathrm{\&sigma;}}_i\mathrm{rect}\left(\frac{t-({2R}_i^{\left(n\right)}/c)}{T_p}\right)\exp [-\mathrm{j\&pi;}{K}_r{(t-({2R}_i^{\left(n\right)}/c))}^2]+u^{\left(n\right)}\left(t\right);---\left(8\right)$

且其中 $s_0^{\left(n\right)}\left(t\right)=\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{\mathrm{\&sigma;}}_i\mathrm{rect}\left(\frac{t-({2R}_i^{\left(n\right)}/c)}{T_p}\right)\exp [-{\mathrm{j\&pi;K}}_r{(t-({2R}_i^{\left(n\right)}/c))}^2]$ 为第n个信号接收点处无噪回波信号，u⁽ⁿ⁾(t)为第n个信号接收点处噪声信号；σ_i为第i个目标的散射系数，

为第i个目标到第n个信号接收点处雷达天线中心的距离，P为场景中点目标的个数，c为光速。 

对比现有技术，本发明的有益效果在于：提出一种基于简化分数阶傅里叶变换的压缩感知合成孔径雷达成像方法，针对稀疏场景，可大幅降低测量数据量，从而有效缓解信号存储 和传输等压力。 

附图说明

图1是本发明提出的基于分数阶傅里叶变换的压缩感知SAR成像方法实施流程图； 

图2是本发明中点目标成像处理示例的距离向压缩感知脉压结果； 

图3是本发明中点目标成像处理示例的距离徙动校正结果； 

图4是本发明中点目标成像处理示例的点目标成像结果； 

图5是本发明中点目标成像处理示例的点目标距离像的位置与反射率图； 

图6是本发明中点目标成像处理示例的以中间点目标为中心的64×64切片。 

具体实施方式

下面将结合附图及具体实施方式对本发明加以详细说明，需要指出的是，所描述的实施例仅旨在便于对本发明的理解，而不起任何限定作用。 

本发明是一种基于分数阶傅里叶变换的压缩感知SAR成像方法，通过简化分数阶傅里叶变换将目标场景稀疏性转换成线性调频回波信号在分数阶傅里叶域的稀疏性，构造了合成孔径雷达成像的压缩感知模型，并采用压缩感知重构算法直接求解出距离脉压结果。其实施流程如图1所示，包括： 

步骤一：对具有稀疏特性的目标场景，发射线性调频(LFM)脉冲信号； 

步骤二：构建测量矩阵，对回波进行测量； 

步骤三：通过压缩感知信号重构方法，优化求解距离脉压结果； 

步骤四：距离徙动校正； 

步骤五：构建方位向脉压矩阵，实施方位向脉压，得到场景最终成像结果。 

SAR成像一般包括三个步骤，一是距离向脉压，二是距离徙动校正，三是方位脉压。 

在本发明中整体思路是：一、我们在距离向脉压时，采用了压缩感知，与已有的采用压缩感知实现脉压不同的是：1.在测量矩阵方面，我们采用的测量矩阵是“随机矩阵与分数阶傅里叶阵的乘积矩阵”实际上是随机矩阵与稀疏基矩阵的乘积矩阵，只不过稀疏基矩阵就是分数阶傅里叶阵，而现有方法采用的测量矩阵仅是“随机矩阵”；2.在重构矩阵方面，我们采用的是前面用到的“随机矩阵”，而现有方法采用的是“随机矩阵与稀疏基矩阵的乘积矩阵”。二、距离徙动校正时所用距离徙动校正矩阵是专门针对利用分数阶傅里叶阵实现距离向脉压的，这一点是无人在公开文献中用到过的。三、方位向脉压我们采用的是传统脉压方式，此外也想表明在距离向采用了我们的压缩感知脉压和距离徙动校正之后，方位向传统方式的脉压方式仍适用。 

下面对本发明的步骤进行解释，具体而言，本发明包含以下几个步骤： 

步骤一：对具有稀疏特性的目标场景，发射线性调频(LFM)脉冲信号； 

具有稀疏特性的目标场景是指诸如只有几艘船的茫茫大海或只有几架飞机的广阔天空等只含少数目标的场景。雷达天线向具有稀疏特性的场景发射以下线性调频(LFM)脉冲信号 

$s_t\left(t\right)=\mathrm{rect}\left(\frac{t}{T_p}\right)\exp \left({j2\mathrm{\&pi;f}}_{c}t\right)\exp \left({-\mathrm{j\&pi;K}}_r{t}^2\right),---\left(4\right)$

其中t为时间变量，T_p、f_c、K_r分别为发射的线性调频脉冲信号的脉冲宽度、载频和调频率，rect(·)为矩形窗函数，定义为

步骤二：构建测量矩阵，对回波进行测量，获取测量数据； 

发射信号被场景中的目标反射形成回波，不失一般性这里假设场景中的目标为点目标，第n个信号接收点处的连续回波s⁽ⁿ⁾(t)具有如下形式： 

$s^{\left(n\right)}\left(t\right)=s_0^{\left(n\right)}\left(t\right)+u^{\left(n\right)}\left(t\right)=\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{\mathrm{\&sigma;}}_i\mathrm{rect}\left(\frac{t-({2R}_i^{\left(n\right)}/c)}{T_p}\right)\exp [-{\mathrm{j\&pi;K}}_r{(t-({2R}_i^{\left(n\right)}/c))}^2]+u^{\left(n\right)}\left(t\right);---\left(8\right)$

且其中 $s_0^{\left(n\right)}\left(t\right)=\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{\mathrm{\&sigma;}}_i\mathrm{rect}\left(\frac{t-({2R}_i^{\left(n\right)}/c)}{T_p}\right)\exp [-{\mathrm{j\&pi;K}}_r{(t-({2R}_i^{\left(n\right)}/c))}^2]$ 为第n个信号接收点处无噪回波信号，u⁽ⁿ⁾(t)为第n个信号接收点处噪声信号；σ_i为第i个目标的散射系数，

为第i个目标到第n个信号接收点处雷达天线中心的距离，P为场景中点目标的个数，c为光速。 

按照压缩感知理论要求的测量形式——被测回波信号向测量矩阵投影（即被测回波信号与测量矩阵各行向量做内积），对回波进行测量，获得比Nyquist采样所需测量数量少的测量数据。测量数据量M由场景的稀疏度K决定，即与场景中目标数目有关。例如：如果以Nyquist采样率对连续回波s⁽ⁿ⁾(t)进行测量所需测量数量(采样点数)为N_r，表征场景中目标数目的稀疏度为K，则用压缩感知方法对回波s⁽ⁿ⁾(t)进行测量所需测量数量M满足K＜M＜＜N_r。 

所用M×N_r维测量矩阵Φ是M×N_r维随机矩阵Θ与N_r×N_r维简化分数阶傅里叶阵

的乘积矩阵，即

随机矩阵Θ可以选取M×N_r维的随机高斯矩阵、二值矩阵、部分傅里叶阵或部分哈达玛阵等随机矩阵，各随机矩阵的构造在D.L.Donoho等人的文章“Extensions of compressed sensing,”(Signal Processing,vol.86,no.3,pp.549-571,2006)中有详细介绍。N_r×N_r维简化分数阶傅里叶阵

的第k行第l列元素为： 

${\left[T_F^{\mathrm{\&alpha;}}\right]}_{\mathrm{kl}}=\exp (-j\frac{2\mathrm{\&pi;}(k-(N_r/2)-1)(l-(N_r/2)-1)}{N_r})\&CenterDot;\exp \left(\frac{j}{2}{(l-\frac{N_r}{2}-1)}^2\left(\cot \mathrm{\&alpha;}\right){\left({\mathrm{\&Delta;}}_t\right)}^2\right),---\left(5\right)$

其中α=arccot(-2πK_r)，k,l＝1,…,N_r，Δ_t为Nyquist采样率下采样间隔，即Nyquist采样率的倒数，N_r为以Nyquist采样率对回波进行测量所需测量数量点数；K表征场景中目标数目的稀疏度； 

采用测量矩阵对每个信号接收点处的回波信号进行测量，第n个信号接收点处的测量信号

是由压缩感知理论要求的测量形式对回波进行测量的结果，第n个信号接收点处的测量结果

为M×1维复向量，

的矩阵表示形式为：即

其中e⁽ⁿ⁾＝Φu⁽ⁿ⁾。 

需要说明的是，按照压缩感知理论要求的测量形式——被测回波信号向测量矩阵投影，即被测回波信号与测量矩阵各行向量做内积，对回波进行测量；本发明实际操作中，用测量矩阵对回波进行测量获得测量结果时，可先以Nyquist采样率对第n个信号接收点处的连续回波信号s⁽ⁿ⁾(t)进行采样，所获得的第n个信号接收点处的离散回波信号s⁽ⁿ⁾的点数为N_r，s⁽ⁿ⁾是N_r×1的行向量，再与测量矩阵Φ做内积获得测量结果或者直接用每个信号接收点处的连续回波信号s⁽ⁿ⁾(t)向测量矩阵投影，即被测回波信号与测量矩阵各行向量做内积，获得测量结果

上面所述先以Nyquist采样率对第n个信号接收点处的回波信号s⁽ⁿ⁾(t)进行采样的这种操作是因为，压缩感知还是一个相当新的领域，将压缩感知理论要求的测量形式——被测模拟信号向测量矩阵投影应用到实用的传感器上还为时尚早。现在只有概念验证模型出现，其中最著名的是RICE大学研制的单像素相机，他们也有了此方面的专利。 

据我们所知，在SAR领域还没有用压缩感知理论要求的测量形式获得的数据，所以现在实际中是先有了N点离散回波数据，再由离散的回波数据向测量矩阵投影得到M点压缩数据，但压缩感知的初衷并不是这样，实际操作上和理论上有出入，理论上的操作步骤是压缩感知要求的，必须按它的要求才能体现压缩感知的优势(即采样时实现数据压缩)。实际中首先对回波信号采样再进行测量操作只是退而求其次的选择。 

RICE大学研制的单像素相机的关键部件是由德州仪器生产的数字微镜芯片（digital micromirror device：DMD）。在该单像素相机结构中DMD扮演的角色即是随机投影矩阵，实现原始图像每个像素的随机投影观测。相信可用于本发明的实施当中用以直接压缩数据。 

步骤三：通过压缩感知信号重构方法，优化求解距离脉压结果； 

将步骤二中测量到的每一个接收点处的测量信号（其中n=1,…,N_a，N_a为信号接收点 数目），通过压缩感知信号重构方法（基于迭代寻优技术的凸优化/线性规划方法，或贪婪算法等，例如压缩感知重构方法包括凸优化算法中的去噪基追踪(BPDN)、梯度追踪算法(GPSR)以及贪婪算法中的匹配追踪(MP)、正交匹配追踪(OMP)和分段正交匹配追踪(StOMP)等），求解优化问题： 

$\arg \min {\left|\right|s_{\mathrm{ofr}}^{\left(n\right)}\left|\right|}_{0}s.t.{\left|\right|s_{\mathrm{cs}}^{\left(n\right)}-{\mathrm{\&Theta;s}}_{\mathrm{ofr}}^{\left(n\right)}\left|\right|}_2\&le;\mathrm{\&epsiv;}---\left(6\right)$

即将满足条件

的

作为距离向在简化分数阶傅里叶域脉压的重构结果，将重构结果

写成列向量形式并依次排列形成二维信号

N_a为信号接收点数目；其中||·||₀为0-范数，即x中非零元素的个数；||·||₂表示2-范数；s.t.表示使得满足条件；ε是预设的噪声存在时优化收敛的门限；本实施例中，所述噪声存在时优化收敛的门限ε≥||e⁽ⁿ⁾||₂，e⁽ⁿ⁾＝Φu⁽ⁿ⁾；u⁽ⁿ⁾为第n个信号接收点处噪声信号。 

下面结合回波信号稀疏形式和数据获取方式推导重构矩阵： 

对第n个信号接收点处的无噪回波信号

进行简化分数阶傅里叶变换，变换的结果 

$s_{\mathrm{ofr}}^{\left(n\right)}\left(f_r\right)={\left(\frac{1}{j2\mathrm{\&pi;}}\right)}^{\frac{1}{2}}\&CenterDot;\{\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{\mathrm{\&sigma;}}_i\&CenterDot;\exp [-j(\frac{4\mathrm{\&pi;}{f}_c}{c}{R}_i^{\left(n\right)}+\frac{4\mathrm{\&pi;}{f}_r}{c}{R}_i^{\left(n\right)}+\frac{{4\mathrm{\&pi;K}}_r}{c^2}{\left(R_i^{\left(n\right)}\right)}^2)]\&CenterDot;T_p\&CenterDot;\sin c[T_p(f_r+K_r\frac{2R_i^{\left(n\right)}}{c})\left]\right\},---\left(9\right)$

其中f_r为回波

的频率，sinc函数定义为

的极值包含了场景中目标散射系数和目标位置信息，因此将

看作场景的像；即经过对无噪回波信号的简化分数阶傅里叶变换，将场景的稀疏性转换为回波信号在简化分数阶傅里叶域的稀疏性； 

在连续信号的情况下，所述简化分数阶傅里叶变换算子

的作用为： 

$T_F^{\mathrm{\&alpha;}}\left[x\left(t\right)\right]\left(f\right)={\left(\frac{1}{j2\mathrm{\&pi;}}\right)}^{\frac{1}{2}}\&times;{\&Integral;}_{-\&infin;}^{\&infin;}\exp (-j2\mathrm{\&pi;ft}+\frac{j}{2}{t}^2\cot \mathrm{\&alpha;})\&CenterDot;x\left(t\right)\mathrm{dt};$

公式(6)获得的场景的像为

无噪回波的简化分数阶傅里叶变换过程即无噪回波信号稀疏形式为

含噪回波信号

结合测量矩阵对步骤二得到的第n个信号接收点处的测量信号进一步展开，得到压缩感知方法下对第n个信号接收点处的测量信号

与稀疏场景像之间的关系满足： 

$s_{\mathrm{cs}}^{\left(n\right)}={\mathrm{\&Phi;s}}^{\left(n\right)}={\mathrm{\&Phi;s}}_0^{\left(n\right)}+{\mathrm{\&Phi;u}}^{\left(n\right)}={\mathrm{\&Theta;T}}_F^{\mathrm{\&alpha;}}{s}_0^{\left(n\right)}+{\mathrm{\&Theta;T}}_F^{\mathrm{\&alpha;}}{u}^{\left(n\right)}={\mathrm{\&Theta;T}}_F^{\mathrm{\&alpha;}}{s}_0^{\left(n\right)}+e^{\left(n\right)}={\mathrm{\&Theta;s}}_{\mathrm{ofr}}^{\left(n\right)}+e^{\left(n\right)}$

其中e⁽ⁿ⁾＝Φu⁽ⁿ⁾；由压缩感知理论得知对稀疏场景像

进行重构的重构矩阵Θ，即步骤二中使用的随机矩阵Θ； 

步骤四：距离徙动校正； 

由于压缩感知重构具有保持复数相位不变的保相性，因此由步骤三得到的压缩感知重构之后的距离脉压信号仍可进行距离徙动校正(RMC)。 

根据雷达参数包括：场景中心到雷达天线中心的最短距离R₀、雷达发射脉冲重复频率PRF、SAR平台飞行速度v、载频f_c，构建距离徙动校正RMC矩阵T_rmc，T_rmc是N_r×N_a维的矩阵，它的第m列为： 

$T_{\mathrm{rmc}}\left(m\right)={(\exp \left(j\frac{{4\mathrm{\&pi;R}}_0{K}_r}{\mathrm{c\&beta;}\left(m\right)}\mathrm{\&tau;}\left(1\right)\right),\&CenterDot;\&CenterDot;\&CenterDot;,\exp \left(j\frac{{4\mathrm{\&pi;R}}_0{K}_r}{\mathrm{c\&beta;}\left(m\right)}\mathrm{\&tau;}\left(N_r\right)\right))}^T,m=1,\&CenterDot;\&CenterDot;\&CenterDot;,N_a---\left(7\right)$

其中

为方位向离散化频率，τ(1),…,τ(N_r)是在Nyquist采样率下的距离向离散化时间， $[\mathrm{\&tau;}\left(1\right),\&CenterDot;\&CenterDot;\&CenterDot;,\mathrm{\&tau;}\left(N_r\right)]=(-\frac{N_r}{2},-\frac{N_r}{2}+1,-\frac{N_r}{2}+2,\&CenterDot;\&CenterDot;\&CenterDot;,\frac{N_r}{2}-1)\&CenterDot;{\mathrm{\&Delta;}}_t,$ 其中Δ_t为Nyquist采样率下距离向采样间隔； 

将步骤三获得的二维信号

做二维傅里叶变换，结果为S_2f；然后通过将距离徙动校正RMC矩阵T_rmc与S_2f做点乘来完成距离徙动校正，即距离徙动校正的结果为： 

S_rmc=T_rmc.*S_2f， 

其中.*表示矩阵点乘； 

步骤五：根据所选用的方位向脉压方式构建相应的方位向脉压矩阵，对步骤四获得的经过距离脉压和距离徙动校正的数据S_rmc实施方位向脉压，得到场景最终成像结果。 

方位向脉压可采用任何传统脉压方式，即方位脉压矩阵H_a可以是频域匹配滤波阵，也可以是时域匹配滤波阵，还可以是Deramp-FFT（去斜-快速傅里叶变换）阵。最常用的如匹配滤波。根据雷达参数（场景中心距离R₀、雷达发射脉冲重复频率PRF、平台飞行速度v、载频f_c等）构建方位向上脉压函数形成方位向脉压矩阵H_a。以频域匹配滤波为例：如果采用频域匹配滤波的方法实现方位向脉压，那么方位向脉压矩阵  $H_a=\mathrm{diag}(\exp \left(j\frac{4\mathrm{\&pi;}}{c}{R}_0{f}_c\mathrm{\&beta;}\left(1\right)\right),\&CenterDot;\&CenterDot;\&CenterDot;,\exp \left(j\frac{4\mathrm{\&pi;}}{c}{R}_0{f}_c\mathrm{\&beta;}\left(N_a\right)\right)).$ 将步骤五中经过距离脉压和距离徙动校正 的数据S_rmc与H_a相乘，得到频域匹配滤波结果为S=S_rmcH_a，最后将S进行二维逆傅里叶变换得到场景最终成像结果。 

以上仅是以频域匹配滤波为例进行说明，实际上任何传统方位向脉压形式在这里都可以应用，即方位脉压矩阵H_a可以是频域匹配滤波阵，也可以是时域匹配滤波阵，还可以是Deramp-FFT（去斜-快速傅里叶变换）阵（Deramp-FFT详见D.C.Munson Jr,andR.L.Visentin,“A signal processing view of strip-mapping synthetic aperture radar,”Acoustics,Speech and Signal Processing,IEEE Transactions on,vol.37,no.12,pp.2131-2147,1989）。 

实施例：为了对本发明进行验证，进行了点目标计算机仿真，仿真参数如表1所示。 

表1：SAR参数表2：点目标参数 

步骤一：按照表1参数向场景发射线性调频脉冲信号 $s_t\left(t\right)=\mathrm{rect}\left(\frac{t}{T_p}\right)\exp \left(j2\mathrm{\&pi;}{f}_{c}t\right)\exp (-\mathrm{j\&pi;}{K}_r{t}^2),$ 其中t为时间变量，K_r＝B／T_p，对位于场景内的三个点目标进行了回波仿真。点目标的参数如表2所示。 

步骤二：将回波向测量矩阵投影得到测量值（其中n=1,…,N_a，N_a为信号接收点数目），即

其中Φ为测量矩阵，本实施例中Φ为随机高斯阵Θ与简化分数阶傅里叶阵

的乘积阵（

具体形式参见具体实施方式中步骤二(5)式）。 

步骤三：距离向脉压问题转化为求解

的优化问题。利用压缩感知重构方法重构距离向在简化分数阶傅里叶域脉压结果

压缩感知重构方法包括 凸优化算法中的去噪基追踪(BPDN)、梯度追踪算法(GPSR)以及贪婪算法中的匹配追踪(MP)、正交匹配追踪(OMP)和分段正交匹配追踪(StOMP)，本实施例中使用的重构方法是OMP算法。OMP算法是一个很成熟的算法，很容易实现，因此不作为重点介绍。将重构的距离向脉压结果写成列向量形式依次排列形成二维信号

N_a为信号接收点数目（即方位向采样点数），仿真结果如图2所示。 

步骤四：对二维信号

做二维傅里叶变换得到S_2f。将距离徙动校正(RMC)矩阵T_rmc与S_2f点乘实现距离徙动校正。RMC矩阵T_rmc的具体形式如具体实施方式步骤四中(7)式。校正结果S_rmc如图3所示。 

步骤五：根据雷达参数（场景中心距离R₀、雷达发射脉冲重复频率PRF、平台飞行速度v、载频f_c等）构建方位向上脉压函数形成方位脉压矩阵，对经过距离脉压和距离徙动校正的数据S_rmc进行方位向脉压。本实施例中采用频域匹配滤波矩阵  $H_a=\mathrm{diag}(\exp \left(j\frac{4\mathrm{\&pi;}}{c}{R}_0{f}_c\mathrm{\&beta;}\left(1\right)\right),\&CenterDot;\&CenterDot;\&CenterDot;,\exp \left(j\frac{4\mathrm{\&pi;}}{c}{R}_0{f}_c\mathrm{\&beta;}\left(N_a\right)\right)),$ 其中 $\mathrm{\&beta;}\left(m\right)=\sqrt{1-\frac{c^2{f}_a^2\left(m\right)}{4v^2{f}_c^2}},$ $f_a\left(m\right)=(m-\frac{N_a}{2})\frac{\mathrm{PRF}}{N_a},$ m＝1,…,N_a为方位向离散化频率，PRF为雷达发射脉冲重复频率，v为SAR平台飞行速度。频域匹配滤波结果为S＝S_rmcH_a，最后将S进行二维逆傅里叶变换得到场景成像结果s，在本实施例中为三个点目标的像，如图4所示。像的位置与预设点目标的位置相同，反射率非常接近，图5显示了点目标距离像的位置和反射率。图6为点目标1为中心的64×64切片图，方位向主瓣宽度为0.82m，距离向主瓣宽度为0.72m。本实施例中最终通过仅使用10%数据，恢复出相当于常规SAR由全数据处理得到的成像结果。 

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉该技术的人在本发明所揭露的技术范围内，可理解想到的变换和替换，都应涵盖在本发明的包含范围之内，因此，本发明的保护范围应该以权利要求书的保护范围为准。 

Claims (7)

1.一种基于分数阶傅里叶变换的压缩感知SAR成像方法，其特征在于，包括以下步骤： 

步骤一：雷达天线向具有稀疏特性的目标场景发射以下线性调频脉冲信号： 

其中

t为时间变量，T_p、f_c、K_r分别为发射的线性调频脉冲信号的脉冲宽度、载频和调频率，rect(·)为矩形窗函数，定义为

步骤二：构建测量矩阵，获取测量数据； 

发射信号被场景中的目标反射形成回波，将场景中的目标假设为点目标，获得每个信号接收点处的回波信号s⁽ⁿ⁾，其中上角标(n)表示第n个信号接收点，n=1,…,N_a，N_a为信号接收点数目，即方位向采样点数；采用测量矩阵对每个信号接收点处的回波信号进行测量，第n个信号接收点处的测量结果

为M×1维复向量： 

对回波进行测量所用M×N_r维测量矩阵Φ是M×N_r维随机矩阵Θ与N_r×N_r维简化分数阶傅里叶阵

的乘积，即

所述N_r×N_r维简化分数阶傅里叶阵

的第k行第l列元素为： 

其中α=arccot(-2πK_r)，k,l＝1,…,N_r，Δ_t为Nyquist采样率下采样间隔，即Nyquist采样率的倒数，N_r为以Nyquist采样率对回波进行测量所需测量数量点数；测量数量M由场景的稀疏度K决定，K表征场景中目标数目的稀疏度；M满足K＜M＜＜N_r； 

步骤三：通过压缩感知信号重构方法，优化求解距离脉压结果： 

对步骤二获得的每一个接收点处的测量信号n=1,…,N_a，求解优化问题： 

即将满足条件

的

作为距离向在简化分数阶傅里叶域脉压的重构结果，将重构结果写成列向量形式并依次排列形成二维信号

N_a为信号接收点数 目；其中||·||₀为0-范数，即x中非零元素的个数；||·||₂表示2-范数；s.t.表示使得满足条件；ε是预设的噪声存在时优化收敛的门限； 

步骤四：距离徙动校正； 

根据雷达参数包括：场景中心到雷达天线中心的最短距离R₀、雷达发射脉冲重复频率PRF、SAR平台飞行速度v、载频f_c，构建距离徙动校正RMC矩阵T_rmc，T_rmc是N_r×N_a维的矩阵，它的第m列为： 

其中

为方位向离散化频率，τ(1),…,τ(N_r)是在Nyquist采样率下的距离向离散化时间； 

将步骤三获得的二维信号

做二维傅里叶变换，结果为S_2f；然后通过将距离徙动校正RMC矩阵T_rmc与S_2f做点乘来完成距离徙动校正，即距离徙动校正的结果为： 

S_rmc=T_rmc.*S_2f， 

其中.*表示矩阵点乘； 

步骤五：根据所选用的方位向脉压方式构建相应的方位向脉压矩阵，对步骤四获得的经过距离脉压和距离徙动校正的数据S_rmc实施方位向脉压，得到场景最终成像结果。 

2.根据权利要求1所述一种基于分数阶傅里叶变换的压缩感知SAR成像方法，其特征在于，步骤二所述第n个信号接收点处的回波s⁽ⁿ⁾(t)具有如下形式： 

且其中

为第n个信号接收点处无噪回波信号，u⁽ⁿ⁾(t)为第n个信号接收点处噪声信号；σ_i为第i个目标的散射系数，为第i个目标到第n个信号接收点处雷达天线中心的距离，P为场景中点目标的个数，c为光速。 

3.根据权利要求1所述一种基于分数阶傅里叶变换的压缩感知SAR成像方法，其特征在于，步骤三中，所述噪声存在时优化收敛的门限ε≥||e⁽ⁿ⁾||₂，

u⁽ⁿ⁾为第n个信号接 收点处噪声信号。 

4.根据权利要求1所述一种基于分数阶傅里叶变换的压缩感知SAR成像方法，其特征在于，步骤二所述的随机矩阵Θ选取M×N_r维的随机高斯矩阵、二值矩阵、部分傅里叶阵或部分哈达玛阵中的一种随机矩阵。 

5.根据权利要求1-4所述任一种基于分数阶傅里叶变换的压缩感知SAR成像方法，其特征在于，步骤五采用频域匹配滤波的方法实现方位向脉压，根据雷达参数包括：场景中心距离R₀、雷达发射脉冲重复频率PRF、SAR平台飞行速度v、载频f_c，构建方位向上脉压函数形成方位向脉压矩阵H_a；将步骤五中经过距离脉压和距离徙动校正的数据S_rmc与H_a相乘，得到频域匹配滤波结果为S=S_rmcH_a，最后将S进行二维逆傅里叶变换得到场景最终成像结果。 

6.根据权利要求5所述一种基于分数阶傅里叶变换的压缩感知SAR成像方法，其特征在于，采用频域匹配滤波的方法实现方位向脉压时，方位向脉压矩阵 

7.根据权利要求1-4所述任一种基于分数阶傅里叶变换的压缩感知SAR成像方法，其特征在于，步骤二所述用测量矩阵对回波进行测量获得测量结果前，先以Nyquist采样率对第n个信号接收点处的连续回波信号进行采样，所获得的第n个信号接收点处的离散回波信号s⁽ⁿ⁾的点数为N_r，s⁽ⁿ⁾是N_r×1的行向量，再与测量矩阵Φ做内积获得测量结果

Images