應用二維經驗模態分解消除水域影像日光返照

國 立 交 通 大 學

土 木 工 程 學 系

碩 士 論 文

應用二維經驗模態分解消除水域影像日光返照

Sun Glint Removal from Water Area Image using

Bidirectional Empirical Mode Decomposition

研 究 生：施伶蒨

指導教授：史天元

應用二維經驗模態分解消除水域影像日光返照

Sun Glint Removal from Water Area Image using Bidirectional

Empirical Mode Decomposition

研 究 生：施伶蒨

Student：Ling-Chien Shih

指導教授： 史天元 Advisor：Tian-Yuan Shih

國 立 交 通 大 學

土 木 工 程 學 系

碩 士 論 文

A Thesis

Submitted to Department of Civil Engineering College of Engineering

National Chiao Tung University in partial Fulfillment of the Requirements

for the Degree of Master of Engineering

in

Civil Engineering July 2012

Hsinchu, Taiwan, Republic of China

I

應用二維經驗模態分解消除水域影像日光返照

學生：施伶蒨 指導教授：史天元

國立交通大學土木工程學系

摘要

水域的航空遙測影像容易因為水面對太陽輻射的鏡面反射而造成日光返照，因而造 成影像資訊損失，影響影像之製圖及解讀，本研究嘗試分析日光返照形成之原理，並使 用能量傳輸解析與影像處理方法測試減低日光返照之影響。

讀取測試影像所得為影像灰度值，無法轉換至輻射值，即使使用Cox and Munk (1954) 模型及 MODO 軟體得到之傳輸量及光源輻射值，亦無法進行輻射值改正；影像拍攝當 時海表面模型未知，所計算之傳輸量無法符合每個點位之差異，難以進行日光返照之改 正。因此先嘗試從影像色彩不均及資訊損失之補償方法著手，使用白平衡及去霧處理兩 種影像處理方法進行測試，測試結果發現日光返照區域有過度增顯情形，範圍增大且亮 度誇張化。但若將日光返照區域視為亮度特徵，把影像RGB 色彩空間轉為 HSV 色彩空 間後，以亮度進行二維經驗模態分解，發現原始影像減去模擬影像之殘餘項與日光返照 之範圍與特徵相當，將其視為雜訊後去除後得到之影像，可以有效改善了日光返照之影 響。 使用平均值、均方根誤差、峰值訊號雜訊比及結構相似性指標對於影像結果進行統 計指標計算，並與肉眼觀察結果進行比較，RGB 波段平均值主要與原始影像色調有關， SSIM 指標值可以看出色彩結構的維持，但主要以人眼判識為主。 關鍵字：日光返照、雙向反射分布函數、MODO 軟體、白平衡、去霧處理、HSV 色彩 空間、二維經驗模態分解、峰值訊號雜訊比、結構相似性指標。

II

Sun Glint Removal from Water Area Image using Bidirectional

Empirical Mode Decomposition

Student：Ling-Chien Shih Advisor：

Tian-Yuan Shih

Department of Civil Engineering

National Chiao Tung University

Abstract

When taking aerial photography images on sea area at daytime, the sea surface reflection of solar radiation cause the sun glint. It will cause the loss of image information, even continuous images stitching. In this paper bring up the definition of sun glint, and use the radiation transfer analysis and image processing methods for testing.

We readed images to get digital numbers, but we can’t transfer it to radiance.Even we use Cox and Munk (1954) Model and MODO to get the transmittance and solar radiance, we still can’t correct the radiance.There was no models of sea surface when shooting images,the transmittance we calculated can’t apply for every point of surface, it can’t be used to correct sun glint. For the uneven image colors and the loss of image information, using white balancing and dehazing for testing, the results found sun glint region increased significantly, increase the range and brightness exaggerated. If sun glint area can be regarded as the brightness characteristics and allowed images from RGB color space into HSV color space, the Value processed with Bidimensional Empirical Mode Decomposition(BEMD), and found that the region and characteristics of residual with sun glint of origin images are the same, and removed the residual. That results can found that BEMD really improve the sun glint.

III

The Mean value, Root Mean Square Error (RMSE) , Peak Signal to Noise Ratio(PSNR), and Structural Similarity (SSIM), which are the statistical indicators calculated for results, and these indicators compared with the naked eyes, and finally make conclusions and suggestions. Keywords: sun glint, BRDF, MODO, white balancing, dehazing, HSV color space, BEMD, PSNR, SSIM.

IV

致謝

轉眼間兩年的碩士生活結束了，曾經在修課的時候遇到很多困難，很多不足都仰賴 老師、學長姊們還有同學們的支持度過，從大學時對於測量只是選修幾門基礎課，這兩 年來對於測量知識一點一滴的累積，都讓我覺得很幸運能成為測量的一份子。 在此感謝兩年的指導教授史天元 博士，老師即使在非常忙碌的時候，依然會抽出 時間與我討論研究的問題，對於我專業知識的不足與論文寫作的困難與錯誤給予指導， 總是不厭其煩，很感謝老師在學生身上所花費的心力，學生在畢業後會更加努力，老師 也要記得保重身體。 在此亦感謝口試委員曾義星 博士、蔡榮得 博士與張智安 博士對於學生論文的指 導與建議，感謝曾老師特別趕來新竹，並且主持學生的口試；感謝蔡老師對於學生論文 的用心，一字一句的建議，想必花了老師不少時間；感謝張老師在兩年期間，指導學生 寫程式所遇到的困難，提供學生軟體的協助，給予許多寶貴的建議，尤其學生在程式的 最後一部遇到關卡時，老師給的建議讓學生卡了幾個禮拜的問題迎刃而解。 兩年中小悶學長、雅信學姊提供的幫助無數，感謝你們容忍我對於軟體與寫作的困 難，不厭其煩的教導我；感謝薇帆、昭儀、郁珊時時刻刻的陪伴與支持，感謝柏溶的幫 助與關心，也感謝世涵、邦和、詠升，兩年的跨年有你們的陪伴好熱鬧，八個人的歌唱 與嬉戲，是我在兩年研究生活中可以繼續的動力，我已經開始期待未來的同學會。也感 謝尚儒和張朔在我忙碌的時候提供協助，以及幫助我口試的準備與紀錄。 感謝我最重要的家人，奶奶、外公、爸媽還有兩個哥哥，還有已經在天堂的外婆。 外婆還在的時候總是不停的跟我說要好好讀書，也總是跟媽媽關心我的近況；爸媽非常 希望我可以完成研究所的學業，但又努力的盡量不給我壓力，謝謝您們相信我；感謝兩 位哥哥在與我談心之餘，告訴我支持我的決定。我親愛的家人們，我很愛你們。 最後感謝許許多多幫助我關心我的長輩與朋友們，有你們的支持才能度過這些年， 祝福每一位未來萬事如意，健康順心。

