內(nèi)容摘要:為解決流動重力觀測儀器數(shù)量多,周期長漂移率復(fù)雜等問題,Geoist團(tuán)隊研發(fā)了貝葉斯平差方法。該方法適用于不同儀器方差未知,絕對重力和多臺相對重力聯(lián)合觀測的數(shù)據(jù)處理問題。根據(jù)觀測時間長短可以分為線性漂移模型和非線性漂移模型兩種。今天我們談?wù)凣eoist的gravity模塊,如何完成具有大量冗余觀測條件下的重力平差任務(wù)。
1、重力平差
地震重力測量核心目標(biāo)是獲取高精度時變重力信號?,F(xiàn)有陸基重力觀測儀器以彈簧型相對型重力儀為主,由于測點(diǎn)之間地理跨度大,一般雙程閉合時間長,如果要獲得高精度重力測量結(jié)果,需要考慮儀器的非線性漂移特性等問題。
Geoist團(tuán)隊基于貝葉斯原理,提出了一種改進(jìn)的重力平差算法,假設(shè)儀器的漂移率光滑是已知的先驗信息,并根據(jù)貝葉斯統(tǒng)計學(xué)方法估計最優(yōu)化參數(shù),提出了一種專用于改進(jìn)現(xiàn)在廣泛應(yīng)用的線性漂移平差模型,并完成了可以用于最優(yōu)估計儀器非線性漂移特性的新平差算法設(shè)計。
重力平差在概念上與水準(zhǔn)平差類似,但是原理上不同,因為相對重力儀測量的段差具有獨(dú)立性,平差最重要的是對儀器漂移率進(jìn)行最優(yōu)估計。
對于水準(zhǔn)測量必須知道每個段差才能計算。舉個例子,假設(shè)水準(zhǔn)測量中間一個數(shù)據(jù)丟失了,后續(xù)平差無法開展。而重力不一樣中間去掉幾個測量數(shù)據(jù),還是可以形成新的段差組合,平差不受影響。
2、經(jīng)典平差
經(jīng)典平差是相對于測網(wǎng)內(nèi)固定的起算點(diǎn),根據(jù)往返閉合差進(jìn)行的點(diǎn)值估計工作。在經(jīng)典重力平差中,每臺儀器的漂移率一般假定為線性。根據(jù)實際測量點(diǎn)的順序和時間,可以構(gòu)成平差方程,采用最小二乘即可估計每個點(diǎn)值。
3、Geoist的平差功能
Geoist的gravity模塊實現(xiàn)了多種平差方法,包括經(jīng)典平差,貝葉斯平差等。graobj里面是平差計算中用到的類定義與接口實現(xiàn)。adjmethods是不同平差方法實現(xiàn)部分。
import geoist.gravity.graobj as gg
import geoist.gravity.adjmethods as adj
在gravity模塊中,平差過程需要儀器、測網(wǎng)、測量和工程四個主要類,分別表示為Meter類,Network類,Survey類和Campaign類。在平差前,需要生成每個類的實例,比如用到兩臺重力儀,就需要用Meter類生成兩個實例。Network是測網(wǎng),通常是我們所說的點(diǎn)之記,里面根據(jù)點(diǎn)的編號,存儲了每個測點(diǎn)的經(jīng)緯度、名稱等信息。Survey是一次測量,實例化后,可以將不同儀器、測網(wǎng)和觀測文件添加進(jìn)去。最后,通過實例化Campaign類,來調(diào)用平差方法,實例化后的Campaign可以添加Survey進(jìn)去,比如一期測量,由多個野外隊完成,先應(yīng)該生成不同的Survey,然后再統(tǒng)一的添加到Campaign中,再開始平差。
注意,要進(jìn)行平差,必須要生成Campaign的實例,然后正確地導(dǎo)入數(shù)據(jù)后,才可以開始計算。下面是一段示例代碼:
m1 = gg.Meter('CG-5','C098')
m1.msf = 1.00009
m2 = gg.Meter('CG-5','C099')
m2.msf = 1.000637
m3 = gg.Meter('CG-5','C097')
m3.msf = 1.000163
m4 = gg.Meter('CG-5','C098')
m4.msf = 1.00009
n1 = gg.Network('NorthChina',1)
n1.read_pnts('./data/SDQSX8.DZJ')
print(n1)
s1 = gg.Survey('DQSW+SX', '201508')
s1.add_meter(m1)
s1.add_meter(m2)
s1.add_meter(m3)
s1.add_meter(m4)
s1.net = n1
s1.read_survey_file('./data/QSCW201508p.098')
s1.read_survey_file('./data/QSCW201508p.099')
s1.read_survey_file('./data/SXCW971508p.ori')
s1.read_survey_file('./data/SXCW981508p.ori')
s1.corr_aux_effect()
print(s1)
#查找一個測量工程中某個測點(diǎn)號對應(yīng)的坐標(biāo)
slon, slat, selev = s1._get_pnt_loc('11003902')
#查找一個測量工程中某臺重力儀的格值
sf_val = s1._get_meter_sf('C097')
ag1 = gg.AGstation('白山洞絕對','11014121','A', 116.169, 40.018, 212.5)
ag1.ref_gra = 1110.54453
ag1.ref_gra_err = 5.0E-3
gravwork = gg.Campaign('IGP201604', 1)
gravwork.add_ag_sta(ag1) #添加絕對點(diǎn)信息 可以添加多次
gravwork.add_surveys(s1) #添加測量到平差任務(wù)
print(gravwork)
#開始平差pre_adj是完成從觀測文件重,生成平差矩陣的
gravwork.adj_method = 1 #1:cls ; 2:Baj; 3:Baj1
gravwork.pre_adj()
gravwork.run_adj('./data/grav_cls.txt',3)
一句話總結(jié):Geoist提供的平差方法對點(diǎn)之記,觀測文件都有一定的格式要求。后面我們會一點(diǎn)點(diǎn)介紹,通過7次課,將主要的平差程序使用方法介紹完,每次課的內(nèi)容也會不定期進(jìn)行修改和完善,力爭做到圖文并茂,通俗易懂。