Volume 23, No.3, 2018, pp. 141-156 DOI 10.6574/JPRS.201809_23(3).0001
1 國立交通大學土木工程學系 博士 收到日期：民國 107 年 04 月 02 日
2 國立交通大學土木工程學系 教授 修改日期：民國 107 年 05 月 23 日
* 通訊作者, 電話: 0912677440, E-mail: elvislikebear@gmail.com 接受日期：民國 107 年 07 月 12 日
應用整合位移偵測演算法於臺灣西部河道變化之研究
吳俊毅
1*史天元
2摘要
本研究利用整合位移偵測方法測試短時期劇烈變化之大安溪中游研究區，將測試成果以實測資料及 前人研究比對，驗證本研究之偵測方法之可行性，探討河道變遷情形。研究成果顯示，大安溪河段地形 隆起後，從 2001 至 2010 年間共可分成 3 個階段，分別為 2001 至 2003 年、2004 至 2006 年及 2007 至 2010 年。河道一開始並無固定流路，行水路徑逐年變化，在氾濫平原區域擺盪，至侵蝕到護甲層流失造成岩 床裸露後，確定河道位置，接著河道開始下切，其下切最大值為 15 m。整合位移偵測方法考量河道水平 擺盪情形加上 DSM 相減成果，確實可擷取出河道之真實三維變化，並詳細描述河道變遷歷程。
關鍵詞： PIV、DEM 相減、河道變遷、多時序分析
1. 前言
臺灣地區板塊運動劇烈，地體抬升和侵蝕作用 造就地形變化劇烈，南北向山脈縱貫臺灣成為主要 分水嶺，河流以東西向為主，特性為流短坡陡、水 流湍急的急流型河川。當河道受到塊體運動和沖刷 作用，使得河道兩旁逐漸開展，河床擴大，逐漸形 成 V 型河谷，河谷地形之橫剖面由於河道側蝕及 下切作用影響，由外圍至河道中心形成谷壁、台地、
氾濫平原、自然堤及主河道等地形特徵。
多時序分析常應用在災害發生前後變遷偵測 (Ren et al., 2013、Anders et al., 2013、Orem and Pelletier, 2015)、河川地形變遷(Croke et al., 2013) 研究上，收集不同時期之測量數據(如：實測 GPS 點位、航空照片、衛星照片及數值地表模型(Digital Surface Model, DSM))以不同演算法分析數據，獲 得地形變化量，有助於了解研究區變化情形。
傳統河道變遷研究利用航空攝影及遙感探測 技術時，經常使用正射影像或遙測影像，以人工方 法數化河道邊界，套疊不同時期之邊界線型以獲得 水平擺盪變化，並利用影像分類技術進行地表分類，
比較不同時期之分類成果，確定研究區變化情形。
另外亦有：使用多時期高程數據相減，獲得高程變
化情況；以剖線擷取高程值，探討高程變化分析；
以及計算地形測計參數，比較不同時期之各參數之 變化量等方法。利用上述方法進行河道變遷之研究，
各有其優缺點，河道邊界人工數化之工作繁複，須 考慮數化人員之專業能力，影像品質及分類演算法 亦影響影像分類之成果。另外在傳統河道變遷方法 中，水平變化和高程變化通常用於分別單獨探討變 化情形，無法在考量研究區之水平變化後，計算高 程變化情形，只能單純顯示相同位置在不同時間之 高程變化。
為了能夠快速獲得河道三維變化，本研究提出 新的整合位移偵測方法，整合位移偵測方法結合影 像匹配技術及 DEM 相減方法，其中影像匹配技術 係使用 PIV (Particle image velocimetry)方法，PIV 技術已有許多用於偵測地表特徵物於二時期影像 間相對之位移量，分析其位移情況特性之研究(謝 宗霈，2007、羅佳明等，2009 及 Tseng et al., 2009)。
本研究選擇臺灣西部河道地形變化劇烈區域作為 方法驗證及演進過程討論，研究區為大安溪中游河 道地區，大安溪主要受自然因素(地震導致地殼隆 起)影響，短時期地形變化劇烈。本研究利用不同 時期航空照片萃取大安溪河道地區短時期(2001 年
~2010 年)水平擺盪變化，分析主河道、自然堤及 氾濫平原之變化情形，接著以相同地區之不同時期 的 DEM 數據，將其相減，獲得河道地區高程變化 情形河道地區著重在主河道及氾濫平原區域之變 化。最後考慮研究區內點位水平變化，計算其隨著 時間經過，相同點位在不同水平位置之高程變化，
獲得研究區三維變化情形，探討研究區行水區區域 之歷年變化，接著以實測資料及前人研究驗證，確 定是否能夠偵測出行水區在河道範圍內(主河道、
自然堤及氾濫平原區域)之變化，並探討其地形演 進過程。
1.1 PIV
PIV 為流場量測技術的一種，中譯為質點影像 測速儀。該技術原多用於流體工程之流場觀測，在 不干擾流場情況下，得到整個流場在一微小時距下 之瞬時速度場(Landreth et al., 1988、Lecordier et al., 1994)。PIV 優點為計算兩張影像不同時間目標物 移動量時，不需要直接相關。套用在地滑移位應用 時，PIV 法是方便且迅速的方法，偵測人力較難到 達區域之變遷情況。PIV 影像測速系統中，兩張影 像經計算後之位移場，為前後期影像待測區域 (interrogation area)中所計算出之平均位移量。單次 曝光之影像對中，前後期影像藉由互相關函數運算 於相同待測區域中計算出位移量，如圖 4 所示，於 互相關計算中所得之相關係數繪製於平面上，最高 峰值處位置與待測區域中心坐標差為計算後待測 區之平均位移量。PIV 計算公式如式(1)所示：
J j
K
k
J
j K
k J
j K
k
v k u j f k j f
v k u j f k j f v
u
1 0 1 0
2 1 2
2
1 0
2 1
) , ( )
, (
) , ( ) , ( )
, C(
... (1) C(u, v)為點位(u, v)偵測之二維水平變化量；f1
和 f2是兩張照片之遮罩區，分別計算 f1(j, k) 和 f2(j-u, k-v)之相關係數，(j, k)為遮罩窗中之點位位置，
J、K 為遮罩窗大小。本研究遮罩窗使用 64 像元×
64 像元大小，每次移動 32 像元繼續計算相關係數，
計算之最大相關係數像元，即可獲得擺盪量。
本 研 究 使 用 PIVTEC GmbH 公 司 研 發 之 PIVview 軟體計算 PIV 成果計算，設定參數後，包 含搜尋視窗、匹配方法、匹配成果之濾波選擇，獲 得研究區水平擺盪量，並將水平擺盪成果輸出。
2. 研究數據及研究方法
本章分別就研究區進行介紹，說明研究區之 地理位置及相關文獻探討，接著描述使用數據及 研究方法。
2.1 研究區介紹
大安溪流域位於臺灣中部(苗栗縣南部及臺中 縣北部)，發源自雪山山脈之大霸尖山西側，東側 (上游地區)地勢較陡峭，愈往西側地勢逐漸平坦。
本河段出露之地層均屬上新世之卓蘭層，岩性主要 包含極厚層塊狀砂岩、極厚層頁岩、砂頁岩互層等，
其中除極厚層塊狀砂岩較具有抗侵蝕能力外，頁岩 及砂岩互層易發生侵蝕。由歷年 GPS 和水準點資 料，大安溪流域因通過多條活動構造，地體屬於抬 升區域，構造運動頻繁(郭基賢，2014)。
1999 年 9 月 21 日發生芮氏規模 7.3 集集大地 震，此次地震為車籠埔斷層活動引起，大安溪蘭勢 大橋附近河道於地震時產生兩道地表破裂面，造成 兩道破裂面之間約 1 km 長的河道向上抬昇，抬昇 之河道破壞原有河道侵淤平衡，使得侵蝕面產生變 化，原河道上覆蓋卵礫石層經過數公尺的降挖疏浚，
造成軟礫石層迅速流失而致使岩床裸露，直接受到 水流侵蝕，加上此河段為軟弱岩床，造成河道快速 向下侵蝕，形成峽谷地形。地震後的構造抬升先使 大安溪河川流路受阻、寬度從大約 40 m 增加至 100 m 以上；2003 年後受到下蝕、側蝕與搬運作用影 響，產生峽谷地形，河道邊界位置無顯著改變，但 河道寬度減少至大約 30 m；至 2008 年時，河道側 蝕、下蝕與搬運作用依然旺盛，河道位置逐漸穩定，
2008 年後變化趨緩(陳宣安，2012)。另有關大安溪 大峽谷地區河道寬度之研究，Cook et al. (2014)使 用 2004 至 2010 年 7 個時期航空攝影測量數據和 2010 至 2013 年 5 個時期光達數據監測大安溪大峽 谷地區之河道寬度，大安溪峽谷地區受到水流作用
的侵蝕而拓寬，拓寬速度一開始較為劇烈(5 m/yr)，
至 2008 年後較為趨緩(1.5 m/yr)，上游地區的水流 流入峽谷地區，導致水流作用對峽谷入口產生強烈 侵蝕，將峽谷地形入口往前推進 90 m 研究區範圍 圖如圖 1 所示，圖中顯示研究區域內有兩條斷層破 裂帶和一條背斜構造，這兩道破裂面為東勢背斜構 造兩翼之褶皺崖，抬升地塊則屬於東勢背斜構造的 一部分。抬升河段岩體突出於河道上，造成砂石淤 積於上游；拱隆河段上由於坡降變陡，且原上覆卵 礫石層進行過數公尺的降挖疏濬，造成軟礫石層迅 速流失而致使岩床裸露，直接受到水流侵蝕，致使 軟弱岩床快速向下侵蝕形成峽谷狀之深槽流路。
由於構造線之分布，將研究區分成 A、B、C 三區，後續分析會以此作為分區探討河道變遷過程 及影響原因。經濟部水利署水利規劃試驗所(2010) 指出 921 地震後，河道流路不明顯，至 2004 年時，
河道均屬於寬廣的多流路型態；2004 至 2006 年河 道逐漸變成窄深流路；2006 年後河道位置大略固 定，有些河道有明顯拓寬跡象。由收集之多時期河 道斷面測量資料，河道高程變化快速侵蝕期間為 2003 年至 2009 年。
圖 1 大安溪研究區範圍圖
表 1 大安溪數據介紹表
年份 拍攝日期 產製方法 儀器 2001 11/12 航空攝影測量 RMKTOP15 2002 9/16 航空攝影測量 RMKTOP15 2003 8/26 航空攝影測量 RMKTOP15 2004 10/3 航空攝影測量 RMKTOP15 2005 10/27 航空攝影測量 RMKTOP15 2006 10/25 航空攝影測量 RMKTOP15 2007 1/31 航空攝影測量 RMKTOP15 2008 6/10 空載光達 Leica ALS50 2009 7/23 空載光達 Optech
ALTM3070 2010 9/12 空載光達 Leica ALS50
在大安溪河道地區，本研究分別使用 2001 至 2010 年，共 10 個年度航攝相片，數據拍攝日期及 相關介紹如表 1 所示。其中 2001 至 2007 年為農林 航空測量所以傳統框標式相機(RMKTOP15)所拍 攝之航照影像，航高 2400 m，相機焦距 152 mm。
2008 年、2009 年及 2010 年使用空載光達點雲之 DSM 成果，空載光達儀器各年度分別為 2008 年使 用 Leica ALS50；2009 年使用 Optech ALTM3070；
2010 年使用 Leica ALS50。利用收集之航空照片，
計算空中三角測量，獲得影像外方位參數，接著使 用影像匹配方法產製 DEM 及正射影像，DEM 及 正射影像解析度為 1 m。另外收集四個時期實測河 道橫剖面數據，以此作為真值和演算法計算成果進 行比較分析，驗證演算法計算出河道變化量之正確 性，河道橫剖面編號分別為剖面 44、44-1 及 45，
實測河道橫剖面數據測量時間為 2000 年 3 月、2007 年 12 月、2008 年 8 月和 2008 年 12 月，其中 2000 年 3 月和 2007 年 12 月河道橫剖面數據由水利署第 三河川局產製；2008 年 8 月和 2008 年 12 月河道 橫剖面數據由水利署水利規劃試驗所產製。
2.2 研究方法
本研究使用多時期正射影像及 DEM，利用影 像匹配演算法計算，獲得地震發生後 10 年期間短 時期河道歷年水平擺盪變化量，並利用 DEM 相減，
獲得河道高程變化。綜合兩者成果，獲得河道地區 地形三維變化，因研究區河道變化劇烈，每年河道 位置皆發生改變，本研究利用整合位移偵測方法建
立短時期河道變化之過程，進而分析河道行水區之 變化，探討臺灣地區行水路徑之變化歷程。
有關多時期航空攝影測量成果，首先使用歷年 拍攝之航照影像，由於航照影像年代較久遠，無法 獲得控制點，本研究由 1995 年農林航空測量所出 版 之 相 片 基 本 圖 掃 描 檔 (1/5000 圖 幅 ； 山 區為 1/10000)點選控制點，另外由於外方位參數已知值 未完全知道，必須倚靠人工量測大量的影像連接點 及地面控制點，以建立像空間及物空間之關係。空 中三角測量平差解算首先以自由網平差進行粗差 偵測，得到觀測精度之估值，其次強制附合控制點 平差。另外由軟體進行影像連接點自動匹配，但自 動匹配的影像連接點常有錯誤匹配情形，需輔以人 工選取重新量測，提高影像連接點網形強度。確定 控制點及影像連接點精度符合條件後，進行空中三 角測量平差解算，求取航照影像外方位元素，還原 地表資訊在物空間座標系統下之型態，進而分析地 面高程起伏變化。
產製數值高程模型精度受空中三角測量解算 成果影響甚鉅，由於歷史航照無法如新拍航照可預 先佈設地面控制標，獲得穩固的控制點位，故必須 檢視歷史資料，尋找適合的點位作為平面控制點和 高程控制點，以求得準確的地面高程。另影像品質 決定影像連接點匹配情形，影像品質較佳可獲得較 好之影像連接情形，提升空中三角解算成果。
本研究使用 BAE system SOCET GXP (BAE System, 2011)進行空中三角解算，並產製數值地表 模型和正射影像。首先設定投影及相機模式，加入 影像並產製影像金字塔，輸入內方位參數。接著自 動匹配影像連接點，檢查連結點是否合宜，若無法 建立影像相對關係，須持續檢查影像連接點連接情 形，並人工新增影像連接點，直至影像間之連接情 形確定。接著輸入控制點，確定控制點分布及平差 精度符合需求，否則將持續確定控制點位置或新增 控制點至符合需求。完成連接點及控制點之量測後，
便可執行空中三角解算。如空中三角解算成果符合 需求精度，即可產製 DSM。若誤差過大則須除錯 及修正，重覆解算直至成果小於需求精度。在有外 方位元素起始值之情況下，自動匹配產生連接點較
可靠。但由於較舊時期之航空照片在拍攝時無 GPS 及 IMU 輔助，因此沒有已知外方位參數，必須依 靠人工作業選取大量影像連接點及控制點，建立影 像物空間關係，需耗大量時間。外方位參數解算完 成後，以自動影像匹配方式產製多時期高程數據。
產製多時期高程及影像數據後，利用多時期正 射影像以 PIV 影像匹配演算法萃取河道擺盪量，
探討大安溪短時間河道變化過程，分析河道水平擺 盪變化情形及探討行水路徑之變化狀況。另外將多 時期 DEM 相減，獲得河道高程變化，DEM 相減 計算公式如下列式(2)所示：
dDSM = DSM(t2) − DSM(t1) ... (2) dDSM 為高程較差值；DSM(t1)為時間 t1 之 高程值；和 DSM(t2)為時間 t2 之高程值。
原本兩者是分開各自計算河道水平及高程變 化，只能分別探討計算成果，但河道地區因為時 間推移後，產生水平擺盪，直接利用相同水平位 置之高程相減只是考慮空間中相同位置之高程變 化，並未考慮相同點位在不同時間時，水平位置 已經產生變動，如此一來無法正確得知同一點位 在不同時間之三維變化。接著整合兩種方法之計 算成果(PIV 技術和 DEM 相減方法，PIV 技術可 獲取水平方向之變化量，DEM 相減可獲得高程方 向之變化量)，提出新的整合方法獲取研究區真實 三維變化，一開始是點位(x0, y0)在初始位置之高 程值(x0, y0, ztime0)，後來經過時間推移，該點位 產生 2 維水平變化(dx, dy)，這個是由 PIV 方法偵 測所得，該點位(x1, y1)位置變成(x0+dx, y0+dy)，
因此該點位之高程值為(x0+dx, y0+dy, ztime1)，高 程較差則是將兩時期之高程值相減，計算公式如 下式(3)所述：
dx = x1 − x0 dy = y1 − y0
dz = DSM(x1, y1) − DSM(x0, y0)... (3) DSM(x0, y0) 是 點 位 (x0, y0) 之 高 程 值 ； DSM(x1, y1) 是點位(x1, y1) 之高程值；和 dz 整合 位移偵測方法之高程較差值。
將計算成果和實測橫斷面予以比較分析，驗證 整合偵測方法之精度。驗證後探討擷取之河道三維 變化量變遷過程，了解研究區短時期之變化，並描 述其變化影響因子。
3. 研究成果
本章以整合位移偵測方法描述短時期河道地 區變化情形，接著利用實測橫剖面數據驗證前述成 果，確定本研究所使用之方法偵測之變化過程和實 際變化之相關性，最後綜合兩者成果，提出整合位 移變化向量，論述此河段在集集大地震後造成河道 地區地殼隆起，河道短時期變遷情況。
(a) 2001 年正射影像 (b) 2009 年正射影像
(c) 2001 年 DSM (d) 2009 年 DSM 圖 2 2001 年正射影像(a)和 DSM(c)及 2009 年正射
影像(b)和 DSM(d)
圖 2 為 2001 年及 2009 年正射影像和 DSM，
由圖可知本研究區河道依地質構造線可分成 3 個 分區，2001 年時，分區 A 有 1 條主要河道；分區
B 沒有主要河道且地勢平坦；分區 C 有 2 條支流流 過。2009 年時，分區 A 依然有 1 條主要河道，但 是水平位置已經改變，不在原始位置上；分區 B 有 1 條固定河道，河道狹窄且陡峭，河道下切作用 劇烈，超過河道側蝕作用；分區 C 由 2 條支流變 成只有 1 條主要河道。另外由高程分布可看出在 2001 年河道地區並無明顯高程變化，至 2009 年時，
由於河道下蝕作用影響，可明顯於高程分布看出河 道位置。由上述描述可知，從 2001 年至 2009 年間，
河道地區短時間產生相當大的變化，為了探究其變 遷過程，本研究使用 PIV 方法計算河道水平擺盪 量，作為河道變遷分析之依據，以下分別探討其計 算成果及成果驗證，並探討河道變遷過程。
3.1 水平擺盪成果
本研究先以 2001 年作為基準，之後年份和 2001 年正射影像和 DSM 進行 PIV 演算法計算，獲 得相較於 2001 年之水平擺盪量，以此探討河道水 平擺盪變化過程，另外將相隔年份之正射影像進行 PIV 演算法計算，獲得每年河道擺盪情形，分析那 個年度造成河道大範圍變化，探討河道變化原因。
如果水平擺盪量逐年增加的話，表示河道位置隨著 時間一直在變化；如果增加到某一年度後，便維持 不變的話，代表河道水平位置已趨於穩定。2007 和 2001 年正射影像及 DSM PIV 成果分布如圖 3 所示，圖 3(a)為 2007 年及 2001 年正射影像 PIV 成 果分布，水平擺盪較大之區域主要分布在河道地區，
顯示河道地區在 2001 年至 2007 年之間，河道位置 並不固定，行水路徑每年皆在改變，水平擺盪以東 北方向為主，分區 A 至 C 皆有偵測到水平擺盪成 果。圖 3(b)為 2007 年及 2001 年 DSM PIV 成果分 布，相較於正射影像，DSM 之 PIV 成果偵測範圍 較小，只有分區 A 和 C 有偵測到水平擺盪成果，
原因推斷為高程變化較劇烈區域(分區 B)，由於高 程灰階值變化太快，PIV 演算法無法偵測其水平擺 盪量，但可偵測之像元，兩者趨勢一致。綜合上述 討論，顯示除了傳統使用影像對計算其擺盪量外，
DSM 也可用來計算擺盪量，只是在高程變化劇烈 地區較難偵測出成果。
圖 3 2001 和 2007 年 PIV 成果，(a)為正射影像成 果、(b)為 DSM 成果
為了驗證 PIV 成果，本研究收集 3 條實測橫 剖面，將 PIV 成果擷取出來和實測數據進行比較，
驗證 PIV 成果。圖 4 至圖 6 為 3 條實測橫剖面(44、
44-1 和 45)之 PIV 成果剖面，以 2001 年為基準，
橫軸為剖面距離，縱軸為水平擺盪值。圖 4 為 2002 至 2010 年橫剖面 44 之 PIV 成果(以 2001 年作為基 準)，河道位置位於剖面距離 300-600 m 處，2003 年時在剖面距離 550 m 處有最大水平擺盪量 2 像元；
2004 年時水平擺盪較大之區域開始擴散，顯示河 道行水路徑開始擴張，剖面距離 300-600 m 皆偵測 出水平擺盪，最大水平擺盪量達到 30 像元；2010 年時在剖面距離 500 m 處有最大水平擺盪量達到 50 像元。圖 5 為 2002 至 2010 年橫剖面 44-1 之 PIV 成果(以 2001 年作為基準)，2002 年時在剖面距離 250 m 處有最大水平擺盪量 15 像元，接著在 2003 年時，水平擺盪較大之區域開始擴散為剖面距離 200-300 m，有最大水平擺盪量 30 像元。2003 年 後河道區域之水平擺盪便呈現穩定狀態，不再有大 幅度變化。圖 6 為 2002 至 2010 年橫剖面 45 之 PIV 成果(以 2001 年作為基準)，可看出此剖面有兩條 支流，分別位於剖面距離 100-200 m 和剖面距離 450-600 m 處，最大水平擺盪量約 30 像元，從 2002 年變化後，之後年份水平擺盪變沒有太大改變，直 至 2010 年時，在剖面距離 100-200 m 有最大水平 擺盪量 50 像元。三條橫剖面偵測之河道位置，和 正射影像相比較，符合河道所在位置，可見 PIV 成 果確實能偵測出河道變化情況，對於行水路徑之偵 測確實有效。
將橫剖面之河道區域 PIV 成果進行統計分析，
探討其河道擺盪變化，其中剖面 45 因為有兩個區 域，故分開進行討論。圖 7 為歷年剖面河道位置 PIV 成果(以 2001 年為基準)平均值及標準差，圖 7(上)可看出剖面 44 隨著時間推移，河道地區擺盪 量逐年增加；剖面 44-1 在 2004 年達到25 像元後，
接下來年份變化程度很小；剖面 45 在一開始擺盪 量較大，之後便逐漸減小，而剖面 45-2 則是擺盪 變化沒有固定趨勢，顯示此河道地區無固定流路。
標準差部分則是各剖面皆在 2005 年達到最大值，
之後標準差值皆較小，顯示河道擺盪程度有穩定 趨勢。
圖 4 2002 至 2010 年橫剖面 44 之 PIV 成果，以 2001 年為基準(橫軸單位：m；縱軸單位：像元)
圖 5 2002 至 2010 年橫剖面 44-1 之 PIV 成果，以 2001 年為基準(橫軸單位：m；縱軸單位：
像元)
圖 6 2002 至 2010 年橫剖面 45 之 PIV 成果，以 2001 年為基準(橫軸單位：m；縱軸單位：像元)
圖 7 歷年剖面河道位置 PIV 成果平均值(上)及標準 差(下)，以 2001 年為基準(縱軸單位：像元) 接著以相隔年份進行 PIV 成果計算，獲得每 個年度河道擺盪量，以此作為河道變化之依據，
計算之 PIV 成果和 3 條實測橫剖面比較。圖 8 至 10 為 3 條實測橫剖面(44、44-1 和 45)之相隔年份 PIV 成果剖面，橫軸為剖面距離，縱軸為水平擺 盪量。圖 8 為 2002 至 2010 年橫剖面 44 之相隔年 份 PIV 成果，河道位置位於剖面距離 300-600 m 處，
2003 年至 2004 年之間發生較大水平擺盪量，約 40 像元，接著 2007 年和 2008 年、2008 年和 2009 年 及 2009 年和 2010 年每年約有 20 像元之水平擺盪 量，顯示河道側蝕能力一直持續在此區作用，由 相隔年份之 PIV 成果也可看出行水區域之變化情 形。圖 9 為 2002 至 2010 年橫剖面 44-1 之相隔年 份 PIV 成果，河道位置在剖面距離 100-200 m，
除了在 2003 年至 2004 年之間有水平擺盪量約 20 像元外，其餘年份皆無水平擺盪，由此可知此區 河道位置只有在 2003 年產生水平擺盪後，河道行 水區域產生變化，之後便不再移動。圖 10 為 2002 至 2010 年橫剖面 45 之相隔年份 PIV 成果，河道 位置在剖面距離 100-200 m 和剖面距離 450-600 m
處，主要有水平擺盪量的年份為 2004 和 2005 年，
最大約 40 像元，其餘 2005 和 2006 年、2007 和 2008 年及 2008 和 2009 年依然有約 30 像元擺盪 量。
接著對相隔年份 PIV 成果之河道位置進行統 計分析，其中剖面 45 因為有兩個區域，故分開進 行討論。圖 11 為相隔年份剖面河道位置 PIV 成果 平均值及標準差，由圖 11(上)可看出剖面 44 和剖 面 44-1 在 2003 年至 2004 年之間發生較大水平擺 盪量，其餘年份之擺盪量則較小；剖面 45-1 在 2004 至 2005 年有較大水平擺盪量；剖面 45-2 在 2003 至 2005 年及 2007 和 2008 年有較大水平擺盪 量，其餘年份之擺盪量則較小。標準差之統計分 布情形和平均值趨勢一致，顯示水平擺盪量大的 時間段，其擺盪量分布情形也較散亂。
圖 8 2002 至 2010 年相隔年份橫剖面 44 之 PIV 成果，
(橫軸單位：m；縱軸單位：像元)
圖 9 2002 至 2010 年相隔年份橫剖面 44-1 之 PIV 成 果，(橫軸單位：m；縱軸單位：像元)
圖 10 2002 至 2010 年相隔年份橫剖面 45 之 PIV 成 果，(橫軸單位：m；縱軸單位：像元)
圖 11 相隔年份剖面河道位置 PIV 成果平均值 (上)及標準差(下)，(縱軸單位：像元)
圖 12 為正射影像 PIV 成果統計分布圖(以 2001 年為基準)，其中正射影像之 PIV 成果由於在不同 年份有不同統計特性，因此將其分成 3 個階段討 論，分別為 2002 年和 2003 年、2004 至 2007 年及 2008 年至 2010 年。圖 12 顯示在 2002 及 2003 年 時，水平擺盪量有 80%及 73%分布在 0-5 像元，最 大水平擺盪量為 40 像元。接著水平擺盪量小於 5 像元的網格數逐漸變少(70%至 50%)，水平擺盪量 往 5-10、10-15 像元之數量開始增加，由此顯示從 2004 年起，河道產生較大之水平變化，且最大水 平擺盪量也增加至 50 像元。河道水平擺盪量另外 一個轉折點是在 2008 年時，水平擺盪量 0-5 像元 只佔 33%，10-15 及 15-20 像元數量增加，顯示河 道水平擺盪量再次增加，且最大水平擺盪量也增 加至 60 像元。
圖 13 為相隔年份之正射影像 PIV 成果統計分 布圖，共可分成 3 種分布情形，第一種是 2002 年 和 2003 年及 2006 年和 2007 年兩個時期，幾乎水 平擺盪量有 98%皆小於 5 像元；第二種是 2004 年
和 2005 年、2005 年和 2006 年、2007 年和 2008 年及 2008 年和 2009 年四個時期，相較於第一種，
水平擺盪量小於 5 像元的網格數逐漸變少(80%)，
水平擺盪量往 5-10、10-15 像元之數量開始增加；
第三種是 2003 年和 2004 年及 2009 年和 2010 年兩 個時期，水平擺盪量 0-5 像元只佔 49%及 50%，
10-15 及 15-20 像元數量增加，顯示這兩個時期河 道擺盪情形特別厲害，河道位置產生劇烈改變，
也呼應之前和 2001 年為基準之 PIV 成果，在 2003 年及 2004 年是兩個時期之區隔。
(a) 2002 年和 2003 年
(b) 2004 至 2007 年
(c) 2008 年至 2010 年
圖 12 正射影像 PIV 成果統計分布圖，以 2001 年為基準
圖 13 相隔年份正射影像之 PIV 成果及 DSM 之 PIV 成果統計分布圖 圖 13 也顯示 DSM PIV 成果在 2002 年有 95%
網格之水平擺盪在 0-5 像元以下；2005 年有 66%
網格之水平擺盪在 0-5 像元以下；2008 年有 61%
網格之水平擺盪量在 0-5 像元以下；2010 年有 63%網格之水平擺盪量在 0-5 像元以下，由此可知 2002 年時河道水平變化相較於其他年份不劇烈，
其餘年份趨勢大致相同，河道水平擺盪量至 2005 年後，便不再持續移動。DSM PIV 成果和正射影 像成果趨勢相同，水平擺盪量級之所佔比例相同，
差別只在於偵測到水平擺盪點位的數量較少。
由於正射影像和 DSM 在 PIV 成果有相同趨 勢，但利用正射影像能偵測出較多之水平擺盪情 況，因此我們以正射影像計算之歷年 PIV 成果作 為分析成果。分區 A 擁有比較明確之河道地區，
PIV 成果可以明確獲取河道變遷情況；分區 B 和 C 由於擁有較複雜之地表分布，使得 PIV 成果受 到混合地表種類影響，無法完全反應河道變遷情 形，但是趨勢仍然可以看出。一開始 2002 年只 有零星水平擺盪較大之區域，隨時時間推移，逐 漸往附近擴散，一直到 2004 年時，水平擺盪較 大之區域才固定下來，接下來年份水平擺盪增加
幅度不大，而相隔年份也獲得相同趨勢，在 2003 年至 2004 年間產生較大河道水平擺盪。
3.2 高程變化分析
圖 14 為 2008 及 2001 年 DSM 相減後，高程 較差分布圖。除了 2002 及 2001 年高程較差幅度較 小外，其餘年份高程較差趨勢和 2008 及 2001 年相 同。從 2001 年至 2010 年間，河道下切最大值為 12 m。另外由分區來看，分區 B 之區域，河道高 程下降之區域較其他分區範圍小，顯示分區 B 之 河道下切能力旺盛，但側蝕能力不顯著，造成該 段河段旁邊邊坡陡峭，而分區 A 及 C 則是河道在 下切作用的同時，側蝕作用依據影響河道地形變 化，因此高程變化區域範圍較大，但下蝕作用區 域仍可由高程變化較大的趨勢明顯看出。
為了驗證歷年 DSM 相減成果，和前述 PIV 成 果一樣，使用 3 條實測橫剖面進行比較，圖 15 為 三條實測橫剖面之歷年 DSM 高程分布圖，圖 15(a) 為橫剖面 44，可看出在橫剖面距離 400 m 處，由 2002 年至 2005 年高程下降 10 m，接著便沒有在持 續下切。圖 15(b)為橫剖面 44-1，2002 年至 2005
年時，在橫剖面距離 320 m 處，高程下降 10 m，
然後至 2007 年時，高程下降至 20 m。圖 15(c)為 橫剖面 45，在橫剖面距離 180 m 處，高程下降 10 m，然後至 2009 年時，高程下降至 15 m，河道邊 坡也隨著時間推移，有往外擴散之趨勢。三條橫 剖面分別對應三個分區，分區 B 河道(剖面 44-1) 下切作用最為強烈，下切量也最大，且河道寬度 也最小，形成獨特峽谷地形。
本研究利用河道寬高比探討不同橫剖面之河 道特性，剖面 44 河道寬高比為 10；剖面 44-1 河道 寬高比為 5；剖面 45 河道寬高比為 13。剖面 44-1 之河道寬高比最小，顯示該地區河道下切作用相 較於側蝕作用劇烈，造成河道坡度陡峭。
圖 14 2008 和 2001 年高程較差分布圖
3.3 成果驗證
本節利用 3 條實測橫剖面數據(剖面 44、44-1 及 45)，將其作為真值，驗證前述 PIV 成果和 DSM 相減之高程較差成果。探討本研究所提之 方法是否能夠描述該區短時期地形變遷情況。圖 16 為實測橫剖面之高程分布，橫軸為剖面距離，
縱軸為高程值，共有 4 個時期數據，虛線代表 2000 年時高程分布；粗實線代表 2008 年 10 月高 程分布。圖 16(a)為橫剖面 44 高程分布，河道地 區為剖面距離 200 至 450 m 處，2000 年至 2008
年高程下降 10 m (位於剖面距離 300 m 處)，和 DSM 相減成果相符。圖 16(b)為橫剖面 44-1 高程 分布，河道地區為剖面距離 150 至 250 m 處，
2000 年至 2008 年高程下降 20 m (位於剖面距離 -400 m 處)，河道寬深比為 5，和 DSM 相減成果 相符。圖 16(c)為橫剖面 45 高程分布，河道地區 為剖面距離 50 至 150 m 處，2000 年至 2008 年高 程下降 15 m (位於剖面距離 100 m 處)，和 PIV 成 果及 DSM 相減成果相符。綜合上述三條實測橫 剖面和本研究所使用之演算法計算之成果驗證比 較，PIV 成果及 DSM 相減成果和實測橫剖面所 比對之結果相符，顯示本研究所用之演算法確實 能夠偵測出河道變遷過程。
(a) 橫剖面 44
(b) 橫剖面 44-1
(c) 橫剖面 45
圖 15 實測橫剖面之歷年 DSM 高程分布圖，橫軸 為剖線距離，縱軸為高程值，(橫軸單位：m；
縱軸單位：m)
(a) 實測橫剖面 44
(b) 實測橫剖面 44-1
(c) 實測橫剖面 45 圖 16 實測橫剖面之高程分布
3.4 整合位移偵測
本研究為了獲得研究區真實三維變化，考慮 PIV 成果作為河道水平擺盪變化，接著和 DSM 相 減方法作結合，考量河道擺盪後之高程值和原始 點位之高程值相減，獲得真實三維變化。PIV 成 果由正射影像及 DSM 獲得，本研究也分別探討這 兩種數據所計算之水平擺盪量，其高程變化之差 異比較，將整合位移偵測成果比較分成全部研究 區及河道地區分別討論，因為全部研究區包括河 道地區及河道地區周圍區域，高程較差最大值及
最小值容易受到影響，使得高程較差值高估，無 法真實反應河道變遷情況。
全部研究區平均水平擺盪量為 3.83 至 16.40 像元(正射影像之 PIV 成果)和 0.92 至 7.85 像元 (DSM 之 PIV 成果)；水平擺盪量標準差為 6.67 至 10.23(正射影像之 PIV 成果)和 3.16 至 10.65(DSM 之 PIV 成果)，正射影像之 PIV 成果較 DSM 之 PIV 成果平均值為大，且水平擺盪量標準差也較大，
顯示正射影像之 PIV 成果有較大的河道擺盪量，
但擺盪幅度也較不一致。水平擺盪變化說明 2002 年至 2005 年期間，河道水平變化劇烈，2005 年後 即維持不變。平均高程變化在 2002 年趨近於 0 m，
接下來的年份皆是負值，說明研究區河道下切作 用旺盛，持續造成研究區高程值下降。高程較差 最小值在 2010 年時，整合位移偵測方法為-18.45 m；DSM 相減方法為- 17.86 m，兩者皆符合相關 文獻所記載之河道下切幅度。整合位移偵測方法 使用 DSM 之 PIV 成果在 2002 年至 2010 年時，不 論在水平擺盪和高程較差部分，能夠確實描述河 道三維變遷情況。
由於受到河道地區周圍影響，在整合位移偵 測方法和 DSM 相減成果都有高程較差值被高估的 情況，因此特別將河道地區以外濾除，探討河道 地區三維變化情形，並和原始 DSM 相減成果進行 比較。河道地區之 PIV 成果平均值變大，標準差 則不變，顯示全部區域確實造成河道擺盪量低估。
利用正射影像之 PIV 成果之整合位移偵測方法在 高程較差最大值、高程較差最小值及高程較差標 準差相較於原始 DSM 相減變化較大，而利用 DSM 之 PIV 成果之整合位移變化偵測方法其統計 成果和原始 DSM 相減相符。圖 17 為兩種整合位 移偵測方法(正射影像之 PIV 成果和 DSM 之 PIV 成果)和 DSM 相減高程較差平均值統計圖，可看 出原始 DSM 相減之平均高程較差值和整合位移 變化偵測方法相符，不論是 PIV 成果由正射影像 獲得或是 PIV 成果由 DSM 獲得，皆有相同之趨 勢，顯示經由整合位移偵測方法獲取之統計成果 和原始 DSM 相減成果相同，但是其有考量河道 位置擺盪變化，真實反應河道三維變化。
圖 17 兩種整合位移變化偵測方法(正射影像之 PIV 成果和 DSM 之 PIV 成果)和 DSM 相減 高程較差統計圖
4. 結論
本研究提出新型整合位移偵測方法，針對臺 灣西部軟弱岩層易遭受侵蝕之河道地區，以短時 期變化之大安溪河道地區為測試區，收集研究區 之多時期數值高程模型及正射影像，獲取河道歷 年真實變化，並以此成果分析河道變遷歷程，剖 析河道的沖蝕現象與變化，建立河道沖蝕演繹過 程，進而說明河道變化影響因子。
綜合 PIV 演算法獲得水平擺盪成果及 DEM 相 減獲得高程變化，本研究提出整合位移偵測方法，
在考量水平擺盪量情況下，計算其高程較差，獲 得河道真實三維變化，剖析河道變化過程，研究 成果顯示此區域在 2001 年至 2010 年間，有平均水 平擺盪幅度為 20 m 和高程下切 15 m，變化幅度相 當劇烈。將整合位移偵測方法所計算出之成果，
和實測橫剖面驗證比較，證明本研究之方法和實 測成果相符，並進一步與相關研究成果比對，說 明河道變遷歷程，確定河道變遷影響因子。
大安溪研究區河道位置在地震造成隆起後，
由於河道回春作用而使原有平衡遭受破壞，重新 形成新的地形風貌，直至恢復原有平衡。由研究 成果可將研究區變化情形分成 3 個階段，分別為 2001 至 2003 年、2004 至 2006 年及 2007 至 2010 年，一開始第一階段河道並無固定流路，經侵蝕 軟弱岩層，河道逐年產生變化，行水路徑逐年變
化，在氾濫平原區域擺盪。至侵蝕到護甲層流失 造成岩床裸露後，始有固定河道位置，擺盪量約 20 像元。接著第二階段，河道水平變化趨緩，水 流之行水路徑開始穩定，加上上游輸砂受阻，河 道下蝕作用旺盛，開始出現河道下切。最後第三 階段時，河道持續下切，最大下切量達 15 m，仍 未有減緩之趨勢。如果河道邊坡岩層能夠承受重 量，河道坡度越來越陡，形成陡峭邊坡地形；不 然的話邊坡上岩層便會受不了重力作用，崩塌至 河道位置，造成河道寬度拓寬。
另外從研究區分區可看出分區 A 和分區 C 變 化趨勢相同，首先河道側蝕作用旺盛，河道以水 平擺盪為主，直到河道位置固定後，河道開始下 切，之後在河道側蝕及下切作用之相互影響下，
保持固定的寬高比；而分區 B 河道下蝕作用十分 旺盛，河道下切劇烈，寬高比值逐年變小，形成 大峽谷地形，此區之高程值逐年下降，且坡度也 越來越陡。
整合位移偵測方法經由本研究測試區驗證後，
研究成果證實可偵測出河道短時期水平擺盪變化 及長時期河道變遷偵測，但仍有可精進之處，例 如加入其他量化指標或其他參考資料(水系、構造 線)、人工編修影像、抑或改用特徵匹配偵測水平 擺盪變化，期望能夠更加詳細的描述河道變遷過 程。另外也可考慮將研究區分成許多小區域，各 自探討其擺盪量變化，詳細描述河道水平擺盪情 形，期望獲取更細部的河道真實三維變化情形。
水流之行水路徑變化也可納入其他相關資料進行 分析，如：颱風資訊、強降雨事件發生日期，以 此探討研究區行水區細部變化情形。對於 PIV 方 法偵測之水平擺盪量，也可更加細部探討，將不 動點及擺盪較大之點位區分出來，詳細評估其成 果，以獲得更詳細之河道水平擺盪情形。
本研究研究區之地形變化劇烈，目前尚未達 成平衡，河道下切侵蝕作用依然在測試區造成地 形變化，未來仍須持續監測大安溪之河道下切侵 蝕變化，納入更多時期之觀測數據，俾以更深入 了解其河道變化，提供給有關單位借鏡，避免後
續災害之發生。
參考文獻
陳宣安，2012，大安溪峽谷河川地形變遷之研究，
國立臺灣大學地理環境資源學研究所碩士論 文。[Chen, S.A., 2012. The study on fluvial morphological change of Daan gorge, Master Thesis, National Taiwan University, Taiwan, ROC.(in Chinese)]
郭基賢，2014，利用大地測量和地形計測資料探 討臺灣地殼垂直變動與構造活動之研究，國 立彰化師範大學地理學系博士論文。[Kuo, C.S., 2014. A study of Taiwan crustal vertical movement and tectonic activity with geodetic surveying and morphometry data, PhD Thesis, National Changhua University of Education, Taiwan, ROC. (in Chinese)]
經濟部水利署水利規劃試驗所，2010，軟弱岩床 劇烈沖蝕河段沖蝕行為之探討-以大安溪為例 總報告。[Water Resources Planning Institute, Water Resources Agency, Ministry of Economic Affairs, 2010. Rapid Erosion Process of Soft Rock Riverbed - The Taan River as an Example.(in Chinese)]
謝宗霈，2007。應用影像計算於地表變位之監測 以紅菜坪地滑為例，國立成功大學地球科學 系碩士論文。[Hsieh, C.P., 2007. Monitoring landslide by image multitemporal analysis geodesy : Application of huangtsaiping landslide, Master Thesis, National Cheng Kung University.(in Chinese)]
羅佳明、林銘郎、董家鈞、張光宗、簡士堯、黃 安斌，2009。應用地形分析、遙測影像判釋 與 PIV 技術於紅菜坪地滑特徵及其分區之研 究，中國土木水利工程學刊，21(2):113-128。
[Lo, C.M., Lin, M.L., Tsung, C.K., Chien, S.Y., and Huang, A.B., 2009. Landslide characterization and zonation of hungtsaiping area based on topography, image of remote
sensing and PIV technology, Journal of the Chinese Institute of Civil and Hydraulic Engineering, 21(2): 113-128. (in Chinese)]
Anders, N.S., Seijmonsbergen, A.C., and Bouten, W., 2013. Geomorphological change detection using object-based feature extraction from multi-temporal LiDAR data, IEEE Geoscience and Remote Sensing Letters, 10(6): 1587 -1591.
BAE System, 2011. Available at:
http://www.geospatialexploitationproducts.com/ , Accessed June 10, 2016.
Cook, K.L., Turowski, J.M., and Hovius, N., 2014.
River gorge eradication by downstream sweep erosion, Nature Geoscience, 7(9): 682-686.
Croke, J., Todd. P., Thompson, C., Watson, F., Denham, R., and Khanal, G., 2013. The use of multi temporal LiDAR to assess basin-scale erosion and deposition following the catastrophic January 2011 Lockyer flood, SE Queensland, Australia, Geomorphology, 184:111-126.
Landreth, C.C., Adrian, R.J., and Yao, C.S., 1988.
Double pulse particle image velocimetry with directional resolution for complex flows, Experiments in Fluids, 6(2):119-128.
Lecordier, B., Mouqallid, M., Vottier, S., Rouland, E., Allano, D., and Trinite, M., 1994. CCD recording method for cross-correlation PIV development in unstationary high speed flow, Experiments in Fluids, 17(3):205-208.
Orem, C.A. and Pelletier, J.D., 2015. Quantifying the time scale of elevated geomorphic response following wildfires using multi-temporal LiDAR data: An example from the Las Conchas fire, Jemez Mountains, New Mexico, Geomorphology, 232: 224-238.
Ren, Z.K., Zhang, Z.Q., Dai, F.C., Yin, J.H., and Zhang, H.P., 2013. Co-seismic landslide topographic analysis based on multi-temporal
DEM-A case study of the Wenchuan earthquake, Springerplus, 2(1), DOI:
10.1186/2193-1801-2-544.
Tseng, C.H., Hu, J.C., Chan, Y.C., Chu, H.T., Lee, J.F., Wei, J.Y., Lu, C.Y., and Lin, M.L., 2009.
Non-catastrophic landslides induced by the Mw
7.6 Chi-Chi earthquake in central Taiwan as revealed by PIV analysis, Tectonophysics, 466(3-4): 427-437.
.
1 Ph.D., Department of Civil Engineering, National Chiao Tung University Received Date: Apr. 02, 2018
2 Professor, Department of Civil Engineering, National Chiao Tung University Revised Date: May. 23, 2018
* Corresponding Author, Tel: 886-912677440, E-mail: elvislikebear@gmail.com Accepted Date: Jul. 12, 2018
Geomorphological Change Detection Using an Integrated Method:
A Case Study on the River Channel of Western Taiwan
Jyun-Yi Wu
1*Tian-Yuan Shih
2Abstract
We propose the integration method that combines the Particle Image Velocimetry (PIV) technique and DEM subtraction, which is described to detect the active kinematics of the river morphology in this thesis. While PIV can provide estimates for the channel change direction and magnitude; vertical analysis such as river incision rate derivation could be performed using the height difference between multi-temporal surface models.
It is observed that the proposed scheme is capable of identifying three river evolution stages from the extracted movements in the Taan case. These stages are: years 2001 to 2003, 2004 to 2006, and 2007 to 2010, respectively.
The channel begins to have no stable path, and have regular path from 2001 to 2010. The maximum value of incision between 2001 and 2010 is about 15 m. This study demonstrates that both PIV and DSM subtraction are effective in river geomorphological change identification. The integration of these two approaches could provide more information when observing the evolution of river morphology.
