遙感圖像污染監測分析的Grab-Cut算法實現

宋占杰 ，湯 楊，王兵振 ，張東亮，田 震

（1．天津大學 數學學院，天津 300354；2．天津大學 電視與圖像信息研究所，天津 300072；3．國家海洋技術中心，天津 300112）

渤海是中國半封閉內海，自凈能力差，周邊又都是重工業發達城市，目前海水污染非常嚴重。據不完全統計，每年排入渤海的污染物達70多萬噸，占全國入海污染物總量的47．7%。渤海沿岸直排排污口攜帶大量懸浮顆粒物、營養鹽、重金屬入海，渤海的海洋生態環境毀壞嚴重、壓力巨大。

天津緊鄰渤海，并且位于渤海的核心區，與渤海其它海域相比，天津周邊海域海洋自凈能力更加脆弱。同時，天津又是重工業較為發達的城市，每年向渤海排放大量污水，使得天津周邊的海洋污染日趨嚴重。圖1是2014年渤海監測的陸源入海排污口狀況。圖中監測點的顏色越紅，污染越嚴重。圖中顯示，天津市周邊污染十分嚴重的排污口主要有3個，對渤海灣造成了嚴重的污染。

圖1 2014年渤海監測的陸源入海排污口狀況

傳統的入海排污口監測考察主要依靠船只人工采樣和樣本分析，這種方法耗時費力，受環境的制約較大。隨著遙感技術的快速發展，以及大量高精度遙感衛星的升空，高分辨率、多頻段、價格低廉的遙感數據的獲取成為現實。這些遙感數據的使用為我們更好地對天津周邊的海洋污染進行追蹤提供了便利，且一些圖像分割算法能夠幫助我們對污染區域進行監測和量化分析。

圖像分割是把圖像分割成一些連續的、非重疊的區域。這種迭代的圖像分割算法目前在計算機視覺方面應用越來越廣泛[1-3]。Boykov[4]第一個在圖像分割領域中提出馬爾科夫隨機場模型，基于這個模型，提出了Graph-Cut圖像分割算法。只需在背景和前景中標記簡單的幾筆，此算法能夠在短時間內分割出想要的圖形。在Boykov和Jolly研究的基礎上，又出現了大量基于Graph-Cut和用戶筆畫的算法模型[4-5]。Slabaugh[6]提出了一種橢圓形狀先驗模型，此算法主要應用于對醫學中共振圖像的血管和淋巴結的提取。Rother[7]提出了Grab-Cut算法，此算法的迭代速度更快，而且所需的用戶筆畫更少，適用范圍廣泛，能夠應用在很多領域。本文主要進一步利用改進的Grab-Cut分割算法對天津沿岸的排污口進行監測。

本研究利用美國USGS Landsat-8衛星的遙感圖像，結合改進的Grab-Cut圖像分割算法，對天津周邊渤海灣的排污口進行了監測，并對污染區域進行了時間和空間上的分析，以期能對相關部門關于河流污染的治理提供參考。

1 Grab-Cut算法

1.1 Graph Cut算法

對于每一幅遙感圖像z=(z1,…,zn,…zN)，zn表示圖像中的每一個像素點的灰度值。定義α=(α1,…,αn,…αN)，其中 αn∈{0,1}，0 代表像素點屬于背景，1 代表像素點屬于要分離出的前景。本文中，在算法進行前會要求使用者在圖形中初步選定一個矩形框，矩形框內αn=1，矩形框外αn=0。

針對αn=1的區域和αn=0的區域，即前景與背景，分別建立像素灰度值的頻率分布直方圖：

前景和背景的頻率分布直方圖中，總會出現一個或兩個所占比例較大的灰度值，這個灰度值就是前景和背景中的主灰度值，前景和背景的區分依賴于這個主灰度值。

利用每個圖像像素點的灰度值向量z=(z1,…,zn,…zN)，和區分前景和背景的向量 α=(α1,…,αn,…αN)，以及前景和背景的灰度頻率分布直方圖θ，可以建立最小能量函數：

最后優化出的向量α?就能夠較好地分離出前景和背景。

其中，E(α,θ,z)被定義為 Gibbs能量：

上式U用來衡量α的正確度，表達式為：

式中：[…]代表一種判斷，滿足括號內的條件為1，不滿足括號內的條件為0。dis(m,n)表示兩點之間的歐幾里德距離。依據 Boykov[4]，β=(2＜(zm-zn)2＞)-1，＜…＞表示均值。(m,n)∈C中的C表示像素點與鄰像素點配對組成的集合，對于某個像素點而言，其鄰像素點為其前后左右4個以及斜對角的4個共8個像素點，對于一個具體的像素點而言，C中的元素共有8個。基于大量的訓練集測試，γ的值選為50較為合適(見文獻[9])。

1.2 Grab-Cut算法模型的建立

