CN103278819B - 基于滑动接收窗的机载高分辨率斜视聚束sar成像方法 - Google Patents

Abstract

本发明公开了一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，包括以下几个步骤：步骤一：读入原始回波数据及相关成像参数；步骤二：方位向解线性调频处理；步骤三：方位向傅里叶变换处理；步骤四：方位向解线性调频残留相位误差补偿处理；步骤五：距离向傅里叶变换处理；步骤六：一致压缩处理；步骤七：斯托尔特（stolt）插值处理；步骤八：方位向傅里叶逆变换处理；步骤九：几何校正处理；步骤十：距离向傅里叶变换处理；本发明提出了一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，解决了目前基于滑动接收窗技术的机载高分辨率斜视聚束SAR原始回波数据没有成像方法的现状。

Description

基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法

技术领域

本发明属于信号处理领域，特别涉及一种基于滑动接收窗的机载高分辨率斜视聚束合成孔径雷达合成孔径雷达(Synthetic Aperture Radar，SAR)成像方法。

背景技术

SAR是一种工作在微波频段的主动遥感器，克服了光学成像受天气和光照条件限制的缺陷，能全天时、全天候、远距离的进行对地遥感观测，并能够穿透天然植被、人工伪装等，大大提高了雷达的信息捕获能力。因此，SAR已经成为雷达技术的热门研究领域，为越来越多的国家所重视。相比于传统的侧视SAR***，机载斜视SAR成像在实际应用中具有很高的灵活性和机动性，通过调整天线波束指向，SAR***可自由灵活地选择观测区域，并且能够快速重访敏感区域，大大提高了SAR的观测能力。此外，随着SAR分辨率的提高，使得SAR目标侦察与识别能力得到显著提升。因此，近些年来，高分辨率斜视成像已经成为一个重要的发展方向。

但是，高分辨率大斜视SAR成像也带来了新的技术难点。一方面，随着斜视角度的增大，SAR原始回波方位向与距离向耦合现象更加严重，导致高精度成像更加困难；另一方面，距离徙动量随斜视角的增加以及分辨率的增加而成几何级数增长，导致SAR原始回波数据量的显著增加，加大了数据存储与实时成像的困难。滑动接收窗技术是指通过改变回波窗开启时刻，消除距离徙动中的一次项，即距离走动。滑动接收窗技术能够减小SAR原始回波距离徙动量，进而减小SAR原始回波数据量。但是，滑动接收窗技术调整了回波数据的录取起始时间，观测区域内每一个目标的多普勒历程也随方位向发生了变化，使得传统成像方法不再适用。

发明内容

本发明的目的是为了解决上述问题，基于滑动接收窗技术特点，结合传统波数域成像方法，提出了一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法。

一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，包括以下几个步骤：

步骤一：读入原始回波数据及相关成像参数；

读入基于滑动接收窗的机载高分辨率斜视聚束SAR二维原始回波仿真复数据S_start以及相应的成像参数，具体包括：方位向采样点数N_a，距离向采样点数N_r，信号采样率f_s，信号带宽Bw，脉冲宽度τ，调频斜率b，脉冲重复频率PRF，参考斜距R_ref，多普勒中心频率fd₀，多普勒调频率f_r0，卫星速度P_v，等效斜视角信号波长λ，信号载频f₀，信号传播速度c；

步骤二：方位向解线性调频处理；

将二维回波仿真复数据S_start进行方位向解线性调频处理：首先复数据S_start同因子H₁相乘，得到复数据S_{1_1}；其次对复数据S_{1_1}做方位向傅里叶变换，即沿每个距离向（按列）进行快速傅里叶变化（FFT），得到复数据S_{1_2}；最后将复数据S_{1_2}同因子H₂相乘，得到最终方位向解线性调频后的复数据S₁；

步骤三：方位向傅里叶变换处理；

将步骤二得到的复数据S₁沿每个距离向（按列）进行快速傅里叶变换（FFT），得到方位向频谱复数据S₂；

步骤四：方位向解线性调频残留相位误差补偿处理；

将步骤三得到的复数据S₂同对应方位时刻的方位向解线性调频残留相位误差补偿因子Ω₁相乘，得到补偿后的复数据S₃；

步骤五：距离向傅里叶变换处理；

将步骤四得到的复数据S₃沿每个方位时刻（按行）进行快速傅里叶变换（FFT），得到二维频谱复数据S₄；

步骤六：一致压缩处理；

