1 引言
合成孔径雷达(Synthetic Aperture Radar,SAR)是一种全天时、全天候、高分辨率的成像雷达,被广泛应用在遥感领域。但是传统的SAR多搭载在卫星和飞机等平台上,这使得雷达的重访时间很长,不能对特定目标持续监测。地基SAR(Ground-Based SAR,GBSAR)是一种小型化、低成本的SAR系统,可实现对特定区域的持续观测[1, 2],结合干涉和差分干涉,GBSAR可以进行形变监测,被广泛应用于塌陷、滑坡等灾害预警[3, 4]。然而,通常GBSAR都安装在直线轨道上,由于轨道长度有限,GBSAR在方位向的观测范围受到了限制。另一种在90年代初提出的旋转式合成孔径雷达(Rotor SAR,ROSAR),通过将雷达置于旋转机械臂上,从而获得了更广的观测范围[5]。ROSAR通常安装在直升机旋翼上,并使用脉冲形式的信号,这使得ROSAR的系统十分复杂。基于ROSAR的运动形式,一种使用调频连续波(Frequency Modulated Continuous Wave,FMCW)体制的改进的旋转式SAR在2012年被提出,并命名为圆弧SAR(Arc-SAR)[6]。
基于FMCW体制的SAR采用两个天线分别收发,进一步降低系统复杂度以及系统对功率的要求[7, 8]。在GBSAR方面,中国科学院电子学研究所(IECAS)成功研制出工作在X波段的FMCW体制的微型SAR(Miniature SAR,Mini-SAR),并用于地基系统。该Mini-SAR具有体积小、质量轻、分辨率高等特点,可以搭载在多个平台上,例如直升机、无人机、汽车等。在2013年,IECAS将Mini-SAR应用到旋转式平台上面,实现了在有限空间内对地面大范围、高精度的持续监测[9]。由于Arc-SAR多用于对矿场边坡的形变监测和高程反演,所以在该场景下,目标的高程不可忽略。文献[10,11,12]的物理模型都假设地面为平坦的,文献[9]的物理模型虽然考虑了地面的高度,但是并没有讨论目标的高程对成像精度的影响。本文主要对FMCW体制的Arc-SAR在高程不可忽略的场景中,给出了一种2维频域的成像算法。该成像算法将场景投影到参考锥面上,可以达到接近BP算法的成像精度。文献[9]虽然也提出了快速成像算法,但是其中做了多次近似,精度有所下降。最后使用实测数据验证了本算法的有效性。
2 Arc-SAR几何结构和物理模型
Arc-SAR把雷达天线放置在一个固定长度的机械臂上,通过机械臂的旋转而实现与目标的相对运动,从而得到多普勒信息。机械臂的长度是相对运动的关键,也是最终图像距离横向分辨率的制约因素。虽然较长的机械臂有助于改善距离横向分辨率,但是机械臂的增长会增加震动,不利于信号的聚焦。
图 1示意了Arc-SAR的几何模型,其中D为俯视锥面;G为地面;H为雷达高度;S为雷达天线的相位中心;θ为S和多普勒中心时刻组成的弧长所对应的圆心角,满足θ=ω(η−τ0)。由Arc-SAR的具体运动模式,可以写出天线与特定点目标的距离关系为:
R(η)=√r20+L2−2r0Lcosβcosω(η−τ0)=r0√1+(Lr0)2−2Lr0cosβcosω(η−τ0)
|
(1) |
其中,r0表示目标与雷达平台旋转中心O的连线;η表示时间,在实际的雷达系统中,η等于方位向慢时间τn与距离向快时间t的和;β表示r0与雷达平台的夹角,以下称之为俯视角;R(η)表示在时间h时雷达天线的相位中心与该目标的距离;L表示旋转机械臂的长度,或是等效旋转半径;ω为平台的旋转角速度;τ0表示最短斜距时刻,也是零多普勒时间。
由于在实际场景中L≪r0,所以对式(1)使用Taylor展开,并保留1阶项可以得到:
R(η)≈r0−Lcosβcosω(η−τ0)
|
(2) |
3 Arc-SAR信号模型及其处理
3.1 信号模型
Arc-SAR使用的信号为FMCW,FMCW体制的信号处理已经有了很多研究[7, 8, 13, 14],这里就不做赘述。接收信号S(t-t0)(t0表示回波延时)会与一个参考信号S(t-tc)(tc为参考时间,或是参考距离rc的延时,满足tc=2rc/c,c为光速)相乘滤波后,回波信号可以写为:
S(τn,t,τ0,r0,β)=σ(τ0,r0,β)exp{−j2πf0(2R(η)c−tc)}⋅exp{−j2πKr(2R(η)c−tc)(t−tc)}⋅exp{jπKr(2R(η)c−tc)2}
|
(3) |
其中,τn表示方位向慢时间;t为距离向快时间;σ表示目标的散射系数;f0表示雷达的载频;Kr表示发射信号的调频率。对式(3)移除残余音频相位(RVP)[15]后,信号为:
S(τn,t,τ0,r0,β)=σ(τ0,r0,β)exp{−j2πf0(2R(η)c−tc)}⋅exp{−j2πKr(2R(η)c−tc)(t−tc)}
|
(4) |
对于FMCW的Arc-SAR,由文献[9]可知,其方位向带宽为:
Ba=4Lωsin(θBeam2)cosβλ
|
(5) |
其中,θBeam为雷达天线方位向的波束宽度;λ为雷达载频信号的波长。对应的方位向的角分辨率和空间分辨率为:
ρag=ωBa=λ4Lsin(θBeam2)cosβρsp=r0ρag=λr04Lsin(θBeam2)cosβ
|
将式(2)代入式(4),并对信号的回波延时作近似(忽略雷达与目标间的距离在信号发射和接收时刻的变化),可以得到:
S(τn,t,τ0,r0,β)σ(τ0,r0,β)⋅exp{−j4πr0λ}⋅exp{j4πLcosβcosω(τn+t−τ0)λ}⋅exp{−j2πKr(2r0c−tc)(t−tc)}⋅exp{j4πLcosβcosω(τn+t−τ0)Krc(t−tc)}
|
(6) |
这就是经过解斜(dechirp)和近似后的回波信号的表达式,该式将用于后续的信号分析和成像处理。
3.2 2维频域处理
由于在文献[9]中使用的BP算法需要做多次的近似,并且需要插值处理,这使得BP算法的计算量非常大。下面给出使用2维频域算法进行成像处理。算法流程图如图 2。
经过dechirp和移除RVP后,式(6)可以写为S(τn,t,τ0,r0,β)=A⋅exp(−jΦ),其中:
Φ(τn,t,τ0,r0,β)=2πf0(2r0−2Lcosβcosω(τn+t−τ0)c−tc)+2πKr(2r0−2Lcosβcosω(τn+t−τ0)c−tc)⋅(t−tc)
|
(7) |
此时的信号在距离时域和方位时域,使用驻定相位原理(Principle Of Stationary Phase,POSP)[16],可得到信号的方位频域表达式。即先求取相位Φ(τn,t,τ0,r0,β)+2πfatn的1阶导数的零点:
∂[Φ(τn,t,τ0,r0,β)+2πfaτn]∂τn=0
|
(8) |
解方程可以得到:
ω(τn+t−τ0)=−arcsin{cfa2Lω[f0+Kr(t−tc)]cosβ}
|
(9) |
使用时-频变换表达:Kr(t−tc)→fr,可以使信号处于2维频域中。则
cfa2Lω[f0+Kr(t−tc)]cosβ→cfa2Lω(f0+fr)cosβ
|
令cfa2Lω(f0+fr)cosβ为W,则式(7)可以写为:
Φ(fa,fr,τ0,r0,β)=4πr0c(f0+fr)+2πfaτ0−4πLcosβ√1−ω2c(f0+fr)−2πfafrKr−2πfaωarcsinω−4πfarcc−4π(f0+fr)rcc
|
(10) |
由于Ω中包含cosβ项,而该项与目标的具体俯视角有关,所以当目标不在参考俯视面上时,2维频域处理的结果不是精确的。
但是在实际中,如果无法得到目标场景的DEM数据,可以使用参考俯视角βref对相位进行补偿,参考信号的形式为:
HRFM(fa,fr,βref)=exp[−ΦR(fa,fr,βref)]
|
(11) |
其中:
ΦR(fa,fr,βref)=4πLcosβref√1−ω2refc(f0+fr)+2πfafrKr+2πfaωarcsinωref+4πfarcc+4π(f0+fr)rcc
|
(12) |
ωref=cfa2Lω(f0+fr)cosβref
|
(13) |
经过补偿后,残余信号的形式为:
G(fa,fr,τ0,r0,β)=GIF(fa,fr,τ0,r0,β)⋅HRFM(fa,fr,βref)
|
(14) |
=σ(τ0,r0,β)⋅exp[−jΦRFM(fa,fr,τ0,r0,β)]
|
(15) |
其中GIF(·)为2维频域中的信号,
ΦRFM(fa,fr,τ0,r0,β)=4πr0λ+4πr0cfr+2πfaτ0−ψ
|
(16) |
ψ=4πL(cosβ√1−ω2−cosβref√1−ω2ref)c⋅(f0+fr)+2πfaω(arcsinω−arcsinωref)
|
(17) |
ψ可以视为相位误差,由于无法知道目标的俯视角,所以这个误差无法补偿。
对式(15)进行2维傅里叶逆变换就可以得到最终的信号聚焦结果。对于参考面上的点目标,聚焦后的点目标形式为:
S(τn,t,β)=σ(τ0,r0,β)⋅exp{−j4πr0λ}⋅pr(t−2r0c)⋅pa(τn−τ0)
|
(18) |
其中pr(·)和pa(·)为sinc类型的函数。
由于式(10)的结果不便于分析,这里给出另一种表达。考虑目标的方位向调频率式(5),由此可以得到对于特定的目标,其方位向带宽与其俯视角有关。由于在Arc-SAR中,天线特有的扫描方式可以保证目标的多普勒中心为零,所以在2维频域中,信号的方位向频谱是对称的。这就表明,不同俯视角下的目标的频谱之间是线性伸缩关系。将式(10)写成形式上不包含 cosb 的式子:
Φ(fa,fr,τ0,r0)=4πr0c(f0+fr)−2π|fa|ω√1ω2−1−2πfafrKr+2πfaτ0−2πfaωarcsinΩ−4πfarcc−4π(f0+fr)rcc=4πr0c(f0+fr)+2πfaτ0−2πfaω(sgn(fa)√1ω2−1+arcsinΩ)−2πfafrKr−4πfarcc−4π(f0+fr)rcc
|
(19) |
其中sgn(·)为符号函数。使用参考俯视角βref对式(19)进行补偿,对应的参考调制函数仍然为式(13)。对于使用参考俯视角的补偿相位函数,相当于使用了参考锥面的一致距离徙动校正,其成像精度和使用相同参考锥面的BP算法一致。相位误差可以认为是由目标俯视角的空变性导致W的方位调频率的不一致而引起的。对特定的参考俯视角,补偿信号在2维频谱的相位为:
ΦR(fa,fr,ωref)=2πfaω(sgn(fa)√1ω2ref−1+arcsinωref)+2πfafrKr+4πfarcc+4π(f0+fr)rcc
|
(20) |
使用式(20)来补偿信号,对应的误差为:
ψ=2πfaω[sgn(fa)(√1ω2−1−√1ω2ref−1)+arcsinω−arcsinωref]
|
(21) |
以验证式(21)和式(17)是等价的。
上述成像平面的参考面的几何形式为锥面,如图 1中的锥面D,非参考面上的目标被投影到参考面上,由于无法知道具体目标的高度,所以无法对非参考平面上的目标准确聚焦。
3.3 算法误差分析
式(17)和式(21)给出了使用参考俯视角对相位作统一补偿所带来的误差,图 3示意了使用雷达转台所在平面为参考面时,俯视角为3°的目标在频域中的相位误差。最大相位误差与俯视角误差之间的关系可以见图 4。从图 4可以看出最大相位误差与俯视角成正相关关系,俯视角差距越大,聚焦效果越差。由式(2)可以得出,目标距离越远,算法精度越高。在实际场景中,由于机械臂长通常较小,所以L≪r0条件容易满足。但对于相同高度的两个目标,距离雷达越远的目标可以得到更小的俯视角,并且关系式(2)也越精确,这就表明,算法对远距目标有更好的成像效果。下一节,通过仿真数据和实测数据验证本算法的适用性。
4 算法验证
4.1 仿真验证与聚焦分析
使用本文提出的2维频域算法对模拟的信号回波进行仿真验证。仿真使用的系统参数为表 1。对仿真点目标成像的图像为图 5,图像未经过几何校正。图中使用了不同高度和不同斜距的点目标。从左到右的点目标对应的最短斜距依次为300,600,900;从上到下对应的高度依次为0,50,100。图 6为图 5经过几何校正的聚焦图像。图 5的部分点目标的放大图和切片图为图 7-图 11。对应的冲激响应宽度(IRW),峰值旁瓣比(PSLR)为表 2。
表 1(Tab. 1)
表 1 仿真参数表Tab. 1Simulation parameters
参数名称 |
参数值 |
参数名称 |
参数值 |
信号载频(GHz) |
9.6 |
平台高度(m) |
100 |
信号带宽(MHz) |
50 |
转动速度(°/s) |
5 |
采样率(MHz) |
1 |
转动半径(m) |
2.5 |
扫描角度(°) |
100 |
参考俯视角(°) |
0 |
脉冲重复频率(Hz) |
12 |
方位向波束宽度(°) |
40 |
|
表 1 仿真参数表
Tab.1 Simulation parameters |
从表 2中可以看出,对于参考面上的目标(图 10,图 11),2维频域算法可以很好地聚焦。对非参考面上的点目标,其散焦程度与目标和天线的相对高度以及目标到雷达的最短斜距有关。并且本文的算法是建立在L≪r0的基础上,即目标越远算法越精确。从图 7和图 8可以看出,虽然两者高度都是0,但是聚焦程度相差很大。雷达所在平面高度为100 m,图 7中的目标的最短斜距为300 m,此时它的俯视角可近似为arcsin(1/3)≈19.5°,而图 8中的目标与雷达相距900 m,它的俯视角近似为arcsin(1/9)≈ 6.4°,所以这验证了聚焦效果与俯视角的依赖关系,即与参考面的俯视角差距越小,聚焦效果越好。图 10和图 11说明了近似条件L≪r0对成像质量的影响,通过对比可以看出,在参考面上的目标,距离雷达越远,其聚焦效果越好。从图 8可以得出,在俯视角差距为6°的情况下,算法的聚焦效果仍然很好。对于距离较近且高度较高的场景,可以使用分层成像来达到较好的聚焦效果。
表 2(Tab. 2)
表 2 点目标的质量参数(已对过采样归一化)Tab. 2Focused target property parameters
最短斜距
r0(m) |
高度
h(m) |
距离向IRW |
距离向PSLR(dB) |
距离横向 IRW |
距离横向PSLR(dB) |
300 |
0 |
0.886 |
-12.961 |
2.827 |
-1.500 |
300 |
100 |
0.885 |
-13.146 |
0.879 |
-12.815 |
600 |
50 |
0.886 |
-13.125 |
0.872 |
-12.772 |
900 |
0 |
0.885 |
-13.123 |
0.883 |
-12.347 |
900 |
100 |
0.884 |
-13.148 |
0.883 |
-12.952 |
|
表 2 点目标的质量参数(已对过采样归一化)
Tab.2 Focused target property parameters |
4.2 实测数据成像对比验证
如图 12显示的是由IECAS研制的安装在室内的Arc-SAR,室内墙壁使用了特殊的对电磁波损耗很小的材料。其机械臂长为2.5 m,扫描角速度为2°~10°/s可调,同时旋转平台还可以在垂直方向上实现1 m的升降,为后期的干涉提供条件。图 14,图 15分别为使用BP算法[9]和本文的2维频域算法对实测数据的成像结果,对应的实际场景为图 13,系统参数为表 3。从图上可以看出实际的目标距离在500 m外,在此区间,本文的2维频域算法可以达到很好的聚焦效果。相比BP算法[9],由于两者都使用了同样的假设条件,所以在没有DEM的辅助下,且都使用斜平面为参考平面时,两者精度一致。依据参考文献[17]的计算算法复杂度的方法可以得出使用BPA和本文的2维频域算法的运算复杂度,表 4给出了BPA和2维频域算法的计算复杂度。其中,用于估算浮点计算的距离向数据为2048点,对应场景最远端约2000 m;方位向点数为8192;对应的观测角度约为80°;插值核长度为8。从表 4可见,BPA的计算量约是2维频域算法的5倍。
表 3(Tab. 3)
表 3 系统参数Tab. 3System parameters
参数名称 |
参数值 |
参数名称 |
参数值 |
信号载频(GHz) |
9.6 |
信号带宽(GHz) |
1 |
采样率(MHz) |
150 |
转动速度(°/s) |
2 |
脉冲持续时间(μs) |
300 |
转动半径(m) |
2.5 |
扫描角度(°) |
100 |
参考俯视角(°) |
0 |
脉冲重复频率(Hz) |
10 |
方位向波束宽度(°) |
40 |
|
表 3 系统参数
Tab.3 System parameters |
表 4(Tab. 4)
表 4 运算复杂度Tab. 4Computation complexity
BPA |
运算复杂度 |
2维频域算法 |
运算复杂度 |
移除RVP |
1.95 |
移除RVP |
1.95 |
距离向IFFT |
0.92 |
方位向FFT |
1.09 |
sinc插值 |
0.50 |
相位补偿 |
0.10 |
相干累加 |
21.19 |
2维IFFT |
1.85 |
全部运算量 |
24.56 |
全部运算量 |
4.99 |
|
表 4 运算复杂度
Tab.4 Computation complexity |
5 结论
调频连续波体制的圆弧式合成孔径雷达是一种新型的用于对地面特定地区持续监测的成像雷达系统,该系统通过两部天线实现了信号的连续发射和接收,降低了宽带脉冲信号对系统功率的要求;同时,旋转械臂带动雷达运动,使其具备方位向合成孔径和360°大范围观测的能力。本文提出了一种在目标斜距远大于旋转机械臂的假设条件下,通过在2维频域补偿相位来聚焦信号的方法。该算法将参考面设定为一个锥面,以此来对回波的频域相位进行补偿。该算法和BP算法在同等假设条件和参考面一致时成像质量相当,并具有计算量小的优点,适合地基系统持续监测的快速成像。