傳統的Graph Cut算法不需要使用像素點的RGB信息，它基于像素點的灰度值直接進行計算。為了更好地使用這些顏色信息，建立了2個高斯混合模型（GMMs）[8-9]，一個對應于背景，一個對應于要分割出的前景。為了更好地應用高斯混合模型，對每個像素 zn賦值 kn，kn∈{1,…,K}，這樣，對于一整幅圖像，有 k={k1,…,kn,…,kN}。

借鑒上述Graph Cut算法，建立如下新的Grab-Cut算法最小能量函數：

V用來衡量整個分割的穩定性，表達式為：

E（α,k,θ,z）定義為如下的 Gibbs能量：

其中 α=(α1,…,αn,…αN)，αn∈{0,1}，0 代表背景，1代表想要分離出的前景。類似于Graph Cut算法，式中U衡量α的正確度，V衡量穩定性。分別定義為如下形式：

其中 p(…)代表高斯分布，π(…)為 αn與 kn的混合加權系數。

參數 γ=50，β=(2＜‖zm-zn‖2＞)-1，其中的＜…＞代表均值；C代表某一像素點與其鄰像素點組合形成的集合[7]。

Grab-Cut算法會迭代計算去優化能量函數，以求得最優解，算法會對公式(6)中的向量α自動逼近，使分割的邊緣部分清晰準確。

1.3 進一步改進的Grab-Cut算法

Grab-Cut算法在以前經常用來處理一些小尺寸照片，很少用在尺寸很大的高分辨率遙感圖像上。

本文利用Grab-Cut算法來處理大尺寸的遙感圖像，并用于海洋污染的監測。因為對整個遙感圖像進行圖像分割的效率低且耗時長，我們基于對污染區域的認識，劃定感興趣的區域，然后再對這些區域進行圖像分割，從而產生更好的分割效果。

2 遙感圖像的獲取

2013年2月11日，NASA成功發射了Landsat-8衛星，衛星上新增兩個主要載荷：運營性陸地成像儀OLI和熱紅外傳感器TIRS，2013年5月30日開始向全球提供免費數據下載。兩個主要載荷所收集的數據包含11個波段的信息，不同波段的組合可以得到不同的圖像，以用于不同的研究。本文的數據只用到其中的紅、綠、藍三個波段。

3 天津周邊渤海灣入海排污口的監測

3.1 天津入海排污口的監測

天津市河流流域眾多，包括薊運河、海河、獨流減河、北排河等等。天津又是重工業較為發達的城市，大量的重工業企業聚集在河流周邊，這些企業排放的污水大部分流入到這些河流當中。每年，這些河流攜帶大量的工業生活污染物排入渤海，給渤海的生態環境造成了巨大的破壞。

圖2是一幅由美國USGS landsat-8衛星于2014年3月21日拍攝的天津市周邊渤海灣的遙感圖像。從圖中可以看出，天津周邊渤海灣的污染非常嚴重，黑色的區域均為河流排向渤海的污染物所形成的異常圖。天津周邊淺海區域的海藻繁殖較為嚴重，整體呈現藍綠色，所以一旦有污水排入會與這些海藻形成鮮明的對比。

天津沿海地區主要分為3個流域，自北向南分別為薊運河、潮白新河、永定新河流域；海河流域；獨流減河、子牙新河、北排河流域。本文對這3個主要流域，利用改進的Grab-Cut算法來分割出各個流域排污口的污染區域，如圖3。圖中綠色方框表示“感興趣區域”，紅色線條表示“標記的前景”，藍色線條表示“標記的背景”。

圖2 天津周邊渤海灣污染情況

圖3 改進的Grab-Cut算法分割排污口污染區域

圖中的結果顯示，算法的分割效果非常好，分割出了污染區域。能夠看出各個流域入海口的污染比較嚴重，污染面積比較大。利用比例尺，估算出了各個污染區域面積的大小。其中，薊運河、潮白新河、永定新河流域和海河流域的污染最嚴重，污染面積分別約為259 km2和217 km2。薊運河、潮白新河、永定新河流域分布區域廣闊，流域內的河流源起河北，流經天津市重工業發達的濱海新區，最后流入渤海；海河流經人口較多的南開區、和平區以及東麗區，最后經過重工業發達的濱海新區注入渤海，河流中攜帶大量的生活污染物和工業污染物。這些原因造成了這兩個流域入海口的污染非常嚴重。獨流減河、子牙新河、北排河流域流經近幾年發展迅速的津南區，區域中的重工業發展迅速，河流也有較大的污染，但相比于上面兩個流域污染相對較輕，總污染面積約為116 km2。

3.2 天津入海排污口污染區域隨時間的變化

針對已有的改進算法作為典型實例，本文分析了2014年到2016年這3年，上述排污口污染區域的變化情況。所選的4幅遙感圖像的拍攝時間分別為2014年2月21日、2015年3月24日、2016年3月26日以及2016年12月23日。

圖4展示了2014-2016年天津渤海灣整體的污染變化情況。圖5展示了各個排污口在2014-2016年局部的污染變化情況。

圖4 2014-2016年天津渤海灣污染整體變化情況

