楊金政， 邱崇濤， 田 宇， 高國林， 祁 程

(1.核工業(yè)航測遙感中心，石家莊 050002;2.中核集團公司 鈾資源地球物理勘查技術(shù)中心(重點實驗室)，石家莊 050002)

0 前言

航空地球物理勘探簡稱航空物探，具有效率高、成本低、便于大面積工作、探測深度較大等優(yōu)點，是基礎(chǔ)性和公益性地質(zhì)調(diào)查、戰(zhàn)略性礦產(chǎn)勘查的重要手段，是地質(zhì)勘查現(xiàn)代化的標志之一[1]。為了降低人員成本與風險，2011年，中國國土資源部試行了航空物探測量無人值守工作方式，2013年由中國核工業(yè)集團公司和中國地質(zhì)研究院又成功進行了無人機航空測量試驗。無人值守的工作方式與傳統(tǒng)方式的一個顯著差異是飛機從起飛到落地，儀器使用同一個測線編號持續(xù)記錄數(shù)據(jù)，首先需要對航跡數(shù)據(jù)進行測線識別、裁剪，目前，許多國內(nèi)、外著名的地球物理數(shù)據(jù)處理軟件(如Montaj Oasis、geoProbe等)，在此環(huán)節(jié)上，處理起來并不十分“舒暢”。筆者利用國內(nèi)地質(zhì)行業(yè)最常用的MapGis地理信息系統(tǒng)，對航測數(shù)據(jù)實現(xiàn)可視化的數(shù)據(jù)編輯與篩選，不僅適用于無人值守飛行數(shù)據(jù)，而且對于常規(guī)飛行數(shù)據(jù)，同樣具有針對性強、工作效率高的特點。

MapGis地理信息系統(tǒng)是武漢中地信息工程有限公司開發(fā)的GIS基礎(chǔ)平臺軟件系統(tǒng)，集地圖輸入、數(shù)據(jù)庫管理及空間分析的綜合系統(tǒng)；

不僅圖形編輯功能強大，而且具有較完善的二次開發(fā)功能[2]；是國內(nèi)地質(zhì)、物探領(lǐng)域中應(yīng)用最廣泛的軟件之一。屬性管理子系統(tǒng)是MapGis的重要組成部分，可實現(xiàn)對屬性數(shù)據(jù)管理[3]。本程序是利用MapGis圖元的屬性結(jié)構(gòu)(包括點、線和區(qū)屬性結(jié)構(gòu))，通過錄入、修改屬性數(shù)據(jù)，借助MapGis強大的二次開發(fā)和圖形編輯功能，交互式地進行數(shù)據(jù)的識別、篩選，最終獲得符合質(zhì)量要求的測量數(shù)據(jù)。

1 設(shè)計背景

1.1 航空物探數(shù)據(jù)簡介

我中心承擔的航測任務(wù)中，大多以航磁、航放同步測量為主。航空數(shù)據(jù)采集是以飛行架次為單位。數(shù)據(jù)記錄為一個二進制文件，先將二進制文件轉(zhuǎn)換為3個ASCII文件(航磁、航放窗數(shù)據(jù)、航放全譜)；每個測量數(shù)據(jù)均包括文件頭和測量數(shù)據(jù)兩個部分；文件頭主要記錄儀器型號、參數(shù)設(shè)置、測量時間、地點、架次等參數(shù)。測量數(shù)據(jù)記錄了該架次所有實測數(shù)據(jù)，包括早晚測試(能譜)、早、晚基線(磁、放)、測線(磁、放)等測量數(shù)據(jù)。圖1為無人機航放窗數(shù)據(jù)ASCII文件的數(shù)據(jù)格式。從圖1中可看出，測量數(shù)據(jù)包括測線號、基點號、經(jīng)緯度、坐標值、高度值、溫度值、放射性元素測量值、日期和時間等內(nèi)容。

1.2 數(shù)據(jù)質(zhì)量控制

