基于地面三維激光點云實現地質體產狀自動提取的方法
【專利摘要】本發明公開了一種基于地面三維激光點云實現地質體產狀自動提取的方法，利用三維點云數據空間位置和顏色信息，結合圖像分割和拐點檢測技術實現地質體出露面和未出露面產狀自動提取。包括如下步驟：點云數據采集和預處理；獲取點云數據的txt文本格式；采用最近鄰法將點云數據重采樣到二維圖像；分割二維圖像；評價圖像分割精度；分別擬合地質體出露面和未出露面，并計算其產狀；其中地質體未出露面在二維圖像上表現為線狀要素，且線條走向的大角度變化代表未出露面產狀的變化，因此，計算未出露面產狀較復雜，涉及拐點檢測、赤平投影等技術。本發明在地質體產狀提取方面，具有十分重要的實際應用價值，大量減少了危險且耗時耗力的實地測量工作。
【專利說明】基于地面三維激光點云實現地質體產狀自動提取的方法

【技術領域】
[0001]本發明涉及三維點云與影像分割技術和構造地質科學領域，特別是一種高精度、高密度、包含空間和顏色信息的三維激光點云實現地質體產狀自動提取的方法。

【背景技術】
[0002]巖體結構面的優勢產狀和幾何信息是評價巖體結構面發育的重要指標，特別是結構面中規模小、數據量大且隨機分布的節理和裂隙，其優勢產狀是評價巖體穩定性的重要參數。野外實地測量不但耗時耗力，而且在一些危險地區，人員無法到達。三維激光掃描是一種遠距離非接觸的全自動立體掃描技術，又被稱為“實景復制技術”。獲取的點云數據精度高且數據量大，涵蓋目標的三維空間位置和真彩色信息，目前已在數字地質露頭、地質工程等領域得到應用。Gigli and Casagli (2011)和劉昌軍等(2011)都提出了基于三維點云數據提取地質體產狀的方法，他們的基本思想是利用結構面上多個點的空間坐標來擬合平面，并根據平面方程計算地質體產狀。
[0003]然而目前研究僅限于出露地質體產狀的人機交互測量，無法實現地質體未出露面產狀的提取。野外，地質體經歷了剝蝕和沉積等過程，其出露面不能真實代表其。例如巖層面產狀的提取，表層的巖層(出露)因為外力作用已經變得破碎不平整，因此要準確計算巖層的產狀，我們需要測量兩地層之間的接觸面，而該接觸面是未出露的；又如作為火山作用與深部地幔性質、演化及殼幔相互作用的紐帶的巖墻，工作者需要測量它與圍巖接觸面來計算巖墻產狀，其產狀有助于計算巖漿超壓，研究巖漿動力學及侵位機制。因此在實際測量工作中，獲取地質體未出露面產狀具有更重要的意義。此外，目前對三維點云數據的研究大多未涉及其顏色信息RGB。鑒于以上情況，充分利用點云數據的空間位置和顏色信息，在出露地質體產狀提取的基礎上，結合影像分割技術，ENVI EX(http://www.esrichina.com.cn/2010/1029/123.html)是目前普遍運用的分割軟件之一，邊界拐點檢測技術，類似的工作包括利用8k (k = 2, 3,4…)鄰域鏈碼和道格拉斯-普克(Douglas-Peuker)抽稀算法(http://www.baike.com/wiki/抽稀)檢測邊界拐點,進一步實現地質體出露面和未出露面產狀的自動提取，具有重要的實用價值。
[0004]Riegl VZ100 (www.riegl.com)是一種目前應用廣泛的三維激光掃描儀，其有效掃描距離最大可達1400m，距工作站點10m內點間距為5mm，水平掃描范圍為0°?360°，垂直掃描范圍為-40°?60°。通過不同站位的多角度掃描，可以獲取覆蓋實驗區的全部點云數據。其配套軟件RISCAN PRO可以對三維點云數據進行賦色、配準、去噪等預處理，從而為后續工作提供較好的實驗數據。


【發明內容】

[0005]本發明的目的在于針對現有技術的不足，利用點云數據的空間位置和顏色信息，結合三維點云數據處理、影像分割和邊界拐點檢測技術提供一種地質體出露面和未出露面產狀自動提取的方法。
[0006]本發明的目的是通過以下技術方案來實現的:一種基于地面三維激光點云實現地質體產狀自動提取的方法，包括如下步驟:
[0007](I)通過三維激光掃描儀采集三維激光點云數據；
[0008](2)對步驟I采集的三維激光點云數據進行預處理，包括賦色、配準、去噪操作；
[0009](3)獲取點云數據的txt文本格式，該文本的每一行包含一個點云數據的空間坐標XYZ和顏色信息RGB;根據其空間坐標，采用最近鄰法將點云數據重采樣到二維圖像，SP將三維點投影到二維圖像上，該圖像包含6個圖層數據，即點云數據的X、Y、Z、R、G、B ;
[0010](4)基于點云數據的顏色信息，對步驟3中獲取的二維圖像進行分割；采用監督評價方法評價圖像分割精度，不斷調整分割參數直到分割精度達到閾值；
[0011](5)地質體出露面和未出露面在二維圖像上分別表現為面狀要素和線狀要素，即分割對象和分割邊界線，當分割精度達到設定閾值時，對于出露地質體，利用組成該對象的三維點，采用最小二乘法擬合結構面并計算其產狀；對于未出露地質體，只能掃描到其邊界線上的點云數據，且邊界線走向的大角度變化代表其產狀的變化，因此其產狀的計算較出露面復雜，具體包括以下步驟:
[0012](5.1)利用8k鄰域鏈碼和道格拉斯-普克抽稀算法檢測邊界拐點，其中，k為大于等于2的整數；
[0013](5.2)以拐點為分隔點，根據邊界線上點云數據的空間坐標XYZ，利用最小二乘法分段擬合地質體未出露結構面并計算其產狀；
[0014](5.3)將步驟5.2中計算得到的地質體產狀進行赤平投影，對產狀相近的地質體進行歸類，并賦予同一標志；合并具有相同標志的相鄰地質體，即相鄰分界線，對于孤立的地質體，如果其相鄰地質體屬于同一類，則將該孤立地質體與其相鄰地質體歸為一類，然后分別根據每類邊界點擬合平面并計算產狀；重復步驟5.3，直到相鄰地質體產狀均不相同。
[0015]進一步地，所述步驟4中所述監督評價方法通過比較分割結果與參考對象之間的相似度或差異度來評價分割精度，具體為:采用基于面積和基于位置的組合度量來評價圖像分割精度。
[0016]本發明的有益效果是:通過結合利用點云數據的空間位置信息和顏色信息，從而實現“三維一二維一三維”的技術路線，將研究較成熟的影像分割技術用于三維點云數據分害I]，并通過精度評價來提高其分割效果，減少對后續操作的影響。本發明基于點云數據的顏色信息獲取了地質體在二維圖像上的面狀和線狀表現，并根據三維點的空間坐標擬合地質體結構面，從而獲取其產狀信息。利用該發明，可以實現出露和未出露地質體產狀的自動提取，從而減少了大量野外實地測量工作，不但節約了時間和資金，還減少了工作的危險性。

