2016年2月6日早晨3时57分,台湾高雄发生6.7级地震,震源深度15千米。震中位于高雄市美浓区,震源深度16.7公里,云林县草岭的最大震度达6级,高雄、屏东、台南、嘉义等地的最大震度均达5级。本节以哨兵1Aslc作为数据源,使用DInSAR的方法,对2016年2月6日台湾南部的地震事件进行干涉测量。
数据情况:
/SARscape/Preferences,设置LoadPreferences为SENTINEL_TOPSAR,这套参数是专门针对哨兵1A数据的TOPSAR(IW)模式做InSAR处理时的系统参数。
图 设置SARscape系统参数
如果计算机支持GPU,在General Parameters面板,设置Opencl PlatformName的类型为相应的设备,如NVIDIA CUDA显卡。
设置Opencl的类型
为了便于后面操作的数据读取和输出,设置ENVI的系统参数中的输入输出路径。ENVI->File->Preferences,在Directories中设置默认的输入路径和输出路径。
设置ENVI系统参数
(1)数据文件输入
网站下载的数据解压,打开数据导入工具/SARscape/Import Data/SAR Spaceborne/SENTINEL1,在Input File List分别输入两景哨兵1A数据的元数据文件manifest.safe。
Optional Input Orbit File List,轨道文件,该文件是可选文件,在这里没有轨道文件,就不输入。
图 数据输入面板
(2)参数设置
切换到Parameters面板,主要参数就是对输出数据命名设置,推荐选择Rename the File UsingParameters:True,可以对输出的数据自动按照数据类型进行命名。如果要输出镶嵌后的SLC数据文件,可以在OtherParameters设置Generate SLCmosaic为True,否则默认是输出各个条带的slc数据。数据类型软件自动识别,不用做设置。
图 数据输入面板参数设置
(3)输出设置
之前设置过默认的输出路径,这里直接按照默认即可,如果要改输出路径,在数据输出路径上点击右键,选择Change OutputDirectries。
图 数据输出面板的输出设置
设置好参数后,点击Exec执行,完成后弹出对话框,点击End。
输出的数据文件包括:整景图像的强度图数据(_pwr)、slc索引文件(.slc_list)带有地理坐标的外边框矢量文件(.shp)、GoogleEarth文件(.kml)。打开矢量文件,ENVI主界面view->ReferenceMap Link,可以自动链接Arcgis Online的在线底图,查看数据的地理范围。
图 数据的地理范围
打开工具/SARscape/General Tools/Digital Elevation ModelExtraction/SRTM-3 Version 4,在InputFile输入两个时相的slc_list文件,DEM/Cartographic System设置GEO GLOBLEWGS84,其他参数按照默认。
注:1)这一步也可以在InSAR界面下进行。
2)如果网速许可,最好下载SRTMVersion4版本的DEM数据。SRTM Version2版本的数据有空洞,会对最终的结果有所影响。
图 根据数据范围下载DEM
自动下载的SRTM DEM数据
打开工具/SARscape/Interferometry/Interferometric Tools/BaselineEstimation,输入主从影像的_slc_list文件,点击Exec执行,输出基线估算的结果。
基线估算面板
基线估算结果:
Normal Baseline (m) =-64.213 Range Shift (pixels) =13.012 Azimuth Shift (pixels) =4.669 Slant Range Distance (m) =879019.561 Absolute Time Baseline (Days) =36 Doppler Centroid diff. (Hz) =17.918 2 PI Ambiguity height (InSAR) (m) =242.615 2 PI Ambiguity displacement(DInSAR) (m) = 0.028 1 Pixel Shift Ambiguity height(Stereo Radargrammetry) (m) = 20379.668 1 Pixel Shift Ambiguitydisplacement (Amplitude Tracking) (m) = 2.330 Master Incidence Angle =39.723 Pair potentially suited forInterferometry, check the precision plot |
基线估算的结果显示,这两景数据的空间基线为64.213米,远小于临界基线6474.081米,时间基线36天,做DInSAR的一个相位变化周期代表的地形变化为0.028米。
打开/SARscape/Interferometry/DInSAR Displacement Workflow工具。
(1)文件输入
注:数据有两种极化方式:VH和VV,选择同极化方式的做差分干涉处理。
图 数据输入面板
输入的数据文件设置好之后,点击Next。弹出自动计算视数的面板,算出来的视数为5:1,点击确定。
图 自动计算的视数
(2)干涉图生成
生成TIFF数据(MakeTIFF):Ture,生成TIFF格式的中间结果,如果需要使用中间结果,如写文章时候作为插图,可以设置为True,其他步骤类似。
图 干涉图生成参数设置面板
处理完成之后,ENVI视窗中自动加载了去平后的干涉图,以及主从影像的强度图,这一步生成的数据结果存放在ENVI默认的临时文件路径下,默认路径为:C:\Users\<计算机用户名>\AppData\Local\Temp,自动生成了一个名为SARsTmpDirXXXXXXXXX的文件夹,这一步生成的结果有:
图 去平后的干涉图_dint
(3)滤波和相干性计算(Adaptive Filter and Coherenace Generation)
提供了三种滤波方法:
这种方法适用于高分辨率的数据(如TerraSAR-X或COSMO-SkyMed)
使用局部干涉条纹的频率来优化滤波器,该方法尽可能的保留了微小的干涉条纹。
这种滤波方法的滤波器是可变的,提高了干涉条纹的清晰度、减少了由空间基线或时间基线引起的失相干的噪声。这种方法是最常用的方法。
图 滤波和相干性生成参数设置面板
点击Next按钮,进行干涉图滤波和相干性生成处理,处理完成之后,自动加载滤波后的干涉图_fint和相干性系数图_cc。这一步处理之后生成的结果有:
图 滤波后的干涉图和相干系数图
(4)相位解缠(PhaseUnwrapping):相位的变化是以2π为周期的,所以只要相位变化超过了2π,相位就会重新开始和循环。相位解缠是对去平和滤波后的相位进行解缠处理,使之与线性变化的地形信息对应,解决2π模糊的问题。
图 相位解缠面板参数设置
说明:
1)解缠方法(Unwrapping Method Type)提供了三种:
DelaunayMCF:和最小费用流法的不同在于,这种方法不是考虑了图像上所有的像元,而是仅考虑了相干性大于阈值的部分,而且不是用正方形的格网而是用了德罗尼三角形格网。只有对相干性高的部分进行解缠,不受低相干像元的影响。对于有大量相干性低的地物存在的时候,如影像上存在大量水体、浓密植被等,推荐使用该方法,最小化相位突变的影响。
2)分解等级(Unwrapping DecompositionLevel):该分解是为了用迭代的方法对数据做多视和疏采样:干涉图以一个较低的分辨率被解缠然后被重采样成原来的分辨率。使用了分解可减少解缠错误,提高处理效率。用户可以指定迭代的次数,每个迭代相当于3次采样过疏。这里可以填的数有:-1、0、1、2、3。其中,-1和0代表不执行分解,用原始的像素采样,当形变很大或是地形很陡峭的情况下(多相位/高密度分布),分解可引起交迭效应,建议设置为-1或0;1代表最小的分解等级。3:最大的分解等级。建议这个次数不要超过3。通常情况设置原始的像素采样(如-1)或者最小的分解等级(1)。
点击Next按钮,进行干涉图滤波和相干性生成处理,处理完成之后,自动加载滤波后的相位解缠结果图_upha。这一步处理之后生成的结果有:
图 相位解缠结果
(5)控制点选择(GCP Selection):输入用于轨道精炼的控制点文件,可以用已有的文件也可在此选择控制点。
在Refinement GCP File(Mandatory)项中,点击
图 选择控制点流程化工具
在控制点生成面板上,点击Next,打开控制点选择工具,鼠标变为选点状态,单击鼠标左键就可以选点。在图像周边远离形变区域的地方、相位好的地方选择若干控制点。
图 控制点选择
单击Cartographic System按钮,查看控制点的参考坐标系统,该坐标系是从参考DEM上自动读取的。
图 设置控制点的坐标系
单击Export按钮,查看控制点的存放路径和文件名。生成的控制点文件为INTERF_out_upha_gcp.xml。
注:如果默认输出的路径有中文字符的文件夹或者文件名,需要改到只有英文字符的文件路径中。
图 控制点文件输出
在控制点生成面板上点击Finish,生成了控制点文件,并自动添加到InSAR流程化处理面板的Refinement GCPFile(Mandatory)项。
图 生成控制点文件
在面板上点击Next按钮,进入下一步。
(6)轨道精炼和重去平(Refinement andReflattening)进行轨道精炼和相位偏移的计算,消除可能的斜坡相位,对卫星轨道和相位偏移进行纠正。这一步对解缠后的相位是否能正确转化为形变值很关键。
图 轨道精炼和重去平参数设置面板
点击Next按钮,进行轨道精炼和重去平处理,处理完成之后,将优化的结果显示在RefinementResults的面板,内容如下:
ESTIMATE A RESIDUAL RAMP Points selected by the user =64 Valid points found = 64 Extra constrains = 2 Polynomial Degree choose =3 Polynomial Type : = k0 + k1*rg +k2*az Polynomial Coefficients (radians): Root Mean Square error (m)=64.1891189212 Mean difference after RemoveResidual refinement (rad) = -0.0431859477 Standard Deviation after RemoveResidual refinement (rad) = 1.9514668017 |
这一步得到的结果有:
(7)相位转形变以及地理编码(Phase to Displacement Conversion andGeocoding):将经过绝对校准和解缠的相位,结合合成相位,转换为形变数据以及地理编码到制图坐标系统,默认得到的是LOS方向的形变。
图 相位转形变主要参数
图 地理编码参数
点击Next按钮,进行相位转形变和地理编码处理。地理编码的坐标系是以参考DEM的坐标系为准。这一步得到的结果有:
图 相位转形变的结果
(8)结果输出(Output):
结果默认输出在ENVI的默认输出路径下,文件名命名为sentinel1_69_20160109_100020274_IW_SIW1_A_VV_slc_list_out_disp。若想保留中间结果便于查看,将DeleteTemporary Files不勾选。
点击Finish,输出结果。结束DInSARDiaplcement处理的工作流。生产的形变数据自动进行密度分割配色展示。
图 形变结果自动配色显示
注:Back和Next按钮,可切换至中间某一步查看参数或调整参数重新处理。
结果显示,本次地震发生的位置,导致了一个区域明显的抬升。如果对整景数据做DInSAR处理,右边的大部分没有形变的区域,都参与到处理中,容易引入一些误差(能看出右侧的稳定区域,结果也有一些抬升的情况)。在这种情况下,可以在_fint的数据上,对地震形变区域手动绘制一个矢量区域,用SARscape的裁剪工具对原始的SLC数据进行裁剪,然后对子区域做DInSAR处理,能得到更精确的结果。
ENVI中手动绘制矢量范围:打开_reflat_fint数据,ENVI主菜单File->New-> VectorLayer,选择_fint数据,绘制矢量区域,如下图所示:
图 绘制子区域的矢量范围
在矢量图层上点击右键,选择Save As,保存为subarea.shp文件。
点击/SARscape/General Tools/Sample Selections/Sample Selection SARGeometry Data工具,打开Sample Selection SAR Geometry Data面板。
图 子区裁剪的数据输入面板
图 输入裁剪矢量范围及参考的强度数据
图 子区裁剪的参数设置面板
图 裁剪结果输出面板
图 裁剪结果强度图
然后对裁剪后的slc数据,进行DInSAR工作流处理(参数和之前相同),得到地震区域的形变结果。在子区域DInSAR处理中,在轨道精炼一步,选择的控制点如下图所示:
子区域DInSAR处理中轨道精炼的GCP
得到的轨道优化的结果报表:
ESTIMATE A RESIDUAL RAMP Points selected by the user =27 Valid points found = 27 Extra constrains = 2 Polynomial Degree choose =3 Polynomial Type : = k0 + k1*rg +k2*az Polynomial Coefficients (radians): Root Mean Square error (m)=16.5580109079 Mean difference after RemoveResidual refinement (rad) = -0.1456119011 Standard Deviation after RemoveResidual refinement (rad) = 0.7900878373 |
图 裁剪后的SLC做DInSAR处理得到的结果
联系客服