CN113406635A - 一种大斜视时变参数sar的多层重叠子孔径成像方法 - Google Patents

Abstract

本发明公开了一种大斜视时变参数SAR的多层重叠子孔径成像方法，采用时变参数体质调整雷达参数，提升了空间分辨性能，省略原本的距离向插值操作，简化了成像处理流程，提高了运算效率。本发明的成像方法通过将沿方位向划分成N层重叠子孔径，提出多层重叠子孔径成像算法，提高了成像幅宽；在此基础上，根据畸变坐标的估计值求解出真实坐标的估计值，用真实坐标的估计值对空变的二次相位误差进行补偿，解决了距离场景中心较远的点不能良好聚焦的问题，进一步提升了成像幅宽。

Description

一种大斜视时变参数SAR的多层重叠子孔径成像方法

技术领域

本发明涉及合成孔径雷达技术领域，具体涉及一种大斜视时变参数SAR的多层重叠子孔径成像方法。

背景技术

合成孔径雷达(Synthetic Aperture Radar,SAR)是一种主动式微波成像雷达，具有全天时、全天候地表信息获取能力，在军事侦察、资源勘探、地形测绘等方面发挥着重要作用，是当前对地观测领域不可或缺的技术手段。大斜视聚束SAR是一种重要的工作模式，其波束视线方向与航迹方向夹角的余角(斜视角)可达八十度甚至八十度以上，可以提前感知场景信息，能显著提高雷达的侦察能力。然而，大斜视SAR成像面临图像分辨差、成像精度低的问题。大斜视构型下，传统的恒定雷达参数体制导致SAR图像中目标点的多层函数(Point Spread Function,PSF)严重畸变，空间分辨恶化。斜视角越大，距离徙动越大，回波时空变耦合严重，且由于二维耦合规律具有时空依赖性，成像处理算法需具备强时变回波解耦合能力。传统的时域成像处理方法具备上述能力，但在大斜视构型下，运算量大，效率太低，不满足SAR实时性需求；传统的距离多普勒域算法处理精度随着斜视角增大而下降，直至不能成像，存在先天缺陷；传统的波数域算法，以极坐标格式算法(Polar FormatAlgorithm,PFA)为代表，通过波数域二维插值实现二维解耦合，但受限于波前平面假设，距离场景中心越远，近似误差越大，成像幅宽受限严重，成像幅宽变小。因此，开展大斜视下的SAR成像体制和算法研究具有相当的必要性。

发明内容

有鉴于此，本发明提供了一种大斜视时变参数SAR的多层重叠子孔径成像方法，能够校正畸变目标点扩展函数，提升成像性能，省略原本的距离向插值操作，简化了成像处理流程，提高了运算效率；并且可以更高精度补偿相位，提升成像幅宽。

本发明具体方案如下：

一种大斜视时变参数SAR的多层重叠子孔径成像方法，包括：

步骤一、采用时变参数体制发射信号，并接收回波信号，进行脉冲压缩处理；

步骤二、对脉压处理后获得的波数域信号，仅进行方位向插值；

步骤三、对步骤二处理后的信号进行距离向压缩，并沿方位向划分成N层重叠子孔径，其中N层重叠子孔径对应N+1维方位向信号；

步骤四、对前N维方位向信号逐个维度进行快速傅立叶变换，得到每一维的畸变坐标的估计值；根据畸变坐标的估计值求解出真实坐标的估计值，用真实坐标的估计值对空变的二次相位误差进行补偿；对第N+1维方位向信号进行快速傅立叶变换，得到第N+1维畸变坐标的估计值；

步骤五、对步骤四处理后的N+1维方位向信号进行重构，得到SAR图像。

进一步地，步骤三中沿方位向划分成N层重叠子孔径，具体划分过程为：

根据单个子孔径长度M₁，以Δ₂为数据抽取比例，抽取M₂个子孔径，构成第一层重叠子孔径；以Δ₂Δ₃为数据抽取比例，抽取M₃个第一层重叠子孔径，构成第二层重叠子孔径；以Δ₂Δ₃Δ₄为数据抽取比例，抽取M₄个第二层重叠子孔径，构成第三层重叠子孔径……以此类推，直至获取N层重叠子孔径；

