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

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

渤海是中國(guó)半封閉內(nèi)海，自凈能力差，周邊又都是重工業(yè)發(fā)達(dá)城市，目前海水污染非常嚴(yán)重。據(jù)不完全統(tǒng)計(jì)，每年排入渤海的污染物達(dá)70多萬(wàn)噸，占全國(guó)入海污染物總量的47．7%。渤海沿岸直排排污口攜帶大量懸浮顆粒物、營(yíng)養(yǎng)鹽、重金屬入海，渤海的海洋生態(tài)環(huán)境毀壞嚴(yán)重、壓力巨大。

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

圖1 2014年渤海監(jiān)測(cè)的陸源入海排污口狀況

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

圖像分割是把圖像分割成一些連續(xù)的、非重疊的區(qū)域。這種迭代的圖像分割算法目前在計(jì)算機(jī)視覺(jué)方面應(yīng)用越來(lái)越廣泛[1-3]。Boykov[4]第一個(gè)在圖像分割領(lǐng)域中提出馬爾科夫隨機(jī)場(chǎng)模型，基于這個(gè)模型，提出了Graph-Cut圖像分割算法。只需在背景和前景中標(biāo)記簡(jiǎn)單的幾筆，此算法能夠在短時(shí)間內(nèi)分割出想要的圖形。在Boykov和Jolly研究的基礎(chǔ)上，又出現(xiàn)了大量基于Graph-Cut和用戶筆畫(huà)的算法模型[4-5]。Slabaugh[6]提出了一種橢圓形狀先驗(yàn)?zāi)Ｐ停怂惴ㄖ饕獞?yīng)用于對(duì)醫(yī)學(xué)中共振圖像的血管和淋巴結(jié)的提取。Rother[7]提出了Grab-Cut算法，此算法的迭代速度更快，而且所需的用戶筆畫(huà)更少，適用范圍廣泛，能夠應(yīng)用在很多領(lǐng)域。本文主要進(jìn)一步利用改進(jìn)的Grab-Cut分割算法對(duì)天津沿岸的排污口進(jìn)行監(jiān)測(cè)。

本研究利用美國(guó)USGS Landsat-8衛(wèi)星的遙感圖像，結(jié)合改進(jìn)的Grab-Cut圖像分割算法，對(duì)天津周邊渤海灣的排污口進(jìn)行了監(jiān)測(cè)，并對(duì)污染區(qū)域進(jìn)行了時(shí)間和空間上的分析，以期能對(duì)相關(guān)部門(mén)關(guān)于河流污染的治理提供參考。

1 Grab-Cut算法

1.1 Graph Cut算法

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

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

前景和背景的頻率分布直方圖中，總會(huì)出現(xiàn)一個(gè)或兩個(gè)所占比例較大的灰度值，這個(gè)灰度值就是前景和背景中的主灰度值，前景和背景的區(qū)分依賴于這個(gè)主灰度值。

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

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

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

上式U用來(lái)衡量α的正確度，表達(dá)式為：

式中：[…]代表一種判斷，滿足括號(hào)內(nèi)的條件為1，不滿足括號(hào)內(nèi)的條件為0。dis(m,n)表示兩點(diǎn)之間的歐幾里德距離。依據(jù) Boykov[4]，β=(2＜(zm-zn)2＞)-1，＜…＞表示均值。(m,n)∈C中的C表示像素點(diǎn)與鄰像素點(diǎn)配對(duì)組成的集合，對(duì)于某個(gè)像素點(diǎn)而言，其鄰像素點(diǎn)為其前后左右4個(gè)以及斜對(duì)角的4個(gè)共8個(gè)像素點(diǎn)，對(duì)于一個(gè)具體的像素點(diǎn)而言，C中的元素共有8個(gè)?；诖罅康挠?xùn)練集測(cè)試，γ的值選為50較為合適(見(jiàn)文獻(xiàn)[9])。

1.2 Grab-Cut算法模型的建立

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

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

V用來(lái)衡量整個(gè)分割的穩(wěn)定性，表達(dá)式為：

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

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

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

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

