一种去除InSAR大气相位的方法与流程

一种去除insar大气相位的方法
技术领域
1.本发明属于合成孔径雷达干涉测量技术（insar）领域，涉及一种去除insar干涉图中的大气相位的方法。


背景技术：

2.insar技术是近三十年来快速发展起来的一项对地观测技术，具有较强的云雾穿透能力，能够全天时、全天候获取大范围地表沉降信息，在滑坡、地震、火山、城市失水等引起的地面沉降监测方面得到了广泛的应用。insar技术利用两次sar影像的相位差来获取目标物的形变信息，其理论精度可达到毫米级。雷达信号在传输过程中会受到大气温度、湿度、气压以及自由电子密度等易变条件的影响，从而造成雷达信号延迟（大气延迟）。由大气延迟引起的大气相位直接影响形变结果的精度，因此需对大气相位进行精确剔除。目前常用的大气相位去除方法主要有多项式拟合、气象信息模型、gps数据和水汽数据校正等，然而以上几种方法并不适用于所有insar干涉图的大气相位去除，以上方法对于部分干涉图的大气相位去除效果较差。


技术实现要素：

3.本发明的目的是提供一种insar干涉图大气相位去除方法。在常规大气相位去除方法效果不佳时，本发明使用bp神经网络模拟insar干涉图的大气相位，然后再将利用bp神经网络模拟的大气相位去除，本发明能够较好的去除insar干涉图的大气相位，进而提高insar技术的监测精度。
4.为实现上述技术目的，本发明采用如下方案：一种去除insar大气相位的方法，包括：获取遥感影像并基于sbas技术获得n个解缠干涉相位图，每个解缠干涉相位图中包括p个永久散射体点；利用二次多项式大气相位校正方法对n个解缠干涉相位图进行大气相位去除，并对大气相位去除效果按照去除效果进行次序排列，取去除效果差的m个结果，查找其对应的未去除大气相位的解缠干涉相位图归为a类干涉相位图；分别以a类干涉相位图中每个干涉相位图的经度和纬度作为输入，对应的相位值作为输出，建立包含输入层、隐含层和输出层的bp神经网络训练模型；将a类干涉相位图中的干涉相位图减去bp神经网络训练得到的大气相位得到去除大气相位的干涉相位图。
5.作为一种优选的实施方式，所述bp神经网络训练模型的建立方式为：对a类干涉相位图中的每个干涉相位图，将p个永久散射体点的经纬度和相位值转换成一个p*3的矩阵i；矩阵i第一列为干涉相位图中每个ps点的经度，第二列为干涉相位图中每个ps点的纬度，第三列为干涉相位图中每个ps点的相位值，对矩阵中前两列数据进行归一化处理得到矩阵i2；对a类干涉相位图中的第一幅干涉相位图，随机选取矩阵i2的部分数据作为bp神
经网络的训练样本，剩余数据作为测试样本，以训练样本作为输入，对应的相位值作为输出建立bp神经网络训练模型；使用训练好的bp神经网络对第一幅图的矩阵i的所有数据进行仿真获得第一幅干涉相位图的大气相位数据，将第一幅干涉相位图减去获取的大气相位得到去除大气相位的干涉相位图；重复步骤直至对a类干涉相位图中所有干涉相位图处理完毕，完成insar大气相位的去除。
6.作为一种优选的实施方式，对矩阵i进行等间距抽稀得到矩阵i1，将矩阵i1的前两列数据归一化到[-1，1]之间得到矩阵i2，用于后续bp神经网络训练/测试。
[0007]
作为一种优选的实施方式，预设参数goal为性能目标进行模型训练，包括：初次训练结束后查看训练目标值是否小于goal，如果目标值不小于goal，则判断训练样本和测试样本存在异常数据，对异常数据进行处理；如果目标值小于goal，利用该次训练的网络对测试样本进行模拟，如果使用该网络模拟的仿真数据与测试样本之间的评价指标满足要求，则训练结束，保存网络参数及结构；否则增加网络的初始隐含层节点个数，利用训练样本重新训练。
[0008]
作为一种优选的实施方式，对异常数据进行处理的方式为：使用当前bp神经网络利用训练样本和测试数据模拟出新的仿真数据，计算仿真数据与矩阵i2第三列数据的相对误差；相对误差超过100%时，将这些数据作为样本中的异常数据点剔除；增加网络的初始隐含层节点个数，再次使用剔除异常数据后的训练样本进行训练，训练结束后查看训练目标值是否小于goal，如果目标值小于goal，利用该次训练的参数对剔除异常数据后的测试样本进行模拟，如果使用该网络模拟的仿真数据与测试样本之间的评价指标满足要求，则训练结束，保存网络参数及结构；否则增加初始隐含层节点个数，去除异常数据，使用去除异常数据后的训练样本进行训练；反复调整网络的初始隐含层节点个数直至bp神经网络对训练样本的训练目标值小于goal和测试样本的评价指标满足要求，结束训练，保存网络参数及结构。
[0009]
作为一种优选的实施方式，当bp神经网络的初始隐含层节点个数达到预设上限值，且使用训练样本的训练目标值不小于goal或测试样本的评价指标不满足要求，则调高目标值goal的取值，再反复调整网络的初始隐含层节点个数直至bp神经网络对训练样本的训练目标值小于新的goal且测试样本的评价指标满足要求，结束训练，保存网络参数及结构。
[0010]
作为一种优选的实施方式，所述评价指标选用误差相关的统计指标，如误差平方和、均方根误差、相对误差百分比等。
[0011]
作为一种优选的实施方式，bp神经网络的隐含层节点数目为2~11。
[0012]
作为一种优选的实施方式，所述方法还包括，将剩余n-m个已去大气相位的干涉相位图归为b类干涉相位图；在完成insar大气相位的去除后，将b类干涉相位图和去除大气相位的a类干涉相位图合并得到新的干涉相位图集合，利用新的干涉相位图集合中的干涉相位图进行形变速率反演，获取整个研究区域高精度的形变速率图。
[0013]
常规去除insar干涉相位图中大气相位的方法往往效果不佳，本发明提供了一种二次多项式结合bp神经网络去除insar干涉相位图中大气相位的方法，先利用二次多项式
大气相位校正方法对n个解缠干涉相位图进行快速的大气相位去除，之后对去除结果按照去除效果进行次序排列，取去除效果差的若干结果，查找其对应的未去除大气相位的解缠干涉相位图，利用所述解缠干涉相位图中ps点的经纬度和相位信息构建出bp神经网络，再利用构建好的bp神经网络模拟出干涉图中的大气相位，进而去除大气相位。使用本发明提供的方法能够较快较好的去除大气相位，进而提高insar技术反演形变速率的精度，为insar技术更好的监测地表形变提供技术支撑。
附图说明
[0014]
图1a是扬中市指南村地理位置图；图1b是裁剪配准后的覆盖扬中市指南村的哨兵影像幅度图（影像获取时间为2016年8月12日）。
[0015]
图2是2016年8月12日获取的哨兵影像与2016年9月29日获取的哨兵影像生成的解缠干涉相位图及使用二次多项式去除大气相位后的解缠干涉相位图。
[0016]
图3中a是最终的训练样本；b是最终的测试样本。
[0017]
图4中a是使用bp神经网络训练样本模拟的大气相位；b是使用bp神经网络测试样本模拟的大气相位。
[0018]
图5中a是原始解缠干涉相位图；b是使用二次多项式去除大气相位后的解缠干涉相位图；c是使用解缠干涉相位图中所有ps点的数据利用bp神经网络模拟的大气相位；d是去除c中大气相位后的解缠干涉图。
[0019]
图6中a是原始解缠干涉相位的统计直方图；b是使用二次多项式去除大气相位后的解缠干涉相位的统计直方图；c是使用本发明提供的方法去除大气相位后的解缠干涉相位的统计直方图。
具体实施方式
[0020]
下面结合附图对本发明的具体实施方式做出详细的说明。
[0021]
1）实施案例选取覆盖江苏省扬中市指南村作为实验对象，其具体范围如图1a所示，选取覆盖研究区域的42景sentinel-1a影像为实验影像，影像获取时间为2015年11月22日至2017年10月30日，90m分辨率的srtm-3 dem为配准及去除地形相位的dem，matlab r2020b及基于matlab开发的sentinel-1影像处理开源软件ctsent_v2为本次实验影像处理软件；ctsent_v2软件自动选取2017年1月27日获取的影像为主影像，使用ctsent_v2软件对42景sentinel-1a影像进行裁剪并配准，配准后的2016年8月12日获取的影像幅度图如图1b所示，将每幅影像按照影像获取时间进行排序，接下来使用sbas-insar对42景配准后的影像进行处理，ctsent_v2软件的sbas-insar时序处理过程中每幅影像仅与排在自己后面的三幅影像生成干涉对，这样可以生成120个影像对，对每个影像对对应的主辅影像进行干涉处理获取120个干涉相位图，使用4：1作为多视比对每个干涉相位图进行多视处理，使用高斯滤波对每个干涉相位图进行滤波降噪处理，将相干性阈值设置为0.35，将幅度振幅离差阈值设置为0.6，对120个干涉相位图进行时序分析选取了p个ps点，对选取ps点后的120个干涉相位图进行解缠处理获取了120个解缠干涉相位图，每个解缠干涉相位图中都包括p个ps点。
[0022]
2）使用二次多项式对120个解缠干涉相位图进行大气相位校正，共获取了120个去
除大气相位的解缠干涉相位图，将未进行大气校正的解缠干涉相位图和已进行大气校正的解缠干涉相位图依次进行对比，找出利用二次多项式大气相位校正方法去除大气相位效果较差的干涉相位图，将此种利用二次多项式方法去大气相位效果较差的干涉相位图对应的未去大气相位的干涉相位图归为a类干涉相位图，a类干涉相位图共有m个，2016年8月12日获取的影像与2016年9月29日获取的影像生成的解缠干涉相位图就属于a类干涉相位图（图2），从图2可以看出，使用二次多项式方法进行大气相位校正后的解缠干涉相位图仍包括大量的大气相位噪声；将利用二次多项式方法去大气相位效果较好的已去大气相位的干涉相位图归为b类干涉相位图，b类干涉相位图共有120-m个。
[0023]
3）构建初步bp神经网络：
①
确认bp神经网络的输入节点个数：本发明的bp神经网络的输入为ps点的经纬度坐标，本发明bp神经网络的输入层神经元个数取值为2；对于输入层数据，由于数据间存在着量级差异，所以本发明分别按照2个输入层神经元将其归一化到[-1,1]之间，以此来消除输入数据间的量级差异问题。
[0024]
②
确认网络的输出节点个数：本发明的bp神经网络的对应输出为ps点的相位值，所以网络的输出层神经元个数为1个。
[0025]
③
确认网络的隐含层节点个数：选择隐含层数的经验公式如下：其中为输入层神经元个数，为隐含层节点个数，为输出层神经元个数，一般为小于或者等于10的正整数。
[0026]
④
构建bp神经网络，按照隐含层数的经验公式，本发明中输入层神经元个数b=2，输出层神经元个数c=1，隐含层节点个数取值在[2，11]之间，网络的初始隐含层节点个数取b1=2；使用matlab中newff函数构建bp神经网络，神经网络隐含层参数设置为2，bp神经网络学习率取0.001，训练步数取1000步，网络的初始训练目标值取goal=0.001，此时已构建好初步的bp神经网络。
[0027]
4）将2016年8月12日获取的影像与2016年9月29日获取的影像生成的解缠干涉相位图中p个ps点的经纬度和相位值转换成一个p*3的矩阵i，第一列为干涉相位图中每个ps点的经度，第二列为干涉相位图中每个ps点的纬度，第三列为干涉相位图中每个ps点的相位值，对矩阵i进行等间距为10的抽稀得到矩阵i1，将矩阵i1的前两列数据归一化到[-1，1]之间得到矩阵i2，随机选取矩阵i2的70%的数据作为bp神经网络的训练样本，将矩阵i2的剩下的30%的数据作为bp神经网络的测试样本，利用已构建好的初步bp神经网络对矩阵i2的训练样本进行训练，然后进行如下步骤：
①
初次训练结束后查看训练目标值是否小于goal，如果目标值不小于goal，则进行第
②
步，如果目标值小于goal，利用该次训练的网络对矩阵i2的30%的测试样本进行模拟，如果使用测试样本模拟的仿真数据与测试样本之间的评价指标满足要求，则训练结束，保存网络参数及结构；如果使用测试样本模拟的仿真数据与测试样本之间的评价指标不满足要求，则取网络的初始隐含层节点个数取b1=3，再次使用训练样本进行训练，如此反复调
整网络的初始隐含层节点个数直至bp神经网络对训练样本的训练目标值小于goal且测试样本的评价指标满足要求，结束训练，保存网络参数及结构。
[0028]
②
初次训练的目标值不小于goal，由于本发明中使用的训练和测试样本存在一定的异常数据，本发明在实施过程中对异常数据做了如下处理：使用当前网络利用训练样本和测试数据模拟出新的仿真数据，仿真数据与样本的第三列数据的相对误差表示如下：其中为相对误差，是使用bp神经网络模拟的仿真数据，为训练和测试样本的第三列数据。
[0029]
当bp神经网络模拟的仿真数据与样本的第三列数据的相对误差超过100%时，将这些数据作为样本中的异常数据点，剔除这些异常数据，取网络的初始隐含层节点个数取b1=3，再次使用剔除异常数据后的训练样本进行训练，训练结束后查看训练目标值是否小于goal，如果目标值小于goal，利用该次训练的参数对剔除异常数据后的测试样本进行模拟，如果使用测试样本模拟的仿真数据与测试样本之间的评价指标满足要求，则训练结束，保存网络参数及结构；如果使用训练样本训练的目标值不小于goal，则取网络的初始隐含层节点个数取b1=4，去除异常数据，使用去除异常数据后的训练样本进行训练，如此反复调整网络的初始隐含层节点个数直至bp神经网络对训练样本的训练目标值小于goal且测试样本的评价指标满足要求，结束训练，保存网络参数及结构。
[0030]
③
当网络的初始隐含层节点个数取b1=11时，使用训练样本的训练目标值和测试样本的评价指标都不满足要求，则将目标值goal的取值调大，再反复调整网络的初始隐含层节点个数直至bp神经网络对训练样本的训练目标值小于新的goal且测试样本的评价指标满足要求，结束训练，保存网络参数及结构。图3a是剔除异常数据点后最终用于训练的训练样本，图3b是剔除异常数据点后最终用于测试的测试样本，图4a是利用bp神经网络使用图3a中的训练样本模拟的大气相位，图4b是利用bp神经网络使用图3b中的训练样本模拟的大气相位。
[0031]
5）上述过程结束后，本发明的bp神经网络的最佳隐含层节点个数取值为6，最佳训练目标值取goal=0.018，至此一个2-6-1的bp神经网络就建立完成，使用该bp神经网络对矩阵i的所有数据进行仿真获得a类干涉相位图中第一个干涉相位图的大气相位数据（图5c），将a类干涉相位图中第一个干涉相位图减去上述利用bp神经网络获取的大气相位即可得到大气相位去除效果较好的干涉相位图(图5d)。三种干涉相位的相位均值和标准差的统计结果如表1所示。通过与图5a和图5b进行对比，可以发现，使用本发明提供的方法去除大气相位的效果明显优于二次多项式方法。图5a中原始相位的均值为3.274rad，相位标准差为1.277rad；图5b中使用二次多项式去除大气相位的干涉相位均值几乎为0rad，相位标准差为1.145rad，从相位标准差上来看，二次多项式几乎无法去除原始相位中的大气相位；图5d中使用本发明提供的方法去除大气相位的干涉相位均值为0.015rad，相位标准差为0.725rad，从相位均值与标准差上来看，本发明提供的方法能够去除大量的大气相位噪声，且去除效果显著。如图5中黑色方框内所示，使用二次多项式去除大气相位后，该处仍含有
大量的大气相位噪声，而使用本发明中的方法去除大气相位后，干涉相位中的大气相位几乎完全去除，效果极其显著。为进一步验证本发明提供的方法的可靠性，将黑色方框内的原始干涉相位、利用二次多项式去除大气相位的干涉相位和利用本发明方法去除大气相位的干涉相位进行直方图统计，统计结果如图6所示，从图中可以看出：原始干涉相位主要分布在[2，9]rad之间，相位均值为5.782rad，相位标准差为1.669rad，说明此处含有大量的大气相位噪声；利用二次多项式去除的干涉相位主要分布在[-2，4]rad之间，相位均值为1.229rad，相位标准差为1.641rad，说明此处仍然含有大量的大气相位噪声；而使用本发明提供的方法去除大气相位后的干涉相位主要均匀分布在[-2，2]rad之间，相位均值为0.050rad，相位标准差为0.799rad，说明此处仅含有极少量的大气相位噪声；从相位标准差上来看，二次多项式方法的相位标准差与原始相位的标准差比几乎没有差别，而本发明提供的方法的相位标准差比二次多项式方法小了2倍，这说明本发明提供的方法去除大气相位的效果相比二次多项式方法提升了2倍。图6的统计结果再次证明本发明提供的方法去除大气相位的效果更加显著。
[0032]
表1 三种干涉相位的统计结果对比 相位均值/rad相位标准差/rad原始干涉相位3.2741.277二次多项式去除大气相位的干涉相位01.145本发明方法去除大气相位的干涉相位0.0150.725原始干涉相位（黑色方框内）5.7821.669二次多项式去除大气相位的干涉相位（黑色方框内）1.2291.641本发明方法去除大气相位的干涉相位（黑色方框内）0.0500.799
6）依次对a类剩余的干涉相位图做步骤3）-步骤5）中的工作共得到m-1个相对应的去除大气相位的干涉相位图，将m个利用bp神经网络去除大气相位的干涉相位图和120-m个b类的干涉相位图进行合并得到新的干涉相位图集合c，利用集合c中的干涉相位图进行形变速率反演，获取整个研究区域高精度的形变速率图。