其中，M₁……M_N为正整数，且M₁……M_N的取值保证每层重叠子孔径的二次相位小于π/4。

进一步地，步骤三中沿方位向划分成N层重叠子孔径，之后获得N+1维方位向信号，N+1维方位向信号表示为：

其中，i′是距离向快速傅立叶变换后的距离向序列，m₁表示-M₁/2,…,M₁/2的序列，m₂表示-M₂/2,…,M₂/2的序列，以此类推；x_p和y_p是成像目标的坐标，U₁和U₂分别为波前假设引起的距离向和方位向几何畸变量，U₃是波前假设引起的方位向二次相位系数，它们只由场景定值参数以及场景目标点坐标x_p,y_p决定，ΔK_y表示方位插值后的距离数波数带宽，Δ₂表示第一个数据抽取比例，Δ₂Δ₃表示第二个数据抽取比例，……Δ₂Δ₃…Δ_N+1表示第N个数据抽取比例；j是虚数符号，u_x是方位向波数间隔。

进一步地，步骤四中根据畸变坐标的估计值求解出真实坐标的估计值为：

根据更高阶的信号估计模型，以解析解的形式求解出真实坐标的估计值，更高阶的信号估计模型表示为：

其中，i和m分别表示快时间序列和慢时间序列，j是虚数符号，x_p和y_p是成像目标的坐标，U₁和U₂分别为波前假设引起的距离向和方位向几何畸变量，U₃是波前假设引起的方位向二次相位系数，K_x(i,m)为方位向波数，K_y(i)距离向波数。

解析解表示为：

其中

其中，n表示第n维方位向信号，n∈[1，2，......N+1]；

Y_c表示合成孔径中心时刻雷达在地面的投影到场景中心的距离，β₀表示合成孔径中心时刻的下视角，φ表示雷达速度的俯冲角，δ表示地面斜视角的余角，H表示雷达合成孔径中心时刻平台飞行高度，y_Q为距离向畸变坐标的估计值，x_Q为方位向畸变坐标的估计值。

进一步地，时变参数体制具体为：沿方位向实时调整雷达参数载频和调频率的方法来抵消由下视角、斜距地面偏角的变化引入的距离向波数的变化；

载频和调频率表示为：

其中c为光速，f_c0和γ₀分别表示初始发射信号载频和线性调频LFM信号调频率，β₀表示合成孔径中心时刻的下视角，β(t)表示平台到场景中心的瞬时下视角，α(t)表示瞬时参考斜距在地平面的投影与OY即Y轴所形成的瞬时夹角。

进一步地，脉冲压缩处理，包括解线频调和去除剩余视频相位操作。

有益效果：

(1)一种大斜视时变参数SAR的多层重叠子孔径成像方法，通过将沿方位向划分成N层重叠子孔径，提出多层重叠子孔径成像算法，提高了成像幅宽；在此基础上，用真实坐标的估计值对空变的二次相位误差进行补偿，以更高精度补偿相位，解决了距离场景中心较远的点不能良好聚焦的问题，进一步提升了成像幅宽。

(2)在获得畸变坐标的估计值之后，根据更高阶的信号估计模型，以解析解的形式求解出真实坐标的估计值，使得补偿相位精度更高，更进一步提升了成像幅宽。

(3)采用时变参数体质调整雷达参数，提升空间分辨性能，省略原本的距离向插值操作，简化了成像处理流程，提高了运算效率。

附图说明

图1为本发明所述的大斜视时变参数SAR的多层重叠子孔径成像方法流程图。

图2为大斜视SAR的几何构型示意图。

图3为方位向N层重叠子孔径划分示意图。

图4为前N维方位向信号的子孔径补偿处理流程图。

图5为实施例中验证时变参数体制优越性的仿真结果示意图。

图6为实施例中PFA、一层OSA和二层MOSA成像幅宽示意图。

图7为实施例中验证MOSA优越性的仿真结果示意图。

图8为大斜视SAR关键参数列表。

具体实施方式

下面结合附图并举实施例，对本发明进行详细描述。

本发明是一种大斜视时变参数SAR的多层重叠子孔径成像方法，如图1所示，其具体步骤包括：

步骤一、采用时变参数体制发射信号，并接收回波信号，进行脉冲压缩处理；进行脉冲压缩处理包括解线频调Dechirp和去除剩余视频相位RVP操作。