将步骤五得到的复数据S₄同对应的一致压缩因子Ω₂相乘，得到粗聚焦复数据S₅；

步骤七：斯托尔特（stolt）插值处理；

对于步骤六得到的复数据S₅，利用辛格（sinc）插值法进行stolt插值处理，得到由二维频域映射到二维波数域的复数据S₆；

步骤八：方位向傅里叶逆变换处理；

将步骤七得到的复数据S₆沿每个距离向（按列）进行快速傅里叶逆变换（IFFT），得到方位时域距离波数域复数据S₇；

步骤九：几何校正处理；

将步骤八得到的复数据S₇同几何校正因子Ω₃相乘，得到经几何校正后的复数据S₈；

步骤十：距离向傅里叶变换处理；

将步骤九得到的复数据S₉沿每个方位时刻（按行）进行快速傅里叶逆变换（IFFT），得到最终的成像结果S_end；

本发明优点在于：

（1）本发明提出了一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，解决了目前基于滑动接收窗技术的机载高分辨率斜视聚束SAR原始回波数据没有成像方法的现状。

（2）本发明提出了一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，具有高精度聚焦成像的特点。由于本发明提出的成像方法是一种改进的波数域成像方法，而波数域成像方法的优势是，只要满足速度恒定这一条件（机载SAR恰好满足这一条件），就能够实现高精度聚焦。因此，利用本发明能够实现场景高精度聚焦成像。

（3）本发明提出了一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，具有适用性强的特点。一方面，由于本发明提出的成像方法是一种改进的波数域成像方法，而波数域成像方法不受斜视角度的限制，因此，在斜视角度很大的条件下，本发明同样能实现场景的精确聚焦，得到高质量的SAR图像。另一方面，本发明能够实现0.1m超高分辨率成像，并得到精确聚焦，因此，本发明适用于目前各种分辨率的成像要求。

附图说明

图1是本发明提出的一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法流程图；

图2是实施例仿真场景示意图；

图3是实施例成像结果；

图4是实施例左上点目标剖面图；

图5是实施例中间点目标剖面图；

图6是实施例右下点目标剖面图。

具体实施方式

下面将结合附图和实施例对本发明作进一步的详细说明。

本发明提出一种基于滑动接收窗的机载分辨率斜视聚束SAR成像方法，处理的对象是基于滑动接收窗的机载高分辨率斜视聚束SAR原始回波数据，得到的结果是一幅高分辨率斜视图像。

滑动接收窗是指SAR在接收回波时，回波接收窗开启时刻随方位时刻而改变，从而减小目标的距离徙动，或者说消除目标的距离走动（距离徙动中的线性部分）。采用滑动接收窗的SAR***在工作中各方位时刻回波窗开启时刻T'(i)为：

$T^{\′}\left(i\right)=T_0-\frac{\mathrm{\λ}\·{\mathrm{fd}}_0\·t\left(i\right)}{c}---\left(1\right)$

其中，T₀是相同条件下传统SAR***固定回波窗开启时刻，fd₀是指多普勒中心频率，λ指信号波长，c指信号传播速度，t(i)是指方位时刻，且,i＝0,1,2,…,N_a-1。滑动接收窗技术减小了回波的距离徙动量，进而减小了回波数据量，缓解了数据存储和实时成像的压力。

本发明是一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，具体流程如图1所示，包括以下步骤：

步骤一：读入基于滑动接收窗的机载高分辨率斜视聚束SAR二维原始回波复数据S_start以及相应的成像参数。其中，S_start是一个二维复数组，大小为N_a×N_r，而成像参数具体包括：方位向采样点数N_a、距离向采样点数N_r、信号采样率f_s、信号带宽Bw、调频斜率b、脉冲重复频率PRF、参考斜距R_ref、多普勒中心频率fd₀、多普勒调频率f_r0、卫星速度P_v、等效斜视角、信号波长λ、信号载频f₀、信号传播速度c；

步骤二：将二维原始回波复数据S_start进行方位向解线性调频处理，具体可以分为以下几个步骤：

（a）构造两个一维序列i,j，其中i代表方位向序列（行），j代表距离向序列（列）；

i＝[1,2,…,N_a]   （2）

j＝[1,2,…,N_r]

（b）获取二维原始回波复数据S_start各行对应的方位时刻t(i)；

$t\left(i\right)=\frac{i-N_a/2}{\mathrm{PRF}}---\left(3\right)$

