使用空載光達點雲求定數值地表高程模型之小波法

(1)

Volume 13, No.2, June 2008, pp. 117-130

使用空載光達點雲求定數值地表高程模型之小波法

鄭邦寧¹ 蔡展榮²

摘要

本文應用小波演算法輔助空載光達點雲過濾，使用小波進行內插生成格網化地形初始面，以網格化資料結構結合Geodesic Dilation法、高程差與其變化率判斷、重新篩選過濾機制，逐步將光達點雲分類為地面點與非地面點，進一步計算產生地表數值高程模型；過程中可參考光學影像特徵，另外針對橋樑區域的光達點雲過濾，應用和橋樑走向正交的掃瞄方式搜尋橋樑種子點，配合區域成長法偵測橋樑位置。從 地形簡單的實驗區一至地形複雜且多斷線的實驗區二的點雲分類成果顯示，分類總錯誤率分別為0.5%和 6%，成果頗佳，構想經實際驗證為可行。

關鍵詞：空載光達、小波、影像、點雲過濾/分類

1.前言

空載光達掃描系統是一個直接定位定向(direct georeferencing)的測量法，能在短時間內透過雷射掃描、GPS 以及 INS/IMU，獲取密集分佈於被探測物表面的大量高精度三維點雲，因此，空載光達成為測製數值(地表)高程模型(Digital Elevation Model DEM)的一個新的測量法。其中，以空載光達點雲製作數值高程模型時，面臨到兩個最大的問題：(一) 如何從複雜的光達點雲覆蓋面上，挑選落於地面的光達點，(二)面對龐大的點雲資料時，如何提升過濾的品質與效率。

為了克服上述二個問題以獲取測製 DEM 之地面點資料，常見的點雲過濾法依不同的操作概念可大致分成四類：(一)斜率為基礎、(二)區塊為基礎、

(三)面為基礎、(四)群集/分塊為基礎。以斜率為基礎的操作概念是建立在非地面點與相鄰地面點之間存有斜率異常變化或斜率驟起之關係，例如：

Vosselman(2000)計算光達點和點之間的斜率，以篩選地面光達點。以區塊為基礎的操作概念是建立在地面點為局部區塊高程值較低的點位，例如：

Petzold et. al. (1999)將搜尋罩窗的尺寸逐步縮小，視

罩窗內與最低點高程值之差值超過門檻值的光達點視為非地面點逐一剔除，以獲取地面點。Zhang and Chen (2003)則逐漸放大搜尋罩窗，透過膨脹 (dilation)和侵蝕(erosion)迭代操作挑選地面點。以面為基礎的操作概念則是建立於非地面點與擬合面之間具有高程異常凸起的特性，例如：Kraus and Pfeifer (1998)以最小二乘法進行曲面擬合，將光達點高程和擬合面高程差超過門檻值之光達點視為非地面點而予以濾除。群集/分塊為基礎的操作則是將光達點高程視為影像灰階值，以影像群集/分塊的方式獲取地面點資訊。例如：邵怡誠及陳良健 (2005)提出以 H-Dome 法進行點雲過濾的演算法，

並於 2006 年改良先前的方法(邵怡誠及陳良健，

2006)，並提出重新篩選地面點的概念，以獲取較佳的地面點分類品質。

為求篩選地面點之工作能以高品質、高效率的方式進行，本文提出以小波重建面(嚴晟瑋，2006) 的方法進行點雲過濾，藉由小波重建面對複雜的空載光達點雲覆蓋面具備良好的描述能力以及過濾點雲雜訊的能力，提供一個和真實環境相近似的擬合面，並根據局部光達覆蓋面的掃描點密度，可將光達所獲取之局部覆蓋面以等間格格點的格式輸

1國立成功大學測量及空間資訊學系碩士

2國立成功大學測量及空間資訊學系副教授

收到日期:民國 97 年 05 月 26 日修改日期:民國 97 年 08 月 02 日接受日期:民國 97 年 08 月 24 日

(2)

出，結合 GD 法、比較基準面、高程差及高程差變化率、及重新篩選地面點等篩選地面點之判斷機制，逐步地提升地面點篩選之成果品質，以獲取測製 DEM 之可靠的光達地面點資料。

此外，本文也嘗試將彩色影像納入，並藉由彩色影像提供的多波段資料，例如：近紅外光(NIR)、

