CN103620438A - 用于盲定位相关源的系统和方法 - Google Patents

Abstract

本发明提供了一种用于到达方向盲估计的系统和方法，其针对J<M个的未知信号的M个接收器的非线性1维阵列，这些未知信号通过未知通道传送并且经受高斯噪声。该方法克服了信号的在时频域中不相交和/或在统计学上独立的条件。似然性函数评估到达方向的概率。加权和正则化方案选择相关的频率分量。所有频率的联合的似然性函数被判读为成本或者目标函数。通过利用在所构造的似然性函数上操作的粒子滤波器和用于方向的简单转移核来将时域包括到估计策略中。

Description

用于盲定位相关源的系统和方法

技术领域

本发明涉及信号的盲分离。更确切而言，本发明涉及由相关的源生成的信号的、也称作信号的到达方向的盲分离。

背景技术

发电主部件（例如涡轮机、发电机、锅炉、变压器，以及相关的辅助设备，诸如锅炉馈送泵、冷却水泵、电扇、阀门、排气净化系统）操作问题通常特征在于其声学特征（在声波和最重要地在超声范围中）中的改变，例如电弧放电、轴承或润滑问题、活动部件、泄露等。基于接触的、对例如对于电厂操作的连续性而言非关键的辅助设备的所有子部件的声学发射监控通常是不经济的。非接触的、基于麦克风的系统可以同时监控带有大量不同类型的设备的大区域。该区域监控的概念在图1中示出。在此，利用包含1、2和3号传感器或麦克风的麦克风阵列以及在处理器上运行的源定位指令确定了来自测试（DUT）号2下的设备的故障。

该场景的挑战是如下。现有技术的、利用时频掩膜的源分离和定位方法将具有由于来自发电设备的声学发射的宽带特性以及继起的在时频域中的信号交叠而引发的问题。尤其，通常使用的信号在时频域中不相交的稀疏假设被违背，并且构建于该假设的算法、例如DUET将失效。被许多算法例如诸如在“M.Jansson,A.L.Swindlehurst,和B.Ottersten所著的Weighted subspacefitting for genral array error models.IEEE Transactions on Signal Processing,46:2484-2498,1997”和“T.Melia和S.Richard所著的Underdetermined blindsource separation in echoic enviroments using desprit,2007”中公开的WeightedSubspace Fitting、ESPRIT或MUSIC所需的另一假设是所有信号是独立的。该假设对于语音信号必须是近似真实的，而对于机器振动不必要。

从而，需要用于并非不相交和/或独立和/或不相关的信号的盲分离的改进的且新颖的方法和系统。

发明内容

在声学电厂监督的情况下，对于通过未知通道传送并且遭受高斯噪声的J<M个未知信号的M个接收器的非线性1维阵列提供了到达方向的盲估计方案。根据本发明的至少一方面提供的方法和系统克服了当前现有技术的“信号在时频域中不相交（W-正交性）和/或在统计上独立”的条件。基于涉及信号和测量的模型，构造了用于到达方向的似然性函数。为了选择相关频率，提供了新的加权和正则化方案且所有频率的相交似然性函数被判读为成本或目标函数。利用在所构造的似然性函数上操作的粒子滤波器和用于方向的简单转移核来将时域包括到估计策略中。还提供了一种用于评估结果的方法。

根据本发明的一个方面，提供了一种用于确定由多个源中的源生成的信号的到达方向的方法，其包括：处理器借助在时频域中操作的粒子滤波器来估计由所述源生成的之前的信号的第一到达方向；由多个接收器接收来自所述多个源的多个信号，该多个信号包括所述信号，至少两个所接收的信号相关联；处理器为所述信号的多个频率分量中每个频率分量确定正则化系数和加权系数；处理器基于该正则化系数和加权系数对于信号的多个频率分量中每个频率分量确定到达方向的似然性；处理器基于对于信号的多个频率分量中每个频率分量得到的到达方向的似然性来确定信号的到达方向的全局似然性；处理器更新粒子滤波器的多个权重；并且处理器通过应用粒子滤波器的更新过的权重来生成对信号的到达方向的估计。

根据本发明的另一方面提供了一种方法，其中，全局似然性L(θ)通过表达式来确定，其中θ表示多个源中所有源的所有到达方向的向量，n表示频率分量，N表示频率分量的总数，

是通过频率分量n的归一化的加权系数λ_n调节的、所有源的所有到达方向的、正则化的似然性。

根据本发明的又另一方面提供该方法，其中，用于频率分量n的加权系数是通过表达式

确定的归一化的加权系数λ_n，其中

而τ_n表示未归一化的加权系数，

表示用于信号子空间上的频率分量n的正则化的射影算子，θ表示多个源中所有源的所有到达方向的向量，X_n表示由多个接收器接收的多个信号的频率分量n的向量，并且

表示非负非递减的加权函数。

根据本发明的又另一方面提供该方法，其还包括：处理器确定所述信号的频率分量的极大后验似然性（MAP）估计。

根据本发明的又另一方面提供该方法，其中，通过

表示的MAP估计通过表达式 ${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n$ 来确定，其中，θ表示多个源中所有源的所有到达方向的向量，A表示导引矩阵（steering matrix），

表示用于频率分量n的导引矩阵的厄密共轭，λ_n表示用于频率分量n的正则化系数，I表示单位矩阵并且X_n表示由多个接收器接收的多个信号的频率分量n的向量。

根据本发明的又另一方面提供该方法，其中，粒子滤波器的更新过的权重

通过表达式 $w_i^k=\frac{w_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}{{\mathrm{\&Sigma;}}_{i=1}^P{w}_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}$ 来确定，其中，k表示更新时间，i表示粒子索引，p表示概率分布函数，P表示粒子总数，X_1:N表示在时刻k上的频率分量的向量，并且

表示在时刻k粒子i的到达方向估计。

根据本发明的又另一方面提供该方法，其中，借助缩放参数调节全局似然性。

根据本发明的又另一方面提供该方法，其中，多个源中的至少一个源是运动的源。

根据本发明的又另一方面提供该方法，其中，信号是声学信号。

根据本发明的又另一方面提供该方法，其中，多个源中的至少一个源与发电机关联。

根据本发明的一个方面提供一种系统，用于确定由多个源中的源生成的信号的到达方向，该系统包括：存储器，其能够存储数据和指令；处理器，其能够检索和执行来自所述存储器的指令，以实施如下步骤：借助在时频域中操作的粒子滤波器估计由所述源生成的之前的信号的第一到达方向；处理由多个接收器接收的来自多个源的多个信号，所述多个信号包括所述信号，至少两个所接收的信号相关；对于所述信号的多个频率分量中每个频率分量确定正则化系数和加权系数；基于所述正则化系数和所述加权系数对于所述信号的多个频率分量中每个频率分量确定到达方向的似然性；基于对于所述信号的多个频率分量中每个频率分量确定的到达方向的似然性来确定所述信号的到达方向的全局似然性；更新所述粒子滤波器的多个权重；以及通过应用所述粒子滤波器的更新过的权重来生成对所述信号的到达方向的估计。

根据本发明的另一方面提供该系统，其中，所述全局似然性L(θ)通过表达式

来确定，其中，θ表示所述多个源中所有源的所有到达方向的向量，n表示频率分量，N表示频率分量的总数，

是通过频率分量n的归一化的加权系数λ_n调节的、对于所有源得出的所有到达方向的、正则化的似然性。

根据本发明的又另一方面提供该系统，其中，用于频率分量n的加权系数是归一化的加权系数λ_n，其通过表达式

来确定，其中而τ_n表示未归一化的加权系数，

表示用于在信号子空间上的频率分量n的正则化的射影算子，θ表示所述多个源中所有源的所有到达方向的向量，X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量，并且

表示非负非递减的加权函数。

根据本发明的又另一方面提供该系统，其还包括：所述处理器确定所述信号的频率分量的极大后验似然性（MAP）估计。

根据本发明的又另一方面提供该系统，其中，由

表示的所述MAP估计由表达式 ${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n$ 来确定，其中，θ表示所述多个源中所有源的所有到达方向的向量，A表示导引矩阵，表示用于频率分量n的导引矩阵的厄密共轭，λ_n表示用于频率分量n的正则化系数，I表示单位矩阵，并且X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量。

根据本发明的又另一方面提供该系统，其中，所述粒子滤波器的更新过的权重

由表达式 $w_i^k=\frac{w_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}{{\mathrm{\&Sigma;}}_{i=1}^P{w}_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}$ 来确定，其中，k表示更新时间，i表示粒子索引，p表示概率分布函数，P表示粒子总数，X_1:N表示在时刻k上频率分量的向量，并且

表示在时刻k对粒子i的到达方向估计。

根据本发明的又另一方面提供该系统，其中，所述全局似然性借助缩放参数来调节。

根据本发明的又另一方面提供该系统，其中，所述多个源中的至少一个源是运动的源。

根据本发明的又另一方面提供该系统，其中，信号是声学信号。

根据本发明的又另一方面提供该系统，其中，所述多个源中的至少一个源与发电机关联。

附图说明

图1和图3A-3F示出了根据本发明的一个方面的系统设置；

图3A-3F示出了应用于根据本发明的一个方面产生信号的设备；

图4A-4D示出了根据本发明的一个方面的信号性质；

图4A-4D-图20示出了根据本发明的一个方面的结果；

图21示出了带有根据本发明的一个方面的处理器的系统；

图22示出了根据本发明的一个方面的应用；以及

图23-28示出了根据本发明的一个或多个方面提供的方法和系统的性能。

具体实施方式

在此，对于宽带源提供了基于贝叶斯范式的到达方向（DOA）估计方法和系统。该根据本发明的一个方面提供的方法以粒子滤波器逼近时频域中DOA的后验概率密度函数。与其它方法不同，不需要时间平均，因此可以追踪运动的源。在此提供的方法将新的低复杂度加权和正则化方案用于融合来自不同频率的信息，以及用于克服当很少传感器可用时的过度拟合问题。

到达方向（DOA）问题是在阵列信号处理中的重要问题并且与盲源分离（BSS）学科紧密关联。两个领域都已经产生丰富的研究和大量方法，毫无疑问地用于在如雷达、无线通信或语音识别那样不同的方面中的许多应用，如在“H.Krim和M.Viberg.所著的Two decades of array signal processingresearch:the parametric approach.Signal Processing Magazine,IEEE,13(4):67-94,1996年7月”中描述那样。DOA估计需要天线阵列并且探索传感器之间的到达时间差。窄带DOA算法将相移用于逼近在不同麦克风的噪声测量之间的到达时间差。这建立了在所接收的和所发射的信号之间的代数统计学关系。存在试图解决该问题的方法。其中大多是如在“R.Roy和T.Kailath所著的Esprit-estimation of signal parameters via rotational invariancetechniques.Acoustics,Speech and Signal Processing,IEEE Transactions on,37(7):984,1989年7月”中描述的ESPRIT或者在“R.Schmidt所著的Multipleemitter location and signal parameter estimation.Antennas and Propagation,IEEE Transactions on,34(3):276,1986年3月”中描述的MUSIC的变型，其利用如在“M.Viberg和B.Ottersten所著的Sensor array processing based onsubspace fitting.Signal Processing,IEEE Transactions on,39(5):1110-1121,1991年5月”中描述的子空间拟合技术并且对于计算解答来说是迅速的。诸如在“P.Stoica和K.C.Sharman所著的Maximum Likelihood Methods forDirection-of-Arrival Estimation.Acoustics,Speech and Signal Processing,IEEETransactions on,38(7):1132,1990年7月”中描述的最大似然性（ML）和在“J.Lasenby和W.J.Fitzgerald所著的A bayesian approach to high-resolutionbeamforming.Radar and Signal Processing,IEE Proceedings F,138(6):539-544,1991年12月”中描述的贝叶斯方法的统计学最优方法长久以来如在“JamesA.Cadzow所著的Multiple source localization–the signal subspace approach.IEEE Transations on Acoustics,Speech,and Signal Processing,38(7):1110-1125,1990年7月”中描述那样被认为是难以驾驭的，然而这些方法近来例如在“C.Andrieu和A.Doucet所著的Joint Bayesian model selection and estimationof noisy sinusoids via reversible jump mcmc.Signal Processing,IEEETransations on,47(10):2667-2676,1999年10月”和“J.Huang,P.Xu,Y.Lu和Y.Sun所著的A novel bayesian high-resolution direction-of-arrival estimator.OCEANS,2001.MTS/IEEE Conference and Exhibition,3:1697-1702,2001”中受到更多的关注。

许多用于宽带DOA的方法在时频域中表达。窄带假说于是对于各个子带或频率点（frequency bin）是有效的。不相干信号子空间方法（ISSM）计算这样的DOA估计，其同时满足在所有子带中信号和噪声子空间正交性条件。另一方面如在“H.Wang和M.Kaveh所著的Coherent signal-subspaceprocessing for the detection and estimation of angles of arrival of multiplewide-band sources.Acoustics,Speech and Signal Processing,IEEE Transactionson,33(4):823,1985年8月”中描述的相干信号子空间方法（CSSM）根据所有数据计算泛空间协方差矩阵（universal spatial covariance matrix）（SCM）。任意窄带信号子空间方法于是可以用于分析泛SCM和估计DOA。

然而，良好的初始DOA估计对于如在“D.N.Swingler和J.Krolik所著的Source location bias in the coherently focused high-resolution broad-bandbeamformer.Acoustics,Speech and Signal Processing,IEEE Transactions on,37(1):143-145,1989年1月”中描述那样正确地使子带SCM与泛SCM一致是必要的。例如在“Ta-Sung Lee所著的Efficient Wideband source localizationusing beamforming invariance technique.Signal Processing,IEEE Transactionson,42(6):1376-1387,1994年6月”中描述的BI-CSSM或者在“Yeo-Sun Yoon,L.M.Kaplan和J.H.McClellan所著的Tops:new doa estimator for widebandsignals.Signal Processing,IEEE Transactions on,54(6):1977,2006年6月”中描述的TOPS的方法被开发以减轻该问题。子空间方法将信号和噪声子空间的正交性用作优化的标准。然而，在数学上更吸引人的方法是使估计过程基于决定理论架构上，从而实现跨层优化。前提是计算DOA的后验概率密度函数（pdf），其可以通过采用粒子滤波器来实现。这种方法在“W.Ng,J.P.Reilly和T.Kirubarajan所著的A bayesian approach to tracking wideband targets usingsensor arrays and particle filters.Statistical Signal Processing,2003IEEEWorkshop on,第510-513页,2003”中采用，其中在时域中表达贝叶斯极大后验似然性（MAP）估计器。

作为本发明的一个方面，贝叶斯MAP估计器利用信号的时频表示来提供。时频分析的优点通过在盲源分离（BSS）领域中使用的、诸如在“ScottRickard,Radu Balan和Justinian Rosca所著的Real-time time-frequency basedblind source separation.In Proc.Of International Conference on IndependentComponent Analysis and Signal Separation(ICA2001),第651-656页,2001”中描述的DUET和在“Thomas Melia和Scott Rickard所著的Underdeterminedblind source separation in echoic environments using desprit,2007”中描述的DESPRIT的技术来示出。这些方法利用不相似的信号特点以将信号分离，并且在第二步骤中估计其DOA。例如，语音信号在时频域中通常不相交，并且可以对于每个频率点单独地计算单个DOA。功率加权的直方图将频率信息组合，并且聚类被用于找到真正的DOA。

作为本发明的一个方面提供的方法将新颖的、试探式的加权方案用于对于各频率将信息组合。粒子滤波器逼近DOA的后验密度，其然后被用于提取极大后验似然性估计。

首先，关于盲定位解释一个方法。关于到达方向还提供一个解释。

图1以俯视图示出了包括在测试（DUT）号码1到3下的设备的空间的平面图，其带有麦克风装置，包括麦克风1、2和3。DUT#2中有故障的子系统的位置通过利用其空气传播的声学发射进行的三角测量而被发现。在此根据本发明的一个方面提供的方法克服了时频正交性的条件并且并不强加任何在信号关联方面的条件。

根据本发明的至少一个方面提供了在定位方法中使用的信号模型。下面将提供源定位方法和该方法的粒子滤波器拓展。随后，还提供了模拟以及现场实验的结果。

问题表达

来自记录J个源信号s₁(t),…,s_j(t)的M个紧密间隔的麦克风的阵列的时间离散混合x₁(t),…,x_M(t)使得

$x_M\left(t\right)=\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{s}_j(t-{\mathrm{\&tau;}}_{\mathrm{mj}})---\left(1\right)$

在此，τ_mj是信号j在麦克风m处相对于麦克风1的延迟。麦克风线性地对齐但是不必以麦克风1与m之间的已知距离d_m相等地间隔。距离被认为相对于源距离而言是小的。尤其，源信号在麦克风之间的衰减被认为是可忽略的并且波前是平面的。延迟τ_mj然后遵循关系式τ_mj=δ_jd_m，其中δ_j是信号j相对于阵列轴线的按照米的延迟。注意到δ_j具有信号j的到达方向信息，并且因此δ_j被称作方向。使用加窗的离散傅里叶变换并且假设

对于每个离散频率ω_n都成立。在该情况下，在信号与记录之间的关系在时频（t-f）域中是线性的，并且可以在每个频率基础上来描述。

在利用傅里叶表达式的大写字母的情况下，该关系读作

$X_m({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)=\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{S}_j({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)e^{-j{\mathrm{\&omega;}}_n{\mathrm{\&delta;}}_j{d}_m}+V_m({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)---\left(2\right)$

其中，τ_k是窗中心，并且独立地和相同地分布（iid.）的零均值循环对称复高斯（CSCG）噪声V_n:=(V₁(ω_n,τ_k),...,V_M(ω_n,τ_k))^T被添加至模型测量噪声、开窗错误和模型不确定性。可以将（2）与“S.Rickard,R.Balan和J.Rosca所著的Real-time time-frequency based blind source separation.In Proc.OfInternational Conference on Independent Component Analysis and SignalSeparation(ICA2001),第651-656页,2001”中的等式（30）-（31）相比较。

一个目标是基于直至时间τ_k的所有可用测量，估计方向δ:＝(δ₁，...，δ_J)^T，其中，源的数目J是已知的。并未在源特性

和噪声方差方面做出假设。

方法概要

首先要说的是对于每个时间点τ_k独立的估计问题和对于方向δ_j构造对数似然性函数，作为所有频率点ω_n的对数似然性的加权和。分数系数将被用作权重来选择带有突出的方向特性的频率点，参见更下方。在第二步骤中，在并非所有信号都活动的频率处的过度拟合问题将通过引入正则化方案来减轻。将集成时间维度并且将提出（非线性）滤波问题，其将利用序贯重要性重采样（SIR）粒子滤波器来解决。

对数似然性

在通过等式（2）限定的模型下，如果以Sn和δ为条件，则观察X_n:＝(X₁(ω_n，τ_k)，...，X_m(ω_n，τ_k))^T是iid.CSCG随机变量，并且所有M个观察的联合概率密度函数（pdf）分解为边缘概率密度函数。（参见“R.G.Gallager所著的Circularly-symmetric gaussian random vectors.Class Notes,2008年1月”，来回顾CSCG随机变量）。

对于对数似然性于是有

$J-\log p\left(X_n\right|S_n,\mathrm{\&delta;})\&Proportional;\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}{|X_M-\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{S}_j{e}^{-i{\mathrm{\&omega;}}_n{\mathrm{\&delta;}}_j{d}_m}|}^2=:L_n(\mathrm{\&delta;},S_n)---\left(3\right)$

其中，t-f索引(ω_n,τ_k)从X_m(ω_n,τ_k)和S_j(ω_n,τ_k)得出。负的对数似然性函数根据下面描述的成本函数架构而使用。如果仅考虑单个频率，则可以通过将和（成本函数L_n）最小化来获得对于δ和S_n的ML估计。如果样本大小M小，则自由度的高数目J可以导致过度拟合。参数δ_j在所有频率中共享，并且其ML估计可以通过解决全局优化问题来找到，其将所有频率组合，或者等效地通过将所有频率的联合似然性最大化。全局对数似然性函数通过在每个频率的对数似然性求和而获得，因为所有噪声项被假设为iid：

$L(\mathrm{\&delta;},S_1,...,S_N)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}{L}_n(\mathrm{\&delta;},S_n)---\left(4\right)$

揭示信号

上面所讨论的方法面对变为双重的问题：并不是所有源在所有频率上都活动。首先，可以存在这样的频率，在这些频率上没有源是活动的却接收大量能量。这可以归因于麦克风特性，或者考虑到小的样本大小M而是由离群值造成的。一个意图是丢弃在其上并未出现信号的频率点。第二，如果信号能量显著不同或者如果并非所有源都在给定的频率上活动，则小的采样尺寸倾向于过度拟合，并且算法可能递送误导性结果。在这种频率上，所有J个可用自由度在估计源信号中的利用通过引入正则化方案而受制约。

似然性的两个修改、频率点的丢弃和正则化建立在系数C^opt上，其被限定为可以被源在方向δ上解释的、所接收的能量的最大量。现在限定的是通过固定的角度/信号对（δ,S）实现的能量减少，

$c(\mathrm{\&omega;},\mathrm{\&delta;},S)=1-\frac{{\mathrm{\&Sigma;}}_{m=1}^M{|X_m-S{e}^{-\mathrm{i\&omega;\&delta;}{d}_m}|}^2}{{\mathrm{\&Sigma;}}_{m=1}^M{\left|X_M\right|}^2}---\left(5\right)$

该系数稍后将被用于评估估计。通过在S上最大化，对于给定方向获得最大可实现的减少

该系数涉及在模型选项中使用的贝叶斯系数（Bayes Factor），在此对带有方向δ的单信号模型和无信号模型的假设加权。（参见“C.P.Robert所著的The Bayesian Choice.Springer,第二版,2001”第7章，用于回顾模型选项和贝叶斯系数。）如果c^opt(ω，δ)，小，则对于方向δ仅能实现小的成本函数减小。因此，该频率点对于估计来自该方向的源而言帮助有限。另一方面，如果c^opt(ω，δ)接近1，则测量值可以被匹配于来自方向δ的信号。换言之，好像信号出现并且在方向δ上，并且该频率可以被用于估计。然而，必须小心，因为即使没有信号出现，一些成本函数减小还是可能的。事实上，在假设没有信号出现、即X_m=V_m(H₀假说)的情况下，则可实现的平均百分数成本函数减小可以被计算为

快速验证的是， $S=\frac{1}{M}{\mathrm{\&Sigma;}}_{m=1}^M{W}_m,$ 其中 $W_m=V_m{e}^{\mathrm{i\&omega;\&delta;}{d}_m},$ 解决该MMSE问题。因此，利用W_m的独立性，

正则化的对数似然性

可实现的百分数成本函数减小系数现在用于创建正则化的对数似然性函数。一个意图是着眼于大信号能量差距的问题：高能频率点对全局似然性函数的影响大，并且将影响所有信号的方向估计，即使仅单个信号出现在该频率上。为了减轻这些所不期望的交叉影响，惩罚λ(ω_n，δ_j)用于信号功率|S_j(ω_n，τ_k)|²

$\mathrm{\&lambda;}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j)=\mathrm{\&alpha;}(\underset{j^{\&prime;}}{\max }{c}^{\mathrm{opt}}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_{j^{\&prime;}})-c^{\mathrm{opt}}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j))---\left(9\right)$

