袁春光，王義剛，黃惠明，蔡 輝，陳 橙

(河海大學(xué)海岸災(zāi)害及防護(hù)教育部重點(diǎn)實(shí)驗(yàn)室，江蘇南京 210098)

海嘯是指由海底地震、火山爆發(fā)和海底滑坡、塌陷所產(chǎn)生的具有超大波長(zhǎng)和周期的大洋行波，當(dāng)其接近近岸淺水區(qū)時(shí)，波速變小，波幅陡漲，有時(shí)可達(dá)30米以上，突然形成“水墻”，瞬時(shí)侵入濱海陸地，沖毀或卷去沿海建筑和人畜［1-2］。其中，海中地震引發(fā)的地殼垂向錯(cuò)動(dòng)(“傾滑型”)可造成大面積水體突然升降，從而形成一系列波長(zhǎng)數(shù)十到數(shù)百公里，周期為2～200 min的長(zhǎng)重力波，一般稱(chēng)之為地震海嘯。

地震海嘯通常分為兩種:一種是橫越大洋或從遠(yuǎn)洋傳播來(lái)的海嘯，成為遠(yuǎn)洋海嘯。例如，1960年智利發(fā)生8.7級(jí)特大地震海嘯，在智利沿岸波高達(dá)20.4 m，海嘯傳到夏威夷時(shí)，波高尚超過(guò)11.0 m，在日本沿岸波高仍有6.1 m。另一種是近海海嘯，海嘯生成源地與其造成危害同處一地，所以海嘯到達(dá)沿岸時(shí)間很短，有時(shí)只有幾分鐘或者幾十分鐘，危害巨大。例如，1983年的印尼6.5級(jí)地震引起的海嘯，遇難人數(shù)為3.6萬(wàn)，克拉客托島三分之一沉入海底。

北京時(shí)間2011年3月11日13時(shí)46分，位于日本宮城縣以東太平洋海域(37.7°N，143.0°E)發(fā)成里氏9.0級(jí)地震(如圖1)［3］，震源深度為24.4 km，14 133人遇難，13 346人失蹤，超過(guò)10萬(wàn)人撤離。在地震發(fā)生后的20多個(gè)小時(shí)里，海嘯波先后到達(dá)俄羅斯、中國(guó)、菲律賓、美國(guó)夏威夷和西海岸智利以及南美其他一些國(guó)家和地區(qū)沿海，各地均監(jiān)測(cè)到了明顯的海嘯波。其中，美國(guó)及南美的智利等國(guó)家監(jiān)測(cè)到的最大海嘯波幅達(dá)到150～200 cm。地震發(fā)生后約4小時(shí)海嘯波到達(dá)我國(guó)臺(tái)灣東部沿海，隨后波及我國(guó)大陸東南沿海，有關(guān)部門(mén)預(yù)計(jì)次日凌晨2點(diǎn)到達(dá)上海和江蘇南通。臺(tái)灣、福建、浙江、廣東、江蘇、上海、海南等省市監(jiān)測(cè)到的海嘯波波幅在3～55 cm之間，其中浙江沈家門(mén)和石浦海洋站分別監(jiān)測(cè)到55 cm和52 cm的海嘯波(如圖2、3)［4］。受其影響，國(guó)家海洋預(yù)報(bào)臺(tái)針對(duì)臺(tái)灣地區(qū)和浙江沿岸發(fā)布2期海嘯藍(lán)色警報(bào)，這是我國(guó)開(kāi)展海嘯預(yù)警業(yè)務(wù)以來(lái)首次發(fā)布海嘯警報(bào)。

圖1 日本地震位置Fig.1 Japan earthquake location

圖2 2011年3月11日日本地震海嘯及最大海嘯波幅圖［4］Fig.2 Japan tsunami maximum amplitude distribution on March 11，2011

1 數(shù)據(jù)來(lái)源和處理方法

將在對(duì)江蘇灌河口處(站點(diǎn)位置如圖4)海嘯期間逐分鐘潮位資料的分析的基礎(chǔ)上，討論“311”日本海嘯對(duì)該地區(qū)的影響，并且討論造成這種現(xiàn)象的原因。

