楊晨，吳銀川，曹亢

西安石油大學(xué) 陜西省油氣井測(cè)控技術(shù)重點(diǎn)實(shí)驗(yàn)室（陜西 西安 710065）

0 引言

感應(yīng)測(cè)井技術(shù)是現(xiàn)階段我國(guó)在石油開采過程中應(yīng)用最廣泛的一項(xiàng)技術(shù)[1]，該技術(shù)通過發(fā)射不同頻率下的激勵(lì)信號(hào)來對(duì)應(yīng)不同地層的探測(cè)深度，當(dāng)采用多種頻率信號(hào)同時(shí)激勵(lì)時(shí)，在加快勘探進(jìn)程的同時(shí)還可以得到地層的多維信息從而提高精細(xì)勘探能力[2]，其中接收傳感器的精密檢測(cè)是感應(yīng)測(cè)井儀器設(shè)計(jì)的關(guān)鍵技術(shù)之一。常見的幅度相位檢測(cè)方法有相敏檢測(cè)和快速傅里葉變換兩種[3]，相敏檢測(cè)技術(shù)是通過鎖相檢測(cè)確定頻率信號(hào)的幅度和相位，在實(shí)時(shí)檢測(cè)系統(tǒng)中可實(shí)現(xiàn)邊采集邊計(jì)算，不需要大量的內(nèi)存空間，實(shí)時(shí)性好。但根據(jù)相敏檢測(cè)原理可知在進(jìn)行多頻信號(hào)檢測(cè)時(shí)，每1次檢測(cè)只能檢出1種頻率信號(hào)的幅度和相位，幾個(gè)頻率的信號(hào)則需要做幾次相敏檢波，因此效率不高?？焖俑道锶~變換（FFT）作為傳統(tǒng)的頻譜檢測(cè)方法，具有原理簡(jiǎn)單易實(shí)現(xiàn)的特點(diǎn)，是數(shù)字信號(hào)處理中頻譜分析最常用的方法[4]，并且可一次檢測(cè)出多種頻率信號(hào)的幅值及相位，提高檢測(cè)效率。

本文基于感應(yīng)測(cè)井技術(shù)的特性，研究設(shè)計(jì)了基于FFT的陣列感應(yīng)信號(hào)測(cè)量分析軟件。闡述了軟件實(shí)現(xiàn)信號(hào)頻譜計(jì)算的原理，推導(dǎo)了FFT算法實(shí)現(xiàn)信號(hào)頻率、幅值、相位公式，測(cè)試并分析了在噪聲環(huán)境下進(jìn)行整周期和非整周期采樣時(shí)，軟件檢測(cè)結(jié)果的誤差范圍；基于檢測(cè)結(jié)果，量化分析了軟件的精確度和穩(wěn)定性。

1 信號(hào)頻譜計(jì)算原理

1.1 信號(hào)頻譜計(jì)算的方法

本文采用FFT方法對(duì)信號(hào)頻譜進(jìn)行檢測(cè)，目前面向?qū)嶋H應(yīng)用層面，基2和基4 FFT算法是使用最廣泛的FFT算法[5]，其中基4時(shí)分FFT算法的復(fù)數(shù)乘法為基2時(shí)分FFT算法的3/4，但基4時(shí)分FFT算法較基2時(shí)分FFT算法可變換的數(shù)據(jù)長(zhǎng)度更少，靈活性較差。為使軟件具有更廣泛的數(shù)據(jù)處理寬度，本文采用DIT-FFT（基2時(shí)間抽取快速傅立葉變換）算法來實(shí)現(xiàn)正弦信號(hào)的頻譜分析要求。基2 FFT算法是利用旋轉(zhuǎn)因子的對(duì)稱性，將有限長(zhǎng)序列不斷進(jìn)行奇偶抽取，直至分解成長(zhǎng)度等于2的短序列，最終只需計(jì)算長(zhǎng)度等于2的短序列的DFT變換[6]。

1.2 頻譜計(jì)算的具體實(shí)現(xiàn)

設(shè)被測(cè)信號(hào)為：

式中：A為信號(hào)幅度，V；f0為信號(hào)頻率，Hz；θ是信號(hào)初相位，（°）。

在滿足奈奎斯特采樣定理的情況下，以采樣頻率fs=1/Ts對(duì)信號(hào)進(jìn)行采樣，可得離散正弦信號(hào)為：

式中：n=0,1,2,…,N-1,N為采樣點(diǎn)數(shù)。根據(jù)基2時(shí)分FFT算法性質(zhì)采樣點(diǎn)數(shù)N應(yīng)滿足2的整數(shù)次冪，對(duì)式中的N個(gè)采樣點(diǎn)進(jìn)行FFT運(yùn)算，其變換結(jié)果見式（3）。

根據(jù)積化和差公式，可得式，令N=lfs，l為整數(shù)，可得式（4）：