圖5 2014-2016年天津渤海灣各個排污口污染局部變化情況

圖4～圖5顯示，各個排污口的污染區域都有變化。利用比例尺，估算出每個區域的污染面積。各個區域的污染面積隨時間的變化情況見表1，折線圖見圖6。

表1 各排污口面積隨時間的變化情況

圖6 各排污口面積隨時間的變化情況

對比污染面積折線圖(圖6)以及各個排污口污染區域變化情況圖(圖5)，薊運河、潮白新河、永定新河流域（圖5左邊一列）和獨流減河、子牙新河、北排河流域(圖5右邊一列)排污狀況有較大改善。但海河流域（圖5中間一列）的排污狀況雖然有所減輕，但相比以上兩個區域，一直非常嚴重，這間接說明海河流經的濱海新區重工業比較發達，大量重工業企業污染物和生活污水排向海河流域，導致海河的污染十分嚴重。另一個原因是天津港位于海河流域下游，港內海水與外界交換很少，自凈能力非常差。因此，對海河周邊重工業企業的排污治理刻不容緩。

總體而言，天津市對于工業排污的治理有一定成效，近3年排污口對渤海的污染明顯減輕，但仍需進一步嚴加管理。

3.3 其它圖像分割算法與“Grab-Cut”算法的比較

本文選擇了一些圖像分割算法與Grab-Cut算法進行比較，從而分析不同的圖像分割算法分割效果的優劣。這里使用“薊運河、潮白新河、永定新河流域入海排污口”圖像，分別使用PhotoShop中的“Quick Search”工具和“Magnetic Lasso”工具，以及“Graph Cut算法”這3種圖像分割算法與“改進的Grab-Cut算法”進行分割效果比對（圖7）。

圖7 幾種圖像分割算法的比較

圖7所列的方法，“Quick Search”和“Magnetic Lasso”是PhotoShop圖像處理軟件中較為常用的圖像分割算法，“Quick Search”在用戶指定前景后，會自動識別與前景相似的區域，但分割效果較差，可能誤將大量無關區域識別為要分割的前景；“Magnetic Lasso”可以較好地分割出前景，但其缺點是需要用戶比較詳細地標記出邊界，這種方法繁瑣、費時且效率低下。有效的圖像分割算法應該是只需要簡單的幾筆標記就能夠完成圖像分割；已有的“Graph Cut”算法雖然標記較少，但是在圖像的邊緣處光滑性差，準確率低，難得到一塊完整的分割結果；改進的“Grab-Cut”算法能夠較好地區別要分割的背景和前景，在較少的用戶標記下，能夠很快地分割出所需區域。綜合比較，改進的“Grab-Cut”算法的圖像分割效果在已知算法中效果最優。

4 結論

對比于傳統的“Quick Search”和“Magnetic Lasso”圖像分割算法，改進的“Grab-Cut”算法的用戶標記要更少；而已有的“Graph Cut”算法雖然標記較少，但圖像邊緣處光滑性差，準確率低，很難得到一塊完整的分割結果。而改進的“Grab-Cut”算法在很少的用戶標記下，能夠很快地分割出所需區域。綜合來看，改進的“Grab-Cut”算法的圖像分割算法相較于前幾種而言效果最好。我們將改進的“Grab-Cut”算法應用于渤海灣的污染監測。發現近3年天津周圍渤海灣的污染非常嚴重，但近1年來污染情況有所改善。

[1]Kolmogorov V,Zabih R.What energy functions can be minimized via graph cuts?[J].Pattern Analysis&Machine Intelligence IEEE Transactions on,2004,26(2):147-159.

[2]Boykov Y,Veksler O,Zabih R.Fast approximate energy minimization via graph cuts[J].Pattern Analysis&Machine Intelligence IEEE Transactions on,2001,23(11):1222-1239.

[3]Boykov Y,Funka-Lea G.Graph cuts and efficient N-D image segmentation[J].International Journal of Computer Vision,2006,70(2):109-131.

[4]Boykov Y Y,Jolly M P.Interactive graph cuts for optimal boundary®ion segmentation of objects in N-D images[J].2001,1:105-112.

[5]Grady L.Random walks for image segmentation.[J].Pattern Analysis&Machine Intelligence IEEE Transactions on,2006,28(11):1768-1783.

[6]Slabaugh G,Unal G.Graph cuts segmentation using an elliptical shape prior[C]//IEEE International Conference on Image Processing.IEEE,2005:II-1222-5.

[7]Rother B C.Kolmogorov V and Blake A(2004)"GrabCut":interactive foreground extraction using iterated graph cuts[C]//ACM Trans Graph.2010.

[8]Ruzon M A,Tomasi C.Alpha estimation in natural images[J].2000,1:18-25.

[9]Blake A,Rother C,Brown M,et al.Interactive Image Segmentation Using an Adaptive GMMRF Model[M]//Computer Vision-ECCV 2004.Springer Berlin Heidelberg,2004:428-441.