（c）将二维原始回波复数据S_start同去旋转因子H₁相乘，得到数据S_{1_1}(i,j)，其中因子H₁(i)是大小为N_a×1的一维复数组，公式为：

H₁(i)＝exp{jπ(fr₀·t²(i)+2·fd₀·t(i))}   （4）

则二维复数组S_{1_1}可由下式得出：

S_{1_1}(i,j)＝S_start(i,j)·H₁(i)   （5）

（d）对复数据S_{1_1}沿每个距离向（按列）进行快速傅里叶变化（FFT），得到复数据S_{1_2}；

S_{1_2}(:,j)＝FFT(S_{1_1}(:,j))   （6）

其中，S_{1_2}(:,j)表示S_{1_2}的第n列，S_{1_1}(:,j)表示S_{1_1}的第n列，FFT(·)表示对一维数组进行快速傅里叶变换。

（e）获取去旋转后的等效脉冲重复频率PRF'；

${\mathrm{PRF}}^{\′}=\frac{N_r\·{\mathrm{fr}}_0}{\mathrm{PRF}}---\left(7\right)$

（g）结合式（7）获取二维复数据每行对应的方位时刻t₁(i)；

$t_1\left(i\right)=\frac{i-N_a/2}{{\mathrm{PRF}}^{\′}}---\left(8\right)$

（h）将复数据S_{1_2}同因子H₂相乘，得到最终方位向解线性调频后的数据S₁。其中因子H₂(i)是大小为N_a×1的一维复数组，其公式为：

$H_2\left(i\right)=\exp \{\mathrm{j\π}\·{\mathrm{fr}}_0\·t_1^2\left(i\right)\}---\left(9\right)$

则二维复数组S₁可由下式得出：

S₁(i,j)＝S_{1_2}(i,j)·H₂(i)   （10）

步骤三：将步骤二得到的复数据S₁(i,j)沿每个距离向（按列）进行快速傅里叶变换（FFT），得到方位向频谱复数据S₂(i,j)；

S₂(:,j)＝FFT(S₁(:,j))   （11）

其中，S₂(:,j)表示S₂的第j列，S₁(:,j)表示S₁的第j列，FFT(·)表示对一维数组进行快速傅里叶变换。

步骤四：将步骤三得到的复数据S₂(i,j)同对应方位时刻的方位向解线性调频残留相位误差补偿因子Ω₁(i)相乘，得到补偿后的复数据S₃(i,j)；

（a）结合式（7）获取二维方位向频域距离向时域复数据S₂(i,j)每行对应的方位频率f_a(i)；

$f_a\left(i\right)=\frac{i-N_a/2}{N_a}\·{\mathrm{PRF}}^{\′}---\left(12\right)$

（b）结合式（12）获取大小为N_a×1的一维补偿因子Ω₁(i)；

${\mathrm{\Ω}}_1\left(i\right)=\exp \{\mathrm{j\π}\·\frac{f_a^2\left(i\right)}{{\mathrm{fr}}_0}\}---\left(13\right)$

（c）获取补偿后的二维复数据S₃(i,j)；

S₃(i,j)＝S₂(i,j)·Ω₁(i)   （14）

步骤五：将步骤四得到的复数据S₃(i,j)沿每个方位时刻（按行）进行快速傅里叶变换（FFT），得到二维频谱复数据S₄(i,j)；

S₄(i,:)＝FFT(S₃(i,:))   （15）

其中，S₃(i,:)表示S₃的第i行，S₄(i,:)表示S₄的第i行。

步骤六：将步骤五得到的复数据S₄(i,j)同对应的一致压缩因子Ω₂(i,j)相乘，得到粗聚焦复数据S₅(i,j)。

（a）根据参考斜距R_ref获取最短斜距R_min；

$R_{\min }=R_{\mathrm{ref}}-\frac{c}{{2f}_s}\·\frac{N_r}{2}---\left(16\right)$

（b）获取二维频域复数据S₄(i,j)每列对应的距离频域f_τ(j)；

$f_{\mathrm{\τ}}\left(j\right)=\frac{j-N_r/2}{\mathrm{Nr}}\·f_s---\left(17\right)$

（c）结合式（12）与式（17）获取二维频域复数据S₄(i,j)每列对应的方位向波数k_x(i)与每行对应的距离向波数k_rc(j)；

$k_x\left(i\right)=\frac{2\mathrm{\π}{f}_a\left(i\right)}{P_v}$ $\left(18\right)$

