大連老虎灘海域潮汐調和分析

許秀娥，張容榕，韋冬妮

（1.國家海洋局東港海洋環境監測站，遼寧 丹東 118000；2.國家海洋局大連海洋環境監測中心站，遼寧 大連 116015）

1 引言

潮汐現象是指海水在天體（主要是月球和太陽）引潮力作用下產生的周期性運動，習慣上把海水垂直方向的漲落稱為潮汐，水平方向的流動稱為潮流。在人類活動比較頻繁的近海，海水養殖、海洋運輸和海洋工程等涉海產業均與潮汐息息相關，因此開展潮汐運動規律的研究對物理海洋學的發展和國民生產有重要意義。此外，潮汐能也是一種可持續利用的清潔再生能源，潮汐能的開發利用有助于實現碳中和目標。中國近海的潮汐能開發利用潛力以東海最佳，黃海次之。在第四次中國近?？稍偕茉凑{查與研究中，遼寧省近海有24個站址入選裝機能量達500 kW 以上的潮汐能站址備選資源，僅次于福建、浙江及上海[1]。

大連灣位于黃海北部遼東半島南端，是我國北方重要的港口。老虎灘驗潮站位于大連灣灣口西南，針對渤海開展的潮汐潮流數值模擬大部分以大連外海為開邊界，或者利用該站的數據進行驗證[2-4]。因此，通過老虎灘驗潮站的實測數據分析掌握潮汐運動規律，提高潮汐預報的精度，對于近海海洋活動的安全保障和防災減災等具有重要意義。

2 數據及方法

2.1 觀測數據

研究數據來源于國家海洋局大連海洋環境監測站老虎灘驗潮站的逐時觀測潮位。潮位觀測資料的時間長短對于其調和分析的結果影響較大，通常只有當觀測時段的長度不短于任意兩個分潮的會合周期時，這些分潮才能彼此分離。由于不同亞群分潮之間的會合周期最長為一個回歸年，所以當觀測時段的長度顯著短于1 a 時，我們就認為記錄是不完整的，因此，潮汐分析時通常選用的時間長度大于1 a[5]。

為此，選取時間長度大于1 a的潮位觀測數據開展潮汐調和分析，時間范圍為2017 年1 月1 日—2018 年1 月4 日共369 d，時間序列中共包含8 857個逐時潮高[6]。同時，為了對比短期（3 M）潮位資料與長期（1 a）潮位資料調和分析結果的差異，又選用了2017 年3 月1 日—5 月31 日3 M 的潮位觀測資料進行分析。對潮位觀測數據進行質量控制，異常值或缺測值進行線性插值，進而獲取逐時無缺測的完整數據樣本，并將時間序列調整為世界時。

本研究采用Foreman程序的自動選取分潮功能進行分潮選取，對于不足1 a 的資料，采用差比法推算隨從分潮，將信噪比大于1的分潮作為顯著分潮。在獲得兩個不同時間長度資料的調和常數后，分別對2019 年的天文潮進行預報，并與2019 年的潮位觀測資料和國家海洋信息中心發布的潮汐表數據進行對比分析和殘差檢驗。

2.2 分析方法

最初的潮汐分析方法是將潮汐漲落與天體運動直接建立相關關系，這類方法稱為非調和法[7]。Newton[8]應用萬有引力定律也對地球的潮汐現象進行了解釋。根據物理學有關原理可知，任何一種周期性的運動都可以由許多簡諧振動組成。潮汐變化是一種非常近似的周期性運動，因而也可以分解為許多固定頻率的分潮波，進而求得分潮的調和常數（振幅和遲角），這種分析潮汐的方法稱為潮汐調和分析。Darwin 等.[9]此后發展了調和方法，通過對引潮力進行調和展開，得到了主要的潮汐分潮頻率。Doodson[10]將引潮勢進一步展開為純調和分潮，引用了月球運動的Brown 系數和Newcomb 表，使引潮勢的展開結果更為精準[7]。近年來，隨著計算機運算性能的迅速提高，潮汐調和分析的能力和預報的速度得到大幅提升，調和分析法成為潮汐分析和預報的主要方法，非調和法退居次要地位。此后，從計算量中解放出來的科學家更加致力于提升分析和預報的準確度。Amin[11]放棄了傳統的交點因子的假設，對19 a 的觀測資料進行了分析。方國洪等[12]采用準調和淺水分潮來表示潮汐的高頻部分。王曉東等[13]對廈門和東山的潮位資料進行譜分析，求出調和常數，利用分析結果回報兩港的潮位，驗證了潮汐譜分析對潮汐預報具有一定的應用價值。

Foreman[14]在1977 年開發了一套潮汐分析與預報的Fortran 程序，2002 年Pawlowicz 等[15]將該程序改寫成Matlab 程序，命名為T_TIDE 工具箱（以下簡稱T_TIDE）。T_TIDE 用復數算式代替了傳統的實數算式，并增加了隨從分潮的推算和非線性誤差分析等一些用戶自定義界面。2013 年，Matte等[16]在研究潮汐受河流入流影響產生河潮（river tides）而出現的非平穩信號時，在T_TIDE 基礎上改進程序從而誕生了NS_TIDE。2018 年，Pan 等[17]利用潮汐調和分析的新方法EHA，開發了工具包S_TIDE。該工具包既可以分析平穩潮，還能分析非平穩潮，進一步完善了潮汐調和分析方法。

老虎灘驗潮站位于遼寧省大連市老虎灘漁人碼頭海港內，北緯38°52′，東經121°41′，四周沒有河流匯入，以規則半日潮平穩信號為主，因此用T_TIDE來分析該驗潮站數據。

T_TIDE 基于調和分析法可以解析分離一組潮位觀測數據中的潮汐和非潮汐信號，最多能分離出45 個天文分潮和101 個淺水分潮，保留了Foreman程序中的自動篩選分潮的方法，用戶可根據需要自行添加淺水分潮。T_TIDE 主要針對1 a 及少于1 a的觀測資料進行分析，對于少于1 a 的觀測資料，用戶還可以根據同一群分潮中主分潮與隨從分潮之間的振幅比和遲角差，選擇利用主分潮來推算隨從分潮，例如利用K1推算P1，S2推算K2。對于1 a 以上的時間序列，T_TIDE 進行的交點因子訂正計算會不夠準確。

T_TIDE 在潮位觀測資料分析[18]和模式結果潮汐調和分析[19-20]等方面有著較為廣泛的應用。本文利用T_TIDE 對大連老虎灘驗潮站一年期和一個季度期兩種時間長度的潮位觀測數據進行了調和分析，并利用得出的調和常數進行天文潮預報，最后將天文潮預報結果與潮位資料及潮汐表數據進行對比。

3 結果與討論

3.1 水位功率譜分析結果

譜分析是時間序列在頻域上進行分析的方法，亦稱頻譜分析或波譜分析?？紤]到海洋運動無論是時間上還是空間上均存在各種尺度的波動現象，所以頻譜分析是分析海洋各種波動現象的常用方法，被廣泛應用在海浪譜分析、潮汐和潮流統計預報、海平面長期變化和氣候變化研究等各個方面。

