要佳敏 莊偉 馮金揚(yáng) 王啟宇 趙陽(yáng) 王少凱 吳書清 李天初

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

1 引言

絕對(duì)重力測(cè)量,即對(duì)重力加速度絕對(duì)值的精確測(cè)量,在輔助導(dǎo)航、資源勘探、地球物理、計(jì)量領(lǐng)域都有著廣泛的應(yīng)用[1-4].現(xiàn)有的高精度絕對(duì)重力儀普遍采用自由落體法,主要包括以角錐棱鏡為下落物體的激光干涉式絕對(duì)重力儀[5-8]和以冷原子團(tuán)為下落物體的原子干涉式絕對(duì)重力儀[9-14],以下分別簡(jiǎn)稱為光學(xué)重力儀和原子重力儀.前者通過(guò)激光干涉儀測(cè)量角錐棱鏡在真空中自由下落的運(yùn)動(dòng)軌跡,擬合求解重力加速度,后者根據(jù)原子團(tuán)在自由下落的過(guò)程中經(jīng)過(guò)3 次拉曼激光脈沖作用形成的原子干涉條紋來(lái)計(jì)算重力加速度.由于二者的實(shí)際測(cè)量對(duì)象是下落物體或下落原子團(tuán)相對(duì)于儀器內(nèi)一個(gè)參考鏡(角錐棱鏡或平面鏡)的運(yùn)動(dòng)加速度,因此參考鏡自身的加速度必然耦合到測(cè)量結(jié)果中.理論上,忽略參考鏡自身振動(dòng)時(shí),現(xiàn)有的絕對(duì)重力儀的測(cè)量精度可以達(dá)到微伽(μGal,1 μGal=1 ×10—8m·s—2)量級(jí);但實(shí)際上,如果直接將參考鏡放置在地面上,地面振動(dòng)引起的測(cè)值離散度將超過(guò)10 μGal,在復(fù)雜振動(dòng)環(huán)境中甚至達(dá)到mGal 或10 mGal 量級(jí).因此,利用測(cè)振儀器獲得參考鏡自身的運(yùn)動(dòng)加速度并在計(jì)算重力加速度時(shí)對(duì)其影響進(jìn)行修正是提高絕對(duì)重力儀測(cè)量精度的有效手段之一,這種修正方法一般稱為振動(dòng)補(bǔ)償.

具體而言,影響原子重力儀的地面振動(dòng)噪聲主要為環(huán)境噪聲,與測(cè)量時(shí)間和測(cè)量點(diǎn)的地基、坐標(biāo)及環(huán)境有關(guān),可以根據(jù)頻率來(lái)粗略區(qū)分.首先是人類活動(dòng)噪聲,一般大于1 Hz,如儀器測(cè)試期間附近的人員走動(dòng),車輛等人工振源產(chǎn)生的振動(dòng),以及建筑物的晃動(dòng)等.其次是來(lái)自地球本身的噪聲,一般在0.1—1.0 Hz 之間,主要包括周期性的脈動(dòng)和地震等地球內(nèi)部運(yùn)動(dòng)引起的振動(dòng),其中地脈動(dòng)噪聲模型的加速度功率譜密度在頻率為0.2 Hz 和3.0 Hz的位置存在峰值[15].最后是氣壓波動(dòng)等大氣運(yùn)動(dòng)噪聲,一般小于0.1 Hz,普遍變化緩慢且幅值較低,對(duì)絕對(duì)重力測(cè)量的影響可以忽略[16,17].由于絕對(duì)重力測(cè)量只依靠垂直方向的位移數(shù)據(jù)進(jìn)行計(jì)算,因此對(duì)原子干涉條紋進(jìn)行振動(dòng)補(bǔ)償時(shí),一般忽略地面振動(dòng)的水平分量,僅需使用測(cè)振儀器輸出的垂直振動(dòng)信號(hào).

目前國(guó)際上較為成熟的原子重力儀是法國(guó)巴黎天文臺(tái)研制的CAG-01 型原子重力儀,曾參加過(guò)多次國(guó)際重力比對(duì),采用的振動(dòng)處理方法為被動(dòng)式垂直隔振系統(tǒng)與振動(dòng)補(bǔ)償?shù)慕Y(jié)合,其中振動(dòng)補(bǔ)償?shù)木唧w原理為:數(shù)據(jù)采集卡同步記錄原子干涉過(guò)程中的地震計(jì)輸出的速度信號(hào)與原子干涉條紋,在控制軟件中利用無(wú)限沖激響應(yīng)(infinite impulse response,IIR)濾波器和非因果低通濾波器對(duì)地震計(jì)信號(hào)進(jìn)行處理,以減小地震計(jì)傳遞函數(shù)幅頻特性不恒定、相頻特性存在非線性的影響,獲得更真實(shí)的參考鏡振動(dòng)信號(hào),最后利用處理后的信號(hào)得到由參考鏡振動(dòng)引入的相位噪聲,修正原子干涉條紋.在該儀器所在的測(cè)量環(huán)境中,此振動(dòng)補(bǔ)償方法本身可以使重力儀靈敏度提升為原來(lái)的3 倍[11].德國(guó)漢諾威大學(xué)提出利用數(shù)字式低通和高通濾波器對(duì)振動(dòng)傳感器的輸出信號(hào)進(jìn)行前處理,以得到更為精確的參考鏡振動(dòng)引入的相位噪聲,從而實(shí)現(xiàn)振動(dòng)補(bǔ)償.其實(shí)驗(yàn)結(jié)果表明,在較安靜的環(huán)境中利用高精度商用地震計(jì)進(jìn)行振動(dòng)補(bǔ)償,或在復(fù)雜振動(dòng)環(huán)境中使用精度較低的商用加速度計(jì)或其自主研制的新型光學(xué)慣性傳感器進(jìn)行補(bǔ)償,均可以有效提高重力儀的靈敏度[18].在我國(guó),浙江大學(xué)自制的原子重力儀同樣采用了振動(dòng)補(bǔ)償方法,其主要特點(diǎn)為利用商用反演軟件修正地震計(jì)實(shí)測(cè)傳遞函數(shù)的影響以獲得更真實(shí)的參考鏡振動(dòng)信號(hào).在其實(shí)驗(yàn)室環(huán)境中,該方法可以將重力儀的靈敏度提升至與使用被動(dòng)隔振系統(tǒng)時(shí)相近的水平[19].國(guó)防科技大學(xué)對(duì)參考鏡垂向振動(dòng)及水平偏轉(zhuǎn)引入的測(cè)量噪聲進(jìn)行了詳細(xì)分析,利用加速度計(jì)測(cè)量原子重力儀工作時(shí)的參考鏡振動(dòng),同時(shí)對(duì)加速度計(jì)與參考鏡之間的傳遞函數(shù)進(jìn)行了較為精確的標(biāo)定.實(shí)測(cè)結(jié)果表明在較復(fù)雜的振動(dòng)環(huán)境中其振動(dòng)補(bǔ)償方法可以實(shí)現(xiàn)顯著的振動(dòng)補(bǔ)償效果[20].軍事科學(xué)院國(guó)防創(chuàng)新院提出的振動(dòng)補(bǔ)償方法則采用具有高精度四通道移相探測(cè)器的邁克耳孫激光干涉儀替代傳統(tǒng)振動(dòng)傳感器.該干涉儀可以直接測(cè)量出原子干涉條紋的相移,并結(jié)合信號(hào)整形、數(shù)字鑒相等方法確定拉曼激光脈沖作用前后原子干涉條紋的相位小數(shù)的變化,得到參考鏡振動(dòng)引起的原子干涉相位偏差.去除該相位偏差,即可對(duì)重力儀測(cè)量到的躍遷概率信號(hào)進(jìn)行修正[21].此外,中國(guó)計(jì)量科學(xué)研究院、中國(guó)科學(xué)院精密測(cè)量研究院、華中科技大學(xué)、中國(guó)科技大學(xué)、國(guó)防科技大學(xué)、美國(guó)斯坦福大學(xué)、德國(guó)洪堡大學(xué)等單位也研制了不同類型的隔振系統(tǒng)[12,13,22-26],以減小地面振動(dòng)對(duì)原子重力儀的影響.

