喬玉，陳凱，陽(yáng)琴

(1.中國(guó)地質(zhì)大學(xué)(北京) 地球物理與信息技術(shù)學(xué)院,北京 100083；2.中國(guó)電子科技集團(tuán)公司 第五十八研究所，江蘇 無(wú)錫 214000)

0 引言

海底大地電磁方法是通過(guò)測(cè)量海底大地電磁信號(hào)，經(jīng)數(shù)據(jù)處理及反演得到海底以下地層的電性信息，進(jìn)而獲取海底深部電性結(jié)構(gòu)。20世紀(jì)末，Steven等[1]已設(shè)計(jì)出可在1000 m深水區(qū)獲得良好的電磁響應(yīng)的海底采集系統(tǒng)。經(jīng)過(guò)近些年的發(fā)展，在設(shè)備研制、資料處理方法上取得了長(zhǎng)足進(jìn)展，但國(guó)外在設(shè)備研制方面一直對(duì)我國(guó)進(jìn)行技術(shù)封鎖。為了打破國(guó)外在OBEM設(shè)備方面的壟斷，國(guó)內(nèi)同行從1998年開(kāi)始進(jìn)行海洋大地電磁相關(guān)的方法理論、儀器設(shè)備、資料處理、反演解釋等研究工作，近年來(lái)獲得了顯著的成功[2-7]。中國(guó)地質(zhì)大學(xué)(北京)在多期“863計(jì)劃”課題及國(guó)家專項(xiàng)的支持下，建立了基于ARM 的海底大地電磁數(shù)據(jù)采集系統(tǒng)——OBEM-Ⅲ，初步實(shí)現(xiàn)了低功耗、小型化以及智能化[8]，在油氣勘探、水合物調(diào)查、深部構(gòu)造研究等領(lǐng)域取得了一定的成果[9-10]，儀器的主要技術(shù)指標(biāo)已經(jīng)與以美國(guó)SIO為代表的國(guó)際先進(jìn)同行齊平。在國(guó)內(nèi)外同行的推進(jìn)下，海底大地電磁方法技術(shù)在海域油氣勘探[11]、地下水調(diào)查[12]等資源環(huán)境領(lǐng)域和洋脊擴(kuò)張、板塊俯沖[13]、海底火山[14]等構(gòu)造地質(zhì)研究領(lǐng)域取得了顯著的成果。

在野外測(cè)量之前[2]，為確認(rèn)儀器性能并校正測(cè)量數(shù)據(jù)中包含的儀器通道響應(yīng)，通常需要對(duì)儀器開(kāi)展通道標(biāo)定。通道標(biāo)定的過(guò)程包括產(chǎn)生標(biāo)準(zhǔn)方波并采集存儲(chǔ)和標(biāo)定計(jì)算兩個(gè)過(guò)程，涉及硬件電路控制以及軟件數(shù)據(jù)處理：① 在硬件電路方面，通過(guò)相關(guān)命令控制集成在采集艙內(nèi)部的標(biāo)準(zhǔn)信號(hào)發(fā)生模塊產(chǎn)生3組不同頻率的方波，并送至采集通道的最前端進(jìn)行采集存儲(chǔ)，生成通道標(biāo)定文件(簡(jiǎn)稱標(biāo)定文件)；② 在數(shù)據(jù)處理程序中，讀取標(biāo)定文件中的數(shù)據(jù)，將該數(shù)據(jù)與原標(biāo)準(zhǔn)方波(標(biāo)定信號(hào)輸入)進(jìn)行比較計(jì)算，得到通道的標(biāo)定結(jié)果，比較計(jì)算的過(guò)程，稱為標(biāo)定計(jì)算。在無(wú)需外接任何信號(hào)發(fā)生器的情況下，由標(biāo)定結(jié)果可判斷儀器工作狀態(tài)，并為后期資料處理提供通道響應(yīng)改正數(shù)據(jù)。