航測數(shù)據(jù)的采集受諸多因素的影響，在確保天氣、儀器性能、地面磁日變站等正常前提下，通過測區(qū)范圍、飛行速度、離地飛行高度、偏航距等指標衡量數(shù)據(jù)質(zhì)量；在物探項目中，項目設(shè)計測線是野外數(shù)據(jù)采集最基本的依據(jù)。離地飛行高度(離地高度)，在航空伽瑪能譜測量中是一項十分重要的參數(shù)。一般情況下，飛機應(yīng)沿地形起伏飛行，高度保持在80 m～ 120 m。由于放射性計數(shù)率隨高度呈指數(shù)衰減，離地高度過高，導致計數(shù)率偏低，將會使地面弱異常被遺漏，同時增大統(tǒng)計漲落誤差；離地高度過低，無法將測線間區(qū)域全覆蓋[4]，同時飛行安全指數(shù)大大降低。偏航距指實際飛行航跡與設(shè)計測線的偏移距離，依據(jù)測量比例尺確定最大偏航距[5]。在實際飛行測量中，一條測線結(jié)束后，飛機轉(zhuǎn)入下一條測線時，測量數(shù)據(jù)雖處于項目設(shè)計之外，但在實際數(shù)據(jù)處理時，應(yīng)盡量多地保留測區(qū)外符合質(zhì)量要求的數(shù)據(jù)。在測區(qū)中，受高大建筑物、地形等因素影響，往往也會出現(xiàn)“超高”、“偏離”的現(xiàn)象，此時必須及時了解具體情況，做出下一步飛行方案。

2 程序?qū)崿F(xiàn)

為了盡可能地保留測區(qū)外的有效數(shù)據(jù)，同時及時準確掌握各架次的飛行質(zhì)量，可視化的數(shù)據(jù)篩選顯得十分必要。

2.1 設(shè)計流程

總體上，利用此程序進行數(shù)據(jù)編輯遵循“數(shù)據(jù)—圖形—數(shù)據(jù)”變換過程，分為5個步驟，除第4步為手工操作外，其余均由程序自動完成(圖2)。

1)設(shè)計測線繪制。根據(jù)設(shè)計坐標，生成具有屬性結(jié)構(gòu)的點、線、區(qū)三個文件，將各自設(shè)計測線信息錄入其中，作為自動識別依據(jù)。

圖2 程序設(shè)計流程圖Fig.2 The flow chart of program

圖3 無人機某一架次物探作業(yè)航跡圖Fig.3 Map of a sortie of flight track of airborne geophysical survey using drone

2)航跡線的繪制。將實測數(shù)據(jù)繪制成圖，線文件為航跡線，包含“測線編號”的屬性。點文件包括：①子圖點文件，表明“測點位置”和“高度等級”；②注釋點文件，包括“測線編號”、“飛行高度”信息(圖3)。

3)依據(jù)設(shè)計測線，將航跡線自動分離、識別并補充錄入“偏航距”信息。

4)在MapGis中，依據(jù)測點信息(點文件)，手工編輯航跡線，剔除不合格數(shù)據(jù)。

5)依據(jù)編輯后的航跡線端點坐標，篩選實測數(shù)據(jù)。

6)換名存儲航跡線，作為下一架次數(shù)據(jù)編輯的參照。

2.2 程序?qū)崿F(xiàn)

此程序基于MapGis軟件提供的COM組件[6]，使用C#語言編寫。有關(guān)MapGis的二次開發(fā)環(huán)境、組件注冊與引用以及點、線、區(qū)的寫入等技術(shù)[7-8]。筆者以2015年在新疆喀什地區(qū)開展的無人機航空物探綜合站測量技術(shù)研發(fā)與應(yīng)用示范項目[9-10]中的航測數(shù)據(jù)為依據(jù)，得到簡述程序?qū)崿F(xiàn)過程。

2. 2.1 設(shè)計測線成圖

