張偉,胡蕾,張釗搏

(1.中國地質(zhì)調(diào)查局 成都地質(zhì)調(diào)查中心，四川 成都 610082; 2.電子科技大學 資源與環(huán)境學院，四川 成都 611731; 3.中科院成都信息技術(shù)股份有限公司，四川 成都 610041)

0 引言

長周期大地電磁法目前是探測地球深部(>150 km)上地幔尺度電性結(jié)構(gòu)特征的主要方法，由烏克蘭Lviv Centre of Institute of Space Research研究所研制的LEMI-417儀器是最早引入到我國用于長時間、穩(wěn)定、連續(xù)觀測超低頻大地電磁場的一套成熟裝備。2007年起，中國科學院地質(zhì)與地球物理研究所[1]、成都理工大學[2]、中國地質(zhì)大學(北京)[3-4]以及中國地震局地質(zhì)研究所[5]等科研院所開始使用該儀器探測青藏高原及鄰區(qū)、龍門山構(gòu)造帶以及川滇黔相鄰區(qū)地震活動帶的深部地質(zhì)結(jié)構(gòu)，用于揭示其關(guān)鍵地質(zhì)帶的深部動力學特征。與該儀器配套的數(shù)據(jù)處理軟件是由俄羅斯Ivan M. Varantsov等編寫的PRC_MTMV商業(yè)軟件，用于實現(xiàn)將大地電磁場從時間域到頻率域的轉(zhuǎn)換，通過傅里葉變換計算得到功率譜，并通過功率譜估計得到大地電磁阻抗要素。

眾所周知，計算機對于二進制格式數(shù)據(jù)的讀取效率遠遠大于文本格式數(shù)據(jù)，因此在進行傅里葉變換時直接讀取二進制數(shù)據(jù)進行計算將大大提高程序的運行效率。但由于烏克蘭儀器生產(chǎn)廠家沒有公開給出其原始數(shù)據(jù)的二進制格式，因此在處理LEMI-417數(shù)據(jù)時，必須首先使用廠家提供的轉(zhuǎn)換程序?qū)x器二進制記錄數(shù)據(jù)轉(zhuǎn)換為文本格式，再使用PRC_MTMV軟件的LEMI2PRC子模塊將ASCII文本格式數(shù)據(jù)再次轉(zhuǎn)換為二進制數(shù)據(jù)后方能進行后續(xù)的傅里葉變換處理[3]，這無疑明顯地增加了一步不必要的處理操作。為此，作者對照V5-2000儀器的二進制數(shù)據(jù)格式的特點[6]，對LEMI-417儀器的原始二進制數(shù)據(jù)格式進行了深入研究，在不斷試錯后，成功剖析出了該儀器二進制記錄數(shù)據(jù)的數(shù)據(jù)結(jié)構(gòu)。

1 LEMI-417型系統(tǒng)的原始數(shù)據(jù)格式

LEMI-417是專門針對長時間采集超低頻段天然電磁場而研制的一款觀測系統(tǒng)，其主要特點是：①采用磁通門傳感器來采集磁場數(shù)據(jù)，觀測范圍為±70 000 nT，可以記錄到幾十萬秒的長周期信號；②電道采集盒子通過專門的電路設(shè)計能夠保證在長時間范圍內(nèi)保持很低的漂移，提供了4組電道接口，可以同時采集4組電極的電場數(shù)據(jù)；③采用高精度GPS同步授時能夠保證在長時間范圍內(nèi)具有很小的時移；④整體功耗很低，能夠保證在1塊12V45ah的鉛酸電池供電下連續(xù)工作一周以上。

該儀器以1 Hz的采樣率來連續(xù)觀測天然電磁場，采集的數(shù)據(jù)以二進制格式實時儲存到儀器內(nèi)部的CF卡中，文件名以L417_×××.B××的規(guī)則來進行命名，其中：L417為儀器的型號名稱；_×××為單點的數(shù)據(jù)存放文件序號，第一個數(shù)據(jù)文件編號為001，當GPS時間為格林威治時間的零點時，儀器系統(tǒng)自動新建一個新文件來存放數(shù)據(jù)，文件序號自動在前一個文件序號的基礎(chǔ)上累加1；B××的首字母B是二進制(binary)的縮寫，××是采集儀器的出廠編號。

1.1 頭文件格式

通常與測點數(shù)據(jù)采集相關(guān)的記錄信息(如儀器型號、采集日期、經(jīng)緯度、高程、電極距等)一般以文件頭的形式存放于二進制文件的首端。LEMI-417I型儀器共分配了32個字節(jié)的大小空間來存放信息(表1)，各字節(jié)位置信息闡述如下。

表1 LEMI-417原始二進制記錄數(shù)據(jù)的頭文件格式(tag)

