EH4數據處理在SCS2D二維反演平臺中的實現

李美英+屈栓柱

摘 要 Matlab高級語言開發平臺的函數資源豐富，對矩陣的計算簡便快捷，文件輸出方便。以這些特點為契機，在Matlab平臺上，轉換EH4儀器采集的原始Z文件為SCS2D二維反演所需的MTD文件，可大大提高EH4數據的處理能力。

關鍵詞 音頻大地電磁法；Matlab平臺；SCS2D反演平臺；代碼

中圖分類號：TP3 文獻標識碼：A 文章編號：1671-7597（2014）05-0041-03

可控源音頻大地電磁法（CSAMT）是在大地電磁法（MT）和音頻大地電磁法（AMT）的基礎上發展起來的一種人工源頻率測深方法。大地電磁法（MT）的觀測頻率很低（n×10-5～n×104 Hz），所以勘探深度很大，可達100千米以上。但由于頻率偏低，所以對淺層的分辨率較差，而且生產效率也比較低。為了更好地研究幾十米到幾千米深度范圍內的地電結構，在MT方法的基礎上，形成了音頻大地電磁法（AMT）。音頻大地電磁法的工作方法、觀測參數與MT方法相同，不過該方法主要觀測由于雷電作用產生的音頻（n×100～n×104 Hz）大地電磁場，由于音頻大地電磁法的觀測頻段內，天然大地電磁場的強度較弱，人文干擾大，所以信噪比很低，野外數據采集十分困難，需要長時間采集并采用多次疊加技術。為了克服上述困難，20世紀70年代初，加拿大多倫多大學的D.W.Strangway教授和他的學生提出沿用AMT的觀測方式，觀測人工供電產生的音頻電磁場。由于所觀測電磁場的頻率、強度和方向可以由人工控制，而其觀測方式又與AMT法相同，所以叫可控源音頻大地電磁法。

EH4連續電導率剖面儀是由美國EMI公司生產的電磁法儀器系統，屬于混合源方法的儀器。它利用的場源可以是天然場，也可以是人工場。數據采集方式和AMT法相同；在頻率（10～1000）Hz范圍內，采集天然電磁場信號；在頻率（750～92K）Hz范圍內，天然場高頻成分信號比較弱，使用人工源信號，這時滿足波區條件的收發距較小，容易實現；并且發送信號的時間短，功率小，裝備和電源輕便。本儀器配備高頻磁探頭，觀測的電磁信號頻率范圍為（10～92k）Hz，可用于測量地下幾米至1千米范圍內介質電阻率的連續分布情況。

近十年來，國際上出現了一批代表性的電磁法數據處理解釋軟件，如美國Zonge公司開發的SCS2D反演軟件、意大利Geosystem公司開發的一個綜合性解釋平臺WinGlink、吉林大學開發的GeoElectro電法數據處理系統和美國某大學實驗室研發的Sinv2d反演軟件等，不同反演軟件產生的反演模型擬合差因數據量及側重深度大小不一，使用一種成熟的處理軟件不僅便于資料對比，同時也為資料的地質解譯提供了較好的物探參考。

隨著電磁法的普遍應用，加拿大鳳凰公司的V8儀器越來越被新疆眾多地質單位青睞，配套的SCS2D反演軟件也逐漸展露其頭角，本文利用第四代計算機語言Matlab為平臺，對EH4儀器采集的數據格式進行轉換，使之能夠通過SCS2D反演，大大提高了數據處理的能力。

1 數據格式

軟件所涉及的文件包括EH4儀器自身采集并初步處理取得的Z文件（文件個數與點位相同）、數據處理所需的參數文件*.par、生成供SCS2D反演的*.mtd文件。

1）參數文件。參數文件（*.par）實際上是一個索引目錄文件，它展示了實際點位、高程與該點位所測原始數據的對應關系。文件共五列，分別為點號、坐標x、坐標y、高程、文件號，可用文本進行編輯。

