牟麗爽, 馮金揚(yáng), 吳書(shū)清, 李春劍, 王啟宇

(中國(guó)計(jì)量科學(xué)研究院,北京 100029)

1 引 言

重力測(cè)量在地質(zhì)研究、資源勘探、地球物理、地震預(yù)報(bào)及計(jì)量科學(xué)領(lǐng)域有著廣泛的應(yīng)用[1～6]。因此重力加速度測(cè)量研究受各國(guó)普遍重視,且已成為世界各國(guó)大力發(fā)展的技術(shù)之一。廣泛的重力測(cè)量應(yīng)用需要重力計(jì)量支撐,其中重力儀溯源是重力計(jì)量的重要組成部分,而絕對(duì)重力儀比對(duì)是重力測(cè)量溯源至SI單位的有效途徑。2011年舉辦的歐洲區(qū)域絕對(duì)重力儀比對(duì),第1次使用可以連續(xù)觀測(cè)重力變化的相對(duì)重力儀—超導(dǎo)重力儀,將各臺(tái)絕對(duì)重力儀不同時(shí)間的測(cè)量值歸結(jié)到相同時(shí)刻。2017年10月,中國(guó)計(jì)量科學(xué)研究院(NIM)作為主導(dǎo)實(shí)驗(yàn)室舉辦了第十屆全球絕對(duì)重力儀國(guó)際比對(duì)[7,8]。在比對(duì)期間,在關(guān)鍵比對(duì)點(diǎn)安置了超導(dǎo)重力儀iGrav-012。由此可以鏈接各臺(tái)絕對(duì)重力儀在不同時(shí)刻的測(cè)量結(jié)果,最終得到關(guān)鍵比對(duì)參考值(KCRV),比對(duì)結(jié)束后,關(guān)鍵比對(duì)點(diǎn)連續(xù)觀測(cè)至今,是目前我國(guó)可以連續(xù)輸出重力變化并溯源至國(guó)際比對(duì)的關(guān)鍵比對(duì)。2019年7月,重力關(guān)鍵比對(duì)參考點(diǎn)及配套超導(dǎo)重力儀等設(shè)備獲批成為重力計(jì)量基準(zhǔn)點(diǎn),此后可以將實(shí)時(shí)的重力觀測(cè)值通過(guò)超導(dǎo)重力儀連續(xù)觀測(cè)的重力變化鏈接至國(guó)際比對(duì)參考值,絕對(duì)重力儀的實(shí)時(shí)測(cè)量結(jié)果可溯源至國(guó)際比對(duì)參考值。超導(dǎo)重力儀的輸出原始信號(hào)為電壓信號(hào)且存在漂移[9～12]。對(duì)其進(jìn)行格值并去除漂移率的影響后可得到連續(xù)重力變化,但此過(guò)程也會(huì)不可避免地引入不確定度。因此需要對(duì)關(guān)鍵比對(duì)點(diǎn)的連續(xù)重力觀測(cè)結(jié)果不確定度進(jìn)行評(píng)估,定期計(jì)算超導(dǎo)重力儀格值及漂移對(duì)關(guān)鍵比對(duì)點(diǎn)觀測(cè)結(jié)果引入的不確定度。

目前格值超導(dǎo)重力儀有4種方法,分別為使用已知的潮汐因子幅值標(biāo)定,使用校準(zhǔn)平臺(tái)標(biāo)定,使用大質(zhì)量塊標(biāo)定及使用絕對(duì)重力儀同址觀測(cè)標(biāo)定[13～17]。其中校準(zhǔn)平臺(tái)標(biāo)定及大質(zhì)量塊標(biāo)定對(duì)實(shí)驗(yàn)室要求較高,因此研究人員通常使用潮汐因子標(biāo)定法與使用絕對(duì)重力儀同址觀測(cè)法,通過(guò)擬合得到格值因子。Imanishi Y等通過(guò)實(shí)驗(yàn)得出結(jié)論,同期同址觀測(cè)法比潮波分析法得到的格值因子精確度高[18];Francis O通過(guò)連續(xù)9 d的同期同址觀測(cè)法得到相對(duì)不確定度為0.08%的格值因子[19];孫和平通過(guò)相同的方法也得到了符合精確度要求的格值因子[16]。格值因子及漂移率的計(jì)算方法已經(jīng)成熟,但重力計(jì)量關(guān)鍵比對(duì)點(diǎn)連續(xù)觀測(cè)結(jié)果不確定度計(jì)算并沒(méi)有相關(guān)研究。本文將使用絕對(duì)重力儀同址觀測(cè)標(biāo)定的方法計(jì)算關(guān)鍵比對(duì)點(diǎn)所使用的的超導(dǎo)重力儀的格值因子,同時(shí)采用較為常用的大時(shí)間跨度與絕對(duì)重力儀同址觀測(cè)的方法,計(jì)算該超導(dǎo)重力儀的漂移率[20～23],最終依據(jù)不確定度傳播律,分析關(guān)鍵比對(duì)點(diǎn)重力連續(xù)觀測(cè)結(jié)果的不確定度。