利用Matlab 提供的pwelch 程序對老虎灘站2017 年的水位做功率譜分析。pwelch 為常用的功率譜分析函數，相關的參數設置如下：窗函數類型選擇hanning 窗，傅里葉變化點數為256 個，窗口重疊數為傅里葉變化點數的一般為128 個，采樣頻率為60 min，功率譜分析結果如圖1所示。

圖1 2017年觀測數據譜分析結果Fig.1 Spectral analysis result of 2017 observation data

功率譜表現出幾個特征：周期介于2 ～48 h 內出現了幾個明顯的譜峰，其中周期在24 h 或12 h 左右的譜密度較高，12 h 左右的譜密度最高，略高于24 h 的譜密度，這說明該海域半日潮占優。周期在8 h、6 h 和4 h 左右也出現了譜峰，但峰值明顯低于12 h和24 h。

3.2 1 a潮位資料分析主要分潮振幅

本文首先對老虎灘驗潮站369 d 的逐時潮位進行調和分析，篩選分析出67 個分潮，其中包含42 個天文分潮和25個淺水分潮。

圖2 為調和分析得到的各分潮的頻譜圖，振幅大于95%顯著線的為主要分潮，振幅小于顯著線的則為次要分潮。結果表明，67個分潮中有52個分潮為主要分潮（藍色），其中包括31 個天文分潮和21個淺水分潮，其他15 個則為不顯著的次要分潮（紅色）。