2）Z文件格式。EH4儀器首先采集的是人工和天然大地電磁場強度變化的時間序列，然后內部進行傅立葉變換得出頻率域幅值響應再通過卡尼亞電阻率公式（見式（1））計算后把結果寫入Z文件，其內容包括頻點、ExHy的相干系數、ExHy標量視電阻率、ExHy標量相位、EyHx的相干系數、EyHx標量視電阻率、EyHx標量相位及各方向阻抗張量的實部及虛部等。

（1）

式中：f—頻率（Hz）；ρ-視電阻率（Ω·m）。

3）*.MTD文件格式。*.MTD是測線所有點位參數的集成文件，用逗號分隔，遺漏數據在行中用雙逗號表示。如圖1，由Z文件轉換而來的數據主要覆蓋5～9列，1～4列根據參數文件（*.par）直接寫入。1～9列表示的含義見表1。

2 Matlab平臺下的數據格式轉換過程

Matlab是矩陣實驗室（Matrix Laboratory）的簡稱，是美國MathWorks公司出品的商業數學軟件，用于算法開發、數據可視化、數據分析以及數值計算的高級技術計算語言和交互式環境。

1）友好的工作平臺和編程環境。Matlab由一系列工具組成。這些工具方便用戶使用Matlab的函數和文件，其中許多工具采用的是圖形用戶界面（圖2）。包括Matlab桌面和命令窗口、歷史命令窗口、編輯器和調試器、路徑搜索和用于用戶瀏覽幫助、工作空間、文件瀏覽器等。隨著Matlab的商業化以及軟件本身的不斷升級，Matlab的用戶界面也越來越精致，更加接近Windows的標準界面，人機交互性更強，操作更簡單。而且新版本的Matlab提供了完整的聯機查詢、幫助系統，極大的方便了用戶的使用。簡單的編程環境提供了比較完備的調試系統，程序不必經過編譯就可以直接運行，而且能夠及時地報告出現的錯誤及進行出錯原因分析。

圖2 程序運行界面

2）簡單通用的程序語言。Matlab使用的第四代計算機語言包含控制語句、函數、數據結構、輸入和輸出語句等，具有面向對象編程的特點。用戶可以在命令窗口中將輸入語句與執行命令同步，也可以先編寫好一個較大的復雜的M文件后再一起運行（圖3為已編好的M文件運行界面）。新版本的Matlab語言是基于最為流行的C++語言基礎上的，因此語法特征與C++語言極為相似，而且更加簡單，更加符合科技人員對數學表達式的書寫格式。使之更利于非計算機專業的科技人員使用。而且這種語言可移植性好、可拓展性極強，這也是Matlab能夠深入到科學研究及工程計算各個領域的重要原因。endprint

圖3 運行程序調用參數文件

3）軟件流程。野外EH4儀器數據采集完畢，測量臺班會記錄并形成一個關于點號、坐標、高程、文件號的索引參數文件，根據該文件讀取原始z文件中xy方向的相關系數、視電阻率、相位、頻率等，Matlab程序會根據提供的各種數據組合成SCS2D反演軟件所需的mtd格式文件，并最終輸出。軟件大致流程如圖4。

圖4 軟件流程圖

4）程序代碼及注釋。

clear all；clc；

[cfilename，cdir]=uigetfile（{'*.par'，'parameter file（*.par）'；'*.*'，'All files（*.*）'}，...

'Choose a parameter file to read'）；%選擇參數文件并讀取。

cfile=[cdir，cfilename]；

fid=fopen（cfile，'r'）；

data=fscanf（fid，'%f'，[5，inf]）；%參數文件分四列，為點號，坐標x，坐標y，高程，文件號。

data=data'；

nsite=size（data，1）；

j=1；

for isite=1：nsite%進入循環語句。

sname=num2str（data（isite，5））；

slength=length（sname）；

if slength==1

sname=['00'，sname]；

elseif slength==2

