鄭 凱 馬 蘭 季 偉 郭博峰

1 武漢大學(xué)測繪學(xué)院，武漢市珞喻路129號，430079

海潮負(fù)荷效應(yīng)會引起測站的位移［1－3］。目前，對海潮負(fù)荷引起的地表位移改正的主要方法是利用全球海潮模型計算，但由于利用衛(wèi)星測高等資料建立的海潮模型在一些地區(qū)（特別是高緯度地區(qū)）模型改正精度不高，一些學(xué)者提出利用GPS技術(shù)測定海潮負(fù)荷位移參數(shù)的研究思路。Zhang等［4］用PPP技術(shù)獲得了南極Amery冰架地區(qū)的海潮波形，采用5d連續(xù)觀測的GPS數(shù)據(jù)反演得到8個分潮在垂直分量上的振幅和相位；Thomas等［5］比較了GPS 和VLBI技術(shù)測定海潮負(fù)荷位移參數(shù)的差異，表明利用GPS和VLBI技術(shù)測定海潮負(fù)荷位移參數(shù)是可行的，且PPP固定解計算得到的參數(shù)并沒有較浮點解得到的參數(shù)在精度上有明顯改善；Yuan等［6］利用香港12個測站3～7 a的GPS觀測資料，采用與Thomas等類似的海潮負(fù)荷位移參數(shù)估算策略，將分潮參數(shù)與PPP模型其余參數(shù)一起估計。該方法雖然理論上比較嚴(yán)密，但增加了PPP 數(shù)據(jù)處理的復(fù)雜度，過度的參數(shù)化還可能造成分潮正余弦系數(shù)和模糊度參數(shù)及鐘差參數(shù)之間的耦合相關(guān)性，也在一定程度上影響PPP定位的收斂速度。本文提出直接利用未顧及海潮負(fù)荷位移改正的動態(tài)PPP 解算得到的長時間（1a）測站坐標(biāo)序列，采用傅里葉變換方法提取海潮周期信息，進(jìn)而反演8個分潮的海潮負(fù)荷位移參數(shù)的方法，將反演結(jié)果與最新的5個全球海潮模型比較，分析評價該方法反演海潮負(fù)荷位移參數(shù)的可行性。

1 海潮參數(shù)反演方法

測站的海潮負(fù)荷位移可由11個主分潮的海潮負(fù)荷位移疊加獲得，其中3 個為長周期分潮。由于反演長周期分潮需要長時間的觀測數(shù)據(jù)，且長周期潮對測站位移的影響較其他分潮要小［7－8］，故在此僅研究8 個主分潮，即4 個半日潮（M2、S2、N2、K2）和4個周日潮（K1、O1、P1、Q1）。

測站的海潮負(fù)荷位移可表示為：

式中，Hk，i和φk，i分別為分潮i在k方向上的振幅和格林尼治相位；t為格林尼治時間；Vi為天文幅角初相（本文天文變量參考J2000.0系統(tǒng)）；ωi為分潮角速度，可查表獲得；fi和μi分別為節(jié)點因數(shù)和天文相角，可由天文變量計算得到；ξk為海潮負(fù)荷引起的測站位移偏差［9］。為便于求解Hk，i和φk，i，將式（1）改為：

式中，A＝fiHk，icos（φk，i－μi），B＝fiHk，isin（φk，i－μi）。由式（2）得：

式中，A、B為待估參數(shù)，ξk為PPP計算得到的E、N、U3個方向的測站坐標(biāo)序列。將式（2）以法方程形式表達(dá)，用經(jīng)典最小二乘準(zhǔn)則求解待估參數(shù)?？紤]到φ變化范圍從－PI到PI，實際解算結(jié)果需要加減180°。

2 PPP處理策略

海潮負(fù)荷對單點定位的影響在沿海地區(qū)可達(dá)mm級甚至cm級［1，4］，因此在對其他誤差精確改正之后，得到的測站位移即由海潮負(fù)荷引起。

