机器人学中的声源定位综述

机器人学中的声源定位：综述

1. 引言

声源定位（SSL）的目标是自动估计声源的位置。在机器人学中， 此功能在多种场景中都非常有用，例如：在服务员类型任务中定位人 类说话者、在无视觉接触的救援场景中，或用于绘制未知声学环境的 地图。其性能对整个机器人听觉系统至关重要，因为其估计结果通常 被用于后续处理阶段，如声源分离、声源分类和自动语音识别。

声源位置有两个分量可以在声源定位中（以极坐标表示）被估计：
• 到达方向估计（可在一维或二维中进行）
• 距离估计。

在真实生活场景中的声源定位需要考虑到环境中可能存在多个活 跃的声源。因此，也有必要估计多个同时发生的声源的位置。此外， 机器人和声源都是移动的，因此通过时间跟踪其位置非常重要。

机器人学界通过改进传统技术，显著推动了声源定位的发展，这 些技术包括：单一声源到达方向（DOA）估计、基于学习的方法（如 神经网络和流形学习）、基于波束形成的方法、子空间方法、时域声 源聚类以及卡尔曼滤波器和粒子滤波等跟踪技术。在将这些技术实现 在机器人平台上的过程中，若干与机器人中声源定位相关的方面变得 尤为突出，包括：所用麦克风的数量和类型、声源的数量和移动性、 对噪声和混响的鲁棒性、所采用的阵列几何结构类型、所构建的机器 人平台类型等。

正如本综述所示，机器人学中的声源定位领域已相当成熟，这一点可从 该主题近期出现的多篇综述中得到证明。例如，[1,2]对双耳机器人听觉进 行了综述，[3]提供了中文版的声源定位综合综述，[4]介绍了一些声源定位 基于双耳技术和多麦克风阵列的研究，[5]全面概述了机器人听觉领 域。本研究旨在综述在各类机器人（如服务型、救援、群体、工业等） 上实现的声源定位文献。我们还回顾了针对机器人平台实施的研究工 作，即使这些工作尚未实际在机器人平台上实现。此外，我们综述了 用于声源定位训练或评估的资源，其中一些资源虽非从机器人角度采 集，但可应用于机器人任务。最后，我们纳入了仅使用单个麦克风进 行声源定位的研究，尽管这些研究尚未应用于机器人平台，但我们认 为其在机器人声源定位领域具有重要的应用潜力。

本文中，我们介绍了：领域的发展（第2节）；声源定位问题的 定义（第3节）；在机器人听觉背景下声源定位技术的分类（第4节）； 用于声源定位的常用跟踪技术概述（第5节）；描述声源定位技术所 应对领域的多个方面（第6节）；当前用于评估声源定位技术性能的 不同评估方法综述（第7节）；以及对机器人听觉领域可能值得关注 的挑战的见解（第8节）。最后，我们强调了机器人听觉领域未来研 究问题的若干动机（第9节）。

2. 声源定位的发展

声源定位在机器人学中的兴起相对较新。据我们所知，它始于 1989年的Squirt机器人，这是首个配备声源定位模块[6,7]的机器人。 Squirt是一个小型机器人，具有两种相互竞争的行为：躲藏在黑暗处 和定位声源。利用声源定位作为驱动机器人交互行为的想法后来被布 鲁克斯自己的研究团队进一步探索，并最终为Cog机器人开发出了一 套SSL系统[8–11]。与此同时，一些日本研究人员也开始研究机器人 中声源定位的潜力。1993年，高梨等人研究了用于机器人的仿人听觉 系统[12,13]（如[10]所述）。此后，该领域取得了一系列重要进展： 千叶机器人[14],、基于RWIB12的机器人[15–18],、Jijo‐2[19,20],、 Robita[21]以及Hadalay[22]。第一代机器人应对了多种复杂场景， 例如人机交互，集成完整的听觉系统（源分离支持语音识别）、主动 定位，处理移动声源和采集系统，并通过探索不同的方法实现鲁棒声 源定位。

20世纪之交，机器人听觉的双耳子领域开始成为一项重要的研究 方向，其中包括声源定位。尽管第一代机器人在技术上具备双耳能力 （例如Squirt、COG、智也、Hadalay），但直到SIG机器人[23] 的出现，双耳机器人听觉领域才开始引起广泛关注。SIG机器人旨在 推动听觉作为机器人的基本技能，并被用作RoboCup人形挑战赛 2000[24]的实验平台。这使得SIG在机器人感知研究中变得流行起 来。此后，其他研究团队也相继跟进双耳机器人听觉的研究，该领域 的进展持续不断[25–36]。

在2000年代，机器人听觉领域在声源定位技术的研究动机方面出 现了重要分歧。双耳听觉基于模仿自然的动机而确立：仅使用两个耳 朵/麦克风。另一方面，提升性能的动机（详见第4.3节）推动了更多 麦克风的使用。这为使用大量传感器的源定位技术（如MUSIC和波 束形成器）在机器人中实现声源定位打开了大门。随后，声源定位问 题的各个方面得以拓展，机器人听觉领域因此产生了多种多样的解决 方案。

在其发展历史中，配备SSL系统的机器人核心目标一直是支持与 人类的互动。在第一代机器人中，一个重要的贡献是能够面向用户， 因为这表明机器人正在关注用户。智也机器人[14]是最早实现基于 注意力交互的机器人之一，其设计理念已融入近期产品如Paro机器人 [37]中。此后，声源定位被应用于更加复杂的场景，其中多种技能 相互结合以实现特定目标，例如：玩马可波罗游戏、充当服务员、提 供辅助并在视觉丢失用户时寻找用户[38]；在与护理人员互动的同 时记录并检测特定声音的来源[39]；玩简化版的捉迷藏游戏，通过手 势检测和声源定位来引导游戏进程[40]；在远程存在场景中，利用 声源提供视觉线索作为补充[41]；以及主持类似问答游戏的活动[42]。 鉴于机器人声源定位的发展，我们确信场景的复杂性将持续增长。事 实上，我们预见未来面临的挑战必将更加严峻（详见第8节的进一步 讨论）。

3. 声源定位问题定义

声源定位（SSL）通过音频数据单独解决估计声源位置的问题。 这通常涉及多个阶段的数据处理。其流程总结于图1。 由于该流程直接从麦克风接收数据并提供声源定位估计，我们考 虑一种端到端的方法来实现这一过程。端到端首先从输入信号中提取 特征，然后进行特征到位置映射，该映射通常依赖于声音传播模型。 这三个阶段在第4节中对每种方法及其相关变体的解释中均以此方式 引用。 本节简要概述了这三个阶段，以便在后续详细说明中便于参考。

3.1. 传播模型

声音传播模型的提出取决于以下因素：麦克风布置，因为其间可 能存在物体；机器人应用，因为用户可能距离麦克风阵列很近或很远； 以及房间特性，因为它们决定了声音在环境中的反射方式。此外，传 播模型通常决定了所使用特征的类型。

最常用的传播模型是自由场/远场模型，该模型假设以下条件：
• 自由场：每个声源发出的声音通过一条直射路径到达每个麦克 风。这意味着声源与麦克风之间以及麦克风之间没有物体，并且假设 环境中没有反射（即无混响）。
• 远场：麦克风间距与声源到麦克风阵列的距离之间的关系使得 声波可被视为平面波。

