1.   引言

合成孔径雷达(Synthetic Aperture Radar, SAR)是一种主动式微波遥感探测技术，可实现全天时、全天候对地观测，并具有一定的地表穿透能力^[1]。通过SAR系统获取目标的3维几何信息可用于对场景的精细化描述，相对于2维平面成像，具有十分重要的军事应用价值。目前较为成熟的3维重建技术是干涉SAR(Interferometric SAR, InSAR)技术^[2,3]，利用两幅相干SAR复图像同名点的相位差反演高程。SAR层析成像(TomoSAR)技术则是在InSAR技术基础上发展起来的一种目标3维重建技术^[4,5]。该技术弥补了InSAR技术在高度向分辨能力缺失的不足，真正实现了对目标的3维分辨。虽然InSAR和TomoSAR技术发展迅速，但其数据的获取需要通过多次航过或多个天线完成，增加了技术实现的难度和成本。

SAR多方位角观测是近年来提出的一种新的观测模式，对同一区域进行不同方位角的观测，不仅能获取场景内更丰富的目标特征信息，还使其具备了3维重建的潜力^[6,7]。圆迹SAR(Circular SAR, CSAR)是多方位角观测的一种特殊形式，通过传感器平台的圆周运动，获取被观测目标多方位乃至360°全向观测信息，突破了常规直线轨迹SAR 2维成像的局限^[8,9]。然而，对于星载SAR系统而言，难以实现单次航过环绕观测场景的多方位角观测。为实现全方位角观测，星载SAR系统需要依靠多航过或多星配合完成。尽管单次航过只能实现部分方位角的观测，但由于卫星轨道自身的弯曲特性，可将其等效为长3维曲线阵列^[9,10]，因此观测数据中依然携带了目标3维几何信息。

实现3维成像的最直接方法是基于回波数据与散射目标系数间的傅里叶变换关系，通过3维FFT实现^[11]。然而，由于多角度SAR数据在高度维是稀疏的，无法通过傅里叶变换实现无模糊3维成像。针对稀疏采样问题，文献[9,12]采用了稀疏重建技术恢复3维几何信息，但其方法基于多航过SAR数据，且运算量大，不适用于实时处理。文献[6]通过合理离散目标函数构造测量矩阵，然后采用稀疏重构算法实现多角度SAR 3维成像，但该方法受限于测量矩阵的设计，必须解决目标稀疏性与发射信息匹配的问题。文献[13]利用了高程相同的目标在多角度SAR图像序列上位置偏移一致的特性，提出了基于低秩和稀疏矩阵分解的目标3维重建方法，但其位置偏移与目标高度的传递关系建立在圆迹SAR观测几何下，不能直接推广至星载SAR系统。文献[14,15]分析了散射中心的分布特性，从非相干处理图像中提取轮廓线，从而估计目标的3维几何信息，但该方法需要建立在已知观测目标模型的基础上，无法重建未知目标。文献[16]通过分析TerraSAR-X凝视聚束模式数据，根据星载SAR轨道特性，给出了方位调频率误差与高程误差的函数关系，提出了一种基于参数估计的高程提取方法，但是其信号模型建立在正侧视几何下，未考虑多角度观测的斜视情况。

针对上述问题，本文提出了一种联合多方位角调频率估计的星载SAR3维成像方法。该方法根据星载SAR多方位角观测几何，首先给出了更精确且适用于斜视几何的多普勒调频率与高程误差间的传递函数，实现联合多方位角SAR子孔径数据的高程估计。然后，根据2维成像的投影几何，实现目标的3维成像。最后，仿真数据验证了方法的有效性。

本文结构如下：第2节介绍多方位角观测星载SAR系统的几何特性；第3节分析星载SAR多方位角观测几何下高程误差与多普勒调频率的关系；第4节，详细地阐述3维成像的全流程；第5节以圆柱点阵目标为例，给出3维成像结果；最后一节对全文进行总结。

2.   星载SAR多方位角观测几何

