Determination of Ice-flow Velocity at the Polar Record Glacierand Dalk Glacier Using DInSAR
-
摘要: 目的 利用1996年 ERS-1/2tandem SAR数据,采用 DInSAR方法提取了极记录和达尔克冰川的冰流速,并利用冰流速转换模型和该地区的实测值进行了对比。结果表明,DInSAR结果很好地体现了整个区域的冰流速。对于内陆地区,DInSAR提取的监测点冰流速和实测值非常吻合;对于达尔克冰川的前缘,受时间差异和潮汐等因素的影响,DInSAR推求的冰流速结果略小于实测值。Abstract: Objective ERS-1/2tandem SAR images acquired in 1996were utilized to map the ice-flow velocity ofthe Polar Record and Dalk Glaciers using DInSAR.The horizontal component was derived from line-of-sight velocity estimates,and compared with the field measurements.The DInSAR results revealedan ice velocity map of the whole area.For the inland areas,the DInSAR results were in good agree-ment with field measurements.As for the front of the Dalk Glacier,the ice velocity derived from DIn-SAR was smaller than the field measurements because of the differences between data acquisitiontimes,tidal motion,etc.
-
Keywords:
- DInSAR /
- ice-flow velocity /
- Polar Record Glacier /
- Dalk Glacier
-
倾斜摄影测量技术已广泛应用于逼真地形景观建模等领域,其基本特点是从多角度影像自动获取密集点云数据,通过简化处理生成不规则三角网格表面模型(triangulated irregalar network,TIN),并自动映射纹理[1, 2]。实际上,受多视影像联合平差精度和密集匹配生成的点云精度的限制,以及点云简化的不足,会导致生成的三角形网格模型表面粗糙,连续区域表面法向量存在噪声,从而造成不一致性。这种不一致性导致对应于同一地物表面的三角面片被投影到不同影像上;同时,由于相机畸变、光照、地物反射特性等原因,不同影像间存在几何和辐射上的差异,最终纹理映射后的结果呈碎片化状态,三维场景真实感不强、不直观。当前主流的倾斜影像处理软件如Street Factory等也都存在类似的问题。
纹理映射碎片化问题的解决方法主要有渐入渐出法[3]、光照连续性重建[4, 5]法以及应用较广泛的多波段融合法[6-9]等。这些方法都是从影像出发,通过匀光或色彩融合等处理,使其在灰度上趋于一致,虽然在一定程度上可以使碎片之间的纹理过渡更加平滑,但计算量大,算法要求高,容易出现模糊和重影现象,并且依然难以消除碎片化的影响。
纹理映射碎片化的根本原因是模型表面上相邻的三角面片由于几何不一致被投影到不同影像上。现有方法只是在影像层面对纹理进行修复,并没有从根本上解决模型表面的几何不一致性问题。针对此,本文提出一种局部区域表面一致性约束的纹理映射方法,采用一致性约束的区域生长方法和分块区域整体的纹理优选方法,解决上述纹理碎片化的问题。
1 算法基本原理
局部区域表面一致性约束的纹理映射方法流程如图 1所示。记三角形网格模型为M[9],构成M的带有法向量的三角形记作有向三角面片t,则M的连续表面S={t1,t2,…,tm},S∈R3。首先,在表面一致性约束的前提下,以有向三角面片t为基元,采用表面一致性约束的区域生长算法将S分割为多个分块连续的区域R1,R2,…,Rn,RiS(1 ≤ i ≤ n),且分割结果应在最大程度上满足人们对建筑物几何表面划分的认知,即连续、光滑建筑物有且仅有一个表面与之对应。为了保证局部区域内纹理映射的连续性,需要根据区域整体的纹理优选方法,从覆盖区域Ri的多张不同角度的影像中,选择一张最佳影像为Ri提供纹理。最后根据投影矩阵计算Ri上的点在影像上的投影,得到纹理,完成纹理映射。
1.1 表面一致性约束的区域生长方法
纹理映射的碎片化是由于三角形网格模型M表面的几何不一致性引起的。产生不一致性的原因是影像平差精度和密集匹配的点云精度的限制以及构网过程中的简化算法不够完善,导致连续区域表面的法向量存在异常噪声。如图 2所示,Γ为建筑物真实的表面,n为其法向量,t1、t2为M表面两个相邻的三角面片,n1、n2为其法向量。
单独来看,t1、t2与Γ之间的夹角θ1、θ2均在噪声容忍范围内,但t1与t2之间的夹角θ=θ1+θ2却超过了噪声容忍范围,形成噪声。表面一致性约束不是剔除噪声点,而是通过约束,使噪声点可以与邻近点形成一个具有一致性的区域,即该区域是连续的,并且区域上任何一点在同一张影像上具有相同的可见性。因此,表面一致性约束包括邻接性约束和平面性约束两方面。
1.1.1 邻接性约束
邻接性约束,即三角面片在空间上是相邻的,能保证区域的连续性。邻接性约束可通过采用树状的索引结构(KD-Tree、R-Tree[10-12]等),进行半径搜索或K邻近搜索来实现。但是由于三角形网格模型密度的不均匀性,进行半径搜索时难以设置合适的阈值以保证每次搜索的有效性,故本文中采用K邻近搜索方法。
K邻近搜索方法实现区域生长的关键是K的取值。平坦的模型表面K一般取较大值,可以加速计算;而对于房屋角点、植被等复杂表面,由于存在多个平面,K一般由大到小渐进取值,这样既能保证分割出每一个平面,又能保证计算效率。本文中,首先令K取一个较大的值(K=20)进行搜索,判断搜索到的点是否满足平面性约束,若不满足则说明是复杂表面,此时令K值减半继续搜索,以此类推。这样既可以保证计算速度,又可以保证区域生长结果在多个平面交界处的准确性。
1.1.2 平面性约束
实际中,许多建筑物的表面由连续曲面构成,但是曲面上的不同点在同一张影像上不具备相同的可见性(不考虑其他面的遮挡)。如图 3(a)所示,Γ 为建筑物真实的表面,I1、I2、I3是不同角度的拍摄的影像,p0 ~p4为Γ上的点,p0 ~p2为I1的可视范围,p1 ~p3为I2的可视范围,p2 ~p4为I3的可视范围,可见Γ上的点不可能全部投影到同一张影像上,因此需要对其进行平面性约束,用多个平面逼近曲面,如图 3(b)所示。
设有向点集P(P∈R3),P的平面度可以通过其协方差矩阵的主成分分析来判断[13],对协方差矩阵进行奇异值分解(singular value decomposition,SVD)后得到特征值σ1 <σ2 <σ3,η1、η2、η3为对应的特征向量,如图 4所示。特征值及其对应的特征向量代表了P在三维空间中的分布情况,σ1越小P越接近平面,极限情况为σ1=0。因此,平面性约束可记为:
(1) 式中,σ2 /σ1 >ξ1可保证点集P近似为平面,而σ3 /σ2 <ξ2则能避免产生过于狭长的条带。本文中采用这种方法作为平面性约束,令ξ1=45,ξ2=15,既有效约束了平面的形状,也较好地包容了噪声点。
若P满足上述约束,则可拟合出平面并记n0为拟合平面的法向量,近似于图 2中建筑物表面Γ的实际法向量n,x为平面上的点。设P中一点pi的法向量ni与n0之间的夹角记为θ(ni,n0),需满足:
(2) 由前面对法向量噪声的分析可知,式(2)避免了法向量噪声对平面的影响。式(1)和式(2)共同构成平面性约束,在区域生长时可将连续的曲面及平滑的过渡地带分割为多个平面逼近。
1.1.3 顾及约束的区域生长方法
区域生长以有向三角面片t为基元,采用迭代的策略由种子点向外生长。为了方便计算,将t的几何重心点及其法向量构成有向点p,则全部有向三角面片t的集合S即对应于有向点集P。区域生长需要解决的三个关键问题为种子点的选取,包含进区域的准则以及生长停止的条件。
以P为数据源创建Kd-Tree,记生长出的区域为R,其边界为B,符号⊕(p)表示添加点p,Θ(p)表示从集合中删除点p,并且添加一个点后都需要将它从KD-Tree中删除,以避免重复计算。首先,从P中任选一点p0作为种子点,R⊕(p0)、B⊕(p0),标记p0为起始点;以B中第j个标记为起始点的点pBj为中心进行K=20的邻近搜索,若搜索到的点{pnew}与R一起满足平面性约束式(1)、式(2),则R⊕{pnew}、B⊕{pnew},标记{pnew}为起始点,并且BΘ(pBj);若不满足平面性约束式(1)、式(2),则K值减半进行搜索,若直到K=1时{pnew}与R一起仍不满足约束条件,则标记pBj为终止点;然后重新令K=20,搜索B中第j+1个标记为起始点的点;逐此迭代,直到B中所有点均标记为终止点,或K邻近搜到的结果为空为止,则停止生长。最后得到的R即为所求的连续区域,它拟合的平面的法向量可视为R的法向量nR,nR是后面选择投影影像时的关键因素。若KD-Tree不为空,则以剩下的点中任意一点为某个区域生长结束后的种子点,进入下一个区域的生长。
1.2 分块连续区域的整体纹理优选方法
倾斜摄影测量具有多角度、多时相的特点,因此区域R会被多张影像I1,I2,…,Im同时覆盖。R与影像间的方位关系不同,在影像上产生投影畸变的程度也不同;同时由于区域间的遮挡关系,R在不同影像上的可见情况不同。因此,需要选择一张最佳的影像来为R提供整体的纹理。
1.2.1 视角关系判断及区域间遮挡检测
倾斜摄影测量中不同角度影像在空间中的方位是固定的,以5头相机为例, 影像范围覆盖同一区域的相机A、B、C、D、E分布如图 5(a)所示。定义ρ为区域与影像之间的视角关系权重,它与两个角度有关:① 区域R的法向量nR与每张影像投影反方向nI之间所成的夹角,记为α;② 以R的中心为起点,连接影像中心,得到向量l,记l与nR之间的夹角为β,如图 5(b)所示,A1、A2、A3为不同时刻相似A在空中的位置。β越小,则R在影像上的投影越靠近影像中心。
则ρ的确定方式如下:
(3) 当ρ大于0时,区域R在影像上可见(不考虑遮挡和超出影像范围)。为了保证投影的质量,令:
(4) 以区域R为整体进行遮挡检测,若R中总三角面片个数为s,被遮挡和超出影像范围的个数为k,记η为R在影像上的遮挡比率:
(5) 1.2.2 最佳影像优选方法
区域R整体的纹理优选方法需同时考虑视角关系ρ和遮挡比率η。首先找出所有ρ大于0对应的影像,并从中选取ηmin所对应的影像作为R的最佳投影影像Ibest。若ηmin大于0,则说明R在Ibest上不完全可见,此时需要将不可见的部分分割出来作为一个单独的区域R′,再以同样的方法为R′选择最佳投影影像。最后,从每一块区域对应的最佳投影影像中提取纹理,进行纹理映射。
2 实验分析
为了验证上述方法的有效性,本文采用SWDC-5的倾斜摄影影像数据来进行实验。整个测区共被15张影像所覆盖,相机A、B、C、D、E分别各有3张影像,影像覆盖区域内包括道路、建筑物(弧形/非弧形)、操场、植被等多种地物。对这些影像进行多视影像联合匹配后生成点云,再利用构网算法得到三角形网格模型,作为纹理映射的模型基础。
对比实验将三角形格网模型M中的有向三角面片t逐个按照§1.2中的方法选择最佳投影影像,然后再将投影到同一张影像上的相邻三角面片连接为一个区域R。这种逐个三角形纹理映射的方法与本文方法之间的区别在于是否进行了一致性约束的区域生长。图 6(a)为逐个三角形纹理映射的结果;图 6(b)为本文方法进行一致性约束的区域生长后的分割结果。将三角形个数小于一定阈值的区域视为碎片,可看出图 6(a)中碎片化现象严重,而图 6(b)中区域连续一致,地面和立面、立面和屋顶被很好地分割开,并且弧形的建筑物表面被分割为6个连续的近似平面,以保证能够正确地选择投影影像。实验中对一个具有60 000个三角面片的三角形网格模型分别在两种方法下形成的连续区域R进行统计分析,其结果见表 1。可看出,本文方法可将碎片化程度减小到逐个三角形纹理映射方法的5%,而这些碎片主要分布在植被不规则的表面、建筑物内凹形成的坑洞表面以及细小地物的表面等处。
表 1 区域生长前后碎片化程度对比Table 1. Comparison of Fragmentation Degrees of the One-by-one Approach and Our Approach方法 总连通区域个数 三角形个数<10的连通区域个数 三角形个数<5的连通区域个数 三角形个数=1的连通区域个数 平均包含的三角形个数 逐个三角形映射方法 11 755 10 766 9 640 5 748 5.1 一致性约束的区域生长方法 909 541 457 216 66.0 图 7为纹理映射后的对比结果。图 7(a)为Street Factory建模并进行纹理映射的结果。从图 7(a)可看出,虽然Street Factory对纹理进行了匀光处理,但仍存在纹理碎片化的现象;图 7(b)、7(c)分别为同一个三角形网格模型在逐个三角面片纹理映射和一致性约束纹理映射的结果,显然逐个三角面片的纹理映射结果受噪声影响较大,尤其是在水面和植被处的纹理碎片化现象严重;而一致性约束纹理映射结果则有显著改善。
3 结 语
本文提出了一种针对三角形网格模型的局部区域表面一致性约束的纹理映射方法,并使用真实数据与Street Factory等商业软件进行了对比实验。表面一致性约束可大大减小三角形网格模型表面噪声(尤其是法向量噪声)的影响,顾及这种约束的区域生长方法能够将三角形网格模型良好地分割为多个连续的区域,最终为每一个区域选择最优的投影影像提供纹理。这种化零为整的思想对消除纹理碎片化现象具有显著效果。进一步的研究将针对更加复杂的建筑物场景的真实感表达进行,以期消除不同建筑物间的纹理接缝,增强场景的真实感。
-
[1] Wen Jiahong,Kang Jiancheng,Sun Bo.RecentProgress in the Study on Antarctic Glaciology andGlobal Change[J].Chinese Journal of Polar Re-search,1997,9(3):232-237(温家洪,康建成,孙波.南极冰川学与全球变化研究新进展[J].极地研究,1997,9(3):232-237)[2] Sun Jiabing,Huo Dongmin,Sun Zhaohui.TheStudy of Remote Sensing on Monitoring Ice Velocityof the Polar Record Glacier and the Dalk Glacier[J].Chinese Journal of Polar Research,2001,13(2):117-128(孙家抦,霍东民,孙朝辉.极地记录冰川和达尔克冰川流速的遥感监测研究[J].极地研究,2001,13(2):117-128)[3] Strozzi T,Luckman A,Murray T,et al.GlacierMotion Estimation Using SAR Offset-tracking Pro-cedures[J].IEEE Transactions on Geoscience andRemote Sensing,2002,40(11):2384-2391[4] Wang Chao,Zhang Hong,Liu Zhi.SpaceborneSAR Interferometry[M].Beijing:Science Press,2002(王超,张红,刘智.星载合成孔径雷达干涉测量[M].北京:科学出版社,2002)[5] Cheng Xiao,Li Xiaowen,Shao Yun,et al.DInSARMeasurement of Glacier Motion in Antarctic GroveMountain[J].Chinese Science Bulletin,2006,51(17):2060-2067(程晓,李小文,邵芸,等.南极格罗夫山地区冰川运动规律DInSAR遥感研究[J].科学通报,2006,51(17):2060-2067)[6] Li Shuanke.On the Glaciation of the LarsemannHills,East Antarctic[J].Antarcitc Research(Chi-nese Edition),1995,7(4):7-16(李栓科.东南极拉斯曼丘陵区的冰川作用[J].南极研究(中文版),1995,7(4):7-16)[7] Bamber J L,Gomez-Dans J L,Griggs J A.A New1km Digital Elevation Model of the Antarctic De-rived from Combined Satellite Radar and Laser Data-Part 1:Data and Methods[J].The Cryosphere,2009,3:101-111[8] Griggs J A,Bamber J L.A New 1km Digital Ele-vation Model of Antarctica Derived from CombinedRadar and Laser Data Part 2:Validation and ErrorEstimates[J].The Cryosphere,2009,3:113-123[9] Scambos T A,Bohlander J A,Shuman C A,et al.Glacier Acceleration and Thinning After Ice ShelfCollapse in the Larsen B Embayment,Antarctica[J].Geophysical Research Letters,2004,31:18 402[10]Nakamura K,Doi K,Shibuya K.Estimation ofSeasonal Changes in the Flow of Shirase Glacier U-sing JERS-1/SAR Image Correlation[J].Polar Sci-ence,2007,1(2):73-83[11]Rott H,Rack W,Nagler T,et al.Climatically In-duced Retreat and Collapse of Northern Larsen IceShelf,Antarctic Peninsula[J].Annals of Glaciolo-gy,1998,27:86-92[12]Rignot E,Vaughan D G,Schmeltz M,et al.Accel-eration of Pine Island and Thwaites Glaciers,WestAntarctica[J].Annals of Glaciology,2002,34(1):189-194[13]Rignot E,Casassa G,Gogineni P,et al.Acceler-ated Ice Discharge from the Antarctic Peninsula Fol-lowing the Collapse of Larsen B Ice Shelf[J].Geo-physical Research Letters,2004,31:18 401[14]Rignot E.Changes in West AntarcticIce Stream Dy-namics Observed with ALOS PALSAR Data[J].Geophysical Research Letters,2008,35:12 505[15]Khazendar A,Rignot E,Larour E.Acceleration943武 汉 大 学 学 报 · 信 息 科 学 版2014年8月and Spatial Rheology of Larsen C Ice Shelf,Antarc-tic Peninsula[J].Geophysical Research Letters,2011,38:9 502[16]King M A,Coleman R,Morgan P J,et al.VelocityChange of the Amery Ice Shelf,East Antarctica,During the Period 1968-1999[J].Journal of Geo-physical Research,2007,112:1 013[17]Pittard M L,Roberts J L,Warner R C,et al.Flowof the Amery Ice Shelf and Its Tributary Glaciers[C].The 18th Australasian Fluid Mechanics Con-ference,Australia,2012[18]Zhou Chunxia,Zhou Yu,E Dongchen,et al.Sea-sonal and Inter-annual Changes Detection in Ice Ve-locity of Antarctic Inland and Outlet Glaciers[C].International Symposium on Glaciers and Ice Sheetsin a Warming Climate,Alaska,USA,2012[19]International Sea Level Indian Ocean Data[OL].ht-tps://www.bodc.ac.uk/data/online_delivery/in-ternational_sea_level/indian_ocean/,2012[20]Wang Qinghua.GPS-derived Horizontal Ice FlowVelocities Along the Traverse Route from Zhongs-han Station to Dome-A,East Antarctica[J].Chi-nese Journal of Polar Science,2003,14(1):35-40[21]Manson R,Coleman R,Morgan P,et al.Ice Veloc-ities of the Lambert Glacier from Static GPS Obser-vations[J].Earth Planets Space,2000,52:1031-1036[22]Zhang Shengkai,E Dongchen,Wang Zemin,et al.Ice Velocity from Static GPS Observations Alongthe Transect from Zhongshan Station to Dome A,East Antarctica[J].Annals of Glaciology,2008,48(1):113-118[23]Michel R,Rignot E.Flow of Glacier Moreno,Ar-gentina,from Repeat-pass Shuttle Imaging RadarImages:Comparison of the Phase Correlation Meth-od with Radar Interferometry[J].Journal of Glac-iology,1999,45(149):93-100[24]Joughin I.Ice-sheet Velocity Mapping:A CombinedInterferometric and Speckle-tracking Approach[J].Annals of Glaciology,2002,34(1):195-201[25]Dong Xiaoyan,Ding Xiaoli,Li Zhiwei,et al.ANew Chain for SAR Pixel Offset Estimation and ItsApplication in Co-seismic Deformation Measurement[J].Geomatics and Information Science of WuhanUniversity,2011,36(7):789-792(董晓燕,丁晓利,李志伟,等.一种新的SAR像素偏移量估计流程及其在同震形变监测中的应用[J].武汉大学学报·信息科学版,2011,36(7):789-792)[26]Bechor N B D,Zebker H A.Measuring Two-di-mensional Movements Using a Single InSAR Pair[J].Geophysical Research Letters,2006 -
期刊类型引用(10)
1. 刘飞鹏. 消费型无人机在河道精细化三维构建中的应用. 测绘与空间地理信息. 2024(03): 29-31+35 . 百度学术
2. 周友华,邢小平. Multigen Creator几何建模中的纹理优化. 电脑编程技巧与维护. 2021(01): 145-146+171 . 百度学术
3. 朱庆,张琳琳,胡翰,翁其强,丁雨淋,李赟,张叶廷. 精细建筑物碎片化纹理优化的二维装箱方法. 西南交通大学学报. 2021(02): 306-313 . 百度学术
4. 王梦威,马秀丽. 基于中介媒介的交互式纹理映射. 电子测量技术. 2020(12): 115-120 . 百度学术
5. 叶震,许强,刘谦,董秀军,王晓晨,宁浩. 无人机倾斜摄影测量在边坡岩体结构面调查中的应用. 武汉大学学报(信息科学版). 2020(11): 1739-1746 . 百度学术
6. 马立华,郑永虎,丛晓明,冷旻航. 无人机倾斜摄影测量三维仿真飞行在龙羊峡水电站高陡边坡数据获取中的应用. 青海科技. 2020(06): 72-76 . 百度学术
7. 乃古色拉,张云生,张明磊,邹峥嵘. 基于图割算法的倾斜影像纹理映射优化方法. 测绘与空间地理信息. 2019(03): 145-147 . 百度学术
8. 朱庆,翁其强,胡翰,王峰,王伟玺,杨卫军,张鹏程. 基于帧缓存的多角度影像精细纹理映射方法. 西南交通大学学报. 2019(02): 269-277 . 百度学术
9. 黄明,贾嘉楠,李闪磊,张建广,龚建辉. 多像位姿估计的全景纹理映射算法. 武汉大学学报(信息科学版). 2019(11): 1622-1632 . 百度学术
10. 梁玉斌,崔铁军. 倾斜摄影测量的研究进展. 天津师范大学学报(自然科学版). 2017(05): 1-6 . 百度学术
其他类型引用(11)
计量
- 文章访问数: 4151
- HTML全文浏览量: 78
- PDF下载量: 769
- 被引次数: 21