波高和周期是刻畫(huà)潮位波動(dòng)的最基本的兩個(gè)特征量。當(dāng)海嘯波趨向近岸淺水時(shí)波形發(fā)生變化，具體表現(xiàn)為波長(zhǎng)變小，波高增加。實(shí)際上，海嘯波高在近岸很大程度上取決于海底地形、坡度和海岸線的形狀與走向。喇叭口或漏斗型的海岸地形有利于波能折射聚集，海嘯將大幅提高。

衛(wèi)星遙感、海上浮標(biāo)和沿海潮位站是獲得海嘯增水?dāng)?shù)據(jù)的主要來(lái)源。由于受到各種天氣和人為等因素的影響，從潮位信號(hào)中精確地提取海嘯波仍然是比較困難的，介紹這方面文章也不多見(jiàn)。于福江認(rèn)為［5］，海嘯影響期間，潮位序列主要表現(xiàn)為較高頻率的海嘯波疊加在長(zhǎng)周期的潮汐信號(hào)上，為了消除高頻噪音，對(duì)潮位信號(hào)進(jìn)行截?cái)嘀芷跒?0～136.53 min的帶通濾波，并且認(rèn)為經(jīng)過(guò)以上處理方法得到的波動(dòng)序列就是“干凈”的海嘯波動(dòng)序列，序列中最大波高即認(rèn)為是海嘯引起的最大增水，以第一個(gè)到達(dá)潮位站的顯著波峰或者波谷為準(zhǔn)作為該站位的海嘯到達(dá)時(shí)間。董非非［6］認(rèn)為，把原始數(shù)據(jù)通過(guò)濾波器，濾掉氣象潮、海洋背景噪音等的影響，得到的潮差(包含潮汐，也可能包含海嘯信息)，扣除潮汐的影響，對(duì)最終得到的數(shù)據(jù)進(jìn)行分析處理，可從周期(100～1 000 s)、波長(zhǎng)等參數(shù)判斷識(shí)別其是否含有海嘯信息。

以上分離海嘯波的方法主要問(wèn)題在于如何確定濾波截?cái)嘀芷诘纳舷藓拖孪?。海嘯波周期范圍一般在2～200 min之內(nèi)，但是海嘯波的波動(dòng)并不會(huì)均勻分布在如此漫長(zhǎng)的周期范圍中，通常只會(huì)聚集在一個(gè)相對(duì)較窄的周期范圍中。如果濾波的截?cái)嘀芷诜秶∮谶@個(gè)周期范圍，則不能分離出全部的海嘯波;如果截?cái)嘀芷诜秶笥谶@個(gè)周期范圍，那么濾波結(jié)果中將摻雜不屬于海嘯波波動(dòng)的“噪聲”，從而造成結(jié)果不準(zhǔn)確。以往對(duì)于海嘯濾波周期范圍的選取，主要依據(jù)個(gè)人經(jīng)驗(yàn)，不同學(xué)者之間取值的差異很大，過(guò)濾出的海嘯波也大小不同，這將影響后續(xù)的研究工作。

圖3 2011年3月11日日本地震海嘯中我國(guó)沿海潮位站最大海嘯波幅［4］Fig.3 Maximum tsunami amplitude in tidal stations of China，March 11，2011

圖4 灌河口測(cè)站位置示意Fig.4 Guanhekou station location

文中將多種信號(hào)分析的方法綜合使用于潮位信號(hào)的分析中，形成一種新的綜合分析海嘯波的模式，更加準(zhǔn)確全面地了解海嘯波的各項(xiàng)物理特征。具體步驟如下:

首先，通過(guò)對(duì)同一潮位站過(guò)往一年369天的資料進(jìn)行潮汐調(diào)和分析，預(yù)報(bào)出該潮位站海嘯發(fā)生時(shí)間段的天文潮位，在不考慮海嘯波與天文潮耦合的情況下，用海嘯期間實(shí)測(cè)潮位減去預(yù)報(bào)潮位，這個(gè)差值即認(rèn)為是由海嘯和其他天氣條件共同引起的增水。這一步主要是為了去除潮位資料中的長(zhǎng)周期天文潮波。

