基于gis的航磁實測切割線網(wǎng)調(diào)平系統(tǒng)的研究

基于gis的航磁實測切割線網(wǎng)調(diào)平系統(tǒng)的研究

0實測切割線網(wǎng)自動調(diào)平技術在航空磁體測量過程中，由于許多因素的影響，測量數(shù)據(jù)不可避免地存在一些誤差。導致這些誤差的因素主要有:(1)測量儀器。由航空磁力儀的動態(tài)噪聲引起的誤差。(2)導航定位。由GPS定位誤差引起的誤差。(3)飛行高度。由于測區(qū)地形起伏,各測線的航高不在同一海拔高度引起的區(qū)域磁場誤差。(4)飛機磁場補償。由磁補償平面八方位的剩余誤差引起的航向磁場誤差。(5)正常場改正。由地磁正常場改正引起的誤差。(6)磁日變校正。由磁日變觀測及校正引起的誤差。(7)磁場水平調(diào)整。在磁場調(diào)平中產(chǎn)生的背景磁場基點誤差。在上述誤差中,有隨機誤差,也有系統(tǒng)誤差。其中系統(tǒng)誤差主要表現(xiàn)為,各測線的背景場不處于同一水平或同一基點上,并具有較強的規(guī)律性。它們在剖面平面圖中,表現(xiàn)為不同飛行架次的測線間,磁場強度曲線總體高低不一致。在平面等值線圖中,表現(xiàn)為串珠、鋸齒、突變、波紋條帶等虛假異常。因此,這種磁場水平的不一致,對成圖質(zhì)量影響較大,給資料推斷解釋會帶來虛假信息。對于上述系統(tǒng)誤差,可利用實測切割線網(wǎng)對測線的磁場水平進行調(diào)整(即調(diào)平),使測區(qū)內(nèi)所有測線背景場的場值統(tǒng)一到同一水平或同一基點上,從而達到消除該系統(tǒng)誤差的目的。在實測切割線網(wǎng)調(diào)平中,必須將測線與切割線在它們交點(即切割點)處的磁場值作為調(diào)平參考值。因此,為了提高調(diào)平精度,要求切割點為平靜磁場點(即平靜切割點)。但在切割線網(wǎng)的實際測量中,切割線難免會穿過異常區(qū),很難保證所測切割點均為平靜點。因此,只有在數(shù)據(jù)后處理中,通過人～機聯(lián)作方式,對切割點是否為平靜點作出自動/人工判定,篩選出較高質(zhì)量的平靜點,以達到提高調(diào)平精度和成圖質(zhì)量的目的。隨著社會的不斷發(fā)展,航磁測量的地域也越來越廣大,觀測數(shù)據(jù)量也越來越巨大。因此,迫切需要一種自動高效的海量數(shù)據(jù)處理手段,來滿足當前社會發(fā)展的要求。鑒于此要求,我們設計開發(fā)了一套基于地理信息系統(tǒng)(GIS)的航磁實測切割線網(wǎng)自動調(diào)平軟件。作者在本文使用的航磁數(shù)據(jù),是中國地質(zhì)科學院地球物理地球化學勘查研究所使用的GB-4A型三光系高精度氦光泵磁力儀采集的數(shù)據(jù)。1原則和算法由于氦光泵磁力儀無零點漂移,所以我們采用算術平均迭代法進行測線數(shù)據(jù)調(diào)平。1.1各測線調(diào)平、調(diào)平先逐條求取各切割線與所有測線在平靜切割點處的磁場差值的算術平均值,并逐條對各切割線進行調(diào)平。在完成各切割線調(diào)平后,再逐條求取各測線與所有切割線在平靜切割點處的磁場差值的算術平均值,并逐條對各測線進行調(diào)平。重復迭代此過程,直到測區(qū)內(nèi)所有測線在滿足給定條件的切割點處的磁場均方差(參見1.2.5測量總精度計算)最小或變化很小為止。1.2該算法1.2.1切割線的優(yōu)化算法因航磁數(shù)據(jù)量巨大,所以在計算測線與切割線的交點時,必須采用優(yōu)化算法,否則將會占用很長的計算時間。在本軟件中,作者采用了如下優(yōu)化算法。(1)確定測量線和切割線是否相交取測線的起始和結(jié)束點形成一條線段;取切割線的起始和結(jié)束點形成另一條線段。判斷這二條線段是否相交。(2)各測點的距離如果方法(1)中判斷為測線與切割線相交,則進一步計算測線上各測點與切割線上各測點的距離的平方值(不是距離本身,因開方運算的計算時間開銷很大),并找出最短距離平方值、測線點號,以及切割線點號。如果此平方值小于給定值(如100.0m2,相當于10.0m距離),則認為此處為切割點。1.2.2安撫消毒點的選擇(1)靜態(tài)切割點的確定在△T正磁區(qū)域,沿測線方向的磁場梯度變化較平緩,且測線與切割線在相交處的飛行高度相差不大的切割點,可視為平靜切割點;在△T負磁區(qū)域,沿測線方向的磁場梯度變化較平緩,測線與切割線的切割點可視為平靜切割點。(2)切割點磁場方差在測線上,以切割點與其前后各n個測點(不包括切割點)磁場的方差,來衡量該切割點附近測點的磁場波動情況。如果該方差小于給定值,則說明該切割點處的磁場梯度變化較平緩,為平靜切割點。切割點磁場方差計算公式見式(1)。σi=√12ni+n∑j=i-n(Τi-Τi+j)2(1)σi=12n∑j=i?ni+n(Ti?Ti+j)2????????????????ue001?ue000ue000(1)式中σi為測線上第i個切割點處的磁場與該切割點前后各n個測點磁場的方差,單位為nT;n為測線上在第i個切割點前后所取的測點個數(shù);Ti+j為測線上第i+j個測點的磁場,單位為nT;Ti為測線上第i個切割點的磁場,單位為nT。1.2.3割斷線程值將一條切割線上各測點的磁場值,減去此切割線與各測線在平靜切割點處的磁場差值的算術平均值。按照此方法,逐步推廣到測區(qū)內(nèi)所有切割線。(1)第三條切割線磁場值的計算ΔˉΤci=1nn∑k=1(Τci-Τj)k(2)ΔTˉˉˉci=1n∑k=1n(Tci?Tj)k(2)式中ΔˉΤTˉˉˉci為第i條切割線與各測線在切割點處的磁場差值的算術平均值,即切割線磁場修正值,單位為nT;n為第i條切割線上參與計算的平靜切割點的個數(shù);Tci為第i條切割線與第j條測線在切割點處的切割線磁場值,單位為nT;Tj為第i條切割線與第j條測線在切割點處的測線磁場值,單位為nT。(2)第i條切割線上第1點的磁場值tc0ijΤcij=Τc0ij-ΔˉΤci(3)Tcij=Tc0ij?ΔTˉˉˉci(3)式中Tcij為第i條切割線上第j點修正后的磁場值,單位為nT;Tc0ij為第i條切割線上第j點修正前的磁場值,單位為nT;ΔˉΤTˉˉˉci為第i條切割線的磁場修正值,單位為nT。1.2.4連接線的值的更改將一條測線上各測點的磁場值,加上此測線與各切割線在平靜切割點處的磁場差值的算術平均值。按照此方法,逐步推廣到測區(qū)內(nèi)所有測線。(1)切割點磁場修正ΔˉΤj=1nn∑k=1(Τci-Τj)k(4)ΔTˉˉˉj=1n∑k=1n(Tci?Tj)k(4)式中ΔˉΤTˉˉˉj為第j條測線與各切割線,在切割點處的磁場差值的算術平均值,即測線磁場修正值,單位為nT;n為第j條測線上參與計算的平靜切割點的個數(shù);Tci為第j條測線與第i條切割線,在切割點處的切割線磁場值,單位為nT;Tj為第j條測線與第i條切割線,在切割點處的測線磁場值,單位為nT。(2)第i條第5項的磁場值Τij=Τ0ij+ΔˉΤi(5)Tij=T0ij+ΔTˉˉˉi(5)式中Tij為第i條測線上第j點修正后的磁場值,單位為nT;T0ij為第i條測線上第j點修正前的磁場值,單位為nT;ΔˉΤTˉˉˉi為第i條測線的磁場修正值,單位為nT。1.2.5切割點精度的計算按照我國航空磁測技術規(guī)范,采用公式(6)可計算出航磁測量總精度。σ=√12nn∑i=1δ2i(6)σ=12n∑i=1nδ2i???????√(6)式中σ為航磁測量總精度(即總誤差),在切割點處的測線磁場與切割線磁場的均方差,單位為nT;δi為第i個切割點處的測線磁場與切割線磁場的差值,單位為nT;n為參與計算的切割點個數(shù)。在計算σ時,允許舍掉磁場梯度較大(在設計的平面定位誤差范圍內(nèi),交叉點處的磁場變化大于3√2σ)32√σ),且磁場差值大于3√2σ的切割點,其余的切割點均應參與計算。2調(diào)整過程2.1數(shù)據(jù)預處理在調(diào)平前,需要對原始數(shù)據(jù)進行一些必要的預處理。(1)數(shù)據(jù)整合將測區(qū)內(nèi)不同架次的測線數(shù)據(jù)合并到一起;將測區(qū)內(nèi)不同架次的切割線數(shù)據(jù)合并到一起。(2)數(shù)據(jù)訪問對測區(qū)內(nèi)的所有測線和切割線進行檢查,刪除壞測線段和重復測線段。(3)日的變化按照同一日變基準值,對測區(qū)內(nèi)的所有測線和切割線進行磁日變改正。(4)正?，F(xiàn)場處理采用IGRF參數(shù),按月對測區(qū)內(nèi)的所有測線和切割線進行正常場改正。(5)投影坐標變換將原始數(shù)據(jù)中的WGS84球面,經(jīng)緯度坐標轉(zhuǎn)換為北京54平面直角坐標。(6)數(shù)據(jù)濾波對測區(qū)內(nèi)的所有測線和切割線的磁場曲線進行遞歸濾波,消除隨機干擾,使曲線保持相對平滑。2.2切割點的顯示在調(diào)平前,必須同時裝載測線和切割線數(shù)據(jù),它們將以GIS的圖層疊加方式顯示在圖形區(qū)內(nèi)(如圖1所示)。測線數(shù)據(jù)以顏色漸變剖面平面圖方式顯示;切割線數(shù)據(jù)以曲線剖面平面圖方式顯示;測線與切割線的切割點以圓點方式顯示。切割點的狀態(tài)分為二種:①平靜切割點;②非平靜切割點(即異常切割點)。平靜切割點在界面上顯示為灰色小圓圈;非平靜切割點在界面上顯示為白色小圓圈。2.3切割點的測量精度計算所有測線與所有切割線的切割點(包括平靜點和異常切割點),此計算過程需要的時間相對較長。在計算完成后,將每個切割點按測線號、切割線號、測線點號、切割線點號的格式保存到文件中,供以后的平靜切割點選定、切割線場值修正、測線場值修正,以及測量總精度計算功能使用。在測線和切割線的條數(shù)及線號、測點數(shù)不改變的情況下,可不再進行重復計算。此時,切割點的默認狀態(tài)為非平靜點,在圖形區(qū)以白色小圓點顯示(如圖1所示)。2.4切割點的選取在切割點計算完成后,可通過以下二種方式,在圖形區(qū)選定切割點的狀態(tài)。選定方法:(1)計算機自動選定。按照公式(1)對切割點的狀態(tài)進行自動判定。(2)人～機聯(lián)作選定。在計算機自動判定切割點的狀態(tài)后,可在圖形區(qū)對點的狀態(tài)進行人工調(diào)整。選定為平靜點:按下鼠標左鍵在圖形區(qū)拉方框,將要作為平靜點的切割點圈在方框內(nèi),松開鼠標左鍵,即可將方框內(nèi)的切割點選定為平靜點。選為平靜點的切割點,在界面變?yōu)榛疑A圈(如圖1所示)。選定為非平靜點:按下鼠標左鍵在圖形區(qū)拉方框,將要作為非平靜點的切割點圈在方框內(nèi),松開鼠標左鍵,即可將方框內(nèi)的切割點選定為非平靜點。選為非平靜點的切割點,在界面變?yōu)榘咨A圈(如圖1所示)。2.5自動提取磁場點在完成平靜點選定后,可按照公式(2)和式(3)進行切割線場值修正。程序?qū)⒆詣痈鶕?jù)平靜點的測線點號和切割線點號,提取測線磁場值和切割線磁場值,并計算出它們的磁場差值。在差值計算后,按切割線逐條求取磁場差值的算術平均值,并以此算術平均值來修正切割線上各測點的磁場值。最后,將修正后的磁場值保存到文件中。2.6切割點磁場差值的計算在完成切割線場值修正后,可按照公式(4)和公式(5)進行測線場值修正。程序?qū)⒆詣痈鶕?jù)平靜點的測線點號和切割線點號,提取測線磁場值和切割線修正后的磁場值,并計算它們的磁場差值。差值計算后,按測線逐條求取磁場差值的算術平均值,并以此算術平均值來修正測線上各測點的磁場值。最后,將修正后的磁場值保存到文件中。2.7測量總精度在完成測線場值修正后,可按照公式(6)進行測量總精度(磁場均方差)計算。切割線和測線場值(本節(jié)2.5和2.6)可進行迭代修正,直到測量總精度最小或變化很小為止。3航空異物綜合站航磁實測切割線網(wǎng)自動調(diào)平軟件在內(nèi)蒙古自治區(qū)大興安嶺中南段1∶50000航空物探綜合站勘查中得到了實際應用。以下是某測區(qū)的航磁數(shù)據(jù)調(diào)平處理結(jié)果。3.1t剖面位置圖2(見下頁)是調(diào)平前、后的航磁△T剖面平面圖。圖2(a)為調(diào)平前的航磁△T剖面平面圖,可以明顯看出,有少數(shù)測線的△T曲線與其周圍測線的△T曲線在總體磁場水平上不一致,顯得稀疏不均。圖2(b)為調(diào)平后的航磁△T剖面平面圖,通過調(diào)平處理,圖2(b)中所有測線的△T曲線的總體磁場水平基本一致,曲線分布較均勻。3.2近南方磁場分布不一致圖3(見后面)為調(diào)平前后的航磁△T等值線平面圖。圖3(a)為調(diào)平前的航磁△T等值線平面圖,可以明顯看出,在圖3(a)中的中部近南東向有串珠狀和鋸齒狀條帶異常存在,此條帶異常為少數(shù)測線的磁場水平不一致所引起,為虛假異常,它們使磁場分布變得雜亂,嚴重扭曲了真實異常。圖3(b)為調(diào)平后的航磁△T等值線平面圖,通過調(diào)平處理,消除了圖3(a)中的虛假條帶異常,顯示出更符合實際的場值分布特征。3.3迭代算法收斂性分析上述調(diào)平處理,共進行了七次迭代計算。該區(qū)的數(shù)據(jù)在計算過程中,總精度(σ)的變化情況如圖4(見下頁)所示。從圖4中可以看出,迭代算法的收斂性較好,在調(diào)平前的總精度為11.82nT,在七次迭代調(diào)平后的總精度為1.75nT,調(diào)平前后的數(shù)據(jù)精度提高了10.07nT。因此,數(shù)據(jù)質(zhì)量得到了明顯改善。由以上對比可以看到,通過該軟件對測線場值的調(diào)平處理,能有效地消除導致背景場水平不一致的系統(tǒng)誤差