全相位與Geiger算法的巖石聲發(fā)射源定位方法*

黃曉紅，孫國慶，張凱月

（1.華北理工大學(xué)信息工程學(xué)院，河北唐山063009；2.華北理工大學(xué)電氣工程學(xué)院，河北唐山063009）

全相位與Geiger算法的巖石聲發(fā)射源定位方法*

黃曉紅1*，孫國慶2，張凱月1

（1.華北理工大學(xué)信息工程學(xué)院，河北唐山063009；2.華北理工大學(xué)電氣工程學(xué)院，河北唐山063009）

Geiger定位算法在巖石定位分析中得到廣泛應(yīng)用，但是對(duì)初始值要求非常嚴(yán)格，若初始值選擇不當(dāng)，則很難進(jìn)入收斂范圍。為了提高定位精度，減少迭代次數(shù)，針對(duì)這些問題，提出全相位與Geiger算法的巖石聲發(fā)射源定位方法。首先進(jìn)行花崗巖（50 mm×100 mm×50 mm）斷鉛實(shí)驗(yàn)，接著充分利用全相位FFT相位不變性的優(yōu)勢(shì)分析斷鉛信號(hào)，通過相位差法求出時(shí)延同時(shí)反演聲發(fā)射源求出Geiger迭代初始值，最后綜合最小二乘法與Geiger算法，迭代求出最優(yōu)解。實(shí)驗(yàn)仿真結(jié)果顯示此算法的平均誤差相比美國PCI-2型聲發(fā)射儀器定位結(jié)果降低了約5 mm，有效解決了Geiger初始值選取的問題，能夠迅速進(jìn)入收斂范圍，提高了收斂速度和定位精度。

信號(hào)分析；聲發(fā)射源定位；全相位相位差；Geiger；最小二乘法；

在一般情況下，由于受力或其它作用下，材料產(chǎn)生瞬間彈性波的現(xiàn)象被稱為聲發(fā)射AE（Acoustic Emission）。

目前對(duì)于巖石類材料破壞研究采用的主要手段仍是聲發(fā)射技術(shù)。最常用的定位技術(shù)為時(shí)差定位，定位的方法主要有：相對(duì)定位法、Geiger定位方法、最小二乘法等。在20世紀(jì)60年代晚期，Mogi［1］最早在巖石彎曲變形條件下，開始了二維聲發(fā)射定位研究。之后Scholz［2］進(jìn)行了聲發(fā)射事件的三維空間定位研究。胡新亮［3］等利用相對(duì)定位法研究了非完整巖體破裂機(jī)理，使聲發(fā)射定位結(jié)果的可靠性與精度有所提高。趙興東［4］等應(yīng)用聲發(fā)射定位方法，研究了不同巖石破裂過程的AE活動(dòng)規(guī)律，應(yīng)用Geiger定位算法對(duì)完整巖石失穩(wěn)破壞過程進(jìn)行了研究，得到較高定位精度?？涤衩罚?］等針對(duì)三維定位中初始值選取影響定位精度的問題，即提出了將最小二乘法與Geiger算法結(jié)合的定位方法，有效提高了算法的收斂速度，但精度有所不足。張會(huì)新［6］等為了獲取傳感器節(jié)點(diǎn)的位置信息，采用TDOA定位算法，并在基礎(chǔ)上進(jìn)行了優(yōu)化，大大簡化了計(jì)算量。楊文鉑［7］等為了提高無線傳感器網(wǎng)絡(luò)的節(jié)點(diǎn)定位精度，通過研究信道噪聲，對(duì)測(cè)量誤差進(jìn)行了修正，提高了定位精度。黃曉紅［8］等針對(duì)巖石力學(xué)聲發(fā)射時(shí)差定位算法中時(shí)延估計(jì)精度的問題，提出了在全相位相位差法基礎(chǔ)上求時(shí)延估計(jì)的新方法，即利用全相位得到時(shí)間差，然后結(jié)合最小二乘法得到定位坐標(biāo)，此方法大大加快了收斂速度，但定位精度有所不足。為了得到更好的定位效果，考慮到Geiger的良好的定位特性，于是在此基礎(chǔ)上提出全相位與Geiger算法的巖石聲發(fā)射源定位方法。

1 全相位與Geiger算法的巖石聲發(fā)射源定位方法

1.1 基于全相位相位差時(shí)延估計(jì)法的巖石聲發(fā)射源定位研究