V

目錄

摘要 ... I  ABSTRACT ... II  致謝 ... IV  目錄 ... V  圖目錄 ... VII  表目錄 ... X  1. 前言 ... 1  1.1 研究背景 ... 1  1.2 研究動機與目的 ... 1  1.3 研究流程 ... 3  1.4 論文架構 ... 4  2. 研究方法 ... 5  2.1  輻射傳輸解析 ... 5  2.1.1 MODTRAN4 模組 ... 5  2.1.2 Cox and Munk (1954)  模型 ... 9  2.1.3 雙向反射分布函數 ... 10  2.2 影像處理 ... 11  2.2.1  白平衡基本原理 ... 12  2.2.2  去霧處理基本原理 ... 13  2.2.3 經驗模態分解原理 ... 16  2.2.3.1 一維經驗模態分解 ... 16  2.2.3.2 二維經驗模態分解 ... 18 

VI 2.2.3.3 實驗過程 ... 18  2.2.3.4 散亂極值點形成包絡面方法‐Delaunay 三角網 ... 20  2.2.3.2  取得 IMF 之過程 ... 21  2.2.3.4 HSV 色彩空間 ... 25  2.3 利用近紅外波段資訊進行校正 ... 26  2.4 影像品質指標 ... 30  3. 研究資料 ... 33  4. 成果與討論 ... 35  4.1 輻射傳輸解析‐MODO ... 35  4.2 白平衡&去霧處理成果 ... 38  4.3 二維經驗模態分解(BEMD)影像成果 ... 43  4.3.1 第一組測試影像 ... 43  4.3.2 第二組測試影像 ... 58  4.4 影像品質分析 ... 63  5. 結論與建議 ... 65  參考文獻 ... 66 

VII

圖目錄

圖 1.1、東沙環礁航空影像拼接(陳雅信，2010) ... 2  圖 2.1、太陽輻射散射路徑圖(李暢游等，2011) ... 6  圖 2.2、MODO 流程圖(MODO, 2010) ... 7  圖 2.3、一般 MODTRAN4 之參數設定檔(tape5) ... 8  圖 2.4、MODO 軟體參數設定介面 ... 9 

圖 2.5、反射幾何示意圖(Cox and Munk, 1954) ... 10 

圖 2.6、BRDF 示意圖(Centennial, 2001) ... 11  圖 2.7、太陽輻射消散示意圖 ... 14  圖 2.8、去霧處理流程圖 ... 15  圖 2.9、去霧處理前後之影像 ... 16  圖 2.10、EMD 架構流程圖 ... 17  圖 2.11、BEMD 分解流程圖 ... 19  圖 2.12、(a)紅點為眾圓之圓心，黑線為圓上任意三點連成之 Delaunay 三角網。(b)將 (a)之紅點連起來為 Voronoi 圖形。(Wikipedia, 2012) ... 20 

圖 2.13、Delaunay 三角網示意圖(胡珺婷, 2012) ... 21  圖 2.14、EMD 分解示意圖(Huang,1998) ... 23  圖 2.15、待分解的原始訊號 ... 24  圖 2.16、尋找局部極大值點 ... 24  圖 2.17、插值得到局部極大值包絡線 ... 24  圖 2.18、插值得到局部極小值包絡線 ... 24  圖 2.19、將局部極大值與極小值包絡線平均後之均值包絡線 ... 25  圖 2.20、原始訊號減去均值後之餘項 ... 25  圖 2.21、HSV 色彩空間示意圖 ... 26  圖 2.22、可見光波段與近紅外波段亮度之線性關係(Hedley, 2005) ... 27 

VIII

圖 2.23、永興島 Worldview2 原始衛星影像(銳悌科技，2010) ... 28 

圖 2.24、(a)、(c)為原始影像圈選之測試區域，(b)、(d)為處理後之影像結果 ... 29 

圖 3.1、IDE UI-2250SE 數位相機(IDE, 2010) ... 33 

圖 4.1、太陽(輻射光源)位置之計算 ... 36  圖 4.2、光譜值之萃取 ... 36  圖 4.3、太陽表面輻射值與傳輸量輻射值之圖示 ... 37  圖 4.4、測試影像 A 原始影像 ... 39  圖 4.5、測試影像 A 去霧處理後 ... 39  圖 4.6、測試影像 A 白平衡後 ... 40  圖 4.7、測試影像 B 原始影像 ... 40  圖 4.8、測試影像 B 去霧處理後 ... 41  圖 4.9、測試影像 B 白平衡後 ... 41  圖 4.10、測試影像 C 原始影像 ... 42  圖 4.11、測試影像 C 去霧處理後 ... 42  圖 4.12、測試影像 C 白平衡後 ... 43  圖 4.13、測試影像 A 原始影像 ... 44  圖 4.14、測試影像 A 之 IMF1 ... 44  圖 4.15、測試影像 A 之 IMF2 ... 45  圖 4.16、測試影像 A 之 IMF3 ... 45  圖 4.17、測試影像 A 之殘餘項 ... 46  圖 4.18、測試影像 A 影像結果(HSV 色彩空間) ... 46  圖 4.19、測試影像 A 二維經驗模態分解結果 ... 47  圖 4.20、測試影像 B 原始影像 ... 48  圖 4.21、測試影像 B 之 IMF1 ... 48  圖 4.22、測試影像 B 之 IMF2 ... 49 

