基于希爾伯特-黃變換的潮汐分析方法研究

張?亮，張佳麗，張學(xué)峰，張安民

基于希爾伯特-黃變換的潮汐分析方法研究

張?亮，張佳麗，張學(xué)峰，張安民

(天津大學(xué)海洋科學(xué)與技術(shù)學(xué)院，天津 300072)

經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)是希爾伯特-黃變換(HHT)的核心組成部分．在對(duì)傳統(tǒng)HHT方法研究與分析的基礎(chǔ)上，發(fā)展了改進(jìn)的集合經(jīng)驗(yàn)?zāi)B(tài)分解的方法(MEEMD)對(duì)海洋潮汐進(jìn)行分析．MEEMD方法根據(jù)水位數(shù)據(jù)序列的特點(diǎn)，通過添加成對(duì)白噪聲和計(jì)算排列熵的方式來剔除分解過程中出現(xiàn)的虛假本征模態(tài)函數(shù)(IMF)，從而抑制了模態(tài)混疊現(xiàn)象．?dāng)?shù)值模擬和實(shí)測(cè)水位數(shù)據(jù)實(shí)驗(yàn)表明，MEEMD方法可以在有效抑制模態(tài)混疊現(xiàn)象的基礎(chǔ)上，對(duì)含有噪聲的原始水位序列實(shí)現(xiàn)自適應(yīng)濾波．MEEMD方法分解得到的高頻IMF1和IMF2分量占原始水位數(shù)據(jù)總能量的90%以上，IMF1和IMF2分量的平均周期和平均振幅分別與調(diào)和分析得到的半日分潮與全日分潮結(jié)果基本一致．通過分析認(rèn)為，IMF1和IMF2之和是主要半日分潮、全日分潮以及淺海分潮的疊加．傳統(tǒng)的潮汐調(diào)和分析方法有賴于長(zhǎng)期的水位記錄，且需要水位數(shù)據(jù)的一些先驗(yàn)知識(shí)，涉及到的天文知識(shí)與數(shù)學(xué)計(jì)算也較為復(fù)雜．準(zhǔn)調(diào)和分析可以處理短期的水位數(shù)據(jù)，但仍然涉及天文知識(shí)與復(fù)雜數(shù)學(xué)計(jì)算，對(duì)數(shù)據(jù)的完整性也要求較高．當(dāng)收集到的水位數(shù)據(jù)因?yàn)榇嬖谠肼暬蛘咝畔⑷笔В蛔阋赃M(jìn)行傳統(tǒng)的調(diào)和或準(zhǔn)調(diào)和分析時(shí)，可以使用HHT方法進(jìn)行有效的補(bǔ)充分析，以獲得主要半日分潮、全日分潮和余水位的幅值信息，因此，該方法可以作為一種很好的潮汐分析的補(bǔ)充方法．

改進(jìn)的集合經(jīng)驗(yàn)?zāi)B(tài)分解；希爾伯特-黃變換；濾波；分潮

海洋潮汐是海水受到月球和太陽的吸引力作用而產(chǎn)生的一種周期性升降運(yùn)動(dòng)．產(chǎn)生潮汐現(xiàn)象的主要原因是由于地球上各點(diǎn)距離月球的相對(duì)位置不同，所受到的引力也不同，這種不同導(dǎo)致地球上的海水做相對(duì)運(yùn)動(dòng)，引起海面升降．而實(shí)際的潮汐受到陸地、海底地形起伏變化，海底摩擦以及地球自轉(zhuǎn)等影響是十分復(fù)雜的，呈現(xiàn)出非平穩(wěn)、非線性的特征．傳統(tǒng)的潮汐分析方法有調(diào)和分析和準(zhǔn)調(diào)和分析．潮汐的調(diào)和分析方法是建立在分潮的概念上，將潮汐看成是以不同頻率傳播的各種分潮疊加產(chǎn)生的現(xiàn)象，目的是求出各個(gè)分潮的振幅和遲角[1-2]．潮汐調(diào)和分析根據(jù)過去一段時(shí)間的潮汐觀測(cè)資料，計(jì)算各個(gè)分潮的調(diào)和常數(shù)，用于潮汐預(yù)報(bào)以及計(jì)算各種基準(zhǔn)面等，其中涉及天文知識(shí)與數(shù)學(xué)計(jì)算，整個(gè)過程較為復(fù)雜．潮汐準(zhǔn)調(diào)和分析適用于手工計(jì)算，通過對(duì)一天或幾天的短期觀測(cè)數(shù)據(jù)計(jì)算獲得準(zhǔn)調(diào)和分潮．準(zhǔn)調(diào)和分潮是以一個(gè)分潮為主，其余分潮的調(diào)和常數(shù)與這個(gè)主要分潮建立一定關(guān)系后合并為一項(xiàng)，此項(xiàng)的振幅和角速率不再是常量，而是隨時(shí)間緩慢變化[3]．希爾伯特-黃變換(Hilbert-Huang transform，HHT)是Haung等于1998年提出的一種分析非線性、非平穩(wěn)信號(hào)的方法[4-6]. HHT作為一種新的信號(hào)分析理論，其關(guān)鍵部分為經(jīng)驗(yàn)?zāi)B(tài)分解(empirical mode decomposition，EMD).大部分通過EMD分解得到的固有模態(tài)函數(shù)(instrinsic mode function，IMF)可以反映真實(shí)的物理過程，同時(shí)也有大量實(shí)驗(yàn)表明，EMD方法是目前提取信號(hào)趨勢(shì)或均值的最好方法．在海洋科學(xué)研究中，EMD方法已應(yīng)用于海平面變化趨勢(shì)、風(fēng)暴潮和海溫異常預(yù)測(cè)等方面的研究[7-10]．Cheng等[11]使用EMD的方法對(duì)潮汐的振幅進(jìn)行分析得到海平面上升和氣候變化的關(guān)系；Pan等[12]應(yīng)用EMD方法對(duì)具有非平穩(wěn)、非線性特征的河流潮汐進(jìn)行分析并驗(yàn)證了方法的有效性．但EMD方法本身存在一些問題，其中一個(gè)主要問題是模態(tài)混疊，針對(duì)這個(gè)問題，很多學(xué)者提出了包括集合經(jīng)驗(yàn)?zāi)B(tài)分解(ensemble empirical mode decomposition，EEMD)[13]、補(bǔ)充的EEMD(comple-mentary EEMD，CEEMD)[14]、改進(jìn)的EEMD(modi-fied EEMD，MEEMD)等解決方法[15]. 另外，黃大吉等[16]針對(duì)HHT存在的端點(diǎn)處理問題，提出了包絡(luò)的極值延拓法和鏡像閉合延拓法，兩種方法都能很好地處理端點(diǎn)問題，使得EMD分解得到合理的各個(gè)IMF模態(tài)．在前述文獻(xiàn)研究的基礎(chǔ)上，本文提出了采用基于MEEMD的希爾伯特-黃變換方法對(duì)潮汐信號(hào)進(jìn)行分析，有效抑制了模態(tài)混疊問題并提高了頻率分辨率，同時(shí)該方法還具有很好的自適應(yīng)濾波效果，是潮汐分析的一種有效方法．