則式（5）可化簡(jiǎn)為式（7）所示

根據(jù)離散正弦信號(hào)的周期性可知當(dāng)lf0+k≠0時(shí)必為0；當(dāng)lf0-k≠0時(shí)必為0。當(dāng)lf0+k=0時(shí)，必有f0=k=0；當(dāng)lf0-k=0時(shí)，有：

其中Ak是根據(jù)復(fù)數(shù)序列第k點(diǎn)所求得的幅值，V，采樣數(shù)據(jù)經(jīng)FFT變換后，生成的復(fù)數(shù)序列xk實(shí)部及虛部組成與被測(cè)信號(hào)相位的關(guān)系，見式：式中：θ為被測(cè)信號(hào)相位（°），信號(hào)頻譜分析中頻率值的計(jì)算與頻譜分辨率有關(guān)，頻譜分辨率是指頻譜分析時(shí)能夠分辨的兩個(gè)頻譜分量的最小間隔，頻譜分辨率F與采樣點(diǎn)數(shù)及采樣頻率有關(guān)，實(shí)際設(shè)計(jì)中應(yīng)根據(jù)被測(cè)信號(hào)頻率大小，設(shè)計(jì)采樣點(diǎn)數(shù)及采樣頻率，防止頻譜分辨率過大，以產(chǎn)生頻譜泄露現(xiàn)象，分辨率公式如式（11）所示，根據(jù)頻譜分辨率的概念，可推出頻率軸第n點(diǎn)所代表的頻率值為：

2 軟件設(shè)計(jì)

2.1 軟件界面和功能

本文基于Visual Studio 2015開發(fā)平臺(tái)C#語言對(duì)軟件進(jìn)行開發(fā)與設(shè)計(jì)，C#語言是一種安全、穩(wěn)定、簡(jiǎn)單易操作同時(shí)不失效率性和魯棒性的面向?qū)ο缶幊陶Z言[7-9]，軟件設(shè)計(jì)框圖及界面運(yùn)行圖如圖1和圖2所示。

圖1 設(shè)計(jì)框圖

圖2 運(yùn)行界面

軟件功能包含數(shù)據(jù)讀取、頻譜分析和報(bào)表生成，數(shù)據(jù)讀取功能支持讀入以TXT格式存儲(chǔ)的正弦信號(hào)數(shù)據(jù)。頻譜分析功能在精確測(cè)量被測(cè)信號(hào)頻率、幅值和相位的同時(shí)生成幅頻圖和相頻圖，多正弦信號(hào)疊加情況下測(cè)量的頻率、幅值和相位以空格分割顯示。報(bào)表生成功能支持將信號(hào)幅頻、相位值以及幅頻圖相位圖存入用戶指定文件中，軟件默認(rèn)以報(bào)表生成時(shí)間及當(dāng)前讀入數(shù)據(jù)點(diǎn)數(shù)及信號(hào)采樣頻率為文件命名。

2.2 數(shù)據(jù)文件讀入方法

軟件數(shù)據(jù)讀取展示及功能流程圖如圖3所示。

圖3 數(shù)據(jù)文件讀取流程圖

軟件利用C#OpenFileDialog類實(shí)現(xiàn)在用戶點(diǎn)擊界面數(shù)據(jù)導(dǎo)入按鈕時(shí)，軟件界面彈出文件選擇對(duì)話框，用戶選定數(shù)據(jù)文件后軟件逐行讀取所選文件內(nèi)容。采用VS中的Textbox控件顯示讀取數(shù)據(jù)的同時(shí)將數(shù)據(jù)存入到字符數(shù)組并計(jì)數(shù)，當(dāng)讀取的數(shù)據(jù)文件行內(nèi)容為空時(shí)停止讀取和計(jì)數(shù)，使用C#中的Convert類將字符型數(shù)組逐個(gè)轉(zhuǎn)換為浮點(diǎn)型數(shù)組，遍歷數(shù)組將數(shù)據(jù)點(diǎn)數(shù)作為chart圖表數(shù)據(jù)源繪出數(shù)據(jù)波形圖，完成數(shù)據(jù)讀取、繪圖和計(jì)數(shù)功能。

2.3 信號(hào)頻譜計(jì)算

根據(jù)奈奎斯特采樣定理，當(dāng)采樣頻率fs大于有限帶寬信號(hào)x(t)帶寬f0的兩倍時(shí)，信號(hào)可以由其采樣值重構(gòu)，否則信號(hào)會(huì)在x(n)中產(chǎn)生混疊，易得采樣數(shù)據(jù)x(n)能表示的頻譜最高模擬頻率為fs2 Hz，因此本文幅頻、相頻圖皆采用半邊頻譜。軟件頻譜計(jì)算整體流程圖如圖4所示。

圖4 頻譜計(jì)算整體流程圖