$k_{\mathrm{rc}}\left(j\right)=\frac{4\mathrm{\π}(f_0+f_{\mathrm{\τ}}\left(j\right))}{c}$

（d）结合式（16）～式（18）获取大小为N_a×N_r的二维一致压缩因子Ω₂(i,j)；

（e）结合式（19）获取一致压缩后的二维复数据S₅(i,j)；

S₅(i,j)＝S₄(i,j)·Ω₂(i,j)   （20）

步骤七：对于步骤六得到的复数据S₅(i,j)，利用sinc插值法进行stolt插值处理，得到由二维频域映射到二维波数域的复数据S₆(i,j)；

（a）根据距离频域到距离波数域的映射关系，结合式（18）获取距离波数域波数k'_rc(i,j)；

（b）遍历获取距离波数域波数k'_rc(i,j)的最大值k'_rc,max与最小值k'_rc,max，并获取距离波数域波数等分间隔Δk'_rc；

$\mathrm{\Δ}{k}_{\mathrm{rc}}^{\′}=\frac{k_{\mathrm{rc},\max }^{\′}-k_{\mathrm{rc},\min }^{\′}}{N_r}---\left(22\right)$

（c）获取二维波数域数据均匀的距离波数域波数

$k_{\mathrm{rc}}^e\left(j\right)=k_{\mathrm{rc}}^{\′}\min +j\·{\mathrm{\Δk}}_{\mathrm{rc}}^{\′}---\left(23\right)$

（d）获取二维波数域复数据各均匀的距离波数域波数在每行对应的不均匀k'_rc(i,:)中的位置p(i,j)。

方法为按行进行下列操作：以获取第一行第一列的位置p(1,1)为例，首先获取绝对差值n＝[1,2,…,N_r]，获取最小的绝对差值Δk_min和对应位置的n，若 $k_{\mathrm{rc}}^e\left(1\right)<k_{\mathrm{rc}}^{\&prime;}(1,n),$ p(1,1)＝n-1，若 $k_{\mathrm{rc}}^e\left(1\right)\&GreaterEqual;k_{\mathrm{rc}}^{\&prime;}(1,n),$ p(1,1)＝n，

以此类推，得到每一个位置p(i,j)。

（e）结合上步得到的位置p(i,j)，获取sinc插值所需采样点间隔q(i,j,n)；

$q(i,j,n)=\frac{k_{\mathrm{rc}}^e\left(j\right)-k_{\mathrm{rc}}^{\&prime;}(i,(p(i,j)+n))}{\frac{4\mathrm{\&pi;}{f}_s}{c}\&CenterDot;\frac{1}{N_r}},$ n＝[-N/2,-N/2+1,…,N/2-1]   （24）

（f）结合式（24）利用sinc插值法，获取出经stolt插值后的二维数据S₆，由于二维数据是复数据，需要分别对S₆(i,j)的实部S_{6_re}(i,j)和虚部S_{6_im}(i,j)分别进行sinc插值法获取得出。

其中，N是插值核长度，sinc(·)是指插值函数S_{5_re}(i,j)是指二维数据S_{5_re}第i行第j列的实部，S_{5_im}(i,j)是指二维数据S_{5_im}第i行第j列的虚部。

步骤八：将步骤七得到的复数据S₆(i,j)沿每个距离向（按列）进行快速傅里叶逆变换（IFFT），得到方位时域距离波数域复数据S₇(i,j)；

S₇(:,j)＝IFFT(S₆(:,j))   （27）

其中，S₆(:,j)表示S₆的第j列，S₇(:,j)表示S₇的第j列，IFFT(·)表示对一维数组进行快速傅里叶逆变换。

步骤九：将步骤八得到的复数据S₇(i,j)同几何校正因子Ω₃(i,j)相乘，得到经几何校正后的复数据S₈(i,j)；

（a）结合式（8）与式（12），获取几何校正因子Ω₄(i,j)；

${\mathrm{\&Omega;}}_4(i,j)=\exp \{-j2\mathrm{\&pi;}\&CenterDot;f_a\left(i\right)\&CenterDot;\frac{\mathrm{\&lambda;}\&CenterDot;{\mathrm{fd}}_0\&CenterDot;t_1\left(i\right)}{c}\}---\left(28\right)$

（b）利用式（30）获取经几何校正后的复数据S₈(i,j)；

S₈(i,j)＝S₇(i,j)·Ω₃(i,j)   （29）

步骤十：将步骤九得到的复数据S₈(i,j)沿每个方位时刻（按行）进行快速傅里叶逆变换（IFFT），得到最终的成像结果S_end(i,j)；