2 重力觀測(cè)及不確定度分析模型

關(guān)鍵比對(duì)點(diǎn)超導(dǎo)重力儀觀測(cè)重力連續(xù)變化結(jié)果可用式(1)表示[24]:

Δg=a·ΔV(t)+k·Δt-gtide(t)-

β·P(t)-gpolar(t)

(1)

式中:Δg為超導(dǎo)重力儀觀測(cè)連續(xù)重力變化;ΔV(t)為超導(dǎo)重力儀輸出電壓值相對(duì)于國(guó)際比對(duì)參考值零點(diǎn)電壓變化量;a為格值因子;k為超導(dǎo)重力儀漂移率;t為連續(xù)觀測(cè)時(shí)間；Δt為相對(duì)參考值零點(diǎn)的時(shí)間變化;gtide(t)為對(duì)重力潮汐模擬理論值;β=-3.0 (nm/s2)/hPa為大氣導(dǎo)納值;P(t)為氣壓觀測(cè)值；gpolar(t)為極移產(chǎn)生的重力變化。

依據(jù)不確定度傳播律,格值因子與漂移率對(duì)關(guān)鍵比對(duì)點(diǎn)連續(xù)重力觀測(cè)結(jié)果引入的不確定度分量分別為:

(2)

(3)

式中:ua為格值因子不確定度;uk為漂移率不確定度;u1為格值因子對(duì)連續(xù)觀測(cè)引入的不確定度;u2為漂移率對(duì)連續(xù)觀測(cè)引入的不確定度。

其合成標(biāo)準(zhǔn)不確定度為:

(4)

超導(dǎo)重力儀格值采用與絕對(duì)重力儀同址觀測(cè)最小二乘加權(quán)擬合的方法[25]。則擬合模型為:

gy=a·gv+b

(5)

式中:gy為FG5絕對(duì)重力儀觀測(cè)值;gv為超導(dǎo)重力儀觀測(cè)值;b為擬合常數(shù)。

分別對(duì)絕對(duì)重力儀和超導(dǎo)重力儀進(jìn)行了n次獨(dú)立觀測(cè),得到n對(duì)觀測(cè)值(gv,i,gy,i),i=1,2，…，n,代入式(5),則有:

(6)

(7)

(8)

超導(dǎo)重力儀漂移率的擬合模型為:

gFG5(t)-giGrav(t)=k·t+b1

(9)

式中:t為觀測(cè)時(shí)間;gFG5(n)為絕對(duì)重力儀FG5觀測(cè)結(jié)果;giGrav(n)為超導(dǎo)重力儀格值后結(jié)果;k為超導(dǎo)重力儀漂移率;b1為擬合常數(shù)。

將觀測(cè)數(shù)據(jù)(gFG5(ti),giGrav(ti))代入式(9),得到:

(10)

(11)

則漂移率的不確定度為:

(12)

由式(2)～式(4)、式(8)、式(12)可得格值因子與漂移率對(duì)連續(xù)觀測(cè)引入的合成標(biāo)準(zhǔn)不確定度為:

(13)

式(13)花括號(hào)中第1項(xiàng)為格值因子對(duì)漂移率引入的不確定度;第2項(xiàng)為漂移率對(duì)連續(xù)觀測(cè)引入的不確定度。由式(13)可以看出,參與擬合時(shí)間越長(zhǎng),連續(xù)觀測(cè)時(shí)間越短,其引入的不確定度越小。因此可通過(guò)定期參與比對(duì)來(lái)減小漂移率對(duì)連續(xù)觀測(cè)結(jié)果引入的不確定度。

3 關(guān)鍵比對(duì)參考點(diǎn)連續(xù)觀測(cè)結(jié)果

NIM的超導(dǎo)重力儀(iGrav-012)于2014年6月開(kāi)始正式運(yùn)行使用,2017年由NIM昌平院區(qū)21號(hào)樓地上一層搬至地下一層比對(duì)大廳。FG5型絕對(duì)重力儀安置于8號(hào)關(guān)鍵比對(duì)點(diǎn),與超導(dǎo)重力儀相距約3 m,可認(rèn)為處于相同重力場(chǎng)(二者的相對(duì)位置如圖1所示)。FG5型絕對(duì)重力儀定期在環(huán)境干擾較少的情況下測(cè)量,每1 h測(cè)量1組,每組落體下落100次,單次落體下落時(shí)間為10 s,每組結(jié)束后暫停0.5 h。