IX 圖 4.23、測試影像 B 之 IMF3 ... 49  圖 4.24 測試影像 B 之殘餘項 ... 50  圖 4.25、測試影像 B 影像結果(HSV 色彩空間) ... 50  圖 4.26、測試影像 B 二維經驗模態分解結果 ... 51  圖 4.27、測試影像 C 原始影像 ... 52  圖 4.28、測試影像 C 之 IMF1 ... 52  圖 4.29、測試影像 C 之 IMF2 ... 53  圖 4.30、測試影像 C 之 IMF3 ... 53  圖 4.31、測試影像 C 之殘餘項 ... 54  圖 4.32、測試影像 C 影像結果(HSV 色彩空間) ... 54  圖 4.33、測試影像 C 二維經驗模態分解結果 ... 55  圖 4.34、測試影像 D 原始影像 ... 56  圖 4.35、測試影像 D 二維經驗模態分解影像結果 ... 57  圖 4.36、測試影像 E 原始影像 ... 59  圖 4.37、測試影像 E 二維經驗模態分解結果 ... 59  圖 4.38、測試影像 F 原始影像 ... 60  圖 4.39、測試影像 F 二維經驗模態分解結果 ... 60  圖 4.40、15 張原始影像鑲嵌結果 ... 61  圖 4.41、15 張測試影像二維經驗模態分解結果之鑲嵌 ... 62 

X

表目錄

表 3.1、東沙地區航線規劃參數 ... 33  表 3.2、相機規格 ... 33  表 3.3、實驗影像資訊 ... 34  表 3.4 飛行及相機相關參數 ... 34  表 3.5、第二組實驗影像資訊 ... 34  表 4.1、影像 A 指標值 ... 64  表 4.2、影像 B 指標值 ... 64  表 4.3、影像 C 指標值 ... 64  表 4.4、影像 D 指標值 ... 64 

1

1. 前言

日光返照(Sun Glint)現象發生在太陽與海表面之入射角及感測器與海表面之反射角 相同時， 受影響區域若為平滑海面會產生鏡面反射，反之，若為粗糙海面則會形成陰 暗色調(Wikipedia, 2012c)。本章簡述研究背景與目的，其後提出解決日光返照問題使用 之方法與背景，對於這些方法使用相同影像進行測試，並在測試後進行成果討論。

1.1 研究背景

當太陽的仰角越大，日光返照之影響區域會越集中，灰度值強度也會更高，對於影 像的資訊損失會更加明顯，在日光返照的去除後，資訊的損失常常無法彌補。 航遙測影像在影像拍攝時，因為海面的波動或是拍攝時間接近正午，都有可能造成 日光返照而使海面下之資訊不明或是顏色亮度過高之情形，因此本研究希望可以藉由能 量傳輸解析以及一般影像處理進行日光返照之改正。

1.2 研究動機與目的

在執行拍攝任務時，拍攝時間雖然可以盡量避免正午或天候不佳，但每日太陽東昇 西降之時間略有差異，甚至拍攝日數與預算之限制等等，在拍攝時仍有無法避免之日光 返照情形，因為返照之影響會造成影像本身資訊損失，光影不均之情況也有造成人為判 斷錯誤的可能性。 受限於像幅大小，拍攝一整個區域之影像必須進行鑲嵌才可得到完整地區影像。將 一批有一定重疊之影像，拼接後去除邊緣部分，組成大範圍之鑲嵌圖，用來鑲嵌之影像 一般具有一定程度的幾何一致性及色調連續性。航空影像即使為連續拍攝，重疊區域仍 然可能因為拍攝角度不同而有影像色調之差異，如圖1.1 所示，當拍攝角度剛好符合表 面的反射路徑時，即會產生日光返照，在拼接時會造成拼接處的光影，破壞影像的連續 性及資訊的遮蔽。

3

1.3 研究流程

一般處理日光返照之文獻或相關工具程式，大多會參考Cox and Munk (1954) 提出 之模型，此模型建立光源之入射與反射的幾何情形，並且考慮大氣或風浪等影響，由影 像所量測之輻射量去除這些影響後，以模型模擬拍攝當時之表面輻射量。使用該模型必 須先得到反射表面之雙向反射分布函數(Bidirectional Reflectance Distribution Function, BRDF)，BRDF 可以由率定後之相機與光源拍攝真實物件而直接測得。

Rosa et al.(2006)使用 Cox and Munk(1954)模型進行 MIVIS 影像之日光返照(Sun Glint)改正，其中太陽輻射量及輻射傳輸量以 MODTRAN(MODerate resolution

atmospheric TRANsmission)進行模擬，因為一般航空拍攝影像為可見光波段(波長範圍約 400~700nm)，因此使用 MODTRAN 進行測試。為了方便使用者輸入參數後來出現了 MODO (MODTRAN Operation)等軟體，本研究測試是否可以使用 MODO 軟體計算之傳 輸量用以改正日光返照。

若將日光返照視為雜訊，參考Huang et al.(1998)提出的經驗模態分解(Empirical Mode Decomposition, EMD)及 Damerval et al. (2005)提出的二維經驗模態分解

(Bidimensional Empirical Mode Decomposition, BEMD)，分解出內建模態函數(Intrinsic Mode Functions, IMF)後，其中一個模態函數會與日光返照有關，將其去除後即可得到改 正之影像。

除了二維經驗模態分解外，使用Tang et al.(2009)提出之單張影像去霧方法及白平衡 這兩種常用於影像改正之處理方法，並且提及處理有近紅外波段資訊影像之日光返照處 理。但本篇論文主要以處理航空拍攝彩色影像之日光返照為目的，並無近紅外波段資訊， 因此在這只做分享不做評論。

4

1.4 論文架構

第一章介紹日光返照之形成背景與研究目的，並提出為了解決日光返照之影響所實 驗之方法種類；第二章對於能量傳輸解析、影像處理及近紅外波段校正三類方法進行文 獻回顧及研究方法介紹，並分享初步的實驗結果；第三章對於實驗影像及使用相機進行 規格介紹；第四章提出實驗所遇到之困難以及成果展示，對於成果進行人眼判識及指標 分析；第五章提出結論與建議。

5

2. 研究方法

遙感探測所拍攝之影像，無論是無論由衛星或航空器，在拍攝時有可能會因為天候、 雲層或拍攝時間接近正午造成影像資訊損失，其中拍攝海域範圍之影像最常見的情況之 一，即為日光返照。 去除日光返照的方法大致可分為輻射傳輸解析、影像處理方法及利用近紅外波段資 訊校正，三種類型。

2.1 輻射傳輸解析

利用輸入環境相關參數及太陽入射角等已知參數與影像輻射值進行計算，可模擬拍 攝當時之輻射值，本研究所使用之影像拍攝區域為海域，表面崎嶇程度隨著風浪變化， 無法得知拍攝當時之反射表面模型，與輻射值無法對應。

