CN105180852B - 基于三重步进的gb‑sar形变监测方法 - Google Patents

Abstract

本发明提供一种基于三重步进的GB‑SAR形变监测方法，雷达系统设计将一种基于三重步进的GB‑SAR用于形变测量：第一重步进指采用步进频连续波信号，合成大的信号带宽，收发天线分离提高隔离度；第二重步进指将雷达安装在线性导轨上，以停‑走‑停方式步进运动，合成较大孔径提高方位分辨率，雷达在停止期发射和接收信号，避免运动期间收发信号带来的相位误差；第三重步进指雷达天线沿方位向依次步进转动，覆盖较宽的监测区域，不同天线位置的数据拼接后用做成像处理。对采集的数据的成像处理是基于GB‑SAR系统的回波信号特点进行的，不需要补偿脉冲间的瞬时斜距差，计算量相对较低，可以满足GB‑SAR系统实时性要求。

Description

基于三重步进的GB-SAR形变监测方法

技术领域

本发明涉及合成孔径雷达SAR成像技术。

背景技术

步进频连续波(Stepped Frequency Continuous Wave，SFCW)GB-SAR(Ground-Based SAR)是一种放置在轨道上的地基合成孔径成像雷达，通过持续发射与接收步进频连续波信号获得大的带宽和高距离向分辨率，通过雷达在导轨上的步进运动合成大的孔径和高方位向分辨率，通过两次运动的成像作干涉处理获得目标的形变或微位移测量值，由于雷达和导轨安装在固定位置，可进行连续高精度测量，实现对特定区域或目标的监测，预警可能发生的自然灾害或基础设施的险情。该雷达系统具有测量精度高、监测区域广、可接触或无接触测量、全天候工作、使用方便的优点，应用领域十分广泛。典型的应用有：滑坡监测,如山体滑坡，滑坡滚石，矿山稳定性等；地形和地表沉降监测；大坝、桥梁监测、大型建筑物、高塔、历史古迹的形变和沉降监测；火山、地震、冰川和雪崩监测；线性物监测，如桥梁斜拉索、“西气东输”管道、高铁等重要的大型线状物。

GB-SAR采集的数据进行成像处理来获得SAR图像，其成像算法主要有：后向投影算法(back-projection algorithm，BPA)和以距离多普勒算法与宽带合成信号算法相结合的算法。通常的形变估计方法是：对两次成像的图像做干涉处理得到干涉图，经相位解缠后对相邻时间的干涉图进行相位估计，再由估计的相位计算出监测目标或区域的形变值，但该方法中相位解缠误差会在相位估计的过程中传播，并且不能控制；估计的相位由三部分构成：形变对应的相位、大气相位和噪声。大气相位部分主要是由于大气的基本性质在不同图像获取期间的变化干扰GB-SAR信号，这一部分通常称为：大气相位屏(Atmospheric PhaseScreen，APS)，主要影响传感器到目标的距离。通常采用的是距离的二次函数模型，仅用两个地面控制点(Ground Control Points，GCPs)来估计，但实际上大气相位不仅是距离的函数，也是方位的函数。

GB-SAR形变监测雷达的难点在于：(1)受外场作业、安装和成本等因素的影响，放置雷达的轨道不能太长，但轨道太短则合成孔进长度小，方位分辨率低，方位向的单元尺寸较大，导致形变图粗略，不能反映具体的形变位置。降低雷达波长也能提高方位向分辨率，但波长太短则对气象条件敏感，降低测量距离且大气相位误差增大。(2)提高距离向分辨率需采用大的信号带宽，常规的脉冲SAR代价较大,其成像算法也不能用于步进频连续波雷达。(3)监测范围是以雷达为中心的扇形区域，通常天线方位向波束宽度仅几度，导致监测区域太小，只能做条带区域监测。降低天线尺寸可增加波束宽度，例如，波束宽度为1弧度或57.3°时，天线尺寸仅与波长相当，天线尺寸太小将导致天线增益降低，性能变差。(4)监测单元的形变值由干涉相位计算出，相位误差将直接影响测量精度，因此，克服相位解缠误差的传播，抑制大气相位噪声十分重要。

发明内容

本发明所要解决的技术问题是，提供一种测量精度更高的形变监测系统。

本发明为解决上述技术问题所采用的技术方案是，基于三重步进的GB-SAR形变监测方法，包括以下步骤：

数据采集步骤：雷达安装在线性导轨上以停-走-停方式步进运动，雷达天线沿方位向步进转动，雷达在进行每一次步进后停下来发射与接收一帧步进频率连续波信号；雷达工作在Ku波段；雷达天线沿方位向完成探测区域中全部位置的步进后，拼接雷达天线各方位向步进转动后获取的步进频率连续波信号完成一次成像数据采集；

