基于計算機視覺的平水橋海灘風暴潮響應(yīng)分析

匡翠萍，劉 旭，夏子龍，朱 磊，叢 新

（1. 同濟大學土木工程學院，上海 200092；2. 深圳市龍華區(qū)水污染治理中心，深圳 518110；3. 河北省地礦局第八地質(zhì)大隊，河北秦皇島 066000）

風暴潮是由強風或氣壓驟變引起的海平面突變現(xiàn)象，其在短期內(nèi)向海灘輸入大量的能量，這種極端天氣條件下產(chǎn)生的強風-浪-流相互作用對海灘地貌的塑造比正常天氣高出1～2個量級［1-2］，進而造成海灘地貌發(fā)生顯著變化［3］。風暴潮過程常伴隨有風暴增水，這往往會使風暴潮對海灘作用的上界向陸推移并造成沿岸侵蝕災(zāi)害；另外，波高增大會造成碎波帶顯著變寬，可能導(dǎo)致風暴潮巨浪的能量在到達海灘之前已經(jīng)耗散殆盡，這種情況不造成強烈的海灘侵蝕［4-5］。風暴潮作為海灘地貌短期演變的重要動力形式，一直是研究的熱點和重點。對海灘進行監(jiān)測是評估風暴潮影響的基本手段之一，監(jiān)測方法包括人工實地測量和攝影測量［6-7］，其中攝影測量可以有效解決風暴潮期間人工測量困難的問題。

隨著計算機視覺技術(shù)的發(fā)展，攝影測量的精度和自動化程度都得到了極大提高。計算機視覺技術(shù)為海岸線的精確提取指明了方向。近些年來，計算機視覺技術(shù)在海岸工程中的應(yīng)用逐漸增多。郭俊麗等［8］基于Argus記錄的舟山東沙海灘圖像，研究了東沙海灘對兩次連續(xù)風暴的響應(yīng)，結(jié)果表明兩次風暴使得海灘濱線較風暴前后退10.91 m。Cao 等［9］將Argus系統(tǒng)監(jiān)測與傳統(tǒng)監(jiān)測相結(jié)合，利用計算機視覺技術(shù)對日照海岸的監(jiān)測圖片進行分析，研究結(jié)果表明養(yǎng)灘工程兩個月后海拔-6 m以上的沉積物體積大致保持穩(wěn)定，通過視頻圖像可以量化波浪的翻越特征。

研究海灘對風暴潮響應(yīng)的方法通常為視頻圖像數(shù)據(jù)或數(shù)值模擬，本文將視頻圖像技術(shù)與數(shù)值模擬相結(jié)合，以提高風暴潮反演精度、更多元化地豐富海灘動力地貌過程研究。首先利用Delft 3D 模型對2016 年7 月20 日的“720”風暴潮過程的動力環(huán)境進行復(fù)演，然后基于計算機視覺技術(shù)處理“720”風暴潮前后的視頻圖像數(shù)據(jù)，反演風暴潮前后地形，結(jié)合模擬的動力環(huán)境與反演的地形，分析極端天氣過程影響下秦皇島平水橋海灘演變短時響應(yīng)特征及恢復(fù)期過程。

1 數(shù)據(jù)來源與方法

平水橋海灘位于秦皇島北戴河西海灘東部，為砂質(zhì)海灘，岸線呈SW～NE走向，海岸沉積以中細沙為主，東側(cè)分布有老虎石天然岬角和人工潛堤岬頭，西側(cè)為平直岸段。平水橋?qū)儆谌醭焙０叮嗄昶骄辈罴s為0.74 m，最大潮差1.50 m。波浪以風浪為主，常浪向為東南向。2016年7月19～21日，受溫帶氣旋影響，渤海沿海出現(xiàn)一次較強的風暴潮過程，最大風暴增水發(fā)生在河北省黃驊站，為115 cm；秦皇島站達到了黃色警戒水位，風暴增水為59 cm［10］。

2014年10月，河北省地礦局第八地質(zhì)大隊在秦皇島平水橋海灘安裝了Argus 視頻監(jiān)測系統(tǒng)，其完整地記錄下了“720”風暴潮的全過程。Argus系統(tǒng)由現(xiàn)場視頻測站、可編程控制器、現(xiàn)場高程控制點系統(tǒng)、專用通訊網(wǎng)、中心服務(wù)器系統(tǒng)等硬件和Deltares公司開發(fā)的后處理軟件2 個部分組成［11］。如圖1 所示，圖1a 為Argus 視頻測站安裝位置，圖1b 為攝像終端，圖1c為中心服務(wù)器系統(tǒng)，圖1d為6個攝像頭監(jiān)測圖像融合后的海灘全貌。Argus 系統(tǒng)對海灘進行24 h連續(xù)監(jiān)測，圖像采集頻率為30 min，記錄2560×1920 pixel2的瞬時圖像和10 min時均圖像。