1～4字節(jié)：4個字節(jié)長度，存放采集儀器的型號信息，對于LEMI-417I型儀器來說，第1個字節(jié)為ASCII字符‘4’(16進制為34)，第2字節(jié)為ASCII字符‘1’(16進制為31)，第3字節(jié)為ASCII字符‘7’(16進制為37)，第4字節(jié)為ASCII字符‘I’(16進制為49)。

5字節(jié)：1個字節(jié)長度，存放采集儀器的出廠編號。

6字節(jié)：1個字節(jié)長度，存放開始采集數(shù)據(jù)時的當前年信息。

7字節(jié)：1個字節(jié)長度，存放開始采集數(shù)據(jù)時的當前月信息。

8字節(jié)：1個字節(jié)長度，存放開始采集數(shù)據(jù)時的當前日信息。

9字節(jié)：1個字節(jié)長度，存放開始采集數(shù)據(jù)時的當前小時信息(24小時制)。

10字節(jié)：1個字節(jié)長度，存放開始采集數(shù)據(jù)時的當前分信息。

11字節(jié)：1個字節(jié)長度，存放開始采集數(shù)據(jù)時的當前秒信息。

12～15字節(jié)：4個字節(jié)長度，存放當前測點GPS同步的緯度。

16字節(jié)：1個字節(jié)長度，存放當前測點坐標是位于北緯還是南緯，以我國位于北半球為例，在我國測量數(shù)據(jù)時應(yīng)為ASCII字符‘N’( 16進制為4E)。

17～21字節(jié)：5個字節(jié)長度，存放當前測點GPS同步記錄的經(jīng)度信息。

22字節(jié)：1個字節(jié)長度，存放當前測點坐標是位于東經(jīng)還是西經(jīng)，以我國位于東半球為例，在我國測量數(shù)據(jù)時應(yīng)為ASCII字符‘E’( 16進制為45)。

23字節(jié)：1個字節(jié)長度，為空字符(對應(yīng)的16進制為00)。

24～25字節(jié)：2個字節(jié)長度，存放當前測點GPS同步的高程(單位：m)。

26字節(jié)：1個字節(jié)長度，存放儀器的采樣頻率。

27字節(jié)：1個字節(jié)長度，存放儀器內(nèi)部的VIN(外部供電電源電壓)信息。

28字節(jié)：1個字節(jié)長度，存放儀器內(nèi)部的VBAT(內(nèi)部紐扣電池電壓)信息。

29字節(jié)：1個字節(jié)長度，存放電道盒子上第1對電極的極距長度(單位：m)。

30字節(jié)：1個字節(jié)長度，存放電道盒子上第2對電極的極距長度(單位：m)。

31字節(jié)：1個字節(jié)長度，存放電道盒子上第3對電極的極距長度(單位：m)。

32字節(jié)：1個字節(jié)長度，存放電道盒子上第4對電極的極距長度(單位：m)。

1.2 時間序列格式

如圖1所示，在頭文件(Tag)后緊接著的是時間序列數(shù)據(jù)，LEMI-417以30個字節(jié)的大小空間來存放一個采樣時間片采集到的電磁場數(shù)據(jù)段(segment，表2所示)，由17個時間片組成一個記錄塊(record)，一個數(shù)據(jù)段結(jié)束后再接下一個文件頭數(shù)據(jù)，依次循環(huán)、直到采集結(jié)束。其中，單個時間片數(shù)據(jù)段中各字節(jié)位置信息闡述如下。

圖1 LEMI-417原始數(shù)據(jù)結(jié)構(gòu)分解示意Fig.1 Schematic diagram of original data structure decomposition

表2 LEMI-417原始二進制記錄數(shù)據(jù)的時間序列格式(segment)

1～3字節(jié)：3個字節(jié)長度，存放磁通門傳感器X方向的測量值(單位：nT)。

4～6字節(jié)：3個字節(jié)長度，存放磁通門傳感器Y方向的測量值(單位：nT)。

7～9字節(jié)：3個字節(jié)長度，存放磁通門傳感器Z方向的測量值(單位：nT)。

10～13字節(jié)：4個字節(jié)長度，存放第1道電極間的測量值(單位：mkV/m)。

14～17字節(jié)：4個字節(jié)長度，存放第2道電極間的測量值(單位：mkV/m)。

18～21字節(jié)：4個字節(jié)長度，存放第3道電極間的測量值(單位：mkV/m)。

22～25字節(jié)：4個字節(jié)長度，存放第4道電極間的測量值(單位：mkV/m)。

26～27字節(jié)：2個字節(jié)長度，存放磁通門傳感器溫度(單位：℃)。