圖2 中頻率大于0.1 cph（cycle per hour）的高頻分潮多為淺水分潮，多數低頻分潮集中在0.08 cph（半日潮）和0.04 cph（日潮），這與水位功率譜分析的結果基本一致。1 a 周期分潮SA 和0.5 a 周期分潮SSA也是顯著分潮，但周期為1 M和0.5 M的低頻分潮并不顯著。老虎灘潮汐振幅最大的前6個主要分潮（按振幅大?。┓謩e為M2、S2、K1、SA、N2和O1。

圖2 老虎灘1 a潮位資料潮汐頻譜圖Fig.2 Frequency spectrum of one-year tidal level data in Laohutan station

潮汐類型以全日分潮和半日分潮的振幅比為量化指標，我國采用以下公式計算[7]：

式中，H表示分潮的振幅。我國判斷潮汐類型的標準為：當0 ＜F≤0.5時，為規則半日潮；當0.5 ＜F≤2.0時，為不規則半日潮；當2.0 ＜F≤4.0 時，為不規則全日潮；當F＞4時，為規則全日潮。經計算，老虎灘潮汐類型值為0.429 1，為規則半日潮海區。

3.3 3 M潮位資料分析主要分潮振幅

對于3 M的短期潮位資料，由于時間長度較短，部分分潮無法分辯，調和分析得出的分潮數會有很大不同。圖3是對老虎灘驗潮站短期觀測潮位數據的頻譜圖，未引入差比關系計算隨從分潮，也未添加任何淺水分潮。對于短期觀測資料，共分析出35個分潮，其中包括19 個天文分潮和16 個淺水分潮，比1 a期潮位資料獲得的67個分潮數少32個。

由圖3 可以看出，短期潮位資料調和分析出的頻譜圖與長期潮位資料的頻譜圖分布較類似，半日潮和全日潮振幅較大，M2、S2、O1和K14 個分潮最大，振幅基本一致，這也說明利用長期或短期資料均能調和出其主要分潮；主要的差異在于分離出的分潮數，短期資料分離出的高于綠色顯著線的分潮有24 個，包括11 個天文分潮和13 個淺水分潮。長周期分潮（比如1 a周期的SA和0.5 a周期的SSA）均未見，得出的周期為1 M和0.5 M的低頻分潮不顯著。

圖3 老虎灘3 M潮位資料潮汐頻譜圖（未引入差比關系和淺水分潮）Fig.3 Frequency spectrum of three-months tidal level data in Laohutan station（difference ratio and shallow-water tidal component are not introduced）

對小于1 a的潮位資料進行分析時，需要通過引入差比關系分析隨從分潮，還可根據需要添加淺水分潮。這里根據給出的主分潮和隨從分潮的差比關系推算出3 M 分潮中所缺的半日和全日隨從分潮，如隨從分潮σ1、P1、π1、K2、2N2和λ2可分別由主分潮O1、K1、K1、S2、N2和M2推出，添加的淺水分潮有SO3、MK4、SK4、2MK6、MSK6和MSN2等。