【專利附圖】

【附圖說明】
[0017]圖1是本發明基于地面三維激光點云實現地質體產狀信息自動提取方法的流程圖；
[0018]圖2是本發明不同站點對同一目標的掃描不意圖；
[0019]圖3是本發明基于面積度量評價分割精度獲取相交面積的流程圖；
[0020]圖4中(a)是具體實施例中對點云數據進行重采樣得到的二維圖像，其中白色線條是手動分割巖墻的邊界，作為圖像分割精度評價中的參考對象；(b)是ENVI EX軟件分割巖墻得到的邊界二值圖；(C)是8k鄰域鏈碼檢測的拐點圖(黑色點為拐點)，其中k為2;(d)是道格拉斯-普克算法檢測的拐點圖(黑色點為拐點)，其中抽稀因子為3.5。

【具體實施方式】
[0021]如圖1所示，本發明基于地面三維激光點云實現地質體產狀自動提取方法，包括以下步驟:
[0022]1.通過三維激光掃描儀采集三維激光點云數據；
[0023]三維激光掃描儀Riegl VZ100在掃描過程中由于受激光掃描儀視場角的限制，單一站點的掃描范圍往往不能覆蓋整個研究區域，因此需要進行多站點掃描，且不同站點之間必須保證存在重疊區，以便實現不同站點間點云數據的配準；具體可參考文獻((Simon Jj Buckley J A，et al，2008.Terrestrial laser scanning in geology: dataacquisit1n, processing and accuracy considerat1ns.Journal of GeologicalSociety,London, 165:625-638》；
[0024]2.對步驟I釆集的三維激光點云數據進行預處理，包括賦色、配準、去噪操作；
[0025]在三維激光掃描儀Riegl VZ100的配套軟件RISCAN PRO中完成三維點云數據的預處理；Riegl VZ1000在掃描的同時會拍照記錄點云數據的顏色信息，在RISCAN PRO中通過Color from images操作，自動完成相應照片對點云數據的賦色，為其增加RGB信息；配準是將不同站點坐標系下的點云數據轉換到一個共同的基準坐標系下，主要包括步驟:1)確定同名點對；2)解算坐標轉換數據函數，通過RISCAN PRO軟件中的Coarseregistrat1n操作完成不同站點點云數據的配準，該步驟需要用戶在兩個站點中選取至少
4組同名點對進行配準，RISCAN PRO軟件會計算配準誤差，根據該誤差確定是否需要重新選點；去噪即是去掉那些受測量儀器或其他環境因素影響產生的噪聲點，該步驟可以通過定義RISCAN PRO軟件中的Filter操作完成，軟件用戶說明書中具體介紹了定義Filter操作的各參數；
[0026]3.在RISCAN PRO軟件中導出點云數據的txt文本格式，該文本的每一行包含一個點云數據的空間坐標XYZ和顏色信息RGB ;根據其空間坐標，采用最近鄰法將點云數據重采樣到二維圖像，即將三維點投影到二維圖像上，該圖像包含6個圖層數據，即點云數據的X、Y、Z、R、G、B,具體的:將三維點云數據投影到xy平面上，并對該區域進行圖像重采樣，最近鄰法可表達為:
[。。27] Disl =+
[0028]其中j和i是正整數，Xj和Y」是點z?數據的X和Y坐標，Xi和Yi是圖像點坐標，Q1和Q2是常數，分別表示點云數據坐標中最小的X值和Y值，Dist是點云數據與重采樣二維圖像中某一確定像素點的距離，對于每一個圖像點，即一個確定的i，我們需要遍歷所有點云數據，即j = I, 2, "^nOi是點云數據的個數)，計算每個點云數據與該像素點的距離，然后選取距離該點位置最近的點云數據(即Distmin)的空間坐標XYZ和顏色信息RGB作為其波段屬性值；
[0029]4.基于點云數據的顏色信息,利用ENVI EX軟件中的Feature Extract1n操作對步驟3中獲取的二維圖像進行分割；采用監督評價方法評價圖像分割精度，該方法是通過比較分割結果與參考對象之間的相似度或差異度來評價分割精度的，在進行該操作前需要準備好參考對象和分割結果，其中參考對象可以通過目視解譯手動分割獲得；通過調整分割參數(ENVI EX軟件中分割參數有Scale level和merge level)以達到分割精度閾值，該閾值由用戶自己設定；利用基于面積和基于重心的組合度量來評價分割精度，組合度量避免了單一度量對分割不足或過度分割的偏好，主要度量公式是:
[0030]基于面積的度量:


