基于非均勻采樣的?-OTDR仿真與信號處理?

（海軍工程大學電子工程學院 武漢 430000）

1 引言

相位敏感光時域反射計（?-OTDR）是一種分布式的光纖傳感技術(shù)。由于其抗電磁干擾、靈敏度高、定位精確、數(shù)據(jù)處理相對簡單等特點，在大型建筑的健康監(jiān)測，海底管線的安全監(jiān)控，國家邊界安全等領(lǐng)域得到了廣泛的應(yīng)用［1］。但是傳統(tǒng)的?-OTDR系統(tǒng)由于受到光波在光纖內(nèi)往返時間的制約，根據(jù)奈奎斯特采樣定理，其最高響應(yīng)頻率只能達到探測脈沖頻率的一半。特別是對于一些長距離的?-OTDR傳感系統(tǒng)，探測頻率的響應(yīng)范圍會出現(xiàn)十分顯著的下降。文獻［2］中的?-OTDR系統(tǒng)傳感距離達到了175km，但是相應(yīng)的最高響應(yīng)頻率只有285Hz。針對這一問題，2017年Zhang提出可以利用非均勻采樣方法，通過隨機調(diào)制探測脈沖的時間間隔，突破奈奎斯特采樣定理的限制，提高?-OTDR系統(tǒng)的最高響應(yīng)頻率。在最高采樣頻率為10kHz的系統(tǒng)中，可以探測到1.153MHz的單頻振動［3］。在文獻［4］中，Zhang進一步對這一方法中各參數(shù)的最優(yōu)取值進行了分析。但是，他僅僅分析了非均勻采樣法如何探測振動頻率，對這一方法如何用于定位并沒有研究。而傳統(tǒng)的?-OTDR系統(tǒng)信號處理技術(shù)，如移動平均并對平均結(jié)果做差分的方法（以下簡稱移動平均法），主要是應(yīng)用于均勻采樣信號的定位［5］，對非均勻采樣信號效果并不理想。本文首先分析了?-OTDR系統(tǒng)探測原理，并借鑒光柵設(shè)計中的離散層剝離思想，對?-OTDR系統(tǒng)探測過程中產(chǎn)生的瑞利散射光信號進行了仿真模擬。其次，針對基于非均勻采樣方法的?-OT?DR系統(tǒng)特性，提出了一種利用頻域信號標準差的改變值來定位振源位置的方法，并對這一方法的抗噪性能，以及振動信號的頻域稀疏度等進行了分析。

2 ?-OTDR系統(tǒng)探測原理

2.1 系統(tǒng)構(gòu)成

?-OTDR的系統(tǒng)框圖如圖1所示，光源發(fā)出穩(wěn)定的超窄線寬激光，經(jīng)由聲光調(diào)制器（AOM）調(diào)制為周期性的光脈沖，在通過摻鉺光纖放大器（ED?FA）并濾除噪聲后經(jīng)環(huán)形器注入到光纖鏈路。光脈沖在傳播過程中，會在光纖沿線不斷產(chǎn)生后向瑞利散射光，并且由于瑞利散射是彈性散射，產(chǎn)生的散射光會發(fā)生多光束干涉，返回的干涉光進入光電探測器（APD）變?yōu)殡娦盘?，然后?jīng)數(shù)字采集器（DAQ）轉(zhuǎn)換為數(shù)字信號，這樣終端處就得到了一條沿時間軸分布的瑞利散射曲線，并且時間點與光纖距離對應(yīng)。當光纖上某點發(fā)生振動時，相應(yīng)位置處瑞利散射光相位會發(fā)生改變，進而造成接收端瑞利散射曲線的變化，通過探測這一變化，就可以定位振動源并求得振動頻率。

圖1 ?-OTDR系統(tǒng)框圖

2.2 探測原理

?-OTDR主要是通過探測光纖中返回到入射端的瑞利散射光的干涉圖樣，從而實現(xiàn)分布式光纖傳感。所謂瑞利散射光，就是指遠小于光波長的粒子所引起的光彈性散射（散射光與入射光頻率相同）。具體到光纖中，則是指光纖中遠小于入射光波長的不均勻性所引起的光彈性散射。這種不均勻性是在光纖拉制過程中固化在光纖中的，可以等效看作是光纖內(nèi)折射率的不均勻起伏［6］。

由于?-OTDR使用的是頻率穩(wěn)定的超窄線寬光源，而且探測脈沖具有一定的寬度，所以在脈沖內(nèi)部不同位置處會產(chǎn)生眾多的瑞利散射光，而且不同時刻產(chǎn)生的散射光在返回時會疊加在一塊，所以接收端接收到的散射光信號實際是多束瑞利散射光干涉的結(jié)果，如圖2所示。

