一种麦克风阵列声源空间实时定位方法
【专利摘要】本发明一种麦克风阵列声源空间实时定位方法，以麦克风阵列作为信号采集和输出设备，通过使用可控功率响应?相位变换法初步给出声源空间位置候选点；通过先验知识进行初步的候选点筛选，并使用可控功率响应?相位变换法计算候选点的可控功率响应输出；用改进的随机区域收缩重新确定搜索边界，提高可控功率响应?相位变换法的效率；最后计算剩余候选点的可控功率响应，选取最大的位置作为最终的声源估计位置；本发明声源定位原理明确，实时性好，实验证明该方法在平面上的定位误差范围可控制在厘米量级上，性能优于基于现有技术的方法；具有较高的运算速度和鲁棒性，可应用于智能家居和智能机器人等需要进行实时声源定位的场景。
【专利说明】
一种麦克风阵列声源空间实时定位方法
技术领域
[0001] 本发明属于语音技术领域，特别涉及一种麦克风阵列声源空间实时定位方法。
【背景技术】
[0002] 语音处理技术中通常需要掌握声源的空间位置信息，为声纹识别、语音内容识别 等交互技术提供帮助。例如，当用户与智能机器人进行语音对话时，通常要确定人的空间位 置并进行转向和靠近，即使在黑暗中也可以仅仅根据声音确定用户的具体方位，并自动靠 近说话人附近，提供必要的服务和帮助;在视频会议中一般需要自动调整会议摄像头，朝向 发言人，并能够对不同发言人进行视频切换。
[0003] 常见的用于声源定位的设备是麦克风阵列，一般定义为：由多个(通常大于三个） 麦克风按照指定的几何规则摆放并完全同步采集声音信号的设备。目前麦克风阵列定位的 常用方法主要包括:基于到达时间差(Time Delay of Arrival，TD0A)的声源定位，定位首 先通过时延估计得到声源信号到达不同阵元的时间差，再通过麦克风阵列的几何构造进行 声源位置判断；基于可控功率响应（Steered Response Power,SRP)的声源定位以及基于高 分辨率谱估计(High-resolution Spectral Estimation)的声源定位等。
[0004] 1、到达时间差（Time Delay of Arrival，TD0A)法：
[0005] 该方法首先通过时延估计得到声源信号到达不同阵元的时间差，再通过麦克风阵 列的几何构造进行声源位置判断。常用的时延估计方法是基于广义互相关(Generalized Cross Correlation,GCC)法，其具体步骤为，设Xi( ω )和Χ」（ω )分别代表第i个麦克风接收 信号Xl(t)和第j个麦克风接收信号&(t)的傅里叶变换，则有如下等式：
[0009] 其中Ψ^(ω)表示两个信号的互功率谱密度，表示两个信号的互相关函数， 为声源到第i个麦克风与第j个麦克风之间的延时估计，W^( ω )是频域加权系数，当此加权 系数为1时，就是基本互相关。常用的频域加权函数还有相位变换(PhAse Transform，ΡΗΑΤ) 加权、ROTH加权、最大似然加权、SCOT加权等，其中PHAT加权较为常用，因为基于广义互相关 的时延估计信息隐藏在互功率谱的相位信息中，与幅度无关，PHAT加权法则去除了幅度信 息的影响，使得相关函数峰值更加尖锐。PHAT加权函数公式为Wij( ω ) = 1/| Ψ^( ω ) |。
[0010] 上述基于到达时间差的定位方法原理简单，计算效率高，但在较大的噪声或混响 干扰下其时延估计性能急剧下降，因此不适用于低信噪比或高混响的声学场景。此外，还有 一些时延估计方法，比如利用信息论中联合熵和互信息，使用迭代方法，或者考虑不同信道 的非对称性等方法，但在强混响和噪声条件下，结果都不太理想，且部分方法计算速度缓 慢，比如迭代法等。
[0011 ] 2、可控功率响应-相位变换（Steered Response Power-PhAse Transform，SRP_ PHAT)法
[0012] 基于可控功率响应-相位变换的定位方法通过对全局空间进行搜索，找到最大可 控功率输出的点，网格越小，搜索分辨率越高。其声源定位原理如下：
[0013] 假定声源位于s处，τι表示声源到第i个麦克风的传输延迟。声源s的可控功率响 应-相位变换输出P PHAT( s)为：
[0015] 其中=不(α>)Χ；(ω)表示两个信号的互功率谱密度，表示s处声源到第j个 麦克风和第i个麦克风的时延差。估计的声源位置:?为：
[0016] s = argmax ρ1'1"1 (Λ)
[0017] 经过以上分析可以知道，SRP-PHAT算法的思路是假定空间中某一位置存在声源， 然后求出该位置的基于相位加权的可控功率输出，最终选取使得P PHAT(S)最大的空间位置8 作为估计出的声源位置。但是SRP-PHAT的思路是先假定空间中某一位置有声源，为了精确 定位出声源信息，一般需要用网格搜索法，遍历空间中所有的网格位置。假设搜索空间c长 宽高分别为X，Y，Z，网格间隔为δ，则网格点数N为：
[0019] 假定X = 3m，Y = 3m，Z = 2m，网格间隔δ = 〇·〇1πι，则网格点数Ν=301*301*201 = 18210801，也就是说要计算接近1800多万个点的⑷，即使只有二维，也要计算301*301 = 90601次，因此效率很低。
[0020] 3、随机区域收缩法
[0021 ]基于SRP-PHAT的随机区域收缩（Stochastic Region Contraction，SRC)方法是Do H，Si lverman H F为了提高SRP-PHAT效率在2007年提出的一种算法，该算法主要思路是通 过随机采样，选取基于PHAT可控功率响应最大的一部分点，根据这部分点重新确定搜索边 界，然后在新的区域内重复上面操作，直至找到全局最大值。该算法鲁棒性的前提是随机采 样时，错失最高峰的概率几乎为零，也就是没有采样到声源所在位置的概率几乎为零。假设 声源体积为Vs，房间体积为Vr。?，则每次随机米样米到峰值的概率Phit为Phit = Vs/Vr。?，错失 峰值概率pmis^phit = Vs/Vr_。假设第一次采样共采集到Μ个点，则错失峰值概率Pmiss(N)为 (从）=Λ!，同时，为了确保第一次采样错失峰值概率低于某个阈值Pthre，采样点个数Μ 要满足pmiss(Μ)彡pthre的条件，所以可以得到Μ 2 logp_ (Λ&?)。在一般情况下，取pthre = 0.005pthre3 = 0.001就能满足鲁棒性的需要。但是在该算法中由于采样点完全是随机的，对 于每个随机采样点都要计算相关向量，因此比计算该点的可控功率响应还要付出更多计算 代价。
[0022] 4、高分辨率谱估计法
[0023] 基于高分辨率谱估计的声源定位方法借鉴了雷达阵列中的定位技术，但是雷达信 号多为窄带信号，相比之下语音信号频带较宽，因此这种方法的效果不尽理想。
[0024] 声源定位在人机交互技术中有广泛的应用，但是传统方法无法同时达到高鲁棒性 和实时的效果。
[0025] 可见，基于麦克风阵列的声源定位方法已被广泛研究，与其他方法相比，基于可控 功率响应的定位方法具有抗环境噪声和混响干扰的优势，但其计算实时性较差，在实际应 用场景中难以发挥作用。
【发明内容】
[0026] 为了克服上述现有技术的缺点，本发明的目的在于提供一种麦克风阵列声源空间 实时定位方法，其声源定位原理明确，实时性好，实验证明该方法在平面上的定位误差范围 可控制在厘米量级上，性能优于基于现有技术的方法;本发明提出的基于麦克风阵列的声 源空间实时定位方法具有较高的运算速度和鲁棒性，可应用于智能家居和智能机器人等需 要进行实时声源定位的场景。
[0027] 为了实现上述目的，本发明采用的技术方案是：
[0028] -种麦克风阵列声源空间实时定位方法，包括以下步骤：
[0029] 首先，以麦克风阵列作为信号采集和输出设备，通过使用可控功率响应-相位变换 (SRP-PHAT)法初步给出声源空间位置候选点；
[0030] 其次，通过先验知识进行初步的候选点筛选，并使用可控功率响应-相位变换 (SRP-PHAT)法计算候选点的可控功率响应输出；
[0031] 随后，用改进的随机区域收缩（Stochastic Region Contraction,SRC)重新确定 搜索边界，提高可控功率响应-相位变换法的效率；
[0032]最后，计算剩余候选点的可控功率响应，选取最大的位置作为最终的声源估计位置。
[0033] 所述初步给出声源空间位置候选点是指确定全部声源的空间位置候选点，方法如下：
[0034] 假设待搜索空间c的长宽高分别为X，Y，Z，网格间隔为δ，则网格点数 即声源空间位置的候选点为Ν个。
[0035] 所述初步的候选点筛选是指初步降低候选点数目，方法如下：
[0036] 假设麦克风阵列包含m个麦克风，则可得到位于s处的声源到麦克风阵列的到达时 延向量TD0As= ，Tus，T2,3,s，T2,4,s，T3,4,s，......]丁，其中 , j, s = (di,s_dj, s )/V 为位于s处的声源到第i个麦克风和第j个麦克风的时延，V表示空气中声音传播速度，心^表 示位于S处的声源到第i个麦克风的物理距离，dp表示位于S处的声源到第j个麦克风的物 理距离；
[0037] 定义采样点个数差(Sample number Difference，SD)向量为：
[0038] SDs= [sdl,2,s,sdl,3,s,sdl,4,s,sd2,3,s,sd2,4,s,sd3,4,s,......]τ
[0039] 当信号采样频率为fs时，有：
[0040] sdi，j，s = round(fs · Ti，j，s)
[0041] 其中round表示对每个元素向最近的方向取整，如果某些候选点所求出的采样点 差向量SD等同，则只保留其中任意一个候选点，删除其他候选点，避免重复计算；
[0042] 进一步求出每两个候选点si，s2之间的采样点个数差：
[0043] SDsi)s2 = abs(SDsi-SDs2)
[0044] 其中abs表示求取绝对值，并选出所有满足max(SDsi,s2Xthreshold的候选点，只 保留其中任意一个候选点，删除其他候选点，避免重复计算，Threshold定义为：
[0046] 其中，λ表示声音波长。在采样频率fs为16000Hz时，设置threshold=l可以有效降 低候选点的个数。
[0047]在实际应用中，上述TD0A向量和采样点差SD向量的计算，以及候选点的挑选需要 预先做好并存储在一个查找表中，这样在声源定位时，只需根据索引查找即可，不用重复计 算。预先建立查找表也是提速的关键。
[0048] 所述候选点的可控功率响应输出计算方法如下：
[0049] 假定筛选得到的候选点位于空间位置s处，Xd ω )和心（ω )分别代表第i个麦克风 接收信号Xi(t)和第j个麦克风接收信号Xj(t)的傅里叶变换表示声源到第i个麦克风的 传输延迟，τ#表示s处声源到第j个麦克风和第i个麦克风的时延差，根据可控功率响应-相 位变换法的定义，可得声源s的可控功率响应输出P PHAT(s)为：
[0051 ]其中'⑷)= ??)<(?)表示两者之间的互功率谱密度，m为麦克风的总个数，对 所有筛选得到的候选点计算P?AT(s)。
[0052]所述重新确定搜索边界是指用快速随机区域收缩算法搜索全局最大值，方法如 下：
[0053]并非通过随机采样，而是对筛选得到的候选点（且其TD0A向量和采样点差SD向量 都可以在第二步的查找表中直接得到)通过计算得到对应的可控功率响应输出，选取前N个 最大值，在其对应的候选点中随机选取Μ个，随后从这些点出发重新确定搜索边界，然后在 新的区域内重复进行这种选取和搜索，直至满足精度要求。
[0054] 所述Ν的值的选择以及之后的随机采样Μ值的确定与具体的麦克风阵列形状及房 间尺寸等有关，Μ和Ν的选择有很多策略，可选择采用如下方式：
[0055] 方式一，选为定值，Μ也选为定值;例如在四元麦克风阵列条件下，一般来说Ν= 100 时效率和准确度达到最好，性能最佳。
[0056] 方式二，每次都根据上一次的Ν个候选点对应的可控功率响应输出，挑选出其中比 均值大的可控功率响应输出对应的候选点，其总数为Ν'，如果N-Ν' <Ν'，则在这些候选点中 随机挑选N-Ν'个作为候选点；如果Ν-Ν'>Ν'，则保留全部Ν'个候选点，这样可以保证每次区 域收缩后，均值不断增加，当这些候选点的可控功率响应输出计算次数多于某个给定的阈 值时停止选取和搜索。
[0057] 通过这样的方法得到的定位结果能够满足一般的精度要求。如果要更精确的定 位，则可以在最终确定的声源结果附近小范围内使用网格法精确搜索;或者根据查找表，首 先得到若干可控功率响应输出最大值的对应区域，找到这些区域附近网格中的所有候选 点，然后计算这些候选点的可控功率响应，选取最大的位置作为最终的声源估计位置。 [0058]与现有技术相比，本发明的有益效果是：
[0059] (1)本发明提出的通过先验信息进行候选点初步筛除的方法，很大程度上降低了 基于候选点上的可控功率响应输出的计算代价，可以适用于多种场景；
[0060] (2)本发明提出的将候选点位置相关向量预先保存在查找表中的方法，原理简单， 计算代价低，可有效提高实时效果；
[0061] (3)本发明提出的两种基于候选点再次精确搜索声源位置的方法，可进一步提高 声源定位的分辨率和精度，且计算复杂度较低，适用于硬件环境配置较低的设备和场景中。
[0062] 本发明提出基于麦克风阵列的声源空间实时定位方法在对声源信号进行空间定 位时，计算代价优于现有技术水平。本发明的声源定位方法具有应用广泛、响应实时等优 点，适用于智能机器人和智能家居等需要使用声音进行准确定位的场景。
【附图说明】
[0063] 图1是基于麦克风阵列的声源空间实时定位方法总体示意图。
【具体实施方式】
[0064]下面结合附图和实施例详细说明本发明的实施方式。
[0065] 如图1所示，本发明实施例的整个计算过程细节构成如下：
[0066] 1、确定全部声源的空间位置候选点
[0067] 假设待搜索空间c长宽高分别为X，Y，Z，网格间隔为δ，则网格点数N为：
[0069] 即声源空间位置的候选点为Ν个。
[0070] 2、计算到达时延向量和采样点差向量
[0071] 假设麦克风阵列包含m多4个麦克风，则可得到位于s处的声源到麦克风阵列的到 达时延向量TD0As为：
[0072] TD0As= [Ti，2，s,T：i，3，s，T：i，4，s，T2,3，s，T2,4，s，T3,4，s，......]丁 （公式 2)
[0073] 其中11^5=((11,5-山,5)八为位于 8处的声源到第1个麦克风和第」个麦克风的时延， 其中ν表示空气中声音传播速度，cU,s表示位于s处的声源到第i个麦克风的物理距离，山, s表 示位于s处的声源到第j个麦克风的物理距离。
[0074]定义采样点个数差向量SDsS:
[0075] SDs= [sdl,2,s,sdl,3,s,sdl,4,s,sd2,3,s,sd2,4,s,sd3,4,s,......]τ (公式 3)
[0076] 当信号采样频率为fs时，有：
[0077] SDi，j，s = round(fs · Ti，j，s) (公式4)
[0078] 其中round表示对每个元素向最近的方向取整，fs为采样频率。进一步求出每两个 候选点si，s2之间的采样点个数差：
[0079] SDsi,s2 = abs(SDsi-SDs2) (公式 5)
[0080] 其中abs表示求取绝对值。
[0081] 3、删除部分候选点
[0082] 在上述第二步计算中，如果某些候选点用公式4所求出的采样点差向量SD等同，则 只保留其中任意一个候选点，删除其他候选点，避免重复计算。并对剩下的候选点通过公式 5再次计算，选出所有满足max(SD si,S2Xthreshold的候选点，只保留其中任意一个候选点， 删除其他候选点，避免重复计算。此处Threshold定义为：
[0084] 其中，λ表示声音波长。在采样频率fs为16000Hz时，设置threshold=l可以有效降 低候选点的个数。
[0085]另外，在实际应用中，上述TD0A向量和采样点差SD向量的计算，以及候选点的挑选 需要预先做好并存储在一个查找表中，这样在声源定位时，只需根据索引查找即可，不用重 复计算。预先建立查找表也是提速的关键。
[0086] 4、计算候选点的可控功率响应输出
[0087] 假定上述步骤完成后得到的候选点位于空间位置s处，&( ω )和心（ω )分别代表第 i个麦克风接收信号Xi(t)和第j个麦克风接收信号Xj(t)的傅里叶变换表示声源到第i个 麦克风的传输延迟，τ#表示s处声源到第j个麦克风和第i个麦克风的时延差。根据可控功率 响应-相位变换法的定义，可得声源s的可控功率输出P PHAT(s)为：
[0089] 其中1/(叻=夂(叫<(?)表示两者之间的互功率谱密度。对所有第二步得到的候 选点计算P PHAT(S)。
[0090] 5、确定Μ和N，重新确定搜索边界
[0091] 对通过上述步骤计算得到对应的PPHAT(s)，选取前Ν个最大值，在其对应的候选点 中随机选取Μ个，随后从这些点出发重新确定搜索边界。
[0092] Ν的值如何选择以及之后的随机采样Μ和Ν值的确定与具体的麦克风阵列形状及房 间尺寸等有关。Μ和Ν的选择有很多策略，最简单的是固定Μ和Ν为确定值。例如在四元麦克风 阵列条件下，一般来说Ν=100时效率和准确度达到最好，性能最佳。另一种选择Ν的策略是， 每次都根据上一次的Ν个候选点对应的P PHAT(s)挑选出其中比均值大的PPHAT(s)对应的候选 点，其总数为N '，如果N-N '彡N'，则在这些候选点中随机挑选N-N '个作为候选点；如果N-N ' > Ν'，则保留全部Ν'个候选点。这样可以保证每次区域收缩后，均值不断增加。当这些候选点 的PPHAT(s)计算次数多于某个给定的阈值时停止选取和搜索。
[0093] 6、重复进行这种选取和搜索，直至满足精度要求后输出候选点
[0094] 在上述步骤完成后，得到新的搜素区域，在此区域内重复进行这种选取和搜索，直 至满足精度要求，随后输出这些满足要求的候选点。通过这样的方法得到的定位结果能够 满足一般的精度要求。
[0095] 7、进一步用网格法或查表法实现更精细的定位
[0096] 如果要更精确的定位，则可以在最终确定的声源结果附近小范围内使用网格法精 确搜索；或者根据查找表，首先得到若干PPHAT(s)最大值的对应区域，找到这些区域附近网 格中的所有候选点，然后计算这些候选点的可控功率响应，选取最大的位置作为最终的声 源估计位置。
【主权项】
1. 一种麦克风阵列声源空间实时定位方法，其特征在于，包括W下步骤： 首先，W麦克风阵列作为信号采集和输出设备，通过使用可控功率响应-相位变换 (SRP-PHAT)法初步给出声源空间位置候选点； 其次，通过先验知识进行初步的候选点筛选，并使用可控功率响应-相位变换（SRP- PHAT)法计算候选点的可控功率响应输出； 随后，用改进的随机区域收缩（Stochastic Region Contraction,SRC)重新确定捜索 边界，提高可控功率响应-相位变换法的效率； 最后，计算剩余候选点的可控功率响应，选取最大的位置作为最终的声源估计位置。2. 根据权利要求1所述麦克风阵列声源空间实时定位方法，其特征在于，所述初步给出 声源空间位置候选点是指确定全部声源的空间位置候选点，方法如下： 假设待捜索空间C的长宽高分别为Χ，Υ，Ζ，网格间隔为δ，则网格点曼即声源空间位置的候选点为Ν个。3. 根据权利要求1所述麦克风阵列声源空间实时定位方法，其特征在于，所述初步的候 选点筛选是指初步降低候选点数目，方法如下： 假设麦克风阵列包含m个麦克风，则可得到位于S处的声源到麦克风阵列的到达时延向 重 TDOAs 二[Tl,2,s,T!,3,s,Tl,4,s,T2,3,s,T2,4,s,T3,4,s,......]了，其中 Ti,j,s = ( cli, ji-dj, s )/v为位 于S处的声源到第i个麦克风和第j个麦克风的时延，V表示空气中声音传播速度，di,s表示位 于S处的声源到第i个麦克风的物理距离，山,S表示位于S处的声源到第j个麦克风的物理距 离； 定义采样点个数差（Sample number Difference，SD)向量为： SDs= [sdl,2,s,sdl,3,s,sdl,4,s,sd2,3,s,sd2,4,s,sd3,4,s,......]了 当信号采样频率为fs时，有： sdi,j,s = round(fs · Ti'j's) 其中round表示对每个元素向最近的方向取整，如果某些候选点所求出的采样点差向 量SD等同，则只保留其中任意一个候选点，删除其他候选点，避免重复计算； 进一步求出每两个候选点si，s2之间的采样点个数差： SDsi, s2 二 abs (SDs 广 SDs2 ) 其中abs表示求取绝对值，并选出所有满足max(SDsl,s2)《threshold的候选点，只保留 其中任意一个候选点，删除其他候选点，避免重复计算，Threshold定义为：其中，λ表示声音波长。4. 根据权利要求3所述麦克风阵列声源空间实时定位方法，其特征在于，所述采样频率 f S为 16000Hz时，设置 threshold = 1。5. 根据权利要求3所述麦克风阵列声源空间实时定位方法，其特征在于，所述TD0A向量 和采样点差SD向量的计算W及候选点的挑选，预先做好并存储在一个查找表中，在声源定 位时，直接根据索引查找，不用重复计算。6. 根据权利要求3所述麦克风阵列声源空间实时定位方法，其特征在于，所述候选点的 可控功率响应输出计算方法如下： 假定筛选得到的候选点位于空间位置S处，Xi( ω )和、（ω )分别代表第i个麦克风接收 信号xi(t)和第j个麦克风接收信号x^t)的傅里叶变换，τι表示声源到第i个麦克风的传输 延迟，表示S处声源到第j个麦克风和第i个麦克风的时延差，根据可控功率响应-相位变 换法的定义，可得声源S的可控功率响应输出pPHAT(s)为：其中宁,,加')=抑蛛1>'")表示两者之间的互功率谱密度，"1为麦克风的总个数，对所有 筛选得到的候选点计算pPHAT(s)。7. 根据权利要求1所述麦克风阵列声源空间实时定位方法，其特征在于，所述重新确定 捜索边界是指用快速随机区域收缩算法捜索全局最大值，方法如下： 对筛选得到的候选点通过计算得到对应的可控功率响应输出，选取前N个最大值，在其 对应的候选点中随机选取Μ个，随后从运些点出发重新确定捜索边界，然后在新的区域内重 复进行运种选取和捜索，直至满足精度要求。8. 根据权利要求7所述麦克风阵列声源空间实时定位方法，其特征在于，所述Ν的值的 选择W及之后的随机采样Μ值的确定与具体的麦克风阵列形状及房间尺寸有关，Ν的值的选 择采用如下方式： 方式一，选为定值，Μ也选为定值； 方式二，每次都根据上一次的Ν个候选点对应的可控功率响应输出，挑选出其中比均值 大的可控功率响应输出对应的候选点，其总数为Ν'，如果Ν-Ν'《Ν'，则在运些候选点中随机 挑选Ν-Ν'个作为候选点；如果Ν-Ν'〉Ν'，则保留全部Ν'个候选点，当运些候选点的可控功率 响应输出计算次数多于某个给定的阔值时停止选取和捜索。9. 根据权利要求8所述麦克风阵列声源空间实时定位方法，其特征在于，所述方式一 中，在四元麦克风阵列条件下，选择Ν=100时。10. 根据权利要求1所述麦克风阵列声源空间实时定位方法，其特征在于，在最终确定 的声源结果附近小范围内使用网格法精确捜索；或者根据查找表，首先得到若干可控功率 响应输出最大值的对应区域，找到运些区域附近网格中的所有候选点，然后计算运些候选 点的可控功率响应，选取最大的位置作为最终的声源估计位置。
【文档编号】G01S5/18GK106093864SQ201610391351
【公开日】2016年11月9日
【申请日】2016年6月3日 公开号201610391351.9, CN 106093864 A, CN 106093864A, CN 201610391351, CN-A-106093864, CN106093864 A, CN106093864A, CN201610391351, CN201610391351.9
【发明人】杨毅, 孙甲松
【申请人】清华大学