OBEM-Ⅲ是基于ARM-Linux平臺(tái)開(kāi)發(fā)的數(shù)據(jù)采集系統(tǒng)，但現(xiàn)有的標(biāo)定計(jì)算需借助Windows平臺(tái)開(kāi)展，每次標(biāo)定計(jì)算前需將標(biāo)定文件下載到本地計(jì)算機(jī)，再通過(guò)Matlab程序進(jìn)行計(jì)算，存在用戶操作繁瑣、作業(yè)效率低等缺點(diǎn)。為提高海上作業(yè)效率，亟需開(kāi)發(fā)一個(gè)標(biāo)定計(jì)算程序，并將其封裝為可在ARM-Linux平臺(tái)上運(yùn)行的APP，可在海底電磁接收機(jī)本地端自主完成標(biāo)定計(jì)算。

1 硬件電路原理

標(biāo)定過(guò)程的硬件電路原理框如圖1所示，主要由選擇開(kāi)關(guān)、采集電路、標(biāo)定信號(hào)產(chǎn)生電路、標(biāo)定計(jì)算電路組成。

選擇開(kāi)關(guān)為單刀雙擲繼電器，可通過(guò)微控制單元(MCU，microcontroller unit)控制采集電路接入信號(hào)，通道標(biāo)定時(shí)切換到標(biāo)準(zhǔn)信號(hào)產(chǎn)生電路，實(shí)際信號(hào)測(cè)量時(shí)切換到外部傳感器信號(hào)輸入。

標(biāo)定信號(hào)產(chǎn)生電路由復(fù)雜可編程邏輯器件(CPLD, complex programmable logic device)、晶振、斬波電路以及衰減電路組成。標(biāo)準(zhǔn)信號(hào)的產(chǎn)生過(guò)程如下：CPLD 對(duì)晶振產(chǎn)生的信號(hào)進(jìn)行分頻，輸出峰峰值為3.3 V的單極性方波，再經(jīng)過(guò)斬波電路和衰減電路，得到峰峰值為100 mVpp 的雙極性方波，即標(biāo)定信號(hào)。斬波電路運(yùn)用比較器的原理，將CPLD輸出的單極性方波轉(zhuǎn)換為與其頻率一致的雙極性方波，通過(guò)改變CPLD對(duì)晶振信號(hào)的分頻系數(shù)，即可輸出不同頻率的標(biāo)定信號(hào)，晶振的頻率準(zhǔn)確性與斬波電路的性能，直接影響產(chǎn)生標(biāo)定信號(hào)的準(zhǔn)確性。

采集電路包括運(yùn)放(amp, amplifier)、模數(shù)轉(zhuǎn)換器(ADC，analogy-to-digital converter)、CPLD、MCU、GPS、SD 卡等，對(duì)輸入信號(hào)進(jìn)行采集、存儲(chǔ)。電路進(jìn)行通道之前，由MCU獲取當(dāng)前GPS時(shí)間并記錄，同時(shí)儀器本地時(shí)間以此為準(zhǔn)繼續(xù)運(yùn)行。通道標(biāo)定開(kāi)始后，由標(biāo)定信號(hào)產(chǎn)生電路生成的信號(hào)通過(guò)運(yùn)放進(jìn)行濾波、放大，經(jīng)ADC進(jìn)行模數(shù)轉(zhuǎn)換、 CPLD 完成數(shù)據(jù)整合，最后由 MCU將采集到的數(shù)據(jù)寫(xiě)入SD卡內(nèi)，生成對(duì)應(yīng)的標(biāo)定文件。后期用Matlab 對(duì)采集到的標(biāo)定文件進(jìn)行復(fù)原，其在時(shí)域上表現(xiàn)為方波且無(wú)明顯畸變，幅值準(zhǔn)確性達(dá)±1%；頻域上基本滿足標(biāo)準(zhǔn)方波的頻譜分布，相位精度小于±0.37152°。

標(biāo)定計(jì)算電路包括MCU 和SD 卡，MCU 讀取SD卡中的標(biāo)定文件，通過(guò)在頻域上與標(biāo)準(zhǔn)信號(hào)進(jìn)行比較計(jì)算，輸出標(biāo)定結(jié)果文件。

