提升聲源位置估測準確度之麥克風陣列幾何配置

(1)

國立交通大學

電控工程研究所

碩士論文

提升聲源位置估測準確度之麥克風陣列幾何配置

Geometrical Arrangement of Microphone Array for

Accuracy Enhancement in Sound Source Localization

研究生：蔡傑名

指導教授：胡竹生博士

(2)

提升聲源位置估測準確度之麥克風陣列幾何配置

Geometrical Arrangement of Microphone Array for

Accuracy Enhancement in Sound Source Localization

研究生：蔡傑名 Student：Chieh-Min Tsai

指導教授：胡竹生博士 Advisor：Dr. Jwu-Sheng Hu

國立交通大學

電控工程研究所

碩士論文

A Thesis

Submitted to Institute of Electrical and Control Engineering College of Electrical Engineering

National Chiao Tung University in partial Fulfillment of the Requirements

for the Degree of Master in

Electrical and Control Engineering

July 2010

Hsinchu，Taiwan，Republic of China

(3)

提升聲源位置估測準確度之麥克風陣列幾何配置

研究生：蔡傑名指導教授：胡竹生博士國立交通大學電控工程研究所碩士班

摘要

本論文提出一套在麥克風數量固定與麥克風間距離有最大限制的情況，在不同訊噪比(SNR)下找出一個對角度估測與位置估測的最佳麥克風陣列排列方法。文中所使用的角度和位置的估測方法乃基於麥克風之間聲源傳遞的時間差(Time Difference of Arrival, TDOA)。因為此時間差的估測在實際狀況下有誤差，因此導致聲源角度和位置的估測不準，而此不準度與陣列的幾何配置有密切關係。本論文提出以估測矩陣的奇異值(singular value)對估測變異量的影響，進而探討麥克風陣列的最佳排列方法。在遠聲場的條件下，本論文得到精確的全域解析解可以使角度估測的 Cramér-Rao Lower Bound (CRLB)最小，而在近聲場的情形下，其結果需要對非線性聯立方程式求解，且與聲源位置有關。因此本論文探討在額外的限制條件下，其數值解的結果與特性。

(4)

Geometrical Arrangement of Microphone Array for Accuracy

Enhancement in Sound Source Localization

Student: Chieh-Min Tsai Advisor: Prof. Juw-Sheng Hu

Institute of Electrical and Control Engineering

ABSTRACT

This work proposes an optimal structure of microphone array for sound source localization under the conditions of fixed microphone numbers, limited distance between each microphone pair and different signal-to-noise ratios (SNR). The location estimation algorithm is based on the time difference of arrival (TDOA) of each microphone pair. The estimation accuracy is influenced by the TDOA errors and is dependent on the microphone array geometry. This work proposes a geometry optimization formulation using the relationship between the singular values of the estimation matrix and the estimation error. Under the far field condition, an analytical solution can be derived to minimize the Cramér-Rao Lower Bound (CRLB). In near field, the result leads to a set of nonlinear equations and is source location dependent. With additional constraints, several numerical solutions are discussed.

(5)

誌謝

兩年的研究碩士生涯，或許也是最後的學生生涯，就要畫下句點了。由衷的感謝實驗室大家長，也是我的指導教授胡竹生博士兩年的教誨，感謝老師在我研究路上給予適時的提點及建議，而老師積極的態度、各領域深厚的研究理論基礎，都是應該追求的典範。感謝實驗室的大家，感謝帶我進入聲音領域，也在我剛進入計畫時給我這小毛頭很多幫助的興哥，祝你今年順利畢業；實驗室的大姐大，也在最近剛成為人妻的鏗元學姊；在我研究生涯完成結婚和生小孩兩件人生大事的永融，祝你們全家幸福；還有跟我同天生日的聯誼達人丸子哥晉源；任何事都難不倒的實驗室強者阿吉；帥氣、認真的明唐，任何事總能解決，在計畫中一起度過了不少時光；同班加同實驗室 6 年的拳霸 judo，祝你在博班不歸路能快速秒殺它；還有做事認真的昌言學長；在我碩一陪我度過一年的強者戰神 dodo，歡迎你回到實驗室；衣著十分有個性、常跟我閒聊的 Lundy 大大；還有熱愛籃球的肉鬆；還有大學是電控也是實驗室同梯的前熱音社社長、衝浪男孩 simon；看似很閒，但還是有在認真做事的沛錡；還有題目換來換去，但都可以應付好的阿 him；還有實驗室的碩一以及碩零學弟；兩位跟我在地下室奮鬥的偉庭、新文；新任管理員&實驗室學妹湘筑；在九樓努力奮鬥的學文、昀軒、建安；兩位未來實驗室台柱昭男、耕維；新進來的學弟們建廷、哲鳴、宗翰。除此之外，更要感謝蛋糕社的夥伴，這我在碩二這年生活更多的充實，阿宅&阿斯豆腐的修明；正在美國旅遊打工，陪我聊天及吃飯玩樂的文文；還有各個 101 級的幹部；正往菸酒生邁進，對人非常好的承閔；同年且打扮超夜店咖的周傑；即將考研究所的 judy 和羅心妤，兩位請多加油。更感謝我的父母和家人，在我有壓力的時候適時給予我鼓勵，讓我可以度過這一切難關。最後感謝所有的人，也感謝這陪伴我六年的交大。

(6)

摘要 ...I ABSTRACT ... II 誌謝 ...III 表列 ... V 圖列 ... VI 第一章緒論... 1 1.1 研究動機 ... 1 1.2 文獻回顧 ... 3 1.3 論文貢獻 ... 6 1.4 論文架構 ... 6 第二章 GCC 和聲源方位估測演算法 ... 7 2.1 GCC 演算法 ... 7 2.2 聲源方位估測演算法 ... 9 2.2.1 聲源方位估測原理 ... 9 2.2.2 有限制的最小平方演算法... 12 第三章麥克風位置擺放最佳化... 14 3.1 奇異數(singular value)... 14 3.1.1 奇異數介紹... 14 3.1.2 奇異數與麥克風位置的關聯 ... 14 3.2 麥克風最佳化排放位置 ... 17 3.2.1 遠聲場的最佳擺放位置... 17 3.2.2 近聲場的最佳擺放位置... 35 3.3 模擬結果與分析 ... 52 第四章結論與未來展望... 59 Reference ... 60

(7)

表列

表一、遠聲場聲源下不同麥克風陣列擺放對不同角度聲源的角度誤差變異量…………...34 表二、近聲場聲源下不同麥克風陣列對90度附近聲源的方位估測誤差(一)………55 表三、近聲場聲源下不同麥克風陣列對90度附近聲源的方位估測誤差(二)………56 表四、近聲場聲源下不同麥克風陣列對45度附近聲源的方位估測誤差(三)………57 表五、近聲場聲源下不同time delay誤差對兩種麥克風擺放對90度來向聲源的方位誤差表………...…...58