紅光 (R) 、以及綠光 (G) ，針對具備多回波特性 (multiple return)的植被區光達點以 NDVI 或 GI 先行篩選植被區光達點，再以本文提出之判斷法剔除植被區非地面點，然後，再透過非植被區光達點及植被區地面點之光達點重新做為輸入資料，以建立較能表示地面資料之小波重建面供後續 GD 法、比較基準面、高程差及高程差變化率及重新篩選等判斷條件，以提升地面點篩選的成果品質。圖 1 為本文提出篩選地面點點雲之流程圖。

2.實驗區資料格式

本文使用的空照彩色影像(Vexcel UltraCamD) 拍攝時間為 2005 年，檔案格式為*.tif，空載光達 (ALTM3070)掃描時間為 2006 年，檔案格式採*.las，

假設光達和影像紀錄的地表內容大致相同。如圖 2 和圖 3 所示，實驗區一是台南一中校園(操場、體育館、籃球場、網球場)，實驗區二是成大校園附近的四維地下道(民族路與前鋒路、北門路之交叉路段)，兩者的空照影像大小都是 551x521 個像元，

包括 R/G/B 彩色影像以及 NIR/R/G 近紅外光影像 (摘錄自 UltraCamD 使用手冊)；另外，兩者的光達 點涵蓋範圍都是 50m x 50m，相鄰掃描線間距約 30cm，相鄰掃描帶間距約 70cm，光達點的平面與 高程之先驗精度分別約±25cm 與±10cm(摘錄自 Optech ALTM 3070 使用手冊)。其中，實驗區一有 14087 個光達點(last return)，平均點間距約 40cm(約 2.4 點/m)，而實驗區二有 18001 個光達點(last return)，平均點間距約 37cm(約 2.7 點/m)。

圖 1 點雲過濾流程圖

圖 2 實驗區一的空照彩色影像(左)及光達點雲三維分布圖(右)，Max(H)-Min(H)= 52.073－39.775 ＝12.298 m

(3)

圖 3 實驗區二之空照彩色影像(左)及光達點雲三維分布圖(右)，Max(H)-Min(H)= 45.549－30.236＝15.313 m

3.點雲資料處理之操作原理

3.1 光達點預處理

光達資料檔*.las 是一個 binary 格式檔案，記錄每一個光達點的三維空間坐標、反射強度值和光達回波序號。以 ALTM 3070 為例，最多可紀錄到四個回波的資料。假設每一個最後回波(last return)有較大的機會落於地表面，因此，本文以光達回波序號做為挑選最後回波值的依據，所處理之光達點資料均為最後回波值的三維坐標。

3.2 小波重建面

本文使用小波法(嚴晟瑋，2006)來進行局部區塊空載光達點雲的最小二乘面擬合，主要考量到小波重建面之真實性及實用性二方面：(一)真實性：

小波重建面能將複雜且具有碎形幾何特性的光達覆蓋面完整的描述。(二)實用性：為配合影像分塊處理之 GD 法，小波重建面能配合點雲資料過濾流程之需求，將某一塊連續面函數以規則網格格點的方式記錄三維坐標。

由於拍攝影像和光達掃描時間不同，使得小波重建面與局部影像內容未能完全一致，例如：行駛於車道上之移動汽車。整體而言，以本文用以做為過濾面的兩個實驗區小波重建面為例，原始光達點 間距約 40cm，小波重建面輸出網格大小 62.5cm 的 DSM，做為點雲過濾判斷的基礎面。

3.3 GD 法

GD 法 (Geodesic Dilation) 是一種形態學方法，常用於影像分塊，其操作概念如圖 4。原始訊號面 I 下降一個高度值 h 後，得到一個操作面 (I-h)，從操作面出發，藉由罩窗反覆迭代滑動，不斷地提高罩窗中心點位高程，直到停止迭代後，可得到 GD 法的重建面，將原始訊號面減去重建面後，可將隆起區分塊出來。根據二原則來提高操作面：(一)罩窗中心點高程不得高於相對應之原始訊號面；(二)罩窗中心點高程不得高於罩窗內相鄰點位之高程最大值。以此原則反覆迭代計算，直至收斂為止，才終止提高操作面之迭代計算。

使用此法分類點雲，需透過此區的地上物高度值給定下降參數值 h，如圖 5 所示，如此一來，

地面分塊的成果較佳，否則提高操作易受錯誤的下降參數值，造成過度提升(無法將高程隆起處準確切割出來)或是高程提升不足(大多數地面點被視為非地面點剔除)，而產生錯誤的點雲分類。

本文測試結果顯示，若下降參數值太小，僅能將吉布斯效應影響較大的建物邊緣及少數植被區挑選出來；若下降過多，往往會產生操作面爬升不足，產生過度過濾的現象。因此，本文使用小波重建面做為 GD 法之起始訊號面，考量高程不連續變化處之吉布斯效應的影響，本文在設定下降參數 h 時，取大於最大地上物高度之整數值為設定標準，

