點云地面點濾波(Cloth Simulation Filter, CSF)

本篇博客參考Wuming Zhang的文章“An Easy-to-Use Airborne LiDAR Data Filtering Method Based on Cloth Simulation”

不方便的小伙伴可以在此：資源鏈接下載。

1. 引言

機載LiDAR可以獲取快速、低成本地獲取大區(qū)域的高精度地形測量值。為了獲取高精度的地形數(shù)據(jù)（厘米級），對機載LiDAR點云數(shù)據(jù)進(jìn)行“濾波”是一個非常重要的步驟。因此近二十年來，國內(nèi)外學(xué)者提出了眾多有效的自動濾波算法，大大降低了人力成本，提高了點云數(shù)據(jù)的應(yīng)用效率。在博客：“點云地面點濾波(Progressive Morphological Filter)算法介紹（PCL庫）”與“點云地面點濾波-progressive TIN densification(PTD)算法介紹”中分別對兩種常用濾波器進(jìn)行了介紹。
以往的眾多算法濾波效果容易受到地形特征的影響（通常在復(fù)雜場景及陡峭地形區(qū)域濾波效果較差）且常常需要用戶對數(shù)據(jù)有較為豐富的先驗知識來進(jìn)行設(shè)置濾波器中的各種參數(shù)。
為了解決這些問題，張吳明教授等人在2015年第三屆全國激光雷達(dá)大會提出的“布料”濾波算法，本篇博客對其進(jìn)行了簡要介紹。

2. 基本思想

傳統(tǒng)的濾波算法大多是考慮在坡度、高程變化之間的不同來進(jìn)行區(qū)分地物點與地面點，而布料”濾波算法從一個完全新的思路來進(jìn)行濾波，首先把點云進(jìn)行翻轉(zhuǎn)，然后假設(shè)有一塊布料受到重力從上方落下，則最終落下的布料就可以代表當(dāng)前地形。如下圖所示，可以幫助我們進(jìn)行理解。

下方的藍(lán)色圖形代表原始測量值，對其進(jìn)行翻轉(zhuǎn)，紅色虛線代表布料，可以反映地形的起伏。

3. CSF算法實現(xiàn)步驟

3.1 “布料”模擬

首先，定義了一個基本公式：

式中，X代表“布料”中的粒子在 t 時刻的位置，F(xiàn)ext(X,t)代表外部驅(qū)動因素（重力，碰撞等），F(xiàn)int(X,t)代表內(nèi)部驅(qū)動因素（粒子間的內(nèi)部聯(lián)系）。

即可以總結(jié)為“布料”粒子的位置受到Fext(X,t)與Fint(X,t)兩方面因素的影響。

隨后如下圖所示，通過四個狀態(tài)對“布料”粒子的移動過程進(jìn)行說明。

（a）為初始狀態(tài)，假設(shè)一個虛擬的布放在翻轉(zhuǎn)的LiDAR測量值上方。黑色點為LiDAR測量值，深紅色點為可移動粒子。
（b）計算重力影響下對粒子產(chǎn)生的位移，由此一些粒子可能會出現(xiàn)在地面的下方。淺紅色點為粒子的舊位置。
（c）相交檢查，檢查粒子新位置是否到達(dá)地面，如果已經(jīng)到達(dá)地面則設(shè)置為不可移動。藍(lán)色點代表粒子不可移動
（d）考慮“布料”的內(nèi)部驅(qū)動，對現(xiàn)有的“可移動”粒子，根據(jù)鄰近粒子所產(chǎn)生的“力”來進(jìn)行移動新的位置。（如上圖小房子處的可移動粒子）

3.2 外部/內(nèi)部因素驅(qū)動

如3.1節(jié)中所述，粒子的位置移動包括外部/內(nèi)部兩方面的驅(qū)動因素，首先假設(shè)只有外部因素Fext(X,t)，并設(shè)置內(nèi)部因素Fint(X,t)為0，可以得到以下公式：

式中，m為粒子的重量，通常設(shè)置為1，△t 為時間步長。且由于G是一個常數(shù)，所以上式非常容易計算，只要給定一個△t 就可以計算出下一次迭代粒子所在的新位置。

為了約束粒子在反轉(zhuǎn)表面空白區(qū)域的反轉(zhuǎn)問題，作者又考慮了使用內(nèi)部因素Fint(X,t) 。任意選取兩個相鄰的粒子，如果兩個粒子都是 可移動 的，則令二者往相反的方向移動同樣的距離；如果一個是不可移動的，則移動另一個；如果兩者具有相同的高度，則不進(jìn)行移動。
位移量可以通過下式進(jìn)行計算：