其中α>0是可以适配于该问题的设计参数。通过这样做，人们偏好如下方向：所述方向相对于不具有高的可实现的百分数成本函数减小的方向来说具有高的可实现的百分数成本函数减小，并且这证明为在信号具有在频域中不同的特点的情况下良好地成立。

将该惩罚添加至对数似然性（3）产生正则化的（负的）对数似然性

$L_n^{\mathrm{reg}}(\mathrm{\&delta;},S_n)=L_n(\mathrm{\&delta;},S_n)+\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}\mathrm{\&lambda;}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j){\left|S_j({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)\right|}^2---\left(10\right)$

从统计学观点来看，正则化还可以被理解为利用对于信号功率的非平坦指数先验（non-flat exponential prior）；即，

$p\left({\left|S_j({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)\right|}^2\right)=\mathrm{\&lambda;}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j)e^{-\mathrm{\&lambda;}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j){\left|S_j({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)\right|}^2}{1}_{[0,\&infin;)}\left({\left|S_j({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)\right|}^2\right)---\left(11\right)$

并且正则化的对数似然性理解为S_n的以δ为条件的后验密度。1_[0，∞)在此是[0，∞)的指标函数。

全局对数似然性

接下来，将提供用于δ和S₁,…S_N的成本函数或者全局对数似然性，其是用于估计器的主要构件模块。等式（4）将被修改为包括来自上面的c^opt系数作为加权系数w(ω_n，τ_k)，

$\hat{w}({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k,\mathrm{\&delta;})=\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{c}^{\mathrm{opt}}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j)---\left(12\right)$

$w({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k,\mathrm{\&delta;})=\frac{\hat{w}({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k,\mathrm{\&delta;})}{{\mathrm{\&Sigma;}}_{n=1}^N\hat{w}({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k,\mathrm{\&delta;})}---\left(13\right)$

将这些权重与用于每个频率的正则化的对数似然性组合产生全局对数似然性

$L(\mathrm{\&delta;},S_1,...,S_N)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}w({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k,\mathrm{\&delta;})L_n^{\mathrm{reg}}(\mathrm{\&delta;},S_n)---\left(14\right)$

如果考虑实时应用，则该J+2NJ维的成本函数的直接最小化在计算上是惊人的。然而，对于固定的δ，每个加数是S_n的凸函数，并且可以明确地找到最小点（见下文）。余数是δ的正弦函数的叠加，这些正弦函数的评估开销大、非凸并且可能具有多个最小点。

粒子滤波器实现方案

为了将时间维度包括到估计策略中，构造了马尔科夫过程随机动态系统；即，用于估计变量δ和S_n的状态等式p(δ_k，S_1，k，...，S_N，k|S_1，k-1，...，S_N，k-1)，其中k索引时间。首先假设对于所有的n、k、l而言δ_k与S_n,l无关，并且然后立即将S_n从该系统丢掉，因为并未作出关于源信号的时间结构的任何假设。将S_n从该系统丢掉在形式上可以通过利用p(S_1，k，...，S_N，k|S_1，k-1，...，S_N，k-1)∝1来实现。尽管该先验不是真正的pdf，但它并不造成问题，因为似然性是表现良好的。

粒子滤波器被用于获得p(δ_k，S_1，k，...，S_N，k|X_1，0∶k，...，X_N，0∶k)的离散逼近，其利用来自时间τ₀到τ_k的所有可用数据。该pdf的估计然后可以用于获得δ_k和S_1，k，...，S_N，k的估计。当新的测量值变得可用时粒子滤波器迭代地将p(δ_k，S_1，k，...，S_N，k|X_1，0∶k，...，X_N，0∶k)更新，并且仅当逼近变差时才使用带有重采样的采样重要性重采样（SIR）变型。粒子的有效数目的估计被用作标准。对于粒子滤波器的指导性纵览可以在【1】S.Arulampalam,S.Maskell,N.Gordon和T.Clapp所著的A tutorial on particle filters for on-linenon-linear/non-gaussian bayesian tracking.IEEE Transactions on SignalProcessing,50:174-188,2001中找到。

转移核和重要性采样

用于源位置的转移概率服务于两个目的：描述源位错和在状态空间中查找可能的源。假设源运动是独立的，并且，因此利用对于各个状态位置没有漂移的1维随机游动；即，

$p\left({\mathrm{\&delta;}}_k\right|{\mathrm{\&delta;}}_{k-1})=\underset{j=1}{\overset{J}{\mathrm{\&Pi;}}}p\left({\mathrm{\&delta;}}_{j,k}\right|{\mathrm{\&delta;}}_{j,k-1})=\underset{j=1}{\overset{J}{\mathrm{\&Pi;}}}N({\mathrm{\&delta;}}_{j,k-1},{\mathrm{\&sigma;}}^2)---\left(15\right)$

其中，N(μ，σ²)用于表示带有均值μ和方差σ²的高斯随机变量的pdf。粒子滤波器理想地从p(δ_k|δ_k-1，X_k)∝p(X_k|δ_k)p(δ_k|δ_k-1)采样以更新其pdf估计。然而，该分部太难以采样。因此，使用重要性采样，其中p(δ_k|δ_k-1)作为重要性密度并且p(X_k|δ_k)用于更新粒子权重

$w_k^i=w_{k-1}^i\exp (-\mathrm{\&gamma;}{L}^{\mathrm{reg}}({\mathrm{\&delta;}}_k^i,S_{1,k}^i,...,S_{N,k}^i))---\left(16\right)$

γ>0是设计参数，其与测量噪声方差关联。如果粒子的有效数目

降至阈值N_th以下，则进行重采样。

信号估计

因为信号转移核被限定为平的（flat），所以在时间步骤之间未传送信息。从而，ML估计器被用于信号（或者如果正则化方案被视为先验，则使用MAP估计器）。尽管对于定位不需要信号估计，还是意欲使用明确的估计而不是边际化，以避免高成本的J维复杂积分。亏得高斯分布的平方指数，所以可以解析地获得ML/MAP估计。所具有的是

$e^{\mathrm{i\&omega;}{\mathrm{\&delta;}}_j{d}_m},\begin{array}{c}{L}_n^{\mathrm{reg}}(\mathrm{\&delta;},S_n)=\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}{|X_m-\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{S}_j{e}^{-\mathrm{i\&omega;}{\mathrm{\&delta;}}_j{d}_m}|}^2+\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}\mathrm{\&lambda;}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j){\left|S_j\right|}^2\\ =\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}({\left|X_m\right|}^2+\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}\underset{j^{\&prime;}=1}{\overset{J}{\mathrm{\&Sigma;}}}{S}_j^{*}{S}_{j^{\&prime;}}{e}^{-\mathrm{i\&omega;}({\mathrm{\&delta;}}_{j^{\&prime;}}-{\mathrm{\&delta;}}_j)d_m}-\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}(X_m^{*}{S}_j{e}^{-\mathrm{i\&omega;}{\mathrm{\&delta;}}_j{d}_m}+X_m{S}_j^{*}{e}^{\mathrm{i\&omega;}{\mathrm{\&delta;}}_j{d}_m}))\\ +\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}\mathrm{\&lambda;}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j){\left|S_j\right|}^2\end{array}---\left(17\right)$

限定了

${\left(H\right)}_{\mathrm{ij}}=\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}{e}^{-\mathrm{i\&omega;}({\mathrm{\&delta;}}_j-{\mathrm{\&delta;}}_i)d_m}---\left(18\right)$

${\left(\mathrm{\&Lambda;}\right)}_{\mathrm{ij}}=\left\{\begin{array}{cc}\mathrm{\&lambda;}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_i)& i=j\\ 0& i\&NotEqual;j\end{array}\right.---\left(19\right)$ $g_i=-\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}{X}_m{e}^{\mathrm{i\&omega;}{\mathrm{\&delta;}}_i{d}_m}---\left(20\right)$

（17）可以按照矩阵向量符号写做

$L_n^{\mathrm{reg}}(\mathrm{\&delta;},S_n)=X^{H}X+S^H(H+\mathrm{\&Lambda;})S+g^{H}S+S^{H}g---\left(21\right)$

并且可以找到如下最小点