(4)

圖 4 GD 法操作展示圖(Vincent, 2003)

圖 5 光達影像剖面圖以分割出高程隆起區。此外，考量實驗區之建物或

高聳樹木等的高度超過 2.5m，因此，本文以 2.5m 做為 GD 法挑選高程隆起區域之高程門檻值，意 即，凡超過 2.5m 之區域即視為高程隆起區，反之，

則視為地面區。

罩窗尺寸較大時，植被區容易受到罩窗內鄰近點高程的影響，使得高程隆起區無法有效挑選出。因此，本文以 3*3 罩窗做為本文用於 GD 法之滑動罩窗，俾計算正確的高程重建面。

3.4 比較基準面

設計比較基準面的主要構想是基於提供一個較為接近地面之高程面做為重新篩選地面點之用，因此，根據不同計算區內建物的一般尺寸、光達掃描線與線間距、帶與帶間距等資訊，定義網格尺寸，做為比較基準面的網格大小。此預先定義的規則網格，高程值以 GD 法獲取的地面候選點為基礎，以網格尺寸做為搜尋半徑，將半徑內地面候選點高程最低點做為該格點之高程值。若格點的搜尋

半徑內無地面候選點，則以該格點做為 3*3 之罩窗中心，取相鄰八個點位中具有地面候選點之格點高程累計取平均，以做為該格點之高程值。若 3*3 罩窗內無搜尋到地面點，則再以 3*3 罩窗滑動計算網格點高程值，直到每一個格點皆有高程值，才終止罩窗滑動。

3.5 高程差及高程差變化率判

由於 GD 法分類時，往往無法將建物邊緣之光達點剔除，使得比較基準面在搜尋局部最低點時仍存有非地面點紀錄於比較基準面上的情況。針對這類情況，本文採用式(1)高程差及高程差變化率的判斷條件，再次檢查比較基準面上是否存有非地面點的資訊並予以剔除。

為配合比較基準面之規則網格結構，且考量比較基準面上之高程值為局部地面點之最低點高程值，本文以 3*3 的罩窗進行高程差和高程差變化率之比較，俾以判斷光達點是地面點或是非地面點。

斷

(5)

高程差變化率門檻值高程差門檻值，

呈對角線，

網格點高程值待處理點

網格中心點

，其中，

非地面點非地面點

=

⎩⎨

⎧

=

−

=

= +

−

=

−

=

−

=

−

=

−

=

>

=

>

thrS thrH

i i

i j i

P P P

H H H

P P P

H H H

H S S dS

H H S H H S end

P

thrS dS

thrH S

elseif P

thrS dS thrH

S if

j i

P P P

j i

P P P

P j i ij

P P j P P i j

ji j

i

ij i

j i

8

~ 5 , 4

4

~ 1 , 4 ,