德國(guó)漢諾威大學(xué)的Richardson 等[18]和我國(guó)國(guó)防科技大學(xué)的研究人員[20]均指出,傳統(tǒng)振動(dòng)傳感器與參考鏡由于安裝位置不同形成的機(jī)械傳遞函數(shù)可能是影響振動(dòng)補(bǔ)償精度的重要因素;在原子重力儀進(jìn)行正式測(cè)量前對(duì)該傳遞函數(shù)中的關(guān)鍵系數(shù)進(jìn)行初始標(biāo)定有助于提高振動(dòng)補(bǔ)償效果.這些分析為提升原子重力儀的振動(dòng)補(bǔ)償精度提供了重要指導(dǎo),但目前尚無(wú)具有普適性的、能夠在任意測(cè)量時(shí)間及地點(diǎn)快速實(shí)現(xiàn)對(duì)傳遞函數(shù)進(jìn)行實(shí)時(shí)修正的具體振動(dòng)補(bǔ)償算法.而在光學(xué)重力儀方面,清華大學(xué)的研究人員給出了對(duì)重力測(cè)量結(jié)果進(jìn)行實(shí)時(shí)修正的基于相同簡(jiǎn)化模型的具體振動(dòng)補(bǔ)償算法,并提供了該算法在不同地點(diǎn)對(duì)延時(shí)系數(shù)和增益系數(shù)進(jìn)行實(shí)時(shí)計(jì)算的結(jié)果及重力測(cè)量的實(shí)測(cè)結(jié)果[27-29].他們采用的傳感器為高精度地震計(jì),利用相關(guān)分析法或黃金分割法搜索出當(dāng)前環(huán)境下這2 個(gè)系數(shù)的具體值,并對(duì)下落物體軌跡信號(hào)進(jìn)行修正.對(duì)同一套振動(dòng)補(bǔ)償算法,輸入在中國(guó)計(jì)量科學(xué)研究院昌平院區(qū)的安靜環(huán)境中測(cè)得的重力數(shù)據(jù)和地震計(jì)數(shù)據(jù),算法的直接輸出結(jié)果表明振動(dòng)補(bǔ)償后光學(xué)重力儀的測(cè)量分辨率可以提升5 倍;輸入在清華大學(xué)所處的復(fù)雜振動(dòng)環(huán)境中的數(shù)據(jù),結(jié)果表明振動(dòng)補(bǔ)償可將重力儀的測(cè)量分辨率提升16 倍,效果更為顯著[27].該算法最大的特點(diǎn)是較強(qiáng)的自適應(yīng)性,在不同時(shí)刻和不同地點(diǎn)的測(cè)量過(guò)程中可以直接使用,以保證真實(shí)的機(jī)械傳遞函數(shù)因環(huán)境變化而發(fā)生改變時(shí)對(duì)其當(dāng)前值進(jìn)行實(shí)時(shí)估算,實(shí)現(xiàn)當(dāng)時(shí)當(dāng)?shù)刈詈玫恼駝?dòng)補(bǔ)償效果.受上述研究結(jié)果的啟發(fā),基于上述簡(jiǎn)化模型,本文給出了一套振動(dòng)補(bǔ)償算法的工作原理和具體流程,該振動(dòng)補(bǔ)償算法可以應(yīng)用于原子重力儀，對(duì)傳遞函數(shù)進(jìn)行實(shí)時(shí)標(biāo)定.我們通過(guò)仿真運(yùn)算對(duì)算法的有效性進(jìn)行了詳細(xì)驗(yàn)證,最后利用實(shí)驗(yàn)室自主研制的原子重力儀實(shí)際評(píng)估了該算法對(duì)原子干涉條紋擬合精度和重力測(cè)值離散度的補(bǔ)償效果,為提高原子重力儀在不同環(huán)境甚至是復(fù)雜振動(dòng)環(huán)境中的測(cè)量精度提供幫助.

2 原 理

2.1 應(yīng)用于原子重力儀的振動(dòng)補(bǔ)償原理

原子重力儀基于冷原子物質(zhì)波干涉原理,將激光脈沖作用到冷原子團(tuán),通過(guò)雙光子受激拉曼躍遷或多光子布拉格衍射過(guò)程來(lái)實(shí)現(xiàn)原子干涉過(guò)程[9].以基于受激拉曼躍遷的原子重力儀為例,冷原子團(tuán)在真空腔中自由下落,真空腔下方的平面參考鏡將垂直向下的拉曼激光原路反射,使兩束對(duì)射的激光同時(shí)作用于原子團(tuán),在此過(guò)程中3 個(gè)連續(xù)的激光脈沖實(shí)現(xiàn)了原子波包的分束、反轉(zhuǎn)和合束,分別為第1 個(gè)π/2 脈沖、π 脈沖和第2 個(gè)π/2 脈沖,相鄰脈沖的時(shí)間間隔為T[30].合束后任一能態(tài)的原子均由2 條路徑上的原子疊加而成,通過(guò)熒光信號(hào)可以探測(cè)脈沖作用后原子團(tuán)的躍遷概率P,即干涉條紋,滿足

式中:ΔФ為原子相位,與當(dāng)?shù)氐闹亓铀俣萭有關(guān);代表?xiàng)l紋對(duì)比度的系數(shù)A和B可以通過(guò)對(duì)干涉條紋進(jìn)行余弦擬合得到.理想情況下,原子相位ΔФ等于其理論值ΔФth,即

式中,keff為有效波矢,α表示激光的掃頻速率.重力儀的控制軟件調(diào)控 (2) 式中的α線性增大,則原子相位ΔФ也將線性變化,使干涉條紋呈現(xiàn)為標(biāo)準(zhǔn)余弦信號(hào)的形式.此時(shí)P應(yīng)與由P擬合出的標(biāo)準(zhǔn)信號(hào)Pfit完全相同,有

但實(shí)際情況下,受到各種相位噪聲的影響,原子相位ΔФ等于其實(shí)際值ΔФm,有

式中Δφvib為參考鏡振動(dòng)引入的相位噪聲,以下簡(jiǎn)稱為振動(dòng)相位;Δφothers為其他噪聲源引入的相位噪聲的總和,以下簡(jiǎn)稱為其他相位噪聲.此時(shí)P由Pfit變?yōu)?/p>

當(dāng)振動(dòng)相位在相位噪聲中占主導(dǎo)時(shí),依據(jù)上述原理可以確定修正干涉條紋P的振動(dòng)補(bǔ)償過(guò)程.第一步,對(duì)原始實(shí)測(cè)條紋P(ΔФth)進(jìn)行余弦擬合,得到擬合出的理論條紋Pfit(ΔФth)和(1)式中的系數(shù)A和B.需要說(shuō)明的是,此時(shí)可以將原子重力儀在正式測(cè)量前的初始測(cè)量過(guò)程中獲得的粗測(cè)值gth代入(2)式以得到擬合條紋Pfit.這是因?yàn)樵又亓x到達(dá)一個(gè)新的測(cè)量地點(diǎn)后一般都需要在啟動(dòng)階段通過(guò)改變拉曼間隔T并掃描掃頻速率α來(lái)獲得正式測(cè)量所用的中心值αcenter,而該測(cè)點(diǎn)的粗測(cè)重力值即為gth=αcenter/keff.此后每一次或每一組躍遷概率數(shù)據(jù)都對(duì)應(yīng)真實(shí)重力值相對(duì)于該粗測(cè)值的偏差Δg,所以實(shí)際上重力儀最終輸出的高精度測(cè)值為g=gth+Δg.因此利用粗測(cè)值gth計(jì)算擬合條紋有利于提高計(jì)算效率,這也符合原子重力儀的一般測(cè)試流程[19,20].第二步,通過(guò)振動(dòng)傳感器獲得一段時(shí)間t內(nèi)參考鏡的真實(shí)振動(dòng)速度vm(t),進(jìn)而根據(jù)

