本研究主要是將快速傅立葉轉換應用於衛星影像圖,繪製出波向線及 推算波浪的入射角並由計算出來的波長估算水深。本文同時收集枋寮海岸 地形資料(民國 99 年度屏東海岸基本資料監測調查計畫報告),此地形資料 可供驗證由影像圖估算水深的準確性。龜山島雖無實測海底地形資料可驗 證,但因其海岸線曲折故有明顯的波浪折射、繞射現象,可說明影像圖推 算波浪之能力。
本文利用直方圖均化法增強影像對比度,使原本不易分辨之波紋線更 明顯,並利用濾波的概念將因波浪反射產生的微小細紋濾除,以提高快速 傅立葉轉換的精度。
3-1 影像前置處理 1. 直方圖均化
直方圖在此表示數位影像圖中統計灰階度像素個數的直方分佈圖,由 於灰階度代表影像的明亮,因此將原先灰階分佈重新均勻等化到新的灰階 值,以達到提高對比的目的。其轉換方程式如下:
) ( f T
g
(3-1)此轉換須滿足下列兩個條件:
(1)T( f)在 fmin f fmax的區間上為單調遞增。
(2)相對於 fmin f fmax,gmin T(f)gmax。
其中 f 為原始灰階度,g 為輸出灰階度。輸入機率密度函數(probability density function)形式和輸出機率密度均勻分佈形式,分別如式(3-2)及(3-3)所示:
max
根據機率密度理論得知其累積機率密度(cumulative density function)相 等,Pg(g)Pf(f),如式(3-4)。將(3-2)及(3-3)式代入(3-4)式得式(3-5)。如圖
min min
)
gray scale
the number of pixel
0 50 100 150 200 250
gray scale
the number of pixel
0 50 100 150 200 250
(a)原影像的直方圖 (b)影像經均化後的直方圖
反應在傅立葉係數中的高頻部分,平滑部分則反應在傅立葉係數中的低頻 部分,依此原理,利用高通濾波即可將低頻的係數衰減,再對處理過的影 像取傅立葉反轉換即可達到高通濾波的效果。
設 f(m,n)為原影像,h(m,n)為濾波器之脈衝響應,g(m,n)為影像濾波後 的結果,原影像與濾波器脈衝響應做迴旋積分(convolution intergal)後可得濾 波後影像表示如式(3-6)。
50 100 150 200 250 50
100
150
200
250
50 100 150 200 250
50
100
150
200
250
圖 3-2 龜山島的部分衛星影像圖
圖 3-3 圖 3-2 經高通濾波後的影像 3. 海陸交界線的擷取
因本文在推算進岸波向時會擷取到陸域影像資料而造成近岸波向推算 不準確,為了增加推算進岸波向的準確性,本文將陸域影像資料取代為海 域影像資料平均值使得方框中取得的資料為海域資料,而海陸分界方式採 常態化水指數標準差(normalized difference water index, NDWI)。
Mcfeeters(1996)提出常態化水指數標準差,用於植物含水率的測定,也 稱為葉面積含水量指標(Leaf area water-absent index)主要是利用綠光波段反
射強度(綠光的波長約 0.50μm)與近紅外光波段反射強度(近紅外光的波長約 0.79μm~0.89μm)的數學和差比來判斷植物葉面的含水量面積,適合用來突顯 含水量高低的差別。NDWI 比值關係式如下式所示:
) (
) (
Green NIR
Green NDWI NIR
(3-10)
其中 NIR 為近紅外光波段反射強度,Green 為綠光波段反射強度。
由於清澈水體對於近紅外光波段的反射率其值接近於零,又反射率最 大值落在綠光波段和紅光波段之間,且土壤對於不同波段的反射率恰與水 體相反,以上所述特徵非常符合常態化水指數標準差(NDWI)的比值關係,
因此衛星影像圖經 NDWI 轉換後,其值大於零可判定為土地,其值小於零 可判定為水體,其值等於零即判定為灘線位置,由此方法即可明顯辨識陸 地及海域範圍。龜山島衛星影像圖經 NDWI 轉換並取二值化後結果如圖 3-4 所示。
圖 3-4 龜山島經 NDWI 及二值化後的灘線位置
3-2 波向推算
間進行轉換,其轉換公式如式(3-11)及式(3-12)所示:
x
(a)龜山島起始區域的位置 圖中,如圖 3-6 所示。
任意選擇一點作為起始中心點,為了滿足傅立葉轉換使用限制第四點,
並且包含足夠的影像資訊,以 512×512 像素的矩陣大小進行二維快速傅立 葉轉換,我們將此矩陣稱為「方框」(frame)。
圖 3-6 龜山島的波向線起始區域選定
步驟 2:方框移動
以往使用此方法皆採用方框兩兩相連的方式,但為了提升推算的精確 度,本文採用每計算一次快速傅立葉轉換後由式(3-13)及式(3-14)推算出波 長與夾角,接著往 x 方向及 y 方向移動 1/4 倍波長的距離並進行第二次計算,
其示意圖如圖 3-7 所示。
此方式雖然會提高計算時間,但明顯表現出波向線的連續性,結果如 圖 3-8 所示。
(b)起始區域的影像圖
(c)起始區域的頻譜等值圖
圖 3-7 方框推進示意圖
圖 3-8 繪製的波向線示意圖 步驟 3:波向入射角的修正
由波浪理論中已知,波浪由遠海傳遞至近岸時波浪入射角並不會有劇 烈的變化,為了使波向線不相互重疊,經測試後每一次推算出的入射角與 前一次推算的入射角相差必須在5以內,若超過範圍,因方框僅移動 1/4 波長,故角度變化小,其修正方式為使用前一次推算出的入射角取代之,
修正前後差異如圖 3-9(a)與圖 3-9(b)所示。
(a)修正前的波向線 (b)修正後的波向線 圖 3-9 波向角修正前與修正後的波向線示意圖 步驟 4:終止計算
為提高近岸推算的精確度,當方框移動至近岸後取海域影像的像素平 均值取代陸地區域的像素,若一個方框中陸地面積占其方框總面積的 1/3 時即停止計算,如圖 3-10 所示,而海陸區隔方式採用常態化差異水指數標 準差。
圖 3-10 龜山島近岸的波向推算方法
步驟 5:加密方框(fine frame)
使用一格為 512*512 的方框僅能計算至離岸約 300 米處,故在方框移 動與陸地接觸時縮減成 256*256,經由此操作即可推算至離岸約 100 米處,
也是港灣結構物施作的範圍,其移動示意圖如圖 3-11 所示。
圖 3-11 近岸波向推算方框加密示意圖
3-3 波向線的繪製
根據上節敘述五步驟並將各點連線即可繪製出一條波向線,如圖 3-12 所示,圖中紅色點表每一格方框的中心點,也作為波浪前進的起始點,黃 色線表波浪的波長,但此方法在波浪繞射區不能有效辨識,故本文在繞射 區部份不做探討。
圖 3-12 繪製龜山島波向線的示意圖 start frame
1200 1400 1600 1800 2000 2200 2400 300
400
500
600
700
800
900
1000
1100
1200
1300
繪製出一條波向線後由起始點往右移動訂出第二個起始點後即進行第 二次波向推算系統,依此方式操作便可推得其餘波向線,如圖 3-13 與圖 3-14 所示。
圖 3-13 龜山島海域的部份波向線
圖 3-14 枋寮海岸的部份波向線
本文推算 96 年 8 月 5 日屏東縣枋寮漁港外側海域的波向並與柯(2012) 使用數值模擬所推算的波向進行角度比較。首先選取 8 條波向線並對應到 數值模擬推算的波向線,如圖 3-15 所示,分別對波向角度進行內插,最後 每條波向線取 323 個點做絕對值誤差及均方根誤差(Root Mean Square Error, RMSE)。RMSE 的定義如式(4-6)所示。
N d d RMSE
N i
c
m
12
(3-15)
其中d 為本模式推算的角度,c d 為數值模擬推算的角度,N 為資料數。此m
值越大表示二者偏差越多,即誤差越大,越不吻合。
本文所提供的方法推算出的波向角度介於 220 度至 232 度之間,而數 值模擬所得的波向角度介於 214 度至 232 度之間,二者的平均絕對值誤差 量為 4.73 度,RMSE 為 5.84 度。
圖 3-15 本文選擇數模推算的波向線(藍色線)的 8 個剖面位置(紅色線)