1?HHT方法簡(jiǎn)介

1.1?MEEMD方法

EMD方法中出現(xiàn)的模態(tài)混疊現(xiàn)象是指在一個(gè)IMF中包含差異極大的特征時(shí)間尺度，或者相近的特征時(shí)間尺度分布在不同的IMF中，造成波形混疊，難以辨認(rèn)[17-18]．為了抑制模態(tài)混疊現(xiàn)象，提出了EEMD方法和CEEMD方法．CEEMD方法是EEMD方法的改進(jìn)，該方法通過在原始信號(hào)中添加成對(duì)幅值有限且均值為零的白噪聲，再通過對(duì)分解結(jié)果進(jìn)行總體平均來消除白噪聲的影響，進(jìn)而達(dá)到抑制模態(tài)混疊現(xiàn)象的目的．但是分解結(jié)果易受到添加白噪聲參數(shù)的影響，若選擇不當(dāng)會(huì)產(chǎn)生虛假IMF分量并破壞分解的自適應(yīng)性．MEEMD方法為了彌補(bǔ)上述CEEMD方法的不足，提出了用計(jì)算排列熵(permutation entropy，PE)[19]的方法來剔除信號(hào)中的虛假分量，篩選出主要IMF分量．排列熵是一種能有效度量信號(hào)復(fù)雜程度的物理量，信號(hào)越復(fù)雜，排列熵值越高，反之，則熵值越低[20]．而含有噪聲的信號(hào)往往具有較高的隨機(jī)性，信號(hào)復(fù)雜程度高，排列熵值相比于不含噪聲的信號(hào)更高．

MEEMD算法具體步驟如下所述．

式中為待分析信號(hào)()的采樣點(diǎn)個(gè)數(shù)．

或

通過MEEMD方法，將原始數(shù)據(jù)序列分解為個(gè)IMF分量及余項(xiàng)之和，余項(xiàng)代表數(shù)據(jù)的平均趨勢(shì)．

1.2?希爾伯特變換

希爾伯特變換是在MEEMD分解的基礎(chǔ)上進(jìn)行的，利用希爾伯特變換求解每階IMF的瞬時(shí)頻率，從而得到信號(hào)的振幅-時(shí)間-頻率表示，即希爾伯特譜．通過計(jì)算希爾伯特譜對(duì)整個(gè)周期的積分得到希爾伯特邊際譜．

將式(5)分解得到的各IMF分量進(jìn)行希爾伯特變換，得到

經(jīng)過計(jì)算獲得信號(hào)的希爾伯特譜可表示為

希爾伯特譜是瞬時(shí)振幅在頻率-時(shí)間平面上的分布，具有良好的時(shí)頻聚集性[5]．

式中為周期．邊際譜描述了信號(hào)在整個(gè)時(shí)頻內(nèi)，同一頻率幅值或能量疊加值在頻域的分布情況．

瞬時(shí)能量表示每個(gè)時(shí)刻在整個(gè)頻域所積累的能量．