圖4 是引入差比關系和淺水分潮的分析結果。與未引入差比關系的結果相比，引入差比關系后分析出的分潮增加了14 個天文分潮和8 個淺水分潮，其中信噪比大于1的顯著分潮有9個（見圖4標注的分潮），其它的13 個都是非顯著的次要分潮。在9個新增的顯著分潮中有3個淺水分潮，其中MSN2分潮的頻率是0.085 cph，是半日潮性質，SO3和2MK6是高頻淺水分潮。9 個分潮中，振幅最大最顯著的分潮是半日潮K2，其次為全日潮P1，淺水分潮SO3的振幅雖然不是很大，但是其信噪比僅次于K2和P1。

圖4 老虎灘3 M潮位資料潮汐頻譜圖（引入差比關系和淺水分潮）Fig.4 Frequency spectrum of three-months tidal level data in Laohutan station（difference ratio and shallow-water tidal component are introduced）

綜上所述，1 a和3 M 潮位觀測資料的調和分析結果有較大區別，3 M 潮位資料得出的分潮數顯然比1 a資料得出的分潮數要少很多，雖然可以通過差比關系推算出其他隨從分潮的調和常數，但是尚有部分主分潮無法從3 M 潮位數據中分析出來。因此，為了保證預報的準確性，本文利用1 a 潮位數據分析出的調和常數進行2019年的天文潮預報，并與老虎灘驗潮站的實測潮位數據和天文潮數據（潮汐表數據）進行對比分析。

3.4 基于T_TIDE 調和分析的潮汐預報結果驗證評估

圖5a 是T_TIDE 得出的2019 年老虎灘天文潮預報值與驗潮站實測值的對比。圖5b 是2019 年國家海洋信息中心發布的老虎灘潮汐預報值與驗潮站實測值的對比圖。圖5c 是潮汐表預報值與T_TIDE 預報值的對比。為了更清晰地進行對比，我們選取了800 h的結果對比繪圖，并統一其起算面。

從對比結果來看，觀測潮位（實測值）中除了含有天文潮外，還包括由溫帶氣旋和熱帶氣旋等引起的余水位。由圖5a 可以看出，觀測水位與T_TIDE預報的天文潮在某些時間差異較大，在第700 h 左右能明顯看出較大的余水位。同樣的現象也出現在潮汐結果上（見圖5b）。利用T_TIDE 預報的天文潮與潮汐表的預測結果接近，差異較小，最大差值僅17.6 cm（見圖5c）。

圖5 T_TIDE預報值、實測值及潮汐表預報值兩兩對比圖Fig.5 Comparison diagrams of T_TIDE predicted value，measured value and predicted value by tidal table

為了更清晰地看出本文的分析結果與實際水位的差異，我們選取了2019 年全年的數據進行對比。圖6 是2019 年老虎灘T_TIDE 預報潮高與實測潮高的散點圖，其中橫軸為驗潮站的實測值，縱軸是T_TIDE 的預報結果，紅線代表預報值與實測值相等。圖7 是預報值與實測值的分位數圖，橫軸為實測分位數，縱軸為預報分位數，圖中的紅線是兩者完全重合的理論線。由兩圖可見，在整個潮高變動區間內，尤其是潮高在100～400 cm 以內時，兩者總體符合良好，而在實測值小于100 cm 和大于400 cm時，T_TIDE的預報結果偏高。

圖6 老虎灘2019年T_TIDE預報潮高與實測潮高散點圖Fig.6 Scatter diagram of T_TIDE prediction and observation of year 2019 in Laohutan station

圖7 老虎灘2019年T_TIDE預報潮高與實測潮高分位Fig.7 Quantile of T_TIDE prediction and observation of year 2019 in Laohutan station

平均絕對誤差（Mean Absolute Error，MAE）和均方根誤差（Root Mean Squared Error，RMSE）是衡量變量精度的兩個常用指標，可以用來衡量預報值與觀測值之間的偏差。兩者定義分別如下：