大斜视SAR的几何构型如图2所示。在图中，以场景中心点为原点O建立坐标系XYZ，XOY表示地平面，OZ即Z轴的指向背离地心且垂直于地表，OY即Y轴的指向则沿着合成孔径中心时刻斜距在地平面的投影方向。点P为场景内任一散射点，其坐标为(xp,yp,0)，A点和S点分别表示平台在合成孔径中心时刻和任一方位向采样时刻的位置。H表示雷达合成孔径中心时刻平台飞行高度，t是方位时间(慢时间)，Rc(t)表示平台到场景中心的瞬时参考斜距，Rp(t)表示平台到散射点P的瞬时斜距。α(t)表示瞬时参考斜距Rc(t)在地平面的投影与OY所形成的瞬时夹角。β(t)表示平台到场景中心的瞬时下视角，β₀表示合成孔径中心时刻的下视角。φ表示雷达速度的俯冲角，δ表示地面斜视角的余角。

假设线性调频信号为发射信号，但与传统的固定雷达平台参数的信号不同，在时变参数体制下，雷达平台参数(载频和调频率)与方位向采样时刻与位置有关，它们表示为

其中，t表示慢时间，f_c0和γ₀分别表示初始发射信号载频和LFM信号调频率，f_c和γ分别表示随着慢时间变化的载频和调频率，α(t)表示瞬时参考斜距Rc(t)在地平面的投影与OY所形成的瞬时夹角。β(t)表示平台到场景中心的瞬时下视角，β₀表示合成孔径中心时刻的下视角。

则发射信号可表示为

其中，j是虚数符号，τ表示距离时间(快时间)，Tp表示脉冲宽度，rect(·)表示窗函数，f_c和γ分别表示随着慢时间变化的载频和调频率。

将信号解调到基带，任意目标点P的回波信号可表示为

其中，c为光速，Rp(t)表示平台到散射点P的瞬时斜距。

接下来对信号进行Dechirp和去除RVP操作，将处理后的信号从模拟域映射到数字域上，可以表示为

其中i和m分别表示快时间序列和慢时间序列(s₁(τ,t)表示模拟信号，s₂(i,m)表示数字信号)，Nr和Na分别为距离向与方位向的采样点数，R_c(m)和R_p(m)分别表示R_c(t)和R_p(t)的数字域表示，F_s是采样率。

根据波前平面假设和如图2所示的几何构型，将式(4)投影到波数域，并在方位向波数K_x和距离向波数K_y进行二维泰勒展开，获得更高阶的信号模型，表示为：

其中，x_p和y_p是成像目标的坐标，U₁和U₂分别为波前假设引起的距离向和方位向几何畸变量，U₃是波前假设引起的方位向二次相位系数，它们只由场景定值参数以及场景目标点坐标x_p,y_p决定；i和m分别表示快时间序列和慢时间序列。方位向波数K_x(i,m)和距离向波数K_y(i)如下式(6)，由于采用时变参数体制，距离向波数只随i变化

其中，α(m)是α(t)的数字域表示。

步骤二、对脉压处理后获得的波数域信号，仅进行方位向插值。

因为采用时变参数体制，公式(5)表示的回波信号在波数域呈梯形分布，无需距离向插值校正，可直接进行方位向插值处理，其通过如下映射实现

方位向插值后的信号表示为

其中，K_x′(m)表示方位向插值后的方位向波数，如下式所示：

其中，u_x是方位向波数间隔。

步骤三、对步骤二处理后的信号进行距离向压缩，并沿方位向划分成N层重叠子孔径，其中N层重叠子孔径对应N+1维方位向信号，具体为第N层重叠子孔径对应第N+1维方位向信号。

沿方位向划分成N层重叠子孔径，具体划分过程为：

根据单个子孔径长度M₁，以Δ₂为数据抽取比例，抽取M₂个子孔径，构成第一层重叠子孔径；以Δ₂Δ₃为数据抽取比例，抽取M₃个第一层重叠子孔径，构成第二层重叠子孔径；以Δ₂Δ₃Δ₄为数据抽取比例，抽取M₄个第二层重叠子孔径，构成第三层重叠子孔径……以此类推，直至获取N层重叠子孔径。