sname=['0'，sname]；%確定z文件名。

end

zfile=['zhs.'，sname]；

zfid=fopen（zfile，'r'）；

zdata=fscanf（zfid，'%f，[15，inf]）；%讀取z文件

zdata=zdata'；

for i=1：39；

if zdata（i，3）<=0；

continue；

else

mdata（j，6）=zdata（i，3）；

mdata（j，5）=zdata（i，1）；

mdata（j，1：4）=data（isite，1：4）；

mdata（j，7）=0.5；

mdata（j，8）=zdata（i，4）*3.1415926*1000/180；

mdata（j，9）=0.5；

mdata（j，10：13）=0；

mdata（j，14）=mdata（j，6）；

mdata（j，15）=mdata（j，8）；

mdata（j，16：17）=0；%生成mtd文件矩陣

j=j+1；

end

end

fclose（zfid）；

end

fclose （fid）；

mtdfile=['line'，num2str（data（1，2）），'.mtd']；

fid=fopen（mtdfile，'w'）；

fprintf（fid，'%s＼n'，'"Stn"，"GridE"，"GridN"，"Elev"，"Freq"，"ARTMobs"，"ARTMerr"，"ZPTMobs"，"ZPTMerr"，"ARTEobs"，"ARTEerr"，"ZPTEobs"，"ZPTEerr"'）；

fprintf（fid，'%s＼n'，'"From S2D2D v2.10l Date：16/12/09 Time：13：34：14"'）；

for i=1：j-1

fprintf（fid，'%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f'，mdata（i，1：9））；

fprintf（fid，'%s'，'，，，，'）；

fprintf（fid，'%6.1f，%6.1f'，mdata（i，14：15））；

fprintf（fid，'%s＼n'，'，，，'）；%輸出mtd文件

end

fclose（fid）；

3 實例分析

1）工區概況。庫普邊防派出所位于奇臺縣八一牧場西北側，距奇臺縣城約197千米，為昌吉邊防支隊的下屬單位，多年來一直處于缺水狀態。為做好雙擁工作，更為解決官兵的實際用水困難，應昌吉國土資源局邀請，物探隊派遣EH4項目組，了解該區地下電阻率分布特征，為尋找儲水構造并布置鉆井提供依據。

該區地處東準噶爾堿性花崗巖體中，尋找裂隙水是本次工作的主攻目標。

2）工作布置。如圖5，限于工區實際條件，布置測線：3線、4線、5線。取點數較多的3線為例，3線307點西北側為學校辦公樓，樓前308點西約6米處為已知水井，該水井為前些年所打，井深10米，基底為花崗巖。由于上游水量降低，該井出水量日益減少，已不能滿足學校最基本生活需求。

圖5 庫普工區點位圖

3）反演成果。如圖6，三種處理方法均能反應測線地下電阻率分布情況：采用美國某大學的REBOCC反演程序側重于地下介質的均勻性，對于地質體整體的分布特征把握較好，處理程度較圓滑，但對于細節反應不明顯；而采用經過本程序格式轉換后的SCS2D反演成果能夠較好的描述整體和局部特征，特別是二維反演成果，非常好的詮釋了地下地質體、斷裂構造的產狀、埋深以及地下水位的分布特征。

圖6 不同處理方法反演成果圖

4）總結。通過本次工作，查清了已知水井缺水的原因，提出了更利于儲水的307點；同時也推斷由于深部下游水源通道的存在，下游不會因為已知水井斷水而沒有開發前景，為本區水源開發利用奠定了較好的基礎。

4 結束語

Matlab工作平臺的使用不僅大大優化了不同數據格式的轉換流程，提高了野外的工作效率，也避免了優秀反演軟件無法通用的弊端，多用于信號和圖像處理、通訊、控制系統設計、測試和測量等眾多應用領域，對地學數據處理解釋起到了推動作用。

參考文獻

[1]鄧明，譚捍東，胡建德，等.超寬頻帶大地電磁信號的采集方法和技術[J].現代地質，1997，11（3）.