2?數(shù)值模擬分析與方法有效性驗(yàn)證

2.1?數(shù)值模擬水位序列的MEEMD與EMD分解

為了驗(yàn)證HHT方法的有效性，采用數(shù)值模擬的方式來進(jìn)行驗(yàn)證．潮汐可以看作是多個(gè)分潮的疊加，調(diào)和分析方法把潮位看成是多個(gè)分潮余弦振動(dòng)之和，根據(jù)最小二乘或波譜分析原理，由實(shí)測(cè)的潮位數(shù)據(jù)計(jì)算出各分潮的平均振幅和遲角，其中分潮的平均振幅和遲角為分潮的潮汐調(diào)和常數(shù)[2]．根據(jù)潮汐調(diào)和分析理論，不考慮非天文潮位，設(shè)潮汐為個(gè)分潮疊加，得到的潮高表達(dá)式為

準(zhǔn)調(diào)和分析是通過引入不同分潮之間的關(guān)系，將不同分潮合并為一項(xiàng)，獲得振幅和角速率隨時(shí)間緩慢變化的準(zhǔn)調(diào)和分潮[3]．這樣，潮位的準(zhǔn)調(diào)和分潮表達(dá)式可以表示為

式中：0代表平均海面高度加上長(zhǎng)周期分潮；和是和天文情況有關(guān)的變量，分別叫作準(zhǔn)調(diào)和分潮的振幅系數(shù)和遲角訂正；代表準(zhǔn)調(diào)和分潮．

這里根據(jù)資料來源處的驗(yàn)潮站數(shù)據(jù)，采用固定的振幅和遲角進(jìn)行潮汐模擬計(jì)算，選取M2、S2、K1和Mf 四個(gè)分潮，其中M2分潮振幅為10m，其余分潮采用相對(duì)振幅．輸入?yún)?shù)如表1所示．

由于水位觀測(cè)序列信息中除了水位信息外還包含了觀測(cè)噪聲等干擾，例如涌浪噪聲會(huì)使水位序列曲線受到高頻振動(dòng)的影響，在水位觀測(cè)時(shí)刻有強(qiáng)烈的上下振動(dòng)，所以在數(shù)值模擬的水位序列中加入了功率為0.05的高斯白噪聲來模擬噪聲等因素的干擾．

對(duì)上述水位模擬信號(hào)分別采用EMD與MEEMD方法進(jìn)行分解，得到的IMF分量個(gè)數(shù)分別為8和5，模擬水位序列和分解結(jié)果如圖1～圖3所示．

表1?水位序列數(shù)值模擬參數(shù)

Tab.1 Numerically simulated parameters of water-level series

從圖2中可以看出EMD分解結(jié)果受噪聲影響很大，出現(xiàn)了明顯的模態(tài)混疊現(xiàn)象．高頻分量IMF1將加入的白噪聲作為一個(gè)模態(tài)與模擬水位序列中的高頻成分疊加分解出來，成為一個(gè)不具備物理意義的虛假分量；而IMF2和IMF3也出現(xiàn)了明顯的模態(tài)混疊現(xiàn)象，其分解結(jié)果也不具備任何物理意義；IMF4～I(xiàn)MF8沒有明顯的模態(tài)混疊現(xiàn)象，但是也沒有周期性規(guī)律．因此，從圖中觀察可知，EMD方法的分解效果易受到噪聲等隨機(jī)因素的干擾，產(chǎn)生模態(tài)混疊現(xiàn)象與虛假分量，使得分解結(jié)果不再具備物理意義，自適應(yīng)性較差．對(duì)于有噪聲的數(shù)據(jù)序列，在進(jìn)行EMD計(jì)算之前還需先進(jìn)行降噪等預(yù)處理．從圖3中可以看出，MEEMD方法分解得到的各分量沒有出現(xiàn)模態(tài)混疊現(xiàn)象，特別是IMF1和IMF2分量的周期性規(guī)律明顯，IMF4和IMF5也有周期性特征．從IMF分量個(gè)數(shù)上來看，MEEMD方法得到了5個(gè)IMF分量，相對(duì)于EMD方法減少了3個(gè)IMF分量，有效減少了虛假分量．其中只有IMF3不對(duì)應(yīng)原始模擬水位序列的任何一個(gè)分潮成分，且沒有明顯的周期特征，是一個(gè)虛假分量．通過分析可知，MEEMD方法相對(duì)于EMD方法有效抑制了模態(tài)混疊現(xiàn)象，減少了虛假分量．除此之外，該方法可以實(shí)現(xiàn)在分解過程中的自適應(yīng)濾波，有效解決了噪聲等因素的干擾．原始的水位模擬序列有4個(gè)分潮，MEEMD方法和EMD方法得到的IMF分量均多于4個(gè)，因此要對(duì)分解結(jié)果做希爾伯特譜分析．