其中，M₁……M_N为正整数，且M₁……M_N的取值应该保证每层重叠子孔径的二次相位小于π/4。如式(25)所示，第一层子孔径的个数M₂，第三个相位里的二次相位

被限制在π/4以内。

对插值后的信号，如式(8)所示，进行距离向压缩即距离向FFT，可以得到

其中，i′是距离向FFT后的距离向序列，ΔK_y表示方位插值后的距离数波数带宽，并得到距离向畸变坐标的估计值

之后，沿着方位向划分成N层重叠子孔径，如图3所示。M₁表示第一层子孔径的长度，Δ₂表示第一个数据抽取比例，M₂表示第一层子孔径的个数，Δ₂Δ₃表示第二个数据抽取比例，以此类推。m₁表示-M₁/2,…,M₁/2的序列，m₂表示-M₂/2,…,M₂/2的序列，以此类推。方位向序列m可以写为

m＝m₁+Δ₂m₂+Δ₂Δ₃m₃+…+Δ₂…Δ_N+1m_N+1 (12)

将式(12)代入式(10)，可得沿方位向划分为N层重叠子孔径的信号，表示为：

步骤四、对前N维方位向信号逐个维度进行快速傅立叶变换，得到每一维畸变坐标的估计值；根据更高阶的信号估计模型，以解析解的形式求解出真实坐标的估计值，对空变的二次相位误差进行补偿；对第N+1维方位向信号进行快速傅立叶变换，得到第N+1维畸变坐标的估计值。

对前N维方位向信号进行子孔径补偿处理，即利用估计值补偿空变的二次相位误差(Quadratic Phase Error,QPE)。

步骤四的流程图如图4所示。首先，为确保m₁维度信号可以被良好聚焦，选取合适的第一层子孔径的长度M₁，使得式(13)第二个相位里的

被限制在π/4以内。对第一维方位向(m₁方向)进行FFT后，第一维方位向FFT后的信号表示为：

其中，m′₁表示m₁方向FFT后的序列，由此得出方位向畸变坐标的第一层的第一次估计即第一个粗估计为：

其中ρ_x1为方位向第一粗分辨率。根据畸变坐标

到真实坐标

的映射，具体推导见公式(17)到公式(24)，得出解析解形式，由式(11)和(15)，求解出距离向真实坐标的估计值和方位向真实坐标的第一个粗估计，如下所示：

其中，f₁和f₂为求解真实坐标的解析解表达式，其具体推导如下：

首先

和

的具体表达式如下

其中，

Y_c表示合成孔径中心时刻雷达在地面的投影到场景中心的距离，

由式(11)和(15)，得到下式

将式(17)里

的表达式代入式(18)里第一个式子，得

R_p0＝R_c0-y_Qsinβ₀ (19)

将上式(19)代入

的表达式里，式(18)变为

因此，

可以表示为

根据R_p0的定义和式(19)，得到如下等式

然后将式(21)代入式(22)，得到下面的方程

然后，通过式(23)可以求出

如下式(24)，之后，

可以通过式(21)求出。至此，由畸变坐标

到真实坐标

的映射推导完毕。

至此，得到解析解的表达式：

其中

其中，

Y_c表示合成孔径中心时刻雷达在地面的投影到场景中心的距离，β₀表示合成孔径中心时刻的下视角，φ表示雷达速度的俯冲角，δ表示地面斜视角的余角，H表示雷达合成孔径中心时刻平台飞行高度，y_Q为距离向畸变坐标的估计值，x_Q为方位向畸变坐标的估计值。

对于不同维度的解析解表达式为：

其中，n表示第n维方位向信号，n小于等于N+1(每一维方位向信号解析解的表达式的形式是相同的，方位向畸变坐标的估计值会随着维数的变化而变化)；

根据式(16)求出

进行相位补偿，补偿后的信号写为

为确保m₂维度信号可以被良好聚焦，选取合适的第一层子孔径的个数M₂，使得式(25)第三个相位里的二次相位

被限制在π/4以内。对第二维方位向(m₂方向)进行FFT后，第二维方位向FFT后的信号表示为：

其中，m′₂表示m₂方向FFT后的序列，得出方位向畸变坐标的第二个粗估计为：

其中，ρ_x2为方位向第二粗分辨率，由式(11)和(27)，得出距离向真实坐标的估计值和方位向真实坐标的第二个粗估计，如下所示