2.1.1 MODTRAN4 模組

根據氣象衛星的觀測結果，除了地物會將能量反射至感測器外，亦包含地物環境散 射與反射以及大氣的散射與反射等。輻射傳輸方程式描述輻射傳播通過介質時在介質間 發生交互作用(吸收、散射、反射等)，而使輻射能夠按照一定的規律傳輸，在地球大氣 條件下之求解非常複雜，只能在一些假定的參數下求解。若Iλ(s)為輻射傳輸量，輻射傳 輸方程式定義如下： I_λ s I e λ _(2.1) 其中大氣消散係數K_λ s λ λ ，表示輻射在大氣中傳播距離ds 時的相對衰減率； 指數項e λ _{為透射率；dIλ為輻射經過距離ds 後的強度變化值。} 如圖2.1 所示，太陽光入射時，感測器所接收之光量由海表面反射光經過散射及水 底反射光經過水中散射所組成，若能由散射後之接收光量計算原始反射量，則可用以估 計日光返照之影響。

7

Reflectance

MODO

Input preparation Tape5

MODTRAN (fortran) Tape7 Solar flux MODO Output evaluation Sensor response Radiance L Convolved L Plots 圖 2.2、MODO 流程圖(ReSe, 2010) MODO 主要有四個限制： (一) MODO 為專業的模擬工具，因此使用者需要有輻射轉換的相關知識。 (二) MODO 暫不支援 MODTRAN4 的 BRDF 功能。 (三) 查詢的輸入檔 tape5 較不適合模擬多維度表面。 (四) 輸入檔內的參數對於懸浮粒子及標準無線電探測剖面沒有有效的定義。 輸入之參數檔tape5 依當時環境選擇參數值，選項如拍攝地區緯度、拍攝季節、拍攝 環境當時氣體含量等等，但是多樣參數是變動值無法確切得之，因此可選擇不填取。一 般在使用MODTRAN4 前必須先製作參數檔 tape5，以代號代替所設定之參數，參數與 參數間須以空格隔開，每個位置代表著特定的參數，在眾多參數的情況下極有可能產生 錯亂之情形如圖2.3。

12

若將日光返照視為雜訊，使用二維經驗模態分解(Bidimensional Empirical Mode Decomposition, BEMD)將影像分成多個內建模態函數，日光返照可能為某個模態函數， 將其扣除後可以得到改正後之影像。

2.2.1 白平衡基本原理

色溫為加熱黑體(Black Body)在該絕對溫度下所顯現之色光，影像拍攝時因為不同 光波所發射的能量不同，會感知不同的色彩變化。當較高色溫與較低色溫之光源照射在 白色物體上，高色溫光源所照射之物體顯色偏藍，低色溫光源所照射之物體顯色偏紅， 因此色溫會影響拍攝目標之色偏情況。(Wikipedia, 2012d) 灰階理論為白平衡處理方法之一，對於色彩均勻之影像處理效果極佳，本實驗使用 色彩均勻之海域影像，因此選擇灰階理論為研究方法。 進行白平衡轉換時需要調整個色彩之乘積值，即為色彩調整係數。「適應白點」為 人眼視覺系統於色溫下感知之恆常性白點，而「採用白點」為相機所拍攝之白點資訊， 在此白平衡轉換利用Von Kries 色適應模式使相機拍攝之採用白點轉換為人眼所感知的 適應白點，即為模擬人眼所見之真實色彩(陳炳旭，2012)。由於適應白點為人所感知， 無法準確量化，一般以純白點即(R,G,B)=(255,255,255)代表適應白點，可得以下方程式： R 255 R , R G 255 G _, G B 255 B , B (2.7) 其中R _, 、G _, 、B _, 為轉換前之白點色彩訊號；R 、G 、B 為轉換前影像各灰度值 色彩訊號；R,G,B 為轉換後影像各灰度值色彩訊號； , 、 , 、 , 為白平衡轉換之「調 整係數」。

13 Chiou (2001)提到灰階理論即為假設人眼所感知的自然界影像色彩成份分布均勻， 三個色彩之強度相同。依此理論，若想在色調偏差下經由色適應模式模擬人眼所感知之 色彩，需將RGB 平均強度調整至相同大小，模式可表示如下： R R K R G G K G B B K B (2.8) 其中K 為影像調整參數；R、G、B 為轉換前各灰度值之強度； R 、G 、B 為整張影像之平均強度；R 、G 、B 為轉換後各灰度值之強 度。

2.2.2 去霧處理基本原理

在室外拍攝影像時，有時因為氣候不佳或大氣的吸收與光源的散射會造成影像有霧 化的情形，造成影像對比度與色彩辨識度有失真的情形，此時會造成資訊的損失，因此 必須進行去霧處理。 霧化影像拍攝當下，相機所接收之輻射值比一般情況低，起因於大氣光量被環境散 射光量所替代，此時的環境散射光量即為霧光(airlight)。霧光的多寡端視於景深，因此 去霧處理的程序與景深有密切的相關，目前廣泛用來描述霧化影像資訊的方程式如下：

15

戶外無霧化影響之影像，非天空範圍的區域會有一些像元在RGB 三個波段至少會有一 個波段之輻射量擁有低強度，這些波段就稱為Dark channel。

首先利用Dark channel prior 方法估計大氣光量，藉由無霧化影像之 J 值趨近於 0 之 假設推導傳輸量圖，最後使用soft matting 方法優化影像。Soft matting 之程序主要為先 限制一小區域，由已知鄰近灰度值或反覆執行非線性α檢定推估前色與背景色，當 matting 方程式導出之α示意圖與霧化模型推導之傳輸量圖相似時，藉由適當的調整可 用來優化傳輸量圖。 其中softmatting 之程序主要為截取影像之一部分，由已知前景深與背景深的像素推 估αmap，以α值替代傳輸量 t(x)。 I(x)=J(x)t(x)+A(1-t(x)) → I(x)=Fα+B(1-α) (2.11) 輸入霧化 影像 Dark Channel Prior 估計大氣光量A 估計傳輸量t(x) Softmatting 優化 去霧影像 圖 2.8、去霧處理流程圖 施伶蒨(2011)使用農林航空測量所於民國 94 年拍攝之東沙環礁數值影像進行測試， 如圖2.9，可以發現在去霧處理後，色彩增顯及視覺效果有所改善，因此選用去霧處理 對於影像日光返照之改正進行測試。