(8)

圖列

圖一、遠聲場情況下二維麥克風最佳排列範例……….…22 圖二、Sensor示意圖………..23 圖三、遠聲場情況下三維麥克風最佳排列範例……….…………25 圖四、三維聲源示意圖……….……25 圖五、遠聲場情況下根據角度變異量設計不同的3種麥克風陣列擺放方式………….…..33 圖六、近聲場聲源對遠聲場假設下之不同麥克風擺放的角度偏差值……….……36 圖七、牛頓法收斂曲線及結果……….………50 圖八、牛頓法步驟的變異量曲線圖……….…………50 圖九、牛頓法收斂曲線及結果(二)……….……….51 圖十、牛頓法步驟的變異量曲線圖(二)……….………….51 圖十一、遠聲場聲源下麥克風擺放對多方位來向聲源的角度與聲速誤差(一) ………….52 圖十二、遠聲場聲源下麥克風擺放對多方位來向聲源的角度與聲速誤差(二) ………….53 圖十三、遠聲場聲源下麥克風擺放對多方位來向聲源的角度與聲速誤差(三)…….…….54 圖十四、近聲場聲源下麥克風擺放對多方位來向聲源的方位誤差(一)……….….55 圖十五、近聲場聲源下麥克風擺放對多方位來向聲源的方位誤差(二)………..…56 圖十六、近聲場聲源下麥克風擺放對多方位來向聲源的方位誤差(三)……….….57 圖十七、近聲場聲源下不同時間延遲誤差不同麥克風擺放對90度來向聲源的方位誤差圖………...58

(9)

第一章緒論

1.1 研究動機

空間中的聲源定位(Sound Source Localization)是聲音訊號處理領域中一項重要的技術，在聲源分離或是降低環境噪音干擾的研究上，目標或是干擾聲源的位置資訊對提升分離或降噪效能上有很大的幫助。此外，在語音處理相關應用中，聲源的位置也是整體系統中的一項重要資訊，例如在視訊會議中確認發言者的位置，或是智慧型機器人辨別交談者的方位等。一般較為準確的聲源定位均需仰賴麥克風陣列，指的是由若干麥克風按照一定的方式佈置在空間不同位置上組成的陣列。麥克風陣列因為其空間選擇性，可以在一定範圍內實現聲源的定位。現有使用麥克風陣列的聲源定位技術大致可分為三類。第一類是基於最大輸出功率的可控波束形成(beamforming)技術。其基本概念是將各陣列元件採集的訊號進行加權求和形成波束，通過搜索聲源的可能位置來引導該波束，並修改權值使得麥克風陣列在聲源方向的輸出信號功率最大。在傳統的波束形成器中，權值取決於各陣列元件上訊號的相位延遲，相位延遲與聲達時間延遲(delay of arrival, DOA)有關，因此稱為延遲求和波束形成器(delay-and-sum beamformer)，第二類是基於訊號頻譜差異的空間特徵估計技術。此方法主要有自我回歸模型(autoregression)、最大熵法(maximum entropy)、最小平方差估計法(minimum variance estimation) 和特徵值分解方法(multiple signal classification，MUSIC ；estimating signal parameters via rotational invariance techniques, ESPIRIT)等方法。第三類是基於聲音到達各陣列元件的時間差(time difference of arrival, TDOA)的定位技術。 TDOA 的估測最常使用的方法就是取兩顆麥克風間的 GCC(Generalized Cross-Correlation)，藉由找出GCC的最大值就等於找到

(10)

了聲音傳遞到這兩顆麥克風的時間差，在經由時間差的資訊便可以推測出聲源相對於這兩顆麥克風的方位。當麥克風個數變多時，即可以利用三角定位法的原理計算聲源的位置。TDOA的方法比之第一、二類方法有兩個優點，第一是可以經過變數轉換，可以間接估測聲音的傳遞速度，第二是其方法可同時適用在近聲場與遠聲場的情況。本論文及以TDOA 的方法來探討定位的準確性。一般很多研究著重在於如何逼近麥克風中的訊號空間轉換使得 TDOA 的誤差縮小，如此對聲源方位或位置估測的誤差就會變小。然而，一旦麥克風的數量增加，麥克風的排列方式也會影響到估測的準確性，例如不恰當的陣列幾何配置可能在某些角度的聲源會造成極大的誤差。其次，對陣列的設計而言，其基本的問題是在無可避免的 TDOA 估測誤差下，如何設計陣列的幾何配置以達到最佳的聲源定位。以往這方面的研究大部分以遠聲場為主，而本論文則嘗試以一統合遠近聲場的定位模型，基於最小方差無偏估計(minimum variance unbiased estimator)，分別探討在近聲場與遠聲場的情況下，麥克風陣列的最佳幾何配置關係。

(11)

1.2 文獻回顧

估測訊號在空間中的方位，或者是入射方向與感測裝置之間的關係，此類研究一般來說稱為Direction of arrival (DOA) estimation，DOA 的技術層面其中一類為基於聲音到達各陣列元件的時間差(time difference of arrival, TDOA)的定位技術，TDOA的估測最常使用的方法就是取兩顆麥克風間的GCC (Generalized Cross-Correlation)，藉由找出這個關係式的最大值就會等效於找到了聲音傳遞到這兩顆麥克風的時間差，而GCC的定義卻延伸出相當多的變化，也是推動此類方法持續演進的關鍵技術，這類方法主要都是轉換到頻率域，並且定義出各自的Weighting function。 Roth[1]所提出的Roth Processor，好處在於可以壓低第一顆麥克風電子雜訊和空間雜訊的所在頻率所產生的影響，估計的正確性。Carter等人[2] 提出的PHAT，如果兩顆麥克風彼此的雜訊的分佈情況是無關的此法可以有相當好的效果，這也是本文中使用的方法。Carter等人[3]又提出的 SCOT，此法是同時去考慮到兩顆麥克風的雜訊所造成的影響，而想較 Roth只考慮到第一顆麥克風的雜訊所造成的影響，可說是SCOT的一個特例。Knapp以及Carter提出的ML (Maximum Likelihood)方法[4]，放鬆了高 SNR 這樣的限制，且假設雜訊以及聲源的統計分佈都是高斯分佈 (Gaussian)，最後推演出新的Weighting function、Brandstein等人的[5]，此方法並不需要求取channel的響應，而是利用新的統計函數(Tukey’s Biweight)來重新定義Weighting function，並且推演出對於反射情況下更穩健的估測方法，以上都是在探討這方面的問題。以上方法都是用於1-D情況，藉由時間差的資訊，便可以推測出聲源相對於這兩顆麥克風的方位，而一旦麥克風數量增加時，麥克風的位置就不一定單單只是一個一維的情況。在多麥克風的情況下，兩兩一對的麥克風對也就增加了，因此可以獲得的資訊也就變多，除了分別利用外，另一方