S_end(i,:)＝IFFT(S₈(i,:))   （30）

其中，S₈(i,:)表示S₈的第i行，S_end(i,:)表示S_end的第i行。

实施例：

本实施例提出一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，仿真场景为3×3点阵，其点与点的间隔为100m，最终分别对仿真场景中左上、中间、右下三个点进行了评估，具体仿真场景如图2所示，其成像过程中涉及的成像参数如表1所示。

表1实施例参数

本实施例具体包括以下步骤：

步骤一：读入基于滑动接收窗的机载高分辨率斜视聚束SAR二维原始回波复数据S_start以及相应的成像参数。其中，S_start是二维复数组，大小为65536×16384，具体成像参数如表1所示；

步骤二：将二维原始回波复数据S_start进行方位向解线性调频处理，具体操作步骤为：

（a）构造一维序列，如式（2）所示，i＝[1,2,…,65536],j＝[1,2,…,16384]；

（b）获取二维原始回波复数据S_start每行对应的方位时刻t(i)，具体过程按式（3）进行；

（c）将数据S_start同去旋转因子H₁(i)相乘，得到数据S_{1_1}(i,j)。其中因子H₁是大小为65536×1的一维复数组，具体获取过程按式（4）进行，而二维复数组S_{1_1}(i,j)获取过程按式（5）进行；

（d）对复数据S_{1_1}(i,j)按式（6）沿每个距离向进行快速傅里叶变化（FFT），得到复数据S_{1_2}(i,j)；

（e）获取去旋转后的等效脉冲重复频率PRF'，获取过程按式（7）进行；

（g）结合式（7）获取二维复数据每行对应的方位时刻t₁(i)，获取过程按式（8）进行；

（h）将数据S_{1_2}(i,j)同因子H₂(i)相乘，得到最终方位向解线性调频后的数据S₁(i,j)。其中因子H₂(i)是大小为65536×1的一维复数组，其获取过程按式（9）进行，而二维复数组S₁获取过程按式（10）进行；

步骤三：将步骤二得到的复数据S₁(i,j)沿每个距离向（按列）进行快速傅里叶变换（FFT），得到方位向频谱复数据S₂(i,j)，具体操作过程按式（11）进行；

步骤四：将步骤三得到的复数据S₂(i,j)同对应方位时刻的方位向解线性调频残留相位误差补偿因子Ω₁(i)相乘，得到补偿后的复数据S₃(i,j)；

（a）获取二维方位向频域距离向时域复数据S₂(i,j)每行对应的方位频率f_a，具体获取过程按式（12）进行；

（b）结合式（7）获取大小为65536×1的一维补偿因子Ω₁(i)，具体获取过程按式（13）进行；

（c）获取补偿后的二维复数据S₃(i,j)，具体获取过程按式（14）进行；

步骤五：将步骤四得到的复数据S₃(i,j)沿每个方位时刻（按行）进行快速傅里叶变换（FFT），得到二维频谱复数据S₄(i,j)，具体操作过程按式（15）进行；

步骤六：将步骤五得到的复数据S₄(i,j)同对应的一致压缩因子Ω₂(i,j)相乘，得到粗聚焦复数据S₅(i,j)。

（a）根据参考斜距R_ref＝19.31km获取最短斜距R_min，具体获取过程按式（16）进行；

（b）获取二维频域复数据S₄(i,j)每列对应的距离频域f_τ(j)，具体获取过程按式（17）进行；

（c）结合式（7）与式（11）获取二维频域复数据S₄(i,j)每列对应的方位向波数k_x(i)与每行对应的距离向波数k_rc(j)，具体获取过程按式（18）进行；

（d）结合式（10）～式（12）获取大小为65536×16384的二维一致压缩因子Ω₂(i,j)，具体获取过程按式（19）进行；

（e）结合式（13）获取一致压缩后的二维复数据S₅(i,j)，具体获取过程按式（20）进行；

步骤七：对于步骤六得到的复数据S₅(i,j)，利用sinc插值法进行stolt插值处理，得到由二维频域映射到二维波数域的复数据S₆(i,j)；

（a）根据距离频域到距离波数域的映射关系，结合式（18）获取距离波数域波数k'_rc(i,j)，具体获取过程按式（21）进行；