其次，對(duì)這個(gè)差值序列信號(hào)以Morlet小波為小波母函數(shù)，應(yīng)用連續(xù)小波分析方法，進(jìn)行小波能量譜分析［7］。小波變換可以在波動(dòng)序列的任意時(shí)刻對(duì)信號(hào)成分進(jìn)行局部化分析，提取差值序列在時(shí)域和頻率域上的局部性特征。小波能量譜分析是基于小波變換理論基礎(chǔ)建立的，具有同時(shí)在時(shí)間域和頻率域內(nèi)識(shí)別振動(dòng)強(qiáng)度的功能，并且可以聚焦到信號(hào)細(xì)節(jié)的分析。通過(guò)這種方法可以找出潮位信號(hào)中和海嘯波動(dòng)周期范圍相似，持續(xù)時(shí)間相當(dāng)，比較可疑的異常強(qiáng)烈波動(dòng)。由于海嘯發(fā)生的概率很低，海嘯引起的水位波動(dòng)只是偶爾發(fā)生，所以通過(guò)小波能量譜分析找出的異常波動(dòng)必須具備嚴(yán)格的罕見(jiàn)性才有可能是海嘯波。同時(shí)還應(yīng)該根據(jù)潮位站位置與海嘯發(fā)生位置之間的距離以及淺水波的傳播速度，初步估算海嘯波動(dòng)到達(dá)潮位站的時(shí)間，并與可疑波動(dòng)發(fā)生時(shí)刻進(jìn)行對(duì)比。如果估算海嘯波到達(dá)時(shí)間與可疑波動(dòng)發(fā)生時(shí)間相近，則大大增加了可疑波動(dòng)就是海嘯波的可能;相反，如果估算海嘯波到達(dá)時(shí)間與可疑波動(dòng)發(fā)生時(shí)間相距甚遠(yuǎn)，那么可疑波動(dòng)很可能是由于其他因素導(dǎo)致。

再次，將滿足以上條件的可疑波動(dòng)找出，確定周期范圍之后，利用2階或者更高階數(shù)的Butter帶通濾波器對(duì)差值序列在該周期范圍內(nèi)進(jìn)行濾波。Butter濾波器的特點(diǎn)是通帶處幅值特性平坦，可以通過(guò)提高階數(shù)的方法提高逼近精度，但是同時(shí)也增大了計(jì)算量。這樣濾出的波動(dòng)即認(rèn)為是海嘯波，序列中最大波高即認(rèn)為是海嘯引起的最大增水。

最后，通過(guò)對(duì)濾波序列進(jìn)行功率譜分析，進(jìn)一步得出各種周期的波動(dòng)對(duì)整體波動(dòng)能量的貢獻(xiàn)，從而更精確和全面的了解各種波動(dòng)成分的運(yùn)動(dòng)特點(diǎn)。文中選用的功率譜分析方法是，先將過(guò)濾出的海嘯波序列由快速Fourier變換獲得能譜密度，再采用Welch算法［8］分成不重疊的6部分，然后取Fourier系數(shù)平均得到比較光滑的譜線。

采用以上海嘯分離方法的優(yōu)點(diǎn)是:1)能夠比較針對(duì)性的對(duì)海嘯波進(jìn)行分離過(guò)濾，尤其是對(duì)受海嘯作用影響不明顯的潮位站，這種方法更加精確可行;2)采用多種方法相結(jié)合的方式對(duì)海嘯波進(jìn)行分析，可以從多個(gè)角度更加全面的了解海嘯波波動(dòng)的物理特性。

對(duì)于灌河口潮位站，海嘯期間實(shí)測(cè)潮位摘取時(shí)間范圍為2011年3月11日13時(shí)46分至3月14日0時(shí)，采樣頻率為1 min。用當(dāng)?shù)?009年8月至2010年8月369天的潮位資料進(jìn)行潮汐調(diào)和分析，預(yù)報(bào)出海嘯摘取期間該地區(qū)天文潮，如圖5。將實(shí)測(cè)潮位值減去天文潮預(yù)報(bào)值，差值表示海嘯或者其他天氣因素引起的增水，如圖6。從圖中可以看出預(yù)報(bào)結(jié)果與實(shí)測(cè)值的誤差大部分在20 cm以下，且波動(dòng)相位吻合良好。

圖5 灌河口海嘯期間天文潮預(yù)報(bào)與實(shí)測(cè)潮位值的比較Fig.5 Comparison of astronomical and measured tide in Guanhekou during tsunami period

圖6 海嘯期間灌河口實(shí)測(cè)與預(yù)報(bào)的差值Fig.6 Difference between measured and predicted tide in Guanhekou during tsunami period