圖2 光纖內(nèi)相干瑞利散射光的生成原理圖

在圖2（a）中，矩形光脈沖在t0時刻前端運動到了a0點，此時aN點是光脈沖的中點。這個脈沖所覆蓋的光纖可以劃分為2N個極窄的光纖微元（寬度為ΔL），設(shè)第i個微元產(chǎn)生的后向瑞利散射光可表示為（圖2中用黑色箭頭表示），則t0時刻一共產(chǎn)生了2N個。設(shè)tN時刻光脈沖的末端運動到aN點（圖2（b）），此時，第2N個微元中除了tN時刻產(chǎn)生的，還包含之前時刻產(chǎn)生的，所以綜上所述，在tN時刻第N個光纖微元中總的后向瑞利散射光Er可以表示為

隨著光脈沖的前進，光纖上的每個ΔL內(nèi)都會產(chǎn)生一個Er，而?-OTDR正是通過探測并記錄各時刻返回的Er，從而實現(xiàn)分布式的光纖傳感。

3 ?-OTDR仿真

3.1 仿真模型設(shè)計

1998年P(guān)ark等將?-OTDR的基本理論模型抽象為將光纖分成N段，在每小段光纖中存在若干個反射率和相位變化隨機取值的反射鏡，用反射鏡的反射光來等效后向瑞利散射光。這一模型很好地解釋了?-OTDR的物理規(guī)律，且符合實驗結(jié)果［7］。但是這一模型沒有將振動造成的折射率波動與散射光的相位變化情況很好地聯(lián)系起來，不利于計算機仿真的實現(xiàn)。因此，本文提出將每小段光纖中的反射鏡替換成一段布拉格光柵，利用光柵設(shè)計的相關(guān)公式，在光纖折射率與瑞利散射光相位之間建立聯(lián)系。

已知均勻布拉格光柵的反射光波長λr=2n0Λ，其中n0是光纖群折射率，Λ是柵距。所以只要控制Λ就可以得到和入射光波長相同的反射光。假設(shè)反射光的頻譜是超窄線寬的，并且滿足光柵相位匹配條件，基于模式耦合理論，可以得到光柵的反射率R和反射光的相位改變值φr［8］：

其中，lF是光柵長度，C是耦合系數(shù)，且C與光纖折射率的波動幅度Δn有如下關(guān)系：

這樣一來，我們只要知道了光纖上每一點的折射率起伏值Δn，就可以算出相應(yīng)位置處后向瑞利散射光的反射率R和相位改變值φr。由式（1）可知，后向瑞利散射光的相干只發(fā)生在長度為NΔL的光纖段內(nèi)，所以借鑒離散層剝離算法的思想［9］，我們可以把對瑞利散射光的求解控制在NΔL層光纖內(nèi)。隨著光脈沖的前進，不斷迭代求得整條光纖上的后向瑞利散射光分布。

3.2 數(shù)值仿真

在光波長為1550nm時，單模光纖中的瑞利散射損耗系數(shù)約為0.16×10-3dB/m，其中后向瑞利散射光占比約為Br≈ 0.082%［10］，所以 R 的均值。已知，只要合理設(shè)置lF的值，通過式（2）、（4）即可求得光纖上折射率的波動幅度均值。假設(shè)光纖中的不均勻性符合瑞利分布，可以利用Matlab生成期望是，滿足瑞利分布的Δn（z）數(shù)組，為每個光纖微元賦值，再代入式（4）求出C（z），進而求得每個光纖微元內(nèi)的R與φr。已知光纖各位置處的R與φr，則只要預設(shè)好入射光脈沖的初始參數(shù)，利用式（1），就可以得到整條光纖上后向瑞利散射光的干涉曲線。

設(shè)仿真數(shù)值如表1所示。

利用表1參數(shù)，仿真得到?-OTDR在一個脈沖周期內(nèi)產(chǎn)生的瑞利散射光信號。如圖3，這是長度為3km和30km光纖上的后向瑞利散射光功率分布曲線。從圖3（b）中可以看到，隨著光纖距離的增長，瑞利散射光強度迅速衰減。因為單模光纖中的后向瑞利散射光信號極其微弱，遠端光纖鏈路上的后向散射光信號通常會淹沒在噪聲中。因此必須使用有效的信號處理手段，從噪聲中提取信號，以實現(xiàn)對光纖鏈路的分布式監(jiān)測。

表1 仿真數(shù)值與說明

圖3 瑞利散射光功率分布曲線