成像处理步骤：对每次采集的成像数据在方位向上快速傅里叶变换后进行距离徒动校正、方位匹配与距离向聚焦，最后在方位向上进行快速傅里叶逆变换得到二维成像图像；

形变估计步骤：对二维成像图像顺序进行配准，再顺序冗余选择配准后的图像两两组成图像对；一组图像对包括一幅主图像与一幅副图像，一组图像对应生成一幅干涉图；对选择的干涉图的低噪声相位进行解缠，再对解缠后的干涉图进行相位估计，通过对所有解缠后干涉图的相位估计得到配准后的相对于第1幅图像的每幅二维成像图像的相位估计值；根据原始相位估计值是监测对象形变相位值与大气相位屏之和，在相位估计值上移除估计的大气相位屏得到监测对象形变相位值。

本发明雷达系统设计将一种基于三重步进的GB-SAR用于形变测量：第一重步进指采用步进频连续波信号，合成大的信号带宽，收发天线分离提高隔离度；第二重步进指将雷达安装在线性导轨上，以“停-走-停”方式步进运动，合成较大孔径提高方位分辨率，雷达在停止期发射和接收信号，避免运动期间收发信号带来的相位误差；第三重步进指雷达天线沿方位向依次步进转动，覆盖较宽的监测区域，不同天线位置的数据拼接后用做成像处理。对采集的数据的成像处理是基于GB-SAR系统的回波信号特点进行的，不需要补偿脉冲间的瞬时斜距差，计算量相对较低，可以满足GB-SAR系统实时性要求。

进一步的，本发明提出一种形变值估计过程中相位解缠后的相位估计采用基于奇异值分解的最小二乘法的相位估计算法，与现有直接相位估计相比，可以检测和校正相位解缠误差，避免误差的传播。

本发明的有益效果是，三重步进实现形变监测能提高距离分辨率与方位分辨率，增加监测的范围。

附图说明

下面结合附图和实施例对本发明进一步说明。

图1是本发明的GB-SAR监测示意图。

图2是本发明的天线转动步进示意图。

图3是GB-SAR发射SFCW信号示意图。

图4是本发明的GB-SAR的成像处理流程图。

图5是本发明的形变监测流程图。

图6(a)是从幅图像中顺序选择两幅图像进行配准的简单网络；(b)从幅配准后图像选择两幅图像的冗余网络。

图7是基于奇异值分解的最小二乘法的相位估计算法。

具体实施方式

雷达工作方式和设计约束：

雷达为步进频连续波体制，采用独立的发射天线和接收天线，天线沿方位向可步进转动，整个雷达安装在长度为L的线性导轨上，如图1所示。监测区域为以雷达为中心点，半径为R_max(雷达最大作用距离)的扇形区域，雷达沿距离向和方位向将监测区域划分为若干个分辨单元。雷达在导轨上按“停-走-停”模式步进运动，每次移动距离ΔL，然后停下来发送和接收一帧SFCW信号，步进频率数目为N，如图3所示。每帧信号的第一个脉冲的载频为f_c，频率步进值为Δ_f，子单脉冲的载频为f_n＝f_c+(n-1)Δf。雷达在导轨上走完单程后返回起点，将收发天线转动角度Δφ，重复上述步进运动和信号收发过程，当天线方位向各位置步进结束后，将各方位向步进转动后获取的数据拼接在一起，完成一次成像数据的采集。共采集P幅图像的数据，分别做合成孔径成像处理，获得P幅SAR图像。对P幅图像做配准后，生成干涉图和相干图，在相位解缠绕后做相位估计，然后去除大气相位噪声，由此计算各监测单元的形变值，所有监测单元的形变值构成监测区域的沿雷达视距的形变图。

雷达系统的设计约束如下。

工作频段：Ku波段；

带宽：其中Δr为距离分辨率，c为光速；

步进频率间隔：其中R_un为最大不模糊距离；

步进频率数目：

轨道长度：其中Δθ为方位向分辨率，f_c为载波中心频率；

雷达在轨道上的步进运动间隔：其中θ_B为天线波束宽度，其中D为天线实孔径；

雷达在轨道上的步进运动数目：

雷达步进运动时间间隔：其中τ为雷达发射N个步进频信号的持续时间；

天线步进转动角度：Δφ≥θ_B；

天线步进转动数目：其中φ为扇形监测区域的方位向角度；