2 數(shù)據(jù)處理結(jié)果分析

通過(guò)小波能量譜圖(如圖7)可以發(fā)現(xiàn)，在整個(gè)3月份的圖譜中，僅在3月11日8時(shí)之后有一股周期范圍大致為25～100 min之間，波動(dòng)能量比較強(qiáng)烈并且持續(xù)時(shí)間相對(duì)較長(zhǎng)的波動(dòng)進(jìn)入灌河口站，而且在其余將近一個(gè)月的時(shí)間內(nèi)，基本上沒(méi)有周期在2～200 min范圍內(nèi)同等連續(xù)并且劇烈的波動(dòng)，這與日本地震海嘯發(fā)生的時(shí)間比較吻合，同時(shí)也在海嘯波的周期范圍之內(nèi)。于是摘取3月11日13時(shí)46分至3月14日0時(shí)這一時(shí)段潮位進(jìn)行進(jìn)一步分析。分析表明，有一列周期范圍在37～98.5 min的水位波動(dòng)從3月12日早上6時(shí)26分左右開(kāi)始影響灌河口站，并于3月13日凌晨2時(shí)26分左右結(jié)束，該周期范圍恰好在海嘯波的波動(dòng)周期范圍之內(nèi)，也是初步懷疑的“可疑波動(dòng)”。同時(shí)，對(duì)全年其他月份的潮位資料進(jìn)行相同的分析，沒(méi)有發(fā)現(xiàn)頻率在海嘯波周期范圍之內(nèi)，能量和持續(xù)時(shí)間和“疑似海嘯波”同等量級(jí)的波動(dòng)，這就愈發(fā)增大了這列波是海嘯波的可能。

從時(shí)空角度上分析，位于江蘇北部的灌河口站與日本地震震源之間的直線距離約為2 131 km，海嘯長(zhǎng)波的傳播速度為，假設(shè)平均水深為120 m，由計(jì)算可得海嘯波大約在地震發(fā)生之后1 035 min(3月12日7時(shí)1分)抵達(dá)灌河口，這與上文分析得到的開(kāi)始影響時(shí)間比較一致。因此，可以基本確定圖7中的波動(dòng)就是由于日本地震引起，傳播到灌河口潮位站的海嘯波。

以帶通周期為37～98.5 min對(duì)潮位信號(hào)進(jìn)行濾波，如圖8所示，橫坐標(biāo)是從海嘯發(fā)生后開(kāi)始計(jì)算的時(shí)間，單位為分鐘。從圖中可以看出，海嘯波最大波高為5.05 cm，這與國(guó)家海洋局《2011年中國(guó)海洋災(zāi)害公報(bào)》中連云港潮位站監(jiān)測(cè)到的海嘯波高5 cm十分吻合(連云港至灌河口距離僅40 km左右)。同時(shí)，分別用2階Butter帶通濾波器和5階Butter帶通濾波器對(duì)潮位信號(hào)進(jìn)行濾波，兩者結(jié)果幾乎沒(méi)有差異，說(shuō)明應(yīng)用2階Butter濾波器可以滿足精度要求。用于福江的海嘯濾波方法對(duì)潮位進(jìn)行濾波，與文中方法濾波的波形進(jìn)行對(duì)比發(fā)現(xiàn)，兩列波動(dòng)總體趨勢(shì)很接近，都是先增大后減小。但是于福江方法的濾波中明顯摻雜了噪聲，波形參差不齊，而且最大波高為6.05 cm，比這里濾波的最大波高增大了20%，可見(jiàn)文中使用的過(guò)濾海嘯波方法更加精確。進(jìn)一步對(duì)濾波序列進(jìn)行功率譜分析(如圖9)可以看出，海嘯波波能主要分布在周期為56～83 min的波動(dòng)中，最大峰值出現(xiàn)在69.4 min，即對(duì)應(yīng)海嘯主周期。

圖7 小波能量譜分析Fig.7 Wavelet power spectrum

圖8 海嘯濾波Fig.8 Tsunami filter

圖9 海嘯波功率譜Fig.9 Tsunami power spectrum

3 結(jié)果討論