Grab-Cut算法會(huì)迭代計(jì)算去優(yōu)化能量函數(shù)，以求得最優(yōu)解，算法會(huì)對(duì)公式(6)中的向量α自動(dòng)逼近，使分割的邊緣部分清晰準(zhǔn)確。

1.3 進(jìn)一步改進(jìn)的Grab-Cut算法

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

本文利用Grab-Cut算法來(lái)處理大尺寸的遙感圖像，并用于海洋污染的監(jiān)測(cè)。因?yàn)閷?duì)整個(gè)遙感圖像進(jìn)行圖像分割的效率低且耗時(shí)長(zhǎng)，我們基于對(duì)污染區(qū)域的認(rèn)識(shí)，劃定感興趣的區(qū)域，然后再對(duì)這些區(qū)域進(jìn)行圖像分割，從而產(chǎn)生更好的分割效果。

2 遙感圖像的獲取

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

3 天津周邊渤海灣入海排污口的監(jiān)測(cè)

3.1 天津入海排污口的監(jiān)測(cè)

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

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

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

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

圖3 改進(jìn)的Grab-Cut算法分割排污口污染區(qū)域

圖中的結(jié)果顯示，算法的分割效果非常好，分割出了污染區(qū)域。能夠看出各個(gè)流域入?？诘奈廴颈容^嚴(yán)重，污染面積比較大。利用比例尺，估算出了各個(gè)污染區(qū)域面積的大小。其中，薊運(yùn)河、潮白新河、永定新河流域和海河流域的污染最嚴(yán)重，污染面積分別約為259 km2和217 km2。薊運(yùn)河、潮白新河、永定新河流域分布區(qū)域廣闊，流域內(nèi)的河流源起河北，流經(jīng)天津市重工業(yè)發(fā)達(dá)的濱海新區(qū)，最后流入渤海；海河流經(jīng)人口較多的南開(kāi)區(qū)、和平區(qū)以及東麗區(qū)，最后經(jīng)過(guò)重工業(yè)發(fā)達(dá)的濱海新區(qū)注入渤海，河流中攜帶大量的生活污染物和工業(yè)污染物。這些原因造成了這兩個(gè)流域入?？诘奈廴痉浅?yán)重。獨(dú)流減河、子牙新河、北排河流域流經(jīng)近幾年發(fā)展迅速的津南區(qū)，區(qū)域中的重工業(yè)發(fā)展迅速，河流也有較大的污染，但相比于上面兩個(gè)流域污染相對(duì)較輕，總污染面積約為116 km2。

3.2 天津入海排污口污染區(qū)域隨時(shí)間的變化

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

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

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

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

圖4～圖5顯示，各個(gè)排污口的污染區(qū)域都有變化。利用比例尺，估算出每個(gè)區(qū)域的污染面積。各個(gè)區(qū)域的污染面積隨時(shí)間的變化情況見(jiàn)表1，折線圖見(jiàn)圖6。

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

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

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

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

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

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

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

圖7所列的方法，“Quick Search”和“Magnetic Lasso”是PhotoShop圖像處理軟件中較為常用的圖像分割算法，“Quick Search”在用戶指定前景后，會(huì)自動(dòng)識(shí)別與前景相似的區(qū)域，但分割效果較差，可能誤將大量無(wú)關(guān)區(qū)域識(shí)別為要分割的前景；“Magnetic Lasso”可以較好地分割出前景，但其缺點(diǎn)是需要用戶比較詳細(xì)地標(biāo)記出邊界，這種方法繁瑣、費(fèi)時(shí)且效率低下。有效的圖像分割算法應(yīng)該是只需要簡(jiǎn)單的幾筆標(biāo)記就能夠完成圖像分割；已有的“Graph Cut”算法雖然標(biāo)記較少，但是在圖像的邊緣處光滑性差，準(zhǔn)確率低，難得到一塊完整的分割結(jié)果；改進(jìn)的“Grab-Cut”算法能夠較好地區(qū)別要分割的背景和前景，在較少的用戶標(biāo)記下，能夠很快地分割出所需區(qū)域。綜合比較，改進(jìn)的“Grab-Cut”算法的圖像分割效果在已知算法中效果最優(yōu)。

4 結(jié)論

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

[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.