全相位相位差時(shí)延估計(jì)法是通過小波分析與全相位頻譜分析相結(jié)合求解相位差，從而得到聲發(fā)射信號(hào)的時(shí)延估計(jì)，進(jìn)而利用時(shí)差定位法反演聲發(fā)射源，實(shí)現(xiàn)精確的定位。其中全相位FFT頻譜［9-10］分析可測(cè)任意含有一個(gè)或多個(gè)頻率的正弦波信號(hào)的取樣點(diǎn)的相位，從而具有很大優(yōu)勢(shì)。通過全相位相位差法得到時(shí)間延遲，其推導(dǎo)過程如下：

設(shè)s1(k)是AE信號(hào)，n1(k)和n2(k)表示隨機(jī)噪聲，Δij為待估計(jì)的時(shí)間延遲，也就是信號(hào)到達(dá)第i個(gè)傳感器和第j個(gè)傳感器的傳播時(shí)間差；α為比例系數(shù)，則兩個(gè)空間分開的傳感器的接收信號(hào)的離散形式分別為：

x1(k)和x2(k)為離散譜線，利用全相位頻譜分析方法［5］求信號(hào)主譜線頻率為k0時(shí)的相位分別為：

時(shí)間延遲為：

式中，fc代表主譜線頻率，fs代表采樣頻率。

定位算法采用最小二乘法［11］，最小二乘法通過構(gòu)造殘差平方和，以模型參數(shù)估量為變量，當(dāng)函數(shù)值達(dá)到最小時(shí)，即可得到參數(shù)估計(jì)值的取值。其最大的優(yōu)勢(shì)為具有良好的估計(jì)特性。公式如式（6）：

式中，x，y，z為試驗(yàn)點(diǎn)坐標(biāo)；t為事件發(fā)生時(shí)刻；ti代表波到達(dá)第i個(gè)傳感器的時(shí)間；xi,yi,zi為第i個(gè)傳感器的位置；要得到最終定位結(jié)果，理論上至少需要四個(gè)傳感器，得到三組方程，進(jìn)行聯(lián)立求解。定位原理及實(shí)現(xiàn)同樣是先在各向同性的介質(zhì)中建立恰當(dāng)?shù)淖鴺?biāo)系，其三維算法模型如下：

分別是四個(gè)傳感器的坐標(biāo)，設(shè)G點(diǎn)為AE源，坐標(biāo)為G（x，y，z）。

式中，xi,yi,zi代表第i個(gè)傳感器的位置坐標(biāo)，v為波速，試驗(yàn)測(cè)得花崗巖樣品中的P波波速為3 128.023 5 m/s（波速獲得途徑通過美國物理聲學(xué)公司（PAC）PCI-2機(jī)器，具體操作方法參考其操作說明，相關(guān)聲速測(cè)試原理見文獻(xiàn)［12］）；Δij為AE信號(hào)從G點(diǎn)到達(dá)第i個(gè)傳感器和第j個(gè)傳感器的時(shí)間差，綜合方程組求解，方可求出AE源G點(diǎn)的坐標(biāo)。

綜合式（7）～式（9），即可求出聲發(fā)射源坐標(biāo)，我們把此結(jié)果作為下面Geiger迭代初始點(diǎn)。

1.2 基于全相位相位差的Geiger優(yōu)化迭代算法

Geiger定位方法［13］在地震學(xué)領(lǐng)域得到廣泛應(yīng)用，它是在給定初始點(diǎn)的基礎(chǔ)上，通過數(shù)次迭代來逼近最終結(jié)果。此處通過最小二乘法計(jì)算修正量，每次迭代都是以此修正量為基礎(chǔ)，把向量加到上次迭代結(jié)果上，就可以得到一個(gè)新的迭代點(diǎn)，判斷新點(diǎn)是否滿足要求，如不滿足，繼續(xù)迭代。每次迭代的結(jié)果由方程式（10）產(chǎn)生。

式中：x，y，z為人為設(shè)定的試驗(yàn)點(diǎn)初始值；t為事件發(fā)生時(shí)刻；ti代表波到達(dá)第i個(gè)傳感器的時(shí)間；xi,yi,zi為第i個(gè)傳感器的位置；v為波速。

設(shè)第i個(gè)傳感器測(cè)到的波到達(dá)時(shí)間為to,i，用試驗(yàn)點(diǎn)坐標(biāo)計(jì)算出的到達(dá)時(shí)間的一階泰勒展開式為：

式中：tc,i為由試驗(yàn)點(diǎn)坐標(biāo)計(jì)算出的波到達(dá)第i個(gè)傳感器的時(shí)間。

對(duì)于N個(gè)傳感器，矩陣形式表示為：

其中：

求解式（13）就得到修正向量：