第二个假设极大地简化了特征与位置之间的建图过程，如第4.1节所 述。在机器人学的声源定位中，还有其他类型的传播模型具有相关性。 伍德沃思‐施洛斯伯格球形头模型[43,pp. 349–361]已被广泛应用于 安装在机器人头部的双耳阵列[23,44]，并在第4.2节中进行了解释。 近场模型[45]假设用户可能靠近麦克风阵列，此时需要将声波视为 圆形波。有一些机器人应用使用了近场模型，例如[46],然而，其使 用频率不如远场模型普遍。事实上，有些方法在近场情况下成功地采 用了改进的远场模型[47]或修改方法设计以考虑近场情况[48]。然 而，如[48],所述，直接将在远场模型应用于近场情况可能会显著降 低声源定位性能。此外，也存在一些传播模型被学习的情况，例如基 于神经网络的方法、[49,50],流形学习[33,51],线性回归[52]以及 作为多模态融合的一部分[11,21]。

3.2. 特征

在所回顾的方法中使用了多种声学特征。在本节中，我们简要概 述最常用的几种：

到达时间差（TDOA） 。它是两个采集信号之间的时间差。在使用外 部耳廓的双麦克风阵列（双耳阵列）中，该特征有时也称为双耳时间 差（ITD）。计算该参数的方法有多种，例如测量信号过零点时刻之 间的时间差[18]或从每个信号计算得到的起始时间之间的时间差 [6,7,14,17]。另一种计算TDOA的方法是假设声源信号为窄带信号。 设频率f处两个信号的相位差为1ϕf。若fm为具有最高能量的频率， 则窄带信号的TDOA（等效于麦克风间的相位差，即IPD）可通过 1ϕ fm 2π fm[23]求得。然而，截至目前最常用的TDOA计算方法是基于互 相关技术，详细内容将在第4.1节中说明。

麦克风间强度差（IID） 。它是给定时刻两个信号之间能量的差异。 当从时域信号中提取该特征时，可用于确定声源位于双麦克风阵列 的右侧、左侧还是前方。为了提供更高的分辨率，需要使用多麦克 风阵列[53]，或采用基于学习的映射方法[10]。IID的频域版本是 麦克风间声级差（ILD），表示为两个短时频变换后采集信号之间的 差谱。该特征也常与基于学习的映射方法结合使用[35]。

与麦克风间声级差（ILD）类似的特征是在频域中对数间隔分布 的一组滤波器（称为滤波器组）的输出差异。这些特征集具有 在使用维度较少的特征向量的同时，表现出比IID[9],更强的抗噪声 鲁棒性，而相较于麦克风间声级差（ILD）则具有更少的尺寸。 在[54],中，ILD是在泛音域中计算的。当频率fo是另一个f 的泛音时，f o= rf（假设r ∈[2, 3, 4,…]），且它们的幅值在 时间上具有高度相关性。这种方法可能对干扰更具鲁棒性，因为 频率之间的相关性意味着它们属于同一声源。

谱凹陷 。当使用外部耳廓1或内耳道时，麦克风信号之间存在轻微的 不对称性。由于这种不对称性，对其信号进行相减后的结果会在某些 频率上出现衰减或放大，这些频率取决于声源方向。通过实验可以将 这些凹陷映射到声源方向[52]。然而，由于外部耳廓的微小变化可 能会阻碍这些观察结果的有效性，因此在使用这些特征时建议采用基 于学习的映射[49]。

双耳/谱线索 。这是一个常用术语，指由IPD和ILD组合而成的特征集。 该特征集通常与基于学习的映射一起使用[50,51]。它们通常在声音 起始时刻被提取，以减少混响的影响[55]。实践中已经表明，对该特 征集进行时间平滑可使最终的建图对中等程度的混响更具鲁棒性[56]。

除了这些特征外，还有其他一些广泛使用的特征，例如 MUSIC伪谱和波束形成器导向响应。然而，它们的应用局限于特定 的端到端方法。因此，其详细说明见第4节。

3.3. 映射方法

声源定位的映射过程旨在将给定的提取特征映射到相应的位置。 实现这一目标的典型方法是直接应用传播模型，例如自由场/远场模 型或伍德沃思‐施洛斯伯格球形头模型，这两种模型在第3.1节中已讨 论。然而，某些类型的特征（特别是用于多声源定位估计的特征）需 要对声源定位解空间进行搜索或优化。一种常见的方法是进行网格搜 索，即在整个声源定位空间中应用映射函数，并记录每个测试声源位 置的函数输出。由此生成一个解谱，其中峰值（或局部最大值）被视 为声源定位解。这是多声源定位估计中最常用的映射过程类型。两个 重要的例子是MUSIC的子空间正交性特征和延迟求和波束形成器的 导向响应。这些将在第4.3节中进一步详细说明。

除了网格搜索之外，还有其他类型的映射方法。这些方法的主要 目的是基于已知位置声源的记录数据来训练映射函数。因此，所学习 到的映射函数隐式地编码了传播模型。在本综述中，这类映射方法被 称为基于学习的映射。它们基于不同的训练方法，例如神经网络 [11,21,49],局部线性回归[57],流形学习[33,51],等。每种映射过 程的更多细节将在第4节所述的方法分类的相关分支中给出。

4. 声源定位端到端方法的分类

如前所述，声源的位置通常被认为由两部分组成：（1）声源的 到达方向（DOA）；（2）声源到麦克风阵列的距离。大多数（即使 不是全部）相关研究采用的参考坐标系如图2所示，其中展示了一个 安装有三麦克风阵列的机器人。

可以看出，麦克风阵列的中心通常被视为原点。方位角方位角平 面（如图2b所示）与物理世界中的地平线平行，而仰角仰角平面（如 图2c所示）则与其垂直。这是本综述中使用的术语和参考坐标系。

在本节中，介绍了机器人学领域用于声源定位的端到端方法的分 类。由于声源位置由两部分组成，一种分而治之策略的声源定位方法 变得流行：分别进行到达方向和距离估计。在大多数情况下，通常仅 报告位置中的到达方向部分。鉴于这种流行趋势，本节主要回顾了对 声源到达方向进行估计的方法。然而，距离估计方面也取得了重要进 展，值得单独作为一个分支进行讨论。

所提出的分类总结如下：
• 1‐维单一声源到达方向估计。在此分支中，描述了在方位平面 上估计单一声源到达方向的技术。
• 2‐维单一声源到达方向估计。在这一分支中，描述了用于估计 单个声源在方位角和仰角平面内到达方向的技术。
• 多到达方向估计。在此分支中，描述了估计多个声源到达方向 的技术。这些技术大多位于方位平面上，但 如何将它们推广到两个平面的过程也进行了描述。该分支进一 步分为三个子分支：
– 基于波束形成的。针对多个波达方向候选进行空间滤波 的方法。– 子空间方法。利用信号子空间与噪声子空间之 间差异的方法。– 时域声源聚类。在不同时间窗口内进行 单个到达方向估计，并通过对这些结果进行聚类来提供多 波达方向解的方法。
• 距离估计。在这一分支中，描述了用于估计声源到麦克风阵列 距离的技术。

已通过一个Excel文件提供了所审查的188项研究的完整列表， 该文件是本文附加外部材料的一部分。

4.1 一维单到达方向估计

许多研究旨在定位和跟踪环境中的单个声源。为实现此目标，最 常用的特征是传感器或麦克风对之间的时延差（TDOA）。截至目前， 估计TDOA的最流行方法是基于两个采集信号之间的互相关向量（ CCV）进行计算。计算CCV的一种最简单方法是基于皮尔逊相关系数， 如公式(1)所示：