28～29字節(jié)：2個字節(jié)長度，存放電道采集盒子的溫度(單位：℃)。

30字節(jié)：1個字節(jié)長度，存放時間片采樣間隔(單位：m)。

2 計算機數(shù)據(jù)類型與程序?qū)崿F(xiàn)

2.1 頭文件數(shù)據(jù)結(jié)構(gòu)與核心代碼

在C++語言中可利用struct類型變量來結(jié)構(gòu)化地讀取和存儲頭文件信息，再按表1所示各信息的字節(jié)位置和字節(jié)長度，定義其內(nèi)部的子變量數(shù)據(jù)類型，其核心代碼如下：

// 頭文件數(shù)據(jù)結(jié)構(gòu)定義

structTagData

{

unsignedcharstr_1; //4

unsignedcharstr_2; //1

unsignedcharstr_3; //7

unsignedcharstr_4; //儀器型號

unsignedcharno; //儀器編號

unsignedcharyear; //年

unsignedcharmonth; //月

unsignedcharday; //日

unsignedcharhour; //時

unsignedcharminute; //分

unsignedcharsecond; //秒

unsignedcharlatitude[4]; //緯度

charN_S; //北緯還是南緯

unsignedcharlongitude[5]; //經(jīng)度

charE_W; //東經(jīng)還是西經(jīng)

charseparator; //空字符

unsignedcharaltitude[2]; //高程

unsignedcharaverage; //采樣率

unsignedcharvin; // 外部電壓

unsignedcharvbat; // 內(nèi)部電壓

unsignedcharL1; //第一對極距長度

unsignedcharL2; //第二對極距長度

unsignedcharL3; //第三對極距長度

unsignedcharL4; //第四對極距長度

}

// … …由于篇幅有限省略部分非關(guān)鍵代碼

// 讀取儀器原始數(shù)據(jù)文件

TagDatatag;

FILE*fp=fopen(strFilePath,"rb");

chardataBuf[32] = "";

fread(dataBuf,sizeof(dataBuf),1,fp);

fclose(fp);

// 結(jié)構(gòu)體各子變量賦值

tag.str_1=dataBuf[0];

tag.str_2=dataBuf[1];

// … …由于篇幅有限省略部分代碼

tag.latitude[0]=dataBuf[11];

tag.latitude[1]=dataBuf[12];

tag.latitude[2]=dataBuf[13];

tag.latitude[3]=dataBuf[14];

// … …由于篇幅有限省略部分代碼

tag.L4=dataBuf[31];

2.2 時間序列數(shù)據(jù)結(jié)構(gòu)與核心代碼

同理，利用struct類型變量來結(jié)構(gòu)化地讀取和存儲時間片數(shù)據(jù)段信息，再按表2所示各信息的字節(jié)位置和字節(jié)長度，定義其內(nèi)部的子變量數(shù)據(jù)類型，其核心代碼如下：

// 單個時間片數(shù)據(jù)段結(jié)構(gòu)定義

structSegmentData

{

charHx[3]; //Hx

charHy[3]; //Hy

charHz[3]; //Hz

charE1[4]; //E1

charE2[4]; //E2

charE3[4]; //E3

charE4[4]; //E4

charTF[2]; //TF

charTE[2]; //TE

charNo; //

}

// … …由于篇幅有限省略部分非關(guān)鍵代碼

// 定義變量

intHx=0,Hy=0,Hz=0,E1=0,E2=0,E3=0,E4=0;shortintTf=0,Te=0;

charNo;

TagDatatag;

SegmentDatasegment;

// 讀取儀器原始數(shù)據(jù)文件

fseek(infile_origin,0L,SEEK_SET);

cur_rc= 0;

intNo= 0;

// 如果當前文件指針讀取位置不是結(jié)尾