希爾伯特譜具有很好的時(shí)頻聚集性，可以在圖4中直觀地看出IMF分量的時(shí)頻特性，頻率的單位為每小時(shí)循環(huán)數(shù)，表示為c/h．從圖4(a)MEEMD得到的Hilbert-Huang譜中可以看到2條清晰的頻率線，這說明分解得到的其中2個(gè)高頻分量頻率集中，分解效果較好．從圖4(b)的邊際譜中可以看到3處頻率集中的位置，其中有1處位于低頻段，這說明分解有效，獲得了低頻分量，但所占能量比例較?。畧D4(c)中EMD得到的Hilbert-Huang譜中頻率線模糊，這說明分量的頻率分散，從圖4(d)的邊際譜中觀察可知分解結(jié)果頻率分散，有較多虛假分量存在，分解效果較差．

上述數(shù)值模擬的結(jié)果初步表明，MEEMD方法能夠?qū)性肼暤乃粩?shù)據(jù)進(jìn)行自適應(yīng)濾波并有效抑制模態(tài)混疊現(xiàn)象，分解效果較好．但是該方法還具有局限性，不能將頻率十分相近的分潮分離出來，而是疊加在了一起，如IMF1是半日分潮M2和S2的疊加，分解結(jié)果對(duì)于長(zhǎng)周期分潮也沒有很好的體現(xiàn)．因此，基于MEEMD的方法可以用于提取水位序列中的半日潮族與全日潮族；對(duì)于淺海海域，頻率較高的淺海分潮也包含在高頻IMF1分量中被提取出來，這樣得到的IMF1和IMF2之和除了全日分潮和半日分潮成分，還有淺海分潮成分．

2.2?數(shù)值模擬水位序列的能量譜

根據(jù)式(12)計(jì)算希爾伯特能量，可以得到水位序列的能量在各IMF分量的分布情況，如圖5所示．從圖中可以看出，高頻的IMF1和IMF2分量占總能量的絕大部分，經(jīng)計(jì)算可知其能量值占總能量的96.39%，是主要的波動(dòng)部分．

（d）EMD分解結(jié)果的邊際譜

圖5?IMF分量能量分布

3?實(shí)驗(yàn)分析

3.1?資料來源

為驗(yàn)證HHT方法的適用性，本文對(duì)實(shí)測(cè)的水位記錄進(jìn)行了實(shí)驗(yàn)與分析．資料來源于美國國家海洋和大氣管理局(NOAA)，本研究選取了位于美國南灘編號(hào)為9435380驗(yàn)潮站水位記錄中，2016—2017年的逐時(shí)水位記錄作為實(shí)測(cè)分析數(shù)據(jù)．NOAA提供的資料顯示，該驗(yàn)潮站建于1967年1月28日，位于北緯44°37.5′N，西經(jīng)124°2.7′W，時(shí)間子午線為120°W，起算基準(zhǔn)面為平均低低潮面(MLLW)．將資料數(shù)據(jù)進(jìn)行MEEMD分解，對(duì)分解得到的各IMF分量進(jìn)行希爾伯特變換后得到希爾伯特譜、邊際譜與能量譜，并計(jì)算各IMF分量的平均振幅與平均周期．將基于HHT方法分析得到的結(jié)果與NOAA提供的該驗(yàn)潮站潮汐調(diào)和分析數(shù)據(jù)以及經(jīng)過準(zhǔn)調(diào)和分析獲得的結(jié)果進(jìn)行對(duì)比，驗(yàn)證了方法的適用性．

3.2?實(shí)測(cè)水位記錄的分析結(jié)果

圖6是2016—2017年的長(zhǎng)期水位序列經(jīng)MEEMD分解后得到的9個(gè)IMF分量與余項(xiàng)．從圖中可以看出高頻IMF分量的分解效果較好，同一IMF分量中沒有出現(xiàn)模態(tài)混疊現(xiàn)象．其中IMF1和IMF2兩個(gè)高頻分量有明顯的周期性波動(dòng)．因?yàn)槌毕牟▌?dòng)是由多個(gè)分潮疊加而成的復(fù)雜振動(dòng)，并且原始數(shù)據(jù)序列中還帶有噪聲等因素的干擾，因此分解的結(jié)果中IMF分量較多，需要進(jìn)一步的分析．

為了更加清晰地了解各IMF分量的時(shí)頻特性，對(duì)IMF1～I(xiàn)MF9進(jìn)行Hilbert譜分析，得到Hilbert-Huang譜以及邊際譜，如圖7所示．圖7(a)是Hilbert-Huang譜，反映了水位序列的時(shí)間-頻率特性．從圖7(a)中可以發(fā)現(xiàn)兩條清晰的頻率線，分別在頻率為0.08c/h與0.04c/h附近，這表示其中兩個(gè)主要高頻IMF分量的頻率聚集性良好，具有明顯的周期性；在低頻0～0.02c/h之間也有頻率聚集現(xiàn)象，這表明在低頻范圍內(nèi)IMF分量的頻率也較為集中，但與兩處高頻聚集區(qū)相比并不明顯．圖7(b)是邊際譜，反映了水位序列的幅值-頻率特性．與Hilbert-Huang譜相比，在低頻0～0.02c/h之間的頻率聚集現(xiàn)象更為明顯．這部分頻率可能有長(zhǎng)周期分潮與之對(duì)應(yīng)，但是相對(duì)于高頻分量所占比例較?。ㄟ^實(shí)驗(yàn)分析可知，實(shí)測(cè)數(shù)據(jù)的分析結(jié)果與數(shù)值模擬的結(jié)果相符合．