圖1 標(biāo)定原理框Fig.1 Block diagram of calibration

2 計(jì)算原理

標(biāo)定計(jì)算程序以二進(jìn)制原始時(shí)間序列為輸入，以十進(jìn)制標(biāo)定計(jì)算結(jié)果為輸出，程序設(shè)計(jì)包括3部分：讀取標(biāo)定文件、對(duì)標(biāo)定數(shù)據(jù)進(jìn)行頻譜分析、比較計(jì)算得到輸出結(jié)果文件。

2.1 讀取標(biāo)定文件

標(biāo)定文件包含TSH、TSM、TSL這3組數(shù)據(jù)，分別對(duì)應(yīng)不同的采樣率以及標(biāo)定信號(hào)頻率，詳細(xì)的參數(shù)見(jiàn)表1，其中基頻代表輸入標(biāo)定信號(hào)的頻率，保留諧波數(shù)是頻譜分析時(shí)保留的最高次諧波。在進(jìn)行數(shù)據(jù)讀取之前，需判斷輸入標(biāo)定文件的類型并設(shè)置對(duì)應(yīng)的參數(shù)。

表1 標(biāo)定文件的部分參數(shù)

每個(gè)標(biāo)定文件包括16 B的塊頭以及1MB的數(shù)據(jù)體。塊頭中包含數(shù)據(jù)塊 ID 和相對(duì)時(shí)間δt，ID 表示當(dāng)前數(shù)據(jù)塊的編號(hào)，δt代表與對(duì)鐘時(shí)間相間隔的時(shí)間。數(shù)據(jù)體中包含8CH×32768個(gè)采樣點(diǎn)的數(shù)據(jù)，每個(gè)采樣點(diǎn)數(shù)據(jù)占4 B，以小端模式存儲(chǔ)，組織形式為：

L+M+H+channel_flag,

(1)

式中：L,M,H分別對(duì)應(yīng)低、中、高位數(shù)據(jù)，channel_flag為通道標(biāo)志位。單通道數(shù)據(jù)的實(shí)際值可表示為式(2)，單位：V：

(2)

2.2 對(duì)標(biāo)定數(shù)據(jù)進(jìn)行頻譜分析

從標(biāo)定文件中讀取出來(lái)的數(shù)據(jù)，均被放置在有限長(zhǎng)的數(shù)組中，等待進(jìn)行下一步的數(shù)據(jù)處理。標(biāo)定計(jì)算是主要完成對(duì)采集信號(hào)的頻譜分析。對(duì)于有限長(zhǎng)時(shí)間序列，通常采用離散傅里葉變換(DFT , discrete fourier transform)進(jìn)行頻譜分析。

對(duì)于N點(diǎn)序列x(n)，其DFT變換的定義為：

(3)

0≤k≤N-1(4)

從式(4)可以看出，N點(diǎn)DFT計(jì)算需要N2次復(fù)數(shù)乘法運(yùn)算、N(N-1)次復(fù)數(shù)加法運(yùn)算。因?yàn)閷?shí)現(xiàn)一次復(fù)數(shù)相乘，需要四次實(shí)數(shù)相乘和兩次實(shí)數(shù)相加，所以復(fù)數(shù)乘法的計(jì)算次數(shù)直接影響程序的運(yùn)行時(shí)間，特別是當(dāng)N很大時(shí)，計(jì)算量呈指數(shù)增加。由于需要處理的數(shù)據(jù)較多，且程序的目標(biāo)運(yùn)行平臺(tái)計(jì)算速度不高，直接套用式(3)，運(yùn)算效率低，耗費(fèi)時(shí)間多，不利于實(shí)時(shí)對(duì)數(shù)據(jù)進(jìn)行處理。為了滿足現(xiàn)場(chǎng)作業(yè)的需求，選用快速傅里葉變換(FFT, fast fourier transform)對(duì)數(shù)據(jù)進(jìn)行處理。

2.2.1 混合基FFT