其中，

是距离向真实坐标的估计值(也可以称为距离向真实坐标的第二次估计值)，

是方位向真实坐标的第二个粗估计，f₃和f₄为求解真实坐标的解析解表达式，进而求出

进行相位补偿，补偿后的信号写为

同理，对进行m₃…m_N维度进行子孔径处理，进行N次补偿后的信号写为

步骤五、对处理后的N+1维方位向信号进行重构，得到SAR图像，即重构完的每一距离门的一维方位向数据，一起组成SAR二维图像。

对第N+1维方位向信号进行快速傅立叶变换(Fast Fourier Transform,FFT)：

在式(30)的信号进行第N+1维方位向(m_N+1方向)进行FFT，得到

此时，所需要的一维方位向数据存储于N+1维矩阵中，如式(31)所示。因此，将FFT后的数据矢量化，重构出一维方位向数据，定义新的输出方位向变量为m′：

其中ρ₁…ρ_xN+1为方位向第一到第N+1个粗分辨率。

根据重构完的一维方位向数据方位向信号，可以获得SAR成像图像。

本发明提出了时变参数SAR的多层重叠子孔径成像方法。在雷达体制层面，针对恒定参数体制目标PSF畸变大的问题，采用时变调整雷达参数(简称时变参数)校正目标PSF，在大斜视构型下实现二维正交PSF，提升空间分辨性能。在成像处理算法层面，在波数域算法PFA的基础上，为解决距离场景中心较远的点不能良好聚焦的问题，本发明采用一种更高阶的信号估计模型，结合重叠子孔径算法(Overlapped Sub-Aperture Algorithm,OSA)，补偿空变的相位误差，并将传统的方位向一层子孔径OSA的多层到N层，提出多层重叠子孔径成像算法(Multi-layer Overlapped Sub-Aperture Algorithm,MOSA)，进一步提升成像幅宽。

实施例

为验证时变参数SAR体制的优越性，使用图(8)中的参数，对比常规后向投影算法(Backprojection Algorithm,BPA)和变参数体制BPA下的点目标成像性能，仿真结果如图5所示。图5中(a)和(b)分别是载频和调频率随时间变化的曲线；(c)和(d)分别是对原点进行仿真的结果图，可以看出，恒定参数下，方位向和距离向不满足正交性，而时变参数体制下，方位向和距离向满足正交性，并且变参数的最差分辨率提升了28％；(e)和(f)分别表示两个相邻点目标的成像结果，可以看出常规BP的成像结果很难分辨两个点目标，但是变参数BP的成像结果可以分辨出两个点目标，证明了时变参数体制的优越性，成像分辨性能更好(提升)。

为验证MOSA在成像幅宽方面的优越性，以方位向二层重叠子孔径的MOSA为例，与一层OSA算法进行对比，仿真结果如图6和图7所示。在图6中，(a)，(b)和(c)的黑色实线以内的区域分别表示PFA、一层OSA和二层MOSA的成像幅宽，由图可知，二层MOSA的成像幅宽与一层OSA相比有明显提升，并且都明显优于PFA。选取在一层OSA和二层MOSA成像幅宽以内的点A(300,-100)，和在一层OSA成像幅宽以外，但在二层MOSA聚集深度以内的点B(290,-400)，进行点目标仿真分析，仿真结果如图7所示。在图7中，(a)，(b)和(c)为两种算法在点A的成像结果，因为在两种算法的成像幅宽以内，成像结果均良好，方位向峰值旁瓣比均能达到-13dB；(e)，(f)和(g)为两种算法在点B的成像结果，二层MOSA的方位向峰值旁瓣比为-13.2dB，而一层OSA的方位向峰值旁瓣比为-10.7dB，二层MOSA在方位向的成像性能有明显提升，证明了MOSA算法的优越性。

以上的具体实施例仅描述了本发明的设计原理，该描述中的部件形状，名称可以不同，不受限制。所以，本发明领域的技术人员可以对前述实施例记载的技术方案进行修改或等同替换；而这些修改和替换未脱离本发明创造宗旨和技术方案，均应属于本发明的保护范围。

Claims (6)

1.一种大斜视时变参数SAR的多层重叠子孔径成像方法，其特征在于，包括：