計(jì)算出真實(shí)的振動(dòng)相位Δφvib,其中S(t)為重力儀的靈敏度函數(shù),T1和T3分別為第1 個(gè)π/2 脈沖和第2 個(gè)π/2 脈沖的作用時(shí)刻,均由重力儀控制軟件的已知參數(shù)決定.第三步,計(jì)算出僅受振動(dòng)影響時(shí)的原子相位ΔФv和躍遷概率Pv,即

此時(shí)重新以ΔФv為橫坐標(biāo)、P為縱坐標(biāo)得到的數(shù)據(jù)對(duì)P(ΔФv)即為修正后的干涉條紋,至此補(bǔ)償過(guò)程結(jié)束.如圖1 所示,與藍(lán)色圓圈表示的原始條紋P(ΔФth)相比,紅色方塊表示的修正條紋P(ΔФv)的橫坐標(biāo)有所移動(dòng),理論上應(yīng)更接近灰色的標(biāo)準(zhǔn)余弦信號(hào),即擬合出的理論條紋Pfit(ΔФth).也就是說(shuō),對(duì)修正后的條紋P(ΔФv)再次進(jìn)行余弦擬合時(shí)的均方根誤差(root mean square error,RMSE)值應(yīng)小于對(duì)原始條紋P(ΔФth)進(jìn)行擬合時(shí)的RMSE 值.

圖1 原子干涉條紋示意圖Fig.1.Schematic diagram of the fringe signal of atomic interferometer.

2.2 基于傳遞函數(shù)簡(jiǎn)化模型的振動(dòng)補(bǔ)償方法

基于傳遞函數(shù)簡(jiǎn)化模型以實(shí)現(xiàn)上述振動(dòng)補(bǔ)償過(guò)程的方法包括硬件和軟件兩方面.

硬件方面,由于目前尚無(wú)以原子重力儀的參考鏡為敏感質(zhì)量的地震計(jì),因此一般而言只能用一臺(tái)獨(dú)立的地震計(jì)或加速度計(jì)放在參考鏡旁邊或下方測(cè)量其振動(dòng).然而,受機(jī)械結(jié)構(gòu)及地震計(jì)自身性能的限制,地震計(jì)輸出信號(hào)Us(t)并不能完全真實(shí)地反映參考鏡的振動(dòng)速度vm(t),二者的關(guān)系如圖2所示,其中Ga表示從地面振動(dòng)速度到參考鏡振動(dòng)速度的傳遞函數(shù),Gb表示從地面振動(dòng)速度到地震計(jì)敏感質(zhì)量振動(dòng)速度的傳遞函數(shù),Gc表示從敏感質(zhì)量運(yùn)動(dòng)速度到地震計(jì)輸出電壓之間的傳遞函數(shù).由于地震計(jì)內(nèi)部的敏感質(zhì)量與參考鏡的實(shí)際位置必然存在水平和垂直差異,因此Ga與Gb不可能相等,且當(dāng)?shù)鼗牧细淖儠r(shí),這兩個(gè)傳遞函數(shù)本身也可能發(fā)生變化.另一方面,地震計(jì)自身的傳遞函數(shù)Gc并非幅頻特性恒定、相頻特性為線性的理想傳遞函數(shù)[31],且有可能隨測(cè)量時(shí)間與環(huán)境的改變而發(fā)生變化(這是因?yàn)闇囟?、濕度、氣壓等環(huán)境參數(shù)會(huì)對(duì)地震計(jì)內(nèi)部的機(jī)電模塊參數(shù)產(chǎn)生影響,且這些參數(shù)可能隨時(shí)間發(fā)生漂移).綜上,參考鏡的真實(shí)振動(dòng)速度與地震計(jì)輸出的電壓信號(hào)之間的傳遞函數(shù)H滿足

圖2 地震計(jì)輸出信號(hào)與參考鏡真實(shí)振動(dòng)的關(guān)系Fig.2.Transfer function between the vibration of reference mirror and the output signal of seismometer.

因此對(duì)傳遞函數(shù)H的求解將直接影響振動(dòng)補(bǔ)償?shù)男Ч?

雖然實(shí)際情況中Ga,Gb,Gc不易測(cè)量且可能改變,導(dǎo)致測(cè)量人員無(wú)法得到總傳遞函數(shù)H的精確值,但由于單次重力測(cè)量任務(wù)中測(cè)量環(huán)境不會(huì)有顯著變化,且測(cè)量時(shí)間一般不超過(guò)24 h,因此可以利用由一個(gè)延時(shí)系數(shù)τ和一個(gè)增益系數(shù)K組成的簡(jiǎn)化模型來(lái)估算H,并認(rèn)為這兩個(gè)系數(shù)的數(shù)值在當(dāng)前測(cè)量過(guò)程中保持不變[27].由此可得圖2 的簡(jiǎn)化模型為

式中,Ks為地震計(jì)的標(biāo)稱靈敏度,由制造商提供.

理論上當(dāng)延時(shí)系數(shù)τ和增益系數(shù)K取最接近真實(shí)值時(shí),利用地震計(jì)輸出電壓Us和(9)式計(jì)算出的參考鏡振動(dòng)速度vm(t)也最接近真實(shí)值,對(duì)條紋的修正效果最好.因此振動(dòng)補(bǔ)償可以選取對(duì)條紋進(jìn)行擬合時(shí)的RMSE 值作為優(yōu)化目標(biāo),基本思路為:給定τ的取值范圍,遍歷搜索使RMSE 值達(dá)到最小值時(shí)對(duì)應(yīng)的最優(yōu)延時(shí)系數(shù)τopt;再給定K的取值范圍,在設(shè)定延時(shí)為τopt的情況下,再次遍歷搜索使RMSE 達(dá)到最小值時(shí)對(duì)應(yīng)的最優(yōu)增益系數(shù)Kopt.最后,在設(shè)定延時(shí)為τopt、增益為Kopt的情況下,得到修正后的條紋P(ΔФv)及其余弦擬合時(shí)的RMSE 值.

由此得到具體的算法流程如圖3 所示,相應(yīng)的處理程序即為振動(dòng)補(bǔ)償方法的軟件部分.上述分析中使用的主要物理符號(hào)匯總在表1 中.

表1 主要物理符號(hào)描述Table 1.Description of main symbols.

圖3 基于系數(shù)搜索的振動(dòng)補(bǔ)償算法流程Fig.3.The algorithm of the vibration correction method based on element searching.

2.3 振動(dòng)補(bǔ)償前后干涉條紋擬合殘差的變化

上述振動(dòng)補(bǔ)償算法的優(yōu)化目標(biāo)為擬合干涉條紋時(shí)的RMSE 值.該值具有簡(jiǎn)單直觀的數(shù)學(xué)含義,但物理意義不明顯.類比基于相同簡(jiǎn)化模型的振動(dòng)補(bǔ)償方法在光學(xué)重力儀中的算法流程[27],可以用另一種更具有物理意義的變量作為振動(dòng)補(bǔ)償效果的最終評(píng)價(jià)指標(biāo),即干涉條紋擬合殘差的變化.

運(yùn)行振動(dòng)補(bǔ)償算法前,先計(jì)算出原始條紋P(ΔФth)相較于理論條紋Pfit(ΔФth)的殘差Rm(ΔФth)的標(biāo)準(zhǔn)差σm,即

運(yùn)行振動(dòng)補(bǔ)償算法后,利用搜索出的最優(yōu)延時(shí)系數(shù)τopt和增益系數(shù)Kopt,根據(jù)(6)式,(7)式和(9)式得到僅受振動(dòng)影響的干涉條紋Pv相較于理論條紋Pfit的殘差Rv.同時(shí),記原始條紋P與僅受振動(dòng)影響的條紋Pv的差值為Rc.顯然Rc也就是原始條紋殘差Rm與僅受振動(dòng)影響的條紋殘差Rv的差值,相當(dāng)于去除振動(dòng)影響后的條紋殘差.理論上振動(dòng)補(bǔ)償效果越好,則振動(dòng)相位Δφvib的計(jì)算越準(zhǔn)確,推算出的僅由振動(dòng)引入的條紋殘差Rv也越接近原始條紋殘差Rm.此時(shí)Rc將趨近于僅由其他相位噪聲Δφothers引入的條紋殘差,其標(biāo)準(zhǔn)差σc也應(yīng)趨于最小值.因此可以以