式中，H表示實測值，h表示預報值，這兩個值均會隨著樣本個數n的增大而增大。當預報值和實測值偏離不大時，這兩個值的比值接近于1。如果預報值中有偏離實測值很大的值存在，那么RMSE 會將偏離的誤差二次平方放大，最終導致RMSE 的值比MAE偏大。

表1是根據式（2）和（3）分別計算的潮汐表預報值及T_TIDE預報值的MAE和RMSE。從表中能看出MAE和RMSE相差不大，兩者的MAE∶RMSE都接近于1，說明潮汐表和T_TIDE 預報值中沒有異常值。T_TIDE 的MAE 和RMSE 比潮汐表的偏小，表明T_TIDE 預報值對實測值的累計誤差少，T_TIDE 的預報結果略優于潮汐表的預報結果。

表1 潮汐表與T_TIDE對實測值的精度指標Tab.1 Precision indexes of Tide Table and T_TIDE predictions to observation

為了進一步對比T_TIDE 和潮汐表的預報效果，將兩者預報的潮差最大值、最小值和實測值進行對比（見表2）。從潮差結果也可以發現T_TIDE的預報效果略好。

表2 潮汐表、T_TIDE及實測值潮差對比（單位：cm）Tab.2 Tidal range comparison between tide table，T_TIDE and observation（unit:cm）

3.5 T_TIDE預報余水位檢驗

實際觀測到的潮位可寫為：

式中，H為實測潮位；為多個調和分潮潮位疊加值；r為余水位或噪聲，它包括由氣象等因素引起的不規則擾動、觀測中存在的誤差、數據處理中的誤差、截斷誤差和被忽略的分潮等，有時也簡單地稱為觀測誤差。如果余水位r的概率呈現正態分布，則表明所得的結果更為可靠[6]。

圖8是T_TIDE預報的余水位的頻率分布，橫軸代表余水位（實測-預報）的預測殘差，縱軸表示各組余水位數據發生的頻率。由圖8 可以看出，T_TIDE 潮高預報的余水位服從正態分布，即高斯分布，這與Powlowicz等[15]的結論一致。

圖8 老虎灘余水位（實測-預報）頻率分布Fig.8 Frequency distribution of residual water level（measured-forecast）in Laohutan station

表3 給出了余水位的均值、方差和95%的置信區間。余水位均值近似為0，可以認為是服從均值為0 的正態分布，方差和區間長度小說明預報的準確度高。

表3 老虎灘余水位統計值Tab.3 Statistical values of residual water level in Laohutan station

4 結論

本文利用T_TIDE 工具箱對老虎灘驗潮站1 a潮位資料和3 M 潮位資料進行調和分析。結論如下：

（1）相比3 M 潮位資料分析結果，1 a 潮位資料能分析出較多的天文分潮和淺水分潮。老虎灘潮汐振幅最大的前6 個主要分潮分別為M2、S2、K1、SA、N2和O1，潮汐類型屬于規則半日潮。

（2）對3 M的潮位資料進行調和分析的結果與利用1 a資料調和分析得出的結果有很大不同，少于1 a 的潮汐資料不能分析出長周期的分潮以及一些主分潮，需要添加差比關系來計算隨從分潮，并且添加必要的淺水分潮。

（3）進一步利用T_TIDE 將1 a 潮位資料計算出的調和分潮常數進行老虎灘2019年天文潮預報，預報殘差符合正態分布，殘差均值小于10-2m，置信區間長度小于10-2，且該預報結果與潮汐表潮汐預報值相差不大。

以上結論表明，調和分析工具箱T_TIDE 對資料長度大于1 a的潮汐資料有較好的分析結果，并能較好地進行潮汐預報。對于資料長度太短的潮汐資料，其不能分析出長周期分潮，以致影響后續潮汐預報。