FFT計(jì)算(又稱蝶形計(jì)算)是式(3)經(jīng)過(guò)變形之后的圖形表達(dá)，其結(jié)構(gòu)可清楚的表現(xiàn)數(shù)據(jù)之間的依賴關(guān)系。每完成一組蝶形運(yùn)算，即可釋放輸入值的存儲(chǔ)空間給輸出值，實(shí)現(xiàn)同址計(jì)算，減少存儲(chǔ)空間消耗。常見(jiàn)的運(yùn)算對(duì)有基2、基3等，表2為長(zhǎng)度為N時(shí)FFT 計(jì)算與DFT 計(jì)算復(fù)數(shù)乘法的運(yùn)算量比較情況，由表可以看出，隨著N的不斷增加，與直接DFT計(jì)算相比，固定基算法計(jì)算量的優(yōu)勢(shì)越來(lái)越明顯，但常見(jiàn)的固定基FFT運(yùn)算只能用于處理長(zhǎng)度為Qk(Q為基底數(shù)，k為整數(shù))的數(shù)據(jù)，由表1可知，3 種文件的FFT計(jì)算長(zhǎng)度并不能滿足這一要求。在Windows系統(tǒng)上，可通過(guò)末位補(bǔ)零法，增加數(shù)據(jù)長(zhǎng)度，使其達(dá)到固定基FFT對(duì)輸入序列的長(zhǎng)度要求但相同的方法在ARM-Linux 平臺(tái)上會(huì)大大增加運(yùn)算時(shí)間。

混合基FFT運(yùn)算是固定基FFT運(yùn)算推廣后的一般情況，受到計(jì)算點(diǎn)數(shù)的限制更少，可以對(duì)點(diǎn)數(shù)為非質(zhì)數(shù)的DFT 進(jìn)行快速計(jì)算，是標(biāo)定計(jì)算的最佳選擇。本程序中按照時(shí)間抽取的方法，先將X(k)逐層展開(kāi)到x(n)，將x(n)進(jìn)行相應(yīng)的排序后，再?gòu)膞(n)逐層計(jì)算回X(k)。

表2 DFT與FFT復(fù)數(shù)乘法運(yùn)算量比較

2.2.2 多基多進(jìn)制排序

由于對(duì)x(n)需進(jìn)行逐次抽取，必須對(duì)輸入序列倒位序，得到的次序相當(dāng)于是原始序列次序的碼位倒置，這一過(guò)程稱為多基多進(jìn)制的排序。以二進(jìn)制數(shù)為例，任意整數(shù)N都可表示為：

(5)

式中:nr,nr-1,…,n2,n1均小于2，等號(hào)左邊為二進(jìn)制表達(dá)方式，括號(hào)外的下角標(biāo)2代表每位都以2為基底；等號(hào)右邊為其對(duì)應(yīng)的十進(jìn)制數(shù)表示，經(jīng)過(guò)倒序后，其新的位序?yàn)?n1n2…nr-1nr)=n1×2r-1+n2×2r-2+…+nr-1×21+nr×20，新的位序也代表了其新的輸入順序。在多基多進(jìn)制的情況下，F(xiàn)FT計(jì)算長(zhǎng)度N可以表示為：

(6)

式中nr,nr-1,…,n2,n1分別小于xr,xr-1,…,x2,x1，等號(hào)左邊為多進(jìn)制表達(dá)方式，括號(hào)外的下角標(biāo)xr,xr-1,…,x2,x1分別代表每一位的基底；等號(hào)右邊為其對(duì)應(yīng)的十進(jìn)制數(shù)表示。同理，其倒敘后的新位序N’表示為:

(7)