area i^X: η Y-)
[0031]OverSegnientat1njj = 1-----


area (Xj)


area (X1 n Yj)
[0032]UnderSeginentat1njj -1----,~~-


areal Yj)
[0033]基于位置的度量:
[0034]Dij = dist (centroid (Xi), centroid (Yj))
[0035]其中X 是參考對象，Y 是分割對象，OverSegmentat1n 和 UnderSegmentat1n 分別表示過度分割和分割不足，Du表示參考對象中心和分割對象中心的距離^rea(Xi)和area(Yj)分別表不參考對象和分割對象的面積，area(Xi Π Yj)是參考對象與分割對象的相交面積，Centroid(Xi)和Centroid(Yj)分別表示參考對象和分割對象的中心；
[0036]OverSegmentat1n> UnderSegmentat1n 和 D 的最優值均為 O,當三者同時取最優值O時即為最優分割，即分割對象與參考對象完全重合；具體可參考文獻《Clinton N, Holt A, Scarborough J, et al, 2010.Accuracy assessment measures forobject-based image segmentat1n goodness.Photogrammetric Engineering andremote sensing, 76(3), 289-299》；
[0037]5.地質體出露面和未出露面在二維圖像上分別表現為面狀要素和線狀要素，即分割對象和分割邊界線，當分割精度達到設定閾值時，對于出露地質體，利用組成該對象的三維點，采用最小二乘法擬合結構面并計算其產狀；
[0038]根據空間幾何原理，結構面的方程可以寫為:
[0039]ax+by+cz+d = O
[0040]其中a，b，c，d為結構面空間常數，x，y，z是三維點空間坐標變量。結構面法線與X軸、Y軸、Z軸的夾角余弦(l，m，n)按以下公式計算:
/ _ a
[0041]/= /.,ο ο
如-+b~ +c~