4 信號處理

4.1 非均勻采樣原理

設(shè)?-OTDR的傳感光纖長度為3km，則光波在光纖內(nèi)往返需要的時間約為0.03ms。所以探測光脈沖的間隔時間需大于0.03ms，若設(shè)系統(tǒng)探測間隔為0.1ms，則相應(yīng)的采樣頻率為10kHz。根據(jù)奈奎斯特采樣定理，當系統(tǒng)等間隔采樣時，其最高響應(yīng)頻率為采樣頻率的一半，即5kHz。當振動信號含有大于5kHz的頻率分量時，振動信號的高頻分量與低頻分量會在頻域發(fā)生混疊。如果對振動信號做離散傅里葉變換（DFT），會在頻域產(chǎn)生偽值。為直觀表述，這里設(shè)存在振動信號s（t），其表達式為

其中f1=3.2kHz，f2=22kHz，f3=68kHz，f4=94kHz，信號圖像如圖4（a）。若采用均勻采樣且采樣頻率為10kHz，則可得采樣后的信號如圖4（b）。對采樣后的信號做DFT可得頻譜圖如圖5所示?？梢园l(fā)現(xiàn)，圖中除了3.2kHz是真實存在的頻譜分量，2kHz和4kHz處的頻譜分量都是由于頻率混疊而出現(xiàn)的偽值。

圖4 信號時域曲線

圖5 均勻采樣信號做DFT得到的頻譜圖

當振動信號在頻域具有稀疏特性時，可以采取非均勻采樣即隨機采樣的方法，突破采樣定理的限制。此時采樣間隔不再是固定值，而是一個時間范圍。依然以式（5）中的信號s(t)為例，若采樣間隔在［0.1ms，0.2ms］之間服從均勻分布，可得采樣后的信號如圖6所示。使用非均勻離散傅里葉變換（NDFT）求解其頻域信號：

其中，x(tn)代表采樣后的離散時域信號，N是總的采樣點數(shù)。另外，NDFT方法的頻域取值間隔fa與探測帶寬B受采樣點數(shù)N的限制，三者的最優(yōu)關(guān)系式為［4］

其中，Nmin為最少采樣點數(shù)。由式（7）可以看出，若想提高頻率分辨率（即減小fa）或增大探測帶寬B，必須要增大Nmin才能做到。但是Nmin增大也會相應(yīng)增加系統(tǒng)的探測時延與運算量。若設(shè)B=100kHz，fa=100Hz，則采樣點數(shù)至少需要N=250。圖7（a）是當N=250時，隨機采樣信號經(jīng)過NDFT后所得的頻譜圖。可以看到，圖中四個較大的頻譜分量正是原始信號s（t）中的四個頻率。但是除了這四個頻率值以外，還存在大量的干擾值，這是由于原始信號在經(jīng)過NDFT時會發(fā)生頻率泄露。由于頻率泄露是均勻分布在整個頻域上的，因而干擾值都比較小。且干擾值可以通過增加采樣點數(shù)來進一步抑制。如圖7（b），這是N=1000時的頻譜圖，可以看到干擾值顯著減小了。

圖6 隨機采樣得到的信號

圖7 隨機采樣信號頻譜圖

4.2 信號處理與分析

為進一步研究隨機采樣方法應(yīng)用到?-OTDR系統(tǒng)時的情況，設(shè)光纖長度為3km，在2km位置處（取前后各5m共10m長光纖）加入一個頻率與式（5）中s（t）相同，且折射率最大改變值為 6.5×10-5的振動。探測光脈沖間隔時間在［0.1ms，0.2ms］區(qū)間上服從均勻分布，仿真產(chǎn)生500條瑞利散射曲線，如圖8所示。

圖8 500條隨機采樣得到的瑞利散射曲線的疊加圖像。

左邊的圖8（a）是無噪聲條件下，振源附近500條瑞利散射曲線的疊加圖，可以在2km處觀察到振動造成的曲線振幅的改變。但是，實際的?-OTDR系統(tǒng)中存在著背景噪聲，而且噪聲會隨著時間累積變得越來越嚴重［11］。這里定義仿真信噪比（SSNR）為散射曲線均值與噪聲功率均值的比值，在圖8（a）曲線中加入SSNR=10dB的噪聲，如圖8（b）?？梢钥吹皆诩尤朐肼暫螅呀?jīng)無法分辨出振動點。