一般地震海嘯波的最大波峰大多出現(xiàn)在第一次波動(dòng)，但是對(duì)比圖8，灌河口海嘯波的最大波峰出現(xiàn)在海嘯持續(xù)時(shí)間的中部，而且波高恰好為前半段平均波高的2倍左右。為了分析這種現(xiàn)象，采用美國(guó)Cornell大學(xué)基于長(zhǎng)波方程淺水理論編譯的COMCOT海嘯模型［9］，進(jìn)行簡(jiǎn)化的海嘯傳播數(shù)值模擬?；A(chǔ)地理信息數(shù)據(jù)采用NGDC的ETOPO1海底地形資料，模擬范圍為118°～150°E ，20°～45°N，球面坐標(biāo)系網(wǎng)格采取5 min的正方形網(wǎng)格。同時(shí)引用GCMT震源機(jī)制解作為震源參數(shù)，具體地震參數(shù)如表1［3］。應(yīng)用Okada理論［10］模型進(jìn)行海嘯初始場(chǎng)計(jì)算，模擬時(shí)間為36 h。

海嘯波數(shù)值模擬傳播結(jié)果見(jiàn)圖10。

圖10 海嘯發(fā)生后不同時(shí)刻的海嘯波模擬結(jié)果Fig.10 Tsunami wave simulation results after breaking out

表1 海嘯地震參數(shù)表Tab.1 Tsunami seismic parameter

從圖10可以看出，在模擬的海嘯波向東海大陸架傳播的過(guò)程中，海嘯波一方面向長(zhǎng)江口及南方的浙江福建沿岸逼近;另一方面也會(huì)向北方的黃海和渤海傳播。向北方黃海傳播的海嘯波會(huì)先抵達(dá)山東半島，然后分裂成兩部分，一部分繼續(xù)北上，向渤海挺近;另一部分則受到山東半島阻擋反射，在山東半島南側(cè)和江蘇沿岸形成第一波海嘯波動(dòng)。由于中國(guó)沿海近似一個(gè)圓弧形狀，向長(zhǎng)江口傳播的海嘯波抵達(dá)長(zhǎng)江口后分裂成兩部分，一部分沿江蘇沿岸北上;另一部分沿浙江沿岸南下。這樣一來(lái)，從長(zhǎng)江口繼續(xù)北上的海嘯波將疊加到江蘇沿岸第一波海嘯波之上，從而最大海嘯波峰值出現(xiàn)在兩股海嘯波疊加的地方，卻不出現(xiàn)在第一波。同時(shí)，由山東半島和江蘇沿岸組成的“喇叭口型”的海岸地形也有利于海嘯波在此疊加。如圖11，與實(shí)測(cè)資料對(duì)比發(fā)現(xiàn)，海嘯波數(shù)值模擬比實(shí)測(cè)數(shù)據(jù)有所提前，而且波形也有所不同。引起這些誤差可能的原因:1)地形資料精度有限，不能完全符合實(shí)際情況，尤其在近岸，這個(gè)問(wèn)題十分突出。2)受水深和岸灘資料不夠準(zhǔn)確影響，網(wǎng)格尺寸比較粗大，只能選擇距離灌河口測(cè)站比較接近的網(wǎng)格節(jié)點(diǎn)水位值代表實(shí)際測(cè)站。3)實(shí)際的海域底部摩阻應(yīng)該用一個(gè)曼寧系數(shù)場(chǎng)來(lái)表示，由于缺少資料這里只用常數(shù)0.013代表。對(duì)于沖繩海槽以東的太平洋海域，曼寧系數(shù)應(yīng)該更小，而黃海海域應(yīng)該更大一些。4)實(shí)際海底地震引發(fā)的水位波動(dòng)比較復(fù)雜，這里用Okada方法進(jìn)行簡(jiǎn)化只是一種近似。同時(shí)，海底地震參數(shù)的測(cè)量也存在技術(shù)障礙，不同研究機(jī)構(gòu)給出的地震參數(shù)各有不同，也直接影響海嘯模擬的精度。其中原因1)和2)主要影響海嘯波傳播速度;3)則直接影響海嘯在傳播中的波高變化。雖然以上數(shù)值模擬存在精度問(wèn)題，仍然可以明顯看出海嘯波動(dòng)先逐漸增大，然后逐漸減小的規(guī)律，這也證明了海嘯波動(dòng)在這里發(fā)生了疊加。

