一种斜视星机双基地合成孔径雷达成像方法

摘要

本发明公开了一种斜视星机双基地合成孔径雷达成像方法，它是利用斜视SA-BSAR系统的特点，通过该系统二维频谱的解析表达式可以得到二维空间波数映射关系，并利用二维非均匀快速傅里叶变换高效实现成像处理，本发明避免了斜视SA-BSAR系统中的二维空间域与二维频域的解耦问题，而且未进行小场景假设，因此该方法适用于斜视角较大、成像场景范围较大的SA-BSAR系统高性能成像。本发明兼具较高的成像性能和运算效率。本发明填补了现有的SA-BSAR频域成像技术无法应用于斜视SA-BSAR高分辨率成像的空白。

说明书

技术领域

本发明属于雷达技术领域，它特别涉及星机双基地合成孔径雷达（简称： SA-BSAR）的高效成像技术。

背景技术

星机双基地SAR（简称SA-BSAR）系统采用星载平台发射信号、机载平台接 收目标回波信号的工作模式实现对目标区域成像。SA-BSAR不但具有良好的隐 蔽性、强的抗干扰能力和战场生存能力，其独特的“远发近收”工作模式还使其具 有独特的优势：系统机动灵活、目标信息获取丰富、抗干扰能力强、空间分辨能 力高、系统造价低。因此，SA-BSAR具有很高的研究价值。德国宇航局（FGAN） 高频物理与雷达技术研究所（FHR）于2009年11月实施了世界上第一次SA-BSAR 系统试验，利用雷达卫星TerraSAR-X和机载雷达PAMIR分别作为系统的发射、 接收系统，验证了SA-BSAR二维成像的可行性。

由于SA-BSAR系统发射机、接收机两平台间速度相差大，为保证发射、接收 波束的空间同步，发射机、接收机往往需要工作在斜视模式。然而，不同于传统 单基地SAR，SA-BSAR系统回波具有二维空变特性，即斜距空变特性（range  variant）和方位空变特性（azimuth variant），致使回波信号的二维空间、二维频 率之间产生了复杂的交互作用，难以对其解耦为两个一维空变，因此传统单基地 SAR斜视频域成像算法无法应用于斜视SA-BSAR。

目前已公开发表的SA-BSAR频域成像文献中，均是通过对二维频域一阶泰勒 展开的方法实现对二维空变的解耦，典型的有文献1：Wang R,Loffeld O,Nies H, et al.Focusing spaceborne/airborne hybrid bistatic SAR data using  wavenumber-domain algorithm.IEEE Transactions on Geoscience and Remote  Sensing,2009,47(7):2275-2283.；文献2：Zhe Liu,Jianyu Yang,Xiaoling Zhang, Yiming Pi Study on spaceborne/airborne hybrid bistatic SAR image formation in  frequency domain.IEEE Geoscience and Remote Sensing,2008,5(4):578-582。这些 成像算法适用于正侧视或斜视角很小的SA-BSAR系统，在斜视角较大（如大于 10度）的情况下，这些成像算法会出现较大的误差，引起成像结果散焦和几何失 真，因而这些成像算法也就无法适用于斜视SA-BSAR系统。此外，文献3：C.Y.Dai, X.L.Zhang,Bistatic Polar Format Algorithm Based on NUFFT Method[J].J.of  Electromagn.Waves and Appl.,Vol.25,2328–2340,2011，提出利用Polar Format  Algorithm结合非均匀快速傅里叶变换（Non-Uniform Fast Fourier Transform,简写 NUFFT）方法实现任意构型BSAR系统的成像，然而该算法只适用于小场景的高 精度成像，因此应用受限。

发明内容

本发明的目的是克服现有频域成像技术无法应用于斜视SA-BSAR的不足，提 供了一种适用于斜视SA-BSAR系统的成像方法，该成像方法充分考虑了斜视 SA-BSAR的系统特点，利用系统回波二维频谱解析式，通过二维非均匀傅里叶 变换（Two-Dimensional NUFF T,简写2D-NUFFT），实现高精度频域成像。该方 法无需通过对二维频域一阶泰勒展开的方法实现二维空变特性进行解耦，而且未 进行小场景假设，因此该方法适用于斜视角较大、成像场景范围较大的SA-BSAR 系统高性能成像。

为了方便描述本发明的内容，首先作以下术语定义：

定义1、SA-BSAR系统相关参数描述

星载平台斜距史$r_S\left(t\right)=\sqrt{r_{0S}^2+v_S^2{(t-t_{0S})}^2}$