本文采用無電離層組合的PPP模型［10－11］：其中，P為無電離層的偽距觀測值，Φ為無電離層載波相位觀測值，ρ為真實幾何距離，dt為衛(wèi)星鐘差，dT為接收機鐘差，Tr為對流層延遲，N為無電離層組合模糊度，λ為組合后的載波波長，ξ為系統(tǒng)誤差（相對論效應(yīng)、固體潮改正、地球自轉(zhuǎn)改正、天線相位中心偏差及變化），ξoc為海潮負(fù)荷引起的位移偏差，εp和εΦ分別為偽距觀測噪聲和相位觀測噪聲。解算過程中，將測站坐標(biāo)、接收機鐘差、無電離層組合模糊度及對流層天頂延遲視為未知數(shù)，采用Kalman濾波作為參數(shù)估計器，衛(wèi)星高度角設(shè)為7°，用TurboEidt方法探測周跳（不作修復(fù)）。同時，對接收機鐘跳進(jìn)行探測修復(fù)［12－13］，采用絕對天線相位中心模型（IGS＿08），用選權(quán)迭代的驗后殘差分析方法進(jìn)行質(zhì)量控制，并將模糊度當(dāng)作實數(shù)處理。為了得到包含海潮信息的測站位移序列，不對海潮負(fù)荷引起的位移進(jìn)行改正。衛(wèi)星鐘差、軌道以及硬件延遲選用IGS的精密產(chǎn)品。

由于相鄰天之間解算數(shù)據(jù)的不連續(xù)以及精密產(chǎn)品不一致等原因，可能導(dǎo)致天與天的結(jié)果文件出現(xiàn)跳變。為解決此問題，將全年觀測文件以7d為1組進(jìn)行合并，如，1～7d，6～12d，…，共73個組合觀測文件。對每個組合文件而言，因上述原因?qū)е聰?shù)據(jù)跳變的點只可能出現(xiàn)在文件的首尾時段。因此，對組合觀測文件進(jìn)行PPP動態(tài)解算之后，截取每組2～6d的解算坐標(biāo)，然后融合成全年的測站坐標(biāo)系列。

3 實例分析

3.1 數(shù)據(jù)

選取高緯度地區(qū)的3 個CORS 站AC67（57.798°N，152.336°W）、AC43（59.521°N，149.709°W）、AC39（58.610°N，152.394°W）以及中緯度地區(qū)的SHAO 測站（31.099°N，121，2°E）2012年全年的觀測數(shù)據(jù)進(jìn)行分析。由于缺乏測站坐標(biāo)真值，故將各測站解算結(jié)果與首歷元計算得到的測站坐標(biāo)求差，并轉(zhuǎn)換為E、N、U3個方向的坐標(biāo)序列。實驗表明，此方法會產(chǎn)生因地殼板塊運動導(dǎo)致的坐標(biāo)漂移。由于海潮負(fù)荷引起的位移僅是余弦波的疊加，故可通過一階線性擬合的方式予以消除?？紤]到數(shù)據(jù)量過大對計算機的性能要求也越高，將數(shù)據(jù)的采樣率由30s稀疏為5min。

3.2 海潮負(fù)荷位移參數(shù)反演

為驗證動態(tài)PPP 解算結(jié)果是否包含分潮信息，利用傅里葉變換對結(jié)果進(jìn)行頻譜分析。以測站AC43為例，將測站單天每個方向的坐標(biāo)序列進(jìn)行傅里葉變換后所得的振幅在頻域內(nèi)進(jìn)行疊乘［14］，可以分離出分潮信號與噪聲。

圖1給出了測站AC43的分潮周期圖。可以很明顯地看出周日潮和半日潮的周期信息，其次是1／3周日潮、1／4周日潮?？梢钥吹?，對測站而言，主要受周日潮和半日潮的影響，而1／3潮、1／4潮的影響僅為兩主潮的1／6～1／7，就目前PPP精度而言可以不予考慮。圖中長周期的分潮信息不明顯，這是因為本文選用的數(shù)據(jù)時間跨度還不夠長，不足以將更長周期的分潮信息提取出來。

圖1 測站AC43海潮頻譜分析周期Fig.1 Ocean tidal periodogram spectrum analysis about station AC43

進(jìn)一步從量化的角度分析。節(jié)點因數(shù)和天文相角近似線性變化如圖2、3所示，相較于其他分潮，K2分潮的節(jié)點因數(shù)和天文相角的變動略為明顯，分別為8%和1%左右。盡管如此，兩者對海潮負(fù)荷位移參數(shù)的影響也在亞mm級［7］，故在單歷元參數(shù)處理中不予考慮，僅對最終結(jié)果取節(jié)點因數(shù)和天文相角的年均值進(jìn)行改正。將測站位移代入式（2）、（3）式進(jìn)行參數(shù)解算，解算結(jié)果與TPXO.7.2海潮模型計算結(jié)果進(jìn)行比較，如表1、2所示。