相較原始條紋殘差Rm的標(biāo)準(zhǔn)差σm的衰減程度作為衡量振動(dòng)補(bǔ)償效果的最終指標(biāo).

3 仿真驗(yàn)證

上述振動(dòng)補(bǔ)償算法的效果將在實(shí)驗(yàn)室自主研發(fā)的NIM-AGRb-1 型原子重力儀[12]上進(jìn)行考察.實(shí)際測(cè)量之前,可以利用美國(guó) MathWorks 公司生產(chǎn)的 MATLAB 軟件驗(yàn)證算法的有效性,具體的仿真流程為:1) 參考NIM-AGRb-1 型重力儀的控制參數(shù)和之前的實(shí)測(cè)數(shù)據(jù),給定合理的計(jì)算參數(shù)值,包括波矢keff、粗測(cè)重力值gth、靈敏度函數(shù)S(t)、拉曼激光的脈沖間隔T=80 ms、1 min 內(nèi)線性增大的掃頻速率序列α(時(shí)間間隔為2 s,相當(dāng)于原子團(tuán)的下落間隔)、決定干涉條紋對(duì)比度的系數(shù)A和B等;2) 給定1 min 內(nèi)代表隨機(jī)振動(dòng)的速度信號(hào)的幅值A(chǔ)rand和代表地脈動(dòng)的速度信號(hào)的頻率序列fs及相應(yīng)的幅值序列As,具體取值參考地脈動(dòng)噪聲功率譜密度的新高噪聲模型 (new high noise model,NHNM)[15],設(shè)定參考鏡真實(shí)振動(dòng)速度信號(hào)vm(t)為這些信號(hào)的總和;3) 根據(jù)(6)式計(jì)算參考鏡振動(dòng)引入的真實(shí)相位噪聲Δφvib,再根據(jù)(4)式和(5)式計(jì)算其他相位噪聲Δφothers為0 時(shí)的躍遷概率Psim,該概率Psim即代表原子重力儀輸出信號(hào)P的仿真值,數(shù)據(jù)量N=30,如圖4(a)所示;4)設(shè)定真實(shí)延時(shí)系數(shù)τset=5 ms,真實(shí)增益系數(shù)Kset=0.9,設(shè)定地震計(jì)靈敏度Ks為將在實(shí)測(cè)中使用的CMG-3 ESP 型地震計(jì)的靈敏度,根據(jù)(9)式反推出此時(shí)地震計(jì)敏感質(zhì)量的振動(dòng)速度和地震計(jì)輸出電壓的仿真值vs(t)和Us(t),注意Us(t)即代表地震計(jì)在原子重力儀工作時(shí)同步輸出的信號(hào),如圖4(b)所示;5)對(duì)上述仿真出的概率Psim和地震計(jì)電壓信號(hào)Us(t)運(yùn)行振動(dòng)補(bǔ)償算法.

圖4 (a) 仿真干涉條紋;(b) 仿真振動(dòng)信號(hào)Fig.4.(a) The simulated fringe signal;(b) the simulated output of seismometer.

記以上仿真過(guò)程為第1 種情況,運(yùn)算結(jié)果如圖5 和圖6 所示.可以看出,不考慮其他相位噪聲時(shí),該算法可以精確地搜索出設(shè)定的延時(shí)系數(shù)τ=5 ms 和增益系數(shù)K=0.9(如圖中黑色虛線所示,下同),并將原始條紋Psim(ΔФth)調(diào)整為更接近理論條紋Pfit(ΔФth)的修正條紋Psim(ΔФv),其中典型的修正位置用橘色圓圈標(biāo)出.補(bǔ)償前對(duì)原始條紋Psim(ΔФth)進(jìn)行余弦擬合時(shí)的RMSE 值為13.1×10—3,而補(bǔ)償后對(duì)修正條紋Psim(ΔФv)進(jìn)行余弦擬合時(shí)的RMSE 值僅為0.02×10—3.同時(shí)該算法計(jì)算出的僅受振動(dòng)影響的條紋相對(duì)于理論條紋的殘差Rv(ΔФth)和原始條紋相對(duì)于理論條紋的殘差Rm(ΔФth)完全相同.依據(jù)2.3 節(jié)提出的評(píng)價(jià)指標(biāo),此時(shí)原始條紋殘差Rm的標(biāo)準(zhǔn)差σm為12.39×10—3,而去除振動(dòng)影響后的條紋殘差Rc的標(biāo)準(zhǔn)差σc僅為0.02×10—3,衰減幅度達(dá)到99.8%,驗(yàn)證了該振動(dòng)補(bǔ)償算法的有效性.

圖5 其他相位噪聲為0 時(shí)仿真得到的擬合RMSE 值與(a)延時(shí)系數(shù)和(b)增益系數(shù)的關(guān)系Fig.5.The RMSE value of fringe fitting as a function of (a) time delay and (b) proportional element,without considering other phase noises,calculated by simulation.

圖6 其他相位噪聲為0 時(shí)仿真得到的(a)修正前后的干涉條紋與理論擬合曲線(典型修正點(diǎn)用橘色圓圈標(biāo)出),(b)原始條紋的擬合殘差與計(jì)算出的僅由振動(dòng)導(dǎo)致的條紋的擬合殘差Fig.6.(a) The simulated fringe signals before and after vibration correction with the fitted curve as the reference (with the typical corrected data marked by orange circles),and (b) the fitting residual of the simulated fringe signals,without considering other phase noises.

進(jìn)一步,修改前述仿真流程的第3 步,設(shè)定其他相位噪聲Δφothers為±10 mrad 范圍內(nèi)的隨機(jī)值.這種設(shè)定可以更準(zhǔn)確地模擬實(shí)際測(cè)量,記該仿真過(guò)程為第2 種情況.此時(shí)雖然振動(dòng)相位Δφvib的范圍約為±400 mrad,是其他相位噪聲Δφothers的40 倍,但仿真運(yùn)算結(jié)果表明其他相位噪聲的干擾仍不可忽視.如圖7 和圖8 所示,此時(shí)算法搜索出的延時(shí)系數(shù)為τ=2.7 ms,與設(shè)定值的區(qū)別達(dá)到毫秒量級(jí);增益系數(shù)為K=0.957,與設(shè)定值相比偏大6.3%.此時(shí)雖然該算法同樣可以將原始條紋Psim(ΔФth)調(diào)整為更接近理論條紋Pfit(ΔФth)的修正條紋Psim(ΔФv)(典型位置如橘色圓圈所示),但偶爾會(huì)因設(shè)定為隨機(jī)值的其他相位噪聲Δφothers而產(chǎn)生隨機(jī)運(yùn)算誤差,如綠色方框所示的數(shù)據(jù)被修正到理論條紋曲線的另一側(cè)而不是曲線上;計(jì)算出的僅由振動(dòng)引入的條紋殘差Rv(ΔФth)和原始條紋殘差Rm(ΔФth)大致相同,但也存在一定差異.最后,此時(shí)原始條紋殘差Rm的標(biāo)準(zhǔn)差σm為17.32×10—3,而去除振動(dòng)影響后的條紋殘差Rc的標(biāo)準(zhǔn)差σc為7.82×10—3,衰減幅度為54.8%,明顯小于第1 種情況的衰減幅度.

圖7 考慮其他相位噪聲時(shí)仿真得到的擬合RMSE 值與(a)延時(shí)系數(shù)和(b)增益系數(shù)的關(guān)系Fig.7.The RMSE value of fringe fitting as a function of (a) time delay and (b) proportional element,with other phase noises considered,calculated by simulation.

