基于相關峰細化時延估計的艙內服務機器人發話人定位研究

王嘯臻，王兆魁，張育林,

(1.國防科技大學 空天科學學院，湖南 長沙 410073； 2.清華大學 航天航空學院，北京 100084)

0 引言

當前空間站的穩定運營主要依賴于航天員的長期在軌工作。航天員需要獨立完成飛船操控、在軌科學實驗操作、科學演示等多項任務，在軌工作效能有限。隨著機器人技術的快速發展，利用機器人輔助完成在軌任務成為可能。國際上已先后提出多個艙內輔助機器人項目，如PSA[1]、Astrobee[2]，以及曾在國際空間站服役的Smart Spheres[3]和Int-Ball等。這些機器人多采用遙控方式運行，其自主運行、人機智能交互能力有限，需要航天員在軌照料。為發展更自主的艙內服務機器人，提升機器人的任務輔助能力，需要解決機器人對航天員的信息感知問題。慣性與視覺導航主要用于獲取艙內的相對位置和姿態，而艙內機器人無法有效實時地獲取航天員的相對位置。利用麥克風陣列進行聲源定位，感知航天員的語音信息，可有效區分航天員與環境，實時獲取航天員角度位置關系，這一方法計算速度快，結構簡單。

目前使用麥克風陣列進行聲源定位的方法一般可分為3類：基于高分辨率譜估計、基于可控功率響應和基于時延估計。而機器人聲源定位系統一般要求麥克風陣列數量少且易搭載，機器人運動特性則要求定位方法高效迅速且具備一定的抗混響和抗噪聲能力[4]。由于基于時延估計的聲源定位方法計算復雜度低，實時性高，且硬件成本較低，因此在機器人系統中得到了廣泛應用。

基于時延估計的定位方法是一種雙步定位方法，它首先估計聲源到達不同麥克風的時間差，再利用得到的時延結合麥克風陣列結構，確定聲源的位置[5]。自1976年以來相關學者提出了許多時延估計的算法，如基于聲學傳遞函數比的算法[6]、基于基音加權的算法[7]、基于自適應特征值分解的算法[8]、基于廣義互相關函數的算法[9]等。廣義互相關函數法由于其計算速度快且符合實時解算要求，因此應用最廣泛，它通過計算兩路信號的互相關函數求取峰值，相關函數峰值對應的時間就是兩路信號的時延差。這類方法原理簡單，但受噪聲和時域采樣頻率的限制，所以時延估計精度不高。

為提高時延估計精度，針對廣義互相關時延估計方法使用快速傅里葉變換(FFT)會在頻域上產生柵欄效應，從而影響時延估計精度的問題，提出了基于相關峰精確插值(FICP)的空間六元陣列發話人定位方法。該方法基于機器人球形結構設計，利用譜減法和二次相關降低噪聲干擾，通過對聲音信號進行頻譜細化提高相關函數的時間分辨率，能有效改善FFT帶來的柵欄效應，提高時延估計及定位精度。改進算法既繼承了相關峰插值算法的優點，又增強了語音處理的抗噪聲能力。本文對傳統廣義互相關算法和改進算法在時延估計準確性方面進行了比較分析，實驗結果表明：改進算法可有效改善發話人定位性能。

1 基于球形結構的空間六元陣列定位原理

考慮到艙內服務機器人為球形結構設計(直徑D=230 mm)，因此系統中的麥克風陣列結構采用空間六元陣。空間六元陣可進行全空域聲源定位，相比平面四元陣在保證分維特性的同時還擁有更優的定位精度[10]。艙內服務機器人結構如圖1所示。

圖1 艙內服務機器人結構Fig.1 Layout of assistant robots

在直角坐標系下對空間六元麥克風陣列進行描述，如圖2所示。6個陣元的坐標分別為M1(d/2，0，0)、M2(0，d/2，0)、M3(-d/2，0，0)、M4(0，-d/2，0)、M5(0，0，d/2)、M6(0，0，-d/2)，其中d為陣列孔徑；r為目標(聲源)S(x，y，z)到陣列中心(坐標原點)的距離；φ為目標方位角；θ為目標俯仰角；C為空氣中聲音傳播速率。假設目標為點聲源，到達各陣元的傳播時間為τi(i=1,2,…,6)，則目標到各陣元的距離為ri=C·τi(i=1,2,…,6)。若到達兩陣元的時延差為τij，則聲源到M1和M3的距離差d1=C·τ31，到M2和M4的距離差d2=C·τ42，到M6和M5的距離差d3=C·τ65。

圖2 空間六元陣列結構Fig.2 Schematic diagram of six-element cone-shaped microphone array

由麥克風陣列和目標的幾何位置關系可聯立方程組

(1)

由直角坐標與球坐標的轉換關系可得目標方位信息r、φ、θ。因r，ri?dj，故可得

(2)

(3)

對式(2)、(3)進行簡化和反三角變換可得方位角φ和俯仰角θ的解算表達式為

(4)

平面四元陣列的方位角和俯仰角的解算表達式[11]為

(5)