基于Argus 監(jiān)測系統(tǒng)采集圖像可以反演得到地形等信息，其工作技術(shù)路線如圖2 所示。首先對所采集的海灘圖像進行處理，結(jié)合潮流模型提供的潮位信息；再根據(jù)攝影測量原理進行坐標轉(zhuǎn)換，即可得到水邊線坐標；最后由不同時刻的水邊線反演地形的變化。

圖2 地形反演流程Fig.2 Flow chart of topography extraction

1.1 數(shù)學模型及參數(shù)設(shè)置

基于圖2 的技術(shù)流程，需要建立精確的潮流模型，以提供拍攝時刻的潮位信息。采用Delft 3DFlow 模塊建立二維潮流模型模擬研究區(qū)域的潮位過程，使用渤海-秦皇島雙重嵌套網(wǎng)格［12］，通過潮位過程來控制模型邊界，其中渤海大模型以大連潮位站至煙臺潮位站的連線為邊界，對大連至煙臺沿線各點的潮位數(shù)據(jù)插值，得到大模型的海域開邊界；秦皇島小模型的邊界由三條海域開邊界和一條岸線閉邊界組成，其中三條海域開邊界的潮位條件由渤海大模型的計算結(jié)果給出，岸線閉邊界流速為零。渤海大模型的曼寧系數(shù)取0.012～0.018，秦皇島小模型的曼寧系數(shù)取均值0.014 5，計算時間步長均設(shè)置為1 min，水平渦黏系數(shù)取15 m2·s-1。

通過布設(shè)在平水橋南面700 m處的波浪監(jiān)測浮標提供每小時的波浪數(shù)據(jù)，浮標處水深約6 m。在“720”風暴潮期間，浮標遭到破壞，造成部分監(jiān)測數(shù)據(jù)的缺失，因而采用第三代SWAN模型模擬研究區(qū)域波浪場，以補充缺失時間段內(nèi)的波浪數(shù)據(jù)。其中風應(yīng)力計算如下：

式中：ρa為空氣的密度；W10表示海洋表面上方10 m處的風速；Cd代表風拖曳力系數(shù)，是精確模擬風暴潮的關(guān)鍵參數(shù)，宋竑霖［13］基于前人建立的風速與風拖曳力系數(shù)之間的關(guān)系，歸納得

采用歐洲中期天氣預(yù)報中心（European Centre for Medium-Range Weather Forecasts，ECMWF）提供的風場數(shù)據(jù)作為模型中的風作用過程。圖3為模型輸入的2016 年7 月風速和風向。可以看出，風暴潮整個生長衰落過程自7 月19 日18 時～7 月21 日7時，歷時約1.5 d。風速從7月19日18時起開始逐漸增大，7月20日9時～7月21日5時風速均在10 m·s-1以上，風向為80°～180°，其后至7月21日7時風速逐漸衰減。

圖3 2016年7月ECMWF風速及風向Fig.3 Wind speed and direction in July 2016 from ECMWF

1.2 數(shù)字高程模型

結(jié)合已建立的潮流模型，對攝得的圖像經(jīng)過預(yù)處理、水邊線提取、坐標轉(zhuǎn)換等步驟，得到拍攝時刻的水邊線實地三維坐標，對多時相的水邊線進行插值，得到研究區(qū)域的數(shù)字高程模型（Digital Elevation Model，DEM）。

1.2.1 預(yù)處理圖像

現(xiàn)場視頻測站測得的圖片包含大量的噪聲，對海灘進行水邊線提取時，需先通過空間濾波圖像處理方法對噪聲進行預(yù)處理。空間濾波用于模糊處理和降低噪聲。典型的隨機噪聲由灰度級的急劇變化產(chǎn)生，通過使用濾波器模板將圖像每個像素的灰度值設(shè)置為其鄰域內(nèi)像素的平均值，從而降低噪聲灰度的急劇變化，有效防止其在邊緣檢測識別中被誤檢。高斯濾波以二維高斯分布為基礎(chǔ)生成濾波模板。濾波分布函數(shù)為

式中：x、y為像素點的坐標；σ為分布函數(shù)標準差；G（x，y）為高斯函數(shù)。

用高斯函數(shù)和輸入圖像的卷積形成一幅平滑后的圖像如下：