除了1年以上的長(zhǎng)期觀測(cè)數(shù)據(jù)，還需要對(duì)短期的水位數(shù)據(jù)進(jìn)行分析．圖8(a)和圖8(b)分別為選取了2017年中某3d和7d的水位數(shù)據(jù)分解結(jié)果．短期觀測(cè)的水位數(shù)據(jù)得到較少的IMF分量，3d和7d兩組水位序列均分解出3個(gè)IMF分量．

圖6?2016—2017年原始水位序列與IMF1～I(xiàn)MF9及余項(xiàng)

圖7?2016—2017年水位序列Hilbert-Huang譜與邊際譜

從圖8中可以看出，短期水位序列的分解結(jié)果沒有出現(xiàn)明顯的模態(tài)混疊現(xiàn)象，各IMF分量均具有周期波動(dòng)特性．將選取的3d水位序列進(jìn)行時(shí)頻分析，得到圖9(a)和圖9(b)的Hilbert-Huang譜和邊際譜圖．從圖中可以看出，短期水位序列也具有類似于長(zhǎng)期水位序列的時(shí)間-頻率特性，頻率主要在0.08c/h和0.04c/h附近聚集，且振動(dòng)能量在0.08c/h附近更為集中；在低頻和更高頻的區(qū)域幾乎沒有能量分布．

圖8?短期原始水位序列及分解結(jié)果

圖9?3 d水位序列Hilbert-Huang譜與邊際譜

圖10是歸一化的各IMF分量的能量以及各IMF分量與原始水位序列相關(guān)系數(shù)的對(duì)比分析．其中相關(guān)系數(shù)表征IMF分量與原始序列的相似程度，相關(guān)系數(shù)較大的IMF分量與原始序列相似程度高，帶有更多能量．從圖中的對(duì)比觀察可知，各IMF的能量分布與相關(guān)系數(shù)分布趨勢(shì)一致，高頻IMF分量占據(jù)原始水位序列的絕大部分振動(dòng)能量，經(jīng)計(jì)算能量占有率為94.98%，這與譜分析的結(jié)果一致．短期水位序列的分析結(jié)果也是如此，高頻IMF1和IMF2分量幾乎占據(jù)了振動(dòng)的全部能量．通過實(shí)驗(yàn)分析可知，不管是長(zhǎng)期還是短期實(shí)測(cè)水位數(shù)據(jù)，分析結(jié)果均與數(shù)值模擬的結(jié)果相符合．

圖10?歸一化能量與相關(guān)系數(shù)對(duì)比

3.3?HHT分析方法與潮汐調(diào)和分析結(jié)果對(duì)比

NOAA提供的資料顯示，該驗(yàn)潮站處由潮汐調(diào)和分析得到共37個(gè)分潮，其中主要半日分潮、全日分潮和淺水分潮的調(diào)和常數(shù)與周期如表2所示，振幅高度以當(dāng)?shù)亻L(zhǎng)期平均海面(MSL)為零點(diǎn)起算．

表2?驗(yàn)潮站分潮參數(shù)

Tab.2?Tidal parameters of tide gauge stations

為了與調(diào)和分析的結(jié)果進(jìn)行對(duì)比，計(jì)算出2016—2017年的長(zhǎng)期水位序列各IMF的平均周期和平均振幅．根據(jù)式(8)與式(9)計(jì)算得到各IMF的瞬時(shí)幅值與瞬時(shí)頻率序列后，再分別求這兩個(gè)序列的均值，可得平均振幅和平均周期如表3所示．

通過對(duì)比表2和表3的分析結(jié)果可知，IMF1和IMF2分量的平均周期與平均振幅分別與半日分潮和全日分潮接近．IMF1分量的平均周期與半日分潮周期一致，特別是與M2分潮的周期最為接近，平均振幅略大于M2分潮的振幅，這說明該分量包含的不是單一分潮，而是多個(gè)半日分潮的疊加，其中可能還包含頻率更高但是比例很小的淺海分潮成分，因?yàn)镸2分潮占比最大，所以更多體現(xiàn)出M2分潮的性質(zhì)．IMF2分量的平均周期與全日分潮周期一致，特別是與K1分潮的周期最為接近，但是平均振幅卻小于K1分潮的振幅，這說明該分量中部分成分被分解到了其他模態(tài)中．因?yàn)樵诜纸膺^程中各IMF分量之間的正交性不能完全保證，MEEMD方法雖然可以有效抑制模態(tài)混疊現(xiàn)象卻不能完整解決，在分解時(shí)會(huì)存在能量損失．

表3?各IMF分量的平均振幅與平均周期

Tab.3?Average amplitude and period of IMF components

3.4?HHT分析方法與潮汐準(zhǔn)調(diào)和分析結(jié)果對(duì)比