根據實際測量經驗，一般得到的方位角和俯仰角誤差較小，而距離的測量誤差較大。考慮到機器人的實際需要，僅得到高精度的角度測量信息即可準確找到人員位置，從而能為視覺導航與人機交互提供輔助，因此本文僅給出角度測量公式。式(5)中測角精度受聲速估計誤差和陣元安裝誤差影響。由式(4)可知，空間六元陣測量方位角φ和俯仰角θ僅與時延τij有關，與聲速C和陣列孔徑d無關，定位精度受時延τij和聲速C的估計誤差影響。對比平面四元陣，空間六元陣的測角受干擾更小，定位精度更高。由此可見，得到高精度的時延估計值即可求得較準確的角度位置信息。

2 基于相關峰精確插值的改進時延估計算法

2.1 相關峰精確插值算法

為提高相關函數的分辨率，對相關峰進行插值計算，一般采用基于信號內插法的時延估計算法，其估計精度受限于內插倍數，復雜度較高。如果要達到0.01采樣間隔估計精度，需要計算100倍插值，不僅計算量大，還會帶來新的誤差源，因此需要計算相關峰的精確算法。相關峰精確插值算法的基本思想是針對兩路長度均為N點的語音信號，計算其互功率密度譜并進行頻域補零，構建完整的頻譜序列以提高頻譜采樣率及時域內相關函數分辨率。FICP算法是在線性調頻Z變換(MCZT)的基礎上，利用互譜的前后兩端計算相關函數的峰值波形，并對相關峰進行細化。FICP算法的具體實現過程如下：

R1(k)=

(6)

2) 對R1(k)補零，將互譜拉長至N2(N2>N1)點，得到新的互譜

R2(k)=

(7)

3) 對R2(k)進行IMCZT逆變換得到相關函數R(n)，由此可突破時域采樣率限制，提高相關函數分辨率，其分辨率可提高N2/N1倍。因陣列信號的時延估計值在有限范圍內，故僅計算相關函數主峰左右各N點，組成長度為2N的新相關函數，降低運算量。

對于n=0,1,…,N-1的時域范圍，相關函數為

(8)

對于n=N2-N,N2-N+1,…,N2-1的時域范圍，取n=N2-N+n1,對應n1=0,1,…,N-1的相關函數為

Rn2(N2-N+n1)=

(9)

由Rn1(n)和Rn2(n)組成的相關函數為

Rn(n)=

(10)

經過復雜公式推導[12]，可得化簡后的相關函數表達式

(11)

(12)

線性調頻Z變換為計算細化頻譜的有效方法，計算精度高且運算速度快。DMCZT和DIMCZT分別代表線性調頻Z變換及其逆變換，其定義式分別為

(13)

x′(n)=DIMCZT[X(k)]=

(14)

2.2 算法分析與比較

選取一段時長0.3 s的語音信號進行分析，采樣頻率fs=44 100 Hz，FICP計算參數N1=44 100，N2=500 000。分別用FICP和FFT計算頻域頻譜，圖3選取了包含譜峰的一段頻譜圖(頻率為400～600 Hz)，可以看出用FICP得到的頻譜信息比FFT更豐富，譜峰高度也更精確。在譜峰附近，用FFT計算得到的頻率依次為470.1、478.2、486.4、502.6、510.7 Hz，用FICP計算得到的頻率依次為484、485、486、487、488 Hz，FICP的頻譜泄漏更少。通過互功率譜計算能提高時延估計值精度，頻譜分辨率從8 Hz提高到1 Hz。用FICP和FFT計算得到的相關函數波形比較如圖4所示。可以看出改進算法相當于對峰值附近的一段相關函數進行插值，使峰值附近的相關函數波形更精細，局部分辨率明顯提高，時延估計值更逼近真實值，精度更高。

圖3 FICP和FFT計算頻譜比較Fig.3 Calculated spectra of FICP and FFT algorithms

圖4 FICP和FFT計算相關函數波形圖比較Fig.4 Signal cross-correlation waveforms calculated by FICP and FFT algorithms

2.3 時延估計算法的改進

FICP算法雖然能提高相關函數的分辨率，但并未改進抗噪聲性能[13-14]。因此，通過語音信號預處理和二次相關進一步提高算法時延估計的魯棒性。

首先利用譜減法對語音信號進行降噪處理，設定前導無話段時長估計環境噪聲能量，將譜減算法與原信號“作差”得到降噪信號，并利用基于倒譜距離的端點檢測篩選有話段，減小互相關函數計算長度，從而提高計算效率。

采用二次相關的目的在于抑制噪聲，每進行一次相關計算，都能提高信號的信噪比[15]。對信號1作自相關，對信號1和2作互相關，將信號1的自相關函數和互相關函數作相關計算，進一步削弱噪聲影響。時延估計改進算法流程如圖5所示。

完整的改進時延估計算法如下：1)利用MCZT取代FFT計算兩路預處理后的語音信號1和2的細化頻譜來提升頻譜精度；2)計算信號1的自功率譜及信號1和2的互功率譜，在此基礎上進行二次相關，提高抗噪性能，得到二次互功率密度譜；3)對上述互功率密度譜進行相關峰精確插值，提高二次相關函數的分辨率，通過峰值檢測計算得到兩路信號的時延差值。

