米 碩
(四川大學(xué)水力學(xué)與山區(qū)河流開發(fā)保護國家重點實驗室,四川 成都 610065)
準確快速地測定海岸線的時空變化,對于海域使用管理具有十分重要的意義[1]。海岸線會在河口泥沙淤積、海平面上升等自然現(xiàn)象的影響下產(chǎn)生緩慢的變化。人類的填海造地工程、港口建設(shè)等人為活動會造成自然岸線的急劇變化。傳統(tǒng)獲取海岸線的時空變化信息的方法多為野外實地調(diào)查。但該方法效率低、周期長,并且數(shù)據(jù)統(tǒng)計處理難度較大[2]。
基于衛(wèi)星圖像,采用圖像處理的方法,可以較為方便的對海岸線進行數(shù)據(jù)統(tǒng)計。與傳統(tǒng)多波段或超光譜遙感圖像相比,Google Earth衛(wèi)星圖像更新快,易獲取且歷史信息全。但是Google Earth衛(wèi)星圖像為普通彩色圖像,可提供波段信息較少,因此需要設(shè)計合適的算法才能對海岸線進行準確有效的提取。
基于Google Earth衛(wèi)星圖像進行海岸線的識別與提取,要求待處理的圖像中包含海洋區(qū)域與岸上區(qū)域,因此待識別海岸線為兩區(qū)域分界線。圖像上海洋與陸地兩區(qū)域顏色、紋理等特性存在較大差異,分界線比較明顯,可針對海陸二者特征的差異,基于圖像處理的方式設(shè)計處理算法,將海陸邊界像素點分離,進而通過搜索方法找尋海岸線像素點坐標。再根據(jù)比例尺換算出實際距離,將海岸線數(shù)據(jù)量化,從而進行后續(xù)的分析。本文海岸線的提取算法流程見圖1。
圖1 海岸線提取算法流程圖
對于單獨的一張Google Earth衛(wèi)星圖像,需要對圖像進行降噪處理[3],否則圖像中無關(guān)噪點的存在有可能會對海岸線數(shù)據(jù)的提取造成干擾。因此本文算法需要在盡可能的保持原始信息即海陸邊緣特征的完整性的同時,又能夠去除信號中無用的信息如海上波浪、船舶等帶來的干擾。
本文采用了雙邊濾波來進行圖像降噪。雙邊濾波是一種可以保留邊緣信息同時降噪的濾波器。它是基于空間分布的高斯濾波函數(shù),可以使得離邊緣較遠的像素不會對邊緣像素產(chǎn)生較大的影響,這樣就避免了邊緣的模糊[4]。
圖像的邊緣檢測[5]是圖像處理的常用算法之一。邊緣檢測強調(diào)的是圖像對比度。檢測對比度也就是檢測亮度上的差別,由此可以增強圖像中的邊界特征,這些邊界出現(xiàn)正是圖像亮度上的差別。因此基于衛(wèi)星圖像中海洋與陸地亮度的顯著差異,邊緣檢測算法可以得到較為良好的效果。常見的邊界算子有Canny算子、Sobel算子、Laplace算子等,本文綜合考慮圖像處理效果及算法復(fù)雜度的影響,采用Sobel算子對Google Earth圖像進行邊緣檢測。Sobel算子主要用作邊緣檢測,它是一種離散型差分算子,用來計算圖像亮度函數(shù)灰度的近似值。
對圖像進行二值化操作可以將圖像上的像素點的灰度值設(shè)置為0或255,因此二值化后的圖像會呈現(xiàn)出明顯的黑白效果。圖像的二值化有利于對圖像進行進一步處理,使得圖像特征變得簡單且數(shù)據(jù)存儲量減少,能強化感興趣的目標的輪廓即海洋與陸地兩部分的圖像特征。將經(jīng)過邊緣檢測的灰度圖進行二值化處理后,海洋部分會以黑色像素點的形式呈現(xiàn),而陸地主體為白色像素點,黑白像素點的交界處即為需要提取的海岸線。因此基于黑白二值圖,就可以從一個種子點開始進行搜索,在遍歷全圖后,可以找到完整的黑白交界處的像素點,其連接起來即為需要提取的海岸線。
圖像形態(tài)學(xué)操作是基于形狀的一系列圖像處理操作的合集,主要有四個操作:膨脹、腐蝕、開運算、閉運算。
膨脹操作可以將圖像的輪廓加以膨脹;腐蝕操作即將圖像外圍的突出點加以腐蝕;開運算可以放大裂縫和低密度區(qū)域,消除小物體,在平滑較大物體的邊界時,不改變其面積,消除物體表面的突起;閉運算可以排除二值化圖像中的小型黑洞,突出比原圖輪廓區(qū)域更暗的區(qū)域,將兩個區(qū)域連接起來,形成連通域。
本文采用的算法先使用腐蝕操作分割出獨立的圖像元素,將不屬于海岸線的像素點獨立成個體,進而可以采取閉運算填充小的孔洞,填補海岸線像素點,使得代表海岸線的像素點平滑而連續(xù)。最后采用膨脹操作填補腐蝕操作時被蠶食的高亮區(qū)域。經(jīng)過這一系列的形態(tài)學(xué)操作,可以在分離出海岸線像素區(qū)域與無關(guān)像素區(qū)域的同時增強海岸線的特征,減少無關(guān)噪點,有利于后續(xù)海岸線數(shù)據(jù)點的提取。
經(jīng)過預(yù)處理及形態(tài)學(xué)操作后的圖像,在計算機中被存儲為一個大小為圖像像素長寬的矩陣,矩陣中的每個數(shù)據(jù)均代表一組像素點坐標。為了提取出海岸線像素點的坐標數(shù)據(jù),本文采用八方向搜索的方式進行海岸線數(shù)據(jù)點的提取,見圖2,首先在圖像矩陣從海洋黑色像素區(qū)域查找到第一個陸地像素即白色像素點作為種子點,坐標為(x,y),然后以(x-1,y),(x-1,y+1),(x,y+1),(x+1,y+1),(x+1,y),(x+1,y-1),(x,y-1),(x-1,y-1)的順序在八個方向查找下一個海岸線像素點位置,直到所有海岸線邊界點被標記。
圖2 八方向搜索提取海岸線數(shù)據(jù)點
圖3展示了本文海岸線提取流程的一個示例,從圖中可以看出經(jīng)過濾波、邊緣檢測、形態(tài)學(xué)操作,二值化等圖像操作后,海岸線能夠準確地被識別并被提取。具有了海岸線的散點數(shù)據(jù)后,可以計算出陸地面積等信息,為海岸線的管理建設(shè)等提供數(shù)據(jù)支持。
圖3 遼寧省菊花島海岸線提取過程
本文選取杭州灣為研究地點,提取了1986年~2016年間22個年份的海岸線,以1996年杭州灣地區(qū)Google Earth圖像為底圖,在1996年~2014年間選取了五個年份進行海岸線的提取,見圖4。
圖4 杭州灣海岸線提取結(jié)果
在獲取到提取的海岸線數(shù)據(jù)后,可以計算出所研究區(qū)域海岸線所包圍陸地面積,進而可以計算出陸地面積的變化趨勢,從而分析港口海岸地區(qū)的時空變化趨勢。提取出的海岸線數(shù)據(jù)為一系列像素點的散點坐標,把圖像邊界像素點添加進散點數(shù)據(jù)集,所有像素點連接即為任意形狀的多邊形,因此可采用任意多邊形面積公式[6]計算陸地的像素面積,公式如下:
其中,(xi,yi)代表海岸線散點的坐標且x1=xN+1,y1=yN+1,計算出陸地的像素面積后,則需要由圖像比例尺換算出真實的陸地面積??梢酝ㄟ^將Google地圖衛(wèi)星圖像的比例尺距離除以其所占像素,進而獲得每像素代表的實際長度,進而可以計算得出真實的陸地面積。
圖5 1986年~2016年杭州灣的面積變化
經(jīng)過圖像處理并提取數(shù)據(jù)后,計算出杭州灣在1986年~2016年間的面積變化,見圖5。圖5數(shù)據(jù)表明,自1996年以來,杭州灣沿岸的陸地面積增長開始加快。從1996年到2014年,杭州灣沿岸陸地面積年均增長率為38.1 km2/a,是1986年~1996年期間杭州灣沿岸陸地面積年均增長面積(17.7 km2/a)的兩倍多。自2014年以來,杭州灣沿岸土地面積趨于穩(wěn)定,并略微呈現(xiàn)出侵蝕的趨勢。經(jīng)過此例可以看出,本文提出的海岸線數(shù)據(jù)提取算法可以基于Google Earth衛(wèi)星圖像對海岸線數(shù)據(jù)進行提取,并進行后續(xù)的計算和分析,避免了傳統(tǒng)野外實地調(diào)查方法效率低、周期長及衛(wèi)星遙感方法時效性不強的弊端,證實了該方法在港口海岸工程方面應(yīng)用的可行性及有效性。