多方位角观测星载SAR系统通过波束方位向扫描实现对同一场景的长合成孔径观测，其数据录取几何如图1所示。图1中所有坐标都是建立在场景局部坐标系 ${{O {{\raisebox{1.5pt}{\tiny{-}}}} \!XYZ}}$ 下；坐标系原点 $O$ 位于观测场景内某处， $X$ 轴所指方向为地距向， $Y\,$ 轴指向方位向， $Z$ 轴指向高度向；卫星沿曲线轨道飞行，从 ${{{P}}\!_1}$ 处开始持续对地面场景进行观测，直至 ${{{P}}\!_2}$ 处停止观测，两处卫星速度矢量分别为 ${{{V}}\!_1}$ 和 ${{{V}}\!_2}$ , ${\theta _1}$ 和 ${\theta _2}$ 为两处斜距矢量与速度矢量的夹角，观测方位角跨度为 $\Delta \varphi$ ；假设场景参考平面对应高程为0；场景内有一个目标位置矢量为 ${{T}}$ ，相对参考平面高程为 $h$ ，且目标到 ${{{P}}\!_1}$ 和 ${{{P}}\!_2}$ 的斜距分别为 ${R_1}$ 和 ${R_2}$ 。若以场景参考平面为成像地平面，对 ${{{P}}\!_1}$ 所在方位角子孔径数据进行成像，该目标在成像平面内的投影位置为 ${{{T}}\!_1}$ ；对 ${{{P}}\!_2}$ 所在方位角子孔径数据进行成像，则该目标在成像平面内的投影位置为 ${{{T}}\!_2}$ 。目标真实位置和两投影位置满足距离-多普勒模型^[17]

图  1  星载SAR多方位角观测几何示意图

Figure  1.  Space-borne SAR multiple azimuth angle observation geometry

下载: 全尺寸图片 幻灯片

式中， $\lambda$ 为雷达载波波长。从式(1)可以看出，投影方向与速度矢量垂直。由于卫星轨道的弯曲，速度矢量 ${{{V}}\!_1}$ 和 ${{{V}}\!_2}$ 不平行，因此两次投影方向不同，即有 ${{{T}}\!_1} \; =\!\!\!\!\!\!/ \;{{{T}}\!_2}$ ，可以根据投影位置偏移量提取目标高程。对于CSAR的圆轨迹，目标在不同方位角SAR图像中的偏移量与高程的关系可以通过几何关系求出^[13]，但这在多方位角观测星载SAR系统中直接推导却十分复杂。因此，可以从信号模型出发，反推偏移量与高程的关系。

3.   目标高程与方位向偏移量的关系

在多方位角星载SAR数据录取过程中，卫星的瞬时速度方向不断变化，对于距卫星同一斜距但高程不同的两个目标而言，多普勒历程是不同的。常见的多普勒频率表达式为

式中， ${{{R}}_{\rm{S}}}$ 为卫星到目标的瞬时斜距矢量， ${{{V}}\!_{\rm{S}}}$ 为卫星的瞬时速度矢量， ${R_{{\rm{st}}}}$ 为卫星到目标的瞬时斜距。因此，多普勒调频率 ${K_{{f_{\rm{d}}}}}$ 可以表示为

式中， ${V_{\rm{s}}}$ 为卫星瞬时速度， ${a_{\rm{s}}}$ 为卫星的瞬时加速度， ${\alpha _{{\rm{nad}}}}$ 为目标瞬时的雷达下视角， ${\alpha _{{\rm{nad}},{\rm{c}}}}$ 为目标在中心方位时刻(正侧视)的雷达下视角， $\theta$ 为斜距矢量与速度矢量的夹角之间的夹角。多方位角观测星载SAR对目标的观测经历了由前斜视到后斜视的过程，不同方位时刻目标的瞬时下视角与斜视情况相关，即 $\cos {\alpha _{{\rm{nad}}}} \approx \sin \theta \cos {\alpha _{{\rm{nad}},{\rm{c}}}}$ 。假设观测方位角跨度为[–20°, 20°]，如果在推导多普勒调频率时忽略斜视角的影响，最大将引入约6%的调频率误差。