(12)

面也可以整合這些資訊來提升估測的效果，但原本的關係式為一個非線性的方程式，所以最常用的方法是固定其中一顆麥克風當做參考點，於是問題就可以轉換成一個線性方程式，Yao[6]中提出利用TDOA和最小平方法 (Least square)，在假設聲速不知道的情況下去估測聲源位置與聲速，更加在[7]中提出了一套有限制的最小平方演算法，其中針對估測值之間的關係加以限制，在最小平方法的基礎和利用已知聲速的理想範圍，加以改良估測求得較合理的估測值，有利於提升估測的準確度。而後[8][9]也在聲速已知的條件下，分別提出限制的最小平方演算法的封閉解。[8]中也是利用估測值之間有個關係式，利用其當成一條限制方程式，使用Lagrange multiplier 來達成限制的最小平方演算法，但因為要求得Lagrange multiplier的係數沒有一個封閉解，對於計算上比較麻煩，所以它提中另外一套近似的最小平方演算法，利用泰勒展開式線性，則可找到對估測的封閉疊代解。[9]中則是深入探討Lagrange multiplier的方法，且針對Range measurement和

Range-difference measurement分別加以討論。

探討陣列形狀的文章，從最簡單的直線分布，[10]探討了在直線不同位置擺放，根據Power spectral dispersion，結論出要利用Stochastic Region Contraction (SRC) algorithm來找到最佳的解， [11]主要對線段上的sensor 做最佳化，不限定在直線上。[11]中，限制麥克風只能分佈在已定的線段上，基於Cramér-Rao Lower Bound(CRBL)去最小化變異量的下界，並分別去探討三種不同的情況，角度估測、距離估測以及位置估測，都有其不同的擺放方式。在[12]中基於Cramér-Rao Lower Bound (CRLB)對聲源位置的變異量做最佳化，在不限制麥克風和聲源位置的情況下，可以得到當2維陣列以聲源為中心呈正多邊形分佈或3維陣列以聲源為中心呈正多面體分布時，可以得到最小的位置變異量，這跟一般的直覺是一樣的。 [13]則延續[12]的工作，給予更多除了[12]提到之外的排法，更特別的探

(13)

討maximum Fisher information matrix, minimum Cramér-Rao Bound,

spherical codes, uniform angular arrays之間的關係。[14]則是基於線性方程式係數矩陣的條件數(condition number)會影響估測誤差的上界，利用使條件數最小化來決定感應器(sensor)的位置，但只探討最簡單的情況，也就是4顆感應器依序分布於原點、以及3軸上，結果可看出來跟其他非條件數最佳化的陣列相比可以有較小的誤差。[15]則是針對遠聲場的聲源，考慮訊號及雜訊的功率譜密度(signal and noise power spectral densities)對時間延遲(time delay)的影響，找出方位角(azimuth)和仰角(elevation)的 CRLB，進而可以推論出麥克風位置對其變異量的影響，為了最小化CRLB 而提出對麥克風位置的一些條件，於這些條件下可以決定麥克風的位置。[16]則是也是針對遠聲場的聲源，加上在聲速已知的條件下，基於線性方程式對azimuth和elevation的估測，分別探討azimuth和elevation的 CRLB，結論出線性方程式係數，也就是麥克風位置對CRLB的影響。

(14)

1.3 論文貢獻

本文假設聲源方位以及聲速未知的情況下，以最小方差無偏估計 (minimum variance unbiased estimator)為基準，結合奇異數(singular value)的觀念，對於麥克風位置加以探討。以遠聲場聲源而言，對於全方位估測，結論出一套可任意擺放和即時調整的方法，該陣列對於所有方向的聲源都有相同的估測準度，並給予其滿足條件的方程式，更佳也證明正多邊形和正多面體都滿足條件；對於已知的小範圍方位估測，可近一步利用縮小角度變異量，可精確其估測準度。對近聲場聲源，對於誤差分布有進一步的分析，分析在方程式係數也有誤差的情況下，對近聲場該如何估測，在估測位置方面，一樣對全方位估測和區域方位估測有所分別探討，並提出一套演算法供區域方位估測設計出精確較佳的麥克風擺放結構。 1.4 論文架構 本論文從聲源估測方法講起，再依序介紹近聲場與遠聲場的差異及分析，最後透過模擬來比較，主要內容如下：第二章：主要GCC演算法與聲源估測演算法及其改良方法。第三章：從聲源估測演算法中探討奇異數與麥克風位置的關聯性，更加深入分析如何找出最好的麥克風位置，分開遠聲場和近聲場的情況，最後為模擬與討論。第四張：結論與未來展望

(15)

第二章 GCC 和聲源方位估測演算法

2.1 GCC 演算法 首先假設空間中有一個聲源，在理想情況下兩個麥克風收到的聲音訊號可以分別表示為： 1( ) 1( ) 1( ) x t s t n t (2.1.1) 2( ) 1( ) 2( ) x t s tD n t (2.1.2) 其中假設s t₁( )，n t₁( )和n t₂( )都是WSS，並且s t₁( )和n t₁( )以及n t₂( )是 uncorrelated。其中D就是真正的延遲(Delay)，而 代表的是改變大小的Scale value，並且其D和的改變相較於s t₁( )是緩慢變化的，而此時麥克風之間的 cross correlation可以表示為：