[2]魏文博.我國大地電磁測深新進展及瞻望[J].地球物理學進展，2002，17（2）.

[3]TAO Shi long，WAN Tian feng.Introduction of Earth Science[M].Beijing： Geological Publishing House， 1999.

[4]http：//baike.baidu.com/view/10598.htm.

作者簡介

李美英（1981-），河北涿州人，助理工程師，主要從事儀器設備維修、維護及數據處理等工作。endprint

圖3 運行程序調用參數文件

3）軟件流程。野外EH4儀器數據采集完畢，測量臺班會記錄并形成一個關于點號、坐標、高程、文件號的索引參數文件，根據該文件讀取原始z文件中xy方向的相關系數、視電阻率、相位、頻率等，Matlab程序會根據提供的各種數據組合成SCS2D反演軟件所需的mtd格式文件，并最終輸出。軟件大致流程如圖4。

圖4 軟件流程圖

4）程序代碼及注釋。

clear all；clc；

[cfilename，cdir]=uigetfile（{'*.par'，'parameter file（*.par）'；'*.*'，'All files（*.*）'}，...

'Choose a parameter file to read'）；%選擇參數文件并讀取。

cfile=[cdir，cfilename]；

fid=fopen（cfile，'r'）；

data=fscanf（fid，'%f'，[5，inf]）；%參數文件分四列，為點號，坐標x，坐標y，高程，文件號。

data=data'；

nsite=size（data，1）；

j=1；

for isite=1：nsite%進入循環語句。

sname=num2str（data（isite，5））；

slength=length（sname）；

if slength==1

sname=['00'，sname]；

elseif slength==2

sname=['0'，sname]；%確定z文件名。

end

zfile=['zhs.'，sname]；

zfid=fopen（zfile，'r'）；

zdata=fscanf（zfid，'%f，[15，inf]）；%讀取z文件

zdata=zdata'；

for i=1：39；

if zdata（i，3）<=0；

continue；

else

mdata（j，6）=zdata（i，3）；

mdata（j，5）=zdata（i，1）；

mdata（j，1：4）=data（isite，1：4）；

mdata（j，7）=0.5；

mdata（j，8）=zdata（i，4）*3.1415926*1000/180；

mdata（j，9）=0.5；

mdata（j，10：13）=0；

mdata（j，14）=mdata（j，6）；

mdata（j，15）=mdata（j，8）；

mdata（j，16：17）=0；%生成mtd文件矩陣

j=j+1；

end

end

fclose（zfid）；

end

fclose （fid）；

mtdfile=['line'，num2str（data（1，2）），'.mtd']；

fid=fopen（mtdfile，'w'）；

fprintf（fid，'%s＼n'，'"Stn"，"GridE"，"GridN"，"Elev"，"Freq"，"ARTMobs"，"ARTMerr"，"ZPTMobs"，"ZPTMerr"，"ARTEobs"，"ARTEerr"，"ZPTEobs"，"ZPTEerr"'）；

fprintf（fid，'%s＼n'，'"From S2D2D v2.10l Date：16/12/09 Time：13：34：14"'）；

for i=1：j-1

fprintf（fid，'%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f'，mdata（i，1：9））；

fprintf（fid，'%s'，'，，，，'）；

fprintf（fid，'%6.1f，%6.1f'，mdata（i，14：15））；

fprintf（fid，'%s＼n'，'，，，'）；%輸出mtd文件

end

fclose（fid）；

3 實例分析

1）工區概況。庫普邊防派出所位于奇臺縣八一牧場西北側，距奇臺縣城約197千米，為昌吉邊防支隊的下屬單位，多年來一直處于缺水狀態。為做好雙擁工作，更為解決官兵的實際用水困難，應昌吉國土資源局邀請，物探隊派遣EH4項目組，了解該區地下電阻率分布特征，為尋找儲水構造并布置鉆井提供依據。