对于匀速直线飞行的SAR平台，式(3)中带有 ${a_{\rm{s}}}$ 的第2项不存在；对于卫星平台， ${a_{\rm{s}}}$ 为卫星的重力加速度。在二体模型下，根据开普勒第三定律，卫星重力加速度为

式中， $G$ 为引力常量(6.674×10^–11 N·m²/kg²), $M$ 为地球质量(5.964×10²⁴ kg), ${R_{{\rm{sat}}}}$ 为卫星轨道半径。

假设目标的高程误差为 $\Delta h$ ，则引入的调频率误差 $\Delta {K_{{f_{\rm{d}}}}}$ 为

式中， ${K_{{f_{\rm{d}}},{\rm{bias}}}}$ 为存在偏差的多普勒调频率， $\Delta {\alpha _{{\rm{nad}}}}$ 为高度误差 $\Delta h$ 引起的雷达下视角误差。

目标高程误差与雷达下视角的几何关系如图2所示。图2中 ${H_{\rm{s}}}$ 为卫星飞行轨道高度， ${R_{\rm{e}}}$ 为该纬度对应的地球半径， ${\alpha _{{\rm{inc}}}}$ 为目标对应的局部入射角， ${{T}}$ 和 ${{{T}}'}$ 分别在目标真实位置和存在高程误差时位置。当 $\Delta {\alpha _{{\rm{nad}}}} \approx 0$ 时， ${{T}}{{{T}}'}$ 可以近似看作直线。结合地球椭球模型^[18]，高程误差引起的下视角误差可以表示为

图  2  高度误差与下视角误差关系几何示意图

Figure  2.  The geometry of height error and nadir angle error

下载: 全尺寸图片 幻灯片

联立式(5)和式(6)可得

式(7)中，多普勒调频率误差与高程误差呈线性关系，与斜距成反比。因此，在方位压缩后，通过对多普勒调频率误差进行估计，即可求出目标真实高程值。

图3为仿真实验与式(7)计算的对比结果。仿真使用的轨道高度为514 km，轨道倾角为97.4°，观测场景选在赤道附近。从图3(a)中可以看出，200 m以内的高程误差引起的调频率误差小于0.2 Hz/s²，调频率误差随下视角减小而增大；图3(b)显示，调频率误差仿真值与式(7)计算值存在约5%的残差。该残差的引入除了推导过程中的两次近似外，还有加速度 ${a_{\rm{s}}}$ 的近似。尽管使用更精确的 ${a_{\rm{s}}}$ 可以提高计算精度，但这却不是必需的。虽然在 ${a_{\rm{s}}}$ 不精确的情况下，单次计算将引入约5%的残差，但可以通过迭代使其快速收敛。假设初始高程残差为200 m，通过2次迭代，由参数不精确引入的高程残余误差即可收敛至2 m以内。

目标高程估计的前提是多普勒调频率估计，常见多普勒调频率估计的方法有最大对比度法和视错位法(Map Drift, MD)^[19]。因为前者无法直接估计调频率的误差方向及幅度，所以此处选用MD法。MD法认为多普勒调频率误差 $\Delta {K_{{f_{\rm{d}}}}}$ 会导致两个方位子孔径图像之间出现方位向偏移。令数据实际的多普勒调频率为 ${K_{{f_{\rm{d}}}}}$ ，前后两子孔径的多普勒中心间隔为 $\Delta {f_{\rm{d}}}$ ，卫星地速为 ${V_{\rm{g}}}$ ，则多普勒调频率误差可以近似为

式中， $\Delta {d_{{\rm{az}}}}$ 为子孔径图像间的方位向偏移量(距离)。结合式(7)和式(8)，可以得出高程误差估计表达式为

