本发明涉及雷达遥感技术领域，尤其涉及一种基于最优参数设置的子带干涉测量的地震形变获取方法，主要应用于测绘、地震行业。
背景技术：
现有传统子带干涉测量技术是一项用于特殊地形测绘的技术，相比于常规雷达干涉测量技术(英文缩写为insar技术)，传统子带干涉测量技术可以在数据处理中避免容易带来误差的相位解缠步骤，从而直接获取到高程梯度较大地区的数字高程模型dem。但是，由于子带干涉技术在子带分割时缩减了带宽，导致了噪声的增大和旁瓣带来的额外干扰，使得结果干涉图质量严重下降或根本无法使用，从而很大程度地影响了测量精度，制约了传统子带干涉测量技术在要求较高的地震形变监测中的应用。
技术实现要素：
本发明的目的在于提供一种基于最优参数设置的子带干涉测量的地震形变获取方法，用于抑制和去除传统子带干涉测量技术中由于子带分割、带宽缩减因素引入的一系列噪声，从而提高地形和形变观测精度，确保子带干涉测量技术在地震形变监测领域的成功应用。
为了实现上述目的，本发明所述基于最优参数设置的子带干涉测量的地震形变获取方法，所述方法包括：
s1，待检测区的主slc、辅slc的精密配准并进行子带分割：其中，待检测区的主slc、辅slc的精密配准具体为：对待检测区的主slc和辅slc进行精配准，并得到主slc和辅slc最高配准度时对应的偏移量多项式；子带分割：设子带带宽为b，待分割slc的全带宽为b，令b＝0.2b，设子带中心频率为f子带中心频率，待分割slc的中心频率为fc，子带分割偏移量为f0，令f子带中心频率＝fc+f0；f0∈[0.3b,0.4b]，在b和f子带中心频率的条件下，分别对主slc、辅slc进行子带分割，得到的上子带主slc和上子带辅slc经过所述偏移量多项式配准后得上子带配准结果，得到的下子带主slc和下子带辅slc经过所述偏移量多项式配准后得下子带配准结果；
s2，进行二次干涉：对上子带配准结果进行干涉，得到上子带干涉相位图；对下子带配准结果进行干涉，得到下子带干涉相位图；再对上子带干涉相位图和下子带干涉相位图进行第二次干涉，去除地形相位保留形变相位，得到子带干涉相位图；
s3，子带干涉相位图的平滑降噪处理：在s2结果的基础上，首先在子带干涉相位图的频谱上精确估计线性条纹频率；然后，以所述线性条纹频率构建极窄陷波器，对频率在[nf-1,nf+1]范围内的频段信号进行过滤，得到平滑降噪处理后的相位图，nf表示所述线性条纹频率；
s4，对平滑降噪处理后的相位图中的子带干涉条纹进行有效性检测：使用insar处理待检测区主slc、辅slc得到干涉图q，利用所述干涉图q获得干涉条纹的有效性范围；在干涉条纹有效性范围的基础上，对平滑降噪处理后的相位图进行掩膜，得到经过有效性检测的相位图；
s5，将经过有效性检测的相位图进行地理编码、地形、形变转换得到待检测区地震后的形变图。
优选地，子带分割过程和二次干涉过程，具体为：s11，设子带带宽为b，待分割slc的全带宽为b，令b＝0.2b，设子带中心频率为f子带中心频率，待分割slc的中心频率为fc，子带分割偏移量为f0，令f子带中心频率＝fc+f0；f0∈[0.3b,0.4b]；s12，在b和f子带中心频率的条件下，分别对主slc、辅slc进行子带分割，得到上子带主slc、下子带主slc、上子带辅slc和下子带辅slc；s13，采用偏移量多项式对上子带辅slc进行配准和重采样，得到重采样后的上子带辅slc，同时，采用偏移量多项式对下子带辅slc进行配准和重采样，得到重采样后的下子带辅slc；s14，将上子带主slc与重采样后上子带辅slc进行干涉，得到上子带干涉相位图，同时将下子带主slc与重采样后下子带辅slc进行干涉，得到下子带干涉相位；s15，将上子带干涉相位与下子带干涉相位进行干涉，得到子带干涉相位。
优选地，s3中，干涉相位图的平滑降噪处理，具体按照下述实现：
s31，利用chirp-z变换精确估计条纹频率；首先，在子带干涉相位图中，沿距离向截取一段长度与距离向相等的区域，且该区域的相位梯度变化平缓、受地形或形变相位干扰少；然后，估计子带干涉相位图中线性趋势条纹的条纹数，并计算该条纹数对应的频率nf；最后，将所述区域所对应的数据在[nf-1,nf+1]频段进行chirp-z变换，获得对该频段加密采样的频谱，找到[nf-1,nf+1]频段出现的振幅峰值，所述振幅峰值所在位置即为精确到小数的条纹数，并将精确到小数的条纹数对应的频率作为条纹频率点ns；
s32，根据条纹频率点ns构建极窄陷波器，对处于[ns-1,ns+1]频段的信号进行过滤，过滤方式是对子带干涉相位图进行逐行过滤；首先，设子带干涉相位图单行的长度为p，在该长度p对应的信号x(n)，n＝1,2,…p后补入一个长度等于延迟量d的均为0值的序列，得到长度为p’＝p+d的新输入信号x’(m)，m＝1,2,…p’；然后，采用构建的极窄陷波器对x’(m)进行滤波处理，得到滤波结果y’(m)；最后，在y’(m)中忽略m＝1,2,…d处的相位延迟，并提取m＝d+1,…p’处的滤波结果，完成干涉相位图的平滑降噪处理。
更优选地，所述极窄陷波器为带阻fir滤波器。
优选地，s4中，使用insar处理待检测区主slc、辅slc得到干涉图q，利用所述干涉图q获得干涉条纹的有效性范围；具体为：
s41，使用insar处理待检测区主slc、辅slc得到的干涉图q，在干涉图q设置尺寸为n×n滑动窗口，并采用公式(1)对任意一个滑动窗口在干涉图q中i位置的窗口内干涉条纹进行二维傅立叶变换；
式中，f(u,v)为在第u行和第v列的频率，u＝0,1,…,n-1；v＝0,1,…,n-1；x和y分别为干涉图q的行号和列号；
s42，由于滑动窗口的尺寸通常远小于干涉图的图幅大小，因此将滑动窗口内的干涉条纹视为线性均匀变化，采用公式(2)计算滑动窗口在干涉图q中i位置的振幅峰值与采样点数比r；
式中，max(|f(u,v)|)表示滑动窗口i内的最大值；
s43，估计阈值rt；s431，在给定的尺寸为n×n的窗口中模拟出不含噪声的干涉相位，并根据公式(2)设置变量rf和rmax，rf表示已知条纹频率所处的振幅与采样点数的比，rmax表示窗口内的振幅峰值与采样点数的比；s432，向窗口内的干涉相位加入多组随机噪声，噪声强度逐步增加；绘制rf和rmax分别对应噪声水平的衰减曲线，噪声作为横轴，r值作为纵轴，获得rf曲线与rmax曲线出现分叉且再未重合时分叉处的r值，将该r值作为rc；s433，重复步骤s431-s432，求出多个rc的算术平均值rt，将rt作为实际操作的阈值，即阈值rt；
s44，将公式(2)计算得到的r与rt对比，如果r≥rt，则滑动窗口i内存在干涉条纹，干涉图q中i位置的相位有效；反之当r<rt，则滑动窗口i内不存在干涉条纹或已被噪声淹没无法被识别，干涉图q中i位置的相位无效；
s45，移动滑动窗口至干涉图q中除i位置外的其他位置，重复s41-s44得到干涉图q的有效相位，根据有效相位进一步获取干涉图q的r分布，最后应用阈值rt获得干涉条纹有效性范围。
本发明的有益效果是：
1、本申请所述方法在子带分割前，设置子带分割的最优参数子带带宽b＝0.2b和子带中心频率f子带中心频率＝fc+f0；f0∈[0.3b,0.4b]，最大程度减少上下子带频谱重叠带来的噪声以及旁瓣效应。2、本申请所述方法增加了干涉像对的常规insar精配准步骤，以获得精确的偏移量多项式，供后续子带干涉对配准处理中使用，避免由分辨率损失带来的配准误差或由频谱变化而引起的线性条纹误差。3、本申请所述方法公开了构建窄带陷波器并配合chirp-z变换的抑制方法，以解决因旁瓣能量泄漏出现的误差。4、本申请所述方法公开了基于采用常规insar的干涉图进行有效性干涉条纹检测，掩膜出子带干涉结果中有效的地形和形变数据区域，以解决利用相干性作为指标掩膜时保留无效区域数据的问题。
附图说明
图1是子带中心频率和带宽配置示意图；图2是图1中方案1～10的均方根误差随b/b变化的曲线示意图；图3是madsen流程与改进流程的结果对比示意图；图3中a图为使用madsen流程获得的子带干涉相位结果示意图；图3中b图为使用改进流程获得的子带干涉相位结果示意图；图4是改进流程与brcic流程的结果对比示意图；图4中a图为使用改进流程获得的子带干涉相位结果示意图；图4中b图为使用brcic流程获得的子带干涉相位结果示意图；图5是本申请中构建窄带陷波器并配合chirp-z的抑制方法对于趋势噪声的抑制效果比较示意图，图5中的a图为子带干涉结果示意图，图5中的b图为利用了趋势噪声抑制方法后的干涉结果示意图；图5中的c图为a图和b图中所示的剖面线示意图；图6是基于条纹频率的子带干涉有效性检测示意图；图6中a图为常规insar差分干涉相位图，图6中b图为a图的相干性检测图，图6中c图为a图的条纹频率检测图，图6中d图为本申请所述方法的子带干涉相位图；图7是基于最优参数设置的子带干涉处理流程示意图；图8是sar数据覆盖范围示意图，sar数据覆盖范围即图中黑色小矩形框，其中背景图是利用landsat方法获取形变场，并且已经转换到地距向；图9是常规insar配准结果对比示意图；图9中a图为滑动窗口为512×512(像素)的配准结果，图9中b图为滑动窗口为64×64(像素)的配准结果；图10是常规insar相位解缠图示意图；图10中a图为从左下端进行相位解缠的结果，图10中b图为从右上端进行相位解缠的结果；图11是常规dinsar地距向形变场示意图；图12是不发生相位缠绕的子带中心频率和带宽配置方案；图13是降噪滤波前后的子带干涉相位示意图；图13中a图为降噪滤波前的子带干涉相位，图13中b图为降噪滤波后的子带干涉相位；图14是子带干涉地表形变场示意图；图15是3种手段获得的同震形变场示意图；图15中a图为模拟得到形变场的示意图，图15中b图为本申请所述方法得到形变场的示意图，图15中c图为landsat-8方法得到形变场的示意图，图15中d图为常规insar方法得到形变场的示意图；图16是4种手段获得形变场的噪声水平示意图；图17是3种手段获得形变场的残差分布示意图；图17中a图为本申请所述方法模拟获得的形变场的残差分布示意图，图17中b图为landsat-8方法模拟获得的形变场的残差分布示意图，图17中c图为常规insar方法模拟获得的形变场的残差分布示意图；图18是rmax和rf随噪声水平变化曲线示意图；图19是形变剖面示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施方式仅仅用以解释本发明，并不用于限定本发明。
为了便于理解本申请，将英文缩写的英文全称和中文全称进行解释说明：sar，syntheticapertureradar，合成孔径雷达。slc，singlelookcomplex，单视复数影像。
关于本申请的几点说明：
如本领域技术人员公知：理想状况下的子带干涉测量应用中，如果sar传感器的轨道和入射角等参数保持一致，利用子带分割出的具备不同中心频率和带宽的子带slc，对同一地物进行干涉测量应获得相同的测量结果。但，在实际的处理过程中往往会受到分辨率降低、二次干涉噪声放大、带通滤波器旁瓣能量泄漏等因素的干扰，给最终结果造成或轻或重的误差影响。为了抑制和去除传统子带干涉技术中由于子带分割、带宽缩减等因素引入的一系列噪声，本申请在几个关键技术环节上提出的处理技术大大提高了地形和形变观测精度，确保把子带干涉技术成功扩展应用到地震形变监测领域，详细说明：
ⅰ、子带干涉的最优化参数选取
子带干涉的本质是通过子带分割(缩减带宽)和两次干涉的方法以增长有效波长，从而减少或消除干涉条纹数，从而能够应用于陡峭地形和形变测量。在子带分割时，有两个参数至关重要：子带中心频率和子带带宽。子带的中心频率设置会影响模拟得到的子带波长大小，而子带带宽的设置则会影响分辨率和噪声水平，进而影响精度水平。因此本申请在子带干涉处理中，合理地设置子带中心频率和子带带宽，以获得更高的测量精度。
基于理论公式，研究学者bamler给出了最优带宽为b/b＝1/3。但在真实的子带干涉测量中，最优带宽除了会受sar影像自身质量、测量地物条件等不可变因素的影响外，最大影响来自带通滤波器的性能，以及子带中心频率的配置，因此子带带宽的选取是否依然以b/b＝1/3为最优，仍需要进行进一步的研究。为了选择出最优子带中心频率和子带带宽，本申请设计以下实验方案(如图1)，以测试不同子带带宽和中心频率对测量精度的影响，以0.05b为间隔，从0.1b开始逐步增大子带带宽，同时中心频率跟着变化，均处于子带带宽中心，fc为原slc的中心频率。
以astergdem数字高程模型作为先验参考，利用逐步参数选取法，进行不同参数下的子带干涉处理，并以干涉图的相干性、高程残差和均方根误差等指标来对各种参数方案配置产生的结果进行评价，从而实现最优化参数的选取。研究发现：子带带宽为0.2b，子带中心频率偏移量f0为0.4b，即中心频率为fc±0.4b时，可认为是最优参数配置方案；当需要进一步增大子带波长时，可维持0.2b的子带带宽，并设置子带中心频率偏移量f0在[0.3,0.4)的区间内移动，也可获得误差较小的结果。
ⅱ流程改进
传统子带干涉技术的处理流程有两种，madsen和brcic处理流程，二者的主要区别集中在子带配准和子带分割的先后次序上。madsen和brcic的处理流程，虽然都能够对slc像对进行子带干涉处理，得到干涉结果，但是均存在一些不足，致使最后结果容易受额外噪声的干扰影响。
madsen是将聚焦后的slc，先做子带分割，然后将上子带像对和下子带像对分别做配准和干涉。采用madsen的处理流程，不但会导致上下子带像对自身出现较大的偏移量误差，也会使上下子带干涉相位之间产生不一致的相位错动。这样的配准误差和相位错动会在二次干涉时对结果造成重大影响，整幅子带干涉相位被一个沿距离向从左至右递减的强趋势面所覆盖，其地形和形变相位值被严重干扰，基本无法分辨。
在brcic的流程中，首先对原主辅slc进行粗配准，然后再做子带分割，这样能够很好地解决上下子带像对单独配准出现的偏移量不一致的问题。但是重采样后的辅slc在频谱上会发生改变，这样的改变不但是幅值在大部分频段上小幅衰减，还有可能使某些频点发生不易观察到的偏移。偏移的频点可能会在最终的子带干涉相位中引入一些周期性的全影像范围覆盖的弱线性条纹误差，给最终结果带来额外的误差影响。
针对上述两种处理流程存在的不足，申请人在配准和子带分割环节提出了改进的流程思路，图7给出了完整的子带干涉处理流程。图7中，本申请增加了干涉像对的常规insar精配准步骤，目的是为了获得精确的偏移量多项式，供后续子带干涉对配准处理中使用。
由上面分析可知，一方面，如果首先进行配准，再进行距离向频谱分割，则会在最终子带干涉相位上出现了一些周期性的线性条纹，这是因为重采样后的辅slc其频谱会产生一定的改变，有可能在某个频段上引入不必要的误差造成的；另一方面，如果对原主辅slc进行子带分割后单独配准，又会因分辨率下降出现配准上的较大误差，以及上下子带干涉相位之间出现一定的错动，最终在二次干涉中出现较大误差。而理想状态下，所有上下子带影像均由同一主辅slc像对分割出来，上子带像对和下子带像对自身的偏移量应与原主辅slc的偏移量相同，且上、下子带干涉之间也不应出现偏移。基于以上分析，本申请提出的改进配准流程是：借助配准精度最高的原主辅slc的偏移量多项式，直接对上子带辅slc和下子带辅slc进行配准和重采样，从而最大化地保证配准精度，避免分辨率损失带来的配准误差和brcic流程中在子带干涉相位中容易出现的线性条纹。
ⅲ子带干涉趋势噪声的抑制方法
子带干涉法获得的最终相位，是由上下子带干涉图进行二次干涉后产生的，而在最终相位中，除了各种噪声会因为二次干涉而被放大，一些在首次干涉中出现的误差也容易因为二次干涉而在最终结果中显露出来。这些在二次干涉中容易突显出的误差包括在处理流程ⅱ中讨论过的配准误差和频点偏移误差，以及因旁瓣能量泄漏有可能出现的误差。配准误差和频点偏移误差由处理流程ⅱ尽可能避免；而旁瓣能量泄漏由带通滤波器的参数设计以及子带中心频率和带宽决定的，并不能在处理流程ⅱ中完全消除。当带通滤波器抑制旁瓣的能力较弱，或上下子带中心频率和带宽过于靠近时，这一误差就有可能显现并影响最终的干涉结果。本申请结合前人的研究成果和子带干涉的特点，提出构建极窄陷波器抑制因旁瓣能量泄漏引起的线性趋势条纹的方法。构建极窄陷波器抑制线性趋势条纹的方法分为2个步骤：首先是在频谱上对线性条纹频率进行精确的估计；然后以该条纹频率点构建极窄的陷波器，对该频点和附近小范围内的频段信号进行过滤。
首先利用chirp-z变换精确估计条纹频率，步骤如下：1.在子带干涉相位图中，沿距离向截取一段相位梯度变化较平缓、受地形或形变相位干扰较少、长度与距离向相等的区域，用以估计条纹频率。2.一般线性趋势条纹在整幅干涉影像上是较为清晰可分辨的，通过目视观察的方法，可粗略估计线性趋势条纹的条纹数nf。3.对该段数据在nf位置及附近的频段(如：nf±1)进行chirp-z变换，获得对该频段加密采样的频谱，找到nf附近出现的振幅峰值，则峰值所在位置即为精确到小数的条纹数。
然后，根据条纹频率点构建极窄的陷波器，对该频点和附近小范围内的频段信号进行过滤。所谓带阻滤波器(bsf)是指能通过大多数频率分量，但将某些范围的频率分量衰减到极低水平的滤波器，本申请要做的就是尽量衰减频率在nf处的能量以达到抑制线性趋势条纹的目的。理想状态下，希望只衰减nf和附近处的频率分量，而尽量保留其他频率分量，这就意味着过渡带是十分陡峭的。公知：使用基于凯塞尔窗的fir滤波器会带来极高的相位延迟。而，本申请采用信号补零的方法可消除相位延迟，具体做法如下：首先，在长度为p的原信号x(n)，n＝1,2,…p后补入一个长度等于延迟量d的0值序列，得到长度为p’＝p+d的新输入信号x’(m)，m＝1,2,…p’；然后，采用构建好的带阻fir滤波器对x’(m)进行滤波处理，得到滤波结果y’(m)；最后，在y’(m)中忽略m＝1,2,…d处的相位延迟，并提取m＝d+1,…p’处的滤波结果。
ⅳ子带干涉条纹的有效性检测方法
在常规insar中通常使用相干性作为衡量干涉图质量和有效性的指标，当相干性低于某个阈值时，会将该相干性对应的区域掩膜去除，不作为最终结果中的有效数据使用。子带干涉相位需要经过多次迭代降噪以消除其中大量的噪声，而一些原本相干性很低的区域，经过多次迭代降噪后，也会受周围较高相干区的影响出现均匀化连续化的现象。这时如果还利用相干性作为指标进行掩膜，容易把无效区域的数据保留下来。因此，本申请采用一种条纹检测的机制，用以判断子带干涉结果在局部地区的有效性。
经过分析发现只要在常规insar中出现干涉条纹的地区，在子带干涉处理中，生成的上下子带干涉相位内该地区仍然保留着相似的条纹分布，在二次干涉后，这一区域的干涉相位仍是较准确可用；而在常规insar中表现为噪声的地区，子带干涉处理后，该地区也仍然为噪声。且由于常规insar中的条纹分辨率更高，不受子带干涉处理带来的额外干扰。因此，本申请提出基于常规insar干涉条纹进行的有效性条纹检测算法，其基本原理是在常规insar干涉条纹图中设置一定尺寸的滑动窗口，寻找窗口内的干涉条纹在频率域中的振幅峰值，并计算该振幅峰值和采样点数的比r，并结合一定的阈值rt来判断该窗口内的条纹是否存在；其中，阈值来源于多次加噪声的模拟干涉实验的平均，当r高于rt时，可判断条纹存在，r越高说明条纹越清晰(条纹的周期性越好)，受噪声干扰的程度越低，反之当r低于rt时，说明条纹可能已经被噪声淹没，或该窗口内可能并不存在干涉条纹。具体步骤如下：
1)对滑动窗口内的干涉图进行二维傅立叶变换：首先在干涉条纹图中设置尺寸为n*n的滑动窗口，并将窗口内的局部干涉条纹进行二维傅立叶变换。公式如下：
式中，u,v＝0,1,…,n-1，x和y分别为干涉图q的行号和列号。
2)计算窗口内的峰值与采样点数比r：由于滑动窗口的尺寸通常远小于干涉图的图幅大小，因此可将窗口内的干涉条纹视为线性均匀变化的。对f(u,v)的振幅取最大值(峰值)，计算该峰值与采样点数的比r，公式如下：
式中max(·)表示取最大值。
3)估计阈值rt：通过多组模拟测试实验给出了常用窗口尺寸下对应rt的平均值，作为实际操作的参考依据。
s431，在给定的尺寸为n×n的窗口中模拟出不含噪声的干涉相位，并根据公式(2)设置变量rf和rmax，rf表示已知条纹频率所处的振幅与采样点数的比，rmax表示窗口内的振幅峰值与采样点数的比；
s432，向窗口内的干涉相位加入多组随机噪声，噪声强度逐步增加；绘制rf和rmax分别对应噪声水平的衰减曲线，噪声作为横轴，r值作为纵轴，获得rf曲线与rmax曲线出现分叉且再未重合时分叉处的r值，将该r值作为rc；
s433，重复步骤s431-s432，求出多个rc的算术平均值rt，将rt作为实际操作的阈值，即阈值rt；
如图18所示，图中rmax和rf两条r曲线分别代表振幅图上的峰值和已知条纹频率所在的幅值，均随噪声增大而持续递减。在噪声水平较小的区间，两条曲线完全重合，说明此时的振幅图上的峰值和已知条纹频率所在的幅值相等。当两条r曲线持续随噪声增大而衰减，会在r为某个阈值rc时出现分叉并再未重合，并且rmax始终大于rf，这说明线性条纹随着噪声逐步增大，在r<rc以后已经完全被噪声淹没，条纹无法识别。
4)应用阈值判断子带干涉相位的有效性：将r与rt做比较，如果r≥rt，则窗口内必然存在干涉条纹，子带干涉相位是有效的；反之当r<rt，则可认为窗口内不存在干涉条纹或已被噪声淹没无法被识别，子带干涉获得的相位也是无效的。
5)移动窗口位置，重复步骤1)～4)计算整幅干涉图的功率谱密度分布，并应用阈值rt获得条纹有效性分布图。
表1rt与窗口尺寸对照表
对本申请提出的四种处理的效果评价：
ⅰ、子带干涉的最优化参数选取结果评价
采用了均方根误差(rootmeansquareerror,rmse)来检测图1中方案1-10计算得到的高程值相对于astergdem高程值的偏离程度，方案中的rmse越小说明该方案的计算精度越高。均方根误差σrmse的计算公式如下：
其中m和n分别为影像的行数和列数，hi是方案1-10计算得到的二维高程图转换成的一维序列，hi是对应的astergdem高程图转换得到的一维序列。由方案1-10计算得到的rmse分别为{47.5844,29.6204,23.0572,24.3553,23.3479,25.7188,23.8768,28.9921,32.4175,37.7373}(米)，与各自对应的子带带宽可绘制出如图2所示的关系曲线。分析图2可知，rmse随子带带宽的增大呈现先下降后上升的趋势。其中最小值出现在0.2b的位置，对应方案3。而在此之前方案1和方案2的rmse均较大，是因为其对应的子带带宽0.1b和0.15b过小，分辨率损失过大所带来的图像模糊造成的，而旁瓣能量泄漏此时造成的影响则非常轻微，基本可以忽略。当子带带宽大于0.2b且继续增大时，对应的rmse在0.25b～0.35b处出现缓慢抬升，而在>0.35b处则出现较快抬升。分析原因，虽然子带带宽的不断增大会对分辨率以及干涉质量(相干性)有所提升，但是也会使谱间干扰问题不断加重，给测量结果带来的负面影响甚至超过分辨率提升带来的增益。
综上，子带带宽处于0.2b～0.3b，且上下上子带带宽相隔最大时对应的rmse均较小，线性趋势条纹也可忽略，是较为理想的带宽选择区间；而bamler等提及的子带带宽在1/3b处，对应的rmse已有所增大，并且在残差分布图中线性趋势条纹也基本清晰可识别，在条纹频率图上峰值也已经出现，说明谱间干扰的负面影响已经开始占据主导，因此该方案并不是最优的带宽配置方案；当子带带宽持续增大，线性趋势条纹带来的干扰会变得十分严重，不应作为子带带宽选择的区间。
ⅱ由技术流程改进产生的效果评价
为了评价改进技术流程的效果，申请人通过实验对改进流程和两种传统流程(madsen和brcic)的结果进行了比较。首先，选取覆盖美国大峡谷部分地区的一对雷达影像进行子带干涉处理，分别采用madsen流程和改进流程进行数据处理来提取美国大峡谷部分地区的dem，图3为利用两种流程提取的dem干涉相位结果，图3amadsen的结果中明显存在一个强趋势相位面，而用本申请提出的改进处理流程提取的dem相位中不存在上下子带干涉相位之间的偏差带来的强趋势相位覆盖，因此地形清晰可识别。
然后，选取覆盖澳大利亚乌鲁汝/艾尔斯岩地区的雷达影像来测试brcic处理流程和本申请提出的改进流程的处理结果，分别用这两种流程进行数据处理来提取该地区的dem，图4为利用这两种流程提取的dem干涉相位结果，图4a为改进流程的处理结果，地形清晰可见，与图4bbrcic流程的处理结果相比，不存在弱线性趋势条纹的影响，相位能更真实地反映地表起伏情况。
ⅲ子带干涉趋势噪声的抑制效果
仍然用覆盖澳大利亚乌鲁汝/艾尔斯岩地区的实际雷达影像数据来测试本申请提取的噪声抑制方法对线性趋势误差的抑制效果。图5a是未使用带阻滤波的子带干涉相位结果，图5b是使用构建的窄带带阻滤波器(bsf)+chirp-z的子带干涉结果，可以明显的看出未进行处理之前的相位图带有条带状背景噪声，而进行了bsf+chirp-z滤波处理后，相位图中这一背景噪声基本消失。图5c为在相同剖面处滤波前的干涉相位剖面，蓝线为滤波前干涉相位，红线为使用本文提出的窄带滤波器处理后的干涉相位，可看到滤波后的干涉相位在趋势边界左右的起伏比抑制前平缓，说明条纹残差得到了一定程度的抑制，线性趋势效应有一定程度上的减轻，但并不能完全消除。
ⅳ子带干涉条纹的有效性检测效果
为了评价本申请提出的基于条纹频率的干涉条纹有效性检测算法的效果，申请人利用2013年巴基斯坦俾路支省(balochistan)mw7.7级地震的真实雷达数据进行测试。图6是本申请提出的条纹频率检测效果比较结果图，图6a是常规insar差分干涉相位图，图6b、c分别是图6a的相干性检测图、条纹频率检测图，对比b、c图可以发现，在分布面积和分布区域上，有效区域(非白色区域)的分布情况大致相同，但是条纹频率检测结果的有效性覆盖区域比相干性检测有了很大的扩增，在有效点数上，本申请提出的方法估计出的点数为4738390个，相干性法估计的点数为4176333个，说明本申请的方法对未降噪的干涉图进行条纹有效识别的能力，与相干性法对经过小幅滤波后的干涉图进行估计的能力相当，本文方法对低相干和噪声的排除能力更高，能够更有效地识别有效的数据区域。
实施例
在上述理论算法和参数选取方法研究的基础上，选取2013年巴基斯坦俾路支省(balochistan)mw7.7级地震作为实例，在国际上首次利用子带干涉技术和最优参数配置，成功提取了该地震同震形变场。通过与常规d-insar的结果进行对比，验证了子带干涉法可有效避免前者在条纹密集区域的解缠错误。通过和offset-tracking、landsat8和模型拟合的结果进行对比，子带干涉法能获得更高测量精度和更丰富的细节展现。具体细节如下：
1震区概况与数据选取
2013年9月24日，在巴基斯坦俾路支省境内的阿瓦兰县北北东66km处发生了一次震级为mw7.7的地震。为了对获取该次地震引起的地表形变，收集相关的数据，用于后面的子带干涉、常规insar、offset-tracking等技术的实验。
外部参考数据方面，选用地面分辨率最接近sar影像的astergdem第二版dem数据，用于生成地形相位。astergdem资料显示，破裂迹线北部山岭分布十分密集，多数山峰在800～1000米以上，地形起伏相较于南部更为复杂。在地表形变场参考上，选用landsat8光学卫星的震前和震后间隔16天的影像获取的地表形变场，以及根据landat8形变场模拟得到的拟形变场。
sar数据方面，表2给出了用于提取此次同震地表形变场的sar影像对的影像参数，影像覆盖了阿瓦兰县部分区域，图8。
表2阿瓦兰地区terrasar-xstripmap影像参数
2常规insar获取地震形变场
首先，采用常规“二轨法”insar技术提取地震的地表形变场，使用专业insar处理软件gamma完成整个处理流程，其基本步骤如下：
1)获取主辅slc影像：所使用的terrasar-x卫星已经取消了提供原始的回波信号数据，即通常所指的0级数据，而是直接提供已经进行成像处理后的slc影像。
2)主辅slc配准：gamma通常使用的配准方法为基于强度信息的最优互相关算法(cross-correlationoptimizationprocedure)。通过多次反复实验发现：当采用进行2次逐步配准的策略，设置滑动窗口分别为512×512(像素)和256×256时，破裂迹线南侧能获得较清晰的干涉条纹，但北侧则完全失相干，图9a，而当配准次数进一步增多，滑动窗口进一步缩小时(滑动窗口分别为：512×512、256×256、128×128、64×64)，发现破裂迹线北侧在原来完全失相干的西北角出现了一些的条纹，图9b，而迹线南侧的干涉条纹仍然存在，说明这一配准方案相较于前者在质量所有了一定的提升。故，将使用这一方案进行后续的干涉处理，并保留该次估计的偏移量多项式，供子带干涉时使用。
3)干涉、去平地相位和去地形相位：由于terrasar-x轨道精度较高，因此可以直接通过主辅像对的轨道参数，而无需借助于条纹频率来估计基线。对步骤2)中产生的两个配准结果，分别进行干涉。并利用基线参数和astergdem模拟出包含地形相位和平地相位的干涉图，与上述两张干涉图进行差分，获得只保留形变相位的子带干涉相位图。
4)降噪滤波：对步骤3)中生成的形变子带干涉相位图进行基于goldstein的自适应降噪滤波处理。由于形变条纹较为密集，为了能尽可能保持条纹的连续性，傅立叶搜索窗口的长度和α值都不应设置过大，将傅立叶搜索窗口的长度和α值分别设为64和0.4。
5)相位解缠：由于受到两条失相干带的阻隔，并且破裂迹线两侧的形变方向不一致，因此选择mcf在南北两端远离形变中心的位置分别进行相位解缠，并设置相干性阈值为0.4，获得如图10a和图10b所示的相位解缠结果。
6)计算形变值和地理编码：首先将图10a和图10b的绝对相位差分别进行地距向形变值的计算；然后将破裂迹线南北两侧的形变场分别进行地理编码，从雷达坐标系转为地理坐标系以便与其他结果进行比较分析。最后沿图11中的虚线，将上述两张地理坐标系下的形变场进行裁剪和合并，得到形变场。
从图11中可以看到，受失相干带的影响，获得的形变场相对于landsat8(图8)有很大的缺失。在破裂迹线南侧，形变值总体为正，北侧形变值总体为负，说明形变方向总体与landsat8观测结果相吻合。
3子带干涉处理及结果分析
3.1最优参数选取。terrasar-x影像覆盖范围内的地距向形变量为–5.64～2.88米，对应的斜距向形变约为–4.69～2.4米。如果需要避免相位发生缠绕，则模拟产生的子带波长λs对应的模糊形变周期damb需覆盖最大形变量-4.69米，damb∈[-4.69,+4.69]，λs＝18.76米，对应的子带中心频率约为16mhz，子带中心频率偏移量f0＝8mhz，详情如图12所示。从图12中可看出上下子带的中心频率十分接近，为了不发生带宽上的重叠，子带带宽b最大值只能设置为16mhz，与表2中全带宽b＝100mhz的比b/b＝0.16。按照对子带中心频率和带宽最优化选取的讨论，一方面，上述方案中b＝0.16b已小于0.2b的最优下限，其保留的分辨率较小，会造成影像的模糊从而带来较严重的噪声误差；另一方面，该方案中上下子带带宽之间已经不存在间隔，此时谱间干扰的影响将最大化，会在子带干涉图中出现严重的线性趋势条纹；而如果为了减小谱间干扰的影响，唯一的途径只能进一步减小子带带宽以提高间隔，这样又会进一步降低分辨率使影像更加模糊。因此，该子带中心频率和带宽的配置方案，并不可取。
另一方面，landsat8提取的形变场显示(图8)紧邻破裂迹线南北两侧，是大形变量集中的区域，而在常规insar相位中(图9b)，这些地区由于形变量过大、地形复杂，以及配准技术的限制等原因，无法获得干涉条纹，出现了完全的失相干。这些在常规insar中失相干的区域，其干涉条纹不可识别，干涉相位被严重的噪声淹没，即使经过子带干涉处理，其结果也是不正确的。因此在设置子带中心频率时可以不考虑这部分地区贡献的形变量，实际有干涉条纹覆盖的区域，其斜距向形变量约为–4.0～2.4米。使用子带干涉法需要模拟产生至少16米的波长，即子带中心频率约为18.7mhz，子带中心频率偏移量f0＝9.35mhz，子带带宽b的最大值为18.7mhz＝0.187b，则不会产生相位缠绕。这一子带参数配制方案，虽然相较于上文有所改进，但上下子带中心频率仍然十分接近，子带带宽仍较小，受到分辨率模糊和谱间干扰的影响仍很大，因此也不应采用该参数方案。
综合上述分析，由于此次地震产生的地表形变量远大于一般子带干涉所能模拟的模糊形变周期，因此在子带中心频率和带宽的配置上，应以减少干涉条纹数为主要目标。所以申请人采用经过实验验证得到的最优参数方案，即子带中心频率为0.8b，子带带宽为0.2b，模拟产生的波长为3.75m，模糊形变周期为1.87米，通过减少干涉条纹数从而降低解缠难度，并保证测量精度。
3.2降噪。在处理流程上，遵循本申请提出的子带干涉流程，计算获得的偏移量多项式，进行子带辅影像的配准和重采样。图13a是依据上述流程进行子带干涉处理后的结果。干涉相位被大量噪声所掩盖，干涉条纹分布十分稀疏，但很难分辨。因此申请人在使用goldstein自适应滤波时，增加滤波器中傅里叶搜索窗口的长度和α值，而不必像常规insar一样过多考虑降噪可能对条纹形态的破坏。申请人使用了3次迭代降噪滤波处理，采用的傅里叶搜索窗口长度分别为128、96、64，α值分别为0.6、0.5、0.4，图13b是图13a经过上述降噪滤波处理后的干涉相位，噪声水平有很明显的下降。可看出干涉条纹分布在破裂迹线南侧十分稀疏，相位在大部分地区处在一个形变周期范围内，只在靠近破裂迹线附近出现了少量的缠绕现象；北侧的干涉条纹出现在影像西北角，并出现了明显的相位缠绕。图14是对破裂迹线南北两侧的子带干涉相位采用mcf方法解缠、计算地距向形变、地理编码并沿虚线剪裁拼接后的子带干涉地距向形变场。
4结果对比
将子带干涉法(图15b)分别与landsat8(图15c)常规insar(图15d)各自提取的地表形变场进行相互对比分析，并将使用模型拟合的形变场(图15a)作为基准参考。为了获得不同结果之间的细节分布，还分别在图9a-d的相同位置从南至北绘制了3条横跨破裂迹线的剖面p1，p2，p3。在噪声水平的评价上，参考了goldstein等人于1998年提出的psd指标法，使用一定尺寸的滑动窗口，对窗口内像素的形变量进行标准偏差的估计，该方法的公式如下：
式中w为窗口内像素坐标的集合，d(x,y)为窗口内像素坐标为(x,y)时对应的形变值，d为窗口内所有形变值的平均值，n为窗口尺寸，本文选用的尺寸为11×11。该公式计算的是整个窗口的形变标准差值σdsd，将其赋于窗口正中的像素，这样可获得与原图像点数相同的σdsd分布，σdsd越小说明图像的噪声水平越低。
landsat8的形变场采用的交叉频谱相关法是一种基于强度信息的子像素精密匹配算法，适用于大形变量的探测，如图15c所示，这种方法相比子带干涉法(图15b)和常规insar(图15d)获得了更完整的地表形变场和更清晰的地表破裂迹线，且迹线南北两侧的运动方向也相互吻合。但这种方法会受到地形、阴影和数据质量等因素的影响，其精度一般要比常规insar低。
1)对比上述3种手段提取的形变场的异常值误差
landsat8提取的形变场基于的像对间隔时间只有16天，在成像时间内基本不存在非地震引发的形变干扰，因此异常值主要表现在影像中出现一些模糊条带状，带有一定线性趋势，并呈南北纵向分布的轨道误差，由于破裂迹线南侧的形变值相比北侧较小，因此在南侧较容易观察得到这一轨道误差，同时噪声也在影像中有一定程度的分布。
子带干涉法和常规insar法提取的形变场，由于在相位解缠阶段已经进行了多次迭代降噪滤波，因此在条纹有效区域保留下来的形变场整体较为光滑，少有数值波动较大的像素点出现。
为了对比分析的一致性，利用公式4对上述3种手段获取的形变场的相同覆盖区域(去掉了区域2和区域4)，分别进行噪声水平估计，得到噪声水平分布(图16)。其中子带干涉法的σdsd(图16a)值域分布为0～1.5，平均值为0.046；landsat8的σdsd(图16b)值域分布为0～0.93，平均值为0.22；常规insar法由于解缠时低相干点被去除，因此参与计算的相位点较少，其σdsd(图16d)值域分布为0～0.027，平均值为0.002，在3种方法中最高。
从上述对比分析中可看出，两种sar干涉方法噪声水平要远低于landsat8。
2)在形变细节的反映上，能从子带干涉法提取的形变场中清晰地观察到位于破裂迹线南侧，区域3中的一处形变量较大且较密集分布的形变中心(黑色矩形框标出)，并能观察到该形变中心所产生的形变场，向南部周边地区呈半圆状扩散。这一形变中心，在landsat8提取的形变场中，也能十分模糊地观察得到，但类似子带干涉法中的半圆状扩散趋势则基本不存在。这说明子带干涉法在细节上拥有比其他两种处理手段更高的精度。在常规insar中，也呈现了这一半圆状的形变中心分布，但十分不明显，说明此处形变量变化速率较为缓慢，形变值整体较小，这一现象将在后续的剖面分析中进行详细地讨论。
3)残差分布
将上述3种手段获取的形变场与模型拟合得到的形变场进行差分，得到各自的残差分布(图17)。landsat8获取的形变场残差，在破裂迹线南侧，由于实际的形变值比北侧偏小，因此受到轨道条带状误差的影响更为严重，如图17中箭头所指能在南部的残差上一些模糊的南北走向的线条状残差，除此之外大部分地区的残差分布较为散乱，残差水平在-1～+1米的范围内，基本符合光学影像基于交叉频谱相关算法估计精度范围内；在破裂迹线北侧，在紧邻迹线的区域(区域4)残差值普遍为正值且较大，在远离迹线的西北区域(区域5)，则表现出总体吻合局部偏大的特点。
子带干涉法和常规insar的残差分布，在区域1中的形变场与模型拟合结果十分吻合，其平均误差分别在0.05～0.08米的范围内，该区域位于影像南部是离形变中心最远的区域，此处形变条纹分布十分稀疏解缠误差较小，较容易获得高质量的测量结果。区域3位于影像中部，破裂迹线的南侧，是形变剧烈的区域。因此在残差图上，此处出现了较大的正值误差，在各自对应的残差图上此处均显示一定程度的正值误差，而常规insar的形变场在此处与模型拟合结果吻合较好。所以认为，上述出现较大误差的2种处理手段(子带干涉、landsat8和)，具有一定的一致性，且是由不同数据源和不同处理手段得到的，因此这一形变值剧烈增大的区域可能客观存在。如果这一假设成立，则模型拟合结果此处形变值变化平缓，可能与四叉树降采样的密度设置有关，也有可能和单条断层模型的设置有关。而insar的结果也未出现形变值剧增的趋势，则可能是因为对密集条纹做解缠时，因噪声和低相干性的影响，造成条纹混叠并丢失导致。在该区域的其他位置，子带干涉测量和常规insar，与模型拟合的结果相比整体偏小残差值在-0.2～-0.5米的范围内。在区域5中，子带干涉法的残差分布由北至南，由形变远场至近场呈现递减的趋势。最远端与模型拟合结果的误差最大，达到了0.7m，随着逐步靠近形变近场，误差逐步控制在0.5～0.2m或更低的水平。而常规insar在该处出现了较大误差，误差范围为3～3.5米，其原因也与该处干涉条纹密集，解缠时条纹丢失有关。对于子带干涉，在该处只发生了一个周期的相位缠绕，可十分容易地.解缠得到更接近真实情况的绝对相位差。注意到，子带干涉法在此处相虽然未出现类似常规insar的解缠错误，而获得较为准确的形变值，但相比影像南部出现了更大的误差，其原因是因为采用了常规insar的多项式全局配准策略，在处理影像局部存在的大偏移、复杂地形情况时出现的局部偏差，这一配准处理方法上的局限，在常规insar中也是存在且难以避免的。
4)剖面线分析
剖面分析可以直观地反映形变场在剖面上分布的细节。在剖面线p1、p2、p3各自对应的形变分布曲线中，随距离由近及远形变趋势均呈现，先缓慢上升至约1.5～2米的最大值后，急剧下降到约-4.3～-4.5米的最小负值，然后又从最小负值逐渐增大趋向于0米的过程。这一过程，印证了破裂迹线南北两侧反向运动的特点，而在形变分布曲线中正负号发生突然转换的区域，正是地表破裂迹线的位置所在。在landsat8的形变分布曲线中，能清晰地观测到破裂迹线所处的位置，但这方法与模型拟合的形变分布曲线在迹线处并不完全重合，而子带干涉法和常规dinsar法在临近破裂迹线的南北两侧都发生了严重的失相干现象，并已在上文的处理中将之掩膜去除，因此在形变分布曲线中，出现了约有10～15公里长的缺失。
在上述形变场的南端，对应图15中区域1，和剖面线p1的0～10公里、剖面线p2的0～17公里处。子带干涉法、insar法的形变曲线围绕模型拟合的曲线上下小幅波动，波动范围在0.01～0.12米之间，并且形变值上升的趋势十分缓慢，上升速率平均只有0.005～0.008米/公里。相比之下，landsat8发生了稍大程度的波动，波动范围在模型拟合的形变曲线上下约0.38～0.49米之间，也存在一定的上升趋势，上升速率约为0.038～0.047米/公里。可看出上述3个形变场数据在该处与模拟结果均较为吻合，且噪声水平均较低。该处由于远离破裂迹线，形变量很小，加之地形由平原和海拔500米左右的山岭组成，因此在常规insar和子带干涉法中该处形变条纹分布稀疏，失相干区主要呈小斑块状分布，干涉质量高易于解缠获得高精度的相位，而在landsat8中，该处受地形山体阴影等因素干扰较小也能获得高质量的形变估计结果。
随着剖面线覆盖范围向北移动，在形变场的区域2位置，对应p1中10～20公里，p2中18～22公里和p3中2.5～10公里处。该处处于子带干涉和常规insar的失相干带上，已经被掩膜去除而没有形变数据。landsat8的形变曲线在此处波动也有一定的增大，波动范围在0.377～1.4米。
在形变场紧邻破裂迹线南侧的区域3位置，对应p1中20～30公里，p2中22～42公里和p3中10～20公里处。其中p1和p3在该处的子带干涉法和常规insar法的形变上升速率与模拟结果均较为接近，约为0.041米/公里，上升速率与区域1处相比并没有很明显的增大变化；两种方法的噪声水平均较低，但形变值始终略低于模拟结果，相差约为0.22米。landsat8在p1的该处位置均出现了较大程度的波动，这可能与landsat8在该处正位于轨道误差范围干扰有关。在p2的38～42公里处，子带干涉和landsat8的形变曲线，相比insar和模拟的形变曲线，出现了幅度较大的突然上升。在4公里的范围内，子带干涉法测量到的形变值从原先的0.5米增大至1.6米的最大值，landsat8的形变曲线上升更为明显，从原先的0.77米增大到了2.44米的最大值。该处正是上文讨论的黑色实线矩形框标出的疑似形变中心所在区域，可看出上述三种测量手段在趋势上均较为吻合，但与常规dinsar和模拟结果相差较大，这也印证了上文分析中常规insar存在的在密集条纹区解缠错误的问题。
在形变场紧邻破裂迹线北侧的区域4位置，对应p1中33～41公里，p2中45～50公里，p3中22～30公里处。该处在子带干涉和常规insar中处于失相干带上，也已被掩膜去除。landsat8的衰减速率有0.11～0.22米/公里这与模拟结果相吻合。造成上述较大误差的原因，与该处地形复杂、山岭林立有关，landsat8出现较多的地形阴影，从而影响了偏移量的估计精度并造成了强噪声的出现。
位于形变场西北部的区域5位置，对应p1中43～50公里，p2中50～62公里，p3中30～35公里处。由模型拟合得到的形变场(图15a)及对应的剖面显示(图19)，该地区形变量约为–1.4～–3.7m，由南至北基本呈线性递减，递减速率(斜率)约为0.15米/功率。landsat8的形变曲线与模拟的形变曲线基本吻合，但存在振幅波动较大的噪声。而常规insar计算得到的形变，与其他方法得到的形变差异很大，这些差异在剖面上同时反映在形变速率以及形变值分布这2个方面：其形变递减速率十分平缓只有0.026米/公里，与0.15米/公里的模拟结果差异很大；而其计算得到的形变值在剖面中约为–0.002～–0.004米，与模型拟合结果以及其他观测处理手段的结果差异也很大(图19中的p2)。这说明在上文分析中，该处解缠出现的严重误差是真实存在的。而子带干涉则很好地克服了这一不足，在该处可十分容易地解缠并计算出较准确的形变结果。在剖面中，由子带干涉获取的形变曲线衰减速率约为0.23米/公里，值域约为–0.87～–3.77m与模型拟合结果相接近，且噪声波动幅度在0.287～0.53米范围内，小于landsat8和offset-tracking0.48～1.67米的噪声幅度。
5子带干涉法的优势
子带干涉法基于sar卫星影像，其分辨率高于landsat8光学影像，因此理论上其测量精度也会优于前者。通过对此次震例进行的研究，可看出子带干涉法在提取大尺度地表形变上具有很大的优势。首先，虽然受失相干的影响其提取的形变场范围相较于landsat8数据采用的交叉频谱相关法法有所缺失，但在共同覆盖的区域，子带干涉法所测量的结果均显示出了更高的精度和更低的噪声水平。其次，在insar存在干涉条纹且十分密集的区域如：影像区域3中的矩形框位置和区域5，子带干涉法不但能更为准确地获得绝对相位差，而且与交叉频谱相关法相比，其反映的形变场细节更丰富。对紧邻破裂迹线的大形变地区的测量，通常是地震研究的重点，而此处也往往会因形变量过大而造成失相干，子带干涉法虽然相比常规insar对相干性的依赖更低，而更依赖于条纹是否存在，但在完全失相干的地区也意味着无法识别出干涉条纹，因此也无法对该处进行测量。
综合上述分析，对于强震的同震形变场提取，其受制约的因素更多，比高程提取更为复杂。子带干涉法可以在大形变梯度产生的密集条纹处获得精度更高细节更丰富的测量结果。
通过采用本发明公开的上述技术方案，得到了如下有益的效果：1、本申请所述方法在子带分割前，设置子带分割的最优参数子带带宽b＝0.2b和子带中心频率f中心频率＝fc+f0；f0∈[0.3b,0.4b]，最大程度减少上下子带频谱重叠带来的噪声以及旁瓣效应。2、本申请所述方法增加了干涉像对的常规insar精配准步骤，以获得精确的偏移量多项式，供后续子带干涉对配准处理中使用，避免由分辨率损失带来的配准误差或由频谱变化而引起的线性条纹误差。3、本申请所述方法公开了构建窄带陷波器并配合chirp-z变换的抑制方法，以解决因旁瓣能量泄漏出现的误差。4、本申请所述方法公开了基于采用常规insar的干涉图进行有效性干涉条纹检测，掩膜出子带干涉结果中有效的地形和形变数据区域，以解决利用相干性作为指标掩膜时保留无效区域数据的问题。
以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视本发明的保护范围。