圖1 比對(duì)大廳平面圖Fig.1 Plane figure of comparison site

超導(dǎo)重力儀連續(xù)輸出觀測(cè)結(jié)果的頻率為1 Hz,在記錄重力變化的同時(shí)也會(huì)記錄環(huán)境中的噪聲。為得到精確度更高的重力觀測(cè)值,可以使用標(biāo)準(zhǔn)的重力數(shù)據(jù)預(yù)處理軟件Tsoft對(duì)數(shù)據(jù)進(jìn)行預(yù)處理[26,27],以減小觀測(cè)數(shù)據(jù)中的噪聲及其對(duì)格值及漂移率計(jì)算的干擾,具體過(guò)程如下:首先,剔除原始數(shù)據(jù)中的干擾信號(hào),并插值修補(bǔ)較小的中斷;其次,使用截止頻率為0.004 166 7 Hz的最小二乘法濾波器對(duì)信號(hào)進(jìn)行濾波,采樣率降為1 min。由此可以得到預(yù)處理后的重力觀測(cè)序列。從2017年經(jīng)過(guò)預(yù)處理的超導(dǎo)重力儀觀測(cè)數(shù)據(jù)與絕對(duì)重力儀同址觀測(cè)數(shù)據(jù)可以看出(如圖2所示),關(guān)鍵比對(duì)點(diǎn)使用的這2種儀器具有很好的一致性。

圖2 同址觀測(cè)數(shù)據(jù)Fig.2 Same site observation data

4 分析與討論

表1列出了2014年至2020年超導(dǎo)重力儀與絕對(duì)重力儀定期同址觀測(cè)的數(shù)據(jù),共10組。

由于大潮前后重力變化較大,可獲得較大的峰峰值以減小格值不確定度,因此絕對(duì)重力儀的數(shù)據(jù)均在此期間測(cè)量得到[6]。2016年實(shí)驗(yàn)樓原址在擴(kuò)建,周圍施工會(huì)對(duì)重力儀測(cè)量重力值產(chǎn)生影響,但儀器仍處于工作中。

格值因子的計(jì)算結(jié)果如圖3所示。10組格值因子的加權(quán)均值為-927.96(nm/s2)/V,加權(quán)均值的標(biāo)準(zhǔn)偏差為0.25(nm/s2)/V。格值因子的標(biāo)準(zhǔn)不確定度為0.80(nm/s2)/V。大潮期間連續(xù)25 h格值的相對(duì)不確定度可達(dá)0.14%,年變化率最大為0.2%。由于超導(dǎo)重力儀觀測(cè)重力變化輸出電壓峰峰值為4 V,依據(jù)式(2)可得格值因子對(duì)超導(dǎo)重力儀連續(xù)觀測(cè)產(chǎn)生的標(biāo)準(zhǔn)不確定度分量為:

=3.20 nm/s2

(14)

表1 超導(dǎo)重力儀與FG5同期同址觀測(cè)數(shù)據(jù)Tab.1 Observation schedule of superconducting gravimeter and FG5 at the same site at the same time

圖3 格值結(jié)果Fig.3 The result of calibration factors

進(jìn)一步,結(jié)合FG5型絕對(duì)重力儀的觀測(cè)結(jié)果來(lái)計(jì)算超導(dǎo)重力儀的漂移率。首先,采用表1所示10組數(shù)據(jù)加權(quán)計(jì)算得到的均值為格值因子,將預(yù)處理后的超導(dǎo)重力儀輸出值轉(zhuǎn)化為重力值；接著，根據(jù)各自的測(cè)量時(shí)間,對(duì)重力值以及絕對(duì)重力儀的測(cè)量值分別進(jìn)行潮汐、大氣壓及極移修正。其結(jié)果如圖4所示。超導(dǎo)重力儀在最開(kāi)始工作的半年時(shí)間內(nèi)漂移較大,這是由于其在此期間還未進(jìn)入平穩(wěn)運(yùn)行狀態(tài)。2017年年中超導(dǎo)重力儀由一層搬運(yùn)至比對(duì)大廳,且周圍設(shè)施不全,FG5仍在一層進(jìn)行實(shí)驗(yàn),因此這期間沒(méi)有參與擬合。隨使用時(shí)間的增加,超導(dǎo)重力儀的漂移也將逐漸變小,最后與絕對(duì)重力儀的觀測(cè)結(jié)果基本一致。可以認(rèn)為需要近半年的時(shí)間才可以保證新啟用的超導(dǎo)重力儀漂移較小;若要使用其工作半年內(nèi)的數(shù)據(jù)進(jìn)行后續(xù)分析,則應(yīng)單獨(dú)計(jì)算在此期間的漂移率。此外,雖然2017年5月超導(dǎo)重力儀進(jìn)行了搬運(yùn),但超導(dǎo)重力儀的漂移率并未因此產(chǎn)生明顯增大,可以在擬合計(jì)算時(shí)認(rèn)為其保持不變。