利用高程误差估计结果修正成像初始参考高程可获得目标真实高程。然后，将目标真实高程和目标在SAR图像中的像素坐标转换成目标3维几何信息，进而完成3维重建。在进行子孔径图像偏移量估计时，噪声的存在将引入偏移量的估计误差，从而影响高程测量精度，可通过切割多个子孔径，并通过对各子孔径的高程估计值求期望将该影响降低。对于多方位角观测SAR系统，由于其观测方位角跨度大，在确保足够子孔径分辨率的同时可以获得足够多的子孔径序列，可以获得更高的测高精度。在实际观测中，背景杂波对目标高程估计的影响相对较大。当杂波纹理特征明显，且杂波对应的场景高程与目标高程相同时，杂波在子孔径图像间的偏移量与目标的偏移量一致，不影响目标高程的估计；当杂波纹理特征明显，且杂波对应的场景高程与目标高程不同时，杂波在子孔径图像中的偏移量与目标不一致，即在计算偏移量时将出现两个峰值，分别属于杂波和目标，峰值大小与信杂比有关，峰值位置与杂波和目标的高程相关，不具有随机性，无法通过多角度观测提高高程估计精度；当杂波纹理特征不明显时，杂波在不同方位角子孔径图像间的配准偏移量具有随机性，可以通过对多方位角估计值求期望降低误差，即本文方法对于该类杂波的影响同样有一定的抑制效果。

4.   处理流程

本文给出了一种多方位角观测星载SAR的3维成像方法，其流程如图4所示。

图  4  多方位角观测星载SAR的3维成像方法流程图

Figure  4.  The schematic diagram of space-borne SAR 3-D imaging method by multiple azimuth angle observation

下载: 全尺寸图片 幻灯片

该方法包含以下步骤：

步骤1　(数据子孔径分割)　根据星载多方位角观测SAR原始数据的孔径长度确定分割子孔径数。为确保各子孔径数据方位聚焦效果对高程的敏感性，子孔径方位角跨度应较大；同时，为了保证子孔径数据内的相干性和划分后子孔径数目足够多，子孔径方位角跨度不宜过大。根据经验，本文以2°作为子孔径分割的参考值，将全孔径数据分割成相同方位向点数的子孔径数据序列 ${{{S}}_1},{{{S}}_2}, ·\!·\!· ,{{{S}}\!_{{N_{{\rm{sub}}}}}}$ ，其中 ${N_{{\rm{sub}}}}$ 为分割后子孔径数。

步骤2　(参考子孔径成像)　从 ${{{S}}_1},{{{S}}_2}, ·\!·\!·,{{{S}}_{{N_{{\rm{sub}}}}}}$ 中选取中间的子孔径数据作为参考子孔径数据 ${{{S}}_{{\rm{ref}}}}$ ，以场景参考平面为参考成像平面进行低分辨率2维成像，获得参考子孔径SAR图像 ${{{I}}_{{\rm{ref}}}}$ 。参考成像平面的高程可以通过先验粗DEM获取，作为高程估计的起始值。当缺少先验DEM时，也可以通过大尺度的高程搜索获取起始值。

步骤3　(强散射点搜索)　考虑到照射场景内地形和地物的空变，先将 ${{{I}}_{{\rm{ref}}}}$ 沿距离向和方位向进行分块，分别统计各块内的像素幅值，并设定相应的幅度门限，利用幅度提取强散射点。然后，合并所有强散射点，得到强点序列 ${{{p}}_1},{{{p}}_2}, ·\!·\!· ,{{{p}}_{{N_{\rm{p}}}}}$ ，其中 ${N_{\rm{p}}}$ 为强点个数。

步骤4　(子孔径序列成像)　在参考成像平面内对子孔径序列 ${{{S}}_1},{{{S}}_2}, ·\!·\!· ,{{{S}}_{{N_{{\rm{sub}}}}}}$ 进行2维成像，获得子孔径图像序列 ${{{I}}_1},{{{I}}_2}, ·\!·\!· ,{{{I}}_{{N_{\rm sub}}}}$ 。