圖8 考慮其他相位噪聲時(shí)仿真得到的(a)修正前后的干涉條紋與理論擬合曲線(典型修正點(diǎn)用橘色圓圈和綠色方框標(biāo)出),(b)原始?xì)埐钆c計(jì)算出的僅由振動(dòng)導(dǎo)致的條紋殘差Fig.8.(a) The simulated fringe signals before and after vibration correction with the fitted curve as the reference (with the typical corrected data marked by orange circles and green square),and (b) the fitting residual of the simulated fringe signals,with other phase noises considered.

上述結(jié)果說(shuō)明其他相位噪聲將明顯限制振動(dòng)補(bǔ)償算法的有效性.當(dāng)振動(dòng)環(huán)境變差,即振動(dòng)相位Δφvib在總相位噪聲中的權(quán)重增加時(shí),算法的有效性將顯著提升.例如,仍取其他相位噪聲Δφothers為±10 mrad 范圍內(nèi)的隨機(jī)值,將參考鏡振動(dòng)信號(hào)vm(t)中隨機(jī)振動(dòng)的幅值A(chǔ)rand增大為原來(lái)的10 倍,并記該仿真過(guò)程為第3 種情況.此時(shí)擬合得到的延時(shí)系數(shù)為τ=4.1 ms,增益系數(shù)為K=0.919,比第2 種情況更接近設(shè)定值;原始條紋殘差Rm的標(biāo)準(zhǔn)差σm為29.28×10—3,去除振動(dòng)影響后的條紋殘差Rc的標(biāo)準(zhǔn)差σc為6.78×10—3,衰減幅度為76.8%,高于第2 種情況.因此,理論上當(dāng)原子重力儀在復(fù)雜振動(dòng)環(huán)境中測(cè)量時(shí),該算法將發(fā)揮更為明顯的作用.

另一方面,增大數(shù)據(jù)量可以在小范圍內(nèi)提升振動(dòng)補(bǔ)償算法對(duì)延時(shí)系數(shù)和增益系數(shù)的計(jì)算精度,但對(duì)條紋修正效果的提升水平有限.例如,設(shè)定其他相位噪聲Δφothers和振動(dòng)相位Δφvib取值等同于第2 種情況,但將原子團(tuán)下落的時(shí)間間隔減小為1 s.在仿真總時(shí)長(zhǎng)仍為1 min 的情況下,此時(shí)躍遷概率信號(hào)的數(shù)據(jù)量擴(kuò)大1 倍,達(dá)到N=60,記該仿真過(guò)程為第4 種情況.此時(shí)擬合得到的延時(shí)系數(shù)為τ=6.8 ms,增益系數(shù)為K=0.876,比第2 種情況更接近設(shè)定值;但原始條紋殘差Rm的標(biāo)準(zhǔn)差σm為14.80×10—3,去除振動(dòng)影響后的條紋殘差Rc的標(biāo)準(zhǔn)差σc為6.56×10—3,衰減幅度為55.7%,和第2 種情況差別不大.

表2 匯總了以上所有運(yùn)算結(jié)果.可以看出,補(bǔ)償前后對(duì)條紋進(jìn)行余弦擬合時(shí)的RMSE 值的衰減幅度和條紋殘差標(biāo)準(zhǔn)差的衰減幅度及變化規(guī)律大致相同.如2.3 節(jié)所述,條紋殘差的變化更具有物理意義,因此在后續(xù)分析實(shí)測(cè)結(jié)果時(shí)不再專門提及RMSE 值的衰減幅度.

表2 仿真運(yùn)算結(jié)果Table 2.Simulated results.

4 實(shí)驗(yàn)結(jié)果

4.1 實(shí)驗(yàn)設(shè)備

該振動(dòng)補(bǔ)償算法的實(shí)際效果在中國(guó)計(jì)量科學(xué)研究院自主研制的NIM-AGRb-1 型原子重力儀上進(jìn)行了驗(yàn)證.NIM-AGRb-1 型重力儀具有兩種工作模式:其一為條紋掃描模式,在此期間拉曼激光波矢keff保持不變,掃頻速率α線性增加,可以得到原子干涉條紋;其二為重力測(cè)量模式,在此期間拉曼激光波矢keff不斷反向,即每一次作用于下落原子團(tuán)時(shí),激光脈沖的波矢keff均與上一次作用時(shí)數(shù)值相同、方向相反.此時(shí),結(jié)合在原子鐘上運(yùn)用的閉環(huán)鎖相技術(shù),可以實(shí)現(xiàn)對(duì)掃頻速率α的交叉鎖定,進(jìn)而保證每一次原子團(tuán)下落后都能得到重力測(cè)值g,而不必等待整個(gè)干涉條紋測(cè)量完畢后擬合出初相位進(jìn)而計(jì)算出重力值,或依據(jù)基于不同拉曼間隔T得到的干涉條紋的交叉點(diǎn)位置來(lái)計(jì)算重力值.配合清華大學(xué)研制的被動(dòng)式零長(zhǎng)彈簧結(jié)構(gòu)的隔振系統(tǒng),該重力儀可實(shí)現(xiàn)44 μGal·Hz—1/2的靈敏度,合成不確定度約為5.2 μGal[12].

驗(yàn)證該振動(dòng)補(bǔ)償算法時(shí),NIM-AGRb-1 型重力儀在第一種條紋掃描模式下工作,此時(shí)重力儀的輸出信號(hào)為實(shí)測(cè)躍遷概率P而非重力測(cè)值g,且最多掃描3 個(gè)整周期的干涉條紋.本實(shí)驗(yàn)的現(xiàn)場(chǎng)情況如圖9(a)所示,測(cè)試地點(diǎn)為中國(guó)計(jì)量科學(xué)研究院昌平院區(qū)的重力實(shí)驗(yàn)樓.實(shí)驗(yàn)過(guò)程中參考鏡直接放置在重力儀探測(cè)腔下方的地面上,旁邊放置英國(guó)Guralp 公司生產(chǎn)的CMG-3 ESP 型地震計(jì),其輸出信號(hào)的靈敏度為Ks=2000 V·m—1·s[31].探測(cè)腔中冷原子團(tuán)下落期間的激光脈沖的時(shí)間間隔為T=80 ms,激光波長(zhǎng)約為780 nm,相鄰原子團(tuán)下落的時(shí)間間隔為2 s.受重力儀自身軟件的控制,采集躍遷概率P的數(shù)據(jù)采集卡的采樣頻率為1 Hz,無(wú)法用于對(duì)地震計(jì)輸出電壓的采樣,因此使用第二臺(tái)美國(guó)National Instruments 公司生產(chǎn)的PXIe-6361 數(shù)據(jù)采集卡[32]來(lái)記錄地震計(jì)的輸出信號(hào)Us(t),采樣頻率為500 kS/s,如圖9(b)所示.該采集卡還同步采集了重力儀的某一路觸發(fā)信號(hào)Utrg(t).重力儀的控制軟件控制采集該信號(hào)總比原子團(tuán)下落時(shí)的第1 個(gè)π/2 脈沖提前40 ms,從而確保兩臺(tái)采集卡分別采樣的數(shù)據(jù)可以在時(shí)間上互相匹配.

圖9 實(shí)驗(yàn)現(xiàn)場(chǎng)圖 (a) NIM-AGRb-1 型原子重力儀;(b)地震計(jì)信號(hào)及觸發(fā)信號(hào)采集系統(tǒng)Fig.9.Photograph of (a) NIM-AGRb-1 atomic-interferometry absolute gravimeter,and (b) the data acquisition device recording the voltage signal of seismometer and trigger.

該測(cè)試地點(diǎn)位于地廣人稀的郊區(qū),周圍無(wú)交通、建筑和大型機(jī)械振源,且具有隔振地基,屬于非常安靜的測(cè)試環(huán)境,其地面振動(dòng)加速度的功率譜密度如圖10 所示.由于NIM-AGRb-1 型重力儀體積較大,尚無(wú)法搬運(yùn)至其他地點(diǎn),因此本次實(shí)驗(yàn)暫且僅在該安靜環(huán)境中進(jìn)行.之后該算法的振動(dòng)補(bǔ)償效果也將在振動(dòng)條件更為復(fù)雜的環(huán)境中驗(yàn)證.