該區地處東準噶爾堿性花崗巖體中，尋找裂隙水是本次工作的主攻目標。

2）工作布置。如圖5，限于工區實際條件，布置測線：3線、4線、5線。取點數較多的3線為例，3線307點西北側為學校辦公樓，樓前308點西約6米處為已知水井，該水井為前些年所打，井深10米，基底為花崗巖。由于上游水量降低，該井出水量日益減少，已不能滿足學校最基本生活需求。

圖5 庫普工區點位圖

3）反演成果。如圖6，三種處理方法均能反應測線地下電阻率分布情況：采用美國某大學的REBOCC反演程序側重于地下介質的均勻性，對于地質體整體的分布特征把握較好，處理程度較圓滑，但對于細節反應不明顯；而采用經過本程序格式轉換后的SCS2D反演成果能夠較好的描述整體和局部特征，特別是二維反演成果，非常好的詮釋了地下地質體、斷裂構造的產狀、埋深以及地下水位的分布特征。

圖6 不同處理方法反演成果圖

4）總結。通過本次工作，查清了已知水井缺水的原因，提出了更利于儲水的307點；同時也推斷由于深部下游水源通道的存在，下游不會因為已知水井斷水而沒有開發前景，為本區水源開發利用奠定了較好的基礎。

4 結束語

Matlab工作平臺的使用不僅大大優化了不同數據格式的轉換流程，提高了野外的工作效率，也避免了優秀反演軟件無法通用的弊端，多用于信號和圖像處理、通訊、控制系統設計、測試和測量等眾多應用領域，對地學數據處理解釋起到了推動作用。

參考文獻

[1]鄧明，譚捍東，胡建德，等.超寬頻帶大地電磁信號的采集方法和技術[J].現代地質，1997，11（3）.

[2]魏文博.我國大地電磁測深新進展及瞻望[J].地球物理學進展，2002，17（2）.

[3]TAO Shi long，WAN Tian feng.Introduction of Earth Science[M].Beijing： Geological Publishing House， 1999.

[4]http：//baike.baidu.com/view/10598.htm.

作者簡介

李美英（1981-），河北涿州人，助理工程師，主要從事儀器設備維修、維護及數據處理等工作。endprint

圖3 運行程序調用參數文件

3）軟件流程。野外EH4儀器數據采集完畢，測量臺班會記錄并形成一個關于點號、坐標、高程、文件號的索引參數文件，根據該文件讀取原始z文件中xy方向的相關系數、視電阻率、相位、頻率等，Matlab程序會根據提供的各種數據組合成SCS2D反演軟件所需的mtd格式文件，并最終輸出。軟件大致流程如圖4。

圖4 軟件流程圖

4）程序代碼及注釋。

clear all；clc；

[cfilename，cdir]=uigetfile（{'*.par'，'parameter file（*.par）'；'*.*'，'All files（*.*）'}，...

'Choose a parameter file to read'）；%選擇參數文件并讀取。

cfile=[cdir，cfilename]；

fid=fopen（cfile，'r'）；

data=fscanf（fid，'%f'，[5，inf]）；%參數文件分四列，為點號，坐標x，坐標y，高程，文件號。

data=data'；

nsite=size（data，1）；

j=1；

for isite=1：nsite%進入循環語句。

sname=num2str（data（isite，5））；

slength=length（sname）；

if slength==1

sname=['00'，sname]；

elseif slength==2

sname=['0'，sname]；%確定z文件名。

end

zfile=['zhs.'，sname]；

zfid=fopen（zfile，'r'）；

zdata=fscanf（zfid，'%f，[15，inf]）；%讀取z文件

zdata=zdata'；

for i=1：39；

if zdata（i，3）<=0；

continue；

else

mdata（j，6）=zdata（i，3）；

mdata（j，5）=zdata（i，1）；

mdata（j，1：4）=data（isite，1：4）；

mdata（j，7）=0.5；