1,2( ) 1( ) 2( ) x x R  E x t x t (2.1.3) 其中E代表期望值，而使得(2.1.3)達到最大值的 即為兩個麥克風之間的時間延遲量，因為實際觀測的時間是有限的，所以cross correlation的估算可以表示為： 1,2 1 2 1 ˆ _{( )} T _{( ) (} ₎ x x R x t x t dt T       



(2.1.4)

其中T代表觀測的時間間隔，而cross correlation 和cross power spectrum 之 間的關係可以用Fourier Transform 表示成下面這個形式： 1 2 1 2 2 , ( ) , ( ) j f x x x x R   G f e  df  

_

(2.1.5) 現在來考慮實際空間的狀態，麥克風收到的聲音訊號是經過空間轉換後的狀態，因此實際麥克風之間的cross power spectrum 可以表示為:

1 2 1 2 * , ( ) 1( ) 2( ) , ( ) y y x x G f H f H f G f (2.1.6) 其中，H₁( )f 以及H₂( )f 分別代表聲源到第一個麥克風以及第二個麥克風的空間轉換函數，因此定義麥克風對之間的Generalized correlation 為: ( ) 2 , ( ) ( ) , ( ) g j f x x g x x R  

_

 f G f e  df (2.1.7)

(16)

其中 * 1 2 ( ) ( ) ( ) g f H f H f   實際上，在有限的時間觀察下只能以估算的 ˆ ₁_,₂_{( )} x x G f 取代 1,2( ) x x G f ，因此改寫為: ( ) 1 2 1 2 2 , ˆ , ˆ g _{( )} _{( )} _{( )} j f x x g x x R   f G f e  df  

_

(2.1.8) 由(2.1.8)就可以估算麥克風對之間的時間延遲量，而_g( )f 的選取對於時間延遲的估計也有一定的影響，而其中一種方法是由Carter 等人提出的名為 PHAT (The Phase Transform)的方法，也就是:

1,2 1 ( ) ( ) g x x f G f   此方法在兩個麥克風彼此的雜訊分佈是無關的時候，有相當好的效果。

(17)

2.2 聲源方位估測演算法 2.2.1 聲源方位估測原理 傳統的聲源定位演算法如MUSIC[17]，必須將聲音傳遞的速度定為已知，雖然關於聲速的數值有理論值可以做為參考，但是實際上，並不知道聲音傳遞的速度，K. Yao 等人[6]提出了一套在聲速為未知的情況下依舊能夠估測聲源位置的方法，此方法以TDE的為基礎去完成的，本論文提出之方法是以此方法為基礎，同時考慮到麥克風陣列的形狀以及聲源距離對估計上的影響，以下將詳細說明及推導此演算法。根據[6]提出的方法([6]中第(13)式)，可以得到以下的式子: 2 2 1 1 1 1 1 1 1 1 ( ) ( ) | | ( ) ( ) | | 2 | | 2 | | i s i i i s s s v t t t t v v             r r r r r r r r r r r r  (2.2.1) 其中，r_s代表聲源的位置向量，r_i [x_i y_i z_i]代表第i個麥克風的位置向量，而t_it_j代表第i個麥克風和第j 個麥克風之間的時間延遲量差(在此都將j 設定為第一個麥克風)，而v 代表聲速。接者定義以下兩個變數: 1 1 ˆ | | s s s    r r r r r 以及 1 1 | | | | i i s     r r r r (2.2.2) ˆ_s r 代表聲源位置向量的單位向量，而_i代表聲源距離和麥克風陣列大小之間的比例，而當聲源屬於遠聲場的時候，_i 1，將(2.2.2)代入(2.2.1)，可以得到以下的式子: 2 2 1 1 1 1 1 ˆ | | 1 ( ) ( ) ( ) 2 | | 2 s i i i i i i i v t t t t v v v            r r r r r r r (2.2.3) (2.2.3)式中，v t( _it₁)代表聲源到第i個麥克風和到第1個麥克風之間的距離，定義這樣的距離差為d_i，也就是: 1 1 ( ) | | | | i i s i s d v t t  r r  r r (2.2.4)

(18)

因此，(2.2.3)可以改寫為: 1 1 ˆ ( ) ( ) 2 s i i i i r r r f t t v       (2.2.5) 其中 1 1 | | | | | | | | i i i i i d d f v v     r r r r (2.2.6) 由(2.2.5)式可以發現，當_i 1的時候，(2.2.5)式可以近似為: 1 1 ˆ ( ) s ( ) i ti t v  r r r   (2.2.7) 將(2.2.7)視為遠聲場下的關係式，所以定義_i為聲場遠近的比例，而 f_i為近聲場的影響因數。由於(2.2.1)式是一個線性的方程式，所以聲源位置的估算可以利用最小平方法解(2.2.1)式而求得，首先定義聲源方向單位向量為: 1 2 3 ˆ [ ] | | T x x x v v   s  s s s r r x r (2.2.8) 以及定義另外兩個變數為: 4 5 1 1 1 , 2 | _s | 2 | _s | v x x v     r r r r (2.2.9) 因此M個麥克風的線性方程式可以寫成: Axb 其中 2 ₂ 2 1 2 1 2 1 2 2 1 1 1 ( ) ( ) ( _M ) _M (_M ) t t t t _ _ _ _ _                  r r r r A r r r r    2 1 1 M t t t t              b  ₄ 5 w w            s x x (2.2.10) 在3D 的環境下，需要求的變數有5個，因此需要至少6個的麥克風才可以求解，利用最小平方法的解為:





ˆ  T -1 T x A A A b (2.2.11) 而聲源位置以及聲速的估算如下: 1 4 4 ˆ ˆ 1 ˆ ˆ ˆ, ˆ ˆ ˆ 2 | | s s s s v or v w x      _  _   x x r r x (2.2.12)

(19)

當聲場接近遠聲場模型的時候(_i 1)，由(2.2.5)式可以得到新的方程式為:  s s A x b 其中 2 1 1 ( ) ( ) s M r r r r       _ _      A  (2.2.13) 由於使用遠聲場模型的關係，無法求得聲源的座標，只能得到聲源的方向:









ˆ | | s s s   -1 T T s s s -1 T T s s s A A A b x r x _{A A} _{A b} (2.2.14) 而聲速的估算為:





1 1 ˆ | _s| v   _-1 T T s s s x _{A A} _{A b} (2.2.15)

矩陣A和A_s必須是full rank，[18]有探討full rank的條件，此條件會跟麥克風

(20)

2.2.2 有限制的最小平方演算法 K. Yao 等人在[7]提出一套有限制的最小平方法解 2 2 1 1 1 1 1 1 1 1 ( ) ( ) | | ( ) ( ) | | 2 | | 2 | | i s i i i s s s v t t t t v v             r r r r r r r r r r r r  經由移項可修改為 2 2 2 1 1 1 1 1 1 ( ) ( ) ( ) | | ( ) | | 2 i s i s i i v t t v t t          r r  r r r r r r (2.2.16) 因此M 個麥克風的線性方程式可以寫成: Byc 其中 2 2 1 2 1 2 1 2 1 1 1 ( ) ( ) ( ) / 2 ( _M ) (_M ) (_M ) / 2 r r t t t t r r t t t t            _ _ _    B    2 2 1 2 1 | | 1 2 | _M | r r r r       _ _  _    c  ₄ 5 x x     _ _     s y y (2.2.17) 其中





2 1 2 3 1, 4 | 1|, 5 T s  y y y  s  y v s  y v y r r r r 利用[7]所提出的方法，在[7](16)式， 2 ' ' v B y c d 2 1 2 1 1 1 ( ) ( ) ' ( _M ) ( _M ) t t t t       _ _       r r B r r   4 ' x        s y y 2 2 1 2 1 | | 1 2 | _M | r r r r       _ _  _    c  2 2 1 2 1 ( ) 1 2 (_M ) t t t t        _ _  _    d  (2.2.18) 對於要求未知的向量，利用有限制的最小平方法解 † ' ' y B c v2 B d'† † ' B  B B( 'T ')1B'T (2.2.19) 利用(2.2.19)，定義 † '  p B c qB d'† 可以把式子表示成： 2 1 1 1 2 2 2 2 2 3 3 3 2 4 4 s y p v q y p v q y p v q v p v q         r (2.2.20) 在假設 r₁ 0的情況下，這些參數之間有個關係式 2 ₂ ₂ ₂ 1 2 3 s  y y y r (2.2.21) (2.2.20)帶入(2.2.21)式

(21)

2 3 2 2 2 ( )v ( )v ( )v 0      (2.2.22) 其中 2 2 2 1 2 3 2 1 1 2 2 3 3 4 2 2 2 1 2 3 4 4 2 4 2( ) 2 q q q p q p q p q q p p p p q p                  對於聲速v來說，這個方程式最多存在3個解，利用聲速的數值是有範圍的，可以選取一個最適合的聲速ˆv，所以對於聲源估測可以寫成 2 ˆ   y p v q 但此方法必須確定所估測的TDOA都是正確的，否則在錯誤的TDOA情況下，他也可能算到合理的聲速，導致估測錯誤。但在TDOA估測無誤的情況下，相較一般的LS可以改善誤差的大小。

(22)

第三章麥克風位置擺放最佳化

3.1 奇異數(singular value) 3.1.1 奇異數介紹 定義一個矩陣A， _ m n ,_m__n A  ，矩陣A的奇異數_i,i1 ~n，滿足



 



1 1 1 for all i i i n n n i           _ _     Av u 0 A v v u u 0    (3.1.1a) , n n m m AVUΣ V U (3.1.1b) where V VT I_{n n}_ ,U UT I_{m m}_ V和U都必須是正交矩陣(orthogonal matrix)。 3.1.2 奇異數與麥克風位置的關聯 M 個麥克風的線性方程式可以寫成: Axb 其中 2 ₂ 2 1 2 1 2 1 2 ₂ 1 1 1 ( ) ( ) ( _M ) _M (_M ) t t t t _ _ _ _ _                  r r r r A r r r r    2 1 1 M t t t t      _ _      b  4 5 x x     _ _     s x x (3.1.2) 首先探討觀測(observation) b的誤差對於估測(estimator) x的影響 假設觀測到的ˆb有一個誤差b ˆ     b b b Ax b (3.1.3) 假設誤差是呈現高斯分佈且兩兩間沒關連 2 _ˆ 2 ( , ) ( , ) N   N  b 0 I b Ax I   (3.1.4)

(23)

利用Cramer-Rao Lower Bound(CRLB)，可以求得我們的最小方差無偏估計 (minimum variance unbiased estimator)和誤差下界(error lower bound)

2 ˆ __N₍ _,_ ₎ b Ax I 為高斯分佈，因此可以得到其機率密度函數(pdf) 2 1 _ˆ _ˆ ( ) ( ) 2 2 / 2 1 ˆ ( ) (2 ) T M p e        b Ax b Ax b (3.1.5) 利用去計算CRLB，(3.1.5)對估測值x做偏微分 2 2 2 ˆ ln ( ; ) 1 _ˆ _ˆ (2 ) ( ) ( ) 2 M T p x x        _    _  _{ }  b x b Ax b Ax (3.1.6a) 2 2 1 _ˆ 1 _ˆ 2 ( ) ( ) 2 T T T    A bAx  A bA Ax (3.1.6b) 1 2 ( ) ˆ T T T         A A A A A b x (3.1.6c) 根據CRLB定義，若可以表示(3.1.7)的形式，最小方差無偏估計則可以達到 CRLB ˆ ln ( ; ) ( )( ( ) ) p g     b x I x x x x (3.1.7) 其中 1 2 ˆ ( ) , ( ) ( ) T T T g    A A  I x x A A A b 由(3.1.7)式，可得到最小方差無偏估計為 1 _ˆ ˆ __g( )_( T ) T x x A A A b (3.1.8) 各估測值的變異量為 1 2 1 ˆ var( ) ( ( ) ) ( ( T ) ) i ii ii x _  _   I x A A (3.1.9) 所以估測的未知變數的分佈為 2 1 ˆ ~ N( , ( T ) ) x x A A (3.1.10) 由(3.1.10)可結論出誤差會正比於 2 1 ( T )  A A 

(24)

跟觀測值ˆb的誤差有關，假設此刻的影響固定，想要讓估測值的變異量和

最小

1 1

ˆ

min var( )_i min tr( ( ) ) min tr(( T ) )

i x    



I x A A (3.1.11a) 因為 1 1 tr( T ) _i, tr(( T ) ) i i i    

_



_

A A A A 1 i i i 1 min tr(( T ) ) min , _   _



A A 是A AT 的eigenvalue (3.1.11b) 根據柯西不等式(Cauchy-Schwarz inequality)









1 1 2 2 1 1 1 1/ 1/ tr( )tr(( ) ) 1/ T D D T T i i i i D D D                  _ _{ } _         



A A A A   (3.1.12a) 2 1 tr(( ) ) , tr( ) T T D   A A  A A D為陣列的維度 (3.1.12b) 當給定一些限制下，tr( T ) A A 可以固定為一個值，也就是說tr(( T ) )1 A A 有一個固定的下界，接下來的目的則是在此限制條件下找到最小的情況，也就是等號成立。底下章節會根據近聲場或遠聲場，分別決定所要的限制條件。很明顯的，等式在 T A A的eigenvalue都相同的情況下成立 cond( ) 1. T _ A A I A  (3.1.13)

所以可以看出來我們估測值的變異數會跟A A的T eigenvalues有關，也就是A

(25)

3.2 麥克風最佳化排放位置 3.2.1 遠聲場的最佳擺放位置 首先先看A矩陣 2 ₂ 2 1 2 1 2 1 2 ₂ 1 1 1 ( ) ( ) ( _M ) _M (_M ) t t t t _ _ _ _ _                  r r r r A r r r r    (3.2.1) 假設為遠聲場的環境下，在此條件下可以簡化(3.2.1)成  s s A x b 2 1 1 ( ) , ( ) s s s s M v         _ _ _{ } _          r r r A x r r r  (3.2.2) 根據最小平方(Least square)方法 1 ˆ_s ( T_s _s) T_s x A A A b (3.2.3) 可另r₁ 0讓我們矩陣更簡單 2 s M             r A r  s A 只跟麥克風的位置有關，也就是麥克風擺放的位置會影響A_s的奇異數根據(3.1.13)，如果要讓估測的變異量最少，所有的singular value要相同。

(26)

先從簡單的2維開始分析





1 1 1 Let (r_ir) a_i_ b_i_ for i2 ~M 1 1 2 2 1 1 1 1 1 1 2 2 1 1 1 1 1 1 1 1 T M M M s s s M M M M M a b a a a b a b a b a b b b a b                       _{ } _   _ _ _ _        A   A A     (3.2.4) 目前要找出(A AT_s _s)det(A AT_s _sI)0 det(A AT_s _s I)0 (3.2.5a) 1 1 1 1 2 2 2 2 2 2 1 1 1 1 ( ) ( )( ) ( ) M M M M i i i i i i i i i i a b a b a b            

_

 

_



_

(3.2.5b) 1 1 1 1 1 2 2 2 2 2 2 2 1 1 1 1 1 ( ) 4 ( ) 2 2 M M M M M i i i i i i i i i i i i i a b a b a b a b                    _  _ _  _        



 



(3.2.5c) 因為 T s s A A 是symmetric和semi-positive definite矩陣，所以_i0, 0 1 1 2 2 cond( ) , , 2 2 s              A (3.2.6) 目標是₁越靠近越好₂ ，同等於要越小越好 ( , )a b _optimum= 1 1 1 1 2 2 2 2 2 , ₁ ₁ ₁ ₁ arg min 4 ( )( ) ( ) M M M M i i i i i i i i i i a b a b a b                        



 a b (3.2.7) 如果我們假設M-1 2 2 i=1 K i i a b 



，也就是第2~M顆麥克風到第一顆的距離平方 2 1 ( _i )  r r 總合為一個定值，又可以進一步簡化(3.2.7)： ( , )a b _optimum= 1 1 1 2 2 2 , 1 1 1 arg min 4 ( )( ) ( ) M M M i i i i i i i a b a b                 



 a b (3.2.8a) = 1 2 1 2 1 2 , ₁ ₁ ₁ arg max ( )( ) ( ) M M M i i i i i i i a b a b            



 a b (3.2.8b) = 1 1 1 1 2 2 , 1 1, 1 1, arg max ( ) ( ) M M M M i j i i j j i j j i i j j i a b a b a b                

 

 a b (3.2.8c) = 2 1 2 1 2 2 2 2 , 1 1 1 1 arg max ( ( )) 2( ) M M M M i j j i i i j j i j i i j i a b a b a b a b                 

 

 a b (3.2.8d)

(27)

= 2 1 2 2 2 2 , ₁ ₁ arg max ( ( 2 )) M M i j i i j j j i i j i a b a b a b a b            

 

 a b (3.2.8e) = 2 1

_

_

2 , ₁ ₁ arg max M M i j j i i j i a b a b           

 

 a b (3.2.8f) i j j i a b a b又可以看成



a b_i, _i



跟



b_j,a_j



的內積，進一步化簡(3.2.8f)













2 1 2 1 ₂ 2 1 1 1 1 , , M M M M i j j i i i j j i j i i j i a b a b a b b a             

 

 (3.2.9a)













2 1 2 , 1 1 , , cos( 90 ) M M i i j j i j i j i a b b a       

_{ }

   (3.2.9b) 若假設



a b_i, _i



對所有i的長度固定=K'： 2 1 2 1 2 2 2 2 , , 1 1 1 1 ' cos ( 90 ) ' sin ( ) M M M M i j i j i j i i j i K  K            

_{ }

  

_{ }

(3.2.10a) 1 2 1 1 2 2 1 1 1 ' sin , 360 j M M M k k i j i k i i K               _ _   

 



 (3.2.10b) M是奇數跟偶數有點差別所以把(3.2.10b)分開表示 1 2 1 2 1 1 sin j M M k i j i k i              

 



( 3) / 2 1 ( 1) ( 1) / 2 ( 1) / 2 1 2 2 ( %( 1)) 1 1 1 ( 2) / 2 1 ( 1) 2 ( %( 1)) 1 1

sin sin , is odd

sin , is even M M j i M i M k n k i j k j i k i M M j i k n i j k j M M                                                    

 



 



(3.2.11) 然而可以證明當所有角度都相同的時候是其中的一個解 2 1 2 1 1 360

ˆ _sin _, _{,is one of the solution.}

1 j M M M M M i j i k i S M            _ _    

 



 (3.2.12) 當M顆麥克風的時候我們的eigenvalue為底下所表示 2 2 ( 1) ( 1) 4 , ( 1) 4 0 2 M M M K n K M S n S          (3.2.13) 代入(3.2.12)的答案

(28)

( 3) / 2 2 1 2 ( 2) / 2 2 1 ( 3) 4 ( 1) 2( 1) 2( 1) cos 2 1 ˆ ( 1) 4 ( 2) 2 ( 1) 2( 1) 2( 1) cos 2 1 M i M _M i M i M M M M M S M i M M M M                  _{ }    _ _ _ _ _  _ 



(3.2.14a) 2 2 ( 3) ( 1) 2( 1) 2( 1) ( 1)( 1 3 2) 0 2 ( 2) ( 1) 2( 1) ( 1) ( 1)( 1 2 1) 0 2 M M M M M M M M M M M M M M                     _ _ _ _ _ _ _ _{ } _{ } _   (3.2.14b) 由(3.2.14b)式可看出滿足最小的情況，也就是兩個eigenvalue相同。可以看出如果是使用硬算eigenvalue的情況下，算式會很複雜，一旦擴展到 3維的空間，會變的更複雜，所以這邊考慮另一種分析方式：

L. Yu. Kolotilina的[19][20]中提到，如果有一個n n Hermitian positive

definite matrix A，把他表示成 11 12

12 * 12 22 2 2 block form = _ _,    A A A A 0 A A ，而且假設 11 A 和A₂₂是正定矩陣，那可以證明A的極端eigenvalue滿足底下的不等式： 12 12 1( ) (1 ), n( ) (1 )  A  A  R  A  A  R R R (3.2.15) 其中 1/ 2 1/ 2 11 12 22    R A A A (3.2.16) 從(3.2.15)可知道，若A₁₂0，則RA₁₁1/ 2A A₁₂ ₂₂1/ 2 0 所以可以看出來 A₁( )一定不會等於 A_n( ) 反過來說，如果想要₁( )A 的下界_n( )A 的上界，A₁₂一定要為0。由式子可以知道要達到cond( ) 1A  有兩個條件以目前2維的例子，22矩陣的情況來說 2 2 11 12 1 1 1 1 1 1 * 2 2 12 22 1 1 1 1 1 1 T M M M s s M M M a a a b a b a b a b b b               _ __ _     _ _   A A A A A A     (3.2.17)

滿足λ (A)'s lower bound=λ (A)'s upper bound₁ _n a b_{1 1}a_M_₁b_M_₁0 (3.2.18) 除此之外，因為a b_{1 1}a_M_₁b_M_₁0

(29)

2 2 1 1 2 2 1 1 0 0 T M s s M a a b b              A A   很容易看出來， 2 2 2 2 1( ) a1 aM 1, n( ) b1 bM 1  A =  _  A =  _ 要滿足₁( )A _n( )A 2 2 2 2 1 M 1 1 M 1 a a _ b b _ (3.2.19) 所以可以結論 1 1 1 1 1 2 2 2 2 1 1 1 1 0 M M n M M a b a b a a b b             _{ }           (3.2.20) 令a_i R_icos

 

_i ,b_i R_isin

 

_i 2 2 1 1 1 1 1 1 2 2 2 2 2 2 2 2 1 1 1 1 1 1 1 1 0 M M M M M M M R a b R a b R a R a R b R b                        (3.2.21a) 2 2 1 1 1 1 1 1 2 2 2 2 2 2 1 1 1 1 1 1

cos sin cos sin 0

M M M M M M R R R R                               (3.2.21b) 2 2 1 1 1 1 2 2 1 1 1 1 sin 2 sin 2 0 cos 2 cos 2 0 M M M M R R R R                     (3.2.21c) 1 1 2 2 2 1 1 0 ( ) 0 i i M M j j i i i i R e  R e     

_

 

_

 (3.2.21d) 也就是說要設計一個麥克風陣列，在決定以一顆麥克風為準，其餘麥克風與它的向量 ji i R e 平方後相加要等於0。可以把(3.2.21d)稍為移項 1 2 1 ( i) 0 M j i i R e   



1 2 2 2 1 1 ( M ) ( i) M j j M i i R e  R e      



(30)

如果對R_i 沒有任何限制，也就是說第1到M 1顆麥克風任意擺放後，可以利用第M 顆來使這個陣列達到最好的位置，位置會由前面M 1顆決定。有M顆麥克風，理論上有 ( 1) 2 M M  個時間延遲資訊，但上面的理論只用到 1 M  個，其中原因也可以用上面來解釋，一旦用了其他的資訊，算出來的變異量有可能很糟，也就是說多用了其他資訊不一定會幫助估測，還有可能讓估測變差。假設R_i 都相同，則可以簡化(3.2.21d)的式子要滿足 1 2 1 0 i M j i e     



極座標單位圓上的點總合為0 明顯的 ( 1) 2 2( 1) i i M      和 2 ( 1) 1 i i M      為其中的一個解但理論上還有其他很多組的解，以上只是提到比較明顯的解。以M=9，也就是9顆的例子為例，圖一兩個情況都會滿足： (a) (b) 圖一、遠聲場情況下二維麥克風最佳排列範例近一步可以拓展，假設我們有K個subarrays，這K個subarrays都使用相同一顆當基準點，而每個都滿足₁_n，則這K個subarrays組成的array必定也滿足  ，所以組合的總類其實是無限多種的。 -10 -5 0 5 10 -10 -5 0 5 10 cond(A)=1 -10 -5 0 5 10 -10 -5 0 5 10 cond(A)=1

(31)

接著拓展到3維分析





1 1 1 1 ( _i ) _i _i _i for 2 ~ Let  r r  a_ b_ c_ i M 1 1 1 1 1 1 s M M M a b c a _ b _ c _     _ _     A    (3.2.22) 2 2 1 1 1 1 1 1 1 1 1 1 2 2 1 1 1 1 1 1 1 1 1 1 2 2 1 1 1 1 1 1 1 1 1 1 M M M M M T s s M M M M M M M M M M a a a b a b a c a c a b a b b b b c b c a c a c b c b c c c                           _       _  _ _ _ _ _ _    A A          (3.2.23)

以L. Yu. Kolotilina的理論來說，(3.2.25)式要滿足condition number=1，需要以下的式子成立： cond( ) 1          A 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 1 1 1 1 1 1 0 0 0 M M M M M M M M M a b a b a c a c b c b c a a b b c c                                 (3.2.24) (3.2.26)式可以分成3區 2 2 2 2 1 1 M 1 M 1 0, 1 M 1 1 M 1 a b a _b _  a a _ b b _ (3.2.25a) 2 2 2 2 1 1 M 1 M 1 0, 1 M 1 1 M 1 a c a _c _  a a _ c c _ (3.2.25b) 2 2 2 2 1 1 M 1 M 1 0, 1 M 1 1 M 1 b c b _c _  b b _ c c _ (3.2.25c) 也就是說在xy x, z y, z平面都要滿足2維條件 i  i  j j  i R Rj 圖二、Sensor 示意圖

(32)

假設a_i R_isin

_{ }

_i cos

_{ }

_i , b_i R_isin

_{ }

_i sin

_{ }

_i and c_i R_icos

_{ }

_i ，如圖二示意

帶入(3.2.25a), (3.2.25b), (3.2.25c)

1

2 2 1

sin cos sin 0

M i i i i i R      



,





1 2 2 2 2 1

sin cos -sin 0

M i i i i i R      



(3.2.26a) 1 2 1

cos sin cos 0

M i i i i i R      



, 1 2 1

cos sin sin 0

M i i i i i R      



(3.2.26b)





1 2 2 2 2 1

sin cos - cos 0

M i i i i i R      



,





1 2 2 2 2 1

sin sin - cos 0

M i i i i i R      



(3.2.26c) 可稍微合併和整理(3.2.26a) ,(3.2.26b) ,(3.2.26c) 1 2 2 2 1 sin i 0 M j i i i R e    



(3.2.27a) 1 2 1 cos 2 i 0 M j i i i R e   



(3.2.27b)





1 2 2 1 1-3cos 0 M i i i R    



(3.2.27c) 簡化式子則假設R_i 都相同，(3.2.27a)可改寫 1 2 2 1 sin i 0 M j i i e     



(3.2.27a-1) (3.2.27b)可改寫成 1 1 cos 2 i 0 M j i i e    



(3.2.27b-1) (3.2.27c)可經過移項寫成





1 2 1 1 cos 3 M i i M     



(3.2.27c-1) 也就是高度平方總合是個定值 1 3 M  ，與麥克風數量有關。

(33)

所以當然當2顆到第M顆是以正多面體的分佈時，是其中的一個特殊解，如

圖三所示。

tetrahed. cube octahed. dodecahed. icosahed.

-10 0 10 -10 0 10 -10 0 10 -10 ₀ 10 -10 0 10 -10 0 10 -10 0 10 -10 ₀ 10_-10 0 10 -10 0 10 -10 ₀ 10 -10 0 10 -10 0 10 -10 0 10 -10 0 10 圖三、遠聲場情況下三維麥克風最佳排列範例同理在相同原點情況下，任意K個正立方體組合也滿足。以上的分析，是針對最小平方法的估測變數讓x_s的變異量最小，但實際上x_s不一定是真正關心的值 1 2 3 cos sin sin sin cos s s s v _x x v v x v                _ _ _ _ _ _          _ _ _ _         r x r (3.2.28) 遠聲場真正關心的是角度(azimuth 和elevation )，如圖四中所表示   圖四、三維聲源示意圖

(34)

其中 2 1 1 atan(x )= ( )g x   x ， 2 2 1 2 2 3 atan( x x )=g ( ) x    x 也就是真正的估測(estimator)是 ( , )  所以CRLB分析的Fish imformation要改寫成 1 1 1 2 1 1 2 ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) T T T g g g g                _ _      _ _      _    x x x x g x g x I I x I x x x x x x x (3.2.29) 其中 1 2 1 ( )  ( T ) I x A A 1 2 1 ( ) ( ) ( ) ( ) ˆ var( ) ( ( ) ) ( ( ) ) T T T               g x g x g x g x I x A A x x x x (3.2.30) 接著要找到一個好的估測方程式由 2 1 _ˆ _ˆ ( ) ( ) 2 2 / 2 1 ˆ ( ) (2 ) T M p e       b Ax b Ax b (3.2.31) 因為現在的估測是角度，所以(3.2.31)式要對角度做偏微分 ˆ ln ( ; )p     b 2 2 2 1 _ˆ _ˆ (2 ) ( ( )) ( ( )) 2 M T          _    _   b Ax b Ax  (3.2.32) 2 1 ( ) _ˆ 2 ( ( )) 2 T T        x A b Ax (3.2.33) 2