式中：f（x，y）為輸入圖像；fs（x，y）為輸出圖像。

圖4為瞬時圖像與高斯濾波后圖像的對比。經(jīng)高斯濾波降噪后，圖像噪點被平滑，更利于后續(xù)水邊線的提取。

圖4 高斯濾波效果Fig.4 Effect of Gaussian filter

1.2.2 提取水邊線

經(jīng)空間濾波圖像預(yù)處理后，對圖像中的水邊線進行提取，主要步驟為梯度計算、非最大抑制、雙閾值處理、像素坐標導(dǎo)出。

計算圖像中每個像素點的梯度幅度M（x，y）和方向a（x，y），則

首先，定義0°、45°、90°和135°四個方向，使用梯度幅度和方向來估計每一點處的邊緣強度和方向，初步檢測出邊緣像素。然后，對邊緣像素應(yīng)用非極大值抑制，以去除邊緣檢測過程中出現(xiàn)的雜散邊緣信息。對每個像素進行非極大值抑制的算法是：①將當前像素的梯度強度與沿正負梯度方向上的兩個像素進行比較；②如果當前像素的梯度強度與另外兩個像素相比最大，則該像素點保留為邊緣點，否則該像素點將被抑制。最后，使用雙閾值處理和連接分析，從邊緣像素中鎖定待提取的水邊線，并使用提取工具導(dǎo)出水邊線邊緣的像素坐標，提取效果見圖5。

圖5 水邊線提取Fig.5 Shoreline extraction

1.2.3 轉(zhuǎn)換坐標

圖像拍攝過程中可能會產(chǎn)生畸變，通過相機標定可消除畸變［14］，并建立起圖像坐標與真實坐標的對應(yīng)關(guān)系。根據(jù)解析攝影測量中的經(jīng)典共線方程，可得到像點坐標與海灘坐標的對應(yīng)關(guān)系如下：

式中：（x，y）為像點的像素坐標；（X，Y，Z）像點對應(yīng)的地面點的實地坐標；（XS，YS，ZS）為拍攝終端的實地坐標；f為拍攝終端主距；a1、a2、a3、b1、b2、b3、c1、c2、c3為外方位元素。

使用6 個高程控制點的實測數(shù)據(jù)，采用最小二乘法平差計算求解得到外方位元素，完成相機標定。水邊線的Z坐標即為圖像拍攝時刻的潮位，將模型計算的潮位數(shù)據(jù)結(jié)合水邊線像素坐標代入已標定的共線方程中，即可輸出水邊線的（X，Y，Z）坐標。

1.2.4 反演地形

將圖片提取出的多個時刻的水邊線視作等高線［15］，然后使用克里金插值算法對多時相的水邊線進行插值，即可得到海灘潮間帶的DEM。

2 模型驗證

采用實測資料對數(shù)學模型進行了驗證，潮位驗證資料為2016 年7 月19～23 日秦皇島海洋站實測數(shù)據(jù)［16］，波浪驗證資料為2016 年7 月19～21 日平水橋波浪監(jiān)測浮標實測波浪數(shù)據(jù)，地形驗證資料為RTKDGPS測量的2016年7月19日地形數(shù)據(jù)，其中波浪和地形資料由河北省地礦局第八地質(zhì)大隊提供。經(jīng)驗證，潮位和波浪的模型計算結(jié)果在大小和相位上與實測值均比較吻合，如圖6所示。7月20日14：00～7月21日3：00期間增水均在0.4m以上，風暴潮增水峰值達0.55m，與實測值0.59m［10］相差0.04m，相對誤差為6.78%；波高誤差均值為0.15 m，波向相對誤差為20.93%［17］。

圖6 潮位及波浪驗證Fig.6 Verification of tide level and wave

在潮間帶內(nèi)測量了36 個高程控制點以驗證反演地形的精度，部分驗證結(jié)果如圖7所示，高程絕對誤差集中在0～0.30 m，平均值為0.12 m，與其他已有研究的誤差相當［8，18］。反演基于同一套轉(zhuǎn)坐標參數(shù)，在評價潮間帶地形變化時候，部分誤差會相互抵消，因此將該地形反演結(jié)果用于評價海灘潮間帶的地貌演變時，具有較好的可靠性。

圖7 反演地形的誤差Fig.7 Deviation of topography extraction

3 結(jié)果與討論

3.1 動力變化

