中華民國一○五年八月 DOI: 10.6574/JPRS
Journal of Photogrammetrty and Remote Sensing
Volume 21 No.1 August 2016
Published by Chinese Society of Photogrammetrt and Remote Sensing
Journal of Photogrammetry and Remote Sensing
發行人：史天元
出版者：中華民國航空測量及遙感探測學會 地址：台北市文山區羅斯福路五段 113 號三樓 信箱：台北市郵政 93-158 號信箱
電話：886-2-8663-3468 886-2-8663-3469 傳真：886-2-2931-7225
電子信件：[email protected] 網址：http://www.csprs.org.tw
PUBLISHER:Peter Tian-Yuan Shih
PUBLISHED BY: Chinese Society of Photogrammetry and Remote Sensing
Address: 3F, No.113, Sec.5, Roosevelt Road, Taipei, Taiwan Mail Address: P. O. Box. 93-158, Taipei, Taiwan
Tel: 886-2-8663-3468 886-2-8663-3469 Fax: 886-2-2931-7225
E-mail：[email protected] Web Site：http://www.csprs.org.tw 總編輯
曾義星
國立成功大學測量及空間資訊學系 電 話：886-6-275-7575 分機 63835 傳 真：886-6-237-5764
電子信件：[email protected]
EDITOR-IN-CHIEF Yi-Hsing Tseng
Department of Geomatics, National Cheng Kung University Tel: 886-6-275-7575 ext. 63835
Fax: 886-6-237-5764
E-Mail: [email protected] 編輯委員(依中文姓氏筆劃排列) EDITORIAL BOARD
王素芬 (彰化師大) 王聖鐸 (師範大學) 江凱偉 (成功大學)
S. F. Wang (National Changhua University of Education) S. Wang (National Taiwan Normal University)
K. W. Chiang (National Cheng Kung University)
何宗儒 (海洋大學) C. R. Ho (National Taiwan Ocean University)
林昭宏 (成功大學) C. H. Lin (National Cheng Kung University)
邱式鴻 (政治大學) S. H. Chio (National Chengchi University)
林老生 (政治大學) L. S. Lin (National Chengchi University)
林唐煌 (中央大學) T. H. Lin (National Central University)
周天穎 (逢甲大學) T. Y. Chou (Feng Chia University)
洪榮宏 (成功大學) J. H. Hong (National Cheng Kung University)
徐百輝 (臺灣大學) 陳朝圳 (屏東科大)
P. H. Hsu (National Taiwan University)
C. T. Chen (National Pingtung University of Science and Technology)
張中白 (中央大學) C. P. Chang (National Central University)
黃金聰 (臺北大學) J. T. Hwang (National Taipei University)
曾義星 (成功大學) Y. H. Tseng (National Cheng Kung University)
詹進發 (政治大學) J. F. Jan (National Chengchi University)
楊明德 (中興大學) M. D. Yang (National Chung Hsing University)
蔡富安 (中央大學) F. Tsai (National Central University)
蔡榮得 (中興大學) J. D. Tsai (National Chung Hsing University)
封面照片說明 About the Cover
歷史航照影像是過去地表景象的忠實紀錄，並且直接地記錄過去某個時間點的地表現象，包括當時
的土地利用、自然環境、道路及聚落等空間分布情形。封面為1947、1948年間美軍航拍任務所拍攝的歷
史航照影像之糾正對位成果，並且套疊於2009年福衛二號2公尺解析度衛星影像，藉以呈現過去真實的地
表樣貌，比較今昔之變遷。
Volume 21, No.1, 2016, pp. 1-12 DOI：10.6574/JPRS.2016.21(1).1
1中國文化大學地質學系 助理教授 收到日期：民國 102 年 08 月 23 日
2中央研究院地球科學所 副研究員 修改日期：民國 102 年 11 月 08 日
3台北科技大學土木與防災研究所 副教授 接受日期：民國 105 年 08 月 02 日
4經濟部中央地質調查所 技士
*通訊作者, 電話, (02)2861-0511#26132, E-mail: [email protected]
應用空載光達數值地形模型於基隆河之 河流地形研究
陳柔妃
1＊詹瑜璋
2張國楨
3謝有忠
4摘 要
本研究利用空載雷射掃描技術(Airborne LiDAR)製作高精度數值地形模型，探討基隆河流域河階分佈 及水系發育，結果顯示河川水系受控於地質及地形條件，主要支流發育在中游的八堵向斜，北側的集水 區較南側平緩，對稱性的水系暗示集水區內並無差異性的地殼變動。基隆河自暖暖到侯硐之間，各次集 水區之面積、形態差異較大，較無法以單一河流階地或構造整體抬升來解釋，而在暖暖附近之坡地地形 具有崩塌特徵，導致古山崩造成河道的堰塞並影響次集水區發育。經野外勘查結果，判定本區可能坡面 破壞形式以岩石傾翻(toppling)或楔型岩石崩落(wedge failure)型態存在。
關鍵字 : 空載雷射掃描技術、數值高程模型、構造地形、古山崩
1. 前言
河流地形常受到不同的岩石特性、地質構造所 影響，透過研究河流水系的幾何分佈與剖面型態特 徵，可瞭解河流地形受到相關地質作用的影響。長 久以來，有關於基隆河流域發育、河流襲奪的相關 議題一直受到學者們的注意(早坂一郎，1930；林 朝棨，1957；Hsu, 1974；周淑文和鄧屬予，1998；
劉明錡，2004)，尤其是在三貂嶺至瑞芳一帶河流 水系呈現出 180 度反向流路、河流兩側遍佈不同比 高的大小河階地、以及存在於河流縱剖面的數個遷 急點等，似乎都隱藏了基隆河發育與地形、地質以 及北台灣地體構造的密切關係。
本研究藉由 LiDAR DTM 去除地表人工建築 物以及植被的功能，完整呈現基隆河流域地區原始 地表起伏的特性。以 2 公尺網格之 LiDAR DTM 進 行本研究區之水系分析、河階分析、河流縱剖面、
以及次集水區盆地分析等。本研究將分析基隆河流
域之水系分佈、河階地分佈、以及次集水區盆地等，
精確建立基隆河的流路變遷，以及進一步瞭解基隆 河階地的形成特性。最後綜合 LiDAR 分析結果與 野外實地查核等相關資料，提出基隆河位於暖暖－
侯硐段，存在的狹長型等高河階地可能形成的原 因。
2. 基隆河流域地形與地質 背景
基隆河位於台北盆地東北為淡水河系三大支 流之一，其源頭於平溪鄉菁桐山，上游流路從菁垌 坑一路流向東北經平溪、六分寮，自三貂嶺、侯硐、
瑞芳至八堵間，呈現近乎直角大轉彎，峽谷、瀑布 與壺穴地形多集中於此河段(圖 1)。瑞芳以西中下 游地區地勢平緩，曲流地形顯著，且因河流下切作 用旺盛，兩岸皆有階地分布，掘鑿曲流和寬廣的河 階地形是基隆河這段地形的主要特徵。八堵為基隆 河中游的起點，流向轉西南，經五堵、汐止至南港，
因河流坡降較緩，向兩側侵蝕較下切作用為強，進 入河谷平原後多處成為低位河階。最後，於北投附 近匯入淡水河，主流約 86 公里，流域面積 501 平 方公里。林朝棨(1957)認為原本古基隆河由東北方 深澳灣出海，因臺北盆地下陷，河蝕作用轉盛，被 西側河流襲奪轉而匯入臺北盆地。基隆河其主要支 流包括東勢坑溪、拔西猴溪、瑪陵坑溪、鹿寮溪、
北港溪、大坑溪及四分溪等大小支流，集水區整體 而言呈現北緩南陡的地勢。
綜合前人在河流縱剖面、水系分析以及河階對
比之研究成果(早坂一郎，1930；林朝棨，1957；
Hsu, 1974；杜友仁，1997；吳麗娟，2000)，認為 基隆河流域的發育與台北盆地的演化息息相關。由 於北臺灣在更新世以來的構造活動，從壓縮型應力 區轉換為伸張型應力區，以致造成台北盆地的陷 落 。在河流地形演育部分，因為侵蝕基準面的改 變，加劇了盆地周圍河流的向源侵蝕，並且改變原 有的河流的流路，反應在基隆河流域的結果是襲奪 了早期向北出海的古基隆河。
圖 1 基隆河流域及其支流水系分佈圖
圖 2 基隆河流域及其支流之地質圖(套疊地調所台北及雙溪圖幅)
深澳港
就地質岩性而言，基隆河流域除東北方地區有 基隆火山群，西北方地區有大屯火山群外，大部分 地區位於臺灣北部第三系沈積岩區，其下游地區已 進入臺北盆地第四系地層中(圖 2)。基隆河的水系 發育主要受控於三個主要的地質構造：石底向斜、
侯硐背斜以及八堵向斜。基隆河流域上游流經石底 向斜的軸部，由侯硐背斜構成本流域中央的嶺線，
在三貂嶺至瑞芳段呈現圓弧狀的流路發育，則是受 到侯硐背斜傾没的影響。基隆河流域中游呈現蜿蜒 的河曲，大致沿八堵向斜的軸部發育。
基隆河流域涵蓋北臺灣褶皺衝斷帶及雪山山 脈帶內之數條逆斷層，各方面之資料指出現今台灣 北部之應力場已由壓縮轉為伸張性，並由山腳斷層 之活動及地表特徵相互印證。Chan (2007)利用地表 LiDAR 資料之分析，顯示基隆河流域內較老之逆 斷層，基本上並無明顯之近期內的正斷層活動跡象，
其結果和 GPS 大地變形測量之結果不謀而合(Yu et al., 1997)。綜合地形及地質相關資料顯示，本區內 之地質構造特徵，其活動年代可能較久遠，並無保 留近期活動之地形特徵。
3. 基隆河流域水系及階地 分析
河流地形常受到不同的岩石特性、地質構造所 影響，透過研究河流水系的幾何分佈與剖面型態特 徵，可瞭解河流地形受到相關地質作用的影響。近 來年，利用數值地形推衍各種地表地形參數，如：
坡度、坡向、曲率、剖面、集水區、水系等，可用 於水文分析與模擬、土壤侵蝕研究、以及其他生態 環境之模擬與分析等。數值地形模型常用者有三類：
包括規則網格(regular grid)、數值等高線(digital contour)、以及不規則三角網(Triangulated Irregular Network, TIN)等(Mark, 1978; Burrough, 1986)。利 用網格式數值地形於水系分析，在運算上是最有效
率的方法(Moore et al.,1991)，因此網格式數值地形 幾乎是數值地形的代名詞(Theobald, 1989)。
水系相關圖層被廣用於輔助水文分析、地形分 析、地質岩性分析、以及地質構造分析等(Bertoldi et al., 2008; Snyder, 2009)，本研究利用 2 公尺網格 式 LiDAR DTM 進行水系分析，利用地表逕流隨著 地形起伏的原理，在萃取坡度、坡向、水系網路及 集水區特性，得到流域內完整的水系網路以及流域 集水區之地形型態與特徵(圖 3)。利用 ArcGIS 地理 資訊系統中地形分析的模組，產製相關之圖層包括:
坡度圖、坡向圖、集水區邊界圖、與水系網路圖等 四種網格圖層以及集水區邊界圖與水系網路圖等 兩種向量圖層，進行萃取河川網系與地形特徵，以 獲得研究區域內之水系分析圖、集水區分佈圖以及 坡度、坡向之空間分析圖等，其目的在於提供 ArcGIS 地理資訊系統中進行各圖層套疊與資料分 析之用。
坡度分析反應出的是地形起伏的劇烈程度，本 區地形起伏平緩處，除了廣布於臺北盆地內，位於 八堵到南港之間也存在著大規模的平緩地形，反應 基隆河中游地區沿岸階地的分布情形。傳統階地判 釋多半藉由航照立體對在立體鏡中的成像，加上人 為的判釋以區別不同的階地。近年來進行河流階地 分析時，利用高解析度數值地面模型資料，配合高 解析度的航空照片，可以清楚地繪製河道與階地分 布，精確完整地呈現原始階地地形面的特性，以及 其兩側沖積平原之範圍(Lane, 2003)。本研究在劃 定河流階地的分佈時，將坡度分析圖中坡度低於五 度視為同一原始地形面，並配合五千分之一高精度 之彩色航空照片，加強地表特徵物的判釋，以區別 該地形面是否為人工開發平坦面或天然的河階地 形，最後藉由地形剖面圖中之河流比高加以區分出 不同時期的階地(圖 4)。
圖 3 基隆河流域及其支流水系分佈圖，底圖為 2 公尺的 LIDAR DTM
圖 4 基隆河流域及其支流之坡度分析圖，紅色區域清楚顯示河流沿線之河階分布
4. 結果與討論
4.1 基隆河流域水系與階地分 析
綜合水系分析的結果，說明了基隆河的發育受 到構造條件以及地形坡度的影響，主要的支流發育 在中游地區較為開闊的八堵向斜軸部附近，位於北
側的集水區也較南側平緩，對稱性的水系代表了近 期並無重大的構造活動造成劇烈的地殼變動。藉由 基隆河流域水系分析的結果，可將基隆河大致分成 三個區段：上游地區菁桐坑－暖暖段，其主流流經 石底向斜在侯硐與瑞芳地區轉向，其支流短小且與 主流正交。中游地區八堵－南港段，其主流發育在 較為開闊的八堵向斜兩側，河谷較為廣寬，曲流發 達、並分佈階地。下游地區南港－關渡段，進入台
北盆地後呈現樹技狀水系。從河流平面的型態，可 以發現基隆河流域之上、中游，其支流方向均大致 與主流方向呈正交，暗示著河流型態與區域構造之 間的關係。
階地分析結果顯示基隆河流域的河階可分為 三群，第一群從南港至八堵之間，存在大規模的河 階，分布於河道兩側，可以視為自由曲流之間的產 物。第二群分布自暖暖到侯硐之間，階地型態以狹 長型為主，其中等高的階地面延伸約 2 公里，我們 發現在暖暖地區存在順向坡地形，經檢視本區為發 生大型順向坡滑動之殘餘地形，而此處之河流剖面 亦顯示出遷急點，推測可能與古山崩崩積物堰塞河 道有關，後續將於 4.4 節進行討論。
第三群分布在上游地區，位於十分寮至菁桐坑 一帶，存在著連續且面積小的階地。中下游階地群 普遍可見二至三階低位河階的分布，而高位河階分 布零星，僅於四腳亭與瑞濱附近可見兩處階地，此 兩階地之高度分別為 95 及 110 公尺，然向下游方 向延伸情況不佳。上游階地群可見兩階河階的分布，
最低的一階在嶺腳瀑布的上游未見其蹤。而高於此 兩河階之上，可見有河階礫石零星分布。
4.2 基隆河流域縱剖面分析
以河道縱剖面來看，基隆河呈階狀分段式均夷 剖面，剖面上明顯轉折點代表新、舊河道剖面的交 會點(稱遷急點 knick point)。其中，十分瀑布為基 隆河上游最明顯的遷急點，透過河階面的對比，推 測為台北盆地變遷引發河流回春作用所造成的遷 急點。十分瀑布遷急點以下，由於河流侵蝕力加大，
加上岩性以厚層砂岩為主，抗蝕力高，因而在大華 及三貂嶺有峽谷與壺穴地形出現。就十分瀑布來說，
後退的機制因瀑壁節理發達，岩石隨著節理不斷崩 落，整體呈平行後退。
本研究利用 2 公尺 DTM 繪製進行基隆河主河 道之縱剖面圖，藉由縱剖面上的河道高度變化(遷 急點)反應出河流侵蝕水準改變後新舊剖面的交界，
結合區域地質圖將主要河床之岩性標示於河流縱 剖面。由河流縱剖面圖(圖 6)中清楚顯示基隆河從 上游至下游依序出現嶺腳瀑布遷急點、十分瀑布遷 急點、侯硐遷急點以及碇內遷急點。尤其是位於南 港層中的十分瀑布遷急點高差將近 40 公尺，明顯 與其下游河段之縱剖面不同。整體來說，這些分佈 於基隆河的遷急點似乎沒有反應在不同的岩性交 接處，亦沒有出現在主要的斷層帶上，說明了這些 遷急點的成因與其河流侵蝕基準面的變化有關，至 於造成侵蝕基準面的改變，可能與近期台北盆地張 裂活動有關。
圖 5 基隆河流域水系及河流階地分佈圖
圖 6 基隆河流域河流縱剖面與其主要遷急點分佈位置圖
4.3 基隆河流域次集水區盆地 分析
如同基隆河階地分析結果，將基隆河流域的河 階分為三群，第一區分佈在南港至八堵之間，受控 於八堵向斜及台北盆地演育的影響，東北岸、西南 岸基本上呈現對稱分布的形式。依據不同的河流級 數亦可區分基隆河流域次集水區盆地，其中河流級 數越高其集水區面積越大，此區兩側之各次集水區 盆地發展大致均等(圖 7)。
第二區分布自暖暖到侯硐之間，本區內各次集 水區之面積、形態差異較大，較無法以單一構造因 素解釋，在基隆暖暖附近之順向坡地形可能導致古 山崩堰塞河道及影響次集水區發育。第三區位於十 分寮至菁桐坑之間，和第一區段之次集水區形貌類 似，受控於石底向斜影響，東北岸、西南岸基本上 亦呈對稱分布。同岸之各次集水區盆地發展大致均 等，唯東北岸位於較平緩之石底向斜東北翼，故次 集水區長度基本上較西南岸長。
4.4 基隆河流域暖暖－侯硐段 之河道堰塞研究
綜合本研究在基隆河流域之河流地形分析，發 現自暖暖到侯硐之間河階分佈、集水區之面積等形 態特徵差異較大，我們進一步針對暖暖－侯硐段進 行塊體滑動地形特徵與河流階地地形之相關研究，
在實際進行野外觀察與調查後推論河道堰塞形成 的機制。
本研究利用 2 公尺 DTM 繪製基隆河流域中游 一帶的河道與階地分佈圖(圖 8(a))，其主要河道兩 側分佈著大小不一的河流階地，特別是在暖暖地區 原有寬廣的河谷，突然成為狹窄的河谷，在地形特 徵上清楚顯示塊體滑移的特徵。基隆河進入八堵之 後，主河道兩側的階地又開始變得成寬廣。利用不 同高度的 DTM 呈現不同的階地，圖 8(b)可以看到 2D 的階地分佈平面圖，配合 3D 河流縱剖面(圖 8(c)) 可以看出河床與階地之相對位置與高度，值得注意 的是古山崩位置其河流比高可達 40 公尺，而上游
地區平均階地比高約為 25 公尺，有別於正常河流 堆積形態(謝孟龍，2007)，因此判斷古山崩形成之 時，可視為一個自然的堤壩，並且造成暫時侵蝕基 準面的改變。
暖暖地區的山坡型態為順向坡地形(圖 9)，順 向坡是指坡地之傾斜與地層之傾斜同一方向，此一 古山崩地形，可能因為受到基隆河侵蝕，坡腳遭切 除致失去下方支撐力，或雨水下滲至地層面上造成 潤滑作用，使上方岩層沿層面下滑產生順向坡滑 動 ，遺留平面狀地形。順向坡所造成的岩層崩落
主要受控地層傾向與傾角，在傾斜、互層狀沉積岩 區此為相當常見之地質災害。
通常順向坡之坡面是由砂岩之組成，坡面有分 布著深層之表土或崩積土等，這種情形最常出現在 本區木山層、大寮層、石底層及南港層等地層出露 之區域，因為前述兩岩層都是由砂頁岩所組成，而 且砂岩層厚且堅硬，加上地層傾角在 20°至 80°之 間，常形成明顯之順向坡地形。基隆河北岸都有大 面積之順向坡區域，坡面傾向東南，坡面傾向變化 很小，使得地形等高線幾乎呈直線狀。
圖 7 基隆河流域及其支流之集水區分析圖
圖 8 (a)基隆河流域暖暖-侯硐地區之河階分佈圖，紅色圓圈地區為推測之古山崩發生區，(b)階地平面分佈 圖可見主河道兩側之河流階地，不同的顏色反應出相異的階地高差，(c)三維型態呈現之基隆河中段 縱剖面圖與階地分佈，垂直放大 30 倍
a
b
)c
)圖 9 基隆河流域暖暖地區古山崩滑動之殘餘地形
圖 10 基隆河流域暖暖地區之大型塊體滑動之殘餘地形及可能堰塞河道的地點，垂直放大 1.5 倍
本研究利用地形特徵分析上的成果，進行野外 考察。經碇內地區野外勘查結果，在暖暖附近之基 隆河道形成一隘口，隘口上游處立即有一老河階地 分佈(圖 11)；經過實際野外調查與資料分析後，初 步 判 定 本 區 之 可 能 坡 面 破 壞 形 式 以 岩 石 傾 翻 (toppling)或楔型岩石崩落(wedge failure)形成，其
坡面穩定分析如圖 12 所示。從坡面與層面資料顯 示其滑動區為順向坡地形，二組節理面所造成的不 穩定區域則反應了其主控崩塌面的發育。綜合暖暖
－侯硐地區所做的相關研究，提供了基隆河流域中 游之大型河階生成的機制。
圖 11 基隆河暖暖附近河谷照片。圖中兩處可能之古地滑區及堰塞形成之河階地
圖 12 坡面穩定分析結果。圖中層面位態為 N20E 傾角為 10 度及三組弱面(含坡面)
5. 結論
本研究利用高精度空載光達資料進行基隆河 河流地形分析，包括集水區盆地、水系分布、河流 發育、地形特徵等分析。綜合各項結果顯示，基隆 河的發育受到區域構造的影響，主要的支流發育在 中游地區的八堵向斜軸部，位於北側的集水區也較
南側平緩，對稱性的水系及次集水區盆地，暗示本 區近期並無重大的構造活動造成劇烈的地殼變動。
利用高精度空載光達數值地形資料，去除地表植被 與建物之餘，有助於呈現原始地形面，增加地形分 析與階地對比之可靠性。另本研究發現基隆河第二 群河階，自暖暖到侯硐之間，各次集水區之面積、
形態差異較大，較無法以單一構造因素解釋。利用 地形特徵分析上的成果，經過實際野外調查與構造
分析後，在暖暖－侯硐地區發現一處順向坡地形可 能導致古山崩堰塞河道及影響次集水區發育，初步 判 定 本 區 之 可 能 坡 面 破 壞 形 式 以 岩 石 傾 翻 (toppling)或楔型岩石崩落(wedge failure)形成，古 山崩的形成造成了一個暫時侵蝕基準，提供了基隆 河流域中游之大型河階生成的機制。
參考文獻
早坂一郎，1930。基隆川溪谷的研究(日文)，臺灣 地學記事，1：60-637。
杜友仁，1997。基隆河之地形研究，國立中央大學 應用地質研究所碩士論文。
吳麗娟，2000。臺灣北部主要河川遷急點之地形學 研究，中國文化大學地學研究所，共 120 頁。
林朝棨，1957。台灣地形，台灣省文獻委員會出版，
共 424 頁。
周淑文，鄧屬予，1998。基隆河襲奪之探討，地質，
18(2)：1-16。
劉明錡，2004。臺灣西北部河階之地形學研究，國 立師範大學地理研究所博士論文。
謝孟龍，2007。臺灣河階地形研究的回顧，檢討與 展望。經濟部中央地質調查所特刊，18：
209-242。
Bertoldi, W., Gurnell, A.M., and Drake, N.A., 2011.
The topographic signature of vegetation development along a braided river: Results of a combined analysis of airborne lidar, color air photographs, and ground measurements. Water Resour. Res., 47, W06525, doi:10.1029/2010WR010319.
Burrough, P.A., 1986. Principles of Geographical Information Systems for Land Resources Assessment. Oxford : Clarendon Press.
Chan, Y.C., 2007. “Acquisition of High-Resolution Digital Elevation Model Using Airborne LIDAR Technique and Its Application to Active Geomorphology in the Taipei Metropolitan Area.(3/3)” Project Final Report, Central Geological Survey in Taiwan. pp. 150.
Hsu, T.L., 1974. Fluvial landforms of northern Taiwan and their neotectonic significance: Bull.
Geol. Surv. Taiwan, 24:109-118.
Lane, S.N., 2000. The Measurement of River
Channel Morphology Using Digital Photogrammetry. The Photogrammetric Record, 16: 937–961. doi: 10.1111/0031-868X.00159 Mark, D.M., 1978. Concepts of Data Structure for
Digital Terrain Model. Proceedings of DTM Symposium, ASP-ASCM, 24-31.
Moore, I.D., Grayson, R.B., and Landson, A.R., 1991.
Digital Terrain Modeling: a review of hydrological, geomorphological, and biological applications. Hydrological Process, 5:3-30.
Theobald, D.M., 1989. Accuracy and bias issues in surface representation. In: Goodchild, M.F., Gopal, S. (Eds.), The Accuracy of Spatial Database, Taylor and Francis, New York, 99–
106.
Snyder, N.P., 2009. Studying Stream Morphology With Airborne Laser Elevation Data. Eos Trans.
AGU, 90(6):45-46, doi:10.1029/2009EO060001.
Yu, S.B., Chen, H.Y., and Kuo, L.C., 1997. Velocity of field of GPS stations in the Taiwan area.
Tectonphysics, 1-3:41-59.
Using airborne LiDAR derived DEMs to analyze river morphology in Keelung stream
Rou-Fei Chen1* Yu-Chang Chan2 Kuo-Jen Chang3 Yu-Chung Hsieh4
ABSTRACT
We organized and performed airborne LIDAR mapping in the Metropolitan Taipei area in order to produce high-resolution and high-precision digital elevation models for geological research. Applying the newly acquired LIDAR DEMs from the MOEA Central Geological Survey, we analyzed geologic and geomorphic features of the Keelung River drainage area. The studied features include topographic scarps and lineaments, river terraces, drainage basins, and landslide scarps. Based on our analysis of the acquired LIDAR DEMs, we arrived at the following findings and conclusions: 1) The river basin analysis based on the LIDAR DTM indicates that the development of the Keelung River was influenced by regional structural patterns and topographic slopes.
Tectonic influences appear to be relatively small for the development of the Keelung River. 2) Along the Keelung River, the LIDAR DTM indicates anomalous terrace morphology and river drainage at the middle section of the river. We interpreted this area to be caused by paleo-landslides, which form a dammed lake in the river and subsequently developed the main terrace of the same level for large area.
Key Words: airborne LIDAR, digital elevation model, river geomorphology, Paleo-landslide
Volume 21, No.1, 2016, pp. 13-30 DOI：10.6574/JPRS.2016.21(1).2
1國立中央大學太空及遙測研究中心 專任助理 收到日期：民國 103 年 11 月 10 日
2國立中央大學太空及遙測研究中心 教授 修改日期：民國 105 年 07 月 01 日
3國立中央大學土木工程學系空間資訊組 博士候選人 接受日期：民國 105 年 07 月 14 日
*通訊作者, 電話：+886-3-4227151 ext. 57619, E-mail：[email protected]
二次微分法於空載全波形光達之高斯波形擬合與 地物分類
盧佑樺
1蔡富安
2*賴哲儇
3摘 要
全波形光達(Full Waveform LiDAR)可完整記錄每條雷射光束在不同時間回傳的反射能量，若能藉 由這些記錄下來的波形衍生出額外資訊，可望增進地形三維重建及地物分類的成效。波形擬合(fitting) 與特徵萃取(feature extraction)是處理和分析全波形光達資料的重要過程。本研究利用高斯函數擬合空 載全波形光達資料之波形，並以二次微分法搜尋迭代計算之初始值。基於波形擬合成果，可以衍生出 振 幅 (amplitude) 、波寬(width) 、背向散射截面 (backscatter cross-section)等波形參數 ，且配合強度 (intensity) 、正規化高程(normalized height)等傳統光達特徵與同步搭載之中像幅正射影像，作為地物分 類之依據。另外，本研究亦比較簡易貝氏(naïve Bayesian)與隨機森林(random forests)等兩種分類器之成 效。研究結果顯示，運用二次微分法決定高斯函數迭代計算之初始值，並搭配隨機森林分類器，能提 供較佳的擬合及分類成果，且波形參數有助於植物類別的辨識。
關鍵字：全波形光達、波形擬合、地物分類、二次微分法、隨機森林分類器
1. 背景
空載光達(Airborne LiDAR)或空載雷射掃描 儀(Airborne Laser Scanner, ALS)為一主動式遙測技 術。利用飛行載具搭載的雷射掃描儀測得載具與地 表之距離及發射角度，配合全球導航衛星系統 (Global Navigation Satellite System, GNSS)及慣性 導航系統(Inertial Navigation System, INS)所提供之 飛 行 載 具 的 外 方 位 參 數 (Exterior Orientation Parameter, EOP) ， 再 以 直 接 地 理 對 位 法 (direct geo-referencing)求得高精度之地表坐標，且高程精 度可達公寸級(Liu, 2008)。由於資料記錄大量三維 地表坐標及其反射強度(intensity)，視覺展示有如 雲霧般，故稱為點雲(point clouds)。常見的空載光 達應用包含製作數值高程模型(Digital Elevation
Model, DEM) (Podobnikar and Vrečko, 2012)、數值 地表模型(Digital Surface Model, DSM) (Stal et al., 2013)、數位城市模型(Sohn and Dowman, 2007)等 等。
然而，過往因硬體儲存設備所限，光達資料 獲取時不得不將資料簡化以節省硬碟空間。現今科 技不斷演進，儲存設備能夠完全保存原始光達點雲 紀錄，使全波形光達資料與技術得以嶄露頭角。相 較於過去傳統光達僅記載地表離散點之三維坐標 及反射強度，全波形光達可額外提供波形資訊，其 為雷射光束所照射路徑上所有回傳能量，其取樣頻 率可達 1Ghz，如圖 1 所示。波形資訊提供了兩個 重要的額外貢獻(Bretar et al., 2008)：其一，波形中 所有的取樣資料皆可反推至實際地表坐標，以增加 點 雲 密 度 ， 尤 其 是 像 森 林 等 具 強 烈 多 重 回 波
(multiple return)特性的區域；其二，藉由波形的模 擬與分析，可進一步反演雷射光束所涵蓋三度空間 範圍內各目標物件的物理特性。因此，全波形光達 可同時提供目標區域豐富的幾何和反射波形資訊，
有很大潛力進行更複雜及先進的應用。
2. 研究動機與目的
全波形光達是空間資訊領域較新興的課題之 一，目前研究方向可概略分成三個主題：一為自行 定義點萃取法，反推波形所有取樣資料的地表坐 標 ，提高點雲密度與三維精度(Chauve et al., 2009;
Qin et al., 2012)；二為森林相關研究，例如產製樹 冠高度模型(Canopy Height Model, CHM)、推估林 木材積以及葉面積指數(Leaf Area Index, LAI)等等 (Reitberger et al., 2009; Yao et al., 2012)，協助森林 資源調查；三為從波形資訊中獲取額外特徵，改善 地物分類(Fieber et al., 2013; Mallet et al., 2011)。後 者常見作法是以高斯函式擬合波形，萃取其波形參 數。然而，高斯函數為非線性函數，必須利用非線
性最小二乘法迭代求解，其初始值會影響波形擬合 成果(Lin et al., 2010)。
因此，本研究嘗試比較兩種初始值偵測法對高 斯函數擬合波形之差異，一為記錄在光達波形資料 檔頭，以儀器自定義方式偵測之點雲於該波形上的 位置(後文以點位稱之，其通常為波峰處)；二為對 波形資料二次微分萃取而得。接續利用高斯波形衍 生出振幅(amplitude)、波寬(width)、背向散射截面 (backscatter cross-section)等波形參數(或稱全波形 特 徵 ) ， 並 配 合 強 度 (intensity) 、 正 規 化 高 程 (normalized height)等傳統光達特徵與同步搭載之 中像幅正射影像，作為地物分類之依據。有關分類 器選擇，本研究採用隨機森林(random forests)演算 法，因其可快速的處理大量資料(Guo et al., 2011)，
並與資料探勘中常用的簡易貝氏(naïve Bayesian) 比較。除此之外，亦探討全波形光達特徵對於整體 分類成果與類別分離之影響。
圖 1 不同目標物反射能量之差別，灰底波形為傳統光達紀錄，黃底波形為全波形光達紀錄(Ullrich and Reichert, 2005)
3. 研究區域與資料
本研究主要實驗資料為利用 Optech ALTM Pegasus 空載全波形光達儀器所獲取之空載全波形 光達資料。表 1 為儀器參數，資料獲取時間為 2012 年 5 月 7 日，地理位置在桃園縣中壢市，飛行高度 約為 2200 m，光跡(footprint)大小約為 0.44 m，總 光達點數有 10,639,828 點，光達點密度約為 0.53 點/m2。另外，載具同時搭載中像幅相機拍攝。圖 2 為研究區域之 DEM、點雲強度影像及中像幅正 射影像。
本研究測試三種案例，各別考慮不同的地物類 別。案例一為評估高斯函數擬合波形成果，以及全 波形光達資料於地物初步分類的可行性。案例二包 含更多地物類別，以測試利用全波形光達資料執行 進階分類任務的能力。案例三僅考慮植物類別，藉 此了解全波形光達對於不同植物判識之潛力。通常 以多光譜遙測影像進行不同植物類別分類時，因頻
譜資訊有限，必須仰賴紋理資訊(例如 Tsai and Chou, 2006)或高光譜資料(例如 Tsai et al., 2007a)，
甚至高光譜的紋理資訊(例如 Tsai et al., 2007b; Tsai and Lai, 2013)，才有機會分離不同植物類別；而案 例三即為測試全波形光達資料是否有辦法達到類 似成果。
表 1 Optech ALTM Pegasus 儀器參數
項目 設置
掃描方式 旋轉鏡
航高(m) 300~2500 波長(nm) 1064 脈衝頻率(kHz) 100~400
發散角(mrad) 0.2 掃描角(deg) 0~60
取樣間隔(ns) 1
資料來源：http://www.chsurvey.com.tw/page04.html
(a) DEM (b) 點雲強度影像 (c) 中像幅正射影像 圖 2 研究區域與使用資料
3.1 案例一
此案例考慮類別皆位於國立中央大學校 區內，地物包含建物、樹木、草地與地面等四 種類別，地真資料的點雲數量如表 2，其與正 射影像套合結果如圖 3。此資料主要是評估高 斯函數擬合波形的成果，並比較儀器提供點位 與二次微分法提供之初始值對擬合成果之影 響，並對這些分離度較高的類別分類，初步測 試全波形光達資料於地物分類之可行性。
表 2 案例一地真資料之光達點數
類別 光達點數
建物 3,865
草地 602
地面 3,004
樹木 3,079
總計 10,550
圖 3 案例一地真資料與正射影像套合結果
3.2 案例二
此案例擴大案例一範圍，並考慮更多地物 類別，總共 7 種類別。地真資料之點雲數量如 表 3，圖 4 為地真資料與正射影像套合結果。
此案例包含更多的地物類別，為測試全波形光 達資料執行進階分類的能力。
表 3 案例二地真資料之光達點數
類別 光達點數
灌木 935
喬木 2,337
草地 4,541
建物 5,362
土壤 2,243
柏油 1,963
PU 1,374
總計 18,755
圖 4 案例二地真資料與正射影像套合結果
3.3 案例三
本案例類別全為植物，計有針葉樹、闊葉 樹、稻作、竹林、果樹及旱作等六類，以了解 全波形光達於植物分類之潛力。地真資料根據 2008 年產製的土地利用圖選取，但與全波形 光達資料獲取時間差了 4 年，因此地真資料已 刻意避開土地利用變遷處。地真資料的點雲數 量及套合正射影像結果如表 4 與圖 5 所示。
表 4 案例三地真資料之光達點數
類別 光達點數
針葉樹 2,144
闊葉樹 4,664
稻作 1,899
竹林 1,425
果樹 1,510
旱作 1,484
總計 13,126
圖 5 案例三地真資料與正射影像套合結果
4. 研究流程與方法
本研究主要包含波形擬合及地物分類 兩個部分，流程及方法詳述如下。
4.1 波形分析
本項工作為比較儀器提供之點位與二次 微分法作為非線性最小二乘平差初始值之高 斯擬合成果。研究流程包含雜訊濾除、二次微 分法偵測初始值、高斯擬合函數計算、檢驗分 析，如圖6所示。各方法理論請見於後續小節。
4.1.1 雜訊濾除
光達感測器接收反射能量時，難免會接 收到雜訊，必須藉由雜訊濾除方法降低影響，
以利後續二次微分萃取初始值及波形擬合等 任務。本研究雜訊濾除過程先消除背景雜訊 (background noise) ， 再 以 均 值 濾 波 器 (mean filter)降低剩餘干擾。
理想而言，沒有接收到反射能量時應無強 度訊號，但受到背景輻射影響，出現微量訊號。
本研究使用單一門檻值(threshold)濾除背景雜 訊，不同光束的門檻值亦不同。而均值濾波器 需給定罩窗尺寸(block size)，此將嚴重影響二 次微分萃取的初始值位置和波形擬合成果。本 研究以 3 ns 及 5 ns 的罩窗尺寸測試所有波形 資 料 ， 並 選 用 擬 合 均 方 根 誤 差 (Root Mean Square Error, RMSE)最小之罩窗尺寸。
圖6 波形分析之流程圖
4.1.2 二次微分法
由於高斯擬合函數為非線性函數，需給定求解 參數之初始值，以進行迭代運算得到最佳解。高斯 函數參數有三：振幅、波形中心與波寬，本研究採 用二次微分法尋找初始值，如式(1)。依波形二次 微分小於零之局部最小值(local minimum)視為波 形中心初始位置，後依該位置記錄之能量及其發射 波形之波寬作為振幅與波寬之初始值。圖 7 為二次
微分尋找波峰之範例，其中圖 7(c)擬合之成果較圖 7(a)儀器能偵測出能量微弱的波峰，顯示波形二次 微分之局部最小值對找出波形波峰有相當良好的 成果(Tsai and Philpot, 1998)。
𝑑𝑦
𝑑𝑡|𝑖 ≈𝑦(𝑡𝑘−∆𝑡)−2𝑦(𝑡𝑘)+𝑦(𝑡𝑘+∆𝑡)
∆𝑡2 ………(1)
其中，𝑦表示波形資料，∆𝑡為時間間距。
(a) 波形與儀器提供之波峰(紅色線條)
(b) 波形二次微分與局部最小值(紅色線條)
(c) 二次微分偵測之波峰與高斯擬合成果 圖 7 儀器提供與二次微分偵測之波峰差異
(c)
4.1.3 高斯擬合
光達波形經雜訊濾除後，其形狀通常可以透過 多組簡單的高斯函數組合加以近似，故高斯函數廣 泛使用於光達波形擬合研究中(Wagner et al., 2006)。
雖然單一高斯分佈未必能完整描述光達波形，但若 使用多個高斯分佈(即使用多組波峰、振幅及波寬 參數加以組合)，則能非常近似地擬合出複雜的波 形 ， 如 圖 7(c) 所 示 。 高 斯 波 形 主 要 包 含 振 幅 (amplitude)及波寬(width)等兩個參數，前者意指波 峰的反射強度，後者表示正負一個標準差的寬度，
如圖8所示。光達波形隱含地物資訊，平坦地的回 波寬度小，且振幅大；傾斜地回波寬度較大且不易 對稱；若表面粗糙或混合多類別時(如樹木)，則波 寬大、振幅較小，並造成多個回波。先前研究(Mallet et al., 2011)指出，比起廣義高斯模型(generalized Gaussian model)或非對稱函數，如Weibull函數、
Burr函數等，使用對稱的高斯模型進行擬合可得到 不錯的成果，且同時兼顧運算速度，且廣義高斯模 型及非對稱函數之參數在分類上的效用不大。因此，
本研究採用高斯函數擬合波形。
圖8 高斯波形與其參數(Heinzel and Koch, 2011)
高斯擬合是將全波形光達每一回波視為高斯 波，將一條光束上的一個或數個回波拆解並擬合出 高斯參數以利後續運用。高斯擬合據雷達方程式改 寫成光達形式(Wagner et al., 2006)如下：
𝑆(𝑡) = 𝐴𝑒−(𝑡−𝑡𝑖)22𝑤2
………...
(2)其中，
𝑆：接收之波形 𝑡：波形接收之時間 𝑡𝑖：高斯波形中心之時間 𝐴：波形之振幅
w：波形之波寬
由於高斯擬合為非線性函數，本研究使用 Levenberg-Marquardt 最佳化為基礎的非線性最小 二乘法求解。此演算法能藉由執行時修改參數達到 結合高斯-牛頓演算法以及梯度下降法的優點，並 改 善 兩 者 缺 點 (Levenberg, 1944; Marquardt, 1963)。
4.1.4 擬合精度評估
本研究波形擬合之成果檢驗使用均方根誤差 (Root Mean Square Error, RMSE)評估，其為殘差的 平方和除以樣本數後開根號，值愈接近 0 誤差愈 小 ，以此決定高斯函數的擬合成果，如式(3)所示。
並比較擬合失敗率及擬合樣本數(波峰數)，擬合失 敗的定義為(1)參數不合理，即以統計的方式濾除 波寬大於三倍標準差之數據，因波寬表目標物於光 束方向上之垂直分布，過大的波寬並不合理；及(2) 高斯函數中心距初始峰值過遠，其大於半高寬 (Full Wave at Half Maximum, FWHM)。
RMSE = √∑𝒏𝒊=𝟏(𝒚𝒏𝒊−𝒚̂)𝒊𝟐
……….
(3)其中，𝑦𝑖：觀測資料
𝑦̂：函數擬合而得之估計值 𝑖 𝑛：觀測資料數目
4.2 地物分類
本研究地物分類之流程如圖 9 所示。使用的 特徵參數包括傳統光達的強度和正規化高程，全波
形光達特徵的振幅、波寬、背向散射截面，以及正 射影像產生的綠指標值(Greenness Index, GI)。分類 器則考量隨機森林及簡單貝氏演算法。精度評估部 分經由比較各檢核點之分類後的土地覆蓋類別與 地真資料的土地類別，可產生一誤差矩陣(Error Matrix)，計算生產者精度(Producer's Accuracy, PA)、
使用者精度(User's Accuracy, UA)、總體準確度 (Overall Accuracy, OA)及 Kappa 值。各項請見以下 小節。
4.2.1 分類特徵
本研究使用之分類特徵參數除了從光達獲取 的強度、正規化高程，全波形光達的振幅、波寬、
背向散射截面外，尚有從正射影像上獲取的綠指標 值。分類成果將會對傳統光達與全波形光達配合其 他特徵之組合進行比較。其中，強度為儀器所提供，
振幅與波寬是從高斯函數擬合出來的結果所取得。
另外，光達反射的能量取決於地表的反射係數與接 觸的面積，背向散射截面即是描述後者，反映目標 物表面的物理特性；背向散射截面根據光達方程式 推導得式(4)。
𝜎𝑖= 𝐶cal𝑅𝑖4𝐴𝑤
……….
(4) 其中，𝜎𝑖：背向散射截面 𝐶cal：校正常數
𝑅𝑖：儀器與目標物之距離
理論上(Bretar et al., 2008; Briese et al., 2008; Hofle and Pfeifer, 2007; Wagner et al., 2006)，接近底點 (nadir)之背向散射截面會如式(5)，其中𝛽𝑡為光束發 散寬度，ρ 為目標物之反射率。
𝜎𝑖= 𝜋𝜌R2i𝛽𝑡2
………...
(5)整合式(4)與(5)可得式(6)，使用人工選取接近 底點(nadir)之柏油路求出 Ccal (假設柏油路的反射
率為 0.25)再將其代回式(4)，計算出所有回波之背 向散射截面(Alexander et al., 2010)。原因為所有地 物的反射率(ρ)不易取得，而柏油路是都市區常見 的地物，故藉由此經驗法則求得 Ccal。
𝐶cal=𝑅𝜌𝜋𝛽𝑡2
𝑖2𝐴𝑤… … … …(6)
正規化高程為光達點至數值高程模型之垂直 距離，在此數值高程模型是用同一組光達地面點內 插而得，其空間解析度為 1 公尺，將其雙線性內 插後得到點雲正下方之地表高，後與光達點雲之 Z 值相減得到正規化高程。
在空載光達獲取資料時，一般都會搭配一部中 像幅相機，加入影像的資訊進行分類與光達資料互 補得到較佳的分類成果。因本研究之全波形光達所 搭配之航照並未配有近紅外光波段，無法計算 NDVI (Normalized Difference Vegetation Index)，故 利用紅光波段與綠光波段計算綠指標。綠指標有助 於 分 辨 植 物 與 非 植 物 ， 再 以 最 鄰 近 法 (nearest neighbor)抓取點雲之綠指標值。綠指標如下：
綠指標 =𝐺𝑟𝑒𝑒𝑛−𝑅𝑒𝑑𝐺𝑟𝑒𝑒𝑛+𝑅𝑒𝑑………...(7)
其中，
Green：綠光波段 Red：紅光波段
4.2.2 隨機森林分類法
隨機森林分類法(Breiman, 2001)是基於決策 樹(decision tree)分類的整合式分類法，此法並不需 要假設任何機率分布，通常是使用經驗法則找出最 佳的決策樹，是資料探勘中常用的方法之一。此分 類方法的優點為可提供分類特徵的重要性，並使用 分類成果最佳的方針建立模型。隨機森林中的每棵 樹都是根據隨機向量所建立，而隨機向量是依據固 定機率分配所產生。其流程如圖 10 所示，詳細演 算步驟如下：
(1) 決定訓練區資料 n 組，分類特徵參數 m 種。
(2) 於 m 個輸入分類特徵參數中，選出 T 個子集合，
T<m。
(3) 運用 bootstrap 抽樣，從訓練資料重複抽樣形成 InBag。
(4) 對每一個 InBag，隨機選擇 m 個分類特徵參數，
根據這些分類特徵參數，計算決策樹，最後進 行投票(vote)選出最佳的成果。
(5) 每棵決策樹都不需要進行修剪、完全生長。
隨機向量的決定方法為：每個節點隨機選取 F 輸入特徵進行分割，所要分割的節點是由選取的 F 特徵中決定，樹木持續生長不需修改。如果原始特 徵個數太少，很難隨機選出互相獨立的輸入特徵，
所以就要增加特徵空間，從這些隨機結合的新特徵 F 中選出最適合分割的節點。在每個節點上隨機從 F 個最好的分割點中選取一個以產生隨機樹，演算 法必須檢查所有節點上的分割特徵。
圖 9 地物分類之流程圖
圖 10 隨機森林法之流程圖(Guo et al., 2011)
5. 研究成果
5.1 擬合成果
本研究使用 Matlab R2010b 的 cftool 配合自定 義的初始值萃取程式進行波形擬合。波形擬合部分 使用案例一分析，比較 RMSE、擬合失敗率及擬合 樣本數(波峰數)。圖 11 所示為擬合失敗之案例，
通常記錄的波形前後會有一段背景雜訊，本研究取 波形前後各五個樣本之最大值為該波形之背景雜 訊門檻，但有極少部分波形因截取波形時，最前頭 有不正常的能量記錄，造成背景雜訊門檻過大而使 後方波形資訊遺失，使擬合時參數不合理。
(a)原始波形與儀器提供之波峰
(b)雜訊消除之波形與儀器提供之波峰
(c) 二次微分與波峰偵測結果
圖 11 以儀器提供之波峰作為初始值擬合之失敗案 例一
以儀器提供之初始值擬合有時會發生如圖 12 之情況，圖 12(b)左邊紅線處於濾除雜訊後並無明 顯波峰，造成最後擬合之參數不合理或高斯函數中 心距初始值過遠。推測可能原因是儀器偵測點位時 未濾除高能雜訊，誤判雜訊造成之微弱波形起伏為 波峰。不過二次微分法則可偵測出兩個局部極小值 (圖 12(c))，提供較接近實際波形情況之波峰位置。
(a) 原始波形與儀器提供之波峰
(b) 雜訊消除之波形與儀器提供之波峰
(c) 二次微分與波峰偵測結果
圖 12 以儀器提供之波峰作為初始值擬合之失敗 案例二
表 5 為案例一以儀器提供波峰與二次微分法 於高斯函數擬合波形之成果比較。如表所示，以儀 器所提供之波峰為初始值，其 RMSE 以及擬合失 敗率都高於二次微分法，而二次微分法偵測出的波 峰數較多，推測可能原因有以下三點：(1)本研究 所使用之儀器在一條光束上僅能記錄 4 個光達點，
若高於此數目將使擬合精度以及成功率下降；(2) 疊合與能量較微弱的回波會影響擬合的精度，若不 能有效偵測初始值會使 RMSE 提升，如圖 7；(3) 儀器誤判波峰造成擬合失敗，如圖 12。
表 5 波形擬合成果比較
儀器點位 二次微分法 RMSE(DN) 1.22 0.43 擬合失敗率 5.05% 0.23%
擬合樣本數 10,550 10,627 擬合失敗樣本 533 24
5.2 分類成果分析
本研究使用 Weka (3.6.8)軟體執行地物分類。
Weka 是以 JAVA 語言為基礎的資料探勘工具，使 用者可以自行開發與加值。本研究使用 Weka 軟體 內建分類模組中的隨機森林以及簡單貝氏演算法 進行實驗。
5.2.1 案例一
本案例比較以儀器點位與二次微分法為初始 值所萃取出的特徵對分類精度的影響，僅以全波形 光達特徵振福、波寬加上正規化高程進行分類。表 6 為分類結果，由成果可看出儀器點位與二次微分 法之初始值所得到的特徵在簡易貝氏分類器的表 現相當接近。但在隨機森林分類器上，以二次微分 法為初始值的分類成果較儀器點位所得之成果佳，
OA 及 Kappa 值可達 89.77%與 0.853，顯示二次微 分法能提供較多可利用的樣本並有較佳的分類成 果。因此後續案例將以二次微分法配合高斯函數擬 合波形，並萃取參數供地物分類。
表 6 案例一之分類結果比較表
初始值 簡易貝氏 隨機森林
OA (%) Kappa OA (%) Kappa 儀器
提供 83.31 0.765 85.89 0.798 二次微
分法 83.52 0.765 89.77 0.853
以表 7 之誤差矩陣進一步探討分類的成果，發 現樹木及建物容易混淆(粗體加底線處)。另外，可 發現表 7 與表 2 的光達點雲數略有不同，原因為地 真資料涵蓋不同地物類別的交界，而這些交界的光 達點雲因反射能量混肴，使得 5.1 節的光達點雲偵 測演算法偵測不出來，因此漏掉 25 個光達點雲 (建 物、地面及樹木個別漏掉 5、2 與 18 個光達點雲)。
圖 13 為分類結果圖，圖 14 為分類正確及錯誤之光 達點位，可觀察到分類錯誤的點主要集中在房屋邊 緣，這些點可能是牆面上的窗台、冷氣等，其波形 會較屋頂破碎，易錯分為樹木。圖 14 右下(黃框) 有一棟建築完全錯分為樹木，因這棟建築的屋頂是 鐵皮斜頂，其波寬較其他單一回波長，且鐵皮屋對 紅外光之反射較強，如樹木一般。
表 7 案例一之誤差矩陣(隨機森林+二次微分法)
建物 草地 地面 樹木 Total PA (%) 建物 3300 14 54 492 3860 85.49
草地 4 584 6 8 602 97.01
地面 40 8 2931 23 3002 97.63
樹木 384 13 31 2633 3061 86.02 Total 3728 619 3022 3156 10525
UA (%) 88.52 94.35 96.99 83.43
圖 13 案例一分類成果圖
圖 14 案例一分類正確(白)與錯誤(紅)之光達點
5.2.2 案例二
此案例擴大案例一範圍，且考慮較多種地物，
包含喬木、灌木、草地、土壤、柏油、PU 舖面以 及建物共 7 種類別。為比較全波形與傳統光達的差
異，並同時比較簡易貝氏與隨機森林分類器的分類 成果。表 8 為使用全波形光達特徵振幅、波寬、正 規化高程與背向散射截面以二次微分法為初始值 與隨機森林分類器的誤差矩陣，其 OA 為 93.39%、
Kappa 為 0.919。觀察其分類結果，草地與土壤以
及喬木與建物錯誤分類的點較多(表 8 粗體加底線 處)，皆屬植物與非植物之區別。因此，加入綠指 標後，前述之錯誤分類已大幅降低(表 9)，其 OA 為 98.23%、Kappa 為 0.978。
表 10 與表 11 分別為使用簡易貝氏與隨機森 林分類器進行分類，所使用的特徵如表格內容所 示 。結果可知，若無加入影像特徵(綠指標)，全 波形光達之表現較傳統光達佳，OA 高出約 4~5%，
加入了綠指標後，兩種光達的分類成果均有所提 升 ，但全波形光達的優勢並不明顯，僅較傳統光
達高出 1.5%左右。比較兩種分類器的成果，隨機 森林分類器的精度較簡易貝氏分類器高約 3~10%，
推測可能原因為光達各種特徵的分布不一定相同，
有 些 甚 至 很 難 以 參 數 化 的 機 率 密 度 函 數 (Probability Density Function, PDF)表示，如高程是 光達分類中相當重要的特徵之一，但在建物高程的 分布上(圖 15)，因房屋樓層數不同，其分布可能集 中在幾個特定的高程值附近。因此，以不需假設特 徵分布的隨機森林分類器所得的成果較為優異。
表 8 案例二誤差矩陣(參數：振幅、波寬、正規化高程、背向散射截面，分類器：隨機森林) PU 土壤 柏油 灌木 喬木 建物 草地 Total PA (%) PU 1366 4 0 1 2 0 1 1374 99.42 土壤 0 1885 6 1 0 0 351 2243 84.04 柏油 0 9 1937 10 2 0 5 1963 98.68 灌木 3 6 15 833 40 1 37 935 89.09 喬木 2 0 8 31 2136 159 1 2337 91.40 建物 0 0 0 4 199 5159 0 5362 96.21 草地 3 308 14 9 7 0 4200 4541 92.49 Total 1374 2212 1980 889 2386 5319 4595 18755
UA (%) 99.42 85.22 97.83 93.70 89.52 96.99 91.40
表 9 案例二誤差矩陣(參數：振幅、波寬、正規化高程、背向散射截面、綠指標，分類器：隨機森林) PU 土壤 柏油 灌木 喬木 建物 草地 Total PA (%) PU 1372 0 0 2 0 0 0 1374 99.85 土壤 0 2242 0 0 0 0 1 2243 99.96 柏油 0 0 1948 7 0 0 8 1963 99.24 灌木 2 0 12 840 36 2 43 935 89.84 喬木 2 0 2 29 2226 77 1 2337 95.25 建物 0 0 0 3 70 5289 0 5362 98.64 草地 1 1 12 12 9 0 4506 4541 99.23 Total 1377 2243 1974 893 2341 5368 4559 18755
UA (%) 99.64 99.96 98.68 94.06 95.09 98.53 98.84
表 10 案例二分類成果比較(簡易貝氏分類器)
使用特徵 OA (%) Kappa
傳統光達 強度、正規化高程、綠指標 93.88 0.924
強度、正規化高程 79.18 0.736
全波形光達 振幅、波寬、正規化高程、背向散射截面、綠指標 95.15 0.940 振幅、波寬、正規化高程、背向散射截面 84.61 0.770
表 11 案例二分類成果比較(隨機森林分類器)
使用特徵 OA (%) Kappa
傳統光達 強度、正規化高程、綠指標 97.69 0.971
強度、正規化高程 89.58 0.871
全波形光達 振幅、波寬、正規化高程、背向散射截面、綠指標 98.23 0.978 振幅、波寬、正規化高程、背向散射截面 93.39 0.919
圖 15 建物之正規化高程直方圖
5.2.3 案例三
案例二顯示加入了綠指標之後，全波形光達與 傳統光達的差異極小，因分類錯誤的部份都屬於植 物與非植物類別的混淆，故綠指標能有效地解決此 問題。案例三的主要目的是比較全波形光達與傳統 光達於植物類別上的辨別力，以及簡易貝式分類器 與隨機森林分類器的差異。表 12 為使用簡易貝氏 分類器之分類結果，傳統與全波形光達特徵的表現
不相上下。表 13 則是使用隨機森林分類器之成果，
在傳統光達上的表現與簡易貝氏分類器差異不大，
但全波形光達特徵的表現明顯提升許多，可觀察到 綠指標在全波形光達的顯著性較案例二小；而傳統 光達在加入綠指標後成果雖有所提升，但相較全波 形光達仍有落差。
表 14 為使用全波形光達特徵以及綠指標分類 成果之誤差矩陣。以隨機森林分類器分類，竹林與 闊葉樹相互混淆的情況相當嚴重(表 14 粗體加底線
處)，而其他類別的 PA 與 UA 至少都有 70%。竹林 與闊葉樹的樹冠層皆十分茂密，葉面較寬，反射面 積大，所以在光達波形的表現上相當類似。表 15
為竹林與闊葉樹各種特徵的分布情形，可發現兩種 類別的特徵分布情形相當類似。
表 12 案例三分類成果比較(簡易貝氏分類器)
使用特徵 OA (%) Kappa
傳統光達 強度、正規化高程、綠指標 68.03 0.582
強度、正規化高程 62.44 0.504
全波形光達 振幅、波寬、正規化高程、背向散射截面、綠指標 67.92 0.585 振幅、波寬、正規化高程、背向散射截面 59.99 0.491
表 13 案例三分類成果比較表(隨機森林分類器)
使用特徵 OA (%) Kappa
傳統光達 強度、正規化高程、綠指標 68.03 0.589
強度、正規化高程 57.24 0.452
全波形光達 振幅、波寬、正規化高程、背向散射截面、綠指標 75.86 0.689 振幅、波寬、正規化高程、背向散射截面 72.70 0.648
表 14 案例三之誤差矩陣，使用振幅、波寬、正規化高程、背向散射截面、綠指標(二次微分法+隨機森 林)
竹林 闊葉樹 針葉樹 旱作 果樹 稻作 Total PA (%) 竹林 405 917 69 0 33 1 1425 28.42 闊葉樹 358 3819 276 0 209 2 4664 81.88 針葉樹 30 480 1521 3 108 2 2144 70.94 旱作 0 0 1 1356 17 110 1484 91.37 果樹 26 166 50 22 1155 91 1510 76.49 稻作 0 0 0 170 28 1701 1899 89.57 Total 819 5382 1917 1551 1550 1907 13126
UA (%) 49.45 70.96 79.34 87.43 74.52 89.20
表 15 竹林與闊葉樹各種特徵分布
強度 振幅 波寬 正規化高
程
背向散射 截面
綠指標
竹林 Mean 535.68 22.65 7.61 8.61 0.30 223.51 Std. 300.69 14.53 1.64 2.64 0.20 19.51 闊葉樹 Mean 606.85 25.87 7.23 8.63 0.35 230.67
Std. 325.03 15.63 1.31 3.14 0.22 18.94
6. 結論與建議
本研究探討不同初始值於全波形光達波形擬 合之分析、全波形光達與傳統光達搭配影像分類之 成效，並比較隨機森林與簡單貝式分類器於全波形 光達分類之適用性。研究成果顯示，以二次微分法 偵測初始值在波形擬合與分類精度方兩方面都較 以儀器提供的波峰為初始值有所提升。在波形擬合 的部分，RMSE 僅有 0.43，擬合失敗率為 0.23%；
分類的部分，以二次微分法配合隨機森林分類器分 類的 OA 與 Kappa 分別為 89.77%與 0.853，以儀器 提供之波峰擬合則分別為 85.89%與 0.798，顯示二 次微分法可提供較有辨別力的特徵以及更多的樣 本數供分類。
比較全波形光達特徵與傳統光達特徵，案例二 的都市區域分類成果顯示：若無加入綠指標，全波 形光達(OA: 93.39%; Kappa: 0.919)優於傳統光達 (OA: 89.58%; Kappa: 0.871)；加入綠指標後，兩 者的分類精度都有所提升，而差異僅剩 0.54%，推 測可能原因為，錯誤分類部分均為植物與非植物 (喬木與建物、草地與土壤)，因此綠指標有相當顯 著的效果。案例三的植物分類成果顯示：若無加入 綠指標，全波形光達 OA 為 72.70%，Kappa 值為 0.648，高於傳統光達的 57.24%與 0.452；在加入綠 指標後，全波形光達的 OA 與 Kappa 值分別為 75.86%及 0.689，傳統光達為 68.03%與 0.589。由 此可知全波形光達特徵在具有多重回波特性的植 物類別之分辨能力有較明顯的優勢，若僅要分辨都 市區域類別時，使用一般光達特徵配合影像資訊即 可獲得不錯的成果並兼顧運算時間。
由隨機森林與簡易貝氏分類器的比較可得知，
於三個測試案例均顯示使用隨機森林法可獲得較 佳的分類結果。推測可能原因為特徵的分布，例如 高程為光達分類中一個相當重要的特徵，但建物從 樓層數不一，高度的範圍從三公尺到數十公尺，分 布的情況會形成多個中心(如圖 15)，很難以參數化 的機率密度函數表示，以此資料代入需假設機率密 度函數的分類器將影響分類成果，因此以決策樹為 基礎的隨機森林分類器較適合使用在光達分類。
案例三顯示闊葉樹與竹林在本研究中使用的 各種特徵分布的情況相當類似，因此分類的成果不 盡理想。未來可嘗試使用影像紋理與光達紋理資訊 分離此二類別，而光達紋理或可考慮使用基於灰度 共生張量場(Gray Level Co-occurrence Tensor Field, GLCTF) (賴哲儇，2009; Tsai and Lai, 2013)所產生 的三維紋理，以全波形光達特徵內插成三維網格進 行分析。
致謝
本 研 究 承 蒙 國 科 會 計 畫 NSC100-2628-E008-014-MY3 經費支助，及資料 提供者國立中央大學太空及遙測研究中心林唐 煌副教授及中興測量有限公司，特此致謝。
參考文獻
賴哲儇，2009。高光譜影像立方體於特徵空間之三 維紋理計算，碩士論文，國立中央大學土木工 程學系。
Alexander, C., Tansey, K., Kaduk, J., Holland, D., and Tate, N.J., 2010. Backscatter coefficient as an attribute for the classification of full-waveform airborne laser scanning data in urban areas. ISPRS Journal of Photogrammetry and Remote Sensing, 65(5): 423-432.
Breiman, L., 2001. Random forests. Machine Learning, 45(1): 5-32.
Bretar, F., Chauve, A., Mallet, C., and Jutzi, B., 2008.
Managing full-waveform lidar data: a challenging task for the forthcoming years. In:
International Archive of Photogrammetry.
Remote Sensing and Spatial Information Sciences, XXXVII (Part B1): 415-420.
Briese, C., Höfle, B., Lehner, H., Wagner, W., Pfennigbauer, M., and Ullrich, A., 2008.
Calibration of full-waveform airborne laser scanning data for object classification. In:
Proceedings of SPIE-The International Society for Optical Engineering, 6950 (Laser Radar Technology and Applications XIII).
Chauve, A., Vega, C., Durrieu, S., Bretar, F., Allouis, T., Deseilligny, M.P., and Puech, W., 2009.
Advanced full-waveform lidar data echo detection: assessing quality of derived terrain and tree height models in an alpine coniferous forest. International Journal of Remote Sensing,
30(19): 5211-5228.
Fieber, K.D., Davenport, I.J., Ferryman, J.M., Gurney, R.J., Walker, J.P., and Hacker, J.M., 2013.
Analysis of full-waveform LiDAR data for classification of an orange orchard scene. ISPRS Journal of Photogrammetry and Remote Sensing, 82: 63-82.
Guo, L., Chehata, N., Mallet, C., and Boukir, S., 2011.
Relevance of airborne lidar and multispectral image data for urban scene classification using random forests. ISPRS Journal of Photogrammetry and Remote Sensing, 66(1):
56-66.
Heinzel, J., and Koch, B., 2011. Exploring full-waveform LiDAR parameters for tree species classification. International Journal of Applied Earth Observation and Geoinformation, 13(1): 152-160.
Hofle, B., and Pfeifer, N., 2007. Correction of laser scanning intensity data: data and model-driven approaches. ISPRS Journal of Photogrammetry and Remote Sensing, 62(6): 415-433.
Levenberg, K., 1944. A method for the solution of certain problems in least squares, Quart. Applied Math., 2: 164-168.
Lin, Y.C., Mills, J.P., and Smith-Voysey, S., 2010.
Rigorous pulse detection from full-waveform airborne laser scanning data. International Journal of Remote Sensing, 31(5): 1303-1324.
Liu, X., 2008. Airborne LiDAR for DEM generation:
some critical issues. Progress in Physical Geography, 32(1): 31-49.
Mallet, C., Bretar, F., Roux, M., Soergel, U., and Heipke, C., 2011. Relevance assessment of full-waveform lidar data for urban area classification. ISPRS Journal of Photogrammetry and Remote Sensing, 66:
S71-S84.
Marquardt, D.W., 1963. An algorithm for least-squares estimation of nonlinear parameters.
Journal of the Society for Industrial and Applied Mathematics, 11(2): 431-441.
Podobnikar, T., and Vrečko, A., 2012. Digital elevation model from the best results of different filtering of a LiDAR point cloud. Transactions in GIS, 16(5): 603–617.
Qin, Y., Vu, T.T., Ban, Y., and Niu, Z., 2012. Range determination for generating point clouds from airborne small footprint LiDAR waveforms.
Optics Express, 20(23): 25935-25947.
Reitberger, J., Schnorr, C., Krzystek, P., and Stilla, U., 2009. 3D segmentation of single trees exploiting
full waveform LIDAR data. ISPRS Journal of Photogrammetry and Remote Sensing, 64(6):
561-574.
Sohn, G., and Downman I., 2007. Data fusion of high-resolution satellite imagery and LiDAR data for automatic building extraction. ISPRS Journal of Photogrammetry and Remote Sensing, 62(1): 43-63.
Stal, C., Tack, F., Maeyer, P., Wulf, A., and Goossens, R., 2013. Airborne photogrammetry and lidar for DSM extraction and 3D change detection over an urban area – a comparative study.
International Journal of Remote Sensing, 34(4):
1087-1110.
Tsai, F., and Philpot, W., 1998. Derivative analysis of hyperspectral data. Remote Sensing of Environment, 66(1): 41-51.
Tsai, F., and Chou, M.J., 2006. Texture augmented analysis of high resolution satellite imagery in detecting invasive plant species. Journal of the Chinese Institute of Engineers, 29(4): 581-592.
Tsai, F., Lin, E.K., and Yoshino, K., 2007a. Spectrally segmented principal component analysis of hyperspectral imagery for mapping invasive plant species. International Journal of Remote Sensing, 28(5-6): 1023-1039.
Tsai, F., Chang, C.K., Rau, J.Y., Lin, T.H., and Liu, G.
R., 2007b. 3D computation of gray level co-occurrence in hyperspectral image cubes.
Lecture Notes in Computer Science, 4679:
429-440.
Tsai, F., and Lai, J.S., 2013. Feature extraction of hyperspectral image cubes using three-dimensional gray-level cooccurrence.
IEEE Trans. Geoscience and Remote Sensing, 51(6): 3504-3513.
Ullrich, A., and Reichert, R., 2005. High resolution laser scanner with waveform digitization for subsequent full waveform analysis. SPIE, 5791:
82-88.
Wagner, W., Ullrich, A., Ducic, V., Melzer, T., and Studnicka, N., 2006. Gaussian decomposition and calibration of a novel small-footprint full-waveform digitising airborne laser scanner.
ISPRS Journal of Photogrammetry and Remote Sensing, 60(2): 100-112.
Yao, W., Krzystek, P., and Heurich, M., 2012. Tree species classification and estimation of stem volume and DBH based on single tree extraction by exploiting airborne full-waveform LiDAR data. Remote Sensing of Environment, 123:
368-380.