17

  圖 2.10、EMD 架構流程圖

Huang et al.(1998)認為 IMFs 須滿足兩個條件：

(1)函數在整個範圍內，每個 IMF 局部的極值點及跨零點數目必須相同或相差一個。 (2)在任意時刻點，局部最大值包絡線(上包絡線)與局部最小值(下包絡線)之平均值需為 0， 並且每個IMF 對稱於局部平均線。

18

2.2.3.2 二維經驗模態分解

EMD 基本原理主要適用於與時間相關之函數，如潮汐之變化及音調之分析等等， 對於處理二維影像時可以保有EMD 能處理非線性及非穩態資料之特性，以二維影像之 行數(或 pixel 數)替代原來的時間作為橫軸。 BEMD 之假設類似於 EMD 如下： (1) 二維平面至少包含一個極大值和一個極小值，或是整個平面沒有極值點但是可以進 行一階或多階求導數運算後得到極值點。 (2) 特徵尺度由極值間的距離差值定義。 Pei (2007)將餘數 rn視為光場分布趨勢，藉由線性轉換將原始影像光影移除，如下式 表示： Gain(x,y)=I(x,y)－rn(x,y) I’(x,y)=Gain(x,y)*255 (2.12) 其中=I(x,y)為原始影像；rn(x,y)為光場分布影像；I’(x,y)為移除光影後之影像。

2.2.3.3 實驗過程

如圖2.11 所示，將原始影像 RGB 空間轉為 HSV 色彩空間，且投影到 XY 平面，取 出影像對應灰度值點的亮度(Value)為 Z 座標值。找出影像亮度的局部極大值與局部極小 值，對找出之散亂極值點進行平面 Delaunay 三角剖分，曲面擬合後得到局部極大值包 絡線與局部極小值包絡線，並求其均值包絡線。 將原始影像亮度減去均值包絡線，檢查是否滿足篩選程序結束條件，若滿足則結束 程序得到IMFs，檢查殘餘項是否滿足影像分解結束條件(殘餘項已為常數無法再分解)， 若滿足則輸出各層結果並結束，若不滿足則繼續篩選過程。由IMFs 或殘餘項可以找到 相應之特徵或趨勢，扣除該趨勢項即可得到改正影像。

19

22 第二次的篩選程序將h1視為原始訊號，重複執行篩選程序k 次後，直到 h1k為IMF 為止。因為篩選過程符合了IMFs 的定義要求，h1理論上應為IMF， 但是實際上原始訊 號常常會有突然遽變的情況，在執行篩選程序後可能會有新的極值的產生，甚至於先前 的極值有過度誇張的情況。在立方曲線插值後有些平均值可能為小數值，無法每次符合 前述IMFs 符合條件(在任意時刻點，局部最大值包絡線與局部最小值之平均值需為 0)， 需對其定義一個誤差範圍。 篩選過程中可以平滑不平整的振幅，但是重複篩選程序到達某一程度時可能會除去 部分具有物理意義的振幅擾動，必須決定讓篩選程序停止的收斂條件，確保IMFs 之振 幅與頻率變動依然符合物理意義。收斂條件簡述為： (1) 跨零點總數與局部極值點總數必須相等。 (2) 標準偏差值(Standard Deviation, SD)設定在 0.2~0.3 之間，

SD

∑

(2.15) 當滿足兩個收斂條件之一則停止篩選程序得到第一個IMF，則停止此篩選程序，繼續進 行下一個IMF 的篩選程序。若 c1為訊號中最好的時間尺度或最短周期的分量，r1可以 從原始訊號中分解出，如下式表示： X(t)-c1=r1 (2.16) 將r1視為新的訊號繼續進行篩選程序， r1-c2=r2 r2-c3=r3   rn-1-cn=rn (2.17) 篩選過程中利用收斂條件停止動作，當rn變成常數或單調函數無法再分解出IMF 就會停止篩選程序，最後的殘餘項會是不對稱的訊號值，也就是訊號的趨勢(trend)。因

30

2.4 影像品質指標

主觀評估以人眼感官判定影像處理後結果之優劣，客觀評估則使用統計方法對於影 像結果之數據進行量化。此篇研究主要以探討日光返照為主，因此只對影像亮度值進行 統計計算。統計指標值之高低不能代表影像處理結果之優劣，只能對於人眼感知外提供 參考，若人眼無法判別兩種方法影像之結果差異，可以參考統計指標值高低比較優劣。 本研究主要討論影像亮度變化，主要希望在處理亮度不均情況時仍然可以保留原始 影像之色彩或結構，以防止失真，因此採用常見比較影像品質優劣的統計指標，以下列 舉均方根誤差(RMSE)、峰值訊號雜訊比(PSNR)及結構相似性指標(SSIM)，不使用標準 化互相關法(NCC)之原因為測試影像為灰度值較為均調之海域，相關係數最大值難以確 認。 (一)平均值(Mean) 利用平均值可以得知在影像處理前後亮度變化情形，其平均值越大表示影像整體亮 度提高，反之，其平均值越小表示整體亮度降低。將影像R、G、B 波段分別取平均， 可用來比較其色偏情形。

μ

∑ ∑ ∑ , (2.21) 其中M,N 為影像行列數；IK(i,j)為影像第 k 波段在位置(i,j)的亮度值。

(二)均方根誤差(Root-Mean-Square Error, RMSE)

均方根誤差(RMSE)常用來量測精度，可以比較原始影像與處理後影像之差異。 RMSE √MSE ∑ ∑ I i, j K i, j (2.22) (三)峰值訊號雜訊比(Peak Signal-to-Noise Ratio, PSNR)

峰值訊號雜訊比(PSNR)一般使用在影像壓縮前後的品質比較，亦常用在檢查影像處 理後品質，確認該處理程序之適當程度。峰值訊號雜訊比之計算方式如下：

MSE ∑ ∑ I i, j K i, j (2.23) PSNR 10 ∙ log (2.24)

31 其中MSE 為兩張影像 Imxn、Kmxn間的均方誤差，MAXI為灰度值最大值(一般設為 255)， PSNR 值之單位為 dB。 PSNR 為目前廣泛使用於比較影像品質的統計指標，影像之 PSNR 值大約在 30~50dB 之間，若在30dB 以下表示為人眼無法接受，一般 PSNR 值越高表示影像品質越好，但 是PSNR 值無法與人眼感知結果一致，有可能 PSNR 較高者在人眼感知反而品質比 PSNR 較低者差，因為人眼感官對於誤差的敏感度受許多因素影響而非絕對。