（b）遍历获取距离波数域波数k'_rc(i,j)的最大值k'_rc,max＝69.12rad/s与最小值k'_rc,max＝-105.06rad/s，并获取距离波数域波数等分间隔Δk'_rc，具体操作按式（22）进行；

（c）获取二维波数域复数据均匀的距离波数域波数具体操作按式（23）进行；

（d）获取二维波数域复数据各均匀的距离波数域波数在每行对应的不均匀k'_rc(i,:)中的位置p(i,j)。方法为按行进行下列操作：先获取第一行第一列的位置p(1,1)，获取绝对差值n＝[1,2,…,16384]，获取最小的绝对差值Δk_min和对应位置的n，若 $k_{\mathrm{rc}}^e\left(1\right)<k_{\mathrm{rc}}^{\&prime;}(1,n),$ p(1,1)＝n-1，若 $k_{\mathrm{rc}}^e\left(1\right)\&GreaterEqual;k_{\mathrm{rc}}^{\&prime;}(1,n),$ p(1,1)＝n，

以此类推，得到每一个位置p(i,j)。

（e）结合上步得到的位置p(i,j)，获取sinc插值所需采样点间隔q(i,j,n)，具体操作按式（24）进行；

（f）结合式（18）利用sinc插值法，选择sinc插值核长度为N＝8，获取出经stolt插值后的二维复数据S₆(i,j)，由于二维数据是复数数据，需要分别对S₆(i,j)的实部S_{6_re}(i,j)和虚部S_{6_im}(i,j)分别进行sinc插值法获取得出，具体操作按式（25）（26）进行。

步骤八：将步骤七得到的复数据S₆(i,j)沿每个距离向（按列）进行快速傅里叶逆变换（IFFT），得到方位时域距离波数域复数据S₇(i,j)，具体操作按式（27）进行；

步骤九：将步骤八得到的复数据S₇(i,j)同几何校正因子Ω₃(i,j)相乘，得到经几何校正后的复数据S₈(i,j)；

（a）结合式（8）与式（12），获取几何校正因子Ω₃(i,j)，具体操作按式（28）进行；

（b）利用式（28）获取经几何校正后的复数据S₈(i,j)，具体操作按式（29）进行；

步骤十：将步骤九得到的复数据S₈(i,j)沿每个方位时刻（按行）进行快速傅里叶逆变换（IFFT），得到最终的成像结果S_end(i,j)，具体操作按式（30）进行；

经过上述步骤的成像处理，得到最终的场景成像结果如图3所示。表2给出了场景左上、中间、右下三个点目标的成像评估结果，图4、图5、图6分别给出了场景左上、中间、右下三个点目标的剖面图。

表二成像评估结果

根据表2评估结果以及图4～图6所示剖面图，可以得出：一方面，该成像方法在斜视角度为70度时仍然能够精确聚焦，说明本发明提出的方法不受斜视角度的限制；另一方面，该成像方法对于0.1m超高分辨仍然能够精确聚焦，说明本发明提出的方法能够对目前各种分辨率实现精确聚焦。因此，本发明所提出的方法可以实现基于滑动接收窗的机载高分辨率斜视角聚束SAR精确成像，得到了高精度的成像结果。

Claims (1)

1.一种基于滑动接收窗的机载高分辨率斜视聚束SAR成像方法，包括以下步骤：

步骤一：读入基于滑动接收窗的机载高分辨率斜视聚束SAR二维原始回波复数据S_start以及相应的成像参数；

S_start是一个二维复数组，大小为N_a×N_r，成像参数包括：方位向采样点数N_a、距离向采样点数N_r、信号采样率f_s、信号带宽Bw、调频斜率b、脉冲重复频率PRF、参考斜距R_ref、多普勒中心频率fd₀、多普勒调频率f_r0、卫星速度P_v、等效斜视角、信号波长λ、信号载频f₀、信号传播速度c；

步骤二：将二维原始回波复数据S_start进行方位向解线性调频处理，具体包括以下几个步骤：

(a)构造两个一维序列i,j，其中i代表方位向序列，j代表距离向序列；

i＝[1,2,…,N_a]

                            (2)

j＝[1,2,…,N_r]

(b)获取二维原始回波复数据S_start各行对应的方位时刻t(i)；

$t\left(i\right)=\frac{i-N_a/2}{\mathrm{PRF}}---\left(3\right)$

(c)将二维原始回波复数据S_start与去旋转因子H₁相乘，得到数据S_{1_1}(i,j)，其中因子H₁(i)是大小为N_a×1的一维复数组，公式为：