机载平台斜距史$r_P\left(t\right)=\sqrt{r_{0P}^2+v_P^2{(t-t_{0P})}^2}$

SA-BSAR系统斜距史$r\left(t\right)=\sqrt{r_{0S}^2+v_S^2{(t-t_{0S})}^2}+\sqrt{r_{0P}^2+r_P^2{(t-t_{0P})}^2}$

星载平台相位史φ_S(t)=k·r_S(t)

机载平台相位史φ_P(t)=k·r_P(t)

SA-BSAR系统相位史φ(t)=k·r(t)+2πf_dt

SA-BSAR系统成像结果坐标系(r,x)，其中x＝v_S·t_0S。

其他参数：τ为快（斜距）时间，t为慢（方位）时间；v_S，v_P分别是星载和 机载平台相对目标的运动速度大小；星载和机载平台分别在t_0S，t_0P时刻距目标 最近，且最近斜距分别为r_0S，r_0P；θ_S，θ_P分别是星载、机载平台的斜视角；f为 对应于快（斜距）时间的频率，ω=2πf为对应于快（斜距）时间的角频率，ω₀为 发射信号中心角频率，c为光速，f_d为对应于慢（方位）时间的多普 勒频率。

定义2、air-phase时间点t_b

air-phase时间点t_b为满足下式成立的时间点：kr′_S(t_b)+2πf_d=0，求解该方程 可以得到t_b的解析解：

$t_b=t_{0S}-\frac{f_d{r}_{0S}}{\sqrt{{\left(\frac{f+f_0}{c}{v}_S\right)}^2-f_d^2}}---\left(1\right)$

在该时间点，系统相位史的变化率φ′(t_b)就等于机载平台的相位变化率kr′_P(t_b)， 因此称该时间点为air-phase时间点。其中r′_S(t_b)和r′_P(t_b)分别是星载、机载平台斜距 史关于时间的一阶导数在t_b时刻的值：

定义3、SA-BSAR系统冲激响应的二维频谱

根据SA-BSAR系统的特点，定义SA-BSAR系统冲激响应的二维频谱 H(f,f_d)为

$H(f,f_d)=\exp \{-\mathrm{j\phi }\left(t_b\right)\}·\exp \left\{j\frac{1}{2}\frac{{\left({\phi }^{\prime }\left(t_b\right)\right)}^2}{{\phi }^{\prime \prime }\left(t_b\right)}\right\}---\left(2\right)$

其中φ′(t_b)，φ″(t_b)分别是系统相位史φ(t)关于时间t的一阶、二阶导数在air-phase 时间点t_b的值：${\phi }^{\prime }\left(t_b\right)={\frac{\mathrm{d\phi }\left(t\right)}{\mathrm{dt}}|}_{t=t_b},$${\phi }^{\prime \prime }\left(t_b\right)={\frac{d^2\mathrm{d\phi }\left(t\right)}{d{t}^2}|}_{t=t_b}.$

定义4、二维非均匀快速傅里叶变换

二维信号Z的离散非均匀傅里叶变换为

$z(p,q)=\underset{m=-M/2}{\overset{M/2}{\Sigma }}\underset{n=-N/2}{\overset{N/2}{\Sigma }}Z(m,n)·\exp \left\{j2\pi \right[\frac{m_{b}p}{M}+\frac{n_{b}q}{N}\left]\right\}---\left(3\right)$

其中，p,q,m,n均为等间隔整数，p=-M/2,...,M/2-1，q=-N/2,...,N/2-1，M,N均 为正整数，m_b,n_b是标识傅里叶变换的非均匀位置，m_b,n_b是非等间隔实数， m_b∈[-M/2,...,M/2-1]，n_b∈[-N/2,...,N/2-1]。

对二维离散非均匀傅里叶变换的快速实现，就是二维非均匀快速傅里叶变换 （2D-NUFFT）。2D-NUFFT是通过对傅里叶指数基底在一个过采样邻域内进行插 值实现的，详情参见参考文献：A.Dutt,V.Rokhlin.“Fast Fourier transforms for  nonequispaced data”.SIAM Journal on Scientific Computing,vol.14,no.6,pp. 1368-1393,1993.

定义5、频域距离压缩方法

设发射信号其中T是发射信号的脉冲宽度，μ是发 射信号的调频斜率，是矩形窗，目标回波信号为 $\upsilon \left(\tau \right)=\mathrm{rect}\left(\frac{\tau -{\tau }_0}{T}\right)·\exp \left(\mathrm{j\pi \mu }{(\tau -{\tau }_0)}^2\right),$其中τ₀是回波时延。