(四)結構相似性指標(Structural similarity, SSIM)

結構相似性指標一般用來衡量兩張影像之相似程度，當以一張影像為無失真影像， 另一張是否為失真影像，則可由結構相似性判斷失真程度，相較於傳統影像品質指標， 例如峰值訊號雜訊比(PSNR)，更能符合人眼對於影像品質之判斷。 Wang(2004)提出之結構相似性指標之性質主要如下： (1)指標值有對稱情形，亦即 SSIM(x,y)=SSIM(y,x)。 (2)指標值有界限，-1≦SSIM≦1；當結構完全相同時，SSIM=1。 由影像結構來看，影像如同由亮度、對比度及結構三個部分組成，結構相似性指標 之計算以平均值(Mean)估計亮度，以標準差(STD)估計對比度，以共變異相關係數估計 結構相似度。計算步驟如下： Step1： 若x 為原始影像，y 為處理後之影像，平均值μx、μy估計了原始及處理後影像之 亮度值。

μ

∑

，μ

∑ (2.25) 其中N 為影像灰度值數目；xi為影像灰度值值；yi為處理後影像灰度值值。 結構相似性指標中亮度估計值L(x,y)計算方式如下：

L x, y

_{μ μ}μ μ (2.26) 其中常數C1之作用為防止μ μ 之值趨近於零，定義 C1=(K1L)2，K1=0.01，L=255。

32 Step2： 結構相似性指標中對比度由標準差σx、σy估計，對比度估計值C(x,y)計算方式如 下： σ ∑ x μ N 1 ，σ ∑ y μ N 1

C x, y

_σσ σ_σ (2.27) 其中常數C2之作用為防止σ σ 之值趨近於零，定義 C2=(K2L)2，K2=0.03，L=255。 Step3： 結構相似性指標中結構相似度由兩張影像之共變異相關係數Cov(x,y)估計，結構相 似性量化後之估計值S(x,y)計算方式如下： Cov x, y σ ∑ x μ y μ N 1

S x, y

, (2.28) 其中常數C3之作用為防止σxσy之值趨近於零，定義C3=C2/2。 Step4：

由亮度估計值L(x,y)、對比度估計值 C(x,y)及結構相似性估計值 S(x,y)將結構相似 性指標(SSIM)定義為：

SSIM x, y L x, y C x, y S x, y (2.29) 當α=β=γ=1 時，上式 SSIM(x,y)值可以表示為：

34 表 3.3、實驗影像資訊 影像 航帶 編號 像素 A 32 36924340165 1600×1200 B 31 32432020013 1600×1200 C 19 62096010166 1600×1200 D 19 62096010164 62096010165 62096010166 1823×1349 第二組測試影像為農航所於民國94 年拍攝之東沙環礁航空影像掃描數值檔，以 RMKTOP15 型 144116 號相機拍攝，影像編號 94R001_16~20、33~37、68~72 ，每張 16858X16806 像素，共計 15 張。 表 3.4 飛行及相機相關參數 相關參數 RMKTOP15型144116號相機 解析度 16858×16806 焦距 152.818mm 平均航高 2370m 像比例尺 1/15500 像元對應地面解析度 21.7cm 表 3.5、第二組實驗影像資訊 影像 編號 像素 E 94R001_ 16858×16806 F 94R001_ 16858×16806 G 94R001_16~20 94R001_33~37 94R001_68~72 35443×43323

35

4. 成果與討論

本章主要就能量傳輸解析及影像處理方法之結果進行展示與討論。MODO 軟體可 以計算太陽輻射量及傳輸量，但沒有雙向反射分布函數無法進一步模擬反射表面輻射量； 白平衡與去霧處理由影像成果可看出對於日光返照改正之效果不佳；二維經驗模態分解 在去除日光返照之雜訊後結果，可以看出日光返照的明顯改善。

4.1 輻射傳輸解析-MODO

MODO 之輸出有三： (一) 太陽(輻射光源)之位置：如圖 4.1，輸入拍攝日期、UTC 標準時間以及拍攝地點經 緯度，即可輸出太陽方位角、太陽仰角及太陽天頂角。 (二) 太陽表面輻射值或傳輸量：如圖 4.2，輸入波長範圍及解析度，選擇輻射值或傳輸 量輸出，即可得到太陽表面輻射值或傳輸量。 (三) 將(二)之輸出量進行繪圖，得到如圖 4.3 之結果，橫坐標為波數(波長之倒數，單位 為cm-1，表示每公分含有的波數)，縱座標為輻射量，黑色與藍色圖示分別為太陽輻 射量及輻射傳輸量。 可見光波段之波數計算後為1.4×104~2.5×104(cm-1)，由圖十可看到可見光範圍內，輻 射消散量(黑色與藍色曲線之差值)接近一個常數，因為在 tape5 介面輸入之參數相同， 因此只能得到單一數值。此時如果能計算影像拍攝區域表面之雙向反射分布函數(BRDF)， 可由表面反射之狀況及輻射消散常數進行表面輻射量之計算，但實驗影像無法求得 BRDF，即無法得知感測器接受輻射增量與表面反射輻射增量之比值。 除此之外，沒有參數可以得知加常數及乘常數將輻射量換算為影像灰度值，即使得 到表面輻射量，若無法換算為影像灰度值，即無法直接對於影像進行日光返照之改正。

38

4.2 白平衡&去霧處理成果