1)準備工作。在程序運行前，建立帶有屬性結(jié)構(gòu)的點、線、區(qū)模板文件。由于該程序設(shè)計僅針對于飛行質(zhì)量，并不對測量值進行修正、處理，所以線、區(qū)文件只添加“測線編號”屬性；點文件添加“基點號”和“測線編號”兩項屬性(圖4)，其他數(shù)值未寫入屬性結(jié)構(gòu)中。上述三個模板文件均為空文件，隨后生成的MapGis文件均基于此模板完成[3]。設(shè)置方法可參見參考文獻[3]。

為提高運算效率，坐標位置也未寫入屬性結(jié)構(gòu)中，而是通過圖中繪制的測線、測點的坐標值換算獲得。

2)設(shè)計測線的繪制。根據(jù)基線、設(shè)計測線的端點坐標，在MapGis中繪制設(shè)計測線圖(圖5)，包括測線編號(點文件)、設(shè)計測線(線文件)、測線識別范圍(區(qū)文件)。測線識別范圍是指單條設(shè)計測線兩側(cè)外擴接近1/2線距，兩端外擴1 km～ 2 km的區(qū)域。由于線、區(qū)文件均是基于上述模板文件生成，所以它們中的圖元均具有“測線編號”的屬性結(jié)構(gòu)。線文件，屬性錄入代碼為：

myLinArea.Load(modelWL);//加載線模板文件

Record att2 = new Record(); //屬性定義并實例化

myLinArea.att.Get(LineCounts, out att2); //獲取原有屬性

att2.Value[2] = LineNo; //獲取屬性(線號)

myLinArea.att.Write(LineCounts, att2);//寫入屬性

……

//另存線文件

myLinArea.Save(WorkArea23 + @"設(shè)計測線.wl");

//區(qū)屬性定義與線屬性類似，(略)

另外，基于測區(qū)范圍外1 km～2 km，小于飛機轉(zhuǎn)彎區(qū)繪制一個范圍框，稱之為“航跡分離框”，用于實測測線的分離(圖4)。

圖4 點、線、區(qū)文件屬性結(jié)構(gòu)參數(shù)Fig.4 Parameters of attribute structure of point, line and region files(a)點文件屬性結(jié)構(gòu);(b)線文件屬性結(jié)構(gòu);(c)區(qū)文件屬性結(jié)

圖5 在MapGis中航跡編輯界面圖Fig.5 The edit interface of flight track on MapGis

2.2.2 實測數(shù)據(jù)成圖

本次測量中，航磁測量采樣率為10 Hz，航放測量為1 Hz[9]。為了提高運算效率，優(yōu)先選用數(shù)據(jù)量較小的能譜數(shù)據(jù)，利用線號、基點號、坐標位置(x，y)、離地高度等數(shù)據(jù)列，繪制具有屬性結(jié)構(gòu)的測點(子圖)、注釋點、航跡(線)(圖3、圖5)。依據(jù)離地高度等級，將測點和注釋點賦予不同的顏色，直觀顯示數(shù)據(jù)質(zhì)量。測點及注釋點的繪制代碼如下：

private void writingCedian(PntArea PntAreaCedian, PntArea PntAreaCDnote, DotInfo5 myDot2, Pnt_Info PntInfoNor, Pnt_Info PntInfoNote) //測點、注釋點子函數(shù)