$\hat{S}=-{(H+\mathrm{\&Lambda;})}^{-1}g---\left(22\right)$

在S中明显为凸的，因此可以通过求解

来找到最小值，其中

是

关于S^H的Wirthinger导数；即 $\&dtri;L_n^{\mathrm{reg}}(\mathrm{\&delta;},S)=(H+\mathrm{\&Lambda;})S+g=0.$

该信号估计可以用于利用（16）来更新粒子权重。

估计评估

为了改进粒子滤波器结果的可判读性和为了精细调节滤波器，根据本发明的一个方面提供了下面的评估方案。各个源位置被分配有分数c(δ_j)，其量化对于在正确位置追踪信号的信任。c(δ_j)的计算遵循（6）中c^opt(ω，δ)的计算，并且替代对于每个方向δj找到最优信号而使用来自（22）的ML估计为了集中于那些对于相应方向相关的频率，根据本发明的一个方面将系数c^opt(ω，δ_j)和所接收的能量E(ω_n，τ_k)用作权重κ。必须对于每个粒子计算c(δ_j)系数。

$E({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)=\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}\left|X_m({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)\right|---\left(23\right)$

${\mathrm{\&kappa;}}_{n,j}^i=\frac{E({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)c^{\mathrm{opt}}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j^i)}{{\mathrm{\&Sigma;}}_{n=1}^{N}E({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k)c^{\mathrm{opt}}({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j^i)}---\left(24\right)$

$c\left({\mathrm{\&delta;}}_j^i\right)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&kappa;}}_{n,j}^{i}c({\mathrm{\&omega;}}_n,{\mathrm{\&delta;}}_j^i,S_j^i({\mathrm{\&omega;}}_n,{\mathrm{\&tau;}}_k))---\left(25\right)$

修改过的转移核

为了协调转移核p(δ_j，k|δ_j，k-1)的两个目的（状态空间查找和源错位），使用取决于上面描述的评估措施的转移方差σ²。当源被追踪时小的方差σ²是优选的，当源没有被追踪时大的方差是优选的。因此，下面的表达式用于取代（15）中的σ²并且σ_max ²作为设计参数：

${\mathrm{\&sigma;}}^2=({\mathrm{\&delta;}}_{k-1},S_{1,k-1},...,S_{N,k-1})=(1-c\left({\mathrm{\&delta;}}_{j,k-1}\right)){\mathrm{\&sigma;}}_{\max }^2---\left(27\right)$

这将在形式上破坏关于S_n和δ所作的独立性假设。

方法总结

根据本发明的一个方面提供的一个定位方法可以通过如下步骤来描述：

方法1：

将粒子状态

j=1,...,J初始化

·for k=1,.,.do:

-for i=1,…,N_粒子do:

*根据（15）生成建议

*根据（6）计算

*根据（9）计算正则化系数

*根据（22）计算信号ML/MAP估计

*根据（3）到（10）计算似然性

*根据（14）计算全局似然性

*根据（23）计算评估系数

*根据（16）更新粒子权重wi。

-end for

-if N_eff<N_th

*重采样粒子

-end if

*end for

模拟结果

如在此根据本发明的一个方面提供的方法的发展受到其在真实数据集上的性能影响。然而，在进入到真实数据之前，首先检验该方法在利用人工生成的数据所控制的环境中的性能。图2a图示了现场中锅炉馈送泵（BFP）室的平面规划图。添加的测试源和麦克风的位置以及功率设备例如BFP#2和BFP#3的活动被示出。图2b示出了带有扬声器源的记录设置#1。设备通过箭头来指出。在空间后方的两个BFP是活动的并且生成高的宽带噪声。

人工数据

为了在理想条件上尝试该方法；即，根据按照本发明的一个方面提供的系统模型，生成了两个测试场景。对于第一场景，将确定性的正弦信号的两个不同的混合用作源。这些信号在频域中部分地交叠。源位置随时间变化并且在模拟的后一半中相交。使用三个不同的信噪比。第二场景由发射相同信号的三个源组成，其又是正弦的确定性混合。在此，源缓慢运动但并不交叉。

结果的讨论

该方法成功分离类似功率的信号，即使将噪声添加给麦克风。然而，在信号与信号的比是2的情况下，弱信号的估计被影响。如果信号与信号的比继续增大，则无法获得有效结果。

真实世界测量

进行测试记录并且该数据被用于对来自现场的应用相关的数据来评估上述方法。在下文中，将讨论所记录的数据和当前定位方法的性能。在带有两个不同的麦克风阵列的锅炉馈送泵室中取得这些记录。第一阵列包括三个非线性间隔的超声麦克风（4Hz–80kHz）。麦克风中心之间的间隔分别是7mm和8mm。第二阵列使用带有在相邻单元的中心之间7mm间隔的4个线性间隔的麦克风（40Hz–16kHz）。不同的声信号在接近功率设备的7个位置处被生成以模拟不正常（out-of-normal）发射。实验的设置在图2a中示出。注意到BFP#12、#13和一些辅助设备在所有实验期间运行并且生成强烈的背景噪声。记录场景由图2b进一步明示。在此测试位置1的记录设置以扬声器源示出。

图3提供了在实验中使用的不同记录设备（上一行）的视图。此外，图3还示出了在记录中使用的不同的源。利用扬声器发射带有限定特性的大量可重复的人造信号。这些测试信号是啁啾声（chirp）、一组频率不同的正弦信号、宽带白噪声、高频噪声和所记录的来自损坏的引擎的噪声。前四个测试信号的时频表示在图4中示出。引擎缺陷信号由于受限制的频率范围而并未在分析中使用。圆锯和锤子在金属管上的冲击用于生成附加于扬声器的高能信号。尽管有可重复的扬声器信号还是增加源的原因是由于一些测试信号的无法听到性而造成的高的背景噪声强度。不同的测试记录在表1中列出。注意到的是圆锯仅与PSC麦克风阵列一起使用，因为其发射在该吵杂环境中仅刚刚能够听见。

总地来说，进行了95个实验而形成1.65GB的声学数据。一个焦点在于这些记录的子集，因为其大的数目和通过穷尽的评估而在信息方面有限的收获。下面对于正弦信号的不同的测试位置示出结果。还对于在测试位置2使用的信号类型提供了结果。

源定位结果

在下文中，示出了带有正弦信号和不同测试位置的实验的结果。选择正弦测试信号，因为其具有最高信噪比（SNR）。实验结果归档在图5-20中并且将在下文中予以概括。这些图中的每个都以曲线示出了粒子滤波器结果。标记曲线表示早先描述的分数函数值c(δ_j)。所接收的混合信号的频谱图也在结果图中示出。还将示出被追踪的信号的所估计的频谱图。注意到的是，不同的测试位置在图2a中示出。图4示出了测试信号的时频绘图：（a）啁啾声测试信号，（b）不同频率处的正弦信号，（c）宽带噪声和（d）窄带高频噪声。

结果的讨论

第一组四个源以近似-60、-50、-5和40度的角度位于机器前方，如通过结果图中虚线所表示那样。位置2和3距离麦克风阵列最近并且实现最高的SNR；因此引起最精确的结果。扬声器在第2和3秒之间产生10KHz处的最清楚的信号。在该时间间隔期间，信号可以在位置1-3正确定位。在2、5、15和20kHz的一些正弦信号在位置2和3被检测到。距麦克风阵列最近的位置3在整个信号持续时间中被正确检测。除了小的距离之外，活动的机器（BFP#12）之一被放置在位置3后方，其辅助信号检测。在所有记录中，BFP#13在大约60度处被持续检测。位置5-7在机器后方并且在麦克风阵列处的SNR对于在此提供的、用于检测合适信号的方法而言过低。

在下面的表1中提供对带有两个麦克风阵列的实验的总结。

下面提供基于不同信号利用所提供的方法的结果的纵览。

以人工数据得到的模拟结果

图5：利用正弦信号的不同组合（底部左方）以两个源进行的模拟。带有指示了真实源位置的虚线的粒子滤波器结果曲线（虚线），以及用于估计结果（顶部）以及所估计的信号（底部右方）的曲线。

图6：利用正弦信号（底部左方）和添加的噪声（其中在所有麦克风处

的不同组合以两个源进行的模拟。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）和所估计的信号（底部右方）。

图7：利用正弦信号（底部左方）和添加的噪声（其中在所有麦克风处

并且

的不同组合以两个源进行的模拟。粒子滤波器结果，其带有真实源位置以及估计结果（顶部）和所估计的信号（底部右方）。

图8：利用相同的正弦信号（顶部左方）以三个源进行的模拟。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部右方）和所估计的信号（底部）。

用于带有正弦信号的记录的结果

图9：测试位置1，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图10：测试位置2，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图11：测试位置3，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图12：测试位置4，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及用于估计结果（顶部）的周期（circles）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图13：测试位置5，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及用于估计结果（顶部）的周期。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图14：测试位置6，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图15：测试位置7，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

用于带有不同信号类型的记录的结果

图16：测试位置2，啁啾声信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图17：测试位置2，正弦信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图18：测试位置2，噪声信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图19：测试位置2，高频噪声信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

图20：测试位置2，锤子在金属上的信号。粒子滤波器结果，其指示真实源位置（虚线）以及估计结果（顶部）。测量值的平均频谱图（底部左方）以及第一和第二源在线性刻度上的所估计的频谱图（底部中间和右方）。

比较来自不同测试信号的结果

为了更好理解和改进在此根据本发明的一个方面提供的方法，采用宽范围的不同源信号。分段地频率恒定的正弦波被用于检查是否可以检测窄带的高SNR信号。此外，使用啁啾声信号，其具有与上面的正弦波相同的特性，但是以缓慢地随时间变化的频率定中心。利用宽带、带有直至24kHz的份额的几乎白噪声信号和带有在20与25kHz之间的份额的高频噪声信号来评估均衡的信号功率和低SNR的情况。高频噪声信号由于扬声器的限制而具有低SNR。最后，为了生成高SNR信号，反复地用钢锤击打金属管。

讨论结果

啁啾声和正弦信号可以由于其高SNR而对于信号持续时间的大部分被正确检测。宽带噪声信号激起对位置估计的可见影响。来自锤的信号对金属管的影响在锤打击能听见的次数期间具有高SNR。在这些时间实例处可以正确定位锤。

总言之，在此根据本发明的一个方面提供的用于盲分离相关源的方法着眼于之前和已知的用于盲分离源的方法的缺点。所提供的方法可以利用非线性的和/或非线性间隔的麦克风阵列来取消在所估计的源位置方面的解疑（disambiguate）。此外，其实现同时定位多个信号和克服已知方法的对于信号在时频域中不相交或在统计学上独立的假设。最后，所实施的方法自动分离预测的源信号。它们和它们的所估计的源位置可以被进一步分析以例如检测异常声和其起源的原因。在此提供的方法的源定位结果示出了关于合成的和真实世界数据两者的改进。

在本发明的另一实施例中，通过假设信号随时间的相干性而增大源定位和分离的鲁棒性。在本发明的又另一实施例中，将对于每个源-接收器对的通道估计用于处理在回声环境中的问题。

如在此提供那样的方法在本发明的一个实施例中在系统或计算机设备上实施。使在图21中示出和如所提供那样的系统能够接收、处理和生成数据。该系统被提供有可以存储在存储器1801上的数据。数据可以从诸如麦克风或麦克风阵列的传感器获得，或者可以从数据源来提供。数据可以在输入端1806上提供。还提供有处理器或者将其以指令集合来编程，或者执行本发明的方法的程序存储在存储器1802上并且提供给处理器1803，其执行1802的指令以处理来自1801的数据。诸如源定位数据或者所估计的源信号或者任意其它得自处理器的数据可以在输出设备1804上输出，其可以是用于显示数据的显示器或者是用于提供声学信号的扬声器。然而，由处理器生成的数据还可以被用于实现信号的进一步处理，包括生成来自源的、所估计的信号的分类。处理器还可以具有通信通道1807用于从通信设备接收外部数据和用于将数据传送给外部设备。该系统在本发明的一个实施例中具有输入设备1805，其可以是传感器、麦克风、键盘、鼠标或可以生成待提供给处理器1803的数据的任意其它设备。该处理器可以是专用硬件。然而，处理器还可以是CPU或者可以执行1802的指令的任意其它计算设备。从而，如在图21中示出那样的系统提供一种系统，其用于得自传感器或其它数据源的数据的处理并且能够执行如在此作为本发明的一个方面提供的方法的步骤。

因此，在此描述了用于相关源的盲定位的系统和方法，例如在用于检测生成信号的机器中的异常操作的系统中。这种系统在图22中示出。系统2200包括接收器2202和处理器2204。接收器2202具有多个麦克风并且接收由与机器相关的不同声源生成的声音信号。该声音信号发送给处理器2204。处理器2204可以具有在前一图中示出的形式。处理器2204处理根据在此描述的步骤接收的信号。根据处理的结果，处理器2204在输出端2206生成指示异常操作模式的信号。

关于到达方向的表达

在前面关于盲定位描述了本发明的方面。还可以关于到达方向问题表达本发明的不同方面。

考虑M个传感器的带有元件间隔d_m-d_m-1的线性阵列。撞击到该阵列上的是来自不同方向θ_j的J个未知波前。波前的传播速度是c。源的数目J被假设为已知并且J≤M。麦克风被假设为在源的远场中。在DFT域中，在第n个子带中第m个传感器处接收的信号可以被建模为：

$X_m\left({\mathrm{\&omega;}}_n\right)=\underset{j=0}{\overset{J-1}{\mathrm{\&Sigma;}}}{S}_j\left({\mathrm{\&omega;}}_n\right)e^{-i{\mathrm{\&omega;}}_n{y}_m\sin \left({\mathrm{\&theta;}}_j\right)}+N_m\left({\mathrm{\&omega;}}_n\right)---\left(28\right)$

其中，S_j(ω_n)是第j个源信号，N_m(ω_n)是噪声，并且

噪声被假设为在各个频率内是循环对称复高斯分布（CSCG）和iid的，即方差会随着ω变化。

可以使用下面的定义：

X_n＝[X₀(ω_n)...X_M-1(ω_n)]^T   (29)

N_n＝[N₀(ω_n)...N_M-1(ω_n)]^T   (30)

S_n＝[S₀(ω_n)...S_J-1(ω_n)]^T    (31)

等式（29）可以以矩阵向量符号的形式再写作：

X_n＝A_n(θ)S_n+N_n    (32)

其中M×J导引矩阵A_n(θ)：

A_n(θ)＝[a(ω_n，θ₀)，...，a(ω_n，θ_J-1)]   (33)

其列是M×1阵列流形：

$a({\mathrm{\&omega;}}_n,{\mathrm{\&theta;}}_j)={[1e^{-i{\mathrm{\&omega;}}_n{v}_1\sin \left({\mathrm{\&theta;}}_j\right)}...e^{-i{\mathrm{\&omega;}}_n{v}_{M-1}\sin \left({\mathrm{\&theta;}}_j\right)}]}^T---\left(34\right)$

和

θ＝[θ₀，...θ_J-1]   （35）

子空间方法

用于解决DOA问题的最常用的方法从所接收数据的样本协方差矩阵计算信号和噪声子空间，并且选择这样的θ_j，其对应的阵列流形a(θ_j)最接近信号子空间，即，

${\widehat{\mathrm{\&theta;}}}_j=\underset{\mathrm{\&theta;}}{\arg \min }a{\left(\mathrm{\&theta;}\right)}^H{E}_N{E}_N^{H}a\left(\mathrm{\&theta;}\right)---\left(36\right)$

其中，E_N的列形成噪声子空间的正交基。刻意地略去了对ω_n的参考。不相干方法对于每个子带计算信号和噪声子空间E_N(ω_n)，并且θ_j被选择为平均起来满足（36）。相反，相干方法通过将数据从所有频率转换到参考频率来计算参考信号和噪声子空间。于是仅对于参考阵列流形校验正交条件（36）。这些其中CSSM和WAVES是两个代表的方法示出了比不相干方法显著更好的性能，尤其是对于高度相关和低SNR信号而言。但是，其需要良好的初始DOA估计以使聚焦过程有效，并且如何可以获得这些估计并不显然。

最大似然性方法

与子空间算法不同，ML方法计算来自A矩阵的信号子空间并且选择在将其在子空间上的投影最大化方面最优地匹配于所观察数据的θ₀，...，θ_J-1：

$\widehat{\mathrm{\&theta;}}=\underset{\mathrm{\&theta;}}{\arg \max }\left|\right|P_n\left(\mathrm{\&theta;}\right)X_n\left|\right|---\left(37\right)$

其中，P_n＝A(A^HA)^-1A^H是在由A_n(θ)的列张开的信号子空间上的投影矩阵。

如果噪声方差对于各频率相等，则可以通过对于各频率求和（37）来获得用于宽带问题的全面对数似然性函数。至今还未处理噪声方差对于各频率改变的问题。

在“C.E.Chen,F.Lorenzelli,R.E.Hudson和K.Yao所著的Maximumlikelihood doa estimation of multiple wideband sources in the presence ofnonuniform sensor noise.EURASIP Journal on Advances in Signal Processing,2008;Article ID835079,第12页,2008.doi:10.1155/2008/835079,2008”中，很详细地研究了如下情况：对于各传感器不一致的噪声却是对于各频率恒定的。与子空间方法相比，ML方法提供了关于阵列形式和信号相关性的更高灵活性，并且通常对于小样本尺寸示出更好性能，然而（37）中的非线性多维优化在计算上是复杂的。近来，对于窄带情况提出了重要性采样方法，以如在“Huigang Wang,S.Kay和S.Saha所著的An importance samplingmaximum likelihood direction of arrival estimator.Signal Processing,IEEETransactions on,56(10):5082-5092,2008”中描述那样有效地解决优化问题。在此处根据本发明的一个方面提供的方法中采用的粒子滤波器以相似方式着手处理非线性优化。

DOA估计

在模型（28）下，如果以S_n和θ为条件，则观察X₀，....，X_M-1是iid CSCG随机变量，并且联合pdf分解为边缘pdf。对于负的对数似然性，存在：

通常计算对于Sn的ML解为：

给定固定的θ并且通过将剩余集中对数似然性（

表示A的Moore-Penrose逆）最大化来对于θ进行求解。等式（39）是简单的最小二乘回归并且必须很注意数据的过度拟合问题，尤其是如果在一些频率上存在功率比其它信号显著更大的信号。

如果噪声方差σ_n ²已知，则可以通过将频率似然性求和来计算全局（负的）对数似然性：

$L(\mathrm{\&theta;},S_{0:N-1})=\underset{n=0}{\overset{N-1}{\mathrm{\&Sigma;}}}\frac{L_n(\mathrm{\&theta;},S_n)}{{\mathrm{\&sigma;}}_n^2}---\left(40\right)$

该条件函数前面已经陈述过，并且如在“James A.Cadzow所著的Multiple source localization–the signal subspace approach.IEEE Transactionson Acoustics,Speech,and Signal Processing,38(7):1110-1125,1990年7月”中描述那样被看做难以处理的。下面根据本发明的一个方面提供了一种粒子滤波器方法以解决对于多个瞬象的滤波问题，其附带地自然解决了优化问题。在实际应用中，良好的加权方案是必要的，如下面作为本发明的一个方面提供那样，以及改进性能的正则化方案是必要的。

加权

等式（40）中的噪声方差σ²无法从单个瞬象来估计。而是，噪声方差被重判读为加权系数τ_n，这是被例如为DUET的BSS算法所采用的观点。在实践中，信号带宽可能未被精确已知并且在一些频率点中J个信号的假设不适用。过度拟合问题在这些频率点中变得严重，并且将这些频率点包括在估计过程中会扭曲结果。提供下面的权重来解决不精确的建模、高噪声频率点和离群值频率点：

${\mathrm{\&tau;}}_n=\frac{{\widehat{\mathrm{\&tau;}}}_n}{{\mathrm{\&Sigma;}}_{n=0}^{N-1}{\widehat{\mathrm{\&tau;}}}_n}---\left(42\right)$

其中，

是非负的、非递减的加权函数。其自变量度量所接收信号的、可以在给定DOA矢量θ的情况下来解释的部分。τ_n是归一化的权重。

过度拟合

如上面指示的，在（37）中计算S_n的ML估计会经受过度拟合，尤其是当传感器的数目小或者J个信号假设不适用时。在岭回归中，对于估计变量引入惩罚项，并且在贝叶斯分析中这些惩罚项转换为用于S_n的先验分布。从而，作为本发明的一个方面的S_n的ML估计由在信号功率的指数先验下的MAP近似来取代。如果选择：

$p\left({\left|S_j\left({\mathrm{\&omega;}}_n\right)\right|}^2\right)=\mathrm{\&lambda;}{e}^{-\mathrm{\&lambda;}{\left|S_j\left({\mathrm{\&omega;}}_n\right)\right|}^2}{I}_{[0,\&infin;)}\left({\left|1S_j\left({\mathrm{\&omega;}}_n\right)\right|}^2\right)11---\left(43\right)$

其中，λ>0并且I是指标函数，在给定θ情况下用于S_n的负的对数后验可以写作：

$L_n^{\mathrm{reg}}(\mathrm{\&theta;},S_n)=L_n(\mathrm{\&theta;},S_n)+\underset{j=0}{\overset{J-1}{\mathrm{\&Sigma;}}}\mathrm{\&lambda;}{\left|S_j\left({\mathrm{\&omega;}}_n\right)\right|}^2---\left(44\right)$

粒子滤波器

用于S_n和θ的似然性函数可以根据等式（40）计算：

$p\left(X_{0:N-1}\right|\mathrm{\&theta;},S_{0:N-1})\&Proportional;e^{-\mathrm{\&gamma;L}(\mathrm{\&theta;},S_{0:N-1})}---\left(45\right)$

只有当在频率n处真实的噪声方差是

时其才是真正的似然性函数。对于下文假设这成立。

借助粒子滤波器将时间维度包括到估计过程中。以旧的后验pdf开始：

$p^{k-1}=p({\mathrm{\&theta;}}^{k-1},S_{0:N-1}^{k-1}|I^{k-1})---\left(46\right)$

其中，I^k-1表示直至时间实例k-1的所有信息（测量值），将贝叶斯规则和基础计算用于当新的测量值变得可用时将后验更新：

p^k∝l^kz^k，k-1p^k-1  (47)

带有似然性：

$l^k=p\left(X_{0:N-1}^k\right|{\mathrm{\&theta;}}^k,S_{0:N-1}^k)---\left(48\right)$

以及用于估计变量θ和S_n的一阶马可夫转移核：

$z^{k,k-1}=p({\mathrm{\&theta;}}^k,S_{0:N-1}^k|{\mathrm{\&theta;}}^{k-1},S_{0:N-1}^{k-1})---\left(49\right)$

假设源运动和信号是独立的，并且信号在时间步骤之间是独立的。还假设源独立地运动。因此转移核可以被分离：

$z^{k,k-1}=p\left(S_{0:N-1}^k\right)\&CenterDot;{\mathrm{\&Pi;}}_{j=0}^{J-1}p\left({\mathrm{\&theta;}}_j^k\right|{\mathrm{\&theta;}}_j^{k-1})---\left(50\right)$

如果粒子的有效数目降至预定阈值以下，则通过重采样着眼于粒子消耗的问题，如例如在“Sanjeev Arulampalam,Simon Maskell,Neil Gordon和Tim Clapp所著的A tutorial on particle filters for on-linenon-linear/non-gaussian Bayesian tracking.IEEE Transactions on SignalProcessing,50-174-188,2001”所述。

对于θ_j，将对于随机游动的如下修改用作转移密度：

$p\left({\mathrm{\&theta;}}_j^k\right|{\mathrm{\&theta;}}_j^{k-1})=\mathrm{\&alpha;}{U}_{[-\frac{\mathrm{\&pi;}}{2},\frac{\mathrm{\&pi;}}{2}]}+(1-\mathrm{\&alpha;})P_{[-\frac{\mathrm{\&pi;}}{2},\frac{\mathrm{\&pi;}}{2}]}N({\mathrm{\&theta;}}_j^{k-1},{\mathrm{\&sigma;}}_{\mathrm{\&theta;}}^2)---\left(51\right)$

其中，

表示在

上的均匀分布的pdf，并且 $P_{[-\frac{\mathrm{\&pi;}}{2},\frac{\mathrm{\&pi;}}{2}]}N({\mathrm{\&theta;}}_j^{k-1},{\mathrm{\&sigma;}}_{\mathrm{\&theta;}}^2)$ 表示投影在 $[-\frac{\mathrm{\&pi;}}{2},\frac{\mathrm{\&pi;}}{2}]{\mathrm{\&theta;}}^{\&prime;}=\mathrm{\&alpha;}\sin \left(\sin \left(\mathrm{\&theta;}\right)\right)$ 上的具有平均值和方差

的正态分布的pdf。使用根据“Yongtao Guan,Roland Fleissner,Paul Joyce和Stephen M.Krone所著的Markov chain monte carlo in smallworlds.Statistics and Computing,16:193-202,2006年6月”的小世界建议密度，因为这可能会加速收敛，尤其在该具有多模态似然性函数的情况中。

计算机模拟

在此根据本发明的一个方面提供的方法以三个不同的计算机模拟的场景来测试。在所有场景中，暂时发射高斯白噪声信号的等功率源位于M个传感器的线性阵列的远场中。任意两个不同的源信号与ρ∈[-1，1]关联；传感器经受iid高斯白噪声。带有矩形窗的L点FFT在传感器输出端上进行，并且在传感器通带f₀±Δf内的N个频率点上进行进一步的处理。SNR值是对于每个传感器的并且用于通带。在前两个场景中，在所有元素（ULA）之间的传感器间间隔是

其中

对于粒子滤波器，将均匀的初始分布用于所有P个粒子的θ_(i)，j。并且

用作加权函数。用于提取MAP估计的柱状图的频率点尺寸是Δ_h。将如下试探用于γ：

${\mathrm{\&gamma;}}^k=\frac{10}{{\mathrm{\&Sigma;}}_{i=0}^{P-1}L({\mathrm{\&theta;}}_{\left(i\right)}^k,S_{\left(i\right),0:N-1}^k}---\left(52\right)$

WAVES、CSSM、IMUSIC和TOPS基于当前的和Q个在前的模块来计算DOA估计。这允许在线动态计算并且在本质上与离线地处理完整的时间系列不同。WAVES和CSSM将如在“H.Hung和M.Kaveh所著的Focusingmatrices for coherent signal-subspace processing.Acoustics,Speech and SignalProcessing,IEEE Transactions on,36(8):1272-1281,1988年8月”描述的RSS聚焦矩阵用于使样本SCM与作为聚焦角度的真实角度符合。这是不现实的假设但是为相干方法提供了性能上边界。WAVES算法如在“E.D.di Claudio和R.Parisi所著的Waves:weighted average of signal subspaces for robustwideband direction finding.Signal Processing,IEEE Transactions on,49(10):2179,2001年10月”中描述那样实施，并且在此将Root-MUSIC用作用于CSSM和WAVES两者的窄带技术。

在下面的表2中列出三个场景的参数值。

表2

将p_all(Δθ)用作评估尺度，其反映其中所有源在真实DOA的±Δθ内被检测到的模块的百分数。每个模拟包括(B+Q)L个样本并且将评估p_all(Δθ)在最后B个模块上平均。所有的结果都是基于对于参数的各个组合的200次Monte Carlo运行的。

第一场景在“H.Wang和M.Kaveh所著的Coherent signal-subspaceprocessing for the detection and estimation of angles of arrival of multiplewide-band sources.Acoustics,Speech and Signal Processing,IEEE Transactionson,33(4):823,1985年8月”和“E.D.di Claudio和R.Parisi所著的Waves:weighted average of signal subspaces for robust wideband direction finding.Signal Processing,IEEE Transactions on,49(10):2179,2001年10月”中使用。

图23示出了场景1：对照用于源关联系数ρ的不同值的SNR示出了其中在2°内检测到所有源的模块的百分数。与CSSM相反，粒子滤波器可以在不同的SNR上区分相关信号。图23中的结果示出了基于粒子滤波器的方法可以以低SNR值和对于任意关联来分辨间隔紧密的信号。相反用于CSSM的分辨率阈值随着关联而增大（并未应用空间平均，因为这约束可能的阵列几何结构）。WAVES结果与CSSM的结果几乎相同并且为了清楚性而并未示出。不相干的方法都未在该低SNR场景中成功分辨所有四个源。该曲线图示出了在对应于约4秒的动态场景平均长度Q中所提出的算法优于CSSM，即使当使用完全聚焦时。

场景2在图24中示出。对照ρ=0（直线）和ρ=0.1（虚线）示出了p_all(2.5°)。对于第二场景（参见图24），使用与声音信号有关的参数。为了实现在双核笔记本电脑上的实时计算，粒子滤波器方法仅使用上10%的高能频率点（N=47）。用于MAP提取的柱状图的频率点尺寸是3°。这些结果示出了，利用仅10%的可用信息的、带有很少粒子的次优选粒子滤波器是CSSM的有效替选方案（再次，WAVES结果近乎相同并且未示出）并且优于IMUSIC。

在图25中示出了场景3，在曲线图中对照时间示出了源位置。在该第三场景中，如在图25中所示，示出了在此提供的方法的用于追踪运动的源的潜力。使用带有坐标d₀，...，d_M-1和距离d_m-d_m-1＝d+Δd(其中

并且

）的M=4个传感器的非均匀的、线性的阵列。信号集中在信号通带 $[f_0-\mathrm{\&Delta;}{f}_{\mathrm{SRC}},,f_0+\mathrm{\&Delta;}{f}_{\mathrm{SRC}}]\&Subset;[f_0-\mathrm{\&Delta;f},f_0+\mathrm{\&Delta;f}]$ 中（其中Δf_SRC=20Hz），并且使用0dB总信号功率对总噪声功率的SNR。当前的方法在运动的源交叉之前和其交叉之后成功估计正确的源位置。

在本发明的一个实施例中，在此提供的方法包括在“C.Andrieu和A.Doucet所著的Joint Bayesian model selection and estimation of noisy sinusoidsvia reversible jump mcmc.Signal Processing,IEEE Transactions on,47(10):2667-2676,1999年10月”的激发中利用可逆的跳跃马可夫模型来估计源的数目。这将阻止当源交叉时的“鬼影”估计。

本发明的方面还在下面的描述中予以解释。

下文提供了对用于宽带源的到达方向估计的进一步描述，其基于根据本发明的一个方面的贝叶斯范式。所提供的方法以粒子滤波器逼近在时频域中的DOA的后验概率密度函数，该方法的一个实施方案命名为缩写“MUST”。与其它早先的方法相反，不需要时间平均，因此可以追踪运动的源。在此提供的方法将新的低复杂度加权和正则化方案用于融合来自不同频率的信息和用于克服当仅很少传感器可用时的过度拟合问题。

如在“H.Krim和M.Viberg所著的Two decades of array signal processingresearch:the parametric approach.Signal Processing Magazine,IEEE,13(4):67-94,1996年7月”中描述那样，数十年的研究形成了解决到达方向（DOA）估计问题的许多方法，并且这些方法在类似雷达、无线通信或语音识别的领域中得到应用。DOA估计需要传感器阵列并且探索传感器之间的到达时间差。窄带方法以相移逼近这些差。用于该问题的现有方法中的大多数是如在“R.Roy和T.Kailath所著的Esprit-estimation of signal parameters viarotational invariance techniques.Acoustics,Speech and Signal Processing,IEEETransactions on,37(7):984,1989年7月”中描述的ESPRIT或者在“R.Schmidt所著的Multiple emitter location and signal parameter estimation.Antennas andPropagation,IEEE Transactions on,34(3):276,1986年3月”中描述的MUSIC的变型,其利用如在“M.Viberg和B.Ottersten所著的Sensor array processingbased on subspace fitting.Signal Processing,IEEE Transactions on,39(5):1110-1121,1991年5月”中描述的子空间拟合技术并且对于计算解答来说是迅速的。一般而言，基于子空间的方法的性能随着信号关联而劣化。诸如在“P.Stoica和K.C.Sharman所著的Maximum Likelihood Methods forDirection-of-Arrival Estimation.Speech and Signal Processing,IEEETransactions on,38(7):1132,1990年7月”中描述的最大似然性（ML）和在“J.Lasenby和W.J.Fitzgerald所著的A bayesian approach to high-resolutionbeamforming.Radar and Signal Processing,IEE Proceedings F,138(6):539-544,1991年12月”中描述的贝叶斯方法的统计学最优方法长久以来被认为是难以驾驭的，如在“James A.Cadzow所著的Multiple source localization–thesignal subspace approach.IEEE Transations on Acoustics,Speech,and SignalProcessing,38(7):1110-1125,1990年7月”中描述那样，然而这些方法近来受到更多的关注，例如在“C.Andrieu和A.Doucet所著的Joint Bayesian modelselection and estimation of noisy sinusoids via reversible jump mcmc.SignalProcessing,IEEE Transations on,47(10):2667-2676,1999年10月”和“J.Huang,P.Xu,Y.Lu和Y.Sun所著的A novel bayesian high-resolution direction-of-arrivalestimator.OCEANS,2001.MTS/IEEE Conference and Exhibition,3:1697-1702,2001”中。

用于宽带DOA的方法通常在时频（t-f）域中表达。窄带假说于是对于各个子带或频率点是有效的。不相干信号子空间方法（ISSM）计算这样的DOA估计，其同时满足在所有子带中信号和噪声子空间正交性条件。另一方面如在“H.Wang和M.Kaveh所著的Coherent signal-subspace processing forthe detection and estimation of angles of arrival of multiple wide-band sources.Acoustics,Speech and Signal Processing,IEEE Transactions on,33(4):823,1985年8月”中描述的相干信号子空间方法（CSSM）根据所有数据计算泛空间协方差矩阵（universal spatial covariance matrix）（SCM）。任意窄带信号子空间方法于是可以用于分析泛SCM。然而，良好的初始估计对于如在“D.N.Swingler和J.Krolik所著的Source location bias in the coherently focusedhigh-resolution broad-band beamformer.Acoustics,Speech and SignalProcessing,IEEE Transactions on,37(1):143-145,1989年1月”中描述那样正确地使子带SCM与泛SCM一致来说是必要的。例如在“Ta-Sung Lee所著的Efficient Wideband source localization using beamforming invariancetechnique.Signal Processing,IEEE Transactions on,42(6):1376-1387,1994年6月”中描述的BI-CSSM或者在“Yeo-Sun Yoon,L.M.Kaplan和J.H.McClellan所著的Tops:new doa estimator for wideband signals.Signal Processing,IEEETransactions on,54(6):1977,2006年6月”中描述的TOPS的方法被开发出以减轻该问题。

子空间算法将信号和噪声子空间的正交性用作优化的标准。然而，在数学上更吸引人的方法是使估计基于决定理论架构。前提是计算DOA的后验概率密度函数（pdf），其可以通过采用粒子滤波器来实现。这种方法在“W.Ng,J.P.Reilly和T.Kirubarajan所著的A bayesian approach to trackingwideband targets using sensor arrays and particle filters.Statistical SignalProcessing,2003IEEE Workshop on,第510-513页,2003”中采用，其中在时域中表达贝叶斯极大后验似然性（MAP）估计器。

在此，贝叶斯MAP估计器利用信号的时频表示来提供。时频分析的优点通过在诸如在“Scott Rickard,Radu Balan和Justinian Rosca所著的Real-time time-frequency based blind source separation.In Proc.Of InternationalConference on Independent Component Analysis and Signal Separation(ICA2001),第651-656页,2001”中描述的DUET和在“Thomas Melia和ScottRickard所著的Underdetermined blind source separation in echoic environmentsusing desprit,2007”中描述的DESPRIT的BSS中使用的技术来示出。这些方法探索不相似的信号特点以将信号分离，并且对于语音信号良好适用。

所示出的多源追踪将新颖的、试探式的加权方案用于对于各频率将信息组合，所述多源追踪的一种实施方式或实施例称为MUST方法。粒子滤波器逼近DOA的后验密度，然后MAP估计被提取。将讨论一些广泛使用的其它方法，并且将提供对MUST的详细描述。将提供利用估计本发明的一个方面在此提供的MAP估计器得到的、和与CSSM、IMUSIC和如在“E.D.di Claudio和R.Parisi所著的Waves:weighted average of signal subspaces forrobust wideband direction finding.Signal Processing,IEEE Transactions on,49(10):2179,2001年10月”中描述WAVES相比较的模拟结果。

进一步的问题表达和相关工作

考虑带有在传感器1和m之间的、记做d_m的距离的M个传感器的线性阵列。撞击到该阵列上的是来自不同方向θ_j的J个未知波前。波前的传播速度是c。源的数目J被假设为已知并且J≤M。麦克风被假设为在源的远场中。在DFT域中，在第n个子带中第m个传感器处接收的信号可以被建模为：

$X_m\left({\mathrm{\&omega;}}_n\right)=\underset{j=0}{\overset{J-1}{\mathrm{\&Sigma;}}}{S}_j\left({\mathrm{\&omega;}}_n\right)e^{-i{\mathrm{\&omega;}}_n{y}_m\sin \left({\mathrm{\&theta;}}_j\right)}+N_m\left({\mathrm{\&omega;}}_n\right)---\left(55\right)$

其中，S_j(ω_n)是第j个源信号，N_m(ω_n)是噪声，并且

噪声被假设为在各个频率内是循环对称复高斯分布（CSCG）以及iid的，即方差会随着ω变化。如果定义：

X_n＝[X₀(ω_n)...X_M-1(ω_n)]^T和类似地定义N_n和S_n，则等式（53）可以以矩阵向量符号的形式再写作：

X_n＝A_n(θ)S_n+N_n   （54）

其中M×J导引矩阵：

A_n(θ)＝[a(ω_n，θ₁)，...，a(ω_n，θ_f)]  （55）

其列是M×1阵列流形：

$a({\mathrm{\&omega;}}_n,{\mathrm{\&theta;}}_j)={[1e^{-i{\mathrm{\&omega;}}_n{v}_1\sin \left({\mathrm{\&theta;}}_j\right)}...e^{-i{\mathrm{\&omega;}}_n{v}_{M-1}\sin \left({\mathrm{\&theta;}}_j\right)}]}^T---\left(56\right)$

子空间方法

用于解决DOA问题的最常用的方法根据所接收数据的样本协方差矩阵计算信号和噪声子空间，并且选择这样的θ_j，其对应的阵列流形a(θ_j)最接近信号子空间，即，

${\widehat{\mathrm{\&theta;}}}_j=\underset{\mathrm{\&theta;}}{\arg \min }a{\left(\mathrm{\&theta;}\right)}^H{E}_N{E}_N^{H}a\left(\mathrm{\&theta;}\right)---\left(57\right)$

其中，E_N的列形成噪声子空间的正交基。刻意地略去了对ω_n的参考。不相干方法对于每个子带计算信号和噪声子空间E_N(ω_n)，并且θ_j被选择为平均起来满足（57）。相反，相干方法通过将数据从所有频率转换到参考频率来计算参考信号和噪声子空间。于是仅对于参考阵列流形校验正交条件（57）。这些方法（其中CSSM和WAVES是两个代表）示出了比不相干方法显著更好的性能，尤其是对于高度相关和低SNR信号而言。但是，其需要良好的初始DOA估计以使聚焦过程有效，并且如何可以获得这些估计并不显然。

最大似然性方法

与子空间算法不同，ML方法计算来自A_n矩阵的信号子空间并且选择在将其在子空间上的投影最大化方面最优地匹配于所观察数据的

$\widehat{\mathrm{\&theta;}}=\underset{\mathrm{\&theta;}}{\arg \max }\left|\right|P_n\left(\mathrm{\&theta;}\right)X_n\left|\right|---\left(58\right)$

其中，P_n＝A(A^HA)^-1A^H是在由A_n(θ)的列张开的信号子空间上的投影算子。决定性的ML估计假设对这些信号并无认识。而如果信号统计已知，则可以计算随机ML估计。

如果噪声方差对于所有各频率相等，则可以通过针对各频率地进行求和（38）来获得用于宽带问题的全面对数似然性函数。至今未处理噪声方差变化的问题。在“C.E.Chen,F.Lorenzelli,R.E.Hudson和K.Yao所著的Maximum likelihood doa estimation of multiple wideband sources in thepresence of nonuniform sensor noise.EURASIP Journal on Advances in SignalProcessing,2008;Article ID835079,第12页,2008.doi:10.1155/2008/835079,2008”中，研究了如下情况：对于传感器不一致的噪声却是对于各频率恒定的。与子空间方法相比，ML方法提供了关于阵列形式和信号相关性的更高灵活性，并且通常对于小样本尺寸示出更好性能，然而（58）中的非线性多维优化在计算上是复杂的。近来，对于窄带情况提出了重要性采样方法，以如在“Huigang Wang,S.Kay和S.Saha所著的An importance samplingmaximum likelihood direction of arrival estimator.Signal Processing,IEEETransactions on,56(10):5082-5092,2008”中描述那样有效地解决优化问题。作为本发明的一个方面采用并且在作为本发明的一个方面提供的方法中、例如在MUST实施方案中采用的粒子滤波器沿着这些线进行优化。

多源追踪（MUST）

在模型（53）下，如果以S_n和θ为条件，则观察X₀，...，X_M-1是iid CSCG随机变量。因此，联合pdf分解为边缘pdf。对于ω_n，通过下式给出负的对数似然性：

-logp(X_n[S_n，θ)∝||X_n-A_n(θ)S_n)||²   (59)

通常将对于S_n的ML解计算为：

表示A_n的Moore-Penrose逆。θ的解通过对于固定的θ求解（60）和通过将剩余集中负对数似然性最小化来迭代地在该估计上进行改进而得到。

如果噪声方差σ_n ²已知，则可以通过将频率似然性求和来计算全局（负的）对数似然性：

$L\left(\mathrm{\&theta;}\right)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}\frac{L_n\left(\mathrm{\&theta;}\right)}{{\mathrm{\&sigma;}}_n^2}---\left(62\right)$

该条件函数前面已经陈述过，并且如在“James A.Cadzow所著的Multiple source localization–the signal subspace approach.IEEE Transactionson Acoustics,Speech,and Signal Processing,38(7):1110-1125,1990年7月”中描述那样被看做难以处理的（1990年）。下面根据本发明的一个方面提供了一种粒子滤波器方法以解决对于多个瞬象的滤波问题，其附带地自然解决了优化问题。在实际应用中，根据本发明的一个方面提供的正则化方案将改进性能。此外，对频率点进行加权是必要的。下面提供根据本发明的一个方面提供的低复杂度方法。

正则化

等式（60）是简单的最小二乘回归并且必须很注意数据的过度拟合问题。如果麦克风数目小或者如果J个信号的假设在一些频率点不适用，该问题变得严重。

在岭回归中，对于估计变量引入惩罚项，并且在贝叶斯分析中这些惩罚项转换为用于S_n的先验分布。

类似于（60），S_n的MAP估计是

${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n---\left(63\right)$

其中，λ_n是用于频率分量ω_n的正则化参数。现在可以消除Sn并且仅以集中对数似然性工作，该集中对数似然性可以写作

$L_n^{\mathrm{reg}}\left(\mathrm{\&theta;}\right)={\left|\right|(I-P_n\left(\mathrm{\&theta;}\right))X_n\left|\right|}^2---\left(64\right)$

其中，

${\hat{P}}_n\left(\mathrm{\&theta;}\right)={A_n(A_n^H{A}_n+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H---\left(65\right)$

加权

等式（62）中的噪声方差σ²无法从单个瞬象来估计。而是，噪声方差被重判读为加权系数λ_n＝σ^-2，这是被例如为DUET的BSS算法所采用的观点。因此，等式（62）可以被重写为其反映了正则化和加权两者。在实践中，信号带宽可能未被精确已知并且在一些频率点中J个信号的假设不适用。过度拟合问题在这些频率点中变得严重，并且将这些频率点包括在估计过程中会扭曲结果。提供下面的权重来应对不精确的建模、高噪声频率点和离群值频率点：

${\mathrm{\&lambda;}}_n=\frac{{\mathrm{\&tau;}}_n}{{\mathrm{\&Sigma;}}_{n=1}^N{\mathrm{\&tau;}}_n}---\left(67\right)$

其中，

是非负的、非递减的加权函数。其自变量度量所接收信号的、可以在给定DOA向量θ的情况下来解释的部分。λ_n是归一化的权重。

粒子滤波器

基于加权和正则化方案，集中似然性函数读作：

p(X_1∶N|θ)∝e^-γL(θ)   (68)

其中，引入缩放参数γ，其确定似然性函数的峰的锐度。下面提供用于γ的试探。然而，其仅当在频率n处的真实噪声方差是时才是真正的似然性函数。在下文中假设是该情况。现在，将时间维度包括到估计过程中。

首先，对于DOA限定马可夫转移核，以将瞬象之间的信息相关：

$p\left({\mathrm{\&theta;}}_j^k\right|{\mathrm{\&theta;}}_j^{k-1})=\mathrm{\&alpha;}{U}_{[-\frac{\mathrm{\&pi;}}{2},\frac{\mathrm{\&pi;}}{2}]}+(1-\mathrm{\&alpha;})P_{[-\frac{\mathrm{\&pi;}}{2},\frac{\mathrm{\&pi;}}{2}]}N({\mathrm{\&theta;}}_j^{k-1},{\mathrm{\&sigma;}}_{\mathrm{\&theta;}}^2)---\left(69\right)$

其中，

表示在

上的均匀分布的pdf，并且

表示具有平均值

和方差

的正态分布的pdf。使用根据“Yongtao Guan,Roland Fleissner,Paul Joyce和Stephen M.Krone所著的Markov chain monte carlo in small worlds.Statistics and Computing,16:193-202,2006年6月”的小世界建议密度，因为这可能会加速收敛，尤其在该具有多模态似然性函数的情况中。

用Ik表示直至瞬象k的所有测量值（信息）。假设旧的后验pdf的离散逼近是可用的：

$p\left({\mathrm{\&theta;}}^{k-1}\right|I^{k-1})=\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{w}_i^{k-1}{\mathrm{\&delta;}}_{{\mathrm{\&theta;}}_i^{k-1}}---\left(70\right)$

其中，

是在的Dirac质量。

称作带有相关权重

的粒子。根据Bayes规则：

$p\left({\mathrm{\&theta;}}^k\right|I^k)\&Proportional;p\left(X_{1:N}^k\right|{\mathrm{\&theta;}}^k)p\left({\mathrm{\&theta;}}^k\right|{\mathrm{\&theta;}}^{k-1})---\left(71\right)$

对于新的后验的逼近可以如在“Sanjeev Arulampalam,Simon Maskell,Neil Gordon和Tim Clapp所著的A tutorial on particle filters for on-linenon-linear/non-gaussian Bayesian tracking.IEEE Transactions on SignalProcessing,50：174-188,2001”中描述那样通过对来自转移核的各个粒子进行重采样来获得：

以似然性来更新权重，并且重归一化

$w_i^k=\frac{w_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}{{\mathrm{\&Sigma;}}_{i=1}^P{w}_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}---\left(73\right)$

参数γ影响粒子滤波器的反应性。小的值将小的置信度给予新的测量，而大的值快速引起粒子消耗，即，所有权重通过很少的粒子来累积。通过实验发现的是，用于γ的、减小重采样粒子的必要性而维持方法的适配速度的良好试探是：

$\mathrm{\&gamma;}=\frac{10}{{\mathrm{\&Sigma;}}_{i=1}^{P}L\left({\mathrm{\&theta;}}_i\right)}---\left(74\right)$

θ的MAP估计可以通过使用基于柱状图的方法来从粒子获得。然而，粒子未能逃脱如在“H. Sawada, R. Mukai, S. Araki和S. Makino所著的Arobust and precise method for solving the permutation problem offrequency-domain blind source separation. Speech and Audio Processing, IEEETransactions on, 12(5):530-538, 2004”中描述的置换不变性问题。如果对于一些粒子来说将θ_i，j和θ_i，j′互换，似然性函数并不改变其值。将简单的类似EM的聚类技术用于应对该问题。

复杂度

在此根据本发明的一个方面提供的方法的MUST实施方案的主要负荷是在等式（63）中计算

其必须对于P个粒子和N个频率点完成。实正定J×J矩阵的逆需要O(J²)运算并且可以利用BLAS例程有效地实施。从而，更新θ的MAP估计的复杂度是O(NPJ²)。请注意源的数目J还确定对于良好逼近所需的粒子数目P。

在图23-25中示出计算机模拟的结果。

上面提供了用于以MAP估计解决宽带DOA估计问题的粒子技术。与广泛使用的子带技术相反，MAP方法对于关联是鲁棒性的并且借助很少的样本工作，这对于追踪运动的源是重要的。在此提供降低了复杂度的加权和正则化方案，以处理模型不精确性和麦克风的小数目，这两个问题经常在实践中碰到。计算能力方面近来的进步使得其能够在如通过（如图24中示出的）场景2的模拟示出那样的真实世界中使用该方法。在此提供的相对于之前的方法的进一步改进是：精调加权和正则化试探法，以及通过包括对源数目的估计。最后，如果不多于M个信号在给定的频率上活动，则在此根据本发明的一个方面提供的方法提供至J>M个信号的延伸。

下面提供对根据本发明的一个方面的过程的总结。（步骤1）以α=1对等式69进行初始化。（步骤2）根据等式60提供对

的计算。（步骤3）根据等式61（和60）提供对L_n(θ)的计算。（步骤4）根据等式66和67进行权重的计算，该权重可以是用于应对噪声σ_n的正则化权重τ_n。（步骤5）根据等式62进行似然性L(θ)的计算。（步骤6）根据等式74进行参数γ的计算。（步骤7）然后根据等式68进行集中似然性的计算。（步骤8）在两个步骤中完成更新前验分布。在给定等式70和带有α<1的转移核等式69的情况下，计算等式71或者 $p\left({\mathrm{\&theta;}}^k\right|I^k)\&Proportional;p\left(X_{1:N}^k\right|{\mathrm{\&theta;}}^k)p\left({\mathrm{\&theta;}}^k\right|{\mathrm{\&theta;}}^{k-1})p\left({\mathrm{\&theta;}}^{k-1}\right|I^{k-1}).$ 利用等式73提供权重更新。

对于每个粒子i更新向量θ。

下面的参考文献提供一般化地与本发明相关的背景信息并且通过引用结合在本申请中：【1】S.Arulampalam,S.Maskell,N.Gordon和T.Clapp所著的A tutorial on particle filters for on-line non-linear/non-gaussian bayesiantracking.IEEE Transactions on Signal Processing,50:174-188,2001。【2】R.G.Gallager所著的Circularly-symmetric gaussian random vectors.Class Notes,2008年1月。【3】M.Jansson,A.L.Swindlehurst,和B.Ottersten所著的Weightedsubspace fitting for genral array error models.IEEE Transactions on SignalProcessing,46:2484-2498,1997。【4】T.Melia和S.Richard所著的Underdetermined blind source separation in echoic enviroments using desprit,2007。【5】S.Rickard,R.Balan和J.Rosca所著的Real-time time-frequencybased blind source separation.In Proc.Of International Conference onIndependent Component Analysis and Signal Separation(ICA2001),第651-656页,2001。【6】C.P.Robert所著的The Bayesian Choice.Springer,第二版,2001。【7】H.Krim和M.Viberg.所著的Two decades of array signal processingresearch:the parametric approach.Signal Processing Magazine,IEEE,13(4):67-94,1996年7月。【8】R.Roy和T.Kailath所著的Esprit-estimation of signalparameters via rotational invariance techniques.Acoustics,Speech and SignalProcessing,IEEE Transactions on,37(7):984,1989年7月。【9】R.Schmidt所著的Multiple emitter location and signal parameter estimation.Antennas andPropagation,IEEE Transactions on,34(3):276,1986年3月。【10】M.Viberg和B.Ottersten所著的Sensor array processing based on subspace fitting.SignalProcessing,IEEE Transactions on,39(5):1110-1121,1991年5月。【11】P.Stoica和K.C.Sharman所著的Maximum Likelihood Methods for Direction-of-ArrivalEstimation.Acoustics,Speech and Signal Processing,IEEE Transactions on,38(7):1132,1990年7月。【12】J.Lasenby和W.J.Fitzgerald所著的A bayesianapproach to high-resolution beamforming.Radar and Signal Processing,IEEProceedings F,138(6):539-544,1991年12月。【13】James A.Cadzow所著的Multiple source localization–the signal subspace approach.IEEE Transations onAcoustics,Speech,and Signal Processing,38(7):1110-1125,1990年7月。【14】C.Andrieu和A.Doucet所著的Joint Bayesian model selection and estimation ofnoisy sinusoids via reversible jump mcmc.Signal Processing,IEEE Transationson,47(10):2667-2676,1999年10月。【15】J.Huang,P.Xu,Y.Lu和Y.Sun所著的A novel bayesian high-resolution direction-of-arrival estimator.OCEANS,2001.MTS/IEEE Conference and Exhibition,3:1697-1702,2001。【16】H.Wang和M.Kaveh所著的Coherent signal-subspace processing for thedetection and estimation of angles of arrival of multiple wide-band sources.Acoustics,Speech and Signal Processing,IEEE Transactions on,33(4):823,1985年8月。【17】D.N.Swingler和J.Krolik所著的Source location bias in thecoherently focused high-resolution broad-band beamformer.Acoustics,Speechand Signal Processing,IEEE Transactions on,37(1):143-145,1989年1月。【18】Ta-Sung Lee所著的Efficient Wideband source localization using beamforminginvariance technique.Signal Processing,IEEE Transactions on,42(6):1376-1387,1994年6月。【19】Yeo-Sun Yoon,L.M.Kaplan和J.H.McClellan所著的Tops:new doa estimator for wideband signals.Signal Processing,IEEE Transactions on,54(6):1977,2006年6月。【20】W.Ng,J.P.Reilly和T.Kirubarajan所著的Abayesian approach to tracking wideband targets using sensor arrays and particlefilters.Statistical Signal Processing,2003IEEE Workshop on,第510-513页,2003。【21】E.D.di Claudio和R.Parisi所著的Waves:weighted average ofsignal subspaces for robust wideband direction finding.Signal Processing,IEEETransactions on,49(10):2179,2001年10月。【22】C.E.Chen,F.Lorenzelli,R.E.Hudson和K.Yao所著的Maximum likelihood doa estimation of multiplewideband sources in the presence of nonuniform sensor noise.EURASIP Journalon Advances in Signal Processing,2008;Article ID835079,第12页,2008.doi:10.1155/2008/835079,2008。【23】Huigang Wang,S.Kay和S.Saha所著的An importance sampling maximum likelihood direction of arrival estimator.Signal Processing,IEEE Transactions on,56(10):5082-5092,2008。【24】Yongtao Guan,Roland Fleissner,Paul Joyce和Stephen M.Krone所著的Markovchain monte carlo in small worlds.Statistics and Computing,16:193-202,2006年6月。【25】H.Hung和M.Kaveh所著的Focusing matrices for coherentsignal-subspace processing.Acoustics,Speech and Signal Processing,IEEETransactions on,36(8):1272-1281,1988年8月。【26】H.Sawada,R.Mukai,S.Araki和S.Makino所著的A robust and precise method for solving thepermutation problem of frequency-domain blind source separation.Speech andAudio Processing,IEEE Transactions on,12(5):530-538,2004。

下面提供对本发明的不同方面的进一步解释。

到达方向估计是一个被良好地研究的主题，并且表示用于较高级别数据判读的重要构建模块。在本论文中提出的贝叶斯算法（MUST）可以估计和追踪可能相关的多个宽带源的方向。MUST以粒子滤波器逼近源方向在时频域中的后验概率密度函数。与之前的其他算法不同，不需要时间平均，因此可以追踪运动的源。MUST将新的低复杂度加权和正则化方案用于融合来自不同频率的信息和用于克服当很少传感器可用时的过度拟合。

数十年的研究形成了解决到达方向（DOA）估计问题的许多算法，并且这些算法在如在“H.Krim和M.Viberg.所著的Two decades of array signalprocessing research:the parametric approach.Signal Processing Magazine,IEEE,13(4):67-94,1996年7月”中描述那样在如雷达、无线通信或语音识别领域中应用。

DOA估计需要传感器阵列并且探索传感器之间到达的时间差。窄带算法以相位差逼近这些差。现有的用于该问题的大多数方法是如在“R.Roy和T.Kailath所著的Esprit-estimation of signal parameters via rotational invariancetechniques.Acoustics,Speech and Signal Processing,IEEE Transactions on,37(7):984,1989年7月”中描述的ESPRIT或者在“R.Schmidt所著的MultipleEmitter Location and Signal Parameter Estimation.Antennas and Propagation,IEEE Transactions on,34(3):276,1986年3月”中描述的MUSIC的变型，其利用如在“M.Viberg和B.Ottersten所著的Sensor Array Processing Based onSubspace Fitting.Signal Processing,IEEE Transactions on,39(5):1110-1121,1991年5月”中描述的子空间拟合技术并且对于计算解答来说是迅速的。

通常，基于子空间的算法的性能随着信号关联而劣化。诸如在“P.Stoica和K.C.Sharman所著的Maximum Likelihood Methods for Direction-of-ArrivalEstimation.Acoustics,Speech and Signal Processing,IEEE Transactions on,38(7):1132,1990年7月”中描述的最大似然性（ML）和在“J.Lasenby和W.J.Fitzgerald所著的A Bayesian approach to high-resolution beamforming.Radar and Signal Processing,IEE Proceedings F,138(6):539-544,1991年12月”中描述的贝叶斯方法的统计学最优方法长久以来被认为是难以驾驭的，如在“J.A.Cadzow所著的Multiple Source Localization–The Signal SubspaceApproach.IEEE Transations on Acoustics,Speech,and Signal Processing,38(7):1110-1125,1990年7月”中描述那样，然而这些方法近来受到更多的关注，例如在“C.Andrieu和A.Doucet所著的Joint Bayesian Model Selectionand Estimation of Noisy Sinusoids via Reversible Jump MCMC.SignalProcessing,IEEE Transations on,47(10):2667-2676,1999年10月”和“J.Huang,P.Xu,Y.Lu和Y.Sun所著的A Novel Bayesian High-ResolutionDirection-of-Arrival Estimator.OCEANS,2001.MTS/IEEE Conference andExhibition,3:1697-1702,2001”中。

用于宽带DOA的算法通常在时频（t-f）域中表达。于是，窄带假设对于每个子带或频率点是有效的。不相干信号子空间方法（ISSM）计算在所有子带中同时满足信号和噪声子空间正交性条件的DOA估计。另一方面，如在H.Wang和M.Kaveh所著的Coherent signal-subspace processing for thedetection and estimation of angles of arrival of multiple wide-band sources.Acoustics,Speech and Signal Processing,IEEE Transactions on,33(4):823,1985年8月”中描述的相干信号子空间方法（CSSM）根据所有数据计算泛空间协方差矩阵（SCM）。任意窄带信号子空间方法于是可以用于分析泛SCM。然而，良好的初始估计对于如在“D.N.Swingler和J.Krolik所著的Sourcelocation bias in the coherently focused high-resolution broad-band beamformer.Acoustics,Speech and Signal Processing,IEEE Transactions on,37(1):143-145,1989年1月”中描述那样正确地使子带SCM与泛SCM一致是必要的。例如在“T.-S.Lee所著的Efficient Wideband Source Localization UsingBeamforming Invariance Technique.Signal Processing,IEEE Transactions on,42(6):1376-1387,1994年6月”中描述的BI-CSSM或者在“Y.-S.Yoon,L.M.Kaplan和J.H.McClellan所著的TOPS:New DOA Estimator for WidebandSignals.Signal Processing,IEEE Transactions on,54(6):1977,2006年6月”中描述的TOPS的方法被开发以减轻该问题。

子空间方法将信号和噪声子空间的正交性用作优化的标准。然而，在数学上更吸引人的方法是使估计基于决定理论架构上。前提是计算DOA的后验概率密度函数（pdf），其可以通过采用粒子滤波器来实现。这种方法在“W.Ng,J.P.Reilly和T.Kirubarajan所著的A bayesian approach to trackingwideband targets using sensor arrays and particle filters.Statistical SignalProcessing,2003IEEE Workshop on,第510-513页,2003”中采用，其中在时域中表达贝叶斯极大后验似然性（MAP）估计器。

贝叶斯MAP估计器利用信号的时频表示来提供。时频分析的优点通过在诸如在“S.Rickard,R.Balan和J.Rosca所著的Real-time time-frequencybased blind source separation.In Proc.Of International Conference onIndependent Component Analysis and Signal Separation(ICA2001),第651-656页,2001”中描述的DUET和在“T.Melia和S.Rickard所著的Underdeterminedblind source separation in echoic environments using DESPRIT.EURASIPJournal on Advances in Signal Processing,2007:Article ID86484,第19页,doi:10.1155/2007/86484,2007”中描述的DESPRIT的盲源分离（BSS）领域中使用的技术来示出。这些方法探索不相似的信号特点以将信号分离，并且对于语音信号良好适用。

所示出的多源追踪（MUST）将新颖的、试探式的加权方案用于对于各频率将信息组合。粒子滤波器逼近DOA的后验密度，然后MAP估计被提取。在本发明的上下文中还将示出一些广泛使用的算法。将提供对MUST的详细描述。示出MUST的模拟结果并且将其与CSSM、IMUSIC和如在“E.D.di Claudio和R.Parisi所著的Waves:weighted average of signal subspaces forrobust wideband direction finding.Signal Processing,IEEE Transactions on,49(10):2179,2001年10月”中描述WAVES方法相比较。

问题表达和相关工作

考虑带有在传感器1和m之间的、记做d_m的距离的M个传感器的线性阵列。撞击到该阵列上的是来自不同方向θ_j的J个未知波前。波前的传播速度是c。源的数目J被假设为已知并且J≤M。通过用于回声路径的附加源来应对回声环境。麦克风被假设为在源的远场中。在DFT域中，在第n个子带中第m个传感器处接收的信号可以被建模为：

$X_m\left({\mathrm{\&omega;}}_n\right)=\underset{j=0}{\overset{J-1}{\mathrm{\&Sigma;}}}{S}_j\left({\mathrm{\&omega;}}_n\right)e^{-i{\mathrm{\&omega;}}_n{y}_m\sin \left({\mathrm{\&theta;}}_j\right)}+N_m\left({\mathrm{\&omega;}}_n\right)---\left(75\right)$

其中，S_j(ω_n)是第j个源信号，N_m(ω_n)是噪声，并且v_m＝d_m/c。噪声被假设为在每个频率内是循环对称复高斯分布（CSCG）并且独立地和相同地分布（iid）的，即

噪声方差会随着ω_n变化。

如果定义：

X_n＝[X₁(ω_n)...X_M(ω_n)]^T   (76)

N_n＝[N₁(ω_n)...N_M(ω_n)]^T   (77)

S_n＝[S₁(ω_n)...S_J(ω_n)]^T   (78)

θ＝[θ₁，...，θ_j]^T   (79)

则等式（75）可以以矩阵向量符号的形式再写作：

X_n＝A_n(θ)S_n+N_n    (80)

其中M×J导引矩阵：

A_n(θ)＝[a(ω_n，θ₁)，...，a(ω_n，θ_J)]   (81)

其列是M×1阵列流形：

$a({\mathrm{\&omega;}}_n,{\mathrm{\&theta;}}_j)={[1e^{-i{\mathrm{\&omega;}}_n{v}_1\sin \left({\mathrm{\&theta;}}_j\right)}...e^{-i{\mathrm{\&omega;}}_n{v}_{M-1}\sin \left({\mathrm{\&theta;}}_j\right)}]}^T---\left(82\right)$

子空间方法

用于解决DOA问题的最常用的算法根据所接收数据的样本协方差矩阵计算信号和噪声子空间，并且选择这样的θ_j，其对应的阵列流形a(θ_j)最接近信号子空间，即，

${\widehat{\mathrm{\&theta;}}}_j=\underset{\mathrm{\&theta;}}{\arg \min }a{\left(\mathrm{\&theta;}\right)}^H{E}_N{E}_N^{H}a\left(\mathrm{\&theta;}\right)---\left(83\right)$

其中，E_N的列形成噪声子空间的正交基。不相干方法对于每个子带计算信号和噪声子空间E_N(ω_n)，并且θ_j被选择为平均起来满足（83）。相干方法通过将所有数据转换到参考频率ω₀来计算参考信号和噪声子空间。于是仅对于参考阵列流形a(ω₀，θ)校验正交条件（83）。这些方法（其中CSSM和WAVES是两个代表）示出了比不相干方法显著更好的性能，尤其是对于高度相关和低SNR信号而言。但是，到参考频率的转换需要良好的初始DOA估计，并且如何可以获得这些估计并不显然。

最大似然性方法

与子空间算法不同，ML方法计算来自An矩阵的信号子空间并且选择在将其在子空间上的投影最大化方面最优地匹配于所观察数据的

其可以被示出为与使似然性最大化等效：

$\widehat{\mathrm{\&theta;}}=\underset{\mathrm{\&theta;}}{\arg \max }\left|\right|P_n\left(\mathrm{\&theta;}\right)X_n\left|\right|---\left(84\right)$

其中，

是在由A_n(θ)的列张开的信号子空间上的投影算子，其中，这些决定性的ML估计假设对这些信号并无认识。而如果信号统计已知，则可以如在“P.Stoica和A.Nehorai所著的On theConcentrated Stochastic Likelihood Function in Array Signal Processing.Circuits,Systems,and Signal Processing,14:669-674,1995.10.1007/BF01213963”描述那样计算随机ML估计。

如果噪声方差对于所有频率相等，则可以通过针对各频率地进行求和（84）来获得用于宽带问题的全面对数似然性函数。至今未处理噪声方差变化的问题。

在“C.E.Chen,F.Lorenzelli,R.E.Hudson和K.Yao所著的Maximumlikelihood doa estimation of multiple wideband sources in the presence ofnonuniform sensor noise.EURASIP Journal on Advances in Signal Processing,2008;Article ID835079,第12页,2008.doi:10.1155/2008/835079,2008”中，研究了如下情况：对于传感器不一致的噪声却是对于各频率恒定的。

与子空间方法相比，ML方法提供了关于阵列形式和信号相关性的更高灵活性，并且通常对于小样本尺寸示出更好性能，然而（84）中的非线性多维优化在计算上是复杂的。近来，对于窄带情况提出了重要性采样方法，以如在“Hu.Wang,S.Kay和S.Saha所著的An importance sampling maximumlikelihood direction of arrival estimator.Signal Processing,IEEE Transactions on,56(10):5082-5092,2008”中描述那样有效地解决优化问题。在MUST实施方案中采用的粒子滤波器沿着这些线进行优化。

多源追踪（MUST）

在等式（75）的模型下，如果以Sn和θ为条件，则观察X₁(ω_n)，...，X_M(ω_n)是iid CSCG随机变量。因此，联合pdf分解为边缘pdf。对于每个频率ωn，通过下式给出负的对数似然性：

-logp(X_n[S_n，θ)∝||X_n-A_n(θ)S_n)||²   (85)

通常将对于Sn的ML解计算为：

表示A_n的Moore-Penrose逆。θ的ML解通过将剩余集中负对数似然性最小化来得到。

如果噪声方差σ_n ²已知，则可以通过将用于所有频率的似然性求和来计算全局（负的）集中对数似然性：

$L\left(\mathrm{\&theta;}\right)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}\frac{L_n\left(\mathrm{\&theta;}\right)}{{\mathrm{\&sigma;}}_n^2}---\left(88\right)$

该条件函数前面已经陈述过，并且如在“J.A.Cadzow所著的Multiplesource localization–the signal subspace approach.IEEE Transactions onAcoustics,Speech,and Signal Processing,38(7):1110-1125,1990年7月”中描述那样被看做难以处理的（1990年）。与子空间方法不同，ML方法和利用对于源信号的ML估计的MUST对于相关源是不敏感的，因为其并不尝试估计级别J信号子空间。

更下面根据本发明的一个方面提供了一种粒子滤波器方法以解决对于多个瞬象的滤波问题，其附带地自然解决了优化问题。在实际应用中，如下面将示出那样，正则化方案将改进性能。此外，对频率点进行加权是必要的。下面提供根据本发明的一个方面提供的低复杂度方法。

正则化

等式（86）是简单的最小二乘回归并且必须很注意数据的过度拟合问题。如果麦克风数目小或者如果J个信号的假设在一些频率点不适用，该问题变得严重。

在岭回归中，对于估计变量引入惩罚项，并且在贝叶斯分析中这些惩罚项转换为用于S_n的先验分布。为了降低复杂度，以用于所有频率和源的单个全局正则化参数λ来使用CSCG先验：

$\left(S_n\right)\&Proportional;\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}\mathrm{\&lambda;}{\left|S_j\left({\mathrm{\&omega;}}_n\right)\right|}^2---\left(89\right)$

类似于（86），S_n的MAP估计是

${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n---\left(90\right)$

现在可以消除S_n并且仅以集中对数似然性工作，该集中对数似然性可以写作

$L_n^{\mathrm{reg}}\left(\mathrm{\&theta;}\right):={\left|\right|(I-{\hat{P}}_n\left(\mathrm{\&theta;}\right))X_n\left|\right|}^2---\left(91\right)$

其中，

${\hat{P}}_n\left(\mathrm{\&theta;}\right)={A_n(A_n^H{A}_n+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H---\left(92\right)$

参数λ是特别选择的。所发现的是，相对于在很少麦克风可用情况下的直至10^-3M个源，在很多麦克风可用情况下的值10^-5M会改进估计。如果关于S_n的信息是可用的，则可以设计更为复杂的正则化模型。

加权

等式（88）中的噪声方差

无法从单个瞬象来估计。而是，噪声方差被重判读为加权系数

这是被例如为DUET的BSS算法所采用的观点。在实践中，信号带宽可能未被精确已知并且在一些频率点中J个信号的假设不适用。过度拟合问题在这些频率点中变得严重，并且将这些频率点包括在估计过程中会扭曲结果。根据本发明的一个方面提供下面的权重来应对不精确的建模、高噪声频率点和离群值频率点：

${\mathrm{\&tau;}}_n=\frac{{\widehat{\mathrm{\&tau;}}}_n}{{\mathrm{\&Sigma;}}_{n=1}^N{\widehat{\mathrm{\&tau;}}}_n}---\left(94\right)$

其中，

是非负的、非递减的加权函数。其自变量度量所接收信号的、可以在给定DOA向量θ的情况下来解释的部分。τ_n是归一化的权重。

粒子滤波器

基于加权和正则化方案，集中似然性函数读作：

p(X_1∶N|θ)∝e^-γL(θ)   (95)

其中，引入缩放参数，其确定似然性函数的峰的锐度。下面提供用于γ的试探。然而，其仅当在频率n处的真实噪声方差是时才是真正的似然性函数。在下文中假设是该情况。现在，将时间维度包括到估计过程中。

首先，对于DOA限定马可夫转移核，以将瞬象k和k-1之间的信息相关：

$p\left({\mathrm{\&theta;}}_j^k\right|{\mathrm{\&theta;}}_j^{k-1})=\mathrm{\&alpha;}{U}_{[-\frac{\mathrm{\&pi;}}{2},\frac{\mathrm{\&pi;}}{2}]}+(1-\mathrm{\&alpha;})N({\mathrm{\&theta;}}_j^{k-1},{\mathrm{\&sigma;}}_{\mathrm{\&theta;}}^2)---\left(96\right)$

其中，表示在

上的均匀分布的pdf，并且表示具有平均值

和方差

的正态分布的pdf。使用在“Y.Guan,R.Fleiβner,P.Joyce和S.M.Krone所著的Markov chain monte carlo in small worlds.Statistics and Computing,16:193-202,2006年6月”中描述的小世界建议密度，因为这可能会加速收敛，尤其在该具有多模态似然性函数的情况中。“Y.Guan,R.Fleiβner,P.Joyce和S.M.Krone所著的Markov chain monte carlo insmall worlds.Statistics and Computing,16:193-202,2006年6月”的作者给出了用于选择α的精确规则，其需要对于后验pdf的精确认识。然而，其还争辩说α∈[10^-4，10^-1]是良好的经验方法。

用I^k表示直至瞬象k的所有测量值（信息）。假设对于I^k-1的特定实现方案而言旧的后验pdf的离散逼近是可用的：

$p\left({\mathrm{\&theta;}}^{k-1}\right|I^{k-1})=\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{w}_i^{k-1}{\mathrm{\&delta;}}_{{\mathrm{\&theta;}}_i^{k-1}}---\left(97\right)$

其中，

是在

的Dirac质量。

和其相关权重都称作粒子。这些粒子包含直至瞬象k-1的所有可用信息。θ的索引i指的是P个粒子之一，并且θ_i＝[θ₁，...，。θ_Ｊ]_i＝θ_i，1，...，θ_i，J]。新的测量值

通过Bayes规则迭代地积分：

$p\left({\mathrm{\&theta;}}^k\right|I^l)\&Proportional;p\left(X_{1:N}^k\right|{\mathrm{\&theta;}}^k)p\left({\mathrm{\&theta;}}^k\right|{\mathrm{\&theta;}}^{k-1})p\left({\mathrm{\&theta;}}^{k-1|}\right|I^{k-1})---\left(98\right)$

对于新的后验的逼近可以如在“S.Arulampalam,S.Maskell,N.Gordon和T.Clapp所著的A tutorial on particle filters for on-linenon-linear/non-gaussian Bayesian tracking.IEEE Transactions on SignalProcessing,50-174-188,2001”中描述那样以两个步骤获得。首先，从转移核对每个粒子进行重采样：

${\mathrm{\&theta;}}_i^k~p\left({\mathrm{\&theta;}}_i^k\right|{\mathrm{\&theta;}}_i^{k-1})---\left(99\right)$

在第二步骤中，以似然性来更新权重，并且重归一化

${\widehat{\mathrm{\&omega;}}}_i^k={\mathrm{\&omega;}}_i^{k-1}p\left(X_{1:N}^k\right|{\mathrm{\&theta;}}_i^j)---\left(100\right)$

${\mathrm{\&omega;}}_i^k=\frac{{\widehat{\mathrm{\&omega;}}}_i^k}{{\mathrm{\&Sigma;}}_{i=1}^P{\widehat{\mathrm{\&omega;}}}_i^k}---\left(101\right)$

参数γ影响粒子滤波器的反应性。小的值将小的置信度给予新的测量，而大的值快速引起粒子消耗，即，所有权重通过很少的粒子来累积。通过实验发现的是，用于γ的、减小重采样粒子的必要性而维持算法的适配速度的良好试探是：

$\mathrm{\&gamma;}=\frac{10}{{\mathrm{\&Sigma;}}_{i=1}^{P}L\left({\mathrm{\&theta;}}_i\right)}---\left(102\right)$

如果粒子的有效数目

$N_{\mathrm{eff}}={\left(\underset{i=1}{\overset{P}{\mathrm{\&Sigma;}}}{\left({\mathrm{\&omega;}}_i^k\right)}^2\right)}^{-1}---\left(103\right)$

落到预定阈值之下，则粒子消耗的问题通过重采样来处理。该粒子滤波器已知为如在“S.Arulampalam,S.Maskell,N.Gordon和T.Clapp所著的Atutorial on particle filters for on-line non-linear/non-gaussian bayesian tracking.IEEE Transactions on Signal Processing, 50:174-188, 2001”中描述的采样重要性重采样（SIR）滤波器。

θ的MAP估计可以通过使用基于柱状图的方法来从粒子获得。然而，粒子未能逃脱如在“H. Sawada, R. Mukai, S. Araki和S. Makino所著的Arobust and precise method for solving the permutation problem offrequency-domain blind source separation. Speech and Audio Processing, IEEETransactions on, 12(5):530-538, 2004”中描述的置换不变性问题。如果对于一些粒子来说将θ_i，j′和θ_i，j″互换，似然性函数并不改变其值。为了应对该问题，使用一些简单的聚类技术，其将θ_i，j′与在之前的时间步骤中从所有粒子计算的的最接近的估计相关联。如果将数个θ_i，j′、θ_i，j″的严分配给相同的源，则该问题通过重分配，如果可能，或者通过忽视MAP估计的计算中的θ_i，j′和θ_i，j″之一来解决。

复杂度

MUST的主要负荷是在等式（90）中计算

其必须对于P个粒子和N个频率点完成。求解J个线性方程组进行需要O(J³)运算并且可以利用BLAS例程有效地实施。从而，更新θ的MAP估计的复杂度是O(NPJ³)。请注意源的数目J还确定对于良好逼近所需的粒子数目P。

计算机模拟

执行三个不同的计算机模拟的场景用于对比。在所有场景中，由M个传感器记录带有关联ρ∈[-1，1]的等功率高斯噪声源。在传感器通带f₀±Δf内的N个频率点上进行处理。WAVES、CSSM和IMUSIC基于当前的和Q个过去的瞬象计算DOA估计。这允许在线动态计算。粒子以均匀分布来初始化。所使用的加权函数是

在前两个场景中，在所有元件之间，传感器间间隔是

其中

参数值在表3中总结。

表3

所有结果都是基于对于参数的每个组合的100次Monte Carlo运行而得出的。

WAVES和CSSM将如在“H.Hung和M.Kaveh所著的Focusing matricesfor coherent signal-subspace processing.Acoustics,Speech and SignalProcessing,IEEE Transactions on,36(8):1272-1281,1988年8月”描述的RSS聚焦矩阵用于使样本SCM与作为聚焦角度的真实角度符合。这是不现实的假设但是为相干方法提供了性能上边界。WAVES算法如在“E.D.di Claudio和R.Parisi所著的Waves:weighted average of signal subspaces for robustwideband direction finding.Signal Processing,IEEE Transactions on,49(10):2179,2001年10月”中描述那样实施，并且将Root-MUSIC用于CSSM和WAVES两者。

第一场景在“H.Wang和M.Kaveh所著的Coherent signal-subspaceprocessing for the detection and estimation of angles of arrival of multiplewide-band sources.Acoustics,Speech and Signal Processing,IEEE Transactionson,33(4):823,1985年8月”和“E.D.di Claudio和R.Parisi所著的Waves:weighted average of signal subspaces for robust wideband direction finding.Signal Processing,IEEE Transactions on,49(10):2179,2001年10月”中用于测试宽带DOA，并且其示出在图26中。图26对照用于源关联系数ρ的不同值的SNR示出了其中在2°内检测到所有源的模块的百分数。标记ρ指的是当所有四个MUST曲线几乎重合（collapse）时的WAVES和CSSM曲线。结果示出了粒子滤波器算法可以在低SNR值处和对于任意关联分辨间隔紧密的信号。相反，CSSM的性能随着关联而降低，IMUSIC未成功分辨所有四个源。

对于第二场景，如在图27中示出那样使用与声音信号有关的参数。对照对于ρ=0（直线）和ρ=0.75（虚线）的SNR示出了其中在2.5度内检测到所有源的模块的百分比。参数选择为示出粒子滤波器的缩减版本的性能，该版本利用包含最多能量和相对小数目的粒子的频率点的仅10%。在这些设置下，实现了在双核笔记本计算机上的实时计算。MUST的性能在IMUSIC和CSSM之间。WAVES结果与CSSM结果几乎相同并且出于易读性原因而未被示出。

在该第三场景中，在图28示出了MUST方法的用于追踪运动的源的潜力。使用带有距离d_m-d_m-1＝d+Δd（其中

并且Δd～U_{[-0.2d，0.2d]}）的M=5个传感器的非均匀的、线性的阵列。信号集中在信号通带 $[f_0-\mathrm{\&Delta;}{f}_{\mathrm{SRC}},f_0+\mathrm{\&Delta;}{f}_{\mathrm{SRC}}]\&Subset;[f_0-\mathrm{\&Delta;f},f_0+\mathrm{\&Delta;f}]$ 中（其中Δf_SRC=20Hz），并且使用0dB总信号功率对总噪声功率的SNR。MUST方法成功估计运动的源的正确的源位置，而该场景提出用于静态子空间方法的问题。

尽管在此示出、描述和指出了本发明的如应用于其优选实施例那样的基本新颖特征，而将理解的是，可以由本领域技术人员做出以所示方法和相同的细节以及其操作中的细节为形式的各种省略、替代和改变，而不偏离本发明的精神。因此，旨在仅如通过权利要求的范围所指示的那样而受限制。

权利要求书(按照条约第19条的修改)

1.一种用于对多个源中的源所生成的信号的到达方向进行估计的方法，包括：

处理器将在时频域中操作的粒子滤波器初始化；

由多个接收器从所述多个源接收多个信号，所述多个信号包括所述信号，至少两个所接收的信号相关；

所述处理器为所述信号的多个频率分量中的每个频率分量确定正则化系数和加权系数；

所述处理器基于所述正则化系数和所述加权系数对于所述信号的所述多个频率分量中的每个频率分量确定到达方向的似然性；

所述处理器基于所述信号的所述多个频率分量中的每个频率分量的到达方向的似然性来确定所述信号的到达方向的全局似然性；

所述处理器基于所述全局似然性更新所述粒子滤波器的多个权重；以及

所述处理器通过应用所述粒子滤波器的更新过的权重来生成对所述信号的到达方向的估计。

2.根据权利要求1所述的方法，其中，所述全局似然性L(θ)通过如下表达式来确定：

$L\left(\mathrm{\&theta;}\right)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}{L}_n^{\mathrm{reg}}\left(\mathrm{\&theta;}\right){\mathrm{\&lambda;}}_n,$

其中，θ表示所述多个源中的所有源的所有到达方向的向量，n表示频率分量，N表示频率分量的总数目，

是所有源的所有到达方向的正则化的似然性，其是通过用于频率分量n的归一化的加权系数λ_n调节过的。

3.根据权利要求1所述的方法，其中，用于频率分量n的加权系数是由如下表达式确定的归一化的加权系数λ_n：

${\mathrm{\&lambda;}}_n=\frac{{\mathrm{\&tau;}}_n}{{\mathrm{\&Sigma;}}_{n=1}^N{\mathrm{\&tau;}}_n},$ 其中其中

τ_n表示未归一化的加权系数，

表示用于信号子空间上的频率分量n的正则化的射影算子，θ表示所述多个源中所有源的所有到达方向的向量，X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量，以及

表示非负非递减的加权函数。

4.根据权利要求1所述的方法，还包括：所述处理器确定所述信号的频率分量的极大后验似然性（MAP）估计。

5.根据权利要求4所述的方法，其中，由表示的MAP估计由下式确定：

${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n$

其中，θ表示在所述多个源中所有源的所有到达方向的向量，A表示导引矩阵，

表示用于频率分量n的导引矩阵的厄密共轭，λ_n表示用于频率分量n的正则化系数，I表示单位矩阵，并且X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量。

6.根据权利要求1所述的方法，其中，所述粒子滤波器的更新过的权重

通过如下表达式来确定：

$w_i^k=\frac{w_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}{{\mathrm{\&Sigma;}}_{i=1}^P{w}_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}$

其中，k表示更新时间，i表示粒子索引，p表示概率分布函数，P表示粒子总数，X_1∶N表示在时刻k的频率分量的向量，并且

表示在时刻k对粒子i的到达方向估计。

7.根据权利要求6所述的方法，其中，所述全局似然性以缩放参数来调节。

8.根据权利要求1所述的方法，其中，所述多个源中的至少一个源是运动的源。

9.根据权利要求1所述的方法，其中，信号是声学信号。

10.根据权利要求1所述的方法，其中，所述多个源中的至少一个源与发电机关联。

11.一种用于对多个源中的源所生成的信号的到达方向进行估计的系统，包括：

存储器，其能够存储数据和指令；

处理器，其能够检索和执行来自所述存储器的指令，以实施如下步骤：

将在时频域中操作的粒子滤波器初始化；

处理由多个接收器接收的来自所述多个源的多个信号，所述多个信号包括所述信号，至少两个所接收的信号相关；

对于所述信号的多个频率分量中每个频率分量确定正则化系数和加权系数；

基于所述正则化系数和所述加权系数对于所述信号的所述多个频率分量中每个频率分量确定到达方向的似然性；

基于对于所述信号的所述多个频率分量中每个频率分量确定的到达方向的似然性来确定所述信号的到达方向的全局似然性；

基于所述全局似然性更新所述粒子滤波器的多个权重；以及

通过应用所述粒子滤波器的更新过的权重来生成对所述信号的到达方向的估计。

12.根据权利要求11所述的系统，其中，所述全局似然性L(θ)通过如下表达式来确定：

$L\left(\mathrm{\&theta;}\right)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}{L}_n^{\mathrm{reg}}\left(\mathrm{\&theta;}\right){\mathrm{\&lambda;}}_n,$

其中，θ表示所述多个源中所有源的所有到达方向的向量，n表示频率分量，N表示频率分量的总数，

是对于所有源得出的所有到达方向的正则化的似然性，其是通过频率分量n的归一化的加权系数λ_n调节过的。

13.根据权利要求11所述的系统，其中，用于频率分量n的加权系数是归一化的加权系数λ_n，其通过如下表达式来确定：

${\mathrm{\&lambda;}}_n=\frac{{\mathrm{\&tau;}}_n}{{\mathrm{\&Sigma;}}_{n=1}^N{\mathrm{\&tau;}}_n},$ 其中

其中

τ_n表示未归一化的加权系数，

表示用于在信号子空间上的频率分量n的正则化的射影算子，θ表示所述多个源中所有源的所有到达方向的向量，X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量，并且

表示非负非递减的加权函数。

14.根据权利要求11所述的系统，还包括所述处理器确定所述信号的频率分量的极大后验似然性（MAP）估计。

15.根据权利要求14所述的系统，其中，由

表示的所述MAP估计由如下表达式来确定：

${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n$

其中，θ表示所述多个源中所有源的所有到达方向的向量，A表示导引矩阵，表示用于频率分量n的导引矩阵的厄密共轭，λ_n表示用于频率分量n的正则化系数，I表示单位矩阵，并且X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量。

16.根据权利要求11所述的系统，其中，所述粒子滤波器的更新过的权重

由如下表达式来确定：

$w_i^k=\frac{w_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}{{\mathrm{\&Sigma;}}_{i=1}^P{w}_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}$

其中，k表示更新时间，i表示粒子索引，p表示概率分布函数，P表示粒子总数，X_1∶N表示在时刻k上频率分量的向量，并且

表示在时刻k对粒子i的到达方向估计。

17.根据权利要求16所述的系统，其中，所述全局似然性借助缩放参数来调节。

18.根据权利要求11所述的系统，其中，所述多个源中的至少一个源是运动的源。

19.根据权利要求11所述的系统，其中，信号是声学信号。

20.根据权利要求11所述的系统，其中，所述多个源中的至少一个源与发电机关联。

Claims (20)

1.一种用于对多个源中的源所生成的信号的到达方向进行估计的方法，包括：

处理器借助在时频域中操作的粒子滤波器来估计由所述源生成的之前的信号的第一到达方向；

由多个接收器从所述多个源接收多个信号，所述多个信号包括所述信号，至少两个所接收的信号相关；

所述处理器为所述信号的多个频率分量中的每个频率分量确定正则化系数和加权系数；

所述处理器基于所述正则化系数和所述加权系数对于所述信号的所述多个频率分量中的每个频率分量确定到达方向的似然性；

所述处理器基于所述信号的所述多个频率分量中的每个频率分量的到达方向的似然性来确定所述信号的到达方向的全局似然性；

所述处理器更新所述粒子滤波器的多个权重；以及

所述处理器通过应用所述粒子滤波器的更新过的权重来生成对所述信号的到达方向的估计。

2.根据权利要求1所述的方法，其中，所述全局似然性L(θ)通过如下表达式来确定：

$L\left(\mathrm{\&theta;}\right)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}{L}_n^{\mathrm{reg}}\left(\mathrm{\&theta;}\right){\mathrm{\&lambda;}}_n,$

其中，θ表示所述多个源中的所有源的所有到达方向的向量，n表示频率分量，N表示频率分量的总数目，

是所有源的所有到达方向的正则化的似然性，其是通过用于频率分量n的归一化的加权系数λ_n调节过的。

3.根据权利要求1所述的方法，其中，用于频率分量n的加权系数是由如下表达式确定的归一化的加权系数λ_n：

${\mathrm{\&lambda;}}_n=\frac{{\mathrm{\&tau;}}_n}{{\mathrm{\&Sigma;}}_{n=1}^N{\mathrm{\&tau;}}_n},$ 其中其中

τ_n表示未归一化的加权系数，

表示用于信号子空间上的频率分量n的正则化的射影算子，θ表示所述多个源中所有源的所有到达方向的向量，X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量，以及

表示非负非递减的加权函数。

4.根据权利要求1所述的方法，还包括：所述处理器确定所述信号的频率分量的极大后验似然性（MAP）估计。

5.根据权利要求4所述的方法，其中，由

表示的MAP估计由下式确定：

${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n$

其中，θ表示在所述多个源中所有源的所有到达方向的向量，A表示导引矩阵，

表示用于频率分量n的导引矩阵的厄密共轭，λ_n表示用于频率分量n的正则化系数，I表示单位矩阵，并且X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量。

6.根据权利要求1所述的方法，其中，所述粒子滤波器的更新过的权重通过如下表达式来确定：

$w_i^k=\frac{w_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}{{\mathrm{\&Sigma;}}_{i=1}^P{w}_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}$

其中，k表示更新时间，i表示粒子索引，p表示概率分布函数，P表示粒子总数，X_1∶N表示在时刻k的频率分量的向量，并且

表示在时刻k对粒子i的到达方向估计。

7.根据权利要求6所述的方法，其中，所述全局似然性以缩放参数来调节。

8.根据权利要求1所述的方法，其中，所述多个源中的至少一个源是运动的源。

9.根据权利要求1所述的方法，其中，信号是声学信号。

10.根据权利要求1所述的方法，其中，所述多个源中的至少一个源与发电机关联。

11.一种用于对多个源中的源所生成的信号的到达方向进行估计的系统，包括：

存储器，其能够存储数据和指令；

处理器，其能够检索和执行来自所述存储器的指令，以实施如下步骤：

借助在时频域中操作的粒子滤波器估计由所述源生成的之前的信号的第一到达方向；

处理由多个接收器接收的来自所述多个源的多个信号，所述多个信号包括所述信号，至少两个所接收的信号相关；

对于所述信号的多个频率分量中每个频率分量确定正则化系数和加权系数；

基于所述正则化系数和所述加权系数对于所述信号的所述多个频率分量中每个频率分量确定到达方向的似然性；

基于对于所述信号的所述多个频率分量中每个频率分量确定的到达方向的似然性来确定所述信号的到达方向的全局似然性；

更新所述粒子滤波器的多个权重；以及

通过应用所述粒子滤波器的更新过的权重来生成对所述信号的到达方向的估计。

12.根据权利要求11所述的系统，其中，所述全局似然性L(θ)通过如下表达式来确定：

$L\left(\mathrm{\&theta;}\right)=\underset{n=1}{\overset{N}{\mathrm{\&Sigma;}}}{L}_n^{\mathrm{reg}}\left(\mathrm{\&theta;}\right){\mathrm{\&lambda;}}_n,$

其中，θ表示所述多个源中所有源的所有到达方向的向量，n表示频率分量，N表示频率分量的总数，

是对于所有源得出的所有到达方向的正则化的似然性，其是通过频率分量n的归一化的加权系数λ_n调节过的。

13.根据权利要求11所述的系统，其中，用于频率分量n的加权系数是归一化的加权系数λ_n，其通过如下表达式来确定：

${\mathrm{\&lambda;}}_n=\frac{{\mathrm{\&tau;}}_n}{{\mathrm{\&Sigma;}}_{n=1}^N{\mathrm{\&tau;}}_n},$ 其中

其中

τ_n表示未归一化的加权系数，

表示用于在信号子空间上的频率分量n的正则化的射影算子，θ表示所述多个源中所有源的所有到达方向的向量，X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量，并且

表示非负非递减的加权函数。

14.根据权利要求11所述的系统，还包括所述处理器确定所述信号的频率分量的极大后验似然性（MAP）估计。

15.根据权利要求14所述的系统，其中，由

表示的所述MAP估计由如下表达式来确定：

${\hat{S}}_n\left(\mathrm{\&theta;}\right)={(A_n^{H}A+{\mathrm{\&lambda;}}_{n}I)}^{-1}{A}_n^H{X}_n$

其中，θ表示所述多个源中所有源的所有到达方向的向量，A表示导引矩阵，表示用于频率分量n的导引矩阵的厄密共轭，λ_n表示用于频率分量n的正则化系数，I表示单位矩阵，并且X_n表示由所述多个接收器接收的所述多个信号的频率分量n的向量。

16.根据权利要求11所述的系统，其中，所述粒子滤波器的更新过的权重

由如下表达式来确定：

$w_i^k=\frac{w_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}{{\mathrm{\&Sigma;}}_{i=1}^P{w}_i^{k-1}p\left(X_{1:N}\right|{\mathrm{\&theta;}}_i^k)}$

其中，k表示更新时间，i表示粒子索引，p表示概率分布函数，P表示粒子总数，X_1∶N表示在时刻k上频率分量的向量，并且

表示在时刻k对粒子i的到达方向估计。

17.根据权利要求16所述的系统，其中，所述全局似然性借助缩放参数来调节。

18.根据权利要求11所述的系统，其中，所述多个源中的至少一个源是运动的源。

19.根据权利要求11所述的系统，其中，信号是声学信号。

20.根据权利要求11所述的系统，其中，所述多个源中的至少一个源与发电机关联。

Images