sin cos cos cos

0

1 _ˆ

( ( ))

cos sin sin sin cos

T T v v v v v                             A b A Ax (3.2.34) 可以看出沒辦法找到一個CRLB形式 ˆ ln ( ; ) ( )( ( ) ) p g          b I 所以說沒有一個對角度有效的估測封閉解。

(35)

但藉由ML(maximum-likelihood) ˆ _{arg max ( ; )} ML p      (3.2.35) 隨著M變大，估測會慢慢趨近CRBL。 ML(maximum-likelihood): 找到 ，使 p( ; )  0 ln ( ; )p  0          ˆ ln ( ; ) 0 p      b

sin cos cos cos

0

ˆ

( ( ))

cos sin sin sin cos

T T v v v v v                           A b A Ax (3.2.36) 可以看出當₍ Tˆ T _{( ))} ₀    A b A Ax 可滿足(3.2.36)式 1 ML( )ˆ ( ) ˆ T T   x  A A A b (3.2.37) 求得x_ML( )ˆ 則可以反推ˆ ˆ    g(x_ML( ))ˆ (3.2.38) 也就是說一樣可以用 1 ML( )ˆ ( ) ˆ T T  _  x A A A b來當作計算估測的式子

(36)

接著來分析Cramer-Rao Lower Bound，也就是誤差下界 11 12 13 1 12 22 23 13 23 33 Let ( T )                      A A 誤差下界 1 2 ( ) 1 ( ) ( ) ( ) T T           g x g x I A A x x 2 2 1 sin cos cos sin sin cos 0 _cos sin sin = v ( ) sin cos sin

