基于分形幾何的BIM地形模型輕量化研究

李勇 段曉峰 王宏遠

蘭州交通大學 土木工程學院，蘭州 730070

鐵路線路地形具有長帶狀特征，現有的BIM 和GIS（Geographic Information Science）軟件使用過程中地形模型數據量過大且冗余過多，導致加載緩慢或長時間無法響應。若大范圍刪減，又可能因刪減不當導致數據點無法較好地表征實際地形特征；同時，構網方式單一化，無法滿足工程多式樣的地形顯示。

在地形建模和輕量化研究方面，王倩［1］研究了三角形網格地形建模的基本方法和理論，包括規則數據及不規則數據三角網的生成算法。文雄志［2］運用層次細節技術和基于網格的模型簡化技術對三維地形建模進行簡化，使其能夠更好地實現建模簡化技術的應用。陸風等［3］對基于地形線的Memoryless 改進算法的簡化效果進行了驗證，顯示水下地形模型簡化后保持了原有特征和細節，但模型會產生跳躍現象。宋偉華等［4］通過Delaunay 逐點插入算法構建TIN 模型，該算法能用更少的數據表達復雜地表形態，但時間復雜度差、運行速度慢。分形理論的出現和發展，為描述自然界的復雜事物提供了解決思路。李鑫［5］利用數字圖像處理技術實現了獲取分形維數的流程，具有較好的測量穩定性，無需占用太多計算機資源即可獲得較高精度。張欣悅等［6］提出利用有理分形插值進行分形曲線建模，驗證了可變參數的有理分形插值曲線在處理實際問題中具有較好的應用價值。彭曉蕾等［7］對蘊含地貌特征的黃土高原地貌數據進行分形設計，將自然景觀特征與公式化的數學思維相結合，為自然景觀的美學研究提供了一個新的方向。

鑒于分形算法由極少部分地形控制點就可表達出原地形的特征，還能在保持地形真實感的狀態下根據實際需求動態控制地形的細化程度，在地形構網上具有優越性，因此本文嘗試基于分形算法解決BIM 地形模型輕量化表達問題。

1 分形理論概述

1.1 分形的概念及維數

若集合F在歐氏空間中的Hausdorff 維數DH恒大于其拓撲維數Dr，即DH>Dr，則該集合是分形集，簡稱分形［8-9］。分形維數是一種非整數維數，用來解決整數維在尺度上出現過大或過小而無法準確評價的問題。如果將分形拓展到表示自然界中的行為或現象，那么分形維數描述的則是零碎的局部特征構成整體系統行為的相關性。

1.2 地形分形維數的計算

地形分形維數計算常用的方法為表面積-體積法和方差圖法。

1.2.1 表面積-體積法

用Hausdorff 面積SH表示地形面積，V表示地形曲面的體積。由量綱分析可得

式中：D為地形曲面的分形維數。

若用S表示地形曲面的歐式面積，則有

式中：k0為系數，又稱形狀因子；l為尺度值。

則可推導出直接計算分形維數的公式為

1.2.2 方差圖法

方差圖法是建立在圖像（分形數據）的高斯統計模型上的。給定一個分形維數D和方差δ，用分形布朗運動（Fractional Brownian Motion，FBM）模型可模擬生成一幅圖像。反之，如果將圖像或地形表面看成是分形布朗曲面，就可以得到其FBM模型的分形維數。

2 基于分形理論的數字地形模型生成算法

用于地形模擬的分形算法主要為基于迭代函數系統（Iterated Function System，IFS）方法和基于分形布朗運動（FBM）方法［9-10］。IFS 方法是對已有數據集進行仿射變換，構造一類插值函數，進而通過迭代過程，得到需要的形狀。FBM 方法通過控制表征值H、δ等參數，加入隨機變量然后構造出所需的分形圖形［11］。本文采用基于FBM的方法對地形進行重構。

FBM 是一種隨機過程，其增量按正態分布［12-13］。將FBM 二維曲面定義為在某概率空間上的隨機過程，且滿足：X（0，0）=0，X（x，y）是（x，y）的連續函數。

對任意變量（x，y），(h，k) ∈R2，其二維增量X(x+h，y+k) -X(x，y)服從期望為0，方差為(h2+k2)2H的正態分布，其概率滿足