圖11 灌河口處海嘯波高數(shù)值模擬結(jié)果Fig.11 Numerical simulation of tsunami wave in Guanhekou

依據(jù)以上觀點(diǎn)，實(shí)際海嘯首波引起灌河口站的波動(dòng)低于5.05 cm，大約只有2.5 cm左右，5.05 cm是后期由于岸線形狀引起波動(dòng)疊加之后的結(jié)果。

4 結(jié)語(yǔ)

1)通過(guò)對(duì)灌河口實(shí)測(cè)潮位資料的分析得出，2011年3月11日發(fā)生在日本的地震海嘯于3月12日早上06∶26左右開(kāi)始影響該站，3月13日凌晨02∶26左右結(jié)束影響。海嘯波周期范圍為37～98.5 min，最大波高為5.05 cm，主周期為69.4 min左右。

2)與其它從潮位中過(guò)濾海嘯波的方法對(duì)比分析，這里使用的濾波方法更加精確，并且更加全面、直觀地表現(xiàn)出海嘯波波動(dòng)的各項(xiàng)物理特性。

3)利用簡(jiǎn)化的海嘯波數(shù)值模型從一定程度上模擬了海嘯波在中國(guó)沿海的傳播。模擬結(jié)果分析得出，由太平洋經(jīng)過(guò)沖繩海槽向中國(guó)沿海傳播的海嘯波會(huì)一方面向東南沿海逼近，另一方面也向黃海海域挺進(jìn)。由于地形的影響，海嘯波產(chǎn)生了反射，這樣最大波峰將出現(xiàn)在兩股波疊加后的時(shí)刻。這就是灌河口站濾出的海嘯波最大波峰不在第一波的主要原因。

4)由數(shù)值模擬可知，灌河口站由海嘯直接引起的海嘯波大約只有2.5 cm左右，但是海嘯波疊加之后，最大波高為5.05 cm，是之前的2倍。雖然海嘯波高量值不大，但是海嘯預(yù)報(bào)中應(yīng)該全面考慮地理位置特點(diǎn)，綜合分析海嘯波運(yùn)動(dòng)的反射疊加等特征，從而提高預(yù)報(bào)質(zhì)量。

［1］ 于福江，葉 琳，王喜年.1994年發(fā)生在臺(tái)灣海峽的一次地震海嘯的數(shù)值模擬［J］.海洋學(xué)報(bào)，2001，6(23):30-39.

［2］ 葉 琳，于福江，吳 瑋.我國(guó)海嘯災(zāi)害及預(yù)警現(xiàn)狀與建議［J］.海洋預(yù)報(bào)，2005，22(增刊):147-157.

［3］ 王培濤，于福江，趙聯(lián)大，等.2011年3月11日日本地震海嘯越洋傳播及對(duì)中國(guó)影像的數(shù)值分析［J］.地球物理學(xué)報(bào)，2012，9(55):3088-3095.

［4］ 國(guó)家海洋局2011年海洋災(zāi)害公報(bào)［ED/OL］.http://www.soa.gov.cn/soa/index.html，2013-05-10-18:00.

［5］ 于福江，原 野，趙聯(lián)大，等.2010年2月27日智利8.8級(jí)地震海嘯對(duì)我國(guó)影響分析［J］.科學(xué)通報(bào)，2010，55(1):1-8.

［6］ 董非非.海嘯在東海大陸架傳播的研究與識(shí)別［D］.蘭州:中國(guó)地震局蘭州地震研究所，2009.

［7］ 謝 利.小波變換及小波能量譜在多尺度分析中的應(yīng)用［J］.中山大學(xué)研究生學(xué)刊:自然科學(xué)版，1998，19(增刊):101-102.

［8］ Welch P D.The use of fast Fourier transform for the estimation of power spectra-A method based on time averaging over short，modified periodograms［J］.IEEE Trans Audio Electroacoust，1967，AU-15:70-73.

［9］ 姚 遠(yuǎn)，蔡樹(shù)群，王盛安.海嘯波數(shù)值模擬的研究現(xiàn)狀［J］.海洋科學(xué)進(jìn)展，2007，25(4):487-494.

［10］ Okada Y.Surface deformation due to shear and tensile faults in a half-space［J］.Bulletin of the Seismological Society of America，1985，75:1135-1154.