倒序之后的序列按照對(duì)應(yīng)的十進(jìn)制大小從小到大進(jìn)行排序。經(jīng)過(guò)排序后的新序列逐層抽取計(jì)算，便可得到最終的FFT計(jì)算結(jié)果。由于旋轉(zhuǎn)因子WNnk具有周期性，計(jì)算量可進(jìn)一步縮小。圖2為N=12=3×2×2時(shí)的逐層抽取的過(guò)程，其中x’(n)為經(jīng)過(guò)倒序后新序列，X11(k)為x’(n)經(jīng)過(guò)第一次基3抽取的結(jié)果，X1(k)序列為X11(k)經(jīng)過(guò)基二抽取的結(jié)果，X(k)為X1(k)經(jīng)過(guò)基二抽取得到的最終結(jié)果。圖3為整個(gè)混合基FFT的流程框，其中index2, index3, index5分別為2,3,5的冪，即混合基計(jì)算的長(zhǎng)度LEN=2index2×3index3×5index5，基2FFT、基3FFT與基5FFT計(jì)算采用式(3)進(jìn)行。

2.3 比較計(jì)算

對(duì)標(biāo)定文件完成FFT計(jì)算之后，需與輸入的標(biāo)定信號(hào)進(jìn)行比較分析得到通道輸入響應(yīng)。在程序中定義幅值分別為1250以及1075000的理想方波來(lái)代替輸入的標(biāo)定信號(hào)，其長(zhǎng)度與FFT計(jì)算長(zhǎng)度相同，也進(jìn)行同樣的頻譜分析。進(jìn)行完FFT計(jì)算后的理想方波與標(biāo)定文件在相同頻點(diǎn)處進(jìn)行幅值相除、相位相減，再存儲(chǔ)到SD卡中，便生成最后的結(jié)果文件。

圖2 12點(diǎn)混合基FFT流程Fig.2 12-point flow chart for Mixed-radix

圖3 混合基FFT程序框Fig.3 Block diagram for Mixed-radix

3 室內(nèi)測(cè)試

編寫(xiě)好的程序通過(guò)交叉編譯后上傳至SD卡，輸入相關(guān)命令，即可在OBEM-Ⅲ進(jìn)行本地標(biāo)定計(jì)算。圖4為本地標(biāo)定計(jì)算的實(shí)際過(guò)程，從計(jì)算開(kāi)始到結(jié)束只需11 s，計(jì)算結(jié)果存儲(chǔ)到SD卡的新建文件中,即圖4中的PM21L.asc文件。

圖5為標(biāo)定計(jì)算結(jié)果的幅值相位圖，由圖可以看出，電道與磁道在通道增益以及相位方面的一致性較好，標(biāo)定計(jì)算的結(jié)果有效可靠。

圖4 標(biāo)定計(jì)算過(guò)程Fig.4 Progress of calibration

4 海上實(shí)測(cè)結(jié)果

2020年7～8月，在南海西南次海盆海底大地電磁探測(cè)科學(xué)任務(wù)中，投入了23臺(tái)海底電磁接收機(jī)，工作水深2 508～4 443 m，完成了46站位海底大地電磁采集任務(wù)。海上任務(wù)開(kāi)展前，均完成了通道標(biāo)定工作，圖6給出了現(xiàn)場(chǎng)經(jīng)預(yù)處理得到的測(cè)深曲線(S33站位)，有效觀測(cè)頻段為30～30 000 s，獲得了高質(zhì)量海底大地電磁測(cè)深數(shù)據(jù)。

5 結(jié)論

針對(duì)海底電磁接收機(jī)標(biāo)定計(jì)算程序依賴于 PC 端的 Matlab 平臺(tái)的不足，開(kāi)展了基于 ARM-Linux 平臺(tái)的計(jì)算程序開(kāi)發(fā)，設(shè)計(jì)最優(yōu)計(jì)算參數(shù)，獲得了通道頻率響應(yīng)，得到了標(biāo)定文件。經(jīng)室內(nèi)測(cè)試驗(yàn)證了實(shí)測(cè)結(jié)果的正確性，計(jì)算耗時(shí)約為11 s。該計(jì)算程序應(yīng)用于海上數(shù)據(jù)采集任務(wù)中，獲得了高質(zhì)量海底大地電磁測(cè)深數(shù)據(jù)，有效提升了海上現(xiàn)場(chǎng)作業(yè)效率。