基于河流水文学使用过程性方法地形建模

'‘’１－擊讀ｈ讓．雜穀耀ｉ露琴，疋＾山．．ｒ；－．顆画壤娘，草動＼；．每＇：調？＾．巧：：ｖ过ＶＶ减ｆ：－巧粗背吉琼却．；拟；；冷书；；心‘、；－．．．学校代码ｔ１０２００，振Ｈ社请争学校名祝东化师范大学韓伊可巧：１跨－、分炎笠型研充生学新繼２１亂：专誦蘇韓Ｖ變：藻眞漏饒与＇．＇：一＇．＇．＇．ｕ＇，；ｒ；挈；，，；Ｖ祥鱗巧点古城巧讚艇較：．＇．議職甚寶遺攀Ｖ；三；ｓ雜＾ｉ．ｉｌｉｉｌ＃＾Ｐ嗎－？－’’、＇吗錢響巧鑛魏纔讀總總：、＇．—＇．公Ｊｔ冬‘＂？；：：．、‘．、；，．柏护知卢：；：：ｒ：，．的巧＞巧料．．茲，：；ｈ，？．没ｆ；齡‘糾某；興＃於作化Ｈ：‘＇．．‘产／觀巧礙巧義揖辦菊若游舅—二丄．．．、；＾，－：．ｖ最省巧痛新．：说巧貧觀硕±学位论文，语雜：．．．＇■．，＇－■一七－＾＇斗的＾１，＾：听＇一；Ｖ今苗Ｔ＇＾与ＶＶ／＇＾．－：…基觀等知占；七攀識韻舉捉設冷：＇遷寶寒和ｒ：．；基于河流水文学使用挂程性方法？建模；、、笔為／成占、Ｂａｓｅｄ０１１ｒｉｖｅｒｈｙｄｒｏｌｏｇｙｕｓｉｎｇｒｏｃ的ｕｒａＬ、占苗ｐ：－．，．．．作Ｖ一户：，导觀如批Ｓ論ｅＳＩ端ｌ為．．＇：’ｍｄｌｉｎｅｄ节ｅｒｒｉｎｏｅｍｔｈｏ．、，ｔａ：‘：．知：．＇ｇ．、拍帯沾、．葬纖讀４妒，＾｜？．．＇山；、．：－－＇■Ｘ：．－Ｙ一、緣巧、、点，接、—‘＇＇＇Ｉ‘—、＾’．寸、．、．？：：卢或：、、幕．絮气－？？．＇．／；％：知－．游：占一．，，．、巧＞：；ｃ證Ｖ．？＇．＇＇－＂去．－＇．妒常Ｖ：山龄．扣．：ｒ＼啦ｆ墙错苗與；．？、？‘、？、‘－：｜一片－—：，：ｒ护矣１伯—．＞端嶺、－’、’Ｖｖ＇ｖ；导师；掩义沪．端教李，，Ｖ公？馨蘇新一‘、－、与，‘；刑＊刪学技术过户去级科．，．：４如六学．、二＿，二、—芭巧一如琪：心．沾茹取学计应用茜结掉二级科＊机沁聚斯捉、＾’爲與＇一’、’．？．，方：：．；ｉｍ巧形学ｒ共托巧巧向计辦准，，巧；＇，薄、＇＇、Ｉ中．＇；义。’；．．．ｖ舞巧片：一学类型；额丈；．－．某位学术Ｊ；沁片：；Ｃ啼．Ｖｆ４聲－－－：：補巧昏雜巧煤站歡馈记ｒ細裙蠢黎掉瓣装爲識海離聲疆喪貧：＇、－＾．．，．纖議議纖操苗實賴．＇賴獨辫哪養驚擊，攀禱？；式：’梧嗦巧寒：Ｖ梯《护巧点梓咕藏泌乾辩驴川听Ｊ纖馨毅、，．’＾－＾：，這藥讓＃觀ｋＶ疑接編挺’‘＇、＇－定‘：於大委员会明、己、东北师龜学学位评又ｖ＼！．誤Ｖ；；、户．．．．＇私．、ｒ；；．、．＾２〇１５甘：Ｖ６月．．鑛觀起嫌篇寅Ｖ益年＾＼纖聲講：护祭争＇’；－砖者；纖辨滿擁？知嫂巧；、，＾？端：嚇巧；捧奇＾强齡羣护＾＇拉、＇｜、醉六：．祐；；巧：Ｌ殘；賊巧粗；．每為 摘要针对传统地形合成方法缺乏可控性，同时合成后的地形不能满足与真实地质地貌相一致等问题，本论文提出了基于使用过程性建模的方法。为了克服过程性建模方法生成的地形模型不符合真实地貌形态的缺点，引入了水文学理论，基于Tokunaga网络的河网形态特征量化结果，进一步给出了基于L系统的河网过程性建模方法，通过给出不同的河网生成器，就能构建出符合地质地貌大体相一致的树状河网或羽状河网。在此基础上提出了构造支持地形山脊的启发式方法。在得到河网和山脊的三维地形骨架后，以此作为特征就可以采用数学的方法拟补出山体所需要的DEM数据，这样就可以模拟出符合真实情况的地形，由于这部分耗时比较长，为了满足用户实时要求的目的，引入了GPU技术来加速山体生成过程。同时利用泊松方程引入离散指导矢量场到参数化网格中，使得虚拟地形DEM可以与真实地形DEM进行无缝融合，并且用户可以根据需求进行局部化的网格编辑，对生成地形的高程产生影响，进而满足用户的需求。最后对25m分辨率的地貌DEM数据进行了实验研究和数据分析对比，结果表明该方法结合了物理建模和过程性建模的优势，具有良好的可控性，支持大规模地形可视化，引入GPU加速后即可达到快速实时目的，同时合成地形具有真实的地形地貌。关键词：过程性建模；河流网络；泊松方程；地形可视化；水文学I AbstractLackofcontrollabilityofthetraditionalterrainsynthesismethods,butcannotmeettheterrainsynthesizedcoincidewiththerealgeologicalfeaturesandotherproblems,thispaperpresentsamethodbasedontheuseofproceduralmodeling.Toovercometheshortcomingsoftheprocessofmodelingterrainmodelgenerationmethoddoesnotmatchtheactualtopography,hydrologytheoryisintroduced,thequantizationofrivernetworkmorphologybasedontheTokunaganetwork,andfurthergivestheriversystembasedonprocessmodelingL-systemthemethod,bygivingadifferentriverBuilder,youcanbuildsubstantiallyinlinewithgeologicalfeaturesconsistenttreeriverorriverplumes.Onthisbasis,theproposedconstructionofterrainridgesheuristics.Afterobtainingthree-dimensionalterrainskeletonriverandridges,asafeatureyoucanusemathematicalmethodsintendedtocomplementthebodycomingouttherequiredDEMdata,sothatyoucansimulatetheterraininlinewiththerealsituation,becausethisparttakeslongerinordertomeettheobjectivesoftheuserreal-timerequirements,theintroductionofGPUtechnologytoacceleratetheprocessofgeneratingtheterrain.WhiletakingadvantageoftheintroductionofthePoissonequationdiscretevectorfieldtoguidetheparametricgrid,SothatthevirtualandtherealterrainDEMseamlessintegration,andtheusercaneditthelocalizedgridaccordingtodemand,thegeneratedterrainelevationgenerationinfluence,andthusmeettheneedsofusers.Finally,the25mresolutionDEMtopographydataofexperimentalresearchanddataanalysisandcomparison,theresultsshowthatthemethodcombinestheadvantagesofphysicalmodelingandprocessmodeling,withgoodcontrollability,supportlarge-scaleterrainvisualization,theintroductionofGPUtoachieveacceleratedfastrealpurpose,whilesyntheticterrainwithrealistictopography.Keywords:proceduralmodeling;rivernet;poissonequation;terrainvisualization;hydrologyII 目录摘要......................................................................IAbstract.................................................................II目录....................................................................III第一章引言................................................................11.1研究意义...........................................................11.2国内外研究现状.....................................................11.3本论文的主要工作...................................................3第二章基于河流水文学的河网构建............................................42.1基于Tokunaga网络的河网形态特征量化................................42.2基于L-System的河网过程性建模......................................52.3基于河流分类的高程计算.............................................6第三章山脊的构造..........................................................83.1Voronoi图.........................................................83.2基于中点位移逆过程的数据扩展.......................................93.3基于Laplace的扩散方程............................................10第四章地形块的嵌入以及变形处理...........................................124.1泊松方程.........................................................124.2地形融合与形变....................................................124.3优化..............................................................14第五章总结与展望.........................................................16参考文献.................................................................19致谢.....................................................................21III 第一章引言1.1研究意义虚拟场景在医药医疗、地质分析、游戏制作、电影创作等领域有着广泛的应用。虚拟地形是核心组成部分之一，也是场景中占主导地位的视觉元素之一。近年来，研究人员已经取得了很好的研究进展，提出了很多有效的地形合成方法。现有技术可大致分为基于过程性的建模方法、基于物理的建模方法以及基于草图或实例的建模方法。然而，过程性方法生成的地形模型通常缺乏可控制性，该类算法通常在大规模地形可视化场景中表现出规则特征，不符合真实地形形态的缺点。基于草图的方法涉及到手工编辑，通常令人感觉枯燥乏味，并且效率较低。基于实例的算法将会受到所提供的输入样本的限制。此外，从地质学的观点，基于物理学算法所产生的虚拟地形是符合真实地貌的，但是与过程性方法存在类似的问题，即缺乏可控性，同时无法完成大规模地形可视化的任务。目前地形生成最重要的问题是缺乏具有可控性，同时能够快速生成包含大规模的地形可视化场景，并且符合地质学物理要求的地形生成算法。数字高程模型(DigitalElevationModel)，简称DEM是大多数地形建模系统的通用数据格式，其优势在于表示方法简单，同时能够快速实现可视化效果，因此大部分方法生成的地形模型通常以DEM形式存储，然而地形通常呈现单一特征，不符合真实地形地貌。我们提出了基于水文学的过程性地形建模方法，实现了虚拟地形生成系统。在本系统的地形生成过程中允许用户交互性控制地形形态，提供良好的可控性，同时支持具有大规模集的地形可视化。我们首先从基础水文网络进行衍生，它保证所构建的地形满足水文学特性。该项工作能帮助地质专家进行辅助决策，同时在虚拟自然环境和虚拟仿真系统中基于大规模集的三维地形也是必备的元素和重要组成部分。1.2国内外研究现状过程性建模是有一种经典的地形建模技术，其原理简单、易于操作，仅需要控制几个参数的变化，便能够生成大范围地形。其中最具代表性的算法是Fournier等人提出的[1]自适应细分方法，它提供了固有的层次细节模型。Perlin等人提出了基于噪声的过程[2]性方法，该方法结合噪声函数在不同规模下提供了不同的细节。Ebert等人对该类方法[3]进行了较为详细的综述。此外，基于分形的方法也通常善于产生大规模地形细节，但这类方法通常缺乏基于地形特征的空间信息以及其他控制信息，很难准确表达河流、山谷等地形地貌。因此他们生成的地形比较生硬、不符合客观实际。实际上，真实的地形通常是通过侵蚀和风化的影响而产生的，一种合理的解决方案是将河网信息引入建模过程。Kelley等人提出了一种过程性方法来生成流域(分水岭，转[4]折点)，采用这种方法可以生成局部草图作为全局地形的输入数据。Prusinkiewicz等人1 [5]结合上下文相关的L-System，引入了中点位移法，将河流的形态使用分形建模系统来[6]嵌入到地形中。2005年，Belhadj和Audibert提出了一种改进的随机细分算法，该方法[7]运用制约分形布朗运动产生了山脊和河流的曲线。2009年，Teoh.等人提出了一种通[8]过生产河网产生地形的相关算法。同样，2011年，Derzapf等人产生一个超大尺度的[9]河流网络。2014年，Tasse等人根据知觉线索，找到了轮廓草图与输入地形的轮廓或者山脊最佳匹配。随后变形算法将被应用到地形中，从给定的透视视角出发使其能够完全匹配的草图，同时确保没有用户定义的轮廓是由地形的另一部分所遮挡的。这些方法虽然已经取得了较好的研究进展，但是存在以下缺点：首先这些算法提供的河网和流域是不相关的；其次河流的路径通常采用随机算法创建，所以在地貌学观察中不符合河流的特点；最后这些算法交互性不好，仅允许用户控制抽象的参数。基于物理的建模方法的基础是建立在地表被暴露于各种自然介质中，例如水的侵蚀、[10]温度的变化以及人类的活动。其中Musgrave等人最早发表了相关论文，引入了简单[11][12][13]的水力和热侵蚀的侵蚀沉积模型。这种方法已经扩展到许多不同的研究领域，他们提供了不同的热力和水力侵蚀算法。虽然这些扩展方法大多数使用了正则高度场，[12]但是能够表示成带有侵蚀模型的分层数据结构，并且在后续研究中推出一个完整的三[13][14]维立体水力侵蚀模型。近年来，Kristof等人引入了结合侵蚀方法的光滑粒子流体[15]力学。Wojtan等人引入了腐蚀模拟形态学算法，其主要缺点是可控性较低。这些方法还不能对具有高分辨率细节的大规模地形建模。当然，面向GPU的物理建模方法能够[16][17]很大程度上提高建模效率，但是大规模地形物理建模的核心效率问题仍然不能解决。为了使构建的地形具有良好的可控性，很多研究者采用交互式编辑的建模方法。[18][19]Rusnell等人采用了基于特征的生成技术。Zhou等人用二维高度场作为样本实例，最终根据特征对地形样本进行合成。该方法的合成效果很好，但是必须以输入地形为基[20]础，这是基于模板的虚拟地形合成算法的缺点。交互式地形编辑和骨架特征提取方法对合成地形具有良好的控制能力，但从地质学角度来说他们可能会导致一些不正确的结果。[21]此外，很多学者尝试了混合建模方法——结合交互式编辑和基于物理学的算法，[22]而该类方法被限于编辑现有的地形，并且只适用于小场景。Hnaidi等人介绍了一种使用基于向量的模型算法直接操作河流轨迹。地形生成是根据给定的一组代表地貌特征，如脊线、悬崖和河床的控制曲线，采用多栅型扩散方程将高程和梯度约束匹配到对应的[23]控制曲线。Yuzhou等人将一种以泊松方程为理论基础的方法用于网格编辑。这种方法的最显着特点是，它通过隐式操纵梯度场修改了原来的网格几何体。该方法的优势在于能够实现全局和局部的编辑操作，如变形、物体合并以及平滑。虽然这两种方法都利用网格的差分性质，并且能够生成大型高分辨率地形，但是它并没有解决河流网络建模中存在的地貌表达问题。由于T样条的出现，使得大量研究者对曲面表示及其细化方面进行了大量研究，其中VibekeSkytt等人使用新颖的曲面表示方法，即LRB-splines曲[24][25]面表示法，并应用这种表示在地理数据的处理中。这些数据集通常非常巨大，并且2 LRB-splines曲面能提供具有局部细节调整机制和整体平滑数据的紧凑表示，但他只是近似逼近方法，对地形的在构建和快速渲染方面提出较好理论基础。在地质河网研究方面，河网形态反映地表水系的形态结构，并影响着河流演变、水[26]文过程、流域地貌、土壤侵蚀等。在Horton-Strahler河网分级的层次上，他们归纳出带有统计平均意义上的Tokunaga分形河网描述方法。1.3本论文的主要工作本文结合水文学理论，利用过程性建模方法的优势，提出了基于具有良好可控性的地形可视化建模方法。本文的主要贡献如下：提出了利用河流作为建模特征实现过程性地形建模方法，优势在于遵循水文学理论，兼顾过程性建模方法的优点，生成的地形具有真实地貌特征，同时满足场景对大规模可视化的需求；给出了生成的河网的规则集，满足景观建模和河流科学观察的精度要求，并且还具有较高的可控性；利用网格本身的差分属性解决了虚拟地形和现有实际DEM地形的无缝融合。研究流程的DFD数据流图如下：3 第二章基于河流水文学的河网构建Horton-Strahler河网层次分级方法是采用以下参数来定量描述河网，包括采用河流[27]长度比、分支比、面积比以及Hack定理。在统计意义上，随着引入自相似的概念，研究者们利用链和分形的概念来描述河网自相似网络的拓扑性质。河网的自相似分形通常有两种方法：一种方法是根据Horton水系定理来确定的分形维数公式，该维数是Horton-Strahler河网分级中的分支和河流长度的关系比，其描述了河网的拓扑结构分形维数特征；第二种方法是探究利用链数和链长来表达河网分形维数。在本文中我们采用了第二种方法，基于河网的平面填充特性，它的分形维数为2，因此在构建河网时必须要满足这一特性。我们结合ZHANG等给出的树状河网的形态特征、自相似河网的基本生成器序列以及生成Tokunaga河网的生成参数，我们提出了一个河网分形规则集，来产生与真实自然景观大体一致的河网。2.1基于Tokunaga网络的河网形态特征量化根据分形的基本定义，我们通过链数和链长来计算分形维数，可以进一步确定自相似河网的拓扑结构和数值特征是可以采用链和分形维数来描述的。分形维数能够反映河网的发育程度，并且给出了能够满足不同程度平面填充性质的一系列基于Tokunaga网络的自相似河网。其优势在于可以构造出与自然景观实际河网几何形态相符的自相似河网形态，对描述河网形态特征的参数ߣ进行量化，并以实际自然河网为模式输入进行计算，发现量化参数ߣ能够较好地描述河网形态。其河网规律符合约束条件：ܶ=ܽܿ௞ିଵ(1)௞其中ܶ௞表示侧支分叉比，是在统计意义上汇入߱级支流中的߱−݇级侧支流的数量。基于Tokunaga的定义，ܽ描述了流入߱级支流的߱−1级侧支流数量，ܿ描述了低级侧支流的数目统计增长趋势。ܽ和ܿ用统计均值表示。同时河网本身具有树状结构的特点，所以一般可以抽象为二叉树，因而我们定义了以二叉树为基础的生成器，为了形成一个有规律性的生成器序列，需要其按照生成器中链数的增加对其排列。表1是4个典型的以二叉树为基础的自相似分形生成器，箭头指向河流出口方向。设λ是生成器序号，ݐ表示迭代次数，当ݐ=1时表示生成器。当前河网的Horton-Strahler级别表示为ߗ，其中当前链或河流的Horton-Strahler级别表示为߱，当前河网中链的总数为ܯ，其中݅表示内链，݋表示外链。以二叉树为基础的自相似分形生成器的数学描述为：ఒߗଵ=2(2)ఒ且ܯଵ=2ߣ+1(3)根据二叉树的基本特点：ܯ௢=ܯ௜+1(4)由(2)式可得4 ఒܯ௜=ߣ，ఒܯ଴=ߣ+1(5)ଵଵ其中参数ߣ是生成序列，构造后生成器序列的外链以生成器ߣ代替、内链以生成器ܲݐ∗(ߣ−1)代替时可以生成分形维数为2的自相似河网。不难发现，这个迭代规则只是[26]在统计意义上有效，在多数情况下ܲݐ不是整数。我们采用了ZHANG提出的在迭代过程中链长的变化规律，如公式(6)所示：ఒƐ=(1/ܲݐ·ߣ)௧ିଵ(6)௧而平面填充河流ܽ和ܿ的值为ܽ=ߣ–1(7)c=Pt·λ=(1+β)/2(8)ଵାఉ其中ܲݐ=，ߚ=√1+4ߣ。ଶఒ表1反映ߣ与ܲݐ的关系Table1therelationshipofߣandܲݐλ1234生成器ࡼ࢚——1.000.7680.640ࡼ࢚·（ࣅ–૚）——1.001.5351.9212.2基于L-System的河网过程性建模过程建模是一类方法的总称，是一定数量的计算机图形学技术创建3D模型和纹理的规则集。分形L-System是程序建模技术，主要是运用算法及产生式规则来生成场景。该规则可能嵌入到算法中，或通过可配置的参数给出。规则的输出我们通常称为过程内容，它可以在计算机游戏、电影、互联网中使用，或者也可以是用户手动编辑的内容。过程模型在表现大规模场景时可以由许多更小的规则来生成。如果所采用的算法产生相同的输出，每次输出不需要被存储。通常，只需用相同的随机种子来实现这个启动算法。大部分建模技术都需要算法来管理和存储数据，过程建模的不同之处在于根据规则集创建一个模型，而不是通过用户输入或者编辑模型。当然，过程建模通常过于繁琐，所以在构建植物、建筑、自然景观时经常采用3D建模工具进行辅助，本文主要实现虚拟地形合成模型，因此建模复杂度适中，不必考虑以上问题。根据相似和分形理论，为了使生成的模型更加符合自然形态，很多学者采用了不同的生成器来生成自相似分形河网。为了可以分析某些极限自相似网络的性质，生成器序列允许λ→∞。使用这个参数的目的是在分形的迭代过程中可以采用简单方便的数学公式进行表达。我们给出了以下规则集对河网进行可视化，则这个规则集对生成器的量化5 为：对于ߣ=2情况外链：G→F[--G]F[-G][+G]内链：F→F[++G]F或F→FF对于λ=3情况外链：G→F[++G]F[--G]F[-G][+G]内链：F→F[--G]F[-G][+G]或F→F[++G]F或F→FF其中F是画一个长度为݈݁݊݃ݐℎ的线，左中括号代表入栈，右中括号代表出栈，“+”和“-”代表旋转角度，对于ߣ>3情况以此类推，如图1所示河网效果图。(a)ߣ=2、݀݁݌ݐh=8的河网生成(b)画四条线矢量线的情况图1生成器产生的河网效果图Fig.1Rivernetworkproducedbygenerator在实际的河网模拟过程中，通过分析用户输入的矢量线，得到起始点和方向，然后设置λ参数就可以得到符合水系的河网：ߣ的值越大，河网呈现羽毛状；ߣ的值越小河网呈现树枝状。为了使得我们生成的自相似河网更具一般性，在我们的迭代实验过程中使用一个完整的生成器来代替外链，而用ܲݐ个生成器来替代内链，进而用不同的生成器和ܲݐ来决定内链和外链的迭代方式，从而对河网的模拟更符合自然特性和河网的平面填充性质。2.3基于河流分类的高程计算[28]我们采用了经典的Rosgen分类河流分类方法，并根据他们的斜坡和轨迹定义了九种河流的类别，来设定河流曲线的类型与坡度。每条河流有一个轨迹类型(A+,A,B,C,D,DA,E,F,G)和河床的坡度。同时包括河床的地质组成：如基岩、岩石、石头、碎石、砂、淤泥或粘土。我们结合河流的斜率和其邻近的海岸来分配每个节点的河流分类。我们采用测地线距离作为邻近海岸的评价标准，它被标记为辫状河。例如，那些由相隔条带所分割的多个通道，记作D或DA，如表2所示。对于每一条河流，边的类型是根据它的两个节点进行计算，并且结点类型是通过一个查找表来确定的。6 表2河流类型与坡度的关系表Table2Rivertypecorrespondingslope坡度>10%4-10%2-4%<2%<4%<0.5%<2%<2%2-4%类型Aa+ABCDDAEFG我们的主要思路是根据河网坡度的斜率以及链的长度，当河流每扩展一步就按照其斜率算出高程（初始高程值为零）。本文描述的过程性模型的交互性体现在用户在空DEM上，进行画简单的有向直线段并指出该段的河流类型，就可以生成逼真的河网。设定Lsystem的具体公式如下：Lsystem.start=Lab.start;Lsystem.direction=Lab.direction;Lsystem.depth=Lab.length/DEM.width*MaxDepth(9)其中Lab.start、Lab.direction、Lab.length分别是线的起始点、方向以及长度，Maxlength是Lsystem能迭代的最大深度，确定Lsystem河流盆地的初始位置(Lsystem.start)，增长方向(Lsystem.direction)以及迭代层数(Lsystem.depth)。7 第三章山脊的构造为了获得更多的形态特征以构建真实地形，我们进一步基于已有的河网产生地形的山脊特征。首先，我们基于河网的节点建立Voronoi图，从而保留并离散化不同流域间的Voronoi图边作为分水岭。接下来，中点位移法的逆过程被用来生成更多的离散高程点，进而更精确地控制山脊形态特征。也能够加速山体的构建。最后，我们采用对已生成的高程点使用基于Laplace的扩散方程进行扩散，从而得到符合地貌特征的全局高度场。3.1Voronoi图在地理信息系统中Voronoi图是解决邻接度问题的常用方法，经常采用它进行快速邻近插值和地理实体影响区域的分析。Voronoi图的基本特性是：每个Voronoi单元ܸ௜内只含有一个数据点ܲ௜；ܸ௜内的点到相应离散点ܲ௜的距离最近；位于ܸ௜边上的点到其两边的离散点ܲ௜和ܲ௝的距离是相等的。很多研究者利用这些特点将他们应用于定性分析、统计分析、邻近分析等。在这里我们将利用这一特性，并根据前面小节介绍的河网结点，来[29]产生地形山脊。在Voronoi图的构建中，我们采用了经典的扫描线（SweepLine）算法。如图2所示。图2根据分形河网中的结点建立VoronoiFig.2BuildVoronoidiagrambasedonfractalriverjunction在建立Voronoi图的基础上，我们将河网与Voronoi图中相交的边去掉，从而形成[30]山脊线。根据河网结点的高程值，我们采用公式(10)计算山脊线上结点的高程值，具体公式如下：ݍ.ݔܽ݉=ݖ(݌ଵ,݌ଶ,݌ଷ)+ߣ·݀(10)其中ݍ是Voronoi结点，݌ଵ,݌ଶ,݌ଷ是距离ݍ最近的三个河网结点，ߣ是可调节高程参数，݀是ݍ到݌ଵ的距离。在我们的实验中，最终形成的三维骨架效果图，如图3所示。8 图3分形产生河网与Voronoi图合并后产生的地形三维骨架Fig.3Three-dimensionalterrainskeletongeneratedbymergingfractalrivernetworkandVoronoidiagram在生成山脊后，我们对其进行离散化，采用了经典的Bresenham方法对平面山脊线[31]段进行离散化参考。不同之处在于，我们在离散化的同时，为了使山脊线高程符合真[2]实自然的情况，离散点的高程值是被随机加上的符合Perlin噪声的随机值，具体公式如下：P-PsP.z=⋅Pe.z+F(P),(11)Pe-Ps其中ܲ点是插值点，ܲ௦和ܲ௘点是山脊线的起始点和终点，ܨ是perlin分形噪声函数。3.2基于中点位移逆过程的数据扩展[6]我们采用的方法与Belhadj等人提出的M.D.I.(MidpointDisplacementInverse)过程相似，该过程是建立在简单的但具有强大事实观测基础上的。当我们在处理一个给定的正方形的中点位移细分时，正方形中点（子节点）使用其边的角点的高程进行插值，所以用来插值为中点的角点（父节点）可被存储在一个列表中。如图4所示，虚线箭头是中点位移法的逆过程，实箭头是对网格实行中点位移法。为了丰富表达山体数据，我们使用基于中点位移法逆过程，生成更多的山脊数据，这为后面的扩散提供更加真实的保障，因此，对于一个给定的正方形面片ܵ（正方形根节点类比为的D.E.M.）我们预先递归地计算对父列表中每个中点。我们称parent(x，z)为点(ݖ,ݔ(到)ݖ݌,ݔ݌(和)ݖ݌,ݔ݌(标坐父是)݀,ݖ݌,ݔ݌(组元三的含包它，表列父的)ݖ,ݔ)的欧氏距离݀。如图4（三角形-边中点位移法，菱形-正方形中点位移法）空心圆代表待插值数据，实心方框代表已有数据。因此，在逆过程中，如果状态为空，我们可以使用其可用的子节点的高程值进行插值。随着层次的升高，误差也随之增大。为了减少突兀的点[6]出现，我们采用了与文献类似的方法，但是不同之处在于我们对算法引入参数ܮ，当݀>ܮ时，我们将停止使用子节点对父节点的插值。最终的效果图如图5所示，原始图的分辨率是1024×1024。9 M.D.I(a)中点位移法的逆过程M.D(b)对网格实行中点位移法图4中点位移法正逆过程示意图Fig.4Schematicdiagramofpositiveandreverseprocessofmidpointdisplacementmethod(a)逆过程前（25900个数据点）(b)逆过程后效果图（57147个数据点）图5中点位移法逆过程前后对比图Fig.5Comparisonaboutreverseprocessesofmidpointdisplacementmethod3.3基于Laplace的扩散方程[32]扩散的方程的物理模型的基本原理是考虑三维空间中材质均匀的、各向同性的物体，这个物体与环境介质具有热传递关系，进一步来研究物体内部温度的分布规律。我们的方法不同之处是原来的方法主要应用在热扩散研究领域，在这里我们引入该方法应用于地形数据的高程差扩散。假定它的内部非空白的高程点作为扩散源，研究最终扩散出地形高程在时刻t的分布规律。其数学模型如下：∂u=DΔu+f(x,y,z,t)(12)∂t10 其中ݐ,ݖ,ݕ,ݔ(݂，程高形地的处)ݖ,ݕ,ݔ(在刻时ݐ示表)ݐ,ݖ,ݕ,ݔ(ݑ)表示单位时间内的非空白高程数据点的数量，ܦ为扩散系数。离散化的拉普拉斯算子为：Δ=4ܷ୧,୨−ܷ୧ାଵ,୨−ܷ୧ିଵ,୨−ܷ୧,୨ିଵ−ܷ୧,୨ାଵ(13)其中ܷ୧,୨是离散网格上的数据高程点。我们将三维骨架作为扩散源，其值是与参考平面的[33]高程差。在扩散之后可能存在突兀的点，我们采用了高斯平滑方法进行去噪。在本文中对于非常平坦的地势，同时可以进行了加噪处理。11 第四章地形块的嵌入以及变形处理4.1泊松方程[23]根据Yu等人提出的离散泊松方程网格处理方法，我们将目标网格的坐标定义为标量场，而每个标量场被定义为一个参数化网格。由于网格是分段线性模型，所以这种标量场实际上是分段线性的，并满足离散场的定义。同时，目标和参数化网格应该具有相同的拓扑结构（顶点连通性），并且它们的顶点应该有一对一的对应关系。应用泊松方程的目的是为了在已知的拓扑结构，但未知顶点坐标的情况下解决未知目标网格。为了获得未知的顶点坐标，泊松方程引入离散指导矢量场到参数化网格中。散度ܦ݅ݒ߱的计算公式如下：(ܦ݅ݒ߱)(݌௜)=∑்ೖఢே(௜)ߘܤ௜௞∙߱|ܶ௞|(14)其中ܦ݅ݒ是散度算子，ܤ௜௞是阶梯函数，ܰ(݅)是݌௜点的邻域，ܶ௞是以݌௜点为邻域的三角形[34]面片。一旦散度被定义，则离散的泊松方程计算公式如下：ܦ݅ݒ(∇Φ)=ܦ݅ݒ߱(15)实际上是一个稀疏线性矩阵：ܣݔ=ܾ(16)在公式15中系数矩阵ܣ是独立于指导域，并且只能使用参数化网格进行获得，ݔ是未知顶点，向量ܾ是从所有顶点散度值的集合中获得的。最后对产生的线性系统进行求解，得到所有顶点的特定坐标，从而获得所有顶点的坐标。整个过程类似于网格克隆，对指导域进行编辑是为了获得目标网格的期望。4.2地形融合与形变本文的实验系统提供了从构建的虚拟地形到现有的实际地形的融合操作。我们选取的由河网构建出来地形为边界约束条件，将利用实际地形数据求得的散度作为输入，从而合成最终的虚拟地形。在网格融合中，原始网格是被给定的，其目标是获得一个融合后的网格。因此，我们很容易将原始网格作为参数化网格。新的梯度矢量定义了参数化网格新的向量场。由于局部变换应用到每个网格可能会有所不同，因此原始网格被分开。同理，新的矢量场不太可能成为任何标量函数的梯度场。为了重建这些向量场网格，我们需要考虑他们作为泊松方程的指导矢量场。给定一个参数化网格和指导矢量场，遵循上一节中的步骤以重构新目标网格。直观上，求解泊松方程，类似于前面将断开的网格拼接在一起的工作。可交互地修改网格梯度场是我们的网格编辑系统的一个关键组成部分。我们的方法优势在于针对当前的主流的DEM格式的地形，交互式选取用户满意的地形，进行快速编辑。12 算法1.地形融合算法输入：源地形块和目标地形块的边界输出：融合结果Step1.选取源地形块。Step2.计算源地形块的散度。Step3.获取目标地形块的边界，作为源地形块边界的初始值。Step4.执行泊松方程，将计算之后的数值，代替目标地形块。Step5.如果没有满足用户的需求：则执行下一次编辑；否则算法结束。(a)实例地形(b)被替换的地形块(c)不同的地形块融合效果图图6融合地形效果图Fig6.Mergedterrain梯度形变的主要思想是通过应用局部变换到网格中，然后重新计算更新后网格的新梯度向量。注意，这里是指一个网格的坐标的梯度向量。局部变换是基于每个网格实现的。由于梯度向量平移不变性，局部变换通常采用在三角形的中心的规范局部标架。进行的变换包括旋转和缩放。可以对山体重新赋予高程值，改变其梯度域。这些高程值可以根据用户的需求对地形进行形变。如果山体的高程不满足用户的目标，我们可以对该地形块执行上述算法，从而得到平滑过度的地形块。当然这样效率会有所下降，根据地形数据块的特有性质，比较适合使用GPU并行计算，算法效率会有所提升。如图6所示，是以一块虚拟地形为例，进行演示是如何操作的。图(a)中左侧线13 框标记的是实例地形，需要被拉普拉斯处理；图(b)中部分框线标记的是将被替换的地形块，我们需要其边界条件；图(c)下半部分是将两个不同的地形块进行融合得到的最终地形结果。4.3优化为了进一步提高算法效率，我们将部分计算移植到GPU上进行执行。尤其是解线性方程组的矩阵运算，我们也利用GPU进行加速运行。我们采用OpenCL进行实验系统的优化。OpenCL是一套针对高性能处理芯片构建并行程序的强大工具集。除了较好的可移植性，OpenCL的优势还在于向量处理和并行编程。在高性能计算中，向量是一种包含相同多个数据类型的数据结构。但是，与其他数据类型不同，我们可以对向量中各个分量并行处理，也就是主机应用程序可以控制多个设备上的数据处理。OpenCL可以帮助我们实现该任务，它将不同的任务发送到不同的设备上，这种处理方式也被称为任务平行编程。此外，OpenCL还能用来加速采用OpenGL（用于3-D渲染的跨平台API）编写的程序。OpenCL应用程序的重要目标是快速处理图像中像素色彩信息相对应的二维矩阵。因此，我们将得到的虚拟地形的数据点，当成像素来处理，将其存储在GPU缓冲区中，这样能够有效地利用计算资源。主机程序设计流程图如下：内核编程设计流程如下:14 关于扩散方程在CPU与GPU上执行次数5500的时间对比，应用两组实验数据扩散方程在CPU上的统计平均时间分别为：数据维度D1024(1024×1024)时间为323.5s，D2048(2048×2048)时间约为21分钟。应用实验数据扩散方程在GPU上时间统计如图7，可以看出在显卡NVIDIAGeForceGT550M拥有2个计算单元上进行8组统计实验，对数据维度D1024执行平均时间达到为6582ms，对数据维度为D2048执行的平均时间为24033ms，是在CPU上的执行效率接近50倍，而在显卡NVIDIAGeForceGTX780Ti拥有15个计算单元进行统计实验，对数据维度D1024执行平均时间为710ms，对数据维度D2048执行的平均时间为1726ms。由此可以满足实时的地形生成的要求。图7GPU优化对比Fig7.GPUOptimization15 第五章总结与展望本文的实验开发环境是：VisualStudio2013作为主要开发平台，结合QT能够提供用户友好的交互性界面，最后利用OpenGL进行地形渲染及其可视化，并使用GLSL着色器语言结合OpenCL利用GPU编程，以保证我们的实验系统具有较好实时性渲染效果。实验通过与Arcgis提取的河网和山脊线进行了对比，如图8所示，可以看出我们的方法创建的地形的山脊线与根据水文学的方法提取的山脊线是基本相符的。因此，我们的方法与传统的过程性建模方法相比，更符合水文学客观要求，同时能够利用产生式规则生成大规模虚拟地形。(a)Arcgis提取出的山谷线(b)Arcgis提取出的山脊线(c)我们实验提取出的山谷线(d)我们实验提取出的山脊线图8我们的方法与Arcgis软件生成的特征线效果图对比Fig8.ContrastbetweenresultsofourmethodandfeaturelinesgeneratedbyArcgis图9给出了我们方法的点云效果图，根据1024×1024的点云，通过扩散方程、插值处理生成的虚拟地形，其满足山脊与河网约束双重条件。(a)以山脊线为视角焦点的地形16 (b)以山谷线为视角焦点的地形图9基于骨架特征的三维地形合成点云效果图Fig9.Diagramofterrainsynthesisbasedonfeaturesofthree-dimensionalskeleton图10给出了应用我们方法的合成的最终效果图，根据开放软件Terragen输入DEM数据生成的可视化虚拟真实地形，展现出本论文的方法产生的地形DEM数据逼近真实地形良好的表现能力。本文给出方法的优势在于能够自动建模具有河流等丰富地理特征的虚拟地形，而实际上地形的形态通常也是围绕河网的，河网构建按照其自然属性本身具有分形特征。我们根据这一思路，首先构造其河网分形产生式规则自适应地构建河网，进而根据河网形态模拟真实地形。我们是根据其Voronoi图的邻近关系来产生山脊，在此基础上形成地形的三维骨架，从而通过扩散方程来产生最终的地形。因此，我们的实验结果具有更好的真实性。本文的结果模型是基于矢量图元的，能够描述山地、峡谷和水道，可以提供多细节层次的几何形状。本文工作的主要基础是水文学中的Tokunaga网络的河网形态特征量化，使我们能够生产重要视觉特征河流，如分形维数和链长。此项工作还可以进一步扩展L-system规则和标记算法，通过实验找出一套良好的观测值，同时也可以在本文提供的方法基础上研究真实河网的复原，进而创建历史地形，这也将是我们未来的主要工作之一。(a)以山谷线为视角焦点的地形17 (b)以山脊线为视角焦点的地形图10基于骨架特征的三维地形合成效果图Fig10.Diagramofterrainsynthesisbasedonfeaturesofthree-dimensionalskeleton18 参考文献[1]FournierA,FussellD,CarpenterL.Computerrenderingofstochasticmodels[J].CommunicationsoftheACM,1982,25(6):371-384[2]PerlinK.Animagesynthesizer[J].ACMSiggraphComputerGraphics,1985,19(3):287-296.[3]EbertDS,MusgraveFK,PeacheyD,etal.TexturingandModeling:AProceduralApproach[J].AcmTransactionsonGraphics,2002,24(3):pags.777-786.[4]KelleyAD,MalinMC,NielsonGM.Terrainsimulationusingamodelofstreamerosion.[J].AcmSiggraphComputerGraphics,1988,22:263-268.[5]PrusinkiewiczP.,HammelM.Afractalmodelofmountainswithrivers.[J].ProceedingsofGraphicsInterface,1993,174–180[6]BelhadjF,AudibertP.Modelinglandscapeswithridgesandrivers:bottomupapproach.[J].Graphite’05ProceedingsoftheRdInternationalConferenceonComputerGraphics&InteractiveTechniquesinAustralasia&SouthEastAsia,2005:447-450.[7]TeohST.RiverLand:AnEfficientProceduralModelingSystemforCreatingRealistic-LookingTerrains[C].Proceedingsofthe5thInternationalSymposiumonAdvancesinVisualComputing:PartI.Springer-Verlag,2009.[8]DerzapfE,GansterB,GutheM,etal.RiverNetworksforInstantProceduralPlanets[C].ComputerGraphicsForum.2011:2031-2040[9]TasseFP,EmilienA,CaniMP,etal.Feature-basedTerrainEditingfromComplexSketches[J].Computers&Graphics,2014:101–115.[10]MusgraveFK,KolbCE,MaceRS.Thesynthesisandrenderingoferodedfractalterrains[J].AcmSiggraphComputerGraphics,1998,23(3):41-50.0[11]ChibaN,MuraokaK,FujitaK.Anerosionmodelbasedonvelocityfieldsforthevisualsimulationofmountainscenery[J].JournalofVisualization&ComputerAnimation,1998,9(4):185--194.[12]BenesB,ForsbachR,ForsbachR.VisualSimulationofHydraulicErosion[J].JournalofWscg,2002,10:79-94.[13]BenesB.,Tesinsk´YV.,HornysJ.,etal.:Hydraulicerosion[J].ComputerAnimation&VirtualWorlds,2006,17(2):99–108.[14]KristofP.,BenesB.,Kriv´AnekJ.,etal.HydraulicErosionUsingSmoothedParticleHydrodynamics[C].ComputerGraphicsForum.2009:219--228.[15]WojtanC,CarlsonM,MuchaPJ,etal.AnimatingCorrosionandErosion.[J].EurographicsConf.nat.phen,2007:15-22.[16]MeiX,DecaudinP,HuBG.FastHydraulicErosionSimulationAndVisualizationOnGPU[C].ComputerGraphicsandApplications,2007.PG"07.15thPacificConferenceon.2007:47-56.[17]VanekJ,BenesB,HeroutA,etal.Large-ScalePhysics-BasedTerrainEditing[J].IEEEComputerGraphics&Applications,2011,31(6):35-44.[18]RusnellB,MouldD,EramianM.Feature-richdistance-basedterrainsynthesis[J].VisualComputerInternationalJournalofComputerGraphics,2009,25(5-7):573-579.[19]ZhouH,SunJ,TurkG,etal.Terrainsynthesisfromdigitalelevationmodels[C].IEEETRANSACTIONSONVISUALIZATIONANDCOMPUTERGRAPHICS.2007:834--848.[20]GainJ,MaraisP,StraerW.Terrainsketching.[J].ProceedingsI3dtheAcmSiggraphSymposiumonInteractivedGraphics&Games,2009:31-38.[21]StavaO.,BenesB.,Brisbin,M.,etal.:Interactiveterrainmodelingusinghydraulicerosion[J].19 EurographicsSymposiumonComputerAnimation,2008.[22]HnaidiH,GuérinE,AkkoucheS,etal.Featurebasedterraingenerationusingdiffusionequation.[J].ComputerGraphicsForum,2010,29(7):2179–2186.[23]YuY,ZhouK,XuD,etal.MeshEditingwithPoisson-BasedGradientFieldManipulation[C].ACMTRANS.GRAPH.2004:644—651[24]JohannessenKA,KvamsdalT,DokkenT.IsogeometricanalysisusingLRB-splines[J].ComputerMethodsinAppliedMechanics&Engineering,2014,269(2):471–514[25]VibekeSkyttn,OliverBarrowclough,DokkenT.Locallyrefinedsplinesurfacesforrepresentationofterraindata[J].Computers&Graphics,2015,49:58–68[26]ZHANGL,DAIBX,WANGGQ,etal.ThequantizationofrivernetworkmorphologybasedontheTokunaganetwork[J].SciChinaSerD-EarthSci,2009,doi:10.1007/s11430-009-0176-y[27]ChorleyR.Horton,R.E.1945:Erosionaldevelopmentofstreamsandtheirdrainagebasins:hydrophysicalapproachtoquantitativemorphology.BulletinoftheGeologicalSocietyofAmerica56,275-370[J].ProgressinPhysicalGeography,1995,(4):533-554.[28]RosgenDL.Aclassificationofnaturalrivers[J].Cantena,1994,22(3):169-199[29]ShamosMI,DanH.Geometricintersectionproblems[C].FoundationsofComputerScience,1975.,16thAnnualSymposiumon.1976:208-215.[30]Jean-DavidGénevaux,EricGalin,EricGuérin,,etal.Terraingenerationusingproceduralmodelsbasedonhydrology[J].AcmTransactionsonGraphics,2013,32(4):96-96.[31]BresenhamJE.Algorithmforcmputercontrolofadigitalplotter[J].IbmSystemsJournal,1965,4(1):25-30.[32]EmilienA,PoulinP,CaniM,etal.InteractiveProceduralModellingofCoherentWaterfallScenes[C].ComputerGraphicsForum.2014.[33]TaubinG.ASignalProcessingApproachToFairSurfaceDesign[J].Siggraph’95ProceedingsofAnnualConferenceonComputerGraphics&InteractiveTechniques,1995:351-358.[34]TongY,LombeydaS,HiraniAN,etal.DiscreteMultiscaleVectorFieldDecomposition[J].AcmTransactionsonGraphics,2003,22(3):445-452.20 致谢时光如梭，如白驹过隙，我的硕士研究生的学习生活即将结束，从刚刚开始的时候对计算机图形学的小白。到经历三年时间的磨砺，坚持对计算机图形学的热爱，最终完成我的硕士论文，回首这三年来的经历路程，自己收获颇丰，也接受到了许多人的关护和真诚的帮助，有的人是在我学术迷茫时候的提纲挈领的指路，有的人是在我生活困难的时候无私无怨不求回报的帮助，感谢我的人生能有你们的出现。首先，我要感谢我的指导教师张慧杰。对待学术严谨，对本方向国内外研究风向的把握，对我们的训练培养，对待生活的智慧，是一个名副其实的指导教师，同时也是我的那个指引人生灯塔的人生导师，至此机会，我谨向张慧杰老师致以深深地谢意。其次，感谢这三年来在同一个实验室下的同学，我们在学术探讨各抒己见，有时思想共鸣，有时针锋相对，但是从那个其中感受到无比的知识，喷发出思想的火花，从中体验的真挚快乐和思想的飞速进步。在生活中，我们努力打造我们“小家”——实验室，温馨的布局，激昂的风格，还有那一排排我们认真收集来的书籍，那一个个给过我们快乐的体育用品，那一列列我们大气磅礴的实验桌，都无不凝集我们每一个人的心血。这种学习氛围，这种生活环境，能在这样一个团队中度过我的研究生生活，是我的人生一大亮彩.我还要感谢我的家人，我的父母和我的爱人，那个我最最坚强的后盾，心灵的港湾。最后，感谢对我的做出论文评阅审读盲审老师，同时也感谢在答辩时，对我提出那些非常中肯建议的老师。使我的三年来学习成果的能有提升水平的机会，让我能够明确今后的思考全面，治学严谨作风，再次感谢所有人，谢谢你们在我的人生的出现，那是我最大人生财富。21'