通過軟件界面選擇數(shù)據(jù)文件讀入采樣數(shù)據(jù)，設(shè)置采樣頻率等參數(shù)，繼而對(duì)讀入的采樣數(shù)據(jù)進(jìn)行DIT-FFT運(yùn)算得出復(fù)數(shù)數(shù)組，根據(jù)公式可得FFT算法變換后的復(fù)數(shù)序列（數(shù)組）xk與幅值序列（數(shù)組）Ak以及信號(hào)相位序列（數(shù)組）關(guān)系，利用公式（12）可得信號(hào)頻率，使用峰值算法提取幅值數(shù)組各個(gè)峰值，將高于閾值的峰值保留下來，通過C#Array.In?dexOf方法峰值所對(duì)應(yīng)的幅值數(shù)組序列號(hào)（相位數(shù)組序列號(hào)）提取出信號(hào)幅值和相位，利用式點(diǎn)數(shù)與頻率分辨率的關(guān)系得出信號(hào)分量頻率，為了消除軟件采用的浮點(diǎn)型數(shù)據(jù)所帶來的測(cè)量幅值的波動(dòng)性，利用用戶在軟件中輸入的信號(hào)測(cè)量閾值篩選掉干擾峰值，精確得出信號(hào)各分量幅頻及相位數(shù)值。

3 實(shí)驗(yàn)結(jié)果及分析

3.1 實(shí)驗(yàn)結(jié)果

實(shí)驗(yàn)選用MIT陣列感應(yīng)測(cè)井儀器[10]的26.325、52.65、105.3 kHz的3種頻率信號(hào)進(jìn)行軟件仿真驗(yàn)證。采用隨機(jī)高斯白噪聲，取信噪比變換范圍為5～30 dB，信噪比定義為信號(hào)功率與噪聲功率之比。取信號(hào)采樣頻率為5 391.36 kHz，采樣點(diǎn)數(shù)為2 048，所選參數(shù)滿足整周期采樣要求，測(cè)50次取其平均值作為最終檢測(cè)值，觀察整周期采樣情況下的測(cè)量精度，幅值及相位檢測(cè)相對(duì)誤差（表1、表2）。其次取采樣點(diǎn)數(shù)2 048，信號(hào)采樣頻率取3 080.777 kHz，測(cè)50次取其平均值，觀察非整周期采樣對(duì)檢測(cè)結(jié)果的影響，見表3及表4。

表1 整周期采樣檢測(cè)的幅值相對(duì)誤差

表2 整周期采樣檢測(cè)的相位相對(duì)誤差

表3 非整周期采樣檢測(cè)的幅值相對(duì)誤差

表4 非整周期采樣檢測(cè)的相位相對(duì)誤差

3.2 結(jié)果分析

由表1及表2可得出,在滿足整周期采樣時(shí)，在軟件信噪比低至5 dB的噪聲環(huán)境下，當(dāng)采樣頻率取5 391.36 kHz，采樣點(diǎn)數(shù)取2 048，滿足整周期采樣的情況下，信號(hào)幅值檢測(cè)相對(duì)誤差不超過0.31%且相位相對(duì)誤差不超過1.14%，驗(yàn)證了本文推導(dǎo)的FFT算法計(jì)算頻率、幅度、相位等參數(shù)的理論計(jì)算公式準(zhǔn)確性。通過表3及表4可得該軟件在信噪比低至5 dB噪聲背景下，當(dāng)信號(hào)采樣頻率取3 080.777 kHz，信號(hào)采樣點(diǎn)數(shù)不變的非整周期采樣情況下，幅值相對(duì)誤差不超過0.49%，相位相對(duì)誤差不超過1.39%，相較整周期采樣檢測(cè)精度有所下降，但檢測(cè)誤差仍保持在1.39%以下，這表明該軟件在噪聲環(huán)境下以及非整周期采樣等非理想情況下仍可實(shí)現(xiàn)信號(hào)頻譜的精確測(cè)量分析功能。

4 結(jié)論

針對(duì)感應(yīng)測(cè)井信號(hào)的處理與分析，設(shè)計(jì)實(shí)現(xiàn)了一種可量化分析信號(hào)幅頻相位參數(shù)值的頻譜分析軟件，引入DIT-FFT、尋峰算法及C#語言構(gòu)建軟件頻譜分析功能及軟件界面。通過在噪聲環(huán)境下以整周期及非整周期兩種采樣方式，對(duì)軟件測(cè)量性能進(jìn)行驗(yàn)證與分析。結(jié)果表明，該軟件在信噪比低至5 dB的高噪聲環(huán)境下，在整周期及非整周期兩種采樣方式下，軟件測(cè)量的幅值及相位誤差皆不超過1.39%，軟件檢測(cè)精度及穩(wěn)定性高，滿足實(shí)際工程使用要求，為測(cè)井信號(hào)數(shù)據(jù)分析提供了一種輔助手段。