H₁(i)＝exp{jπ(fr₀·t²(i)+2·fd₀·t(i))}               (4)

则二维复数组S_{1_1}由下式得出：

S_{1_1}(i,j)＝S_start(i,j)·H₁(i)                   (5)

(d)对复数据S_{1_1}沿每个距离向进行快速傅里叶变化，得到复数据S_{1_2}；

S_{1_2}(:,j)＝FFT(S_{1_1}(:,j))                       (6)

其中，S_{1_2}(:,j)表示S_{1_2}的第j列，S_{1_1}(:,j)表示S_{1_1}的第j列，FFT(·)表示对一维数组进行快速傅里叶变换；

(e)获取去旋转后的等效脉冲重复频率PRF'；

$\mathrm{PR}{F}^{\&prime;}=\frac{N_r\&CenterDot;{\mathrm{fr}}_0}{\mathrm{PRF}}---\left(7\right)$

(g)结合式(7)获取二维复数据每行对应的方位时刻t₁(i)；

$t_1\left(i\right)=\frac{i-N_a/2}{{\mathrm{PRF}}^{\&prime;}}---\left(8\right)$

(h)将复数据S_{1_2}同因子H₂相乘，得到最终方位向解线性调频后的数据S₁；其中因子H₂(i)是大小为N_a×1的一维复数组，其公式为：

$H_2\left(i\right)=\exp \{\mathrm{j\&pi;}\&CenterDot;{\mathrm{fr}}_0\&CenterDot;t_1^2\left(i\right)\}---\left(9\right)$

则二维复数组S₁由下式得出：

S₁(i,j)＝S_{1_2}(i,j)·H₂(i)                   (10)

步骤三：将步骤二得到的复数据S₁(i,j)沿每个距离向进行快速傅里叶变换，得到方位向频谱复数据S₂(i,j)；

S₂(:,j)＝FFT(S₁(:,j))                        (11)

其中，S₂(:,j)表示S₂的第j列，S₁(:,j)表示S₁的第j列，FFT(·)表示对一维数组进行快速傅里叶变换；

步骤四：将步骤三得到的复数据S₂(i,j)同对应方位时刻的方位向解线性调频残留相位误差补偿因子Ω₁(i)相乘，得到补偿后的复数据S₃(i,j)；

(a)结合式(7)获取二维方位向频域距离向时域复数据S₂(i,j)每行对应的方位频率f_a(i)；

$f_a\left(i\right)=\frac{i-N_a/2}{N_a}\&CenterDot;{\mathrm{PRF}}^{\&prime;}---\left(12\right)$

(b)结合式(12)获取大小为N_a×1的一维补偿因子Ω₁(i)；

${\mathrm{\&Omega;}}_1\left(i\right)=\exp \{\mathrm{j\&pi;}\&CenterDot;\frac{f_a^2\left(i\right)}{{\mathrm{fr}}_0}\}---\left(13\right)$

(c)获取补偿后的二维复数据S₃(i,j)；

S₃(i,j)＝S₂(i,j)·Ω₁(i)                       (14)

步骤五：将步骤四得到的复数据S₃(i,j)沿每个方位时刻进行快速傅里叶变换，得到二维频谱复数据S₄(i,j)；

S₄(i,:)＝FFT(S₃(i,:))                       (15)

其中，S₃(i,:)表示S₃的第i行，S₄(i,:)表示S₄的第i行；

步骤六：将步骤五得到的复数据S₄(i,j)同对应的一致压缩因子Ω₂(i,j)相乘，得到粗聚焦复数据S₅(i,j)；

(a)根据参考斜距R_ref获取最短斜距R_min；

$R_{\min }=R_{\mathrm{ref}}-\frac{c}{2f_s}\&CenterDot;\frac{N_r}{2}---\left(16\right)$

(b)获取二维频域复数据S₄(i,j)每列对应的距离频域f_τ(j)；

$f_{\mathrm{\&tau;}}\left(j\right)=\frac{j-N_r/2}{\mathrm{Nr}}\&CenterDot;f_s---\left(17\right)$

(c)结合式(12)与式(17)获取二维频域复数据S₄(i,j)每列对应的方位向波数k_x(i)与每行对应的距离向波数k_rc(j)；

$\begin{array}{c}{k}_x\left(i\right)=\frac{2\mathrm{\&pi;}{f}_a\left(i\right)}{P_v}\\ k_{\mathrm{rc}}\left(j\right)=\frac{4\mathrm{\&pi;}(f_0+f_{\mathrm{\&tau;}}\left(j\right))}{c}\end{array}---\left(18\right)$

