基于MapGis平臺可視化篩選無人值守航空物探數據

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

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

0 前言

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

MapGis地理信息系統是武漢中地信息工程有限公司開發的GIS基礎平臺軟件系統，集地圖輸入、數據庫管理及空間分析的綜合系統；

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

1 設計背景

1.1 航空物探數據簡介

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

1.2 數據質量控制

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

2 程序實現

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

2.1 設計流程

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

1)設計測線繪制。根據設計坐標，生成具有屬性結構的點、線、區三個文件，將各自設計測線信息錄入其中，作為自動識別依據。

圖2 程序設計流程圖Fig.2 The flow chart of program

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

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

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

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

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

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

2.2 程序實現

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

2. 2.1 設計測線成圖

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

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

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

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 + @"設計測線.wl");

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

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

圖4 點、線、區文件屬性結構參數Fig.4 Parameters of attribute structure of point, line and region files(a)點文件屬性結構;(b)線文件屬性結構;(c)區文件屬性結

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

2.2.2 實測數據成圖

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

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

{

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分別為測點、注釋點的參數值，以不同的點參數值標示飛行質量等級。

2.2.3 測線識別

首先，利用航跡分離框將本架次航跡線分離為獨立的航跡線段(圖5)。代碼如下：

double m_Radiu = 0.0001;//模糊半徑

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

//myFrame-分離框rawLinArea-原始線文件desLinArea-結果線文件

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

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

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++)

{ //遍歷所有設計測線范圍

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

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

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

{ //最大距離為0時，表明測線在設計測線范圍內

//獲取設計測線范圍的屬性，即線號

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

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

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

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

break;

}

}

}

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

2.2.4 編輯航跡

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

2.2.5 數據篩選與輸出

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

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

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++)

{ //遍歷線圖元節點

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為數據篩選后的結果數據。從圖6中可以看出，數據格式未變，測線編號已將原來的“005”變為了實際測線號，并刪除了不符合質量要求的數據。

另外，該程序還進行了本架次的數據質量、工作量等信息計算、統計。

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

2.3 有人值守航測數據的處理

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

3 結語

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