{

double fontSize = Convert.ToDouble(textBox11.Text);

Record att3 = new Record();//屬性定義并實例化

D_Dot pntXY3cd = new D_Dot();

pntXY3cd.x = myDot2.x; //測點X坐標

pntXY3cd.y = myDot2.y; //測點Y坐標

D_Dot pntXY3cdAlt = new D_Dot();

pntXY3cdAlt.x = myDot2.x+1; //注釋點X坐標

pntXY3cdAlt.y = myDot2.y - fontSize*0.75; //注釋點X坐標

PntAreaCedian.Append(pntXY3cd, "a", PntInfoNor);

int DotCounts = PntAreaCedian.count - 1;

PntAreaCedian.att.Get(DotCounts, out att3); //取點屬性

att3.Value[1] = myDot2.fn;//基點號賦值

att3.Value[2] = myDot2.lineNo;//測線號賦值

att3.Value[3] = myDot2.ralt;//雷達高度號賦值

PntAreaCedian.att.Write(LineCounts, att3);//寫入屬性值

//寫入測點

PntAreaCDnote.Append(pntXY3cd,myDot2.fn,PntInfoNote);

//寫入注釋點,高度值,偏航距為零,測線識別后,再寫入

PntAreaCDnote.Append(pntXY3cdAlt,myDot2.ralt.ToString() + "-" + myDot2.dis, PntInfoNoteNor);

}

程序中myDot2、PntInfoNor和PntInfoNote分別為測點、注釋點的參數(shù)值，以不同的點參數(shù)值標示飛行質(zhì)量等級。

2.2.3 測線識別

首先，利用航跡分離框?qū)⒈炯艽魏桔E線分離為獨立的航跡線段(圖5)。代碼如下：

double m_Radiu = 0.0001;//模糊半徑

MAPGISBASCOM1Lib.IAnalysis pClip = new MAPGISBASCOM1Lib.Analysis();

//myFrame-分離框rawLinArea-原始線文件desLinArea-結(jié)果線文件

pClip.ClipLin(myFrame, rawLinArea, desLinArea, m_Radiu, Enum_Clip_Type.gisOVLY_INCLIP);

然后，遍歷各航跡線段，利用測線識別范圍(區(qū)文件)，逐一賦予“測線編號”屬性。傳遞屬性關(guān)鍵代碼如下：

for (int intK = 1; intK < myLinArea.count; intK++)

{ //遍歷所有實測測線

//獲取實測測線上的點集合

int i = myLinArea.Get(intK, out dataSet, out inf, out dima);

…….(在測線尋找3個合適的點坐標，代碼較長，略)

myLinArea.att.Get(intK, out att2); //獲取實測測線屬性

for (int intN = 1; intN < myRegArea.count; intN++)

{ //遍歷所有設(shè)計測線范圍

//計算實測測點與設(shè)計測線范圍的最大距離

double minDis1 = myRegArea.MaxDistOfPntToReg(dot1,intN); double minDis2 =myRegArea.MaxDistOfPntToReg(dot2,intN);

if (minDis1 == 0 && minDis2 == 0)

{ //最大距離為0時，表明測線在設(shè)計測線范圍內(nèi)

//獲取設(shè)計測線范圍的屬性，即線號

myRegArea.RegAtt.Get(intN, out myRecord);

string LineNo = myRecord.Value[3]; //賦值

att2.Value[2] = LineNo; // 傳遞給線屬性(線號)

myLinArea.att.Write(intK, att2); // 錄入線屬性

break;

}

}

}

最后，根據(jù)設(shè)計測線位置，計算每個航跡點的偏航距并繪制更新偏航距的標示(圖5)。

2.2.4 編輯航跡

借助MapGis軟件平臺，使用不同的顏色、符號標示出每個測點數(shù)據(jù)質(zhì)量等級，將實測數(shù)據(jù)以圖形展現(xiàn)，可以很清晰了解本架次的數(shù)據(jù)質(zhì)量。參照以往架次的航跡線和測點標識，利用MapGis圖形編輯功能，逐一修剪本架次航跡線，確定需保留數(shù)據(jù)范圍。

2.2.5 數(shù)據(jù)篩選與輸出

數(shù)據(jù)編輯基于MapGis中確定的航跡線端點，從源數(shù)據(jù)文件中提取質(zhì)量合格的數(shù)據(jù)。