cos cos sin cos sin

0 sin T                                                     A A (3.2.39) 於是可以歸納出方位角跟仰角的變異量









2 2 1 2 v ˆ

var( ) sin cos 0 ( ) sin cos 0

sin T T           A A  (3.2.40a)





2 2 2 2 11 12 22 2 v ˆ

var( ) sin 2 sin cos cos

sin              (3.2.40b)





2 2 1 cos cos ˆ

var( ) v cos cos sin cos sin ( ) sin cos

sin T                    _ _      A A (3.2.41a) 2 2 2 2 2 2 2 11 12 22 2 13 23 33 ˆ

var( ) v ( cos cos 2 cos sin cos sin cos

2 cos cos sin 2 sin cos sin sin )

                            (3.2.41b) 如果 T c  A A I 1 2 2 2 2 1 2 11 2 2 1 v v ˆ var( ) sin sin M i i a              _ _ 



 (3.2.42a) 1 1 2 2 2 2 2 11 1 ˆ var( ) v v M i i a            _ _ 



 (3.2.42b) 看的出來_ˆ_{只跟目前聲源的}__有關，_ˆ_{則跟聲源在哪都沒有關係} 對於全方位估測可說是最好的擺法。

(37)

如果 T __c A A I 那聲源的跟對相同陣列就會造成影響