监测区域面积

GB-SAR的成像处理过程：

第一步：方位快速傅里叶变换FFT(Fast Fourier Transform)。回波信号在方位向上是具有线性调频特性的时域信号，根据驻定相位原理，对回波信号进行方位FFT将数据变换至多普勒域。设脉冲重复周期为T_r，即脉冲宽度T_P。以导轨的方向为X轴，其垂直的方向即距离向为Y轴，建立二维的坐标系。对于位置坐标为(0,R_ob)的目标，雷达收到单窄带脉冲的回波信号经过解调到基带信号，忽略幅度信息，为：

其中，m表示雷达发射的第m帧SFCW信号，m＝0，1，…M-1。n表示在第m帧中发射的第n个子脉冲，n＝0，1，…N-1。R(m)表示雷达与目标之间的距离，由图1和图2可知，

由于R_0b>>L/2，R(m)可以近似为：

将其代入式(1)中，信号为：

从第二项指数项可以看出方位相位调制，相位是[(M-1)/2-m]²的函数，故回波信号在方位向具有线性调频特性，调频率为：

利用驻定相位原理，得到方位向上的时频关系为：

f_a＝K_a[(M-1)/2-m]

方位FFT后信号为：

第二步：距离徒动校正。雷达在导轨上步进运动，其与监测区域目标之间的瞬时距离随着其运动而改变，在雷达照射时间内目标的轨迹经过不同的距离单元，需要校正目标的轨迹在相同的距离单元，即距离徒动校正RCMC(Range Cell Migration Correction)。

距离多普勒域中，式(2)变为：

式(7)中第二项距离徒动，需要校正。可以通过sinc函数插值进行校正。距离徒动校正是精确的，信号变为：

第三步：方位匹配函数相乘。方位匹配函数相乘，即通过匹配滤波器进行数据的方位聚焦。由于RCMC后数据处于多普勒域，因而在该域中进行方位匹配滤波效率较高。

方位的聚焦通过频域的匹配滤波，匹配滤波器为式(10)第三项指数项的复共轭，

为了进行方位压缩，将式(8)与式(9)相乘，得到

式子中，第一项与第三项表示的是距离向信息，第二项与第四项表示的是方位向信息。距离与方位没有耦合，后面分开处理，不影响结果。

第四步：距离向IFFT。通过距离向的快速傅里叶逆变换IFFT(Inverse FastFourier Transform)进行数据的距离向聚焦。雷达对停止期接收的N个子脉冲回波信号进行采样，提取N个离散采样值进行IFFT后的得到sinc函数，所以雷达的回波信号在距离向是频域信号，进行IFFT就可以完成距离压缩。

对信号进行距离IFFT，得到sinc函数，就完成了距离压缩。

其中，l＝0,…N-1。

第五步：方位IFFT。通过方位向的IFFT将数据变换到方位时域，得到二维压缩后的图像。

其中，h＝0,…M-1。

最后，信号的经过距离和方位压缩后，为：

成像算法如图4所示。

监测对象的形变测量

本发明的形变监测方法的主要步骤如图5所示，S1对成像处理得到SAR图像的配准；S2干涉图和相干图生成；S3相位解缠和相位估计；S4APS估计；S5形变值计算。

S1图像配准：

雷达在导轨上重复运动，经过P次数据采集，得到P组回波数据,对获得的P组数据分别进行成像处理。从P幅SAR图像选择图像对，选择方式按照图6(a)所示的方式顺序选择，总共有P-1对图像，对每对图像分别进行配准，得到P幅配准图像。

配准算法采用以获得高质量的干涉条纹为出发点的相关配准方法。要生成高质量的干涉条纹，需要对SAR图像先进行内插处理。GB-SAR处理中采用高精度的三次卷积内插法，利用多项式S(ω)来逼近理论上的最佳插值函数sinc(x)，其中

利用周围的16个点灰度进行内插。由两幅SAR主副图像的幅度信息的图像按逐个像素点方式二维显示出来，从中选出两幅图像中最具特点、最相近的两个区域的中心点(x₁,y₁)和(x₂,y₂)。根据选出的这两个坐标点，分别在第一幅和第二幅图像中取出两块子图像M₁和M₂，以M₁为“基模板”，以M₂为“搜索模板”，计算两者之间的相似度。搜索模板搜索整个基模板后，从中选出最大的相似度值，其所对应的基模板和搜索模板之间的位移坐标Δx和Δy。根据Δx和Δy，对两幅SAR图像进行匹配，完成了两幅图像的配准处理。

S2干涉图和相干图的生成：