式中，d為粒子的位移量；當(dāng)粒子可移動時，b等于1，不可移動時b等于0；pi為p0的相鄰粒子，n是把點進(jìn)行標(biāo)準(zhǔn)化到垂直方向上的單位向量（0，0，1）T。

移動過程如下圖所示：

同時還引入了一個新的參數(shù)rigidness, RI，用來描述粒子的移動次數(shù)，RI為1時，則移動一次，且移動的距離是兩個粒子之間高差的1/2；RI為2時，則移動兩次，且移動3/4。類似的，RI為3時，移動三次，移動距離為兩個粒子高差的7/8。RI的值越大，則“布料”就越硬，一般來說，平地需要設(shè)置RI的值較大，陡坡、山區(qū)需要設(shè)置的RI值較小。

3.3 具體實現(xiàn)

綜上所述，布料模擬的簡要步驟可以總結(jié)如下：
1）通過第三方軟件去除原始LiDAR點云數(shù)據(jù)中的“孤立點”。
2）反轉(zhuǎn)LiDAR點云。
3）初始布料格網(wǎng)，這時用戶需要設(shè)置格網(wǎng)的大小（grid resolution, GR）。初始“布料”的位置通常在最高點的上方。
4）把所有LiDAR點與格網(wǎng)粒子投影到同一個水平面，并找到每一個粒子的最近鄰點（corresponding point, cp），記錄其投影前的高程（intersection height value, IHV）
5）對于每一個可移動的格網(wǎng)“粒子”，計算其受到重力影響產(chǎn)生的位移，并與當(dāng)前粒子對應(yīng)cp點的IHV進(jìn)行比較，如果粒子的高度低于或者等于IHV，則把粒子的高度設(shè)置為IHV并設(shè)置為不可移動點。
6）對于每個格網(wǎng)“粒子”，計算其受到內(nèi)部驅(qū)動因素影響所產(chǎn)生的的位移。
7）重復(fù)上述5）、6）步驟，直到所有粒子的最大高承諾變化足夠小或者迭代次數(shù)到達(dá)用戶的預(yù)先設(shè)置值，則停止模擬過程。
8）計算LiDAR點云與格網(wǎng)粒子之間的高度差異。
9）區(qū)分地面點與非地面點，如果LiDAR點與模擬粒子之間的距離小于預(yù)先設(shè)置的閾值hcc，則認(rèn)為其是地面點，反之則認(rèn)為其為非地面點。

3.4 后處理

上述方法處理過后對于平坦區(qū)域可以取得較好的效果，但是對于“陡坡”區(qū)域仍然會產(chǎn)生較大的誤差，因為“布料”會由于其粒子間的內(nèi)部約束而落在地面測量值的上方。如下圖所示：

為了消除這些誤差，可以通過一個后處理方法來解決。在每一個“可移動”粒子的四個鄰接點中找到“不可移動”粒子，然后比較“最近鄰點”之間的高程差，如果高程差值小于閾值hcp，則把“可移動”粒子移動到地面，并設(shè)置為不可移動。
例如，上圖中的D點，找到了一個不可移動的粒子A，然后比較C,B之間的高程差值（對于D點、A點之間的最近鄰點CP）。如果高程差值小于閾值hcp則待判斷點D就被移動到C點，且被設(shè)置為不可移動。重復(fù)上述操作，直到所有的“可移動”粒子被正確處理。
如果逐行掃描“粒子”來進(jìn)行實現(xiàn)此后處理操作，結(jié)果可能會受到特定掃描方向的影響。因此，作者采用了先建立一個“強連接元素集”（strongly connected components, SCC）包含一個關(guān)聯(lián)的可移動粒子集。

如上圖所示，每個SCC中包括兩種粒子，第一種是至少包含一個近鄰“不可移動”粒子（上圖中的黃色點M1），第二種是近鄰點中不包含“不可移動”粒子（上圖中的紅色點M2）。執(zhí)行過程中使用外部的黃色點為“初始種子點”，通過“廣度優(yōu)先”的方式逐步向內(nèi)遍歷SCC。如上圖中的從1-18點，實現(xiàn)了從邊緣逐步向內(nèi)部靠近處理過程，而不受掃描方向的影響。

4. 算法使用

Cloth Simulation Filter已經(jīng)在各種點云處理軟件，以及MATLAB中都可以直接進(jìn)行使用，對用戶操作友好，非常方便。

各種使用方式的超鏈接如下：
MATLAB中mex插件
cloud compare 中CSF濾波
dll與lib下載
以及源碼都在github上開源了

到此這篇關(guān)于點云地面點濾波(Cloth Simulation Filter, CSF)"布料"濾波算法詳解的文章就介紹到這了,更多相關(guān)python 點云地面點濾波內(nèi)容請搜索腳本之家以前的文章或繼續(xù)瀏覽下面的相關(guān)文章希望大家以后多多支持腳本之家！