, ,

) (

&

0 0

0

0 0

由式(1)可以看出，若待處理點高程和網格中心點高程的差值小於或等於 thrH，則待處理點有機會為地面點，若單純以 thrH 做為約制條件，則高程變化一致之上坡面很可能被視為非地面點剔除，為降低上述錯誤過濾的情況發生，本文另外納入高程差變化率(thrS)做為挑選非地面點的判斷條件。由 thrS 之計算可以得知當 thrS 值愈接近 0，則操作視窗內之剖線方向為地勢起伏平緩處。因此，透過 thrS 沿著操作視窗剖線方向的高程差約制條件，可降低上坡面被視為非地面點剔除的情形。如圖 6 所示。

圖 6 之 x 軸為水平方向，z 軸為高程，＋符號代表高程上升，－表示高程下降。以圖 6(a)可以看出若為坡地地形，以式(1)之處理流程，P_i_將被視為地面點保留，而圖 6(b)則P_i _{視為非地面點剔除。}

圖 6 高程差判斷說明圖

3.6 重新篩選

經由比較基準面、高程差及高程差變化率的判

斷後，比較基準面以較貼近地面的姿態，提供重新篩選光達地面點之用。首先，根據光達的平面坐標透過雙線性內插計算各個光達點在比較基準面上的高程值，再將各光達點高程和內插高程之差值統計於直方圖內，根據實驗區之光達點覆蓋面定義出高程差門檻值，若差值滿足門檻值的條件，視為地面點保留，否則，視為非地面點剔除。最後，將重新篩選後的地面點輸入小波重建面之計算來求定 DEM。

3.7 影像介入輔助篩選植被區光達點

為降低植被區多重回波效應對點雲分類品質之影響，本文嘗試以空照彩色影像所提供的波段資料，以遙測常用以偵測植被區的二種指標—NDVI (Normalized Difference Vegetation Index) (Lillesand and Kiefer, 2000)、GI (Greenness Index) (Niederöst, 2001)，挑選影像植被區，其中，以GI 搜尋植被區是建立在影像上綠色像元處多為植被區發生處的概念上。計算公式分別如下：

^NDVI ⁼ _NIR^NIR₊⁻_R^R (2)

R G

R GI G

+

= − (3)

式中，NIR、R、G 分別為物體對近紅外光、紅光、

綠光之反射率，NDVI 及 GI 的值域皆為[-1,1]。以 NDVI 或是以 GI 判斷影像植被區時，主要藉由植物對於近紅外光(0.7~1.3 mμ )或綠光(0.5~0.6 mμ )有較佳的反射率，且與紅光(0.6~0.7 mμ )反射率差值相比後，植被區會因較大的反射率反差挑選出來。

在空載光達點雲與影像之間的套合關係上，經空三計算而獲取影像外方位元素，透過共線式求得光達點雲和影像之間的套合偏差量約為 1~4 個像 元，對應於地面約 12~48cm，以光達點平均點間距 40 cm 而言，此套合結果足以提供影像篩選植被區 的光達點。

如圖 7(b)所示，Y1 到 Y2 實線為小波重建面之剖面，黑點則表示落於圖 7(a)所示的 Y1、Y2 剖線 (1)

(6)

的兩側各 1m 範圍內之光達點，從剖面圖可看出植 被區的非地面光達點與小波重建面之間的關係，受到多重回波效應的影響，使得植被區非地面光達點高程與小波重建面高程之差值多大於或等於三倍光達高程先驗精度，因此，當影像資料介入時，本文使用 NDVI 或是 GI 指標搜尋植被區，再以共線式(collinearity condition equation)倒投影的方式將落於影像植被區之光達點挑選出來，並以三倍的光達先驗高程精度做為門檻用以剔除植被區非地面光達點。

關於 NDVI 或是 GI 門檻值的給定，本文是透過直方圖輔助人工肉眼判讀的方式來決定，以實驗區一為例，將 NDVI 值分做六個值域區間並以灰階表示，如圖 8(b)所示，配合 NDVI 直方圖，如圖 9(a)，

從圖 8(a)可以看出實驗區一之植被區面積約佔影像的 1/3，從直方圖可知，總像元數之 1/3 約落在 NDVI 為 0.2-0.4 的區間內，透過人工的調整以及肉眼的檢查，以 NDVI 大於或等於 0.2 界定影像植被區最為恰當，其成果如圖 9(b)所示，白色表示 NDVI 值大於或等於 0.2 之植被區，黑色為非植被區。

圖 7 小波重建面之 Y1、Y2 剖線位置圖(a)與 Y1、Y2 高程剖面圖(b)

圖 8 實驗區一 NIR/R/G 影像(a)及 NDVI 值域分布圖(b)

B1 B2

A1 A2

Y1 Y2

(a) (b)

(7)

圖 9 實驗區一 NDVI 統計直方圖(a)和植被區二元影像(b) 同樣地，以 GI 指標篩選影像上綠色像元時，

本文採相同方式，直方圖輔助人工檢視得知，設定 0 為綠色像元門檻值，影像上綠色像元的搜尋成果為最佳。

式(4)為本文使用 NDVI 或 GI 進行植被區點雲分類之判斷法則，式中 thrV 即以上述方式給定的影像植被區門檻值。

小波重建面高程光達點高程值

門檻值或

用來判斷植被區影像之

值或於影像上的對應像點之

光達點其中，

植被區地面點植被區非地面點

=

−

=

≥

DSM wavelet P

i P

i P P P

H H

H H dH

GI NDVI thrV

GI NDVI P

V P end

end V else V

dH if

thrV V

for

i i i

_

_ 0

|

3σ

3.8 橋樑偵測

如圖 10 所示，由於橋樑兩端橋樑面和其外接路面的高程常相等，造成橋樑面的光達點常被誤判為地面點，卻未找出橋樑下方的地面點(如圖中虛線處)。因此，為了增加橋樑區篩選地面點之正確率，本文提出了一個橋樑偵測的方法。

圖 10 橋樑剖面示意圖

此方法需要給定橋樑面的離地高度並且給定和橋樑走向近似正交的掃瞄搜尋方向，俾搜尋橋樑種子點，從種子點出發，往上下左右四個方向以區域成長法的方式，定義橋樑所在區。詳細處理步驟如圖 11 所示。

圖 11 橋樑偵測演算流程圖小波重建面

1. 給定掃描方向與橋樑高度 2. 搜尋小波重建面上可能為橋樑

之種子點

給定用以偵測橋樑之區域成長法的高程差門檻值

以區域成長法判斷小波重建面之橋樑位置 (4)

(8)

3.8.1 篩選橋樑種子點

根據橋樑和鄰近地面具有高程驟起之特性，本文用來偵測橋樑之操作視窗需和橋樑的走向相互正交，以來回掃描的方式獲取橋樑種子點，若橋樑走向為東西向，則以 5*1 縱視窗南北向的方式掃描；若橋樑為南北向，則以 1*5 橫視窗東西向掃描，以獲取橋樑種子點。判斷橋樑種子點的法則如式(5)所示。

以實驗區二為例，由於實驗區二之橋樑為南北向，所以採用 1*5 的橫視窗來搜尋橋樑種子點；透過實地目測及光達三維分布剖面圖輔助判斷，給定 式(5)的橋樑高度值(Bh)為 6m。此外，由於橋樑南 側有高聳植被區，為避免此處的吉布斯效應之影響，因此以較靠近橋樑中心之格點做為區域生長法之種子點，如圖 12 所示。

透過式(5)之演算內容可得知，若建物高度與橋樑和地面之間的高程差相同，則建物頂面上的光達點會被誤判為橋樑種子點。實驗區二內因橋樑所佔的面積較大，且其周邊建物高度未必與橋樑高度值 (Bh)相同，因此，圖 12 可以看出種子點落於橋面上確實佔多數。將所有橋樑種子點位(X,Y)上的小波重建面高程 Z(X,Y)值以直方圖表示，如圖 13 所示，挑出橋樑種子點之高程值有高一致性，對於後續用於偵測橋樑之區域成長法能夠提供一個較佳的起始生長環境。

預設之橋樑高度值

網格點高程值其中，

橋樑種子點之高程橋樑種子點之高程

=

>

−

>

=

>

−

<

+

−

− +

−

+ +

− +

−

Bh

bH bH

end bH

Bh bH

bH if

bH bH

bH bH elseif

bH

Bh bH

bH if

bH bH

if

i i

i i i

i i

i i i

i i

2 2

2 2 2

2 1

1 2

2 2

1 1

2

~

&

圖 12 實驗區二之小波重建面高程影像(81*81 個像元，對應的地元大小為 0.625m X 0.625m)及橋樑種子點(白點)

圖 13 橋樑種子點位(X,Y)上的小波重建面高程 Z(X,Y)之直方圖

3.8.2 區域成長法偵測橋樑位置

獲取橋樑種子點後，以式(6)進行區域成長法，

找出橋樑所在位置。圖 14 為種子點(B⁰)操作視窗。

由於橋樑表面又為一個近似平坦的區域，因此，從種子點出發，以圖 14 為操作視窗來搜尋鄰近點，

將滿足高程差門檻值 dB 之鄰近點 Bⁱ視為下次搜尋橋樑之新種子點，如此反覆迭代計算，直至無新的橋樑種子點產生為止，所生長出之種子點區即視為橋樑區。

(5)

(9)

橋樑區高程差門檻值種子點高程

待處理點高程待處理點

種子點其中，

橋樑種子點

=

≤

−

dB H H

i B

B end

end B

dB H

H if

B all For

B B i i

B B

i i

0 0

4

~ 1 ,

|

0 0

圖 14 種子點操作視窗

以實驗區二為例，以第 3.8.1 節搜尋到的橋樑種子點，透過式(6)的生長條件，以人工肉眼檢查之 試誤法發現，當橋樑區高程差之門檻值為 0.4m 時，能將大多數屬於橋樑的區域挑選出來，如圖 15 所示。

圖 15 大小圓圈內的範圍為非橋樑區，圖 16(上) 顯示實驗區二 X1、X2 剖線位置。由於大圈勾勒的範圍對應圖 16(下)水平刻畫 2.5~5 之範圍，可以知道實驗區二落於大圈內的區域其高程與橋樑面高程近似相等，因此在式(6)的判斷下，最大圈內的區域歸類為橋樑區。第二大圈的位置則對應圖 16 (下) 紅圈圈選的下坡區。由於式(6)只根據鄰近點與種子點兩者的高程差來判斷，所以也將上(下)坡區誤判為橋樑區。

另外，對偵測橋樑區來說，圖 15 三個較小圓圈內的地區為偵測錯誤區，但以式(6)的判斷下，小圓圈內種子點所生長出的區域為高起之樓頂，對於篩選地面點而言，此法能提供部分濾除非地面點之功能。

圖 15 區域成長法偵測得到的橋樑區位置圖(小波重建面之高程影像：81*81 個像元，對應的地元大小為 0.625m X 0.625m)

圖16 實驗區二 X1、X2 剖線位置圖(上)和此剖線區的小波重建面剖面圖(下)

(6)

(10)

4.實驗成果與分析

本文使用 MLidarView 軟體(Wang and Tseng，

2005)用以展繪不同視角的光達點雲三維分布圖並配合影像之輔助判別來定義地真參考點，根據分類成果來製作表 1 並計算分類成果之漏授誤差、誤授誤差及總錯誤率，據以分析點雲分類成果的品質。

漏授誤差(第一類型錯誤，Type I Error，或稱為 Omission Error)之定義乃光達點地真資料認定為地面點，卻被分類法遺漏而未被分類為地面點(被分類法視為是非地面點)的錯誤百分率。誤授誤差(第二類型錯誤，Type II Error，或稱為 Commission Error) 則定義為光達點地真資料為非地面點，卻被分類法誤判為地面點之錯誤百分率。總錯誤率(Total Error) 則定義為光達點分類結果不屬於地真資料之錯誤百分率。

表 1 點雲分類精度統計表分類結果

參考點地面點非地面點

地面點 a b

非地面點 c d

※漏授誤差= b/(a+b)

※誤授誤差= c/(c+d)

※總錯誤率= (b+c)/(a+b+c+d)

兩個實驗區之點雲分類精度依照有無影像介入輔助判斷分別統計，表 2 和表 3 的 PART A 顯示以式(4)對植被區光達點雲初步分類的成果。再以小波重建面為基礎，透過 GD 法之操作，將大多數屬於非地面點之點位剔除，其分類成果記錄於 PART

B，而 PART C 則紀錄重新篩選地面點之成果；圖 17 及圖 18 為無影像介入時，透過本文提出之地面點篩選法所得之點雲過濾成果。

實驗成果顯示，本文提出小波演算法來進行點雲分類，不僅能過濾光達點位雜訊，且由於點雲密度高，所以也能同時以規則網格的方式逼真的表現出光達掃描的複雜覆蓋面之細節特徵，快速獲取初步過濾的地面光達點。

然後，配合影像分塊常用的 GD 法，以人工半自動化的方式給予隆起地物的高度做為下降參數值，從實驗的結果中發現，以 3*3 的罩窗迭代滑動，能降低 GD 法在提昇操作時受到相鄰點高程驟起的影響，而減少隆起地物處無法剔除的機會，再 佐以地物高度 2.5m 做為判斷隆起的非地面區之門 檻值，獲得較佳的地面點分類成果。

此外，實驗結果顯示，GD 法能濾除大多數的非地面點，但卻無法有效的濾除高程不連續的隆起物邊緣之非地面點，針對此現象，以 GD 法分類為地面點的點雲為準，以局部地面點之最低點建置比較基準面的方法，能有效改善建物邊緣非地面點分類錯誤的現象，進一步佐以高程差及高程差變化率的條件過濾比較基準面上的非地面點，藉由雙線性內插計算光達點在比較基準面上之高程、以及它與 本身高程之差值，發現以 1m 做為高程差門檻值篩 選地面光達點的成果最好。經過光達點重新篩選 後，兩實驗區的總錯誤率由原先的 3.6%和 14.8%分 別降低至 0.5%及 5.9%，顯示本文提出之比較基準 面供重新篩選地面點的方法有助於提升點雲過濾的品質。

圖 17 實驗區一重新篩選後之地面點(左)與非地面點(右)

(11)

圖 18 實驗區二重新篩選後之地面點(左)與非地面點(右) 表 2 實驗區一之點雲分類結果

影像輔助無影像輔助

(%) NDVI GI

Type I Err. 2.3 1.8

Type II Err. 20.0 21.0 A

Total Err. 4.9 3.9

************ GD 法過濾 ************

Type I Err. 0.9 1.0 1.8

Type II Err. 13.9 14.7 21.8

B

Total Err. 2.1 2.2 3.6

************ 重新篩選 ************

Type I Err. 0.2 0.2 0.2

Type II Err. 3.8 3.8 3.8

C

Total Err. 0.5 0.5 0.5

表 3 實驗區二之點雲分類結果

影像輔助無影像輔助

(%) NDVI GI

Type I Err. 8.6 16.7

Type II Err. 64.9 70.4

A

Total Err. 36.8 42.7

************ 橋樑偵測和 GD 法過濾 ************

Type I Err. 11.6 13.5 11.7

Type II Err. 18.6 17.6 18.6

B

Total Err. 14.2 15.4 14.8

************ 重新篩選 ************

Type I Err. 10.0 10.6 10.3

Type II Err. 0.8 1.1 0.8

C

Total Err. 5.7 6.2 5.9

(12)

針對橋樑區域的光達點雲過濾處理，本文提出一個與橋樑走向垂直的掃描方向進行偵測，經由人工輸入橋樑高度，找出橋樑的種子點，根據局部橋樑面近似是一個等高面之幾何特性，以區域成長法找出橋樑區的位置。由實驗區二偵測橋樑的成果可知，本文透過橋樑邊緣之吉布斯效應能偵測出位在橋樑上的種子點，再以區域成長法找出與種子點高程相近似的區域，而此區域即為欲搜尋之橋樑區。

實務使用上，仍須以人工檢查的方式來給定生長種子點之門檻值，否則，本文使用的區域成長法將會在地勢緩慢變化之坡面產生誤判。

當影像納入時，本文的實驗影像以 NDVI 門檻值為 0.2 或 GI 門檻值為 0 時較能挑選出影像植被區或綠色像元，且在找出植被區或綠色像元區後，

從實驗區一的成果發現，以光達高程先驗精度的三倍為高程門檻值，能將大多數落於植被區的非地面點剔除，但在實驗區二的分類結果卻不如實驗區一，此原因歸咎於實驗區二之覆蓋面 DSM 複雜且多斷線、植被生長於橋面上，使得實驗區二的分類結果不佳。

但是，光達點雲經由影像 NDVI 或 GI 介入輔助、再以 GD 法和橋樑偵測暨橋樑區點雲過濾分類、及最後重新篩選等三個步驟的過濾成果來看，

實驗區一從總錯誤率 4.9%、2.1%降低到 0.5%，實 驗區二從 36.8%、14.2%，再降低到 5.7%，顯示出 本文提出的小波法具有層層過濾點雲且逐步提升點雲分類精度的良好效果。此外，經由重新篩選 後，從地形單調的實驗區一分類總錯誤率 0.5%和 地形複雜多變化的實驗區二分類總錯誤率 5.7%，

配合圖 19 提供實驗區二細節的成果資訊來看，小波法能篩選出一定品質程度的地面點雲。

若以 NDVI、GI 與無影像介入重新篩選的成果來判斷優劣，在不考慮人工誤判地真點資料的前提 下，從實驗區一 0.5%、0.5%、0.5%和實驗區二 5.7%、6.2%、5.9%的成果來看，此成果說明影像介 入輔助植被區光達點的成效不彰顯，此乃比較基準面造成之結果，由於比較基準面之建立是透過 GD 法地面點的局部最低點而來，使得在搜尋半徑內的所有光達點均和該格點最低點有關，因此，無論有

無影像輔助挑選植被區光達點，所有光達點均需和比較基準面相比，因此，有影像和無影像介入時之點雲分類成果具有高一致性。

圖 19 實驗區二非地面點俯視圖(上)及細部展示圖

5.結論與建議

總體說來，本文研究的具體結論包括：

1. 從地形簡單的實驗區一至地形複雜且多斷線的實驗區二的點雲分類成果品質來看，本文應用小波法過濾空載光達點雲，確實具有良好的過濾分類效果。

2. 從兩個實驗區的分類總錯誤率 0.5% 和 6%，顯示出本文以 GD 法、高程差判斷以及重新篩選等過濾步驟，能獲得不錯的篩選地面點效果，

同時能兼具電腦的處理效率。

3. 由實驗區二偵測橋樑的成果可知本文提出橋樑偵測的方式，可將多數屬於橋樑的非地面點濾除，而獲得橋樑下真正屬於地面的地表高程模型。

對於未來的工作與建議，由於本法於實務上操 (下)

(13)

作時，涉及的資料繁多，在影像上，需透過人工給定植被區或綠色像元之門檻值，於光達資料處理上，需透過光達剖面圖的輔助給定 GD 法下降參數值，且在偵測橋樑上，針對東西向或南北向的橋樑，演算法需給定橋樑的高度值以偵測位在橋樑上的種子點，以及區域成長法生長的高程差門檻值等，在未來研究上，可將重點放在如何以小波重建面逐塊過濾光達點雲的方式，自動化給定植被門檻值、下降參數值或是橋樑高度值等，或是在橋樑偵測上納入橋樑寬度的判斷條件，以高自動化、高成功率之演算流程，獲取製作 DEM 之地面光達點雲資料。

6. 致謝

本研究承蒙國家科學委員會支助，專案編號 NSC 95-2221-E-006-339，謹此致謝。承蒙審查委員們提供寶貴的修正意見，讓本文更臻完善，亦謹此致謝。

參考文獻

邵怡誠、陳良健，2005。利用光達資料於 DEM 生產及房屋偵測之研究，航測及遙測學刊，第十卷，第一期，民國 94 年 3 月，第 47-58 頁。

邵怡誠、陳良健，2006。空載光達點雲於 DEM 自動生產與精度評估--使用 ISPRS 測試資料為例，航測及遙測學刊，第十一卷，第一期，

民國 95 年 3 月，第 1-12 頁。

嚴晟瑋，2006。以小波為基礎的最小二乘物面重建演算法，國立成功大學測量及空間資訊學系碩士論文。

Kraus, K. and N., Pfeifer, 1998. Determination of terrain models in wooded areas with airborne laser scanner data, ISPRS JPRS, Vol. 53, pp. 193-203.

Lillesand, T. and R. Kiefer, 2000. REMOTE SENSING AND IMAGE INTERPRETATION, 4th edition, John Wiley & Sons, Inc., ISBN 0471255157. page 448.

Niederöst, M., 2001. Automated update of Building Information in Maps Using Medium-scale Imagery

(1:15,000), Automatic Extraction of Man-Made Objects from Arial and Space Images, Vol. 3, pp.161-170.

Optech Incorporated, 2007. ALTM 3100 System Specifications, downloaded on May 15, 2007, from http://www.optech.ca/pdf/Spacs_altm_3100.pdf Petzold B., P. Reiss and W. Stossel, 1999. Laser scanning – surveying and mapping agencies are using a new technique for the derivation of digital terrain models, ISPRS Journal of Photogrammetry and Remote Sensing, Vol. 54, pp.95-104.

Vexcel Corporation, 2005. UltraCam Technical Specifications, downloaded on May 15, 2007, from http://www.vexcel.com/products/photogram/

ultracam/downloads.html

Vincent, L. 2003. Morphological Grayscale Reconstruction in Image Analysis: Applications and Efficient Algorithms, IEEE Trans. On Image Processing, Vol. 2, No. 2, pp.176-201.

Vosselman, G., 2000. Slope Based Filtering of Laser Altimeter Data, IAPRS, Vol. XXXIII, Part B3, Amsterdam, The Netherlands, pp.935-942.

Wang, M. and Y.-H. Tseng, 2005. Automatic 3D Feature Extraction from Structuralized LIDAR Data, ACRS 2005, 26th Asian Conference on Remote Sensing, 7-11 Nov 2005, Hanoi, Vietnam.

(CD-ROM)

Zhang, K. and S-C. Chen, 2003. A Progressive Morphological Filter for Removing Non-Ground Measurements from airborne LIDAR data, IEEE Transactions on Geoscience and Remote Sensing, Vol. 41, No. 4, April 2003, pp.871-882.

(14)

Wavelets-Based Method for Determining DEM by Using Airborne LiDAR Data

Pang-Ning Cheng¹ Jaan-Rong Tsay²

ABSTRACT

In this study, the wavelet approach is used for airborne LiDAR data filtering. An initial digital surface model (DSM) is generated by wavelet approach. The method utilizes then geodesic dilation, height difference and second order difference of heights of neighboring grid points, and then re-selects LiDAR points by generating a reference surface to select more reliable terrain points. Finally, those selected terrain points are used to generate the digital elevation model (DEM) in that area. In the process, digital aerial color images and their derived feature figures such as NDVI and GI can be also applied. In addition, a method is proposed to exclude those non-terrain points in a bridge region. A point scanning method is presented to sort the seeds located on the bridge surface. These seeds can be adopted for region growing process. Then all non-terrain points within the bridge area can be found. Both test area I with simpler surface topography and test area II with more complicated surface topography as well as break lines show that the total error of LiDAR point classification is 0.5% and 6%, respectively. They verify the method is applicable for LiDAR point filtering and DEM generation.

Keywords:

Airborne LiDAR System, Wavelet, Image, Poiut cloud Filtering /Classification

1 Master Student, Department of Geomatics, National Cheng Kung University

2 Associate Professor, Department of Geomatics, National Cheng Kung University

Received Date: May. 26, 2008 Revised Date: Aug. 02, 2008 Accepted Date: Aug. 24 2008

使用空載光達點雲求定數值地表高程模型之小波法