圖10 測(cè)試地點(diǎn)的地面振動(dòng)加速度功率譜密度Fig.10.The power spectrum density of the ground vibration acceleration at the measurement site.

4.2 實(shí)驗(yàn)結(jié)果

本次實(shí)驗(yàn)首先采集了10 組數(shù)據(jù),每組數(shù)據(jù)包含30 次原子團(tuán)下落過(guò)程,其他控制參數(shù)如前所述.對(duì)每組數(shù)據(jù)分別運(yùn)行振動(dòng)補(bǔ)償算法,結(jié)果如表3 所示.

表3 10 組數(shù)據(jù)的實(shí)際測(cè)量結(jié)果Table 3.Measurement results obtained from 10 data sets.

可以看出該振動(dòng)補(bǔ)償算法在安靜環(huán)境中已具有良好的條紋修正效果,去除振動(dòng)影響后的條紋殘差Rc的標(biāo)準(zhǔn)差σc相比原始條紋殘差Rm的標(biāo)準(zhǔn)差σm的衰減幅度的均值達(dá)到44.4%,最大可達(dá)58.2%,與表1 中第2 種仿真情況的運(yùn)算結(jié)果基本一致;相比對(duì)原始條紋進(jìn)行余弦擬合時(shí)的RMSE值,擬合修正條紋時(shí)的RMSE 值也顯著減小.另一方面,各組數(shù)據(jù)分別計(jì)算出的最優(yōu)增益系數(shù)Kopt的標(biāo)準(zhǔn)差為0.115(均值為1.119),計(jì)算誤差較小;而最優(yōu)延時(shí)系數(shù)τopt的標(biāo)準(zhǔn)差為2.47 ms(均值為3.56 ms),計(jì)算誤差較大.前述表2 中第2 種情況的仿真結(jié)果已說(shuō)明,即使Δφothers比Δφvib小1 個(gè)量級(jí),τopt的計(jì)算誤差也將達(dá)到毫秒量級(jí).因此最有可能導(dǎo)致實(shí)測(cè)結(jié)果中τopt的標(biāo)準(zhǔn)差達(dá)到毫秒量級(jí)的原因仍是其他相位噪聲Δφothers.這再次說(shuō)明在后續(xù)工作中有必要在復(fù)雜振動(dòng)環(huán)境中驗(yàn)證該算法的振動(dòng)補(bǔ)償效果,即考察振動(dòng)相位Δφothers在總相位噪聲中所占的權(quán)重更大時(shí),σc相比σm的衰減幅度以及延時(shí)系數(shù)和增益系數(shù)的計(jì)算精度能否如仿真預(yù)測(cè)一樣提升.

可以以第1 組數(shù)據(jù)的計(jì)算結(jié)果為例觀察振動(dòng)補(bǔ)償算法的工作過(guò)程.如圖11 所示,該算法可以有效搜索出擬合條紋的RMSE 達(dá)到最小值時(shí)對(duì)應(yīng)的唯一最優(yōu)延時(shí)系數(shù)τopt和最優(yōu)增益系數(shù)Kopt.從圖12(a)可以看出,相比原始條紋P(ΔФth),對(duì)橫坐標(biāo)進(jìn)行振動(dòng)補(bǔ)償后的修正條紋P(ΔФv)更接近理論條紋Pfit(ΔФth)(典型位置如橘色圓圈所示).補(bǔ)償前對(duì)原始條紋P(ΔФth)進(jìn)行余弦擬合時(shí)的RMSE 值為11.2 × 10—3,而補(bǔ)償后對(duì)修正條紋Psim(ΔФv)進(jìn)行余弦擬合時(shí)的RMSE 值減小為4.8 × 10—3.原始條紋P相對(duì)于理論條紋的殘差Rm(ΔФth)和從地震計(jì)實(shí)測(cè)信號(hào)推算出的僅受振動(dòng)影響時(shí)的條紋Pv相對(duì)于理論條紋的殘差Rv(ΔФth)也基本吻合,如圖12(b)所示.

圖11 第1 組實(shí)測(cè)數(shù)據(jù)的擬合RMSE 值與(a)延時(shí)系數(shù)和(b)增益系數(shù)的關(guān)系Fig.11.The RMSE value of fringe fitting as a function of (a) time delay and (b) proportional element,calculated from the 1st dataset of measurement.

圖12 第1 組實(shí)測(cè)數(shù)據(jù)(a)修正前后的干涉條紋與理論擬合曲線(典型修正點(diǎn)用橘色圓圈標(biāo)出),(b)原始?xì)埐钆c計(jì)算出的僅由振動(dòng)導(dǎo)致的條紋殘差Fig.12.(a) The fringe signals before and after vibration correction with the fitted curve as the reference (with the typical corrected data marked by orange circles),and (b) the fitting residual of the fringe signals,calculated from the 1st dataset of measurement.

實(shí)際上還可以采用基于黃金分割原理的二維搜索方法重新編寫振動(dòng)補(bǔ)償算法流程,以實(shí)現(xiàn)對(duì)最優(yōu)延時(shí)系數(shù)τopt和最優(yōu)增益系數(shù)Kopt的全局搜索.計(jì)算結(jié)果表明,此時(shí)振動(dòng)補(bǔ)償算法的效果與表3 所示的分別搜索τopt和Kopt時(shí)的效果相當(dāng).這也再次證明(9)式所示的傳遞函數(shù)簡(jiǎn)化模型具有合理性.

基于表3 所列的10 組干涉條紋及后續(xù)測(cè)量得到的50 組干涉條紋進(jìn)行計(jì)算的結(jié)果如圖13 所示,其中圖13(a)為這60 組條紋各自修正前的殘差標(biāo)準(zhǔn)差σm以及修正后的殘差標(biāo)準(zhǔn)差σc,圖13(b)為由此得到的σc相對(duì)于σm的衰減比.可以看出每組條紋殘差標(biāo)準(zhǔn)差的衰減比的均值仍為45%左右,且分布特征具有隨機(jī)性,導(dǎo)致該隨機(jī)性分布最可能的原因仍是原始數(shù)據(jù)中其他相位噪聲造成的隨機(jī)干擾.

圖13 修正前后60 組實(shí)測(cè)干涉條紋的(a)殘差標(biāo)準(zhǔn)差和(b)殘差標(biāo)準(zhǔn)差的衰減比Fig.13.(a) The standard deviations of the fitting residuals of 60 sets of fringe signals with and without the vibration correction,and (b) the reducing factor between them.

如前所述,進(jìn)行本實(shí)驗(yàn)時(shí)NIM-AGRb-1 型重力儀工作在條紋掃描模式,重力儀本身不輸出依據(jù)干涉條紋得到的重力測(cè)值g.此時(shí),為了分析該振動(dòng)補(bǔ)償算法對(duì)重力測(cè)值的影響,可以參考浙江大學(xué)提出的方法[19]計(jì)算重力測(cè)值的離散度:1)對(duì)每一組條紋進(jìn)行余弦擬合時(shí),進(jìn)一步計(jì)算擬合曲線的初相位φinitial(n)(n為條紋序號(hào),n=1,2,···);2)去除φinitial(n)的均值,得到初相位的偏差序列Δφinitial(n);3)根據(jù)(4)式,可得重力偏差序列

也就是說(shuō),基于修正前后的干涉條紋P(ΔФth)和P(ΔФvib),可以依據(jù)(12)式得到振動(dòng)補(bǔ)償前后的重力偏差序列Δgm(n)和Δgc(n).由于本實(shí)驗(yàn)所在的測(cè)量環(huán)境無(wú)固定振源,僅有與測(cè)量過(guò)程在時(shí)間上不相關(guān)的地脈動(dòng)、輕微的人員走動(dòng)或?qū)嶒?yàn)樓外的車輛移動(dòng)等隨機(jī)振源,理論上參考鏡振動(dòng)只影響重力測(cè)值的離散度,不會(huì)產(chǎn)生恒定的系統(tǒng)誤差.因此補(bǔ)償前后重力偏差序列Δgm(n)和Δgc(n)的對(duì)比可以合理地說(shuō)明該振動(dòng)補(bǔ)償算法的有效性.