步骤5　(相邻子孔径分块配准)　从每个子孔径图像中取出 ${N_{\rm{p}}}$ 图像块，图像块的中心分别对应 ${N_{\rm{p}}}$ 个强散射点，图像块大小为 $5 \times 5$ 个像素，整个孔径图像序列共可以获得 ${N_{{\rm{sub}}}} \times {N_{\rm{p}}}$ 个图像块。然后，对 ${N_{{\rm{sub}}}} \times {N_{\rm{p}}}$ 个图像块进行升采样，并配准插值后相邻子孔径图像同一强点所在的图像块，得到 $\left( {{N_{{\rm{sub}}}} - 1} \right) \times {N_{\rm{p}}}$ 组偏移量。

步骤6　(基于配准偏移量的高程估计)　将偏移量代入式(9)，估计出每个强散射点的 ${N_{{\rm{sub}}}} - 1$ 个高程值。

步骤7　(多方位角高程估计值平均)　对于每个强点，将 ${N_{{\rm{sub}}}} - 1$ 高程估计值取平均，认为是该强点的最终高程估计结果。

步骤8　(3维信息重建和点云显示)　根据强点在 ${{{I}}_{{\rm{ref}}}}$ 中的位置和高程估计值，计算各点在场景中的3维坐标，最终以3维点云的形式显示。

5.   仿真结果

仿真1　单点目标仿真　仿真参数如表1所示，表中的雷达波束中心下视角特指雷达波束正侧视照射场景时对应的波束中心下视角，仿真的点目标处于波束中心，设置不同信噪比和信杂比进行仿真实验。在不同信噪比的仿真实验中，SAR图像中加入的噪声为高斯白噪声，信噪比设有5 dB, 10 dB, 15 dB, 20 dB, 25 dB 5个梯度值，分别进行50次高程提取蒙特卡洛模拟实验，以高程估计的均方根误差(Root Means Square Error, RMSE)作为评价指标。在进行杂波仿真实验时，仿真方法与噪声的仿真方法不同，不直接向SAR图像中添加杂波。 $\kappa$ 分布是目前应用最广的分布模型之一^[20]，在高分辨情况下可以在很宽范围内匹配杂波数据的幅度分布。因此仿真实验中使用服从 $\kappa$ 分布的杂波模型生成场景的散射场。然后，用该散射场与目标进行多角度回波仿真。仿真数据的信杂比设有10 dB, 15 dB, 20 dB, 25 dB 4个梯度值(5 dB信杂比将导致目标淹没在杂波中)，各进行50次蒙特卡洛实验，同样以高程估计的均方根误差作为评价指标。本文中信杂比的定义为信号功率与杂波功率(杂波幅度的均方根)的比值。

在不同信噪比仿真实验中，为验证本文方法在斜视处理时相对于文献[16]方法的优势，截取了中心方位角为16°，方位角跨度为4°对应的斜视回波数据，分别用文献[16]方法和本文方法进行处理，对比结果如图5(a)所示。由于文献[16]中方法的提出未考虑斜视几何，因此相比于本文方法在高程估计上存在更大的偏差。为进一步验证联合多方位角观测对提高测高精度的有效性，开展了单方位角观测和多方位角观测数据处理对比实验。其中，单方位角观测的方位角范围为[–2°, 2°]，处理时等分为两个子孔径，每个子孔径宽度为2°；多方位角观测的方位角范围为[–16°, 16°]，等分为16个子孔径，每个孔径宽度为2°，处理结果如图5(b)所示。可以看出，联合多方位角观测估计高程的RMSE比单方位角观测估计结果小，证明多方位角观测具有更好的抗噪性能；随着信噪比的升高，噪声的影响减弱，联合多方位角数据处理结果与单方位角数据处理结果逐渐逼近，且两者的RMSE不断下降，最终趋近于1 m，该残差与图3的分析一致。当信噪比大于10 dB时，多方位角观测高程估计精度优于2 m。

图  3  不同下视角对应的多普勒调频率误差与高程误差

Figure  3.  The relationship between Doppler frequency error and height in different nadir angles

下载: 全尺寸图片 幻灯片