根據(jù)潮汐的準(zhǔn)調(diào)和分析方法，分別選取3組長(zhǎng)度為3d的水位記錄，進(jìn)行準(zhǔn)調(diào)和分析，得到K1、O1、M2、S2、M4、MS4準(zhǔn)調(diào)和分潮的振幅(單位：m)，如表4所示．

表4?準(zhǔn)調(diào)和分潮振幅

Tab.4?Quasi-harmonic tidal amplitude m

對(duì)3組實(shí)驗(yàn)數(shù)據(jù)進(jìn)行HHT分析，得到各IMF分量的平均周期和平均振幅如表5所示．

表5?各次實(shí)驗(yàn)IMF分量的平均振幅和平均周期

Tab.5 Average amplitude and period of IMF components in each experiment

通過表5中的數(shù)據(jù)可知，實(shí)驗(yàn)1和實(shí)驗(yàn)3分別獲得了3個(gè)IMF分量，而實(shí)驗(yàn)2獲得了2個(gè)IMF分量，但是3次實(shí)驗(yàn)獲得的IMF1和IMF2分量的平均振幅和平均周期十分相近，平均振幅可分別對(duì)應(yīng)于表4中準(zhǔn)調(diào)和分潮結(jié)果中M2和K1分潮的振幅．而通過表3和表5的對(duì)比可知，在平均周期上IMF1和IMF2分量的平均周期分別介于M2、S2和K1、O1之間，所以可以認(rèn)為這兩個(gè)分量分別為半日準(zhǔn)調(diào)和分潮和全日準(zhǔn)調(diào)和分潮的疊加；因?yàn)樵谡穹螹2和K1分潮占主導(dǎo)地位，且在分解過程中存在能量損失，因此IMF1和IMF2分量分別主要對(duì)應(yīng)于M2和K1準(zhǔn)調(diào)和分潮．IMF3分量的平均振幅和平均周期與表4中的準(zhǔn)調(diào)和分析結(jié)果沒有明顯的對(duì)應(yīng)關(guān)系，且實(shí)驗(yàn)1和實(shí)驗(yàn)3的結(jié)果差別很大，所以IMF3分量不具備明顯的分潮特性，可能為其他長(zhǎng)周期分潮的疊加．因?yàn)槊看芜x取的實(shí)驗(yàn)數(shù)據(jù)長(zhǎng)度只有3d，對(duì)應(yīng)于一個(gè)朔望月(29.5306d)的不同時(shí)期，所以獲得的平均振幅結(jié)果會(huì)有一定差異．例如當(dāng)朔或望時(shí)為大潮期，獲得的分量平均振幅較大；而當(dāng)上、下弦月時(shí)為小潮期，獲得的分量平均振幅較小，表3和表4的振幅分析結(jié)果與此相一致．3d水位數(shù)據(jù)不能很好地體現(xiàn)長(zhǎng)周期分潮的性質(zhì)，因此雖然IMF3具有明顯的波動(dòng)規(guī)律，但選取不同時(shí)間的水位數(shù)據(jù)得到的IMF3分量的平均周期結(jié)果可能不同．準(zhǔn)調(diào)和分潮獲得了淺水分潮成分，在HHT分析結(jié)果中并沒有體現(xiàn)，可能包含在IMF1分量中或是損失于分解過程中．

4?結(jié)?語

希爾伯特-黃變換(HHT)是一種新的處理非平穩(wěn)、非線性信號(hào)或數(shù)據(jù)的方法，本文采用基于改進(jìn)的總體平均經(jīng)驗(yàn)?zāi)B(tài)分解(MEEMD)的HHT方法對(duì)原始水位數(shù)據(jù)進(jìn)行分析．通過MEEMD方法，把具有潮汐波動(dòng)與其他噪聲干擾的復(fù)雜驗(yàn)潮站水位數(shù)據(jù)記錄，分解為有限數(shù)目且具有一定物理意義的本征模態(tài)函數(shù)(IMF)．通過希爾伯特變換計(jì)算各IMF的瞬時(shí)頻率、幅值等得到數(shù)據(jù)序列的時(shí)頻表示與能量分布情況．通過數(shù)值模擬與實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)分析，HHT方法可以對(duì)原始數(shù)據(jù)進(jìn)行自適應(yīng)濾波分解，得到的高頻IMF分量的時(shí)頻聚集性良好，但是包括長(zhǎng)周期分潮在內(nèi)的低頻分量的分解效果并不理想．將分解結(jié)果與傳統(tǒng)的潮汐調(diào)和分析得到的分潮進(jìn)行對(duì)比可知，IMF1和IMF2分量分別與主要半日分潮和全日分潮一致，其和可以認(rèn)為是主要半日分潮、全日分潮以及高頻淺海分潮的疊加．雖然MEEMD方法有效抑制了模態(tài)混疊現(xiàn)象，但是對(duì)于結(jié)構(gòu)復(fù)雜的信號(hào)或數(shù)據(jù)仍有模態(tài)混疊現(xiàn)象存在，需要對(duì)各分量做進(jìn)一步細(xì)致的分析．對(duì)于3～7d短期的實(shí)測(cè)數(shù)據(jù)，HHT方法的分解結(jié)果沒有出現(xiàn)明顯的模態(tài)混疊現(xiàn)象，高頻IMF1和IMF2分量的平均振幅分別與準(zhǔn)調(diào)和分析獲得的M2和K1相一致，且分析結(jié)果具有較高精度．雖然本研究沒有使用HHT方法對(duì)淺水分潮進(jìn)行進(jìn)一步的分析，但該方法相比于傳統(tǒng)調(diào)和分析方法來說，仍然具有一定優(yōu)勢(shì)，尤其當(dāng)收集到的水位數(shù)據(jù)沒有經(jīng)過預(yù)處理(如含有涌浪等高頻噪聲)或水位數(shù)據(jù)殘缺不全、數(shù)據(jù)來源信息不明而不足以進(jìn)行傳統(tǒng)的調(diào)和或準(zhǔn)調(diào)和分析時(shí)，可以使用HHT方法進(jìn)行有效的補(bǔ)充分析，以獲得主要半日分潮、全日分潮和余水位的幅值信息．該方法目前還不能夠用于精確地預(yù)報(bào)潮位，但是在采用該方法提取有效水位特征信息基礎(chǔ)上，可以借助神經(jīng)網(wǎng)絡(luò)等方法實(shí)現(xiàn)對(duì)潮位的預(yù)測(cè)和預(yù)報(bào)，這部分內(nèi)容需要做進(jìn)一步的研究．