根據(jù)上述分析,由修正前后的條紋數(shù)據(jù)得到的修正前后的重力偏差如圖14 所示.對(duì)全部60 組數(shù)據(jù),振動(dòng)補(bǔ)償后重力測(cè)值的離散度相比補(bǔ)償前降低了一半左右;由修正前的實(shí)測(cè)條紋P(ΔФth)得到的重力偏差Δgm的標(biāo)準(zhǔn)差為σgm=51.5 μGal,而由修正后的條紋P(ΔФvib)得到的重力偏差Δgc的標(biāo)準(zhǔn)差σgc=22.0 μGal.σgc相對(duì)于σgm的衰減幅度達(dá)到57.3%,說(shuō)明通過(guò)該振動(dòng)補(bǔ)償算法提高干涉條紋的擬合精度可以顯著提升重力測(cè)量的精度.

圖14 由60 組實(shí)測(cè)干涉條紋得到的修正前后的重力偏差Fig.14.The variation of the measured g value deduced from 60 sets of fringe signals with and without the vibration correction.

目前NIM-AGRb-1 型原子重力儀作為計(jì)量基準(zhǔn)裝置已被設(shè)定為長(zhǎng)期工作在第二種重力測(cè)量模式,即掃頻速率α處于鎖定狀態(tài),無(wú)法直接采用本文提出的振動(dòng)補(bǔ)償算法進(jìn)行測(cè)量.因此本文仍以考察該算法對(duì)原子干涉條紋殘差的衰減效果為主.目前實(shí)驗(yàn)室正在研制采用條紋掃描模式來(lái)長(zhǎng)期工作的可移動(dòng)式原子重力儀,該振動(dòng)補(bǔ)償算法也將配合此重力儀在不同的測(cè)試環(huán)境中進(jìn)行長(zhǎng)期測(cè)量,驗(yàn)證其對(duì)包括阿倫方差、重力測(cè)量靈敏度在內(nèi)的關(guān)鍵重力參數(shù)的提升效果.

5 討論

在前述仿真和實(shí)驗(yàn)過(guò)程中,振動(dòng)補(bǔ)償算法查找最優(yōu)延時(shí)系數(shù)時(shí)的優(yōu)化目標(biāo)是擬合條紋時(shí)的RMSE 值,即取RMSE 達(dá)到其最小值時(shí)對(duì)應(yīng)的延時(shí)系數(shù)為最優(yōu)延時(shí)系數(shù)τopt.另一方面,根據(jù)理論分析及圖6(b)、圖8(b)和圖11(b)所示的仿真、實(shí)測(cè)結(jié)果,理論上振動(dòng)補(bǔ)償效果好時(shí),從地震計(jì)信號(hào)推算出的僅受振動(dòng)影響時(shí)的條紋相對(duì)于理論條紋的殘差Rv(ΔФth)非常接近原始條紋相對(duì)于理論條紋的殘差Rm(ΔФth).因此,參考基于相同簡(jiǎn)化模型的振動(dòng)補(bǔ)償方法在光學(xué)重力儀中的算法流程[27],實(shí)際上還可以嘗試以Rv和Rm的相關(guān)系數(shù)CR作為優(yōu)化目標(biāo)搜索最優(yōu)延時(shí)系數(shù)τopt.

相關(guān)系數(shù)表示2 個(gè)列向量之間線性相關(guān)的程度,取值范圍為—1 至1,絕對(duì)值越大表示相關(guān)性越強(qiáng)[33].因此,相關(guān)系數(shù)CR的絕對(duì)值越大,說(shuō)明Rv(ΔФth)和Rm(ΔФth)越接近,即推算出的振動(dòng)相位Δφvib越準(zhǔn)確,則修正后的條紋越接近擬合曲線,其RMSE值應(yīng)趨于最小值.需要說(shuō)明的是,計(jì)算相關(guān)系數(shù)CR時(shí)需要對(duì)Rc和Rv進(jìn)行歸一化處理,因此無(wú)法以其為優(yōu)化目標(biāo)來(lái)搜索最優(yōu)增益系數(shù)Kopt.

以CR達(dá)到最大值時(shí)對(duì)應(yīng)的延時(shí)系數(shù)為最優(yōu)延時(shí)系數(shù)τopt,相應(yīng)的仿真和實(shí)測(cè)數(shù)據(jù)的計(jì)算結(jié)果分別如表4 和表5 所示.總體而言,采用兩種優(yōu)化目標(biāo)時(shí)分別計(jì)算出的最優(yōu)延時(shí)系數(shù)τopt非常接近.如圖15 所示,RMSE 值隨延時(shí)系數(shù)的變化曲線與相關(guān)系數(shù)CR隨延時(shí)系數(shù)的變化曲線具有幾乎完全相同的變化規(guī)律,但變化方向相反,符合理論預(yù)期.個(gè)別情況下采用兩種優(yōu)化目標(biāo)時(shí)分別計(jì)算出的最優(yōu)延時(shí)系數(shù)τopt存在毫秒量級(jí)的差異,例如表4中第4 種情況和表4 中第2 組數(shù)據(jù),其中前者對(duì)延時(shí)系數(shù)的遍歷搜索結(jié)果如圖15(d)所示,RMSE 值達(dá)到最小值時(shí)的延時(shí)系數(shù)和相關(guān)系數(shù)CR達(dá)到最大值時(shí)的延時(shí)系數(shù)略有區(qū)別(分別用黑色和粉色點(diǎn)劃線標(biāo)出).導(dǎo)致這種情況最可能的原因仍是原始數(shù)據(jù)中其他相位噪聲造成的隨機(jī)干擾.

表4 以相關(guān)系數(shù)CR 為優(yōu)化目標(biāo)時(shí)的仿真運(yùn)算結(jié)果Table 4.Simulated results with the correlation value CR as the objective of the algorithm.

圖15 擬合RMSE 值和相關(guān)系數(shù)CR 與延時(shí)系數(shù)的關(guān)系 (a) 仿真時(shí)的第1 種情況;(b) 仿真時(shí)的第2 種情況;(c) 實(shí)測(cè)時(shí)的第1 組數(shù)據(jù);(d) 仿真時(shí)的第4 種情況Fig.15.The RMSE value of fringe fitting and the correlation value CR between the fitting residual of fringes before and after the vibration correction as a function of time delay of (a) 1st simulated dataset,(b) 2nd simulated dataset,(c) 1st measured dataset and(d) 4th simulated dataset.

表5 以相關(guān)系數(shù)CR 為優(yōu)化目標(biāo)時(shí)的10 組數(shù)據(jù)的實(shí)際測(cè)量結(jié)果Table 5.Measurement results obtained from 10 data sets with the correlation value CR as the objective of the algorithm.

對(duì)每一組數(shù)據(jù),依據(jù)振動(dòng)補(bǔ)償算法計(jì)算出最優(yōu)延時(shí)系數(shù)τopt和增益系數(shù)Kopt后,可以得到最終確定的僅受振動(dòng)影響時(shí)的條紋相對(duì)于理論條紋的殘差Rv,進(jìn)而得到最終確定的相關(guān)系數(shù)CR.圖16展示了60 組數(shù)據(jù)中每一組數(shù)據(jù)的相關(guān)系數(shù)CR以及修正后的條紋殘差標(biāo)準(zhǔn)差σc相對(duì)于修正前的條紋殘差標(biāo)準(zhǔn)差σm的衰減比的對(duì)應(yīng)情況.可以看出,相關(guān)系數(shù)CR越大,σc相對(duì)于σm的衰減比也越大,即振動(dòng)補(bǔ)償?shù)男Ч胶?二者的變化趨勢(shì)完全吻合,與理論預(yù)期相符,再次驗(yàn)證了該振動(dòng)補(bǔ)償算法原理的合理性.

圖16 最終確定的Rv 和Rm 的相關(guān)系數(shù)CR 與修正前后條紋殘差標(biāo)準(zhǔn)差的衰減比Fig.16.The correlation value CR between the fitting residual of fringes before and after the vibration correction,and the corresponding reducing factor of the standard deviation of the fitting residuals of fringe signals after the vibration correction.