(d)结合式(16)～式(18)获取大小为N_a×N_r的二维一致压缩因子Ω₂(i,j)；

(e)结合式(19)获取一致压缩后的二维复数据S₅(i,j)；

S₅(i,j)＝S₄(i,j)·Ω₂(i,j)                       (20)

步骤七：对于步骤六得到的复数据S₅(i,j)，利用sinc插值法进行stolt插值处理，得到由二维频域映射到二维波数域的复数据S₆(i,j)；

(a)根据距离频域到距离波数域的映射关系，结合式(18)获取距离波数域波数k'_rc(i,j)；

(b)遍历获取距离波数域波数k'_rc(i,j)的最大值k'_rc,max与最小值k'_rc,max，并获取距离波数域波数等分间隔Δk'_rc；

${\mathrm{\&Delta;k}}_{\mathrm{rc}}^{\&prime;}=\frac{k_{\mathrm{rc},\max }^{\&prime;}-k_{\mathrm{rc},\min }^{\&prime;}}{N_r}---\left(22\right)$

(c)获取二维波数域数据均匀的距离波数域波数

$k_{\mathrm{rc}}^e\left(j\right)=k_{\mathrm{rc}}^{\&prime;}\min +j\&CenterDot;{\mathrm{\&Delta;k}}_{\mathrm{rc}}^{\&prime;}---\left(23\right)$

(d)获取二维波数域复数据各均匀的距离波数域波数在每行对应的不均匀k'_rc(i,:)中的位置p(i,j)；

具体为：以获取第一行第一列的位置p(1,1)为例，首先获取绝对差值n＝[1,2,…,N_r]，获取最小的绝对差值Δk_min和对应位置的n，若p(1,1)＝n-1，若p(1,1)＝n，同理，以此类推，得到每一个位置p(i,j)；

(e)结合上步得到的位置p(i,j)，获取sinc插值所需采样点间隔q(i,j,n)；

$q(i,j,n)=\frac{k_{\mathrm{rc}}^e\left(j\right)-k_{\mathrm{rc}}^{\&prime;}(i,(p(i,j)+n))}{\frac{4\mathrm{\&pi;}{f}_s}{c}\&CenterDot;\frac{1}{N_r}},n=[-N/2,-N/2+1,...,N/2-1]---\left(24\right)$

(f)结合式(24)利用sinc插值法，获取出经stolt插值后的二维数据S₆，由于二维数据是复数据，需要分别对S₆(i,j)的实部S_{6_re}(i,j)和虚部S_{6_im}(i,j)分别进行sinc插值法获取得出；

其中，N是插值核长度，sinc(·)是指插值函数S_{5_re}(i,j)是指二维数据S_{5_re}第i行第j列的实部，S_{5_im}(i,j)是指二维数据S_{5_im}第i行第j列的虚部；

步骤八：将步骤七得到的复数据S₆(i,j)沿每个距离向进行快速傅里叶逆变换，得到方位时域距离波数域复数据S₇(i,j)；

S₇(:,j)＝IFFT(S₆(:,j))                      (27)

其中，S₆(:,j)表示S₆的第j列，S₇(:,j)表示S₇的第j列，IFFT(·)表示对一维数组进行快速傅里叶逆变换；

步骤九：将步骤八得到的复数据S₇(i,j)同几何校正因子Ω₃(i,j)相乘，得到经几何校正后的复数据S₈(i,j)；

(a)结合式(8)与式(12)，获取几何校正因子Ω₄(i,j)；

${\mathrm{\&Omega;}}_4(i,j)=\exp \{-j2\mathrm{\&pi;}\&CenterDot;f_a\left(i\right)\&CenterDot;\frac{\mathrm{\&lambda;}\&CenterDot;{\mathrm{fd}}_0\&CenterDot;t_1\left(i\right)}{c}\}---\left(28\right)$

(b)利用式(30)获取经几何校正后的复数据S₈(i,j)；

S₈(i,j)＝S₇(i,j)·Ω₃(i,j)                     (29)

步骤十：将步骤九得到的复数据S₈(i,j)沿每个方位时刻进行快速傅里叶逆变换，得到最终的成像结果S_end(i,j)；

S_end(i,:)＝IFFT(S₈(i,:))                     (30)

其中，S₈(i,:)表示S₈的第i行，S_end(i,:)表示S_end的第i行。