［1］ 張鳳燁，魏澤勛，王新怡，等. 潮汐調(diào)和分析方法的探討[J]. 海洋科學(xué)，2011，35(6)：68-75.

Zhang Fengye，Wei Zexun，Wang Xinyi，et al. Tidal harmonic analysis[J]. Marine Sciences，2011，35(6)：68-75(in Chinese).

［2］ 宋艷朋. 潮汐調(diào)和分析預(yù)報(bào)與基準(zhǔn)面計(jì)算軟件實(shí)現(xiàn)及南海應(yīng)用研究[D]. 青島：山東科技大學(xué)，2017.

Song Yanpeng. The Implementation of Tide Harmonic Analysis Forecasting and Datum Calculation Software and Application in South China Sea[D]. Qingdao：Shandong University Of Science and Technology，2017(in Chinese).

［3］ 方國洪. 潮汐和潮流的分析和預(yù)報(bào)[M]. 北京：海洋出版社，1986.

Fang Guohong. Analysis and Forecast of Tides and Tidal [M]. Beijing：China Ocean Press，1986(in Chinese).

［4］ Huang N E，Shen Z，Long S R，et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[J]. Proceedings：Mathematical，Physical and Engineering Sciences，1998，454(1971)：903-995.

［5］ Huang N E，Shen Z，Long S R. A new view of nonlinear water waves：The Hilbert spectrum[J]. Ann Rrev Fluid Mech，1999，31(1)：417-457.

［6］ 黃誠惕. 希爾伯特-黃變換及其應(yīng)用研究[D]. 成都：西南交通大學(xué)，2006.

Huang Chengti. Study on Hilbert-Huang Transform and Its Application[D]. Chengdu：Southwest Jiaotong University，2006(in Chinese)

［7］ Chen X，F(xiàn)eng Y，Huang N E. Global sea level trend during 1993—2012[J]. Global and Planetary Change，2014，112：26-32.

［8］ Chambers D P. Evaluation of empirical mode decomposition for quantifying multi-decadal variations in sea level records[J]. Nonlinear Processes in Geophysics，2015，1(2)：157-166.

［9］ 楊金湘，袁方超，李郅明，等. EMD方法在中國沿岸風(fēng)暴潮增水分析中的應(yīng)用[J]. 海洋通報(bào)，2016，35(2)：157-162，169.

Yang Jinxiang，Yuan Fangchao，Li Zhiming，et al. Application of the EMD method in the storm surges at the coast of China[J]. Marine Science Bulletin，2016，35(2)：157-162，169(in Chinese).

［10］ 李嘉康，李其杰，趙?穎，等. 基于CEEMD-BP神經(jīng)網(wǎng)絡(luò)的海溫異常預(yù)測(cè)研究[J]. 數(shù)學(xué)的實(shí)踐與認(rèn)識(shí)，2017，47(24)：163-171.

Li Jiakang，Li Qijie，Zhao Ying，et al. Sea surface temperature anomaly forecast research based on CEEMD-BP neural network[J]. Mathematics in Practice and Theory，2017，47(24)：163-171(in Chinese).

［11］ Cheng Y，Ezer T，Atkinson L P，et al. Analysis of tidal amplitude changes using the EMD method[J]. Continental Shelf Research，2017，148：44-52.

［12］ Pan H，Guo Z，Wang Y，et al. Application of the EMD method to river tides[J]. Journal of Atmospheric & Oceanic Technology，2018，35：809-819.

［13］ Wu Z H，Huang N E. Ensemble empirical mode decomposition：A noise-assisted data analysis method[J]. Advances in Adaptive Data Analysis，2009，1(1)：1-41.