還可以直接以2.3 節(jié)提出的條紋殘差標(biāo)準(zhǔn)差這一指標(biāo)作為振動(dòng)補(bǔ)償?shù)膬?yōu)化目標(biāo),即取(11)式所示的標(biāo)準(zhǔn)差σc達(dá)到最小值時(shí)的延時(shí)系數(shù)和增益系數(shù)為最優(yōu)值.相應(yīng)的仿真及實(shí)測(cè)數(shù)據(jù)的計(jì)算結(jié)果表明,在絕大多數(shù)情況下以RMSE 值或σc為優(yōu)化目標(biāo)時(shí)搜索出的τopt和Kopt均完全相同.個(gè)別情況下搜索出的延時(shí)系數(shù)存在毫秒量級(jí)的差異,增益系數(shù)也可能存在微小差異,如表6 所示.導(dǎo)致這種結(jié)果最可能的原因仍是原始數(shù)據(jù)中其他相位噪聲造成的隨機(jī)干擾.

表6 以殘差標(biāo)準(zhǔn)差為優(yōu)化目標(biāo)時(shí)的仿真及實(shí)測(cè)數(shù)據(jù)運(yùn)算結(jié)果Table 6.Simulated and measurement results with the standard deviation of the fitting residuals as the objective of the algorithm.

綜上,當(dāng)原子重力儀在振動(dòng)條件更復(fù)雜的環(huán)境中測(cè)量時(shí),即振動(dòng)相位Δφvib的影響遠(yuǎn)大于其他相位噪聲Δφothers的影響時(shí),該振動(dòng)補(bǔ)償算法的計(jì)算

精度和對(duì)條紋的修正效果有望得到進(jìn)一步提升.在后續(xù)工作中,該振動(dòng)補(bǔ)償算法將配合實(shí)驗(yàn)室正在研制的可移動(dòng)式原子重力儀在復(fù)雜環(huán)境中進(jìn)行測(cè)試.

總體而言,本文在第2 節(jié)定性分析了原子重力儀所用參考鏡與振動(dòng)傳感器之間的傳遞函數(shù)簡(jiǎn)化模型中的延時(shí)系數(shù)τ和增益系數(shù)K發(fā)生變化的原因,描述了對(duì)這兩項(xiàng)系數(shù)進(jìn)行搜索求解的具體振動(dòng)補(bǔ)償算法流程;在第3 節(jié)采用不同的仿真輸入對(duì)算法的有效性進(jìn)行了詳細(xì)的定量論證;在第4 節(jié)利用重力儀在非常安靜的環(huán)境中的實(shí)測(cè)數(shù)據(jù)驗(yàn)證了該算法對(duì)原子干涉條紋以及重力測(cè)值的補(bǔ)償效果.該算法最大的優(yōu)點(diǎn)是自適應(yīng)性較強(qiáng):相比現(xiàn)有的其他振動(dòng)補(bǔ)償研究成果,它具有較強(qiáng)的技術(shù)意義,可以直接應(yīng)用于任何能夠掃描出干涉條紋的原子重力儀在任意時(shí)間和地點(diǎn)的測(cè)量結(jié)果,以實(shí)現(xiàn)對(duì)傳遞函數(shù)中的兩項(xiàng)關(guān)鍵系數(shù)的實(shí)時(shí)標(biāo)定,實(shí)現(xiàn)當(dāng)前最好的振動(dòng)補(bǔ)償效果.另一方面,受限于NIM-AGRb-1型原子重力儀的體積、工作模式和測(cè)試任務(wù),該算法的有效性還未能在較為復(fù)雜的振動(dòng)環(huán)境中進(jìn)行驗(yàn)證;而目前也尚無(wú)現(xiàn)有的其他振動(dòng)補(bǔ)償方法在與本次實(shí)驗(yàn)相同的安靜環(huán)境中應(yīng)用于與該重力儀精度水平相當(dāng)?shù)脑又亓x上的實(shí)驗(yàn)數(shù)據(jù).從該算法與其他振動(dòng)補(bǔ)償方法在不同測(cè)量環(huán)境及具有不同噪聲水平的原子重力儀上的實(shí)驗(yàn)結(jié)果來(lái)看,該算法尚未展現(xiàn)出顯著優(yōu)于其他方法的補(bǔ)償效果,有待更進(jìn)一步的測(cè)試驗(yàn)證,其工作原理和具體流程也將在后續(xù)工作中不斷優(yōu)化.

6 結(jié)論

應(yīng)用于原子干涉式絕對(duì)重力儀的振動(dòng)補(bǔ)償方法的主要作用是修正參考鏡振動(dòng)引入的相位噪聲,從而提高干涉條紋的擬合精度,最終提高重力加速度g的測(cè)量精度.現(xiàn)有的振動(dòng)補(bǔ)償方法在確定的靜態(tài)或動(dòng)態(tài)環(huán)境中可以實(shí)現(xiàn)較高的測(cè)量精度,但大部分公開成果沒有介紹對(duì)于因測(cè)量環(huán)境改變而產(chǎn)生的原子重力儀中某些部件的機(jī)械傳遞函數(shù)的變化進(jìn)行實(shí)時(shí)修正的具體應(yīng)對(duì)流程.從長(zhǎng)期穩(wěn)定性的角度來(lái)看,這可能在一定程度上限制重力儀在不同環(huán)境尤其是復(fù)雜振動(dòng)環(huán)境的測(cè)量精度.本文介紹了一種應(yīng)用于原子重力儀的基于搜索傳遞函數(shù)簡(jiǎn)化模型系數(shù)的振動(dòng)補(bǔ)償方法,著重描述了其具體的算法流程,并利用仿真及實(shí)測(cè)結(jié)果對(duì)其合理性和有效性進(jìn)行了詳細(xì)論證.具體而言,該方法將參考鏡的振動(dòng)信號(hào)與測(cè)量該振動(dòng)的傳感器的實(shí)際輸出信號(hào)之間的真實(shí)傳遞函數(shù)簡(jiǎn)化為一個(gè)延時(shí)系數(shù)和一個(gè)增益系數(shù)的結(jié)合,通過(guò)計(jì)算軟件對(duì)重力儀的干涉條紋和傳感器的輸出信號(hào)進(jìn)行分析,搜索出當(dāng)前測(cè)量環(huán)境下這兩個(gè)系數(shù)的最優(yōu)估算值,由此可以更加精確地推算出當(dāng)前環(huán)境下參考鏡的真實(shí)振動(dòng)信號(hào).該方法具有自適應(yīng)性,可以在測(cè)量環(huán)境改變時(shí)對(duì)結(jié)構(gòu)的真實(shí)傳遞函數(shù)進(jìn)行合理估算,以實(shí)現(xiàn)在當(dāng)時(shí)當(dāng)?shù)貙?duì)原子干涉條紋最好的修正效果.實(shí)測(cè)結(jié)果表明,在極為安靜的環(huán)境中,本文提出的振動(dòng)補(bǔ)償算法可以將原子重力儀干涉條紋的余弦擬合殘差的標(biāo)準(zhǔn)差大幅衰減,最大衰減幅度達(dá)到58%.利用修正后的干涉條紋計(jì)算出的重力值的標(biāo)準(zhǔn)差相比修正前衰減57%.根據(jù)理論分析和仿真運(yùn)算,在復(fù)雜振動(dòng)環(huán)境中進(jìn)行重力測(cè)量時(shí),該算法將表現(xiàn)出更好的干涉條紋修正效果.因此,該振動(dòng)補(bǔ)償算法有望直接應(yīng)用在不同時(shí)刻和不同地點(diǎn)對(duì)原子重力儀測(cè)量結(jié)果的實(shí)時(shí)修正,提高原子重力儀在不同環(huán)境尤其是復(fù)雜振動(dòng)環(huán)境中的測(cè)量精度,其補(bǔ)償效果也將在之后的工作中進(jìn)行更充分的驗(yàn)證及優(yōu)化.