mdata（j，8）=zdata（i，4）*3.1415926*1000/180；

mdata（j，9）=0.5；

mdata（j，10：13）=0；

mdata（j，14）=mdata（j，6）；

mdata（j，15）=mdata（j，8）；

mdata（j，16：17）=0；%生成mtd文件矩陣

j=j+1；

end

end

fclose（zfid）；

end

fclose （fid）；

mtdfile=['line'，num2str（data（1，2）），'.mtd']；

fid=fopen（mtdfile，'w'）；

fprintf（fid，'%s＼n'，'"Stn"，"GridE"，"GridN"，"Elev"，"Freq"，"ARTMobs"，"ARTMerr"，"ZPTMobs"，"ZPTMerr"，"ARTEobs"，"ARTEerr"，"ZPTEobs"，"ZPTEerr"'）；

fprintf（fid，'%s＼n'，'"From S2D2D v2.10l Date：16/12/09 Time：13：34：14"'）；

for i=1：j-1

fprintf（fid，'%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f，%6.1f'，mdata（i，1：9））；

fprintf（fid，'%s'，'，，，，'）；

fprintf（fid，'%6.1f，%6.1f'，mdata（i，14：15））；

fprintf（fid，'%s＼n'，'，，，'）；%輸出mtd文件

end

fclose（fid）；

3 實例分析

1）工區概況。庫普邊防派出所位于奇臺縣八一牧場西北側，距奇臺縣城約197千米，為昌吉邊防支隊的下屬單位，多年來一直處于缺水狀態。為做好雙擁工作，更為解決官兵的實際用水困難，應昌吉國土資源局邀請，物探隊派遣EH4項目組，了解該區地下電阻率分布特征，為尋找儲水構造并布置鉆井提供依據。

該區地處東準噶爾堿性花崗巖體中，尋找裂隙水是本次工作的主攻目標。

2）工作布置。如圖5，限于工區實際條件，布置測線：3線、4線、5線。取點數較多的3線為例，3線307點西北側為學校辦公樓，樓前308點西約6米處為已知水井，該水井為前些年所打，井深10米，基底為花崗巖。由于上游水量降低，該井出水量日益減少，已不能滿足學校最基本生活需求。

圖5 庫普工區點位圖

3）反演成果。如圖6，三種處理方法均能反應測線地下電阻率分布情況：采用美國某大學的REBOCC反演程序側重于地下介質的均勻性，對于地質體整體的分布特征把握較好，處理程度較圓滑，但對于細節反應不明顯；而采用經過本程序格式轉換后的SCS2D反演成果能夠較好的描述整體和局部特征，特別是二維反演成果，非常好的詮釋了地下地質體、斷裂構造的產狀、埋深以及地下水位的分布特征。

圖6 不同處理方法反演成果圖

4）總結。通過本次工作，查清了已知水井缺水的原因，提出了更利于儲水的307點；同時也推斷由于深部下游水源通道的存在，下游不會因為已知水井斷水而沒有開發前景，為本區水源開發利用奠定了較好的基礎。

4 結束語

Matlab工作平臺的使用不僅大大優化了不同數據格式的轉換流程，提高了野外的工作效率，也避免了優秀反演軟件無法通用的弊端，多用于信號和圖像處理、通訊、控制系統設計、測試和測量等眾多應用領域，對地學數據處理解釋起到了推動作用。

參考文獻

[1]鄧明，譚捍東，胡建德，等.超寬頻帶大地電磁信號的采集方法和技術[J].現代地質，1997，11（3）.

[2]魏文博.我國大地電磁測深新進展及瞻望[J].地球物理學進展，2002，17（2）.

[3]TAO Shi long，WAN Tian feng.Introduction of Earth Science[M].Beijing： Geological Publishing House， 1999.

[4]http：//baike.baidu.com/view/10598.htm.

作者簡介

李美英（1981-），河北涿州人，助理工程師，主要從事儀器設備維修、維護及數據處理等工作。endprint