对图像进行干涉处理。从P幅已经配准图像中，按照一定方式冗余选择图像对，生成Q幅干涉图和相干图。

从P幅配准图像中按照图6(b)所示方式选择图像对，得到幅干涉图。选择其中一对配准图像，U和S，干涉图为：

I＝U·S^*

干涉相位为：

其中，为干涉缠绕相位U为主图像，S^*为副图像的共轭。干涉相位的具体表达式为：

其中，是获取两幅图像期间发生的形变对应的相位，是大气相位；是噪声相位；2kπ是相位解缠项，k是未知整数。

干涉图中给定像素点的干涉相位是否可以利用取决于低的干涉相位的像素点可以选择。像素点的选择标准是基于相干图的像素的幅度选择，即选择所谓的永久散射系数。像素选择的关键参数是幅度分散(dispersion of amplitude，DA)，低DA值表明了“好”像素点，即具有低

对于每个干涉图，通过下面的式子生成相干图，计算相干性通常是使用n×m像素的窗，相干性为：

对于Q幅相干图，DA可以表达为：

其中，Stdev表示标准差(standard deviation)，A_k表示第k幅相干图图像像素(i,j)的幅度。MA(i,j)为给定像素(i,j)的均值幅度(mean amplitude,MA)

像素的选择，需要确定DA值门限，也就是需要对相位质量和可利用形变测量的密度进行权衡。

S3相位解缠和相位估计：

基于幅度进行一定的像素选择，然后采用最小费用流法在选择的像素点基础上进行相位解缠。对解缠后的干涉图进行相位估计，本发明采用一种以奇异值分解的最小二乘法为基础的相位估计算法。

(1)相位解缠。对Q幅相干图的缠绕的干涉相位依次进行解缠，也就是对式子(17)常数k的估计，在选择像素点基础上采用最小费用流算法进行相位解缠操作。

(2)相位估计采用的是基于奇异值分解的最小二乘法的相位估计算法。以式子(17)的观测方程为基础，其中k是相位解缠求出的整数。相位估计是假设第1幅图像的相位值为0，估计的相位是之后第i幅与第1幅相位的真实的差值。假如配准后有a，b，c，d四副图像，估计的相位可以是b相对a的干涉相位，也可以是c相对于a的干涉相位，或者是d相对于a的干涉相位。

假设第一步得到P幅配准图像中的第一幅的相位值为零，剩下的配准图像中第i幅的估计相位是是形变相位与大气相位屏的和为：

其中，i＝2,…,P-1；为：

式中，表示的第u幅图像和第s幅图像形成的干涉图相位解缠后的干涉相位。

第i幅图像的估计相位为可以写出i-1个系统方程，例如i＝4， 因而，对于每个像素点，可以写出Q个系统方程和P-1个未知数，Q是相位解缠后干涉图数，P是配准图像数，应用奇异值分解的最小二乘法来求解未知数，得到P-1幅配准图像每个像素点的估计相位。

相位估计整个推导过程中的关键参数：残差(residuals)，也就是输入的相位值和后验估计相位值：

其中，代表着最初的干涉相位的差值(解缠之前)。与残差点相联系的观测值(解缠之后的干涉相位)称为异常值。

算法的主要步骤如图7所示：

(1)应用奇异值分解的最小二乘法解系统方程，同时计算对应的残差点；

(2)分析残差点，找大于固定门限的残差点，选择绝对值大的残差点，找出异常值；

(3)从观察值中移除异常值，进行新的相位估计，重新计算和

(4)应用奇异值分解的最小二乘法解系统方程，计算异常值对应的残差点；

(5)分析异常值的旧残差点和和新残差点：校正、移除、保留。异常值的残差如果是2π的倍数，这个观测值可以校正和保留；通过对异常值加或减2π的倍数对异常值进行校正。

另外，对于异常值的保留还是移除，根据旧残差点和新残差点的对比来决定：新的残差点与旧的残差点如果都是2π的倍数，这个观测值就可以保留。

(5)新的相位估计，应用奇异值分解的最小二乘法解系统方程，同时，计算对应的残差点。重新开始步骤(2)，过程迭代是从(2)到(5)，直到没有剩余的异常值，也就是没有残差值大于门限。

S4大气相位屏APS(Atmospheric Phase Screen)估计：APS主要影响传感器到目标的距离，也在不同方式上影响方位向的目标。采用的方法是以现场获得气象数据和湿度函数的初始校正为基础，建立二维模型。

在大气相位屏估计过程中，使用的是最小二乘法估计。主要假设是APS是空间相关，APS的模型是变换程度g的二维多项式