讀取每條航跡線的端點坐標和航跡點(點文件)屬性結(jié)構(gòu)，從航跡點中獲取“測線編號”和“基點號”。以“基點號”為依據(jù)，比對、篩選航磁、航放源數(shù)據(jù)，更換“測線編號”，獲得最終數(shù)據(jù)。從航跡線中提取信息的代碼如下：

for (int intL = 1; intL < myLinArea3.count; intL++)

{

//獲取圖元位置

myLinArea3.Get(intL, out lineSet, out linA, out dim);

myLinArea3.att.Get(intL, out attLine); //獲取圖元屬性

for (int intI = 0; intI < lineSet.count; intI++)

{ //遍歷線圖元節(jié)點

tempXY.x = Convert.ToDouble (lineSet[intI].x.ToString ("0.00"));//提取X坐標

tempXY.y = Convert.ToDouble (lineSet[intI].y.ToString ("0.00"));//提取Y坐標

tempXY.lineNo =attLine.Value[2]; //提取線號

lineLinInfo.Add(tempXY); //寫入泛型集合

}

}

從航跡點中提取信息的代碼如下：

for (int intK = 1; intK < myPntArea3.count; intK++)

{ //遍歷點圖元

myPntArea3.att.Get(intK, out myAtt);//獲取圖元屬性

myPntArea3.GetPos(intK, out myXY); //獲取圖元位置

//提取X坐標

myDot1.x = Convert.ToDouble(myXY.x.ToString("0.00"));

//提取Y坐標

myDot1.y = Convert.ToDouble(myXY.y.ToString("0.00"));

myDot1.lineNo =myAtt.Value[2]; //提取線號

myDot1.Fn =myAtt.Value[1]; //提取基點號

//提取離地高度值

myDot1.ralt = Convert.ToDouble(myAtt.Value[3]);

myDot2.Add(myDot1); //寫入泛型集合

}

圖6為數(shù)據(jù)篩選后的結(jié)果數(shù)據(jù)。從圖6中可以看出，數(shù)據(jù)格式未變，測線編號已將原來的“005”變?yōu)榱藢嶋H測線號，并刪除了不符合質(zhì)量要求的數(shù)據(jù)。

另外，該程序還進行了本架次的數(shù)據(jù)質(zhì)量、工作量等信息計算、統(tǒng)計。

將本架次的航跡線更換不同顏色，以飛行日期和架次命名保存此次航跡線，供下一架次的數(shù)據(jù)裁截做參照對比。

2.3 有人值守航測數(shù)據(jù)的處理

該程序同樣適用于有人值守航測數(shù)據(jù)篩選。在正常情況下，由于每個架次的測線編號已事先錄入，可直接手工編輯航跡線，再根據(jù)線號、基點號比對，篩選實測數(shù)據(jù)。如若在儀器發(fā)生故障、操作員記錄錯誤等情況下，容易圈定錯誤范圍，利用航跡線屬性修改、圖形編輯篩選數(shù)據(jù)，可明顯提高工作效率。

3 結(jié)語

MapGis軟件是國內(nèi)地質(zhì)行業(yè)使用范圍最廣的軟件之一，基于其屬性結(jié)構(gòu)進行二次開發(fā)，具有開發(fā)效率高，操作簡單快捷等優(yōu)點。利用程序提供的圖形界面，交互式地編輯圖形，獲取合格的實測數(shù)據(jù)，不僅提高了數(shù)據(jù)裁切的針對性，而且可以熟知飛行數(shù)據(jù)質(zhì)量和工作進度。同時，在數(shù)據(jù)截取后，對飛行數(shù)據(jù)質(zhì)量進行統(tǒng)計，避免了重復(fù)工作，提高了工效。此程序已應(yīng)用于新疆喀什地區(qū)無人機試驗項目(2015年)、黑龍江省完達山—太平嶺地區(qū)航測項目(2016年)和黑龍江省龍江—嫩江地區(qū)航測項目(2016年)中，經(jīng)不斷完善，實用效果較好。