本發(fā)明涉及三維空間聲源識別領域
背景技術：
：憑借旋轉對稱性好、聲場記錄全面、衍射作用強、記錄信號信噪比高等優(yōu)勢，實心球形傳聲器陣列廣泛流行于三維聲源識別領域。球諧函數波束形成(SphericalHarmonicsBeamforming,SHB)是其普遍采用的聲信號后處理方法。但是受球諧函數階次截斷、傳聲器離散采樣等因素影響，該方法輸出的聲源成像圖被大量高水平旁瓣污染，同時空間分辨率低，這些缺陷會使聲源識別結果的分析具有嚴重不確定性。目前三維空間已有的反卷積(DAMAS)方法在SHB方法輸出結果、陣列點傳播函數、聲源分布之間建立線性方程組，通過在反復迭代過程中引入正約束來定解該方程組，從而提取真實聲源分布，可獲得超高空間分辨率，但存在嚴重耗時問題；實心球的函數型延遲求和方法(FDAS)在實心球延遲求和(DAS)方法的基礎上發(fā)展而來，能快速衰減DAS方法產生的旁瓣并能準確探測弱源，但該方法僅能輕微改善DAS方法空間分辨率，無法準確分辨彼此靠近的聲源，并且實際應用中對聲源的量化偏差也較大。技術實現(xiàn)要素：本發(fā)明的目的是針對實心球陣列三維聲源識別低旁瓣、超高分辨率聲學圖像的快速獲取問題提供一種方法，同時克服FDAS方法對聲源量化偏差大的不足。為實現(xiàn)本發(fā)明目的而采用的技術方案如下，實心球聲源識別低旁瓣超高分辨率聲學圖像快速獲取方法：1)：構建測量系統(tǒng)：參見圖1，一個處于三維空間內的實心球陣列；所述實心球陣列包括一個半徑為a的實心球，以及分布在所述實心球表面的Q個傳聲器；這些傳聲器的編號為q，q＝1,2,…,Q；所述實心球的球心作為原點，實心球陣列所處的三維空間內任意一點的位置坐標用(r,Ω)描述，r表示所描述位置與原點間的距離，表示所描述位置的方向，θ、分別為其仰角、方位角；(a,Ωq)為q號傳聲器的位置坐標；2)獲取聲壓各個傳聲器接收聲壓信號p(ka,Ωq)；其中，處于所述三維空間內的單極子點聲源，輻射聲波為球面波；波數k＝2πf/c，聲波頻率為f，聲速為c；3)通過各個傳聲器接收聲壓信號P，得到互譜矩陣C：C=ppH‾,]]>P為各個傳聲器接收的復聲壓信號，p=Σ(r0,Ω0)∈Bs(kr0,Ω0)t(kr0,Ω0),]]>B表示所有聲源位置坐標組成的集合，(r0,Ω0)和(r0',Ω0')均表示聲源的位置坐標，s(kr0,Ω0)表示聲源的強度，t(kr0,Ω0)＝[t1(kr0,Ω0),t2(kr0,Ω0),…,tQ(kr0,Ω0)]T，tq(kr0,Ω0)表示聲源到q號傳聲器的聲場傳遞函數，tq(kr0,Ω0)≡Σn=0∞Σm=-nnRn(kr0,ka)Ynm*(Ω0)Ynm(Ωq),]]>分別為Ω0、Ωq方向的球諧函數，n、m為球諧函數的階次，Rn(kr0,ka)為聲信號傳播徑向函數，Rn(kr0,ka)=-4πihn(2)(kr0)(jn(ka)-jn′(ka)hn(2)′(ka)hn(2)(ka)),]]>其中，被減數項為球面波在自由場中傳播的徑向函數，減數項為實心球引起的球面波散射徑向函數，為虛數單位，jn(ka)為n階第一類球貝塞爾函數，為n階第二類球漢克爾函數，j'n(ka)、分別為jn(ka)、的一階導數。上標“-”、“H”、“*”和“T”分別表示平均運算、轉置共軛運算、共軛運算和轉置運算；4)確定各個頻率對應的球諧函數階次截斷長度：表示將數值向正無窮方向圓整到最近的整數；5)將互譜矩陣C和球諧函數階次截斷長度N代入FDAS方法輸出式：WF(krf,Ωf)=1Q(vH(krf,Ωf)C1ξv(krf,Ωf)||v(krf,Ωf)||22)ξ---(I)]]>其中：(rf,Ωf)為聚焦點的位置坐標，v(krf,Ωf)＝[v1(krf,Ωf),v2(krf,Ωf),…,vQ(krf,Ωf)]T為聚焦列向量，vq(krf,Ωf)表示q號傳聲器的聚焦分量，Rn(krf,ka)為聚焦徑向函數。||||2表示2范數，指數參數ξ推薦取值為16。6)對(Ⅰ)構造脊探測(RD)對象，L(Ωf)＝|WF(krf,Ωf)|“||”表示取模運算，對每個聚焦方向，可按下列步驟檢測其是否為極大值脊：步驟1.構造L(Ωf)的二階偏導數矩陣，即黑塞矩陣：H=∂2L∂φf2∂2L∂φf∂θf∂2L∂θf∂φf∂2L∂θf2]]>步驟2.特征值分解矩陣Η，找出絕對值最大的特征值λ及相應的特征向量uλ，其中uλ指示最大表面曲率方向，輸出值在uλ方向上若為極大值則λ為負，若為極小值則λ為正。步驟3.計算標量ρ(Ωf)為：ρ(Ωf)=-12sign(λ)|sign(▿L(Ωf+χuλ)·uλ)-sign(▿L(Ωf-χuλ)·uλ)|]]>其中，“sign”為符號函數，為梯度算子，“．”表示求向量內積，χ為常數，用以控制檢測空間精度，本發(fā)明取為0.2。當ρ(Ωf)＝1時，聚焦方向為極大值脊，保留FDAS的輸出值；反之，用零取代FDAS在該聚焦方向的輸出值來剔除非脊點。7)對DAS方法結果進行反卷積(DAMAS)求取各聲源的平均聲壓貢獻：定義聲源平均聲壓貢獻PAC(kr0,Ω0)為聲源在陣列各傳聲器處產生聲壓信號自譜的均值，即為互譜矩陣C的跡，并帶入DAS輸出量表達式，W(krf,Ωf)=Σ(r0,Ω0)∈BPAC(kr0,Ω0)vH(krf,Ωf)t(kr0,Ω0)tH(kr0,Ω0)v(krf,Ωf)||v(krf,Ωf)||22||t(kr0,Ω0)||22]]>定義DAS的點傳播函數為：psf((krf,Ωf)|(kr0,Ω0))=vH(krf,Ωf)t(kr0,Ω0)tH(kr0,Ω0)v(krf,Ωf)||v(krf,Ωf)||22||t(kr0,Ω0)||22]]>其表示(r0,Ω0)位置處的單位平均聲壓貢獻點聲源在聚焦點(rf,Ωf)處的波束形成貢獻量。DAS輸出量為各聲源平均聲壓貢獻與對應點傳播函數的乘積的和。其中，t(kr0,Ω0)＝[t1(kr0,Ω0),t2(kr0,Ω0),…,tQ(kr0,Ω0)]T表示聲場傳遞函數列向量。6)中RD方法獲取的脊中包含聲源點，且聲源點處的輸出值相對很大，可默認一定動態(tài)范圍外的脊和非脊處無聲源，輸出值為零，僅對動態(tài)范圍內的脊進行DAMAS運算來獲得超高空間分辨率，此時G為參與運算的脊的數目，該數值很小，可使反卷積具有極高的效率。W為G×1維列向量，對應元素為參與運算各脊處的DAS輸出量。A為G×G維的點傳播函數矩陣，PAC為G×1維列向量，對應元素為參與運算各脊處聲源的平均聲壓貢獻，為未知非負數。如下線性方程組成立：W＝APAC該線性方程組采用高斯-賽德爾迭代方案求解PAC，移除點傳播函數的影響，能夠顯著衰減旁瓣并獲得超高空間分辨率。初始化基于第l次迭代計算的結果進行第l+1次迭代的具體步驟為：步驟1.計算殘差re：re(l)(krf,Ωf)=Σ(r0,Ω0)∈Dpsf((krf,Ωf)|(kr0,Ω0))PAC(l+1)(kr0,Ω0)+Σ(r0,Ω0)∈Epsf((krf,Ωf)|(kr0,Ω0))PAC(l)(kr0,Ω0)-W(krf,Ωf)]]>其中，D為已完成l+1次迭代計算的聲源點組成的集合，E為未完成l+1次迭代計算的脊點組成的集合，二者的并集包含所有參與運算脊點聲源。步驟2.確定PAC(l+1)(kr0,Ω0)=max(PAC(l)(kr0,Ω0)-re(l)(krf,Ωf)psf((krf,Ωf)|(kr0,Ω0)),0)rf=r0,Ωf=Ω0]]>按照該方案依次計算每個運算脊點，即可得出8)根據迭代所得PAC即能定位聲源位置并且能繪制聲源成像圖實現(xiàn)三維聲場可視化。本發(fā)明的方法推導過程如下：步驟1：獲取各傳聲器接收復聲壓信號的互譜矩陣圖1為球陣列波束形成坐標系，原點位于陣列中心，三維空間內任意位置可用(r,Ω)描述，r表示所描述位置與原點間的距離，表示所描述位置的方向，θ、分別為其仰角、方位角。符號“●”代表嵌在實心球表面的傳聲器，記a為球半徑，Q為傳聲器總數，q＝1,2,…,Q為傳聲器索引，(a,Ωq)為q號傳聲器的位置坐標。假設聲源為單極子點聲源，輻射聲波為球面波。聲波頻率為f，聲速為c，波數k＝2πf/c，(r0,Ω0)表示聲源的位置坐標，tq(kr0,Ω0)表示聲源到q號傳聲器的聲場傳遞函數，根據散射理論，tq(kr0,Ω0)≡Σn=0∞Σm=-nnRn(kr0,ka)Ynm*(Ω0)Ynm(Ωq)---(1)]]>其中，分別為Ω0、Ωq方向的球諧函數，n、m為球諧函數的階次，上標“*”表示共軛運算，Rn(kr0,ka)為聲信號傳播徑向函數：Rn(kr0,ka)=-4πihn(2)(kr0)(jn(ka)-jn′(ka)hn(2)′(ka)hn(2)(ka))---(2)]]>其中，被減數項為球面波在自由場中傳播的徑向函數，減數項為實心球引起的球面波散射徑向函數，為虛數單位，jn(ka)為n階第一類球貝塞爾函數，為n階第二類球漢克爾函數，j'n(ka)、分別為jn(ka)、的一階導數。用s(kr0,Ω0)表示聲源的強度，B表示所有聲源位置坐標組成的集合，則q號傳聲器接收的聲壓信號p(ka,Ωq)可寫為：p(ka,Ωq)=Σ(r0,Ω0)∈Bs(kr0,Ω0)tq(kr0,Ω0)---(3)]]>令p＝[p(ka,Ω1),p(ka,Ω2),…,p(ka,ΩQ)]Tt(kr0,Ω0)＝[t1(kr0,Ω0),t2(kr0,Ω0),…,tQ(kr0,Ω0)]T分別為傳聲器聲壓信號列向量和聲場傳遞函數列向量，上標“T”表示轉置運算，對應式(3)，p=Σ(r0,Ω0)∈Bs(kr0,Ω0)t(kr0,Ω0)---(4)]]>各傳聲器接收聲壓信號的互譜矩陣C寫為：C=ppH‾=Σ(r0,Ω0)∈BΣ(r0′,Ω0′)∈Bs(kr0,Ω0)s*(kr0′,Ω0′)‾t(kr0,Ω0)tH(kr0′,Ω0′)---(5)]]>其中，上標“-”、“H”分別表示平均運算和轉置共軛運算，(r0',Ω0')也為聲源位置坐標。步驟2：獲取實心球延遲求和(DAS)輸出函數過程201定義實心球DAS輸出函數DAS輸出函數為：W(krf,Ωf)=1QvH(krf,Ωf)Cv(krf,Ωf)||v(krf,Ωf)||22---(6)]]>其中，(rf,Ωf)為聚焦點的位置坐標，v(krf,Ωf)＝[v1(krf,Ωf),v2(krf,Ωf),…,vQ(krf,Ωf)]T為聚焦列向量，||||2表示2范數。v(krf,Ωf)中，元素vq(krf,Ωf)表示q號傳聲器的聚焦分量，被定義為：vq(krf,Ωf)=Σn=0NΣm=-nnRn(krf,ka)Ynm*(Ωf)Ynm(Ωq)---(7)]]>其中，為Ωf方向的球諧函數，Rn(krf,ka)為聚焦徑向函數，其表達式與式(2)類同，只需將式(2)中的r0替換為rf即可，N為球諧函數階次截斷長度，需人為設定。式(6)也可寫為：W(krf,Ωf)=Σq,q′=1QCqq′vq*(krf,Ωf)vq′(krf,Ωf)Q||v(krf,Ωf)||22---(8)]]>其中，q'也為傳聲器索引，Cqq'為q號傳聲器接收聲壓信號相對于q'號傳聲器接收聲壓信號的互譜。過程202變形實心球DAS輸出函數定義平均聲壓貢獻為聲源在陣列各傳聲器處產生聲壓信號自譜的均值，記為PAC(kr0,Ω0)，等于傳聲器接收聲壓信號互譜矩陣的跡平均。對于單聲源，互譜矩陣C可寫為：C＝S(kr0,Ω0)t(kr0,Ω0)tH(kr0,Ω0)(9)其中，被定義為具有功率意義的聲源強度。相應地，PAC(kr0,Ω0)=1Qtr(C)=1QS(kr0,Ω0)||t(kr0,Ω0)||22---(10)]]>其中，“tr”表示求矩陣的跡。根據式(10)，用PAC(kr0,Ω0)表示S(kr0,Ω0)后代入式(9)可得：C=PAC(kr0,Ω0)Q||t(kr0,Ω0)||22t(kr0,Ω0)tH(kr0,Ω0)---(11)]]>將式(11)代入式(6)得：W(krf,Ωf)=PAC(kr0,Ω0)vH(krf,Ωf)t(kr0,Ω0)tH(kr0,Ω0)v(krf,Ωf)||v(krf,Ωf)||22||t(kr0,Ω0)||22---(12)]]>psf被定義為實心球DAS聲源識別方法的點傳播函數，表示單位平均聲壓貢獻單極子點聲源的響應，根據式(12)該DAS方法的psf表達式為：psf((krf,Ωf)|(kr0,Ω0))=vH(krf,Ωf)t(kr0,Ω0)tH(kr0,Ω0)v(krf,Ωf)||v(krf,Ωf)||22||t(kr0,Ω0)||22---(13)]]>對單聲源，該DAS方法的輸出結果等于聲源的平均聲壓貢獻與相應點傳播函數的乘積。過程203確定球諧函數階次截斷長度球諧函數階次截斷長度按下式確定；其中，符號表示將數值向正無窮方向圓整到最近的整數。步驟3獲取實心球函數型延遲求和(FDAS)輸出函數函數波束形成(FB)方法首先特征值分解陣列傳聲器接收聲壓信號的互譜矩陣C，得：C＝UΣUH(15)其中，U＝[u1,u2,…,uq,…,uQ]為酉矩陣，Σ＝diag(σ1,σ2,…,σq,…,σQ)為對角矩陣，σq(q＝1,2,…,Q)為矩陣C的特征值，滿足σ1≤σ2≤…≤σQ，uq為σq對應的特征向量。然后，引入指數參數ξ，令：C1ξ≡Udiag(σ11ξ,σ21ξ,...,σq1ξ,...,σQ1ξ)UH=Σq=1Qσq1ξuquqH---(16)]]>最后，采用已有聲源識別方法后處理并計算所得結果的ξ次方。根據步驟2，構建的FDAS方法的輸出函數為：WF(krf,Ωf)=1Q(vH(krf,Ωf)C1ξv(krf,Ωf)||v(krf,Ωf)||22)ξ(ξ≥2)---(17)]]>顯然，當ξ＝1時，式(17)即為步驟2中DAS方法輸出函數的表達式。假設三維空間內僅存在一個單極子點聲源，平均聲壓貢獻為PAC(kr0,Ω0)，理論上，下列關系成立：σ1＝σ2＝…＝σQ-1＝0(18)σQ＝tr(C)＝QPAC(kr0,Ω0)(19)uQuQH=t(kr0,Ω0)tH(kr0,Ω0)||t(kr0,Ω0)||22---(20)]]>聯(lián)立式(16)～(20)可得：WF(krf,Ωf)=PAC(kr0,Ω0)(vH(krf,Ωf)t(kr0,Ω0)tH(kr0,Ω0)v(krf,Ωf)||v(krf,Ωf)||22||t(kr0,Ω0)||22)ξ=PAC(kr0,Ω0)psfξ((krf,Ωf)|(kr0,Ω0))---(21)]]>表明：FDAS方法對單極子點聲源的響應等于該聲源的平均聲壓貢獻與其點傳播函數ξ次方的乘積。步驟4對函數型延遲求和(FDAS)結果進行脊探測(RD)若二元函數在某位置的輸出值在其最大表面曲率方向上為極值，則該位置被稱為脊。特定頻率及聚焦距離下，F(xiàn)DAS算法的輸出函數WF(krf,Ωf)是關于聚焦方向Ωf＝(θf,φf)的二元函數，其在聲源方向的輸出值遠大于在其它聚焦方向的輸出值，聲源方向為明顯的極大值脊，可通過脊檢測剔除主瓣內的非脊點，以縮減主瓣寬度、提高空間分辨率。根據定義，脊也是最大表面曲率方向上一階導數符號改變的轉折點。記L(Ωf)＝|WF(krf,Ωf)|為脊檢測對象，“||”表示取模運算，對每個聚焦方向，可按下列過程檢測其是否為極大值脊：過程401.構造L(Ωf)的二階偏導數矩陣，即黑塞矩陣：H=∂2L∂φf2∂2L∂φf∂θf∂2L∂θf∂φf∂2L∂θf2---(22)]]>過程402.特征值分解矩陣Η，找出絕對值最大的特征值λ及相應的特征向量uλ，其中uλ指示最大表面曲率方向，輸出值在uλ方向上若為極大值則λ為負，若為極小值則λ為正。過程403.計算標量ρ(Ωf)為：ρ(Ωf)=-12sign(λ)|sign(▿L(Ωf+χuλ)·uλ)-sign(▿L(Ωf-χuλ)·uλ)|---(23)]]>其中，“sign”為符號函數，為梯度算子，“．”表示求向量內積，χ為常數，用以控制檢測空間精度，本發(fā)明取為0.2。當ρ(Ωf)＝1時，聚焦方向為極大值脊，保留FDAS的輸出值；反之，用零取代FDAS在該聚焦方向的輸出值來剔除非脊點。步驟5對延遲求和(DAS)結果進行反卷積(DAMAS)記G為聚焦點總數，W為G×1維列向量，對應元素為各聚焦點的延遲求和輸出量，由式(12)計算得出，A為G×G維的點傳播函數矩陣，由式(13)計算得出，PAC為G×1維列向量，對應元素為各聚焦點處聲源的平均聲壓貢獻，為未知非負數，對于不相干聲源如下線性方程組成立：W＝APAC(24)(24)式采用高斯-賽德爾迭代方案求解PAC，移除點傳播函數的影響，能夠顯著衰減旁瓣并獲得超高空間分辨率。初始化基于第l次迭代計算的結果進行第l+1次迭代的具體過程為：過程501.計算殘差re：re(l)(krf,Ωf)=Σ(r0,Ω0)∈Dpsf((krf,Ωf)|(kr0,Ω0))PAC(l+1)(kr0,Ω0)+Σ(r0,Ω0)∈Epsf((krf,Ωf)|(kr0,Ω0))PAC(l)(kr0,Ω0)-W(krf,Ωf)---(25)]]>其中，D為已完成l+1次迭代計算的聲源點組成的集合，E為未完成l+1次迭代計算的聲源點組成的集合，二者的并集包含所有聲源。過程502.確定PAC(l+1)(kr0,Ω0)=max(PAC(l)(kr0,Ω0)-re(l)(krf,Ωf)psf((krf,Ωf)|(kr0,Ω0)),0)rf=r0,Ωf=Ω0]]>按照該方案依次計算每個聚焦點，即可得出直接用DAMAS清晰化DAS的成像結果需要考慮所有聚焦點，龐大的聚焦點數會使反卷積承受嚴重耗時問題。步驟4中對FDAS結果進行脊探測獲取的脊中包含聲源點，且聲源點處的輸出值相對很大，可默認一定動態(tài)范圍外的脊和非脊處無聲源，輸出值為零，僅對動態(tài)范圍內的脊進行DAMAS運算來獲得超高空間分辨率，此時，G不再為聚焦點總數，而是參與運算的脊的數目，該數值很小，可使反卷積具有極高的效率。步驟6定位聲源位置并且能繪制聲源成像圖實現(xiàn)三維聲場可視化。根據步驟5求得的PAC即能定位聲源并能繪制聲源成像圖實現(xiàn)三維聲場可視化。附圖說明圖1是本發(fā)明球陣列波束形成坐標系圖；圖2是本發(fā)明測量布置示意圖；圖3是DAMAS方法3000Hz單聲源識別聲源成像圖；圖4是FDAS方法3000Hz單聲源識別聲源成像圖；圖5是本發(fā)明方法3000Hz單聲源識別成像圖；圖6是DAMAS方法3000Hz不相干雙聲源識別聲源成像圖；圖7是FDAS方法3000Hz不相干雙聲源識別聲源成像圖；圖8是本發(fā)明方法3000Hz不相干雙聲源識別成像圖；具體實施方式下面結合實施例對本發(fā)明作進一步說明，但不應該理解為本發(fā)明上述主題范圍僅限于下述實施例。在不脫離本發(fā)明上述技術思想的情況下，根據本領域普通技術知識和慣用手段，做出各種替換和變更，均應包括在本發(fā)明的保護范圍內。實心球聲源識別低旁瓣超高分辨率聲學圖像快速獲取方法:1)：構建測量系統(tǒng)：參見圖1，一個處于三維空間內的實心球陣列；所述實心球陣列包括一個半徑為a的實心球，以及分布在所述實心球表面的Q個傳聲器；這些傳聲器的編號為q，q＝1,2,…,Q；實施例中，Q＝36；所述實心球的球心作為原點，實心球陣列所處的三維空間內任意一點的位置坐標用(r,Ω)描述，r表示所描述位置與原點間的距離，表示所描述位置的方向，θ、分別為其仰角、方位角；(a,Ωq)為q號傳聲器的位置坐標；作為驗證，本實施例中，聲源為穩(wěn)態(tài)白噪聲信號激勵的揚聲器，采用公司、半徑為0.0975m、集成4958型傳聲器的36通道實心球陣列采樣聲壓信號，揚聲器聲源到陣列中心的距離均為1m。參見圖2，為本實施例測量布置示意圖。2)獲取聲壓通過硬件采集系統(tǒng)采集各個傳聲器接收聲壓信號p(ka,Ωq)；其中，波數k＝2πf/c，聲波頻率為f，聲速為c，(a,Ωq)為q號傳聲器的位置坐標；3)通過各個傳聲器接收聲壓信號，得到互譜矩陣C(可采用B&K公司的pulse軟件獲得)：C=ppH‾,]]>p=Σ(r0,Ω0)∈Bs(kr0,Ω0)t(kr0,Ω0),]]>其中：B表示所有聲源位置坐標組成的集合，(r0,Ω0)表示聲源的位置坐標，s(kr0,Ω0)表示聲源的強度，t(kr0,Ω0)＝[t1(kr0,Ω0),t2(kr0,Ω0),…,tQ(kr0,Ω0)]T，tq(kr0,Ω0)表示聲源到q號傳聲器的聲場傳遞函數，tq(kr0,Ω0)≡Σn=0∞Σm=-nnRn(kr0,ka)Ynm*(Ω0)Ynm(Ωq),]]>分別為Ω0、Ωq方向的球諧函數，n、m為球諧函數的階次，上標“-”、“H”、“*”和“T”分別表示平均運算、轉置共軛運算、共軛運算和轉置運算；4)確定各個頻率對應的球諧函數階次截斷長度：符號表示將數值向正無窮方向圓整到最近的整數。5)將互譜矩陣C和球諧函數階次截斷長度N分別代入FDAS和DAS輸出式：WF(krf,Ωf)=PAC(kr0,Ω0)(vH(krf,Ωf)t(kr0,Ω0)tH(kr0,Ω0)v(krf,Ωf)||v(krf,Ωf)||22||t(kr0,Ω0)||22)ξ]]>和W(krf,Ωf)=1QvH(krf,Ωf)Cv(krf,Ωf)||v(krf,Ωf)||22]]>其中：PAC(kr0,Ω0)為互譜矩陣C的跡，(rf,Ωf)為聚焦點的位置坐標，t(kr0,Ω0)＝[t1(kr0,Ω0),t2(kr0,Ω0),…,tQ(kr0,Ω0)]T為傳聲器聲壓信號聲場傳遞函數列向量，v(krf,Ωf)＝[v1(krf,Ωf),v2(krf,Ωf),…,vQ(krf,Ωf)]T，vq(krf,Ωf)表示q號傳聲器的聚焦分量，ξ為指數參數，推薦取值為16，||||2表示2范數，上標“T”表示轉置運算；6)對WF(krf,Ωf)進行脊探測(RD)并剔除非脊點和30dB動態(tài)范圍之外的脊點，然后針對剩余脊點對延遲求和結果W(krf,Ωf)進行反卷積(DAMAS)求得PAC。7)根據PAC繪制聲源成像圖，實現(xiàn)聲場可視化，準確定位聲源位置。用DAMAS方法、FDAS方法和本發(fā)明方法分別識別3000Hz單聲源，成像動態(tài)范圍均為30dB，聚焦聲源面被設定為與陣列同心的球面，仰角從0°到180°，方位角從0°到360°，間隔均為5°，共37×73個聚焦點，聚焦距離為1m，反卷積的迭代次數均取為5000。該揚聲器聲源平均聲壓貢獻的實測值為49.1649dB。參見圖3，為DAMAS方法3000Hz單聲源識別聲源成像圖，輸出主瓣各聚焦點處輸出值的線性疊加為49.1866dB，與實測值相差0.0217dB，具有超高分辨率并且旁瓣很??；參見圖4，為FDAS方法3000Hz單聲源識別聲源成像圖,輸出主瓣峰值為48.2983dB,與實測值相差較大為0.8666dB且分辨率較低；參見圖5，為本發(fā)明方法3000Hz單聲源識別成像圖,輸出主瓣各聚焦點處輸出值的線性疊加為49.1649dB,與實測值相差很小為0.0275dB，同時也具有超高分辨率并且旁瓣衰減能力極好；參見圖6，為DAMAS方法3000Hz不相干雙聲源識別聲源成像圖；參見圖7，為FDAS方法3000Hz不相干雙聲源識別聲源成像圖；參見圖8，為本發(fā)明方法3000Hz不相干雙聲源識別成像圖；顯見，本發(fā)明方法對單聲源及不相干聲源均能快速準確地定位聲源，成像圖清晰干凈。相對于直接對延遲求和結果進行反卷積(DAMAS)，本發(fā)明方法的計算時間僅為前者用時的千分之四，極大地提高了效率，同時與DAMAS方法具有同等的超高空間分辨及更優(yōu)的旁瓣衰減能力；相對于函數型延遲求和(DAS)方法，本發(fā)明延續(xù)了前者的強旁瓣衰減能力同時分辨率獲得極大提高，并且能更準確地量化聲源。當前第1頁1 2 3