［14］ Yeh J R，Shinh J S，Huang N E. Complementary ensemble empirical mode decomposition：A novel noise enhanced data analysis method[J]. Advances in Adaptive Data Analysis，2010，2(2)：135-156.

［15］ 鄭近德，程軍圣，楊?宇. 改進(jìn)的EEMD算法及其應(yīng)用研究[J]. 振動(dòng)與沖擊，2013，32(21)：21-26，46.

Zheng Jinde，Cheng Junsheng，Yang Yu. Modified EEMD algorithm and its applications[J]. Journal of Vi-

bration and Shock，2013，32(21)：21-26，46(in Chinese).

［16］ 黃大吉，趙進(jìn)平，蘇紀(jì)蘭. 希爾伯特-黃變換的端點(diǎn)延拓[J]. 海洋學(xué)報(bào)，2003，25(1)：1-11.

Huang Daji，Zhao Jinping，Su Jilan. Practical implementation of the Hilbert-Huang transform algorithm[J]. Acta Oceanologica Sinica，2003，25(1)：1-11(in Chinese).

［17］ 胡愛軍，孫敬敬，向?玲. 經(jīng)驗(yàn)?zāi)B(tài)分解中的模態(tài)混疊問題[J]. 振動(dòng)、測(cè)試與診斷，2011，31(4)：429-434，532-533.

Hu Aijun，Sun Jingjing，Xiang Ling. Mode mixing in empirical mode decomposition[J]. Journal of Vibration，Measurement & Diagnosis，2011，31(4)：429-434，532-533(in Chinese).

［18］ 練繼建，榮欽彪，董霄峰，等. 抑制模態(tài)混疊的HHT結(jié)構(gòu)模態(tài)參數(shù)識(shí)別方法研究[J]. 振動(dòng)與沖擊，2018，37(18)：1-8.

Lian Jijian，Rong Qinbiao，Dong Xiaofeng，et al. Structural model parameter identification method based on an improved HHT for suppressing mode mixing[J]. Journal of Vibration and Shock，2018，37(18)：1-8(in Chinese).

［19］ Bandt C，Pompe B. Permutation entropy：A natural complexity measure for time series[J]. Physical Review Letters，The American Physiological Society，2002，88(17)：174102-1-174102-4.

［20］ 陳?哲，李亞安. 基于多尺度排列熵的艦船輻射噪聲復(fù)雜度特征提取研究[J]. 振動(dòng)與沖擊，2019，38(12)：225-230.

Chen Zhe，Li Ya’an. Research on complexity feature extraction of ship radiated noise based on multi-scale permutation entropy[J]. Journal of Vibration and Shock，2019，38(12)：225-230(in Chinese).

Study of the Analysis of Tides Based on the HHT Method

Zhang Liang，Zhang Jiali，Zhang Xuefeng，Zhang Anmin

(School of Marine Science and Technology，Tianjin University，Tianjin 300072，China)

Empirical mode decomposition(EMD)is the core component of the Hilbert-Huang transform(HHT).Based on research and analysis of the traditional HHT method，we developed a method known as modified ensemble empirical mode decomposition(MEEMD)to analyze ocean tides．With respect to the characteristics of a water-level data series，we used MEEMD method to eliminate the false intrinsic mode function(IMF)in the decomposition process by adding paired white noise and calculating the permutation entropy，thus suppressing the occurrence of mode aliasing．A comparison of the numerical simulation results with measured water-level data shows that MEEMD method can effectively suppress modal aliasing and realize adaptive filtering for original water-level sequences with noise．In this study，the high-frequency IMF1and IMF2components decomposed by MEEMD method accounted for more than 90% of the total energy of the original water-level data．The average period and amplitude of the IMF1and IMF2components were basically consistent with the results of the semidiurnal-constituent and diurnal-constituent tidal components obtained by harmonic analysis，respectively．The results show that the sum of IMF1and IMF2is the superposition of the main semidiurnal，diurnal and shallow water tidal components，respectively. Traditional tidal harmonic analysis methods rely on long-term water-level records and require some prior knowledge of the water-level data．The astronomical knowledge and mathematical calculations involved are also more complex．Quasi-harmonic analysis can process short-term water-level data，but still requires astronomical knowledge and complex mathematical calculations，as well as high data integrity．When the collected water-level data are inadequate for traditional harmonic or quasi-harmonic analysis due to noise or lack of information，the HHT method can be used for effective supplementary analysis to obtain the magnitude information of the main semidiurnal constituent，diurnal constituent，and residual water levels．Therefore，MEEMD-HHT method can be used as a good supplementary method for tidal analysis．

modified ensemble empirical mode decomposition(MEEMD)；Hilbert-Huang transform(HHT)；filtering；tidal constituent

P731.23

A

0493-2137(2020)07-0725-11

10.11784/tdxbz201908071

2019-08-30；

2019-10-07.

張?亮（1983—??），男，博士，講師，liangzhang@tju.edu.cn．

張學(xué)峰，xfz_nmdis@126.com.

國家自然科學(xué)基金資助項(xiàng)目(41876014，11801402).

Supported by the National Natural Science Foundation of China(No. 41876014，No. 11801402).

(責(zé)任編輯：金順愛)