频域距离压缩方法为

$V_{\mathrm{RC}}\left(f\right)=V\left(f\right)·V_0^{*}\left(f\right)---\left(4\right)$

其中V_RC(f)是距离压缩后信号傅里叶变换的结果，V(f)是目标回波信号υ(τ)傅 里叶变换后的结果，是发射信号υ₀(τ)经过傅里叶变换后的复共轭,详情参 见参考文献：保铮，邢孟道，王彤等，雷达成像技术，电子工业出版社，2004 年。

定义6、最近邻域插值方法

设信号γ(u,v)，其中u,v是整数。对信号γ(u,v)的最近邻域插值为

γ(u′,v′)=γ(uⁿ,vⁿ)                （5）

其中u′,v′是非整数，uⁿ,vⁿ是对非整数u′,v′取整后的结果。

本发明提供一种斜视星机双基地合成孔径雷达频域成像方法，如图1所示， 包含如下步骤：

步骤1、回波信号距离压缩

将斜视星机双基地合成孔径雷达系统的原始回波数据s(τ,t)存放在一个M行 N列的回波数据矩阵中，其中τ和t分别是快时间和慢时间，M和N均为正整数， 回波数据矩阵的每列数据是存放慢时间回波信号的采样，回波数据矩阵的每行数 据是存放快时间回波信号的采样；

采用传统频域距离压缩的方法，对回波数据矩阵s(τ，t)进行距离压缩，得到 经过距离压缩后的信号其中*表示复共轭， S₀(f)是发射信号s₀(τ)经过傅里叶变换后的结果，f是对应于快时间的频率；

步骤2、方位向傅立叶变换

针对步骤一中得到的距离压缩后的回波信号数据矩阵中的每一列， 做傅立叶变换，得到经过距离压缩后信号的二维频谱其中f_d是对应于 慢时间的频率；

步骤3、参考点相位补偿

将步骤二中得到距离压缩后信号的二维频谱数据矩阵，与参考点目 标二维频谱的复共轭H₀(f,f_d)逐点相乘，得到经过参考点相位补偿后信号的二维 频谱其中H₀(f,f_d)是参考点目标二维 频谱，H₀(f,f_d)是利用公式（2）：$H(f,f_d)=\exp \{-\mathrm{j\phi }\left(t_b\right)\}·\exp \left\{j\frac{1}{2}\frac{{\left({\phi }^{\prime }\left(t_b\right)\right)}^2}{{\phi }^{\prime \prime }\left(t_b\right)}\right\}$得到的， 其中公式（2）中，t_b是参考点目标P₀(r₀,x₀)的air-phase时间点： 其中t_0S0是卫星平台距离距参考点目标P₀(r₀,x₀)最近的 时刻，r_0S0表示卫星平台距离参考点目标的最近斜距，φ(t)是SA-BSAR系统的系 统相位史，φ′(t_b)和φ″(t_b)分别是系统相位史φ(t)关于时间t的一阶、二阶导数在 air-phase时间点t_b的值，f₀为发射信号中心频率，v_S是卫星平台相对目标的运动 速度大小，f_d为对应于慢时间的多普勒频率，f为对应于快时间的频率，c是光 速大小；

步骤4、二维非均匀快速傅里叶变换

对步骤三中得到的经过参考点相位补偿后信号的二维频谱数据矩 阵，做二维非均匀快速傅里叶变换，得到二维空域图像d(p,q)： $d(p,q)=\underset{m=-M/2}{\overset{M/2}{\Sigma }}\underset{n=-N/2}{\overset{N/2}{\Sigma }}{S}_{\mathrm{RF}}^{\tau ,t}(m,n)·\exp \left\{j2\pi \right[\frac{m_{b}p}{M}+\frac{n_{b}q}{N}\left]\right\},$其中m_b,n_b是傅里叶变换的非均匀 位置，m_b,n_b分别由下式得到

$m_b=\frac{f_{x}M}{f_{\mathrm{xs}}},$$n_b=\frac{f_{y}N}{f_{\mathrm{ys}}}---\left(6\right)$

其中f_xs=max{f_x}-min{f_x},f_ys=max{f_y}-min{f_y}， p=-M/2,...,M/2-1，q=-N/2,...,N/2-1，P₀表示参考点目标的位置；

步骤5、空间重采样