3 算法仿真與實驗測試

3.1 FICP算法仿真

實驗中分別選取頻率fs為20 000、40 000、80 000 Hz的單頻正弦模擬信號，預設兩接收信號時延差Td=1.7×10-4s，采樣點數N=5 000，FICP算法參數N1=10 000，N2=40 000，預設信噪比RSN=8 dB，仿真次數為10次，如圖6所示。

圖6 10次計算結果比較Fig.6 Comparison of calculation results for ten times

由圖6可知，當真實時延值為采樣間隔的整數倍時，基于FFT計算的GCC法計算的時延值始終解算為采樣間隔的整數倍，而FICP法計算的時延值則可準確求取真值，有

Ts=(1.7±0.02)×10-4s

(15)

式中：Ts為時延估計值。

3.2 真實語音信號處理測試

實驗中選取采樣頻率fs=16 000 Hz，采樣點數N=16 347的兩路真實語音信號，FICP算法參數N1=16 000，N2=100 000。經過譜減法降噪處理、倒譜距離端點檢測篩選有話段后計算時延差。端點檢測后得到有話段為58～137幀，計6 321個采樣點。語音信號處理如圖7所示。

圖7 語音信號處理Fig.7 Voice signal processing

從圖7可看出，用FICP法得到的計算頻譜比GCC法更精細，頻譜分辨率得到提升，峰值附近的相關函數波形更精細，局部分辨率明顯提高。用GCC法計算的時延估計值為3.75×10-4s，即6個采樣時長；用FICP法得到的時延估計值為3.94×10-4s，位于3.75×10-4s(6個采樣時長)和4.375×10-4s(7個采樣時長)之間。實驗結果表明：用FICP法計算的時延值更接近時延差的真實值。

3.3 麥克風陣列定位測試實驗

為驗證空間六元麥克風陣列的定位效果，進行了靜止目標的定位實驗，如圖8所示。實驗條件如下：六元陣列及信號處理模塊預設在簡易結構的球形模型中，模型直徑23 cm，放置在高約20 cm的空心基座上。采用全指向性的PDM數字MEMS麥克風，頻率響應20～20 kHz，信號采樣頻率為16 000 Hz。實際測量中，聲源(音響)固定在距離麥克風陣列的一定范圍內，通過改變目標聲源的相對位置，完成對方位角和俯仰角的實時測量。

圖8 空間六元麥克風陣列Fig.8 Six-element cone-shaped microphone array

第1組實驗中，聲源相對陣列距離分別為50、60、70 cm，保持同一水平高度，方位角φ間隔30°進行音頻采集，定位求解結果如圖9所示。

圖9 方位角測定值Fig.9 Results of azimuth angle measurement

第2組實驗中，聲源保持同一水平高度，方位角φ設為40°和70°，聲源相對距離為40、50、60、70、80、90 cm的6組，采集12次音頻，定位求解結果如圖10所示。

圖10 方位角測定值Fig.10 Results of azimuth angle measurement

圖9、10對比分析了基于廣義互相關與相關峰精確插值這2種聲源定位方法的定位結果。發話人定位旨在為機器人視覺導航提供與人員的相對位置關系，由于機器人服務工作在極其有限的艙內空間，優于5°的穩定定位精度即可滿足機器人對人員的定位需求。隨著相對距離的增大，基于GCC法的方位角測量誤差逐漸增大，而基于改進算法的定位準確率則能保持在90%以上。可以看到，基于改進算法的聲源定位結果明顯優于基于GCC法的聲源定位結果。

第3組實驗中，聲源保持同一方位角和水平距離，方位角φ設為70°，俯仰角θ為70°、80°、90°、100°、110°的5組，采集10次音頻，定位求解結果如圖11所示。

圖11 俯仰角測定值Fig.11 Results of pitching angle measurement

圖11對比分析了空間六元陣列與平面四元陣列的俯仰角測量結果。可以看到，平面四元陣的俯仰角測量誤差很大，而空間六元陣列的俯仰角測量精度相對較高。綜上所述，基于改進時延估計算法的空間六元陣列定位方法可有效改善發話人定位性能。

4 結束語

發話人定位是空間站艙內服務機器人實現人員相對位置關系確定及人機交互任務輔助的重要前提。本文提出了基于FICP和二次相關的發話人定位算法，并對算法的流程進行了推導說明。該方法基于機器人球形結構設計，利用信號預處理和二次相關降低噪聲干擾，通過對聲音信號進行頻譜細化以提高相關函數的時間分辨率，能有效改善FFT變換帶來的柵欄效應，提高時延估計及發話人定位精度。實驗結果表明：相對基于廣義互相關的定位方法，基于FICP的空間六元陣列發話人定位方法能得到更高的時延估計精度，從而對發話人進行更準確的定位。未來將圍繞發話人定位與視覺導航方法的融合問題展開進一步研究，實現艙內輔助機器人的穩定跟蹤。