使用三張影像A、B、C 進行去霧處理及白平衡，影像 A 原始影像及結果如圖 4.4~4.6， 影像B 原始影像及結果如圖 4.7~4.9，影像 C 原始影像及結果如圖 4.10~4.12。 影像A 為海域，海面下有不明顯之礁塊，去霧處理後影像(圖 4.5)色彩在強度較強 區域過於增顯，呈現整張影像色彩失真之情況，右半部有嚴重的資訊損失，但是隱藏之 礁塊範圍比較可以分辨；影像A 白平衡後影像(圖 4.6)色調被均勻化，可以看出右下角 區選區域日光返照區域過於強化，海底下之原本不明顯之礁塊範圍比較可以分辨，但是 無論是去霧處理或是白平衡只能增顯影像之資訊，對於影像之日光返照無法處理。 影像B 為較單純海域，因為風浪的走向而使日光返照呈現明顯的波紋，去霧處理後 影像(圖 4.8)因為過度誇張化，整張影像偏亮，有嚴重色偏；白平衡後影像(圖 4.9)將色 調均化後，右上角之色彩資訊較為明顯，但也看出了日光返照區域反而有強化的反效 果。 影像C 有明顯之礁塊，去霧處理後影像(圖 4.11)，整張影像過度增顯，亮度過亮以 致整張影像失真；白平衡後影像(圖 4.12)色彩均勻化，但日光返照區域沒有明顯改善。 由影像A、B、C 之測試，可以發現去霧及白平衡處理對於海域影像直接進行處理 是不可行的，甚至於去霧處理後整張影像過亮，可能是整張影像範圍皆有日光返照之影 響，因此在增顯後過度誇張。

58

4.3.2 第二組測試影像

圖4.36 為測試影像 E 原始影像，涵蓋區域為東沙島環礁的一部分，經過二維經驗 模態分解扣除殘餘項後，得到如圖4.37 之結果，右下角日光返照區域與附近色調相近， 改善了亮度過亮的情況，但是殘餘項內可能含有部分環礁紋理，扣除後雖然改善日光返 照影響，卻也讓紋理辨識程度下降。 圖4.38 為測試影像 F 原始影像，涵蓋整個島嶼，經過二維經驗模態分解扣除殘餘 項後，得到如圖4.39 之結果，影像處理後日光返照區域與周圍顏色較為接近，但是島嶼 邊界以及海域特徵都有模糊之情形。 圖4.40 為 15 張原始影像鑲嵌結果，在日光返照區域有光暈情形，使得接合處有不 連續的情況，經過二維經驗模態分解扣除殘餘項後(圖 4.41)，色彩較為均化，但日光返 照區域無法完全消除，並且每張影像單獨進行處理時，因為海底底質不同造成單張影像 主要組成色彩不同，在影像處理過後顏色略有深淺之差異。 比起第一組影像，因為影像解析度較差，經過二維經驗模態分解之訊號無法單獨分 出影響日光返照之亮度，因此在移除殘餘項後模糊情況較第一組嚴重，顯示使用此方法 進行亮度分解時會受到影像解析度之限制。

63

4.4 影像品質分析

除了人眼判識外，以指標值對於去霧處理、白平衡及 BEMD 之結果進行指標品質 分析，分析對象為影像灰度值，主要分析項目有平均值(Mean)、標準差(STD)、均方根 誤差(RMSE)、標準化互相關法(NCC)、峰值訊號雜訊比(PSNR)及結構相似性指標，由統 計指標值探討影像處理前後之灰度值變化情形。 均調區域影像A、B 經過白平衡處理後，R、G、B 平均值皆增加；礁塊區域影像 C 原始色彩為藍綠色調為主，因此在白平衡後G、B 值降低，而 R 值增加。由 RMSE、PSNR 及SSIM 值皆可看出影像 C 因為有礁塊特徵，因此結構與原始影像較為接近。加上人眼 判識後，白平衡處理雖然可以使影像色調均勻，但無法改正日光返照之影響。 影像A、B、C 經過去霧處理後，對於日光返照區域過度增顯，因此 R、G、B 平均 值皆誇張化，由PSNR 與 SSIM 值亦可看出其視覺效果與結構已被破壞，可知整張影像 進行去霧處理產生的增顯效果對於日光返照影響無法去除。 影像A 及 B 為日光返照影響區域較小之影像，經過 BEMD 處理後，因為返照區外 亮度增加，因此整體灰度平均值反而有增加的情形；影像 C 及 D 為日光返照影響幾乎 整個範圍之影像，在改正日光返照區域之亮度後，灰度平均值皆有降低之情形。以均方 根誤差(RMSE)來看，測試影像 D 與原始影像之灰度值偏差較小，測試影像 C 次之，與 肉眼觀察結果符合；測試影像B 之灰度值偏差最高，以肉眼觀察影像色偏情形，的確有 一定的色偏情形。 BEMD 處理後之影像，PSNR 值皆 40 以上，較無法提供視覺判斷；以結構相似性 指標(SSIM)來看，測試影像 C 與 D 保留原始影像灰度值結構的情況較好，測試影像 A 與B 在影像灰度值結構有較多的差異，可能會有色彩失真，與肉眼觀察結果相符。但是 礁塊區域之細節在 BEMD 處理後有平滑趨勢，表示對於特徵特徵可能會有資訊損失之 情形。

64 表 4.1、影像 A 指標值 (A)處理 方法 RMean (前/後) GMean (前/後) BMean (前/後) RMSE PSNR(dB) SSIM 白平衡 43.122/68.875 59.503/88.615 79.914/107.936 22.128 28.775 0.632 去霧化 43.122/161.909 59.503/206.382 79.914/242.150 122.015 19.554 0.462 BEMD 43.122/53.822 59.503/76.792 79.914/104.308 17.320 49.449 0.752 表 4.2、影像 B 指標值 (B)處理方 法 RMean (前/後) GMean (前/後) BMean (前/後) RMSE PSNR(dB) SSIM 白平衡 29.855/70.500 40.535/91.187 69.961/130.354 62.764 21.632 0.588 去霧化 29.855/92.118 40.535/142.43 69.961/228.577 75.939 17.331 0.377 BEMD 29.855/39.560 40.535/55.322 69.961/100.820 24.778 50.855 0.774 表 4.3、影像 C 指標值 (C)處 理方法 RMean (前/後) GMean (前/後) BMean (前/後) RMSE PSNR(dB) SSIM 白平衡 80.699/104.458 95.495/90.645 108.855/107.479 7.201 45.311 0.833 去霧化 80.699/238.381 95.495/252.323 108.855/254.106 147.377 11.232 0.293 BEMD 80.699/74.696 95.495/88.237 108.855/102.577 28.230 50.560 0.697 表 4.4、影像 D 指標值 (C)處 理方法 RMean (前/後) GMean (前/後) BMean (前/後) RMSE PSNR(dB) SSIM BEMD 84.610/63.220 102.329/77.008 113.777/86.698 1.874 43.882 0.922

65

5. 結論與建議