对步骤四中得到的二维空域图像d(p,q)，利用传统最近邻域插值方法进行二 维空间重采样，得到成像结果g(p,q)： g(p,q)=d(p+a₁p²+a₂q²+a₃pq,q+b₁p²+b₂q²+b₃pq)，其中，空间重采样参数 a₁□a₃,b₁□b₃分别由下式确定

$a_1=\frac{{\partial f}_{x2}}{{\partial f}_x},$$a_2=\frac{{\partial f}_{y2}}{{\partial f}_x},$$a_3=\frac{{\partial f}_{\mathrm{xy}}}{{\partial f}_x},$$b_1=\frac{{\partial f}_{x2}}{{\partial f}_y},$$b_2=\frac{{\partial f}_{y2}}{{\partial f}_y},$$b_3=\frac{{\partial f}_{\mathrm{xy}}}{{\partial f}_y}---\left(7\right)$

其中，$f_{x2}={\frac{{\partial }^2\phi \left(t_b\right)}{{\partial x}^2}|}_{P_0},$$f_{y2}={\frac{{\partial }^2\phi \left(t_b\right)}{{\partial r}^2}|}_{P_0},$$f_{\mathrm{xy}}={\frac{{\partial }^2\phi \left(t_b\right)}{\partial r\partial x}|}_{P_0},$P₀表示参考点目标的位置；

对于步骤五得到的成像结果g(p,q)作坐标变换：得到具有 较高分辨率的目标成像结果g(r,x)。

需要说明的是：

经过步骤一至步骤四的处理，成像场景内非参考点目标回波中仍然存在未被 补偿的部分，会导致成像结果出现几何畸变，无法满足高分辨率成像的要求，因 此需要步骤五进一步的处理。

本发明的实质是利用斜视SA-BSAR回波二维频谱解析式，以及二维非均匀 傅里叶快速实现方法，推导得到一种适用于该系统的高效、高精度成像方法。根 据回波二维频谱的解析表述，获取SA-BSAR的二维空间波数映射关系；最后通 过研究该二维空间波数映射关系，利用二维非均匀快速傅里叶变换，实现对斜视 SA-BSAR回波的聚焦成像。

本发明的创新点在于从SA-BSAR二维频谱解析式中，合理提取能反映 SA-BSAR回波二维空变特性的二维空间波数映射关系，并利用这种映射关系， 通过二维非均匀快速傅里叶变换，高效实现斜视SA-BSAR的成像处理。

本发明的基本原理是针对斜视SA-BSAR系统回波二维空变的特点，利用 SA-BSAR系统的几何模型及回波二维频谱解析模型，分析系统回波二维空域波 数与二维频率之间的耦合关系，利用二维非均匀快速傅里叶变换，有效补偿成像 场景内非参考点目标二维空变。

本发明解决的技术问题：现有SA-BSAR频域成像方法无法解决斜视工作模 式下二维空间域与二维频域的准确解耦，因此在斜视角较大情况下（如斜视角大 于10度），成像结果会散焦；本发明利用斜视SA-BSAR系统的特点，通过该系统 二维频谱的解析表达式可以得到二维空间波数映射关系，并利用二维非均匀快速 傅里叶变换高效实现成像处理，本发明避免了斜视SA-BSAR系统中的二维空间 域与二维频域的解耦问题，而且未进行小场景假设，因此该方法适用于斜视角较 大、成像场景范围较大的SA-BSAR系统高性能成像。

本发明的有益效果：充分利用SA-BSAR系统的特点，提取能反映回波二维 空间波数映射关系，避免了因斜视SA-BSAR系统中的二维空间域与二维频域的 解耦带来的误差，此外，本发明通过二维非均匀快速傅里叶变换实现成像处理， 因此本发明适用于斜视工作模式下的SA-BSAR高精度成像，本发明兼具较高的 成像性能和运算效率。本发明填补了现有的SA-BSAR频域成像技术无法应用于 斜视SA-BSAR高分辨率成像的空白。

附图说明

图1为本发明的工作流程框图。

图2为仿真成像场景内九个点目标原始相对位置关系图，

其中P₀是参考点目标，P₁□P₈为八个非参考点目标；

图3为斜视星机双基地SAR系统回波数据s(τ,t)经过本发明步骤一至步骤三处 理后的成像结果示意图；其中，P₀是参考点目标，P₁□P₈为八个非参考点目标， 横轴r表示斜距向，纵轴x表示方位向，m是长度单位米；