之后，以(θ+Δθ)為試驗(yàn)點(diǎn)繼續(xù)迭代，直至滿足Δθ＜φ（這里Δθ為修正量，φ為一小量）。試驗(yàn)中設(shè)定本次迭代結(jié)果中Δt（式（13）中）值與上次Δt值的差值，即 φ(t)=Δtn-Δtn-1≤10-4，當(dāng) φ(t)＞10-4數(shù)量級(jí)時(shí)，迭代結(jié)束。由于Geiger初始值一般選取幾個(gè)傳感器包圍的中心點(diǎn)坐標(biāo)，這樣選取極容易受傳感器位置的影響，所以結(jié)合相位差法，選擇式（7）～式（9）中結(jié)果作為初始迭代值會(huì)使最終定位精度更精準(zhǔn)。黃曉紅等提出的基于相位差時(shí)延估計(jì)法是對(duì)同一主頻的子帶信號(hào)做全相位頻譜分析，從而得到每個(gè)子帶信號(hào)的中間點(diǎn)相位，再依據(jù)相位差［14］與時(shí)間差的關(guān)系，分別求出AE信號(hào)到達(dá)不同傳感器之間的時(shí)間差，如式（5），進(jìn)而結(jié)合最小二乘法算法求解，實(shí)現(xiàn)了聲發(fā)射源更精確的定位?；谌辔幌辔徊罘ǖ腉eiger優(yōu)化迭代算法則是利用全相位相位差法結(jié)合最小二乘法求解，從而得到最優(yōu)初始迭代點(diǎn)，然后利用Geiger迭代法迭代計(jì)算，將全相位相位差法和Geiger算法結(jié)合使用，能夠迅速進(jìn)入收斂范圍，僅需幾次迭代就能夠迅速收斂。此算法最大優(yōu)勢(shì)即省去尋找最優(yōu)步長的大量時(shí)間，加快求解速度，縮短了計(jì)算時(shí)間，提高了定位準(zhǔn)確度。

2 實(shí)驗(yàn)分析

試驗(yàn)樣品為花崗巖，試件尺寸加工為50 mm× 100 mm×50 mm，加工機(jī)器采用DJ-4型自動(dòng)巖石切割機(jī)。聲發(fā)射系統(tǒng)采用PCI-2型聲發(fā)射測(cè)試分析系統(tǒng)（美國物理聲學(xué)公司PAC）。傳感器型號(hào)為R6a型，諧振頻率是90 kHz，頻率范圍35 kHz～100 kHz。聲發(fā)射系統(tǒng)的門檻值設(shè)定為40 dB；前置放大器型號(hào)為2/4/6，其頻帶范圍為20 kHz～1 200 kHz；采樣頻率為2 MHz。為減小傳感器自身不可避免的誤差，耦合劑用凡士林，將探頭耦合在巖樣上。本試驗(yàn)中放置的8個(gè)傳感器的坐標(biāo)分別為A1（25，80，0）、A2（25，20，0）、A3（50，80，25）、A4（50，20，25）、A5（25，80，50）、A6（25，20，50）、A7（0，80，25）、A8（0，20，25），單位為mm。8個(gè)傳感器的坐標(biāo)如表1所示。

表1 傳感器布置坐標(biāo)

試驗(yàn)中，首先確定傳感器的安裝位置，然后以此點(diǎn)為圓心以傳感器為半徑，畫圓，布置傳感器使其和圓形區(qū)域高度重合，保證傳感器中心點(diǎn)和安放位置中心點(diǎn)吻合。試驗(yàn)共做了10個(gè)坐標(biāo)點(diǎn)（25，60，50）、（25，45，50）、（50，60，25）、（50，45，25）、（25，60，0）、（25，45，0）、（0，60，25）、（0，45，25）、（25，100，35）、（25，100，15）的斷鉛試驗(yàn)。對(duì)試驗(yàn)中8個(gè)傳感器獲取的AE信號(hào)分別進(jìn)行小波去噪處理和小波分解，比較分解后的子帶信號(hào)，找到每個(gè)傳感器中主頻最接近的4個(gè)子帶信號(hào)，然后進(jìn)行全相位頻譜分析，分別找到主頻對(duì)應(yīng)的相位（對(duì)于全相位頻譜分析而言，主頻位置對(duì)應(yīng)的相位即截取信號(hào)中間樣點(diǎn)的相位），利用相位差法進(jìn)行時(shí)間差測(cè)量，結(jié)合最小二乘法進(jìn)行定位，確定初始迭代試驗(yàn)點(diǎn)的坐標(biāo)，最后用Geiger算法迭代求解，求出最優(yōu)解。圖1列舉在坐標(biāo)點(diǎn)（50，45，25）處進(jìn)行斷鉛實(shí)驗(yàn)，8個(gè)傳感器采集到的信號(hào)和小波分解第5層細(xì)節(jié)信號(hào)頻譜圖。表2為8個(gè)傳感器接收到的信號(hào)以及小波分解后子帶的信號(hào)主頻。表3為美國生產(chǎn)PCI-2型聲發(fā)射儀器Geiger定位結(jié)果，表4為基于全相位相位差的Geiger優(yōu)化迭代算法的定位結(jié)果。