利用這500條瑞利散射曲線（SSNR=10dB），我們可以在光纖鏈路的每一點上都得到一個N=500的隨機采樣信號，代入式（6）即可求出各個位置處的頻譜圖。已知噪聲信號在頻域不具有稀疏性，如果僅僅是去掉幾個比較大的頻率分量，并不會改變信號的標準差。但是振動信號在頻域具有稀疏特性，即除了幾個比較大的頻率分量外，其余頻率分量都比較小。所以如果去除掉這幾個比較大的頻率分量，則信號的標準差會發(fā)生較大的改變。利用這一變化，就可以定位振動源的位置。在這里我們提出一種閾值法來得到各位置處的較大的頻率分量［12］，閾值定義如下：

其中，α是頻域信號的均值，Nf是頻域采樣點數(shù)。若所求頻域帶寬為100kHz，頻率分辨率取fa=100Hz，則Nf=1000。利用閾值T，我們可以將頻域信號中大于T的頻率分量數(shù)目設(shè)為KT，KT?Nf而且各位置處KT值會隨信號頻率大小自適應(yīng)調(diào)整，如圖9（a）。

圖9 （a）光纖各位置處KT的取值；（b）利用信號標準差改變值得到的定位圖像；（c）信號經(jīng)移動平均法處理后的定位圖像；（d）兩種方法的信噪比對比圖

將各位置處的KT個較大的頻率分量去除，并計算去除前后該位置處頻域信號標準差的改變值，將求得的改變值歸一化后繪成曲線如圖9（b）所示。可以看到，在2km位置處存在峰值。

將2km位置處的峰值與其他位置處的最大值之比定義為定位信噪比（LSNR），則可得圖9（b）中的LSNR=14.3dB。這一結(jié)果是在SSNR=10dB的情況下得到的。將這一結(jié)果與傳統(tǒng)?-OTDR系統(tǒng)的信號處理方法進行對比，繪制出相同噪聲條件下，傳統(tǒng)的移動平均法［5］得到的定位圖像，如圖9（c），圖中的LSNR=3.4dB。對比兩種方法的抗噪性能，得到它們的LSNR與SSNR的關(guān)系圖，如圖9（d）。圖中散點是五次結(jié)果取平均后得到的，相應(yīng)的擬合曲線為5階擬合曲線。其中實心黑點是利用標準差改變值定位時得到的LSNR值，而空心圓圈則是利用移動平均法定位時得到的LSNR值。由圖9（d）可得，本文所提定位方法應(yīng)用于隨機采樣信號時優(yōu)于傳統(tǒng)的移動平均法。而且當SSNR=0時，這一方法的LSNR值依然大于0，說明這一方法具有較好的抗噪性能。

定位振源位置之后，就可以求出相應(yīng)的振動頻率。為了求解準確，將圖9（b）中峰值處求得的頻域信號相加取平均，并利用式（8）找出其中大于閾值的頻率分量，可得振動信號頻譜圖如圖10（a）所示，其中明顯突出的四個頻率分量就是我們要求的頻率分量。

圖10 （a）振動位置處的信號頻譜圖；（b）FSNR與信號頻域稀疏度K的關(guān)系圖

最后需要指出，隨機采樣求振動信號的方法也有其局限性，那就是只能求頻域符合稀疏特性的振動信號。設(shè)振動信號的頻域峰值與干擾值峰值之比為頻率信噪比（FSNR），且設(shè)振動信號的頻域稀疏度為K。FSNR與K之間的關(guān)系如圖10（b）所示，圖中散點是十次結(jié)果求平均得到的均值，曲線是對應(yīng)散點的5階擬合曲線。由圖可知，若想得到FS?NR較高的頻譜圖，K應(yīng)該越小越好。當振動信號中的頻率分量大于10時，信號處理結(jié)果存在較大干擾。

5 結(jié)語

基于非均勻采樣方法的?-OTDR系統(tǒng)，可突破奈奎斯特采樣定理對系統(tǒng)最高響應(yīng)頻率的限制，但是由于這一方法需要隨機調(diào)制系統(tǒng)探測脈沖的間隔時間，所以傳統(tǒng)的信號處理方法，如移動平均法應(yīng)用于此處效果并不理想。針對這一問題，本文提出了一種利用離散頻域信號標準差的改變值來定位振源的方法，通過設(shè)計仿真模型，對其抗噪性、信號頻域稀疏度等方面進行了研究，并與傳統(tǒng)移動平均法進行了對比。從仿真結(jié)果來看，對于具有頻域稀疏特性的振動信號，本文方法的定位效果與抗噪性能都要優(yōu)于傳統(tǒng)的移動平均信號處理方法。這對于基于非均勻采樣方法的?-OTDR系統(tǒng)的進一步優(yōu)化與應(yīng)用有重要意義。