图4为斜视星机双基地SAR系统回波数据s(τ,t)经过本发明步骤一至步骤五处 理后的成像结果示意图；其中，P₀是参考点目标，P₁□P₈为八个非参考点目标， 横轴r表示斜距向，纵轴x表示方位向，m是长度单位米；

图5为本发明实施例中采用的斜视星机双基地SAR系统平台参数。

具体实施方式

本发明主要采用仿真实验的方法进行验证，所有步骤、结论都在 MATLAB7.11（R2010b）上验证正确。

本实施例分别采用TerraSAR-X卫星和机载PAMIR作为发射、接收平台，两 平台平行同向飞行，发射、接收天线斜视角分别为29.1°和10°，发射、接收 天线波束速度分别为2394.8m/s和400m/s，发射信号的中心频率为9.65GHz，发 射信号带宽为100MHz，脉冲重复频率为3800Hz，回波信号距离向采样频率为 200MHz，其他系统平台仿真参数如图5所示，仿真成像场景内包含9个点目标， 其相对位置关系如图2所示，其中参考点目标为P₀位于场景中心。

步骤一、回波信号距离压缩

将斜视星机双基地合成孔径雷达系统回波数据s(τ,t)以一个3665行4096列 的数据矩阵存放，其中每列数据是存放慢时间回波信号的采样；每行的数据是存 放快时间单脉冲回波信号的采样。

把距离压缩参考信号s₀(τ)作傅立叶变换得到参考信号频谱S₀(f)，把回波信 号s(τ,t)逐行做傅立叶变换得到S^τ(f,t)，将S^τ(f,t)逐行与S₀(f)共轭相乘得到 实现距离压缩。傅立叶变换可以通过快速傅立叶变换（Fast Fourier  Transform,简称FFT）实现；

步骤二、方位向傅立叶变换

针对距离压缩后的回波数据矩阵的每一列做FFT得到距离压缩后的 二维频谱

步骤三、参考点相位补偿

选择成像场景中心点P₀为参考点目标，卫星、飞机平台距参考点目标P₀的最 近斜距分别为635.34km，3.71km，卫星平台距参考点目标P₀最近的时刻为48.73 s。利用图5所示的系统参数及公式（1）可以得到该系统

关于参考点目标的air-phase点t_b，根据公式（2） $H(f,f_d)=\exp \{-\mathrm{j\phi }\left(t_b\right)\}·\exp \left\{j\frac{1}{2}\frac{{\left({\phi }^{\prime }\left(t_b\right)\right)}^2}{{\phi }^{\prime \prime }\left(t_b\right)}\right\}$可以获取该系统关于参考点目标的冲激响 应二维频谱H₀(f,f_d)，将其共轭矩阵与数据矩阵逐点相乘得到 经过参考点相位补偿后回波信号的二维频谱数据矩阵

步骤四、二维非均匀快速傅里叶变换

根据公式（6）得到非均匀频域位置m_b,n_b，对经过参考点相位补偿的回波信 号二维频谱矩阵做二维非均匀快速傅里叶变换，得到二维空域图像 d(p,q)；

步骤五、空间重采样

利用公式（7）得到二维空间重采样参数a₁□a₃,b₁□b₃，并对二维空域图像d(p,q) 进行二维空间重采样g(p,q)=d(p+a₁p²+a₂q²+a₃pq,q+b₁p²+b₂q²+b₃pq)，得到成像结 果g(p,q)；

经过上述步骤处理，作坐标变换：就可以从SA-BSAR系统 接收到的目标回波数据s(τ,t)中获取具有较高分辨率的目标成像结果g(r,x)。

图3为回波数据s(τ,t)经过上述步骤一至步骤三处理后的成像结果。从图3 中可以看出除参考点目标之外的其余点目标，其成像结果位置不正确，出现了严 重的二维非线性几何畸变，因此，无法满足高分辨率成像的要求，需要步骤四至 步骤五进一步的处理。

图4是回波数据s(τ,t)经过上述步骤一至步骤五处理后的最终成像结果。其 中，横轴表示斜距向，纵轴表示方位向，坐标单位均为米，坐标原点为参考点目 标所在位置。从图4中可以看出，采用本发明提供的成像方法处理后，各点目标 均被很好的聚焦，且各点目标分别位于各自正确的空间位置且聚焦效果良好。因 此，本发明提供频域成像方法适用于斜视SA-BSAR系统，能有效校正成像结果 几何畸变、散焦等严重失真问题，可以用于其实现其高分辨率成像处理。