基于FBM 的插值方法主要有泊松階躍法、隨機中點位移法（Random Midpoint Displacement，RMD）、傅立葉濾波法、逐次隨機增加法等［13］。本文采用基于FBM的改進隨機中點位移法，通過對OpenRail 軟件二次開發以改進優化原有地形生成時的臃腫情況。

中點位移法是標準的分形幾何法，其在細分過程中，通過對兩個或多個頂點之間插值進行地形建模，具有高速度以及為已有形狀增加細節的能力，但在實際生成地形過程中地形曲面會出現褶皺現象。為此采用一種加入與分形特征表征值H、σ相關的補償項的改進隨機中點位移法。該方法既能滿足用分形布朗運動特征的分形表面來表達自然地形，同時又能很好解決此問題［14-16］。

設01，02，…，0n為初始柵格數據點，格網間距為d0。每一級細分包括兩步內插過程。

1）第一步：diamond步

由最鄰近的四個格網點高程h0i（i=1，2，3，4），內插其中心點h1i的高程。高程值等于四個格網點高程的平均值加上一個隨機位移Δ ～N(0，V2?)。h1i計算式為

此時，原始的正方形格網旋轉45°成為菱形格網，格網間距也由d0變為d1，即

2）第二步：Square步

由第一步得到的已知菱形四個角點高程內插其中心點高程，即

經過第二步，格網間距又變為d2，即

如此，隨機中點位移法的一個完整步驟就完成了，如圖1所示，接著按照所需迭代次數反復進行迭代計算就可得到所需的分形曲面。但存在隨機變量的取值問題。通常，對Δ 的選擇要使待插點的高程增量滿足FBM 所需的冪指數規律。令格網間距的變化尺度為γn（n為迭代水平），那么隨機位移的方差為

圖1 隨機中點位移法計算示意

由d1和d2可知γ=2-n/2，因此演變出常用的隨機位移公式，即

式中：sc為位移大小調整因子；G為高斯隨機分量，G～N(0，1)。

式（10）一般用于虛擬分形地形的生成，不適用于自然世界中實測數據的插值。為使待插點的分形特征符合實際數據的分形特征，Yokaya 提出了式（11）的隨機位移增量。

本文采用文獻［9］提出的改進Yokaya 隨機增量，即

3 基于分形算法的地形模型重構技術

為了直接將分形地形與實際線路BIM 設計相結合，選用OpenRail Designer進行二次開發。

3.1 OpenRail二次開發環境搭建

本文所使用的軟件及配置為：①操作系統為Windows10，內存16G，顯卡GTX1650，顯存8G；②編譯平臺為Visul Studio2015；③開發環境為.NET4.6.2 及以上版本的框架，與OpenRail版本匹配的SDK。

3.2 基于分形算法重構地形模型

3.2.1 原始數據重采樣

本次地形數據采用無人機航測數據，在原有數據的基礎上進行重采樣，以滿足算法需要。

由分形布朗運動相關知識可知，要計算分形維數D，需要不斷改變其尺度。因此本文在提取分形特征值時，考慮實際地形數據特征，將采樣間隔設置為20、40、50、60 m，且為判斷不同軸采樣對分形特征值的影響，分別沿數據的x軸、y軸和對角線采樣，獲得不同的采樣結果，見圖2、圖3。

圖2 高程點重采樣

圖3 重采樣柵格圖

3.2.2 分形特征提取

假設f（x）為一個隨機分形布朗運動，則滿足式（13）的概率分布函數［14-15］。

式中：f（x）為采樣值；Δx為采樣間隔；H為頻譜指數；F（y）為累積概率分布函數，服從標準正態分布（0，δ2）。而分形維數D=N+1-H，N為拓撲維數。對地形曲面進行分形模擬時，H表示地形復雜情況的一種抽象和概括，H越大分形布朗運動變化趨于平緩，H越小變化趨于精細；δ反映地形曲面的起伏特征。

式（13）又可改寫為

最小二乘法的直線回歸函數為

由此求得系數H、C，便可推得D和δ。

為比較采樣方向對地形特征的影響，將整體數據劃分為四個大小一致的區塊依次計算各分塊分形特征值，用來比較區塊分形特征值與整體分形特征值的差別，見圖4、圖5。

圖4 整體數據分形特征擬合