图  5  不同信噪比下的高程估计均方根误差

Figure  5.  Height estimation RMSE in various signal noise ratio

下载: 全尺寸图片 幻灯片

在不同信杂比仿真实验中，用上述方法仿真的匀质杂波如图6(a)所示。图6(b)给出了杂波像素幅值统计结果和 $\kappa$ 分布概率密度曲线的对比，在统计时将杂波幅值范围等分成80个区间，对比结果显示仿真的杂波像素幅值服从 $\kappa$ 分布。实验数据后续的处理细节与信噪比仿真实验的数据处理细节一致，不赘述。图7(a)和图7(b)记录了50次单方位角数据处理和多方位角数据联合处理的蒙特卡洛实验结果，可以看出基于多方位角数据高程提取结果的波动范围要小于单方位角数据的高程提取结果的波动范围，即有更稳定的估计结果。图8给出了不同信杂比下50次蒙特卡洛实验的高程估计均方根误差，从图中的曲线可以看出基于多方位角数据联合高程估计结果精度优于单方位角数据高程估计精度，证明了本文所提方法在目标背景中含有纹理特征不明显杂波的条件下依旧有效。

图  6  杂波仿真结果

Figure  6.  Simulation result of clutter

下载: 全尺寸图片 幻灯片

图  7  不同信杂比下的蒙特卡洛高程估计结果记录

Figure  7.  Monte Carlo height estimation record in various signal clutter ratio

下载: 全尺寸图片 幻灯片

图  8  不同信杂比下的高程估计均方根误差

Figure  8.  Height estimation RMSE in various signal clutter ratio

下载: 全尺寸图片 幻灯片

仿真2　圆柱点阵仿真　主要仿真参数同表1，仿真的SAR图像信噪比为15 dB，在添加噪声时，参考的信号功率是点阵目标所在像素幅值的均方根值。图9显示了仿真的半径40 m、高20 m圆柱点阵模型，颜色对应高度。图10(a)是参考子孔径成像的幅度图，该图用于提取强散射点位置；图10(b)显示了强点提取结果，由于噪声存在和各点聚焦质量的差异，漏检了部分点目标，并且少量噪点被检测为目标；图10(c)和图10(d)分别是文献[16]方法和本文方法提取的强散射点最终的3维成像结果，两种方法都能在一定程度上还原点阵模型，其中本文所提方法的3维成像结果轮廓更清晰，证明本文方法的高程提取精度优于文献[16]。

图  9  圆柱点阵仿真模型

Figure  9.  The model of cylider scatter array

下载: 全尺寸图片 幻灯片

下载: 全尺寸图片 幻灯片

表  1  仿真参数

Table  1.  Simulation parameters

参数	数值
轨道高度(km)	514
轨道倾角(°)	97.4
雷达载频(GHz)	9.70
雷达带宽(MHz)	600
雷达波束中心下视角(°)	35
数据方位角跨度(°)	[–16, +16]
观测场景纬度(°N)	0

下载: 导出CSV 

| 显示表格

6.   结论

针对星载SAR多方位角观测在高度维采样的稀疏性，无法直接通过3维FFT实现无模糊成像，且目标在不同方位角的投影与高程间的传递相对复杂等问题，本文提出了一种联合多方位角调频率估计的星载 SAR 3维成像方法。该方法首先给出了多普勒调频误差与高程误差进行关联，然后通过MD法估计多普勒调频率，接着通过联合多方位角高程估计值提高高程估计精度，最终根据高程估计结果和参考子孔径SAR图像恢复目标3维信息，实现3维成像。该方法利用了卫星轨道几何特点，建立了多普勒调频率与高程误差的关系，避免了求解目标在不同方位向投影位置的困难。同时，发挥了多方位角数据观测方位角大的优势，提高了高程测量精度，仿真实验验证其高程测量精度可达米级。本文3维成像的对象是强点目标，在实际场景中大多是扩展点目标(面目标)，如何实现对这些目标的3维成像将是下一步工作的重点。