圖1 斷鉛實(shí)驗(yàn)8個(gè)傳感器接收到的信號(hào)及其小波分解信號(hào)

表2 8個(gè)傳感器接收到的信號(hào)及小波分解后子帶信號(hào)主頻 單位：kHz

表3 美國生產(chǎn)PCI-2型聲發(fā)射儀器Geiger定位結(jié)果

表4 基于全相位相位差法的Geiger定位結(jié)果

圖2為4種方法定位結(jié)果平均絕對(duì)誤差比較。

圖2 4種方法定位結(jié)果平均絕對(duì)誤差比較

通過比較，基于全相位相位差的Geiger優(yōu)化迭代算法的定位結(jié)果相比美國PCI-2型聲發(fā)射儀器的定位結(jié)果的平均誤差降低了約5 mm；相比全相位相位差法聲發(fā)射定位結(jié)果，平均誤差降低了約2 mm；相比美國PCI-2型聲發(fā)射儀器Geiger定位結(jié)果，平均誤差降低了約2 mm。

3 結(jié)論

①在小波分解基礎(chǔ)上，全相位頻譜分析與相位差法結(jié)合進(jìn)行AE信號(hào)時(shí)間差估計(jì)，有效地減小了誤差。

②與單獨(dú)使用最小二乘法確定初始迭代試驗(yàn)點(diǎn)相比，全相位相位差法能夠更快進(jìn)入收斂范圍，避免了對(duì)傳感器坐標(biāo)位置選取嚴(yán)格的缺點(diǎn)。

③全相位相位差法的Geiger優(yōu)化迭代定位算法定位精度比美國PCI-2型聲發(fā)射儀器的定位結(jié)果精度明顯提高了。從圖1中發(fā)現(xiàn)，在定位實(shí)驗(yàn)點(diǎn)最后兩位的效果不如前幾個(gè)點(diǎn)效果，出現(xiàn)原因?yàn)镚eiger定位對(duì)初始迭代值要求嚴(yán)格，所以結(jié)果好壞很大程度由全相位求初始值決定，所以第一步是關(guān)鍵，初始值選擇需進(jìn)一步優(yōu)化，信號(hào)分析能力需進(jìn)一步提高。特別指出各試驗(yàn)點(diǎn)相互獨(dú)立，所以個(gè)別點(diǎn)不代表整體效果。因?yàn)槊總€(gè)試驗(yàn)點(diǎn)的傳感器接收的信號(hào)不同，所以每個(gè)點(diǎn)的的初始迭代值是不同的，為了避免出現(xiàn)不理想的情況，在小波分解細(xì)節(jié)信號(hào)層數(shù)時(shí)需要選出效果最好的一層，所以為了得到好的效果，需要盡可能多的層數(shù)分解，通過一個(gè)點(diǎn)的實(shí)驗(yàn)比較，找出效果最好的一層，再往下進(jìn)行別的點(diǎn)的定位實(shí)驗(yàn)，以達(dá)到最優(yōu)定位效果。

④基于全相位相位差的Geiger優(yōu)化迭代定位新算法為巖石聲發(fā)射定位研究提供了一種途徑。

［1］Mogi K.Study on Elastic Shocks Caused by the Fracture of Hetero?geneous Material and Its Relation to Earthquake Phenomena［J］.Bull of the Earthquake Research Institute，1962，40（6）：831-853.

［2］Scholz C.Experimental Study of the Fracturing Process in Brittle Rock［J］.Journal of Geophysical Research，1968，73（5）：1447-1454.

［3］胡新亮，馬勝利，高景春，等.相對(duì)定位方法在非完整巖體聲發(fā)射定位中的應(yīng)用［J］.巖石力學(xué)與工程學(xué)報(bào)，2004，23（2）：277-283.

［4］趙興東，劉建坡，李元輝，等.巖石聲發(fā)射定位技術(shù)及其實(shí)驗(yàn)驗(yàn)證［J］.巖土工程學(xué)報(bào)，2008，30（10）：1472-1476.