圖5 各區塊分形特征擬合

由圖4（a）可知，沿x軸采樣數據整體只呈現一個無標度區，采樣點數據分布趨勢基本一致，擬合可得該地形沿x軸的分形特征值為H=0.432 8，δ=3.845。

由圖4（b）可知，沿y軸采樣數據整體只呈現一個無標度區，無錯誤點，擬合可得該地形沿y軸的分形特征值為H=0.476 5，δ=1.906。

由圖4（c）可知，沿對角線采樣數據整體只呈現一個無標度區，無錯誤點，擬合可得該地形沿對角線的分形特征值為H=0.659 1，δ=0.845。

由圖5（a）可知，對區塊1 進行采樣獲得的數據以lg Δr=8 為界限呈現兩個無標度區，在lg Δr=7.8 后出現另一個無標度區，與區塊1 多數數據呈現的特征不符，因此需要剔除后面的采樣點，最終擬合可得該區塊地形分形特征值為H=0.500 7，δ=1.998。

由圖5（b）可知，對區塊2 進行采樣獲得的數據整體呈現一個無標度區，最終擬合可得該區塊地形分形特征值為H=0.430 3，δ=2.228。

由圖5（c）可知，對區塊3 進行采樣獲得的數據以lg Δr=7.8 為界限呈現兩個無標度區，在lg Δr=7.8后出現另一個無標度區，與區塊3 多數數據呈現的特征不符，因此擬合時只對第一個無標度區進行擬合，最終擬合可得該區塊地形分形特征值為H=0.529 1，δ=1.935。

由圖5（d）可知，對區塊4 進行采樣獲得的數據呈現兩個無標度區，在lg Δr=8.2 后出現另一個無標度區，與區塊4多數數據呈現的特征不符，因此擬合時只對第一個無標度區進行擬合，最終擬合可得該區塊地形分形特征值為H=0.543 6，δ=2.202。

對所有的數據結果進行整理，結果見表1?？芍?，無人機獲取的試驗區域地形數據，無論采用三個方向對整體數據采樣，或者是對數據劃分，提取出的分形特征值基本一致，分形維數D始終保持在2.5 左右。因此在進行分形插值構網時，可以采用全局使用一個分形特征值進行插值算法的編譯。

表1 地形分形特征值統計

3.2.3 曲面生成

根據分形布朗運動原理，基于地形數據、重采樣數據及擬合出的分形特征值，利用改進隨機中點位移法，通過編制分形程序并加載到BIM 軟件OpenRail Designer，可實現不同特征值的地形分形模擬。

對同一地形（圖6），采用區塊分形特征值（H=0.501，δ=2.091）和整體分形特征值（H=0.523，δ=2.199）進行插值構網，迭代一次得到圖7（a）和圖7（b），視覺上兩個地形的表達差異極小。

圖6 原始地形

當使用整體分形特征值迭代兩次以后，得到的地形［圖7（c）］顯示效果優秀，已經完全能夠表達原地形的起伏特征，從表2的平均坡率也能驗證這一點。

圖7 三種工況的分形地形

表2 原始地形與分形地形指標對比

由表2 可知：第一次迭代點數為8 320，僅為原始地形點數的0.831%，網格數也僅為原始網格數的0.821%；經第二次迭代，點數與網格數分別是原始數據的2.270%和2.286%。結合對比運算時間，采用分形算法對地形進行輕量化是切實可靠的，可用極少地形控制點完美實現地形特征表達，在節約內存、提高地形構網效率方面也有較大優勢。

此方法能更好地解決鐵路BIM 與GIS設計過程中長帶狀地形數據量龐大、加載緩慢的問題。

4 結論

1）通過比較區塊與整體分形特征值，無論是采用三個方向對整體數據采樣，或是對數據劃分，對每個區塊采樣進行分形特征值的擬合，提取出的分形特征值基本一致，分形維數D始終保持在2.5左右。因此，在進行分形插值構網時，可以采用全局使用一個分形特征值進行插值算法的編譯。

2）基于分形算法的地形模型重構可以由極少部分地形控制點較好地表達原地形的特征，并能根據實際需求動態控制地形的細化程度，解決軟件構網方式單一化問題，滿足工程實際中多式樣的地形顯示需求，解決了地形輕量化表達問題。