圖8為Argus系統(tǒng)3號攝像機記錄下的2016年7月20 日海灘圖像，此區(qū)域為海灘的平直岸段，且坡度較緩，水位增加對岸灘影響較大?？梢钥闯?，在風暴潮的襲擊下，3號攝像機監(jiān)視的區(qū)域顯示出不同的狀態(tài)。13：30水位明顯較平時高；17：00海水漫上岸灘；18：30 海水進一步蔓延，對附近房屋造成破壞；19：00，水位繼續(xù)增高。

圖8 Argus3號攝像頭下的2016年“720”風暴潮過程Fig.8 Process of 2016 "720" storm surge captured by the Argus camera No.3

圖9為風暴潮期間波浪和潮位過程，有效波高從7 月19 日18 時的0.63 m 開始增加，至20 日17 時達到峰值2.22 m，21 日7 時衰減到0.70 m，整個風暴潮過程的平均有效波高為1.47 m。岸線與正北方向順時針成65°～80°角，故波向在65°～280°之間的波浪均可以直接入射海灘。由圖9a可知，風暴潮期間的波浪均可以直接入射海灘，風暴潮對海灘存在持續(xù)且較強的影響。由圖9b可知，風暴潮發(fā)生時平水橋海灘正處于大潮潮情，最大潮差為1.45 m，高于全年平均潮差0.74 m，最大潮位出現(xiàn)在20 日21 時，與有效波高峰值出現(xiàn)時刻接近，增強了對海灘的影響，最大風暴增水達0.55 m。

圖9 2016年“720”風暴潮期間波況、潮位及風暴增水Fig.9 Wave, tidal level and storm surge elevation during the"720"in 2016

3.2 海灘地形演變

3.2.1 風暴期海灘沖淤演變

圖10a、圖10b分別為視頻反演的2016年7月19日、7 月22 日的潮間帶數(shù)字高程模型DEM，風暴潮對潮間帶地形產(chǎn)生了重大影響，促使潮間帶坡度趨于平緩，泥沙自灘肩向水下輸移。進一步地，將兩天的地形相減，得到風暴潮前后海灘地形沖淤變化，其中正值為淤積，負值為沖刷，如圖10c所示。風暴潮期間，海灘沉積物整體呈離岸輸運的趨勢。海岸線（平均大潮高潮線）高程以上，海灘處于沖刷狀態(tài)，潮間帶處于淤積狀態(tài)。大風浪導(dǎo)致岸灘及高潮帶泥沙掀揚侵蝕，而在中、低潮帶區(qū)域落淤，且自岸向海方向淤積厚度增加，低潮帶最大淤積厚度達0.6 m。

圖10 2016 年7 月19、22 日潮間帶DEM 及7 月19~22 日地形變化Fig.10 Topography of intertidal zone and its differ?ence between July 19 and 22,2016

潮間帶沖淤分布存在明顯的空間差異。海灘中部受風暴潮影響最為劇烈，其坡度比兩側(cè)坡度更為陡峭，且無岬角等地形掩護，根據(jù)地形動力學的觀點，更陡的坡度對波能耗散性更弱，海灘受到風暴潮的影響也更為強烈［19］。海灘西側(cè)處于平直岸段，地形變化最為劇烈；海岸線以上出現(xiàn)強烈的侵蝕，最大沖刷深度達0.6 m，潮間帶淤積明顯。風暴潮期間波浪破碎，波生流方向自東向西，主導(dǎo)泥沙向西輸移。平水橋海灘東側(cè)，由于天然岬角及NE-SW向的人工潛堤岬頭的遮蔽，海灘對風暴潮響應(yīng)較弱，潮間帶地形變化不明顯，但沖淤變化情況較為復(fù)雜，整體以淤積為主，自東向西淤積厚度增加，靠近岬角側(cè)存在部分區(qū)域侵蝕，天然岬角和人工潛堤岬頭對岸段起到一定保護作用，提高了抵御風暴能力。

3.2.2 風暴期岸線變化

“720”風暴潮在短期內(nèi)對海灘輸入了大量能量，造成海灘地形的劇烈變化。圖11為平水橋海灘2016年7月18日、22日的同一潮高下4、5號攝像頭的實景照片。通過對比風暴潮前后的實景照片，發(fā)現(xiàn)灘肩上存在較為明顯的沖刷，砂質(zhì)粗化，如方框區(qū)內(nèi)的灘面沉積物明顯變粗；在同一潮高下的海灘出露面積明顯減小，表明潮間帶淤積，印證了風暴潮造成的灘肩及高潮帶侵蝕，海灘整體坡度趨于平緩，岸灘及高潮帶泥沙向海輸移。對圖10海灘2016年7月19日和22日地形進行處理，得到對應(yīng)的海岸線位置，如圖12所示。風暴潮發(fā)生前（7月19日），灘肩寬度為130.80 m，風暴潮之后（7月22日）為135.14 m，風暴潮導(dǎo)致灘肩沖刷，灘肩高程降低，沖刷下的泥沙落淤在潮間帶，導(dǎo)致海岸線向海推進4.34m。

