[发明专利]一种宽带合成孔径雷达子带干涉数据处理方法

摘要

一种宽带合成孔径雷达子带干涉数据处理方法，该方法有六大步骤：步骤一：主副图像的自适应精密配准及副图像的重采样；步骤二：子带频谱位置和带宽的优化；步骤三：主副图像的带通滤波；步骤四：子带干涉图和子带差分干涉图的生成；步骤五：参考面和地形高程(DEM)的相位补偿；步骤六：子带差分干涉图的堆栈处理。本发明在大多数情况下无需相位解缠，即可获取绝对雷达干涉相位，从而可以通过转换直接得到数字地形图或者地表形变图，消除了由于相位解缠错误引起的误差扩散；它特别适用于高山地区的数字地形测绘和大形变速率的滑坡及矿山塌陷监测。

主权项

1.一种宽带合成孔径雷达子带干涉数据处理方法，其特征在于：该方法具体步骤如下：步骤一：主副图像的自适应精密配准及副图像的重采样由于子带干涉技术主要是应用在高山地区的高程测量或者是具有较大形变速率的滑坡和矿山塌陷测量，在主副图像进行配准时，不仅要考虑雷达干涉系统引起的地面同一目标在主副图像上的位置的偏移，还要考虑地面目标因较大的高程或较大的形变位移引起的主副图像上的位置偏移；在这里，副图像利用Sinc函数进行重采样，Sinc函数的位移量由下式决定：Offset(i,j)=(△i,△j)=((△i,△j)_i,j)_System+((△i,△j)_i,j)_Object      (3)上式中，offset(i,j)表示了副图像中位于第i行第j列的像元相对于主图像上同一像元位置的偏移量，它是系统偏移和目标本身因高程或形变引起的偏移量的和；offset(i,j)由副图像s中位于(i,j)的像元的邻域与相应的主图像p中像元的邻域的协方差函数取最大值时的位移量决定，见下面的(4)和(5)式cps(n,m)i,j=Σk=-N/2k=N/2Σl=-M/2l=M/2(p(k,l)-p~)·(s(k+n,l+m)-s~)---(4)]]>Δi=nmaxΔj=mmax---(5)]]>上两式中，p和s分别表示主副图像，c表示在主副图像中以第i行第j列为中心的N行M列的两个窗口的协方差函数，n和m分别表示行和列偏移量，k和l表示窗口中参与相关计算像元的行列位置，(△i，△j)表示当c取得最大值时的行列号，即（nmax，mmax）；步骤二：子带频谱位置和带宽的优化在理想的情况下，希望子带带宽尽可能的宽，这样能使雷达图像分辨率损失不要太大；另一方面，希望上下子带图像没有相同的频率成分，否则这些相同的频率成分会以噪声的形式出现在二次干涉图中；而根据信号处理理论，理想带通滤波器是不可实现的，也就是说，输入信号在带通滤波器的阻带内只能受到有限的衰减，因此，希望上下子带的位置相距越远越好；显然，两个理想状态是不可能同时得到满足的，一个折中的选择，是上下子带的带宽bw都等于雷达原始信号带宽B的三分之一，上下子带位于f_c±f₀的位置上，即bw=B3f0=B3---(6)]]>在这样的优化选择下，模拟波长与载波波长的比值为：λ′λ=3fc2B---(7)]]>一般情况下：f0=B-bw2λ′=fcB-bw·λ---(8)]]>在上面三式中，B是雷达脉冲带宽，fc表示载波频率，λ表示载波波长，bw表示上下子带带宽，f₀表示上下子带中心频率对于载波频率的偏移量，λ′表示子带差分干涉图的模拟波长；步骤三：主副图像的带通滤波主副雷达SLC图像通过上下子带带通滤波器滤波，生成各自的上下子带复数图像；上子带带通滤波器的归一化中心频率为0.333，归一化带宽为0.333；下子带带通滤波器归一化中心频率为0.666，即-0.333，归一化带宽为0.333；为了减少频谱能量的泄漏，带通滤波器的频谱需加窗，窗口函数选用凯撒窗即Kaiser,其频谱表达式如下：W(f;β,F)=I0(β1-(2f/F)2)I0(β)-F2<f<F2---(9)]]>上式中，F为通带带宽，I₀是零阶贝塞尔函数，β是控制窗口通带宽和阻带衰减折中关系的一个常数，通常取3～6之间；步骤四：子带干涉图和子带差分干涉图的生成主图像的上下子带图与对应的副图像的上下子带图进行干涉，产生上子带干涉图和下子带干涉图；将上下子带干涉图再次进行干涉，生成子带差分干涉图；记P和S分别为主副图像，PS为它们的子带差分干涉图，则PS=(Pup·Sup*)·(Plow·Slow*)*---(10)]]>上式中，P_up，P_low，S_up，和S_low分别表示由主副图像分裂出来的上下子带图，PS表示由这四副子带图像构成的子带差分干涉图；步骤五：参考面和地形高程(DEM)的相位补偿地球的WGS84椭球体模型作为干涉测量的参考面，如果要进行差分干涉测量，地球的WGS84椭球体表面还需叠加所测区域的DEM，参考面和地形高程(DEM)引起的模拟干涉相位要从子带差分干涉图中去除，为此，首先需计算主副图像中各像元与雷达的距离R，它满足下面的非线性方程组：R=(Rs→-R→t)·(Rs→-R→t)2λR(Vs→-V→t)·(Rs→-R→t)=0xt2(a+h)2+yt2(a+h)2+zt2b2=1---(11)]]>上式中，R_s为雷达位置向量，R_t为像元位置向量，即(x_t,y_t,z_t)，V_s和V_t分别为雷达和像元目标运动速度，h为像元高程，若仅补偿参考面相位，取h=0，a和b是地球椭球体模型的长短轴；计算出每一个像元分别在主副图像中与雷达的距离R以后，通过像元在主副中的距离的差△R，便可以计算出子带差分干涉图中的参考面和地形高程的模拟相位：φ=-4πλ′·ΔR---(12)]]>上式中的波长应为子带干涉后的模拟波长，利用(1)式和(2)式所表达的高程和形变模糊度便可将子带差分干涉图转化为地表高程或地表形变场；步骤六：子带差分干涉图的堆栈处理a)数字高程图即DEM的生成子带差分干涉图的高程模糊度虽然被提高了很多倍，但同时也将干涉相位的误差的方差亦扩大了几近相同的倍数，单张子带差分干涉图的相位会含有更大的噪声；因为干涉图的相位通过相位模糊度的因子正比于高程，故可以将多张子带差分干涉图组成一个堆栈，然后用最小二乘法估计每一个像元的高程：MinhΣ(h-Δφi2πH2π,i)2---(13)]]>上式中，△φi表示像元i在堆栈中某一干涉图里的相位，H_2π,i表示该点的高程模糊度，h表示该点高程的估计值；b)地表形变场的测量由于子带差分干涉图的形变位移模糊度很大，故由差分相位转换为形变量时几乎不用进行相位解缠，但为了消除大气影响、高程误差及其它的各种噪声，还是有必要利用一组子带差分干涉图的时序分析，特别是短基线集法，来求得精确的地表形变场及其动态变化；子带差分干涉图的时序分析方法同传统的全带宽数据干涉时序分析法基本是相同的，唯一不同的是高相干点的选取方法；这里利用了两个判据来确定候选的高相干点：σμ<cSCR>d---(14)]]>上式中，μ和σ是像元在堆栈时间范围内强度值的均值和方差，c是对方差与均值的比值设定的一个阈值，它基本上能衡量该像元在堆栈时间范围内幅值的稳定与否，一般取值为0.2左右；SCR是像元的信杂比即Signal to Clutter Ratio，在子带差分干涉图中，进行时序分析的像元的信杂比必须大于一个阈值d,一般取值25dB；选取了候选点后，利用一般的周期图法或者短基线集法获得地表形变场的速率，线性与非线性形变量；如果借助人工反射器，大尺度的突变和非线性形变的测量都不会受到限制。