圖4 長(zhǎng)期同址觀測(cè)結(jié)果Fig.4 Long term observation results at the same site

假定搬運(yùn)前后超導(dǎo)重力儀具有相同的漂移率,把超導(dǎo)重力儀搬運(yùn)前后數(shù)據(jù)分別代入式(9)進(jìn)行擬合,可得超導(dǎo)重力儀漂移率為-14.66(nm/s2)/年,擬合標(biāo)準(zhǔn)不確定度為uk1=0.03 (nm/s2)/年。參與擬合的數(shù)據(jù)時(shí)間長(zhǎng)度為5.5 年,因此漂移率的標(biāo)準(zhǔn)不確定度為:

(15)

重力關(guān)鍵比對(duì)點(diǎn)連續(xù)觀測(cè)起點(diǎn)為2017年10月14日,連續(xù)觀測(cè)時(shí)間長(zhǎng)度為t=3年。依據(jù)式(3)和式(13),可以得到超導(dǎo)重力儀的漂移率對(duì)關(guān)鍵比對(duì)點(diǎn)連續(xù)觀測(cè)引入的標(biāo)準(zhǔn)不確定度為:

u2=|t|·uk=6.30 nm/s2

(16)

最后,關(guān)鍵比對(duì)點(diǎn)超導(dǎo)重力儀連續(xù)重力觀測(cè)結(jié)果合成標(biāo)準(zhǔn)不確定度為:

(17)

關(guān)鍵比對(duì)點(diǎn)重力連續(xù)觀測(cè)結(jié)果標(biāo)準(zhǔn)不確定度為7.07 nm/s2,主要來(lái)源于超導(dǎo)重力儀的漂移率,其引入的標(biāo)準(zhǔn)不確定度6.30 nm/s2。

根據(jù)式(15)～式(17)可知,漂移率不確定度主要與參與擬合的數(shù)據(jù)時(shí)間長(zhǎng)度和基準(zhǔn)連續(xù)觀測(cè)的時(shí)間長(zhǎng)度有關(guān)。因此,可以通過(guò)定期參加比對(duì)、得到更多的比對(duì)參考值、減小連續(xù)觀測(cè)的時(shí)間長(zhǎng)度等方式減小漂移率引入的不確定度,從而減小關(guān)鍵比對(duì)點(diǎn)重力連續(xù)觀測(cè)值的不確定度。

完成格值因子及漂移率的計(jì)算可得到關(guān)鍵比對(duì)點(diǎn)的重力連續(xù)觀測(cè)結(jié)果,如圖5所示。由圖5去除漂移前后對(duì)比可知,處理后的超導(dǎo)重力儀觀測(cè)結(jié)果與FG5觀測(cè)結(jié)果幾乎吻合,說(shuō)明格值計(jì)算及漂移率計(jì)算結(jié)果可靠準(zhǔn)確。

圖5 超導(dǎo)重力儀連續(xù)觀測(cè)結(jié)果Fig.5 Continuous observation results of superconducting gravimeter

5 結(jié) 論

對(duì)我國(guó)重力計(jì)量基準(zhǔn)中用于關(guān)鍵比對(duì)點(diǎn)連續(xù)觀測(cè)的超導(dǎo)重力儀的測(cè)量結(jié)果不確定度進(jìn)行了計(jì)算分析。首先,根據(jù)7年共10組絕對(duì)重力儀與超導(dǎo)重力儀同期同址的觀測(cè)結(jié)果,使用最小二乘加權(quán)擬合的方法,得到該超導(dǎo)重力儀的格值因子為(-928.00±0.222) (nm/s2)/V。其次,通過(guò)線性擬合對(duì)超導(dǎo)重力儀與絕對(duì)重力儀長(zhǎng)期同址觀測(cè)的結(jié)果進(jìn)行分析,得到超導(dǎo)重力儀的漂移率為(-14.66±2.1) (nm/s2)/年。由此可得格值因子與漂移率對(duì)重力連續(xù)觀測(cè)結(jié)果引入的合成標(biāo)準(zhǔn)不確定度為7.07 nm/s2,主要來(lái)源于漂移率。在之后的工作中,可以通過(guò)定期參與比對(duì)來(lái)減小關(guān)鍵比對(duì)點(diǎn)重力連續(xù)觀測(cè)結(jié)果的不確定度。