b
[0042]m =卜,

^Jcr + Ir +c~

c
[0043]11 = I 2 ,2 2.V? +b +c
[0044]根據結構面的法向向量計算其傾向α和傾角β (傾向是傾斜平面上走向線的垂線在水平面上的投影所指的沿平面向下傾斜的方位，其中走向線是傾斜平面與水平面的交線，傾角是傾斜平面與水平面所成的角度)；

? ?ηλ
[0045]cc - arctan — Λ O
? \
^η
[0046]β = arctan.， ，

w/2+^2 J
[0047]其中Q 是常數，當 1>0，m>0 時，Q = 0° ;當 1>0，m〈0 時，Q = 360° ;當1〈0 時，Q = 180。；具體可參考文獻《Gigli G，Casagli N, 2011.Sem1-automatic extract1n ofrock mass structural data from high resolut1n LIDAR point clouds.1nternat1nalJournal of Rock Mechanics and Mining Sciences，48(2)， 187-198》
[0048]對于未出露地質體，只能掃描到其邊界線上的點云數據，因此利用其邊界線上的三維點(只要有3個不在同一直線上的空間點便能擬合平面，三維激光儀的多角度掃描和野外露頭的不平整可以保證邊界線上的點不在同一直線上)來擬合結構面，且邊界線走向的大角度變化代表其產狀的變化，因此其產狀的計算較出露地質體復雜，具體包括以下步驟:
[0049](5.1)利用8k(k = 2,3,4…)鄰域鏈碼和道格拉斯-普克(Douglas-Peuker)抽稀算法檢測邊界拐點；8k鄰域鏈碼和道格拉斯-普克抽稀算法分別適用于柵格圖像和矢量圖像，針對不同的目標物，需要通過多次實驗來確定最優k值或抽稀因子。
[0050]8k(k = 2，3，4…)鄰域鏈碼適用于具有閉合邊界的二值柵格圖，該算法首先利用8k鄰域鏈碼求得二值圖像中物體邊界個點的曲率，然后通過尋找曲率的局部峰值檢測出拐點；具體可參考文獻《尚振宏，劉明亞.二值圖像中拐點的實時檢測算法.中國圖像圖形學報，2005，10 (3)))；
[0051]道格拉斯-普克(Douglas-Peuker)抽稀算法適用于矢量圖，其基本思路為:
[0052](a)對曲線的首末點虛連一條直線，求曲線上所有點與該直線的距離，并找出最大距離值dmax,用dmax與事先給定的閾值D相比:
[0053](b)若dmax〈D，則將這條曲線上的中間點全部舍去；若dmax彡D，保留dmax對應的坐標點，并以該點為界，把曲線分為兩部分，對這兩部分重復步驟I和2，直到所有dmax均〈D，即完成對曲線的抽稀；
[0054](5.2)以拐點為分隔點，根據邊界線上點云數據的空間坐標XYZ，利用最小二乘法分段擬合地質體未出露結構面，然后按照與計算出露面產狀相同的方法來計算地質體未出露面的產狀；
[0055](5.3)將步驟5.2中計算得到的地質體未出露面產狀在軟件DIPS中進行赤平投影，對產狀相近的地質體進行歸類，并賦予同一標志，合并具有相同標志的相鄰地質體；對于孤立的地質體，如果其相鄰地質體屬于同一類，則將該孤立地質體與其相鄰地質體歸為一類，然后根據每一類地質體邊界上的三維點重新擬合結構面并計算產狀；并重復步驟(5.3)，直到相鄰地質體產狀均不同。該步驟通過合并產狀相近的相鄰地質體，可以增加擬合未出露地質體的三維點，從而避免因范圍過小而擬合出錯誤的地質體未出露面。
[0056]下面以巖墻產狀提取為例進行說明。
[0057]巖墻作為一種厚度比較穩定的板狀侵入體，巖墻面幾乎不出露，出露地表的僅是巖墻的厚度，利用出露面上多個點的空間坐標構建平面方程來獲取結構面產狀的方法不適用于巖墻。因此，利用本發明提出的基于地面三維激光點云實現地質體產狀自動提取的方法來提取巖墻產狀信息。
[0058]圖1是基于地面三維激光點云實現地質體產狀自動提取方法的流程圖。使用者按步驟I采集原始點云數據；按步驟2對點云數據進行預處理，包括賦色、配準、過濾操作；按步驟3導出點云數據的txt文本格式，其中該文本的每一行包含一個點云數據的空間坐標XYZ和顏色信息RGB，并根據最近鄰法將其重采樣為二維圖像，實現三維到二維的轉換；按步驟4分割研究對象巖墻，利用監督評價方法評價分割精度，圖3是基于面積度量評價分割精度獲取相交面積的流程圖，通過不斷調整分割參數使分割精度達到設定閾值(80%);鑒于巖墻屬于未出露地質體，在二維圖像上表示為線狀要素，因此按步驟5.1檢測巖墻邊界拐點；按步驟5.2分段計算巖墻產狀，實現二維到三維的轉換；按步驟5.3不斷對分段巖墻進行歸類，并計算產狀。
[0059]參閱圖4所示，圖4a是根據最近鄰法將點云數據重采樣得到的二維圖像，圖4a中的白色線條是手動分割巖墻的邊界，作為圖像分割精度評價中的參考對象。圖4b是利用ENVI EX軟件對巖墻進行分割后獲取的邊界二值圖，以圖4a中的分割巖墻為參考對象，采用基于面積和位置的組合度量評價圖4b巖墻的分割精度，總體精度達80%以上，但比較圖4a和圖4b可知:ENVI EX軟件對巖墻的分割忽視了左側巖墻分支。利用8k(k = 2,3,4...)鄰域鏈碼和道格拉斯-普克(Douglas-Peuker)抽稀算法分別檢測巖墻邊界拐點(黑色點)得圖4c (k取2)和圖4d(抽稀因子取3.5)，其中圖4c和圖4d分別包含35和49個拐點。對于該實驗數據，D-P算法抽稀檢測出的拐點更多，對巖墻的分段更細。最后通過分割線上的空間點分段擬合巖墻未出露面，并不斷合并產狀相近的相鄰面，最終得到兩個大平面，其產狀(傾向、傾角)分別是Ν55.4° W，76.4° (圖4b中近豎直方向的巖墻①)和Ν25.8° E，60.6° (圖4b中右側巖墻分支②)。
【權利要求】
1.一種基于地面三維激光點云實現地質體產狀自動提取的方法，其特征在于，該方法包括如下步驟: (1)通過三維激光掃描儀采集三維激光點云數據； (2)對步驟I采集的三維激光點云數據進行預處理，包括賦色、配準、去噪操作； (3)獲取點云數據的txt文本格式，該文本的每一行包含一個點云數據的空間坐標XYZ和顏色信息RGB ;根據其空間坐標，采用最近鄰法將點云數據重采樣到二維圖像，即將三維點投影到二維圖像上，該圖像包含6個圖層數據，即點云數據的X、Y、Z、R、G、B ; (4)基于點云數據的顏色信息，對步驟3中獲取的二維圖像進行分割；采用監督評價方法評價圖像分割精度，不斷調整分割參數直到分割精度達到閾值； (5)地質體出露面和未出露面在二維圖像上分別表現為面狀要素和線狀要素，即分割對象和分割邊界線，當分割精度達到設定閾值時，對于出露地質體，利用組成該對象的三維點，采用最小二乘法擬合結構面并計算其產狀；對于未出露地質體，只能掃描到其邊界線上的點云數據，且邊界線走向的大角度變化代表其產狀的變化，因此其產狀的計算較出露面復雜，具體包括以下步驟: (5.1)利用8k鄰域鏈碼和道格拉斯-普克抽稀算法檢測邊界拐點，其中，k為大于等于2的整數； (5.2)以拐點為分隔點，根據邊界線上點云數據的空間坐標XYZ，利用最小二乘法分段擬合地質體未出露結構面并計算其產狀； (5.3)將步驟5.2中計算得到的地質體產狀進行赤平投影，對產狀相近的地質體進行歸類，并賦予同一標志；合并具有相同標志的相鄰地質體，即相鄰分界線，對于孤立的地質體，如果其相鄰地質體屬于同一類，則將該孤立地質體與其相鄰地質體歸為一類，然后分別根據每類邊界點擬合平面并計算產狀；重復步驟5.3，直到相鄰地質體產狀均不相同。
2.根據權利要求1所述的一種基于地面三維激光點云實現地質體產狀自動提取的方法，其特征在于，所述步驟4中所述監督評價方法通過比較分割結果與參考對象之間的相似度或差異度來評價分割精度，具體為:采用基于面積和基于位置的組合度量來評價圖像分割精度。
【文檔編號】G06T17/05GK104183017SQ201410367744
【公開日】2014年12月3日 申請日期:2014年7月29日 優先權日:2014年7月29日
【發明者】陳寧華, 倪妮娜, 陳建裕 申請人:浙江大學