步骤一、采用时变参数体制发射信号，并接收回波信号，进行脉冲压缩处理；

步骤二、对脉压处理后获得的波数域信号，仅进行方位向插值；

步骤三、对步骤二处理后的信号进行距离向压缩，并沿方位向划分成N层重叠子孔径，其中N层重叠子孔径对应N+1维方位向信号；

步骤四、对前N维方位向信号逐个维度进行快速傅立叶变换，得到每一维的畸变坐标的估计值；根据畸变坐标的估计值求解出真实坐标的估计值，用真实坐标的估计值对空变的二次相位误差进行补偿；对第N+1维方位向信号进行快速傅立叶变换，得到第N+1维畸变坐标的估计值；

步骤五、对步骤四处理后的N+1维方位向信号进行重构，得到SAR图像。

2.如权利要求1所述的成像方法，其特征在于，步骤三中所述沿方位向划分成N层重叠子孔径，具体划分过程为：

根据单个子孔径长度M₁，以Δ₂为数据抽取比例，抽取M₂个子孔径，构成第一层重叠子孔径；以Δ₂Δ₃为数据抽取比例，抽取M₃个第一层重叠子孔径，构成第二层重叠子孔径；以Δ₂Δ₃Δ₄为数据抽取比例，抽取M₄个第二层重叠子孔径，构成第三层重叠子孔径……以此类推，直至获取N层重叠子孔径；

其中，M₁……M_N为正整数，且M₁……M_N的取值保证每层重叠子孔径的二次相位小于π/4。

3.如权利要求1所述的成像方法，其特征在于，所述步骤三中沿方位向划分成N层重叠子孔径，之后获得N+1维方位向信号，所述N+1维方位向信号表示为：

其中，i′是距离向快速傅立叶变换后的距离向序列，m₁表示-M₁/2,…,M₁/2的序列，m₂表示-M₂/2,…,M₂/2的序列，以此类推；x_p和y_p是成像目标的坐标，U₁和U₂分别为波前假设引起的距离向和方位向几何畸变量，U₃是波前假设引起的方位向二次相位系数，它们只由场景定值参数以及场景目标点坐标x_p,y_p决定，ΔK_y表示方位插值后的距离数波数带宽，Δ₂表示第一个数据抽取比例，Δ₂Δ₃表示第二个数据抽取比例，……Δ₂Δ₃…Δ_N+1表示第N个数据抽取比例；j是虚数符号，u_x是方位向波数间隔。

4.如权利要求1所述的成像方法，其特征在于，步骤四中，所述根据畸变坐标的估计值求解出真实坐标的估计值为：

根据更高阶的信号估计模型，以解析解的形式求解出真实坐标的估计值，所述更高阶的信号估计模型表示为：

其中，i和m分别表示快时间序列和慢时间序列，j是虚数符号，x_p和y_p是成像目标的坐标，U₁和U₂分别为波前假设引起的距离向和方位向几何畸变量，U₃是波前假设引起的方位向二次相位系数，K_x(i,m)为方位向波数，K_y(i)距离向波数；

所述解析解表示为：

其中

其中，n表示第n维方位向信号，n∈[1，2，......N+1]；

Yc表示合成孔径中心时刻雷达在地面的投影到场景中心的距离，β₀表示合成孔径中心时刻的下视角，φ表示雷达速度的俯冲角，δ表示地面斜视角的余角，H表示雷达合成孔径中心时刻平台飞行高度，y_Q为距离向畸变坐标的估计值，x_Q为方位向畸变坐标的估计值。

5.如权利要求1所述的成像方法，其特征在于，所述时变参数体制具体为：沿方位向实时调整雷达参数载频和调频率的方法来抵消由下视角、斜距地面偏角的变化引入的距离向波数的变化；

所述载频和调频率表示为：

其中c为光速，f_c0和γ₀分别表示初始发射信号载频和线性调频信号调频率，β₀表示合成孔径中心时刻的下视角，β(t)表示平台到场景中心的瞬时下视角，α(t)表示瞬时参考斜距在地平面的投影与OY即Y轴所形成的瞬时夹角。

6.如权利要求1所述的成像方法，其特征在于，所述脉冲压缩处理，包括解线频调和去除剩余视频相位操作。

Images