while(cur_rc

{

if(i% 17 == 0)

{

fread(tag,sizeof(unsignedchar),32L,infile_origin);

cur_rc=ftell(infile_origin);

i++;

continue;

}

fread(segment,sizeof(unsignedchar),30L,infile_origin);

cur_rc=ftell(infile_origin);

i++;

j++;

// 將二進制數(shù)值轉(zhuǎn)換為整型變量，詳見下文

convertSegmentData(segment,&Hx,&Hy,&Hz,&E1,&E2,&E3,&E4,&Tf,&Te,&No);

pre_rc=cur_rc;

}

2.3 二進制數(shù)值到整型變量的轉(zhuǎn)換

LEMI-417原始數(shù)據(jù)是在儀器采集時存儲的二進制數(shù)據(jù)，將它賦值到計算程序語言的整型變量時必須進行二進制位操作轉(zhuǎn)換。以磁通門傳感器的x方向磁場數(shù)據(jù)(Hx)為例，由于整型變量(int)無論是在32位還是64位操作系統(tǒng)其字節(jié)長度均為32位，幸運地避開了不同操作系統(tǒng)程序代碼的特殊處理問題，位轉(zhuǎn)換的實質(zhì)就是將原始數(shù)據(jù)中各字節(jié)位置的數(shù)值完整地賦值到整型變量的各字節(jié)位置上，其核心代碼如下：

inttemp1=0,temp2=0,temp3=0;

temp1 =segment.Hx[2];

temp1=temp1 ? 16;

temp1 =temp1 & 0x00FF0000;

temp2 =segment.Hx[1];

temp2 =temp2 ? 8;

temp2 =temp2 & 0x0000FF00;

temp3 =segment.Hx[0];

temp3 =temp3 & 0x000000FF;

Hx=temp1|temp2|temp3;

if(Hx>=0)

{

Hx=Hx& 0x00FFFFFF;

}

else

{

Hx=Hx| 0xFF000000;

}

3 應(yīng)用開發(fā)實例

在對原始數(shù)據(jù)格式成功剖析的基礎(chǔ)上，筆者針對該LEMI-417儀器在實際工作中存在的一些問題，補充性地開發(fā)出一套基于MFC的、可以在任何Windows系統(tǒng)運行的界面程序，輔助一線野外工作人員對個別偶然會出現(xiàn)的數(shù)據(jù)問題進行人為糾正。

實際應(yīng)用中，LEMI-417型儀器在無GPS信號時會出現(xiàn)記錄紊亂的現(xiàn)象[4]。在峽谷地貌或林木茂密區(qū)開展工作時可能會出現(xiàn)個別時段GPS衛(wèi)星信號偶然解鎖的現(xiàn)象，使得儀器無法同步當前的GPS時鐘和測點坐標，在對應(yīng)的文件頭信息中記錄下多組差別迥異的信息。例如，野外工作中儀器部署在深切山谷地貌區(qū)時，在部分時段會出現(xiàn)儀器GPS無法鎖定到≥3顆衛(wèi)星的情況，此時儀器記錄的測點經(jīng)緯度和高程均變?yōu)?，個別情況甚至會出現(xiàn)1顆衛(wèi)星都無法鎖定的情況，此時便會出現(xiàn)衛(wèi)星失鎖問題，無法同步到正確的GPS時鐘，使采集數(shù)據(jù)時序出現(xiàn)混亂，從而導(dǎo)致在后期數(shù)據(jù)導(dǎo)出時無法順利形成連續(xù)的時間序列數(shù)據(jù)，需要進行2次返工，造成不必要的成本支出。

針對該常見問題，開發(fā)出了Lemi-417數(shù)據(jù)修改器(圖2)，可以將當前測點采集的連續(xù)多天的數(shù)據(jù)加載、合并到一起，遍歷每個文件頭中的GPS時鐘、測點經(jīng)度、緯度、高程信息，判斷每段數(shù)據(jù)對應(yīng)的GPS時鐘是否是以設(shè)置的采樣率進行連續(xù)觀測，會缺失的、不連續(xù)的時間序列進行三次樣條插值修復(fù)，此外判斷每段數(shù)據(jù)對應(yīng)的經(jīng)度、緯度和高程信息是否出現(xiàn)較大的數(shù)據(jù)偏差，通過取眾數(shù)方法，識別出正常狀態(tài)下和故障狀態(tài)下的記錄，對故障狀態(tài)下的記錄進行自動修復(fù)。最后通過最小二乘法對正常狀態(tài)下的GPS經(jīng)度、緯度和高程進行線性回歸擬合得到高精度的測點定位信息。由此，通過計算機程序解決了以往由于GPS偶然失鎖而需要重復(fù)返工的問題。

圖2 窗口程序界面截圖Fig.2 Screenshot of the window program interface

4 結(jié)論

隨著我國廣大學者對大地電磁時頻轉(zhuǎn)換算法研究的不斷深入[7-9]，國產(chǎn)化數(shù)據(jù)處理軟件[10]、儀器設(shè)備[11]的研制工作也已取得長足進步。本文深入剖析了LEMI-417I型地球深部探測主流觀測系統(tǒng)的原始數(shù)據(jù)格式，較詳細地給出了對應(yīng)的核心程序代碼，為廣大從事天然電磁場現(xiàn)代數(shù)字信號處理方面研究的學者提供一定的前期基礎(chǔ)支撐，對今后同類型國產(chǎn)儀器、軟件的研發(fā)工作具有一定的借鑒意義。此外，也針對LEMI-417I型儀器存在的GPS解鎖時的記錄紊亂問題，提出了一個經(jīng)濟有效的計算機輔助糾正技術(shù)方案。