$$ CCV[\tau]= \frac{\sum_t(x_1[t] -\bar{x}_1)(x_2[t - \tau] -\bar{x}_2)} {\sqrt{\sum_t(x_1[t] -\bar{x}_1)^2}\sqrt{\sum_t(x_2[t - \tau] -\bar{x}_2)^2} $$

其中x1和x2是被比较的两个离散信号； τ是x2进行线性移动并计算 相关性的点；而x1和x2分别是x1和x2,的均值。声源的时延差（τo） 是使 τ互相关向量CCV取得最大值的值。如前所述， 自由场/远场传播模型是最常用的（见公式(2)），它提供了 τo与时延 差之间的简单特征到位置映射（θo）：

$$ \theta_o= \arcsin\left( \frac{V_{\text{sound}} \cdot \tau_o}{f_{\text{sample}} \cdot d} \right) $$

其中，Vsound为声速（∼343 m/s）；fsample为采样频率（单位： Hz）；d为麦克风间距（单位：米）； τo为声源的时延差（以采样 点数表示）。为了简化实际环境中的声源定位，麦克风被布置成它们 之间的假想连线与方位平面平行，从而使得 θo为该平面内的到达方 向。如果需要估计仰角，则可将麦克风对布置为跨越仰角平面。然而， 该方法仅能估计麦克风阵列所在平面内的到达方向。

通过互相关向量估计时延差对混响和其他噪声源非常敏感[59, 第213页–215]。在这些情况下，相关性值会“扩散”到其他时延差中 [60],导致相关性出现宽广的峰以及到达方向估计误差[61]。为了应 对这一问题，可以通过另一种方法计算类似的互相关向量。基于频域 的互相关器（CCF）[62]在公式（3）中给出：

$$ CCF[f]= X_1[f]X_2[f]^* $$

其中X 1 和X 2 分别是x 1 和x 2 的傅里叶变换，见公式(1)；{.} ∗表示复共 轭运算符 操作；且CCF是基于频域的互相关器。需要指出的是，所得 F−1(CCF) 呈现的相关性信息方式与互相关向量不同。然而， F−1(CCF)中的峰 值位置与互相关向量在 −τmax ≤ τ ≤ τmax范围内的高相关性位置 相同，其中 τmax是 τ在物理上可达到的最大值。

通过在频域中执行此操作，可以应用如公式(4)所示的加权函数， 该方法称为广义互相关（GCC）[60]：

$$ GCCF[f]= \psi[f]X_1[f]X_2[f]^* $$

ψ[f]根据相关向量的目标不同而变化。如果 ψ[f]= 1，则得 到的 F−1(GCCF)等同于 F−1(CCF)，其对混响敏感，类似于通过公式 互相关向量由公式(1)计算得出的情况。因此，值得关注的是，狄拉 克δ函数仅在相关性较高的位置出现在相关向量中。由于狄拉克δ函数 的傅里叶变换在所有频率上的幅值均为1，通过对GCCF中的幅值进行 归一化，可迫使在高相关性位置出现狄拉克δ函数的近似。为了归一 化这些幅值，ψ[f]被设为信号乘积幅值的倒数，如公式(5)所示：

$$ \psi[f]= \frac{1}{|X_1[f]X_2[f]^*|} $$

通过在ψ[f]应用公式(4)，相位信息在GCCF中得以保留，因此 公式(5)中的ψ[f]被称为相位变换（PHAT）。带相位变换的广义互相关[60]（GCC‐PHAT）如公式(6)所示：

$$ PHAT_F[f]= \frac{X_1[f]X_2[f]^ }{|X_1[f]X_2[f]^ |} $$

执行此归一化等效于对输入信号进行“白化”，因为所有频率的 幅值均为1。这已被证明会产生更尖锐的互功率谱[63]。因此，由实 际声源或环境反射（即混响）产生的干扰源也倾向于在相关向量中表 现为其他峰值。这会削弱它们在其他时延差上的相关计算中的影响， 而基于皮尔逊的方法（见公式(1)）则不具备这一特性。如[61]中的 仿真所示，这使得广义互相关‐相位变换对混响具有鲁棒性。此外， 在多源场景中[64]所示，其在高信干比(SIR)情况下也对干扰源具有 鲁棒性。

PHAT加权函数通常在整个频带内均匀应用，这会引入对宽带噪 声源的灵敏度。为解决这一问题，可以修改PHAT加权函数，使其根 据频带的“噪声程度”设置附加权重。[47,65,66],中给出了此类方 法的良好示例，其中根据频带信噪比在PHAT加权函数中增加了一个 附加加权项，从而提高了抗噪声鲁棒性。然而，在频带上应用非零权 重可能导致噪声泄漏到其他频率中。为避免此问题，[67],提出了一 种改进方法，采用硬二值掩膜，即仅向PHAT加权函数添加二值权重： 当信噪比超过某一阈值时为1，否则为0。遗憾的是，使用这种硬掩膜 会导致到达方向估计中出现泄漏，并产生不希望出现的主导峰。该问 题在[68],中通过在噪声和语音窗之间使用过渡掩膜得以解决。

需要指出的是，这些掩码方法需要在线噪声估计，以计算每个频 率的窄带信噪比。另一种方法是创建一个二值掩码，仅将应用程序所 针对的声源类型所使用的频率范围之外的频带置零。在[69],的情况 下，作者旨在仅跟踪语音声源，因此将语音频率范围之外的所有频带 均置零。

广义互相关‐相位变换（GCC‐PHAT）可能是机器人听觉中用于 单到达方向估计的最常用的时延差估计技术，因其具有良好的鲁棒性 和易于实现的特点。例如，在[70,71],中，GCC‐PHAT被用于通过 机器人探索来构建环境的声学地图。其他将GCC‐PHAT作为服务机 器人声源定位系统组成部分的研究可在[28,72–76]中找到。

有趣的是，在代表目标信号的频段之外的PHAT频段中出现的峰 值可能构成其他声源，这些声源被视为干扰。因此，一些研究提出使 用PHAT作为一种简单的方法来估计多个到达方向[64,77]。然而， 即使通过修改ψ[f]使其适用，峰值的出现仍取决于多个声源之间的 功率比[64]。据我们所知，这种GCC‐PHAT技术的变体尚未在机器 人听觉背景下得到应用。

如前所述，自由场/远场传播模型是单波达方向估计中最常用的 模型，但也可使用其他声音传播模型进行一维单到达方向估计。在 [78],假设球形头部位于麦克风之间。为此，作者使用了伍德沃思‐施 洛斯伯格头部模型[43,第349‐361],见公式(7)：

$$ \tau(\theta)= \frac{d}{2V_{\text{sound}}} (\theta+ \sin(\theta)) $$

其中d/2表示头部的半径。此时可观察到两条传播路径：一条通过头 部前方传播，另一条通过头部后方传播。尽管使用前向传播路径已足 以进行时延差估计，但后向传播路径会对此估计造成干扰。为抵消此 影响，引入了多径干扰补偿因子。所得的传播模型见公式(8)：

$$ \tau(\theta)= \frac{d}{2V_{\text{sound}}} (\theta+ \sin(\theta))+ \frac{d}{2V_{\text{sound}}} (\text{sign}(\theta)\pi -2\theta)|\sin(\theta)| $$

其中sign(θ的定义见公式(9)：

$$ \text{sign}(\theta)= \begin{cases} -1, & \theta< 0 \ 1, & \theta \ge 0 \end{cases} $$

在[65],对公式(8)中的传播模型进行了补充，以考虑一个衰减 因子 βm（通常按作者建议设为0.1），如公式(10)所示：

$$ \tau(\theta)= \frac{d}{2V_{\text{sound}}} (\theta+ \sin(\theta)) + \frac{d}{2V_{\text{sound}}} (\text{sign}(\theta)\pi -2\theta)|\beta_m \sin(\theta)| $$

在[44],中，作者从听觉极线几何（AEG）的角度得出了与公 式(2)和(7)中所提出的相同模型。[23]极线几何在立体计算机视觉中 被广泛使用，用于物理定位从两个已知位置的摄像头同时捕获的两幅 图像中提取的特征。[79]在[44]中对AEG进行的修订（RAEG）类 似于公式(7)所提出模型的傅里叶变换。作者应用了以下网格搜索映 射来执行一维声源定位：（1）使用RAEG，针对每个可能的f和到达 方向计算一组麦克风间的相位差（IPDf）；（2）根据选定f值下的输 入信号估计一个IPDf ；以及（3）

最终的 θo是与IPDf相关联的那个，该IPD来自与估计的IPDf最相似的集 合。

或者，可以使用机器学习方法来解决声源定位问题。这些方法不 定义声音传播模型，而是从声音信号的特征空间 κ到声源位置 θ学习 一个映射函数 Ψ。

$$ \Psi_\varphi: \kappa \mapsto \theta $$

该映射基于用于识别参数 φ的录制样本（训练数据集）。传统上， 会记录额外的样本来评估训练好的模型（测试数据集），并期望其在 真实世界环境中表现良好。一些已用于声源定位的训练技术包括：循 环神经网络[80–83]；仿生脉冲神经网络[84]；以及最近正在探索的 深度学习架构[85–87]。

然而，为了训练模型，训练数据多样性对于泛化至关重要[88]。 这是因为训练过程需要“观察”一组代表性样本集，从中尝试推广到 所有可能的情景。如果模型在真实世界环境中部署的情况与其训练时 的情况大不相同，则大多数机器学习技术的一个基本假设就被忽略了 [89]。例如，如果训练数据是在安静的办公室中录制的，而机器人被 部署在繁忙的餐厅中，噪声水平的差异将导致较差的声源定位性能。 为了避免这种情况，需要一个高度多样化的训练数据集，其变化因素 包括：房间响应、声源数量和位置、麦克风布置、外部耳廓以及/或 耳道的使用、头部模型等。在[90],提出了一种既能定位又能识别用 户的方法，作者建议将评估其在训练与测试条件不匹配情况下的声源 定位性能作为未来工作。有关此现象如何影响SSL系统性能的示例可 参见[57],详见第4.2节。

4.2. 二维单到达方向估计

在三维环境中，需要估计声源在方位角和仰角平面上的到达方向。 如果麦克风阵列具有三维几何结构，则可以使用第4.1节中讨论的相 同基于时延差的技术来计算仰角，并结合方位平面上的到达方向，实 现二维到达方向估计。然而，这需要使用额外的麦克风。

在双耳听觉中，麦克风通常被放置在两者之间的身体（打破自由 场假设）和/或配有外耳或耳道。这样做是为了使声源信号在被麦克 风捕获之前以依赖于其在两个平面上方向的方式进行‘滤波’。尽管 自由场假设被打破，这种滤波效应仍可与非自由场传播模型结合使用， 以实现二维到达方向估计。

例如，人类即使使用双麦克风阵列[91]，也能够估计声源的方位 角和仰角。当位于人头一侧的声源发出声波时，离声源较远的耳朵接 收到的声音是较近的耳朵所接收声音经过修改后的版本。这种修改是 在声波通过和/或被3人头、躯干、外部耳廓以及内耳道[92]反射时， 由多个可测量的物理现象共同作用的结果。这些现象可以通过放置

当音频信号穿过身体时所观察到的物理现象通常出现在低频率中，而反射的物理 现象则出现在高频率中。

在专用人工头的两侧各放置一个麦克风，并在消声室中捕捉其脉冲响 应（以避免环境影响）。这些测量可用于计算麦克风之间的物体如何 根据声源位置（特别是其方位角和仰角）改变所接收信号的特性。这 一组特性可用于创建所谓的头相关传递函数（HRTF），这是一种滤 波器，旨在模拟给定预设二维到达方向[93]时改变声源音频波的物理 现象。HRTF可用于多种应用，例如空间音频重放，其中声源被虚拟 地“定位”在佩戴耳机的听者周围[94]。为了实现单到达方向估计， 提出一组候选二维到达方向，并针对这些方向测量一组头相关传输函 数（HRTFs）。然后，可以从每个HRTF到到达方向的对应关系中计 算出执行HRTF逆函数的一组滤波器数据库（此处称为IHRTF）。随 后，这些IHRTF可通过网格搜索映射应用于输入信号。每次应用后可 测量一个提出的度量，例如输出信号之间的相关性。使该度量最大化 的IHRTF所对应的到达方向即被作为估计的到达方向[95]。需要注意 的是，外部耳廓在HRTF估计中也可能起到重要作用，如[96]所示， 其中使用螺旋耳进行二维到达方向估计。

然而，这类方法需要测量数据库‐IHRTF（通常基于在低噪声消 声室中通过测量估计得到的HRTF）在实际环境条件下采集信号时的 影响（实际情况通常包含混响和噪声）。为应对这一问题，可以预先 测量房间特性，并将其纳入数据库滤波器的计算中（在概念上不同于 前述的IHRTF），从而使滤波器与真实世界环境保持一致。遗憾的是， 将机器人部署到真实世界环境中意味着声学特性无法预先获知[97]。 此外，房间响应依赖于麦克风阵列在房间内的位置。这意味着当机器 人移动（无论是线性移动还是旋转）时，必须重新测量房间响应，并 重新计算滤波器数据库以适应这种变化[98]。可以通过仿真来模拟此 类真实世界环境[99,100],，从中获取房间特性以实现滤波器数据库 修改的自动化。然而，进行此类仿真可能不可行，因为环境的尺寸和 材料可能无法预先获知。在已知这些信息的情况下，此类仿真也十分 耗时[44]。

改进的听觉极线几何（RAEG）[44],如第4.1节所述，可用于克 服真实环境与训练条件之间差异的问题。基于RAEG的声源定位无需 在消声室中估计一组IHRTF，因为机器人头部的形状已被考虑在内。 然而，它仅提供方位平面上的到达方向。通过散射理论，可将其推广 至二维定位。

散射理论是物理学的一个领域，用于模拟波或粒子在传播路径中 因遇到物体而产生的扰动[101,102]。对于声学或电磁波而言，这种 建模过程可用于源定位[103]。如[97],所示，通过假设麦克风之间的 物体形状为球形，可利用给定的IPD和IID，结合散射理论方程来计 算方位角和仰角到达方向。在[97]中采用了一种网格搜索映射方法 （类似于前述基于IHRTF和RAEG的方法），该方法使用接近球形的 机器人头部置于麦克风之间。该方法的步骤概括如下：（1）针对每 个可能的方位角到达方向、仰角到达方向和f组合，利用ˆ ˆ 散射理论方程计算一组[IPDf, IID f ]元组；（2）估计出一个IPDf和 IIDf元组

来自所选f中的输入信号；（3）最终的二维到达方向是与估计元组 最相似的集合中对应的那个元组。

在[104],中提出了一种在二维空间中估计单个到达方向的另一种 方法，称为空间域距离（SDD）方法。该方法依赖于一种距离度量， 将其应用于捕获的时频（TF）单元与计算出的时频（TF）单元之间， 后者用于估计如果信号位于某一给定方向时所捕获的信号。该度量基 于球面傅里叶变换（SFT）[105],，其目的是以球坐标系来研究具 有旋转对称性的物体，例如来自球形机器人头部的音频信号扰动。 SFT系数通过测量或仿真进行
来自所选f中的输入信号；（3）最终的二维到达方向是与估计元组 最相似的集合中对应的那个元组。

在[104],中提出了一种在二维空间中估计单个到达方向的另一种 方法，称为空间域距离（SDD）方法。该方法依赖于一种距离度量， 将其应用于捕获的时频（TF）单元与计算出的时频（TF）单元之间， 后者用于估计如果信号位于某一给定方向时所捕获的信号。该度量基 于球面傅里叶变换（SFT）[105],，其目的是以球坐标系来研究具 有旋转对称性的物体，例如来自球形机器人头部的音频信号扰动。 SFT系数通过测量或仿真进行估计，并可用于选取包含直达路径信号 信息的时频（TF）单元，前提是已知预定义的到达方向。在[104]中 提出的SDD度量方法，用于衡量由输入信号计算出的直达路径时频 （TF）单元与利用给定预定义到达方向的SFT系数得到的相同时频 （TF）单元之间的距离。基于该度量，通过网格搜索寻找使该距离 最小化的方向。该技术对混响具有鲁棒性，并且对声源的频率范围不 变。在[106]中测试旋转头部对单个到达方向估计的影响时也使用了 该方法。SDD的最小化假设仅存在一个声源，但作者认为该方法有可 能扩展至多到达方向估计。

先前描述的仰角估计方法使用了网格搜索映射。然而，也可以采 用其他类型的基于建模的方法将特征空间映射到仰角估计。例如，在 [57]中提出的方法通过应用一种称为概率分段仿射映射（PPAM） 的逆回归方法来隐式编码HRTF。在这种情况下，PPAM的目标是将 一系列特征映射到一个位置。该特征集的一部分包括整个频谱的 ILDs和IPDs。噪声估计通过时间平均进行，由此计算出活动性度量。 该活动性度量也被用作训练特征，以提高映射过程对自噪声、背景噪 声和低混响的鲁棒性。当麦克风移动时，该技术的方位角精度不受影 响，但其仰角估计的性能会显著下降。作者指出，为了使该方法对这 些移动具有鲁棒性，“需要结合来自不同真实和/或模拟房间的训练 数据”。值得一提的是，[57]中描述的方法能够定位两个声源，使其 成为一种多二维波达方向估计器。然而，这种情况需要相应的训练数 据。这意味着，若要使其能够通用地定位多个声源，则必须使用大量 数据进行训练。不过，这是使基于学习的二维声源定位在现实世界和 训练条件变化下更具鲁棒性的重要努力。

基于学习的二维自监督学习的更多示例包括：在[49]中，采用 受人耳启发的耳廓来提取谱凹陷，并将其输入神经网络；在双耳系统 中使用流形学习来估计方位角和仰角[33]，并已扩展以包含视觉信 息[107]；在[52]中采用线性回归方法作为映射过程，利用从带有螺 旋耳的iCub机器人头部提取的谱凹陷来学习传播模型。

4.3. 多方向到达估计

对于考虑听觉场景中存在多个声源的场景，可以采用多种多声源 到达方向估计技术，这些技术可分为三类：
基于波束形成的

子空间方法 时域声源聚类 为了便于描述，这些方法的一维版本被呈现出来。因此，估计的 到达方向被限定在方位平面内。然而，需要注意的是，当使用三维阵 列几何结构时，可以将这些方法扩展到二维（方位角和仰角）。

4.3.1. 基于波束成形

波束成形是一种基于传感器阵列所捕获信号的空间加权的滤波技 术，其输出为从预指定到达方向（DOA）接近传感器阵列的信号。使 用该技术进行到达方向估计可能看起来违反直觉，因为它假设声源的 到达方向是已知的。然而，该方法的整体步骤如下：
1. 提出一个包含L个候选到达方向的集合，用于搜索声源。2. 创 建一个波束形成器，并引导其指向每个候选到达方向。3. 测量每 个波束形成器的响应，通常以其输出的能量来衡量。4. 生成导向 响应谱。5. 在谱中寻找峰值，并将其位置作为声源的到达方向。

导向响应谱是一个由所有按到达方向排序的波束形成器响应构成 的 1 × L向量。如果以波束形成器的输出来衡量其响应，则导向响应 谱可被视为一种能量谱，用于显示麦克风阵列从各个到达方向接收到 的能量大小。显然，网格搜索映射是基于波束形成的多到达方向估计 中最常用的方法。

最简单的波束成形形式被称为延迟求和波束成形（DAS），如图3所 示。 假设采用自由场/远场声传播模型，声源的到达方向（DOA）与 传感器之间采集信号的时延差（TDOA）具有直接关系，如公式(2)所 示。因此，DAS波束形成器旨在对信号进行人工移位以补偿这种时间 差，然后将移位后的信号相加以获得其输出。该方法的基本思想是， 来自预定义到达方向（或导向方向）的采集信号分量在各个移位后的 信号中被对齐。这使得相对于其他方向上的其他声源，该分量在波束 成形器的输出中被增强（如图3上部所示）。由此可见，这种增强程 度与所使用的传感器数量成正比。

假设传感器为全向麦克风，则DAS波束成形器的输出可由公式(12) 表示：

$$ \hat{s} \theta[t]= \sum {n=1}^{N} x_n[t - \tau_n(\theta)] $$

其中， ˆsθ 是波束成形器的输出，表示从导向方向θ到达阵列的地面 真实信号sθ 的估计值；t 是时间帧；xn 是第n个麦克风接收到的信 号；N 是麦克风数量；而 τn(θ)是声源相对于导向方向 θ在第n 个 麦克风中的时延差。需要指出的是，时延差 τn 通常基于参考麦克风 进行计算。只要这些参数与阵列的几何结构和传播模型一致，该方法 就能适用于任意阵列维度。此外，它还能适应多种非自由场/非远场 传播模型，例如近场模型[45]和球形头模型[43, 第349‐361]页。

指向 θ的波束成形器输出的能量可以通过公式(13)计算： (Eθ)

$$ E_\theta= \sum_{t=1}^{T} \hat{s}_\theta[t]^2 $$

在公式(12)中，时间偏移是在时域中进行的。也可以通过对频域 中的采集信号进行操作来实现时间偏移，如公式(14)所示：

$$ X_{n\tau_n(\theta)}[f]= X_n[f]e^{-2\pi f \tau_n(\theta)} $$

其中Xn是xn的f为频率帧；且Xnτn(θ)是xn(t − τn(θ)的与导向方向相关 的、应用于采集信号的时移排列可表示为复数值N × F矩阵Wθ，如 公式(15)(15)所示：

$$ W_\theta= \begin{bmatrix} 1 & 1 & \cdots & 1 \ e^{-2\pi f_1\tau_2(\theta)} & e^{-2\pi f_2\tau_2(\theta)} & \cdots & e^{-2\pi f_F \tau_2(\theta)} \ e^{-2\pi f_1\tau_3(\theta)} & e^{-2\pi f_2\tau_3(\theta)} & \cdots & e^{-2\pi f_F \tau_3(\theta)} \ \vdots & \vdots & \ddots & \vdots \ e^{-2\pi f_1\tau_N(\theta)} & e^{-2\pi f_2\tau_N(\theta)} & \cdots & e^{-2\pi f_F \tau_N(\theta)} \end{bmatrix} $$

其中，f为频带；N 为麦克风数量；F 为频率窗大小；Wθ 为宽带导 向矩阵，其每一列代表一个窄带导向矢量。

波束形成器输出的傅里叶变换（S）可通过公式(16)构造：

$$ \hat{S} \theta[f]= W \theta[f]^H X[f] $$

其中，{.} H算子表示共轭转置4；X是一个N × F复值矩阵，其行 包含所有Xn；Wθ[f] 是W θ 的N × 1 复值列，包含频率f的所有 波束成形权重；X[f] 是N × 1 列，包含所有采集信号在频率f的 频域信息；S θ[f] 是指向 θ方向、频率f下的波束形成器复值

ˆ 输出。 S θ 是 ˆsθ 的傅里叶变换，使得 ˆsθ = F−1(

S θ)。 在频域中工作为W θ 提供了进一步优化以提升波束成形器性能的 可能性。例如，一组权重（Aθ MVDR），其在保持与Wθ 相同方向的同时 最小化波束成形器能量，如公式(17)所示：

$$ A_{\theta \text{MVDR}} [f]= \frac{R[f]^{-1}W_\theta[f]}{W_\theta[f]^H R[f]^{-1}W_\theta[f]} $$

共轭转置用于消除W θ 中的时延差，从而在将信号导向声源的到达方向时产生对 齐的信号。

其中，R[f] 是X[f]的N × N协方差矩阵。这就是著名的卡彭波束 形成器，也称为最小方差无失真响应（MVDR）[108]。在实现 MVDR时，一个重要考虑因素是R[f]的计算。一种常用的R[f] 估计 方法是样本协方差矩阵R[f],，它通过在T个时间窗口上对X[f,t − T: t]X[f,t −T: t]H 求平均值得到，其中t是当前时间窗口。这 是MVDR的局限性，因为在R[f] ≈ 收敛到R[f], 之前，它无法提供 可靠的结果，而这可能需要多个时间窗口才能完成。

能量Eθ来自公式(13)基于估计的 ˆsθ，这意味着声源定位的性 能取决于导向方向上的声源被增强的程度。因此，当波束形成器指向 声源时，其输出的信号与干扰比（SIR）会影响导向响应谱中对应峰 值的相对高度。如果SIR过低（即声源未被充分增强），则在谱中可 能难以找到该峰值。由于这种增强与所使用的麦克风数量有关，通常 采用麦克风数量较多的麦克风阵列（8或以上）（详见第6.3.1节）。

使用DAS波束成形进行多到达方向估计时出现的一个问题是，导 向响应谱中出现的峰值通常非常宽，导致网格搜索的分辨率很差。解 决此问题的一种方法如ManyEars项目所提出的[109],是将波束形成 器的输出首先用互相关向量（CCV）表示，如公式(18)所示：

$$ E_\theta= \sum_{t=1}^{T} \hat{s} \theta[t]^2= \sum {t=1}^{T} \left( \sum_{n=1}^{N} x_n[t - \tau_n(\theta)] \right)^2 = \sum_{n=1}^{N} \sum_{t=1}^{T} x_n[t - \tau_n(\theta)]^2 + 2 \sum_{n_1=1}^{N} \sum_{n_2=1}^{n_1-1} \sum_{t=1}^{T} x_{n_1}[t - \tau_{n_1}(\theta)]x_{n_2}[t - \tau_{n_2}(\theta)] = \sum_{n=1}^{N} \sum_{t=1}^{T} x_n[t - \tau_n(\theta)]^2 + 2 \sum_{n_1=1}^{N} \sum_{n_2=1}^{n_1-1} \text{CCV} {x {n_1},x_{n_2}} [\tau_{n_1}(\theta)- \tau_{n_2}(\theta)]. $$

最后一步是通过使用公式(1)中的互相关计算来完成的。如[60], 所示，广义互相关向量（如公式(4)所示）可以用来替代公式互相关 向量中的(18)。因此，[63]的作者转而使用公式相位变换中的(6)向 量来锐化所得导向响应谱中的峰值。

然而，由于相位变换执行了幅度归一化，每个频率帧在相关计算 中具有相同的贡献。这使得该过程对噪声敏感。为应对这一问题，如 [109],可应用频谱加权以降低窄带信噪比低的频率的贡献。为此， 在每个时间窗口内计算所有麦克风的平均功率谱密度，称为1 × F复 数值向量Xmean。然后，通过时间平均Xmean来估计频域噪声信号（称 为 1 × F复数值向量Y）。接着，将式（19）中给出的加权函数ψ[f ]应用于式（4）以进行相关计算：

$$ \psi[f]= \begin{cases} \psi_{\text{PHAT}}[f], & X_{\text{mean}}[f] \le Y[f] \ \psi_{\text{PHAT}}[f]\left( \frac{X_{\text{mean}}[f]}{Y[f]} \right)^\gamma, & X_{\text{mean}}[f]> Y[f] \end{cases} $$

其中ψPHAT[f]= $\frac{1}{| X_1[f]X_2[f]^*|}$是相位变换；而 0< γ<2是一个标 量，可通过校准来定义对具有高窄带信噪比的频率赋予多大权重。这 提高了抗噪声鲁棒性。

利用波束成形进行多到达方向估计的一个显著优势是，由于波束 形成器执行了一种源分离，因此可以从分离出的声源中提取特征，并 用作网格搜索映射的指标。例如，在[110]感兴趣的声源为人类的情 况下，采用语音相似性度量来平滑导向响应谱，使得最高的峰值可被 视为人声源而非噪声。人类元音声音具有特定的谱峰组合，即峰值特 征（J，一个 1×F维复数值向量）。该特征与其他类型的声音不同， 可用于语音活动检测[111]。通过采集多位说话人的元音录音语料库， 可以获得一组峰值特征模板[110]。这些模板仅在谱峰对应的频带处 取值为1，其余频带均为0。通过计算波束形成器在假设到达方向上的 输出（Sθ[f]）与峰值特征（J）之间的峰谷差（PVD），如公式（ 20）所示：

$$ \text{PVD}= \frac{\sum_{f=1}^{F} \hat{S} \theta J[f]}{\sum {f=1}^{F}J[f]} - \frac{\sum_{f=1}^{F} \hat{S} \theta(1 -J[f])}{\sum {f=1}^{F}J[f]} $$

PVD通过测量两个频谱在共同峰值周围的平均能量与其它频带的 平均能量之差，来衡量两个频谱之间的距离。

此外，在[110],另一种称为广义旁瓣相消器（GSC）的波束形 成器被应用于来自平滑导向响应谱的n个最佳波达方向候选的方向上。 这样做是为了增强“最佳”方向上的信号。GSC被用于[110]中， 因为它设计用于适应噪声/干扰变化[112]。图4展示了GSC的概览。 如图所示，GSC在其上支路中采用了延迟求和波束形成器（ DAS）。在下支路中，它进行噪声估计，然后从DAS输出中减去该噪 声估计值。噪声估计是通过应用一个阻塞矩阵来实现的，其目的是通 过减去延迟后的信号（本质上是一个反向波束形成器），估计出除指 定到达方向以外其他方向存在的噪声。得到的信号随后经过滤波，使 得它们在相加后能够提供一个可以直接从DAS输出中减去的噪声估计 值。这些滤波器通过最小均方方法，利用当前估计输出和过去的噪声 估计，随时间不断优化，从而实现对环境变化的适应。所提出的工作 在[110]中，是我们发现的唯一一项将GSC用作其多到达方向估计方法 一部分的机器人听觉工作。

另一个例子是[113], ，其中将基于频率的选择方法应用于波束 成形器的输出，以去除衰减噪声。计算基于延迟求和的导向响应谱， 并从中获得最高能量和第二高能量声源的到达方向。将最高能量波束 成形器输出中低于第二高能量对应值的频率滤除。该滤波后的输出从 波束形成器的输出中减去，并重复执行此过程，直到仅剩背景噪声。 这使得对低能量声源具有高灵敏度，而这些声源在存在高能量干扰时 难以检测。然而，隐含假设是这些低能量声源的能量水平高于背景噪 声水平。

使用波束成形进行多到达方向估计的一个重要缺点是，每个提议 的到达方向都需要一个波束形成器。根据所采用的波束成形技术，执 行该操作可能需要大量的计算资源[69]。

4.3.2. 子空间方法

该分支中最流行的方法之一是多重信号分类（MUSIC）[114]。 可以说，它是该分支所有其他方法的基础，因此本节主要重点描述 MUSIC。

MUSIC的核心思想在于搜索与代表感兴趣信号的子空间相交的 到达方向。考虑如公式(21)所示的采集信号模型：

$$ X= W_s S+ V $$

其中，X是一个N × F复值矩阵，其行包含所有的Xn；Xn 是频域 变换后的输入离散信号（xn），每个表示传感器n处的采集信号；S 是一个D × F复值矩阵（其中D为源信号数量），其行包含频域中 的所有源信号sm；V是一个N × F复值矩阵，其行包含各传感器在 频域中存在的噪声信号集合。Ws[f] 是一个N × D复值矩阵，用于 建模每个源信号（d）在给定频率f下相对于其到达方向的各麦克风 （n）的时延差（τn:d），如公式(22)所示：

$$ W_s[f]= \begin{bmatrix} 1 & 1 & \cdots & 1 \ e^{-2\pi f \tau_{2:1}} & e^{-2\pi f \tau_{2:2}} & \cdots & e^{-2\pi f \tau_{2:D}} \ e^{-2\pi f \tau_{3:1}} & e^{-2\pi f \tau_{3:2}} & \cdots & e^{-2\pi f \tau_{3:D}} \ \vdots & \vdots & \ddots & \vdots \ e^{-2\pi f \tau_{N:1}} & e^{-2\pi f \tau_{N:2}} & \cdots & e^{-2\pi f \tau_{N:D}} \end{bmatrix} $$

需要指出的是Ws[f]与W θ 矩阵类似，见公式(15)，但在本例中， 它建模了每个声源接收到的信号的时延差。可以看出，该方法假设信 号是窄带的，中心频率为f，而语音信号并非如此。我们目前将假设 信号是窄带的，但稍后会给出MUSIC算法的宽带变化示例。

MUSIC的目标是估计由Ws的列所张成的子空间。为此，首先对 采集信号在频率ff处的样本协方差矩阵R[f]（其计算方式与第4.3.1 节中详细描述的最小方差无失真响应波束形成器相同）进行

特征分 解，如公式(23)所示：

$$ \hat{R}[f]= Q[f]\Lambda[f]Q[f]^{-1} $$

其中Λ[f]是N×N复值对角矩阵，包含频率f的协方差的特征值， 并按降序排列，如公式(24)所示：

$$ \Lambda[f]= \begin{bmatrix} \lambda_1[f] & 0 & \cdots & 0 \ 0 & \lambda_2[f] & \cdots & 0 \ \vdots & \vdots & \ddots & \vdots \ 0 & 0 & \cdots & \lambda_N[f] \end{bmatrix} $$

Q[f]是一个N×N复值矩阵，其列包含频率f的协方差的特征向 量集合，按Λ中确定的顺序排列。如公式(25)所示：

$$ Q[f]= \begin{bmatrix} q_{1:1}[f] & q_{2:1}[f] & \cdots & q_{N:1}[f] \ q_{1:2}[f] & q_{2:2}[f] & \cdots & q_{N:1}[f] \ \vdots & \vdots & \ddots & \vdots \ q_{1:N}[f] & q_{2:N}[f] & \cdots & q_{N:N}[f] \end{bmatrix} $$

Q[f]可以划分为两个子空间，在索引 λl处，如公式(26)所示：

$$ Q_s[f]= Q[f][1: \lambda_l] \quad Q_v[f]= Q[f][\lambda_l+ 1: N] $$

其中，Qs[f] 是一个N × λl 矩阵，其列向量为表示信号子空间的 特征向量；Qv[f] 是一个N ×(N − λl)矩阵，其列向量为表示噪 声子空间的特征向量，两者均针对频率f。 λl通常设为Ws[f]的秩。 如果难以估计（例如在噪声环境中），则可将 λl设为高特征值集合 中最小特征值的索引，即那些不接近零的特征值。这基于一个假设： 噪声子空间的特征值明显小于信号子空间的特征值。另一种方法是在 对协方差矩阵进行特征分解之前先白化噪声，使得噪声子空间的特征 值接近1，而信号子空间的特征值大于 1[115,116]。该白化过程首 先通过在声源静音时捕获噪声V[f]来估计N ×N噪声协方差矩阵 K[f]（见公式(21)）。然后计算K[f] =V[f]=V[f]V[f]H，并执 行公式(27)中所示的特征分解。

$$ K[f]^{-1} \hat{R}[f]= Q[f]\Lambda[f]Q[f]^{-1} $$

由于K[f]^{-1}可以是任意类型的矩阵，只要能使K[f]^{-1} ˆR[f]矩阵 变为方形的，这就形成了所谓的广义特征分解(GEVD)。当该分解作 为MUSIC的一部分使用时（即GEVD‐MUSIC），已证明其对环境中 的非相关噪声具有较强的鲁棒性[115]。

完成子空间划分后，可以启动波达方向搜索，在该过程中会提出 多个波达方向候选。为此，针对每个候选N × 1向量b[θ,f]按如下 方式计算，如公式 θ所示，见式(28)：

$$ b[\theta,f]= \begin{bmatrix} 1 \ e^{-2\pi f \tau_1(\theta)} \ e^{-2\pi f \tau_2(\theta)} \ \vdots \ e^{-2\pi f \tau_N(\theta)} \end{bmatrix} $$

其中N为传感器数量；而 τ n( θ)是麦克风n相对于候选 θ的时延差。 这些 τ n( θ)值根据给定的传播模型计算得出。如前所述，自由场/远 场是最常用的传播模型，在本方法版本中也是如此。但与第4.3.1节 所述的基于波束形成的方法一样，也可采用其他传播模型，例如近场 模型[45]和球形头模型[43, pp. 349–361]。

由于特征向量之间的正交性，当b[θ,f]与Qv[f],中的特征向量正 交时，它表示信号子空间（Qs[f]）中的一个源信号，因此指向源信号 的方向。为了检验这种正交性，应用公式(29)：

$$ P_{\text{MUSIC}}[\theta,f]= \frac{1}{b[\theta,f]^H Q_v[f] Q_v[f]^H b[\theta,f]} $$

其中向量PMUSIC[f]表示频率f处的MUSIC伪谱。其峰值的位置代表声 源的到达方向，类似于基于波束成形的方法的导向响应谱。PMUSIC[θ, f]在公式(29)中分母为0时无定义，即当b[θ,f]与Qv[f]正交时。然而， 噪声与源信号完全不相关的情况极为罕见，因此通常只会遇到接近0 的值。这意味着PMUSIC[f]中的峰值往往具有非常大的数值。

需要指出的是，MUSIC要求至少有一个特征向量张成噪声子空 间，否则Qv[f]将为空。因此，MUSIC的一个基本要求是存在的声源 数量最多为N − 1个。为了克服这一限制，通常采用具有较多麦克 风（8或以上）的麦克风阵列（详见第6.3.1节）。

由于MUSIC假设源信号为窄带，因此不能直接适用于宽带的语 音信号。解决此问题的一种方法是针对每个频率f计算窄带MUSIC谱， 并将其按到达方向的平均值作为宽带MUSIC谱[117,118]。然而，当 使用抗噪声能力强的GEVD‐MUSIC变体进行实时处理时，这会导致 需要大量的计算资源。为应对这一问题，采用广义奇异值分解 (GSVD)替代广义特征分解(GEVD)，可提供MUSIC实现所必需的相 同正交性，如公式(30)所示：

$$ K[f]^{-1} \hat{R}[f]= Q_l[f]\Lambda[f]Q_r[f]^H $$

其中Ql[f]和Qr[f]是频率f的左奇异向量和右奇异向量，它们彼此 正交且为酉矩阵。GSVD‐MUSIC使用Ql[f]代替从公式Q[f]起的 (27)。关于计算成本，[116]表明，使用GSVD‐MUSIC每个时间窗 口的处理时间明显少于GEVD‐MUSIC。此外，在同一项工作中，不 再在整个可能的到达方向集合上进行网格搜索映射，而是首先进行粗 略搜索，然后在粗略搜索的峰值处进行精细搜索。作者将此方法称为 分层声源定位(H‐SSL)，并且可以实现实时下低于 1◦的分辨率。

在[119], MUSIC对噪声环境的鲁棒性通过以下方式增强：(1) 通过对频率帧之间的噪声相位进行线性回归来校正样本协方差矩阵以 去除污染；(2)通过执行基于子空间的信噪比估计来忽略噪声频率帧。 需要记住的是，相位数据 ˆ 存在于样本协方差矩阵R[f]中。所提出的方 法使用双麦克风阵列，这导致X[f]仅具有两行，因此只能定位一个声 源。然而，如后文所示，该方法可推广为多波达方向估计器。

第X1[行]被视为参考麦克风接收到的信号，因此，X2[f]=X1[ f]e − i2 πf τ。在此场景中，可以从协方差矩阵中提取相位数据，如公式 (31)所示：

$$ \hat{R}[f]= X[f]X[f]^H = \begin{bmatrix} X_1[f]X_1[f]^H & X_1[f]X_2[f]^H \ X_2[f]X_1[f]^H & X_2[f]X_2[f]^H \end{bmatrix} = \begin{bmatrix} |X_1[f]|^2 & |X_1[f]|^2 e^{i2\pi f \tau} \ |X_1[f]|^2 e^{-i2\pi f \tau} & |X_1[f]|^2 \end{bmatrix} $$

其中 ζ(f τ)=e−i2πf τ表示相位数据，这是一种f 与 τ之间的映 射关系。传统上，在无噪声环境中，这两者之间的关系接近线性。然 而，作者发现在线性噪声环境中这种线性关系会丧失，从而导致到达 方向估计出现误差。为此，作者提出采用以下步骤来校正相位数据： ˆ1.计算样本协方差矩阵R[f]。ˆ2.在每个频率f 处从R[f]中 提取相位数据 ζ(f τ)。3.回归阶段（重复执行直到误差达到预 定义的最小值）： (a)对相位数据进行一阶回归，使其接近作者所称的相位线。 这些是通过计算来自不同方向的观测信号的相位数据获得的， 在无噪声环境中应为直线。(b)更新低于其在估计相位线上 对应值的相位数据。(c)计算更新后的相位数据与估计相位 线之间的误差。 4.使用更新后的ζ(f τ)重构样本协方差矩阵。

此外，基于公式(32)中所示的指标进行频带选择：

$$ \epsilon[f]= \log\left( \frac{ \lambda_1 \hat{R}[f] }{ \lambda_2 \hat{R}[f] } \right) $$

ˆ 其中 λ1[f]是R[f]的第一个特征向量对应的特征值，而λ2[f]是 第二个特征值。由于在此场景下只能估计一个源，因此 λ1[f]表示信 号子空间的特征值， λ2[f]表示噪声子空间的特征值。因此， ϵ[f]表 示从计算R[f]的频率f处的窄带子空间信噪比（子空间信噪比）。可 通过预设阈值来忽略子空间信噪比较低的频带[119]。

尽管作者未明确说明，但该方法或许可以推广到使用两个以上麦 克风的阵列（并进行多到达方向估计）。子空间信噪比可以通过计算 信号子空间中特征值之和与噪声子空间中特征值之和的比值得到，如 公式(33)所示：

$$ \epsilon[f]= \log\left( \frac{ \sum_{n=1}^{\lambda_l} \lambda_n[f] }{ \sum_{n=\lambda_l+1}^{N} \lambda_n[f] } \right) $$

在[120–123]中可以找到将MUSIC作为服务机器人声源定位的一部 分应用的更多示例。

4.3.3. 基于时间的到达方向聚类

由于这类多到达方向估计器的特性，首先定义“时间窗口”的含 义非常重要，因为它是基于时间上的到达方向聚类（DOA聚类）的 声源定位方法的关键部分。前面各节所述的技术均假设提供给估计器 的所有信息都来自由过去样本数量（即时间窗口）定义的某一时间段。 基于此，前述估计器仅使用该数据提供一个多到达方向结果。

此处描述的技术通过多个时间窗口估计单个5到达方向，并将不同的 到达方向分组 簇，每个簇代表一个声源。一个很好的例子是[15],，其中如果某个到达方向 接近其上一次的到达方向，则将其分配给该到达方向的簇（称为“峰值轨迹”）。 [69]的作者通过使用广义互相关‐相位变换方法的一种变体，计 算三角形阵列[124]中三个麦克风对之间的时延差，从而执行冗余的 单个到达方向估计。利用其冗余度量以及同时语音的非重叠特性，该 单个到达方向估计器即使在多说话人场景下也能从一个源估计出单个 到达方向。在时间窗口t计算出的到达方向可能与在时间窗口t+1 估计出的声源不同。为此，采用一种基于雷达跟踪技术的聚类方法6 [69]，将连续估计出的单个到达方向聚类形成轨迹，每条轨迹代表 一个声源。在每个时间窗口，将每条轨迹的所有到达方向输入到卡尔 曼滤波器中。所提出的系统通过对所有轨迹进行到达方向预测，实现 多到达方向估计。在[122],中也使用了类似的基于卡尔曼的方法， 但仅定位了一个声源。

另一种由[126],提出的方法通过基于高斯混合模型的跟踪器进 行多到达方向估计。该方法使用概率假设密度滤波器（GM‐PHD） 持续估计所有当前说话人