能量傳輸模式可由太陽幾何模式計算入射輻射量，再大氣相關參數計算傳輸量，此 傳輸量即為日光返照影響量，將影像讀取得到之傳輸量扣除影響量得到改正後影像。但 是測試影像沒有拍攝當時反射表面模型，計算之傳輸量無法反映表面起伏；另外傳輸量 轉為灰度值時，因為缺乏乘常數與加常數，無法轉成灰度值對影像進行改正。 影像有色偏或資訊損失時，為了使色調均勻或增顯，使用了白平衡與去霧處理兩種 方法。由測試影像結果可以發現白平衡改善了色偏情形，去霧處理對於不明顯之特徵有 增顯之效果，但對於日光返照之影響無益。 由第一組測試影像結果發現殘餘項之特徵和範圍與原始影像日光反影之影響相當， 將原始影像扣除殘餘項後即為改正影像。影像結果發現在海面灰度均質區域可改善日光 反影之影響，但在礁塊區域對於礁塊會有資訊損失。鑲嵌之影像在邊界區域拼接良好， 無日光反影造成之光影不均情形。 由第二組測試影像結果發現，此組影像可能受解析度影響，相較於第一組測試影像 清晰度較低，因此較難細分訊號，殘餘項含有部分日光返照區域外之亮度，導致扣除殘 餘項後的影像結果過於模糊，顯示二維經驗模態分解用在日光返照之改正亦受到解析度 的限制。 使用平均值(Mean)、均方根誤差(RMSE)、峰值訊號雜訊比(PSNR)及結構相似性指標 (SSIM)對於影像結果進行統計指標計算。白平衡可以改善色調不均，RMSE 相對較小， SSIM 值亦顯現其保留原始結構程度較高，但對於日光反影之改正以視覺而言無益；去 霧化對於影像過度增顯，可以看出灰度值平均及RMSE 都相當高，PSNR 值在 30 以下， 為人眼無法接受，SSIM 值極低，表示其結構已被破壞。 均質海域影像經由 BEMD 處理，影像品質指標可幫助判斷影像品質。 BEMD 在礁 塊區域有資訊損失，導致影像處理後灰度值平均降低，並且SSIM 值較海域均質區域為 低，建議在均質區域可以使用BEMD 方法，若有礁塊則需要進行遮罩。

66

參考文獻

李暢游、孫標、高占義、張生，2011，呼倫湖多波段遙感水深反演模型研究，水利學報 42(12)：1423-1431。 施伶蒨、史天元、陳雅信、陳杰宗、陳佳勳，2011，以白平衡及去霧處理增顯東沙航空 影像，測量及空間資訊研討會。 胡珺婷，2012，二維經驗模態分解應用於航遙測影像增強之研究，碩士論文，國立中 山大學海洋環境及工程學系研究所。 陳炳旭，2002，數位相機自動白平衡：結合灰界理論與色域對應之研究，私立世新大 學平面傳播科技學系碩士論文。 陳雅信、史天元、陳杰宗、陳佳勳，2010，海域空照影像鑲嵌：東沙環礁案例探討，測 量及空間資訊研討會。

Cox, C. and W. Munk, 1954. Measurement of the Roughness of the Sea Surface from Photographs of the Sun’s Glitter. Journal of the Optical Society of America, 44(11): 838-850.

Centennial, N., 2001. Bidirectional Reflection Distribution Function, http://math.nist.gov/~ FH unt/appearance/brdf.html. (last access:2012/06/20)

Chiou, T.S. and C.S. Fuh, 2001. Automatic White Balance for Digital still Camera, MSC thesis, National Taiwan University.

Delaunay, B., 1934. Sur la sphere vide . Izvestia Akademii Nauk SSSR, Otdelenie Mate matic heskikhi Estestvennykh Nauk, 7:793-800.

Damerval, C., Sylvain and Valérie Perrier, 2005. A Fast Algorithm for Bidimensional EMD. IEEE Processing Letters, 12(10): 701-704.

Huang, N. E., Z. Shen, and S. R. Long et al, 1998. The Empirical Mode Decomposition Method and The Hilbert Spectrum for Non-stationary Time Series Analysis. Proc. R. Soc. Lond. A, 454:903-995.

Hochberg, E. J., Serge Andréfouët and M. R.Tyler, 2003. Sea Surface Correction of High Spatial Resolution Ikonos Images to Improve Bottom Mapping in Near-Shore Environments. IEEE Transactions on Geoscience and Remote Sensing, 41(7): 1724-1729.

Hedley, J.D., A. R. Harborne and P.J. Mumby, 2005. Simple and robust removal of sun glint for mapping shallow-water benthos. International Journal of Remote Sensing, 26(10): 2107-2112.

He, K., J. Sun and X. Tang, 2009.Single Image Haze Removal Using Dark Channel Prior. IEEE Conference on Computer Vision and Pattern Recognition.1957–1963.

Kay, S., J. D. Hedley, and S. Lavender, 2009. Sun Glint Correction of High and Low Spatial Resolution. Remote Sensing, 1: 697-730.

67

ed-ssim-m- code. (last access: 2012/06/20)

Pei, S.C. and Mary Tzeng, 2007. Uneven Illumination Removal and Image Enhancement Using Empirical Mode Decomposition. 20th Conference in Computer Vision, Graphics and Image Processing(CVGIP).

Rosa, M. C., P. Stefano, and Z. Elisabetta, 2006. Correction of Sun Glint Effect on MIVIS Data of the Sicily Campaign in July 2000. Annals of Geophysics, 49(1): 277-286.

ReSe, 2010. MODO, http://www.rese.ch/products/modo/index.html. (last access: 2012/09/17) Schläpfer, D., Dr. sc. Nat. and D. Odermatt, 2010. MODO User Manual. Rese Application

Schläpfer under support of the Remote Sensing Laboratories(RSL) of the University of Zurich.

Wikipedia, 2012a. Peak signal to noise ratio, http://en.wikipedia.org/wiki/Peak_signal-to-noi se _ ratio. (last access: 2012 /06/20)

Wikipedia, 2012b. Delaunay triangulation, http://en.wikipedia.org/wiki/Delaunaytriangulati on. (last access: 2012/06/20)

Wikipedia, 2012c. Sun glint, http://en.wikipedia.org/wiki/Sunglint. (last access: 2012/09/17) Wikipedia, 2012d.色溫, http://zh.wikipedia.org/wiki/色溫(last access: 2012/09/17)