圖2 全年節(jié)點因數(shù)變化趨勢Fig.2 Yearly change trend about nodal corrections which account for the modulating effect of the lunar node

圖3 全年天文相角變化趨勢Fig.3 Yearly change trend about astronomical phase angle

3.3 結(jié)果評價

選用5個最新的全球海潮模型作為參考（表3），并采用均方根（RMS）［7］來評定PPP反演的參數(shù)和海潮模型計算值的差異。計算方法如下：

其中，A、φ分別代表振幅和相位，n（1，2，…，N）為測站號，j為分潮號，k為坐標(biāo)分量，i為虛數(shù)符號。

圖4給出了4個測站反演結(jié)果與TPXO.7.2模型計算結(jié)果的均方根。造成測站間差異的主要原因與測站所處的地理環(huán)境有關(guān)。例如，不同地區(qū)海潮負(fù)荷對地殼形變的影響程度不同，其中SHAO 測站在U方向上RMS 較小，尤其是K2分潮RMS約為其他3個站的1／3，是因為利用海潮模型計算的K2分潮在SHAO 測站的振幅僅為亞mm級，而其他幾個測站K2分潮的振幅達(dá)到了幾個mm，說明SHAO 測站受K2分潮的影響（尤其在垂直方向）較其他測站小很多。此外，在淺海地區(qū)由于復(fù)雜的地勢，海潮模型并不能十分準(zhǔn)確地反映分潮信息［15］。4 個測站在K1、K2分潮的均方根最為顯著，主要是受軌道誤差和多路徑效應(yīng)的影響［14］。

為進(jìn)一步說明分潮的反演情況以及各模型間的差異，圖5給出了PPP反演參數(shù)與5個海潮模型計算值之間的均方根誤差統(tǒng)計結(jié)果?？梢钥吹?，除U方向Q1分潮各模型計算值有較大差異外，所選的5個全球海潮模型差異不大，說明在研究區(qū)域5個模型具有一致的精度。圖中還發(fā)現(xiàn)，利用PPP反演出與太陽有關(guān)的分潮參數(shù)（S2、K2、K1、P1）與模型間的偏差較大，其中在水平方向上K2、K1略為明顯，而在垂直方向上P1分潮的偏差則要高于K1。N方向的均方根誤差整體優(yōu)于2.5mm，E方向次之優(yōu)于3.5mm，U方向相對偏差最大但仍優(yōu)于6 mm，主要是垂直方向的定位精度低于水平精度所致。而影響均方根誤差的原因除了和PPP 反演結(jié)果與模型間的系統(tǒng)偏差之外，數(shù)據(jù)長度也是一個重要原因。

表1 測站半日潮負(fù)荷位移參數(shù)反演結(jié)果與TPXO.7.2模型半日潮計算結(jié)果比較Tab.1 Compared PPP－derived results with TPXO.7.2about semi－diurnal constituents

表2 測站周日潮負(fù)荷位移參數(shù)反演結(jié)果與TPXO.7.2模型周日潮計算結(jié)果比較Tab.2 Compared PPP－derived results with TPXO.7.2about ter－diurnal constituents

表3 全球海潮模型Tab.3 Global ocean tide model

圖4 PPP反演結(jié)果與TPXO.7.2模型之間的比較Fig.4 RMS of the PPP－derived OTL displacement estimates to the TPX.O.7.2model estimates

圖5 PPP反演結(jié)果與模型的均方根誤差Fig.5 RMS of the PPP－derived OTL displacement estimates to all ocean tide model estimates

在不考慮系統(tǒng)差的前提下，利用1a的GPS觀測數(shù)據(jù)反演得到的參數(shù)與模型間的均方根誤差與Yuan等［6］計算結(jié)果基本相當(dāng)，這也說明利用本文方法反演海潮負(fù)荷位移參數(shù)可以較快地得到收斂。

4 結(jié) 語

1）本文方法計算的海潮負(fù)荷位移參數(shù)與模型間的均方根誤差為mm級，與Yuan等［6］方法解算結(jié)果基本相當(dāng)，說明利用本文方法計算海潮負(fù)荷位移參數(shù)是可行的。