圖11 2016年7月18日和7月22日同一潮高下海灘圖像Fig.11 Images of beach at the same tidal level on July 18 and July 22,2016

圖12 2016年7月19日和7月22日高潮線位置Fig.12 Location of high tide line on July 19 and July 22,2016

3.2.3 恢復(fù)期海灘演變

風暴潮過后的一段時間里，海灘處于恢復(fù)階段。為更詳細地了解恢復(fù)期的海灘地形演變規(guī)律，對風暴潮后14d（8 月3 日）的20 張Argus 圖片進行解譯，提取水邊線，進而插值反演得到風暴過后8月3日的地形，如圖13a 所示；將8 月3 日與7 月22 日的地形相減，得到7 月22 日～8 月3 日地形變化，即風暴潮過后14d的海灘地形沖淤變化，如圖13b所示。

由圖13 可知，自7 月22 日～8 月3 日，潮間帶整體以淤積狀態(tài)為主，在空間上將海灘大致分為東、中、西三個區(qū)域，海灘東部受人工潛堤岬頭和天然岬角影響，地形變化較為復(fù)雜，高潮帶輕微沖刷，中低潮帶淤積，淤積厚度小于0.2 m；海灘中部為平直岸段，整體呈上淤下沖的狀態(tài)，高潮帶區(qū)域淤積厚度小于0.1 m，中低潮帶沖刷，沖刷深度小于0.2 m；海灘西側(cè)也為平直岸段，地形演變規(guī)律與海灘中部較為一致，即上淤下沖，但變化程度更為劇烈，高潮帶最大淤積厚度達0.3 m，低潮帶最大沖刷深度達0.5 m。

綜觀海灘整體，發(fā)現(xiàn)在7月22日～8月3日的潮間帶演變中，平直岸段與岬角掩護岸段變化呈現(xiàn)出不同的規(guī)律，岬角掩護岸段高潮帶微沖而中低潮帶淤積，泥沙離岸輸移；平直岸段演變主要受中小恢復(fù)浪的作用，呈現(xiàn)出上淤下沖的態(tài)勢，泥沙由中低潮帶向高潮帶輸移。

4 結(jié)論

本文首先建立了“720”風暴潮期間的水動力數(shù)學模型，復(fù)演了風暴潮過程的動力環(huán)境，然后結(jié)合視頻監(jiān)測，運用計算機視覺技術(shù)，反演出平水橋海灘2016年“720”風暴潮發(fā)生前后的潮間帶地形，結(jié)合動力數(shù)據(jù)分析風暴潮對海灘的影響，得出以下結(jié)論：

（1）“720”風暴潮發(fā)生時平水橋海灘正處于大潮潮情，整個過程自7 月19 日18 時起～7 月21 日7 時止，共歷時約1.5 d，風暴潮期間平均有效波高為1.47 m，最大有效波高為2.22 m，波浪均為入射海灘，最大風暴增水0.55 m。

（2）平水橋海灘對風暴潮過程短時響應(yīng)劇烈。在強浪的影響下，泥沙向西、向海輸移，造成灘肩明顯沖刷，砂質(zhì)粗化，潮間帶整體淤積，從而導(dǎo)致海灘坡度變緩，岸線向海推進4.34 m。受波向、地形等因素影響，海灘對風暴潮的響應(yīng)存在明顯的空間差異，其中以平直岸段的海灘中部對風暴潮響應(yīng)最為劇烈，海灘東部受人工潛堤岬頭和天然岬角影響，地形變化幅度最小。因此可認為人工潛堤岬頭對其影響岸段起到了一定的保護作用，提高了岸段的抗風暴潮能力。

（3）風暴潮后的恢復(fù)期，平水橋海灘平直岸段與岬角影響岸段呈現(xiàn)出不同的變化規(guī)律。岬角掩護岸段高潮帶微沖而中低潮帶淤積，泥沙離岸輸移；平直岸段演變主要受中小恢復(fù)浪的作用，呈現(xiàn)出上淤下沖的態(tài)勢，泥沙由中低潮帶向高潮帶輸移。

作者貢獻聲明：

匡翠萍：論文撰寫及修改。

劉旭：數(shù)據(jù)處理、論文撰寫。

夏子龍：數(shù)據(jù)處理、圖像繪制。

朱磊：現(xiàn)場工作。

叢新：論文修改。