其中，(i,j)是给定像素的坐标，是APS的相位，a_r,(g-r)是多项式参数，(i₀,j₀)是参考像素点的坐标。使用监测区域的地理信息，鉴别稳定区域，位于其中的像素点用来估计式(25)所表示的APS模型的参数。对于稳定区域内的像素点，其相位为零，则

APS估计算法包含下面的步骤：

(1)使用监测区域的地理信息，鉴别稳定区域，理想的监测情形是稳定区域包围着形变区域。

(2)使用式(25)在选择的稳定区域上估计APS。通过LS调整异常值来进行估计。估计的APS相位项如下，

其中，是估计参数，(i,j)是选择点的图像坐标，(i₀,j₀)是参考点的坐标。

(3)最后用式(26)估计形变区域的APS相位项。

S5形变估计：经过大气相位估计后，可以得到由形变引起的真实相位值，并进一步计算出监测区域的累积形变。

对所有的选择的像素点，从第三步得到的估计相位移除形变区域的APS来获得形变对应的相位雷达视线(Line Of Sight，LOS)方向的形变值为

在GB-SAR形变监测过程中，由于观测距离较近，通常为几十米到几公里之间。形变前雷达与目标之间的距离为R，雷达入射角，雷达LOS形变值ΔR_Los与形变真值ΔR_Defo的关系近似为：

利用余弦定理，LOS形变值与形变真值的关系为：

由式(30)和式(29)相减可得，两种计算所产生的误差为：

根据式(31)可以计算出与形变真值近似计算误差δ_R之间的关系，最后可以得到累积形变图。

Claims (2)

1.基于三重步进的GB-SAR形变监测方法，其特征在于，包括以下步骤：

数据采集步骤：雷达安装在线性导轨上以停-走-停方式步进运动，雷达天线沿方位向步进转动，雷达在进行每一次步进后停下来发射与接收一帧步进频率连续波信号；雷达工作在Ku波段；雷达天线沿方位向完成探测区域中全部位置的步进后，拼接雷达天线各方位向步进转动后获取的步进频率连续波信号完成一次成像数据采集；

成像处理步骤：对每次采集的成像数据在方位向上快速傅里叶变换后进行距离徒动校正、方位匹配与距离向聚焦，最后在方位向上进行快速傅里叶逆变换得到二维成像图像；

形变估计步骤：对二维成像图像顺序进行配准，再顺序冗余选择配准后的图像两两组成图像对；一组图像对包括一幅主图像与一幅副图像，一组图像对应生成一幅干涉图；对选择的干涉图的低噪声相位进行解缠，再对解缠后的干涉图进行相位估计，通过对所有解缠后干涉图的相位估计得到配准后的相对于第1幅图像的每幅二维成像图像的相位估计值；

根据原始相位估计值是监测对象形变相位值与大气相位屏之和，在相位估计值上移除估计的大气相位屏得到监测对象形变相位值；

进行相位估计的具体方法是：

1)基于奇异值分解的最小二乘法求解：

为第1幅二维成像图像配准后的相位估计值，为第i幅二维成像图像配准后的相位估计值，i＝2,…,P-1，P为二维成像图像总数，表示第i幅配准后的二维成像图像作为主图像和第i-1幅配准后的二维成像图像作为副图像形成的干涉图在解缠后的相位估计值；

同时计算各干涉图在解缠后的相位估计对应的残差点其中，为解缠之前主图像U和副图像S配准后形成的干涉图的干涉相位，和是通过对解缠之后干涉图的进行的相位估计得到主图像U和副图像S配准后的相位估计值；

2)将大于固定门限的残差点所对应的干涉图解缠后的相位估计值作为异常值，判断是存在异常值，如是进入步骤3)，如否将最新得到的相位估计结果作为最终的相位估计值；

3)当首次处理异常值，则直接移除异常值后返回步骤1)进行新的相位估计；

否则，分析异常值并对异常值进行校正、移除或保留处理后返回步骤1)进行新的相位估计；当异常值上一次计算的残差点与本次计算的残差点均是2π的倍数，则对该异常值进行校正或保留，否则移除该异常值；通过对异常值加或减2π的倍数对异常值进行校正。

2.如权利要求1所述基于三重步进的GB-SAR形变监测方法，其特征在于，雷达信号频率步进间隔由最大不模糊距离限制；轨道长度由方位分辨率和雷达波长限制；雷达在轨道上的步进运动间隔由波束宽度和波长限制；步进运动的时间间隔应由收发信号的最大时间限制；天线步进转动间隔由天线波束宽度限制。