2）除Q1分潮TPXO.7.2 和FES2004 與其他海潮模型有差異外，在中高緯地區(qū)本文所選5個海潮模型差異不明顯。

［1］周江存，孫和平.海潮負(fù)荷對GPS 基線的影響［J］.大地測量與地球動力學(xué)，2005，25（4）：27－32（Zhou Jiangcun，Sun Heping.Effect of Ocean Tide Loading on GPS Baseline Measurement［J］.Journal of Geodesy and Geodynamics，2005，25（4）：27－32）

［2］鄭祎，伍吉倉，王解先，等.GPS精密定位中的海潮位移改正［J］.武 漢 大 學(xué) 學(xué) 報：信 息 科 學(xué) 版，2003，28（4）：405－408（Zheng Yi，Wu Jicang，Wang Jiexian，et al.Ocean Tidal Displacement Corrections in GPS Precision Positioning［J］.Geomatics and Information Science of Wuhan University，2003，28（4）：405－408）

［3］Dragert H，James T S.Ocean Loading Corrections for Continuous GPS：A Case Study at the Canadian Coastal Site Holberg［J］.Geophysical Research Letters，2005，27（14）：2 045－2 048

［4］Zhang X H，Andersen O B.Surface Ice Flow Velocity and Tide Retrieval of the Amery Ice Shelf Using Precise Point Positioning［J］.J Geod，2006，80：171－176

［5］Thomas I D，Peter M A K，Clarke J.A Comparison of GPS，VLBI and Model Estimates of Ocean Tide Loading Displacement［J］.J Geod，2007，81（5）：359－368

［6］Yuan L G，Ding X L，Zhong P.Estimates of Ocean Tide Loading Displacement and Its Impact on Position Time Series in Hong Kong Using a Dense Continuous GPS Network［J］.J Geod，2009，83：999－1 015

［7］郗欽文.精密引潮位展開的精度評定［J］.地球物理學(xué)報，1992，35（2）：150－153（Xi Qinwen.The Evaluation on the Precision of the Development of the Tidal Generating Potential［J］.Acta Geophysica Sinica，1992，35（2）：150－153）

［8］周旭華，吳斌，李軍.高精度大地測量中的海潮位移改正［J］.測 繪 學(xué) 報，2001，30（4）：327－330（Zhou Xuhua，Wu Bin，Li Jun.Ocean Tide Displacement Corrections in High Precision Geodesy［J］.Acta Geodaetica et Cartographica Sinica，2001，30（4）：327－330）

［9］黃祖珂，黃磊.潮汐原理與計算［M］.青島：中國海洋大學(xué)出版社，2005（Huang Zuke，Huang Lei.Tidal Theory and Calculation［M］.Qindao：Publishing House of Ocean University of China，2005）

［10］Héroux P，Kouba J.Precise Point Positioning Using IGS Orbit and Clock Products［J］.GPS Solution，2001，5（2）：12－28

［11］Kouba J.A Guide to Using International GNSS Service（IGS）Products，Natural Resources Canada［EB／OL］.http：／／igscb.jpl.nasa.gov／components／usage.html，2012

［12］Zhang X H，Guo B F，Guo F.Influence of Clock Jump on the Velocity and Acceleration Estimation with a Single GPS Receiver Based on Carrier－Phase－Derived Doppler［J］.GPS Solutions，2013，17（4）：549－559

［13］Ruymbeke M V，Zhu P，CadicheanuNat N.Very Weak Signals（VWS）Detected by Stacking Method According to Different Astronomical Periodicities（HiCum）［J］.Hazards Earth Syst Sci，2007（7）：651－656

［14］King M.Kinematic and Static GPS Techniques for Estimating Tidal Displacements with Application to Antarctica［J］.J Geodyn，2006，41：77－86

［15］李大煒，李建成，金濤勇.利用驗潮站資料評估全球海潮模型的精度［J］.大地測量與地球動力學(xué)，2012，32（4）：106－110（Li Dawei，Li Jiancheng，Jin Taoyong.Accuracy Estimation of Recent Global Ocean Tide Models Using Tide Gauge Data［J］.Journal of Geodesy and Geodynamics，2012，32（4）：106－110）