2 2 1 2 v ˆ

var( ) sin cos 0 ( ) sin cos 0

sin T T           A A  (3.2.40)





2 2 1 2 v ˆ

var( ) ( ) , sin cos 0 1

sin T T         d A A d d   





2 2 1 cos cos ˆ

var( ) v cos cos sin cos sin ( ) sin cos

sin T                           A A (3.2.41)





2 2 1 ˆ

var( )  v T( T ) , cos cos  sin cos  sin 1

  e A A e e    可以看出來，一但決定一個陣列形狀， T A A就決定了要 1 min xT(A AT ) x，很明顯可以知道是要找到(A AT )1最小eigenvalue的 eigenvector，由此eigenvector可反推出跟。 1 ( T ) A A 最小eigenvalue也就是A AT 的最大eigenvalue。但根據d和e來看，可以明顯看出來它們算出來的會差個90度，也就是說 ˆ var( ) 小，則var( )ˆ 會大，反之亦然。或是也可以在已知跟下去找出最好的陣列位置舉個例子：假設已知聲源 90 , 固定





2 2 11 2 v ˆ var( ) sin      (3.2.42) 2 2 2 2 22 23 33 ˆ

var( )  v ( cos 2 cos sin   sin ) (3.2.43)

提升聲源位置估測準確度之麥克風陣列幾何配置

國 立 交 通 大 學

電控工程研究所

碩 士 論 文