［5］康玉梅，劉建坡，李海濱，等.一類基于最小二乘法的聲發(fā)射源組合定位算法［J］.東北大學(xué)學(xué)報(bào)，2010，31（11）：1648-1651.

［6］張會(huì)新，陳德沅，彭晴晴，等.一種改進(jìn)的TDOA無線傳感器網(wǎng)絡(luò)節(jié)點(diǎn)定位算法［J］.傳感技術(shù)學(xué)報(bào)，2015，28（3）：412-415.

［7］楊文鉑，邢鵬康，劉彥華.一種基于自適應(yīng)RSSI測(cè)距模型的無線傳感器網(wǎng)絡(luò)定位方法［J］.傳感技術(shù)學(xué)報(bào)，2015，28（1）：137-141.

［8］黃曉紅，張艷博，田寶柱，等.基于相位差時(shí)延估計(jì)法的巖石聲發(fā)射源定位研究［J］.巖土力學(xué)學(xué)報(bào)，2015，36（2）：382-383.

［9］盧新寧.基于全相位FFT的鐵道信號(hào)頻率檢測(cè)算法研究［D］.海南：海南大學(xué)，2013.

［10］龐前娟，孫山林，魏麗娜.基于全相位FFT的故障信號(hào)檢測(cè)方法［J］.制造業(yè)自動(dòng)化，2010，32（7）：96-99.

［11］馬學(xué)敏.巖石聲發(fā)射信號(hào)定位研究［D］.河北：河北聯(lián)合大學(xué)，2014.

［12］馬榮田.巖體聲波測(cè)試技術(shù)及其在巷道圍巖巖體波速測(cè)試中的應(yīng)用［J］.鐵道建筑，2006，10：37-40.

［13］劉建坡，王洪勇，楊宇江，等.不同巖石聲發(fā)射定位算法及其實(shí)驗(yàn)研究［J］.東北大學(xué)學(xué)報(bào)，2009，30（8）：1194-1195.

［14］黃曉紅，鄧振淼.改進(jìn)的相位展開算法及其在瞬時(shí)頻率估計(jì)中的應(yīng)用［J］.電子學(xué)報(bào)，2009，37（10）：224-230.

黃曉紅（1973-），女，漢族，河北省唐山人，教授，博士，華北理工大學(xué)信息工程學(xué)院碩士生導(dǎo)師，主要從事信號(hào)與信息處理方面的教學(xué)和研究工作，tshxh@163.com；

孫國慶（1988-），男，漢族，河北省張家口人，華北理工大學(xué)電氣工程學(xué)院碩士，主要從事信號(hào)與信息處理方面研究工作，457297857@qq.com。

張凱月（1992-），女，漢族，河北省滄州人，華北理工大學(xué)信息工程學(xué)院碩士，主要從事信號(hào)與信息處理方面研究工作，10322872757@qq.com。

Rock Acoustic Emission Source Localization Method Based on All Phase and Geiger Algorithm*

HUANG Xiaohong1*，SUN Guoqing2，ZHANG Kaiyue1
（1.College of information engineering，North China University of science and technology，Tangshan Hebei 063009，China；2.College of electrical engineering，North China University of science and technology，Tangshan Hebei 063009，China）

Geiger localization algorithm is widely used in the analysis of rock localization，however，it is very strict on the initial value，if the initial value is not properly selected，it is difficult to enter the convergence range.In order to improve the positioning precision and reduce the number of iterations，rock acoustic emission source localization method based on all phase and Geiger algorithm is put forward.First of all，the experiment was carried out on the granite（50 mm×100 mm×50 mm）.Secondly，making full use of the advantage of apfft phase invariance to analyse the lead signal，and the delay time is obtained by the phase difference method，at the same time，the initial value of the Geiger iteration is obtained by the inversion of the acoustic emission source.In the end，combining least square method and Geiger algorithm to obtain optimum solution.Experimental simulation results are shown：Comparing the positioning results of the United States PCI-2 type acoustic emission instrument，the average error decreased by about 5 mm.This method can effectively solve the problem of selecting the initial value of Geiger，so that it can quickly enter the convergence range，improving the convergence speed and positioning accuracy.

signal analysis；acoustic emission source location；all-phase phase difference；Geiger；least square method EEACC：7230；7310N

10.3969/j.issn.1004-1699.2016.11.013

TU45

A

1004-1699（2016）11-1705-06

項(xiàng)目來源：國家自然科學(xué)基金面上項(xiàng)目（51174071；51374088）；河北省鋼鐵聯(lián)合基金項(xiàng)目（E2016209357）；華北理工大學(xué)培育基金項(xiàng)目（GP201506）

2016-05-14 修改日期：2016-06-27