技术领域
[0001] 本
发明涉及音频
信号处理领域，具体涉及一种基于麦克风阵列的回声消除方法与系统。
背景技术
[0002] 实时通信设备进行麦克风拾音时，会受环境噪声、混响、
电路噪声的干扰，同时还受回声的干扰。当远端话者的讲话通过受话器播放时，麦克风在拾取话者声音的同时也拾取受话器信号。受话器信号会直接传达到麦克风，同样也会经过物体反射，衍射，共振等方式传递到麦克风。因此，远端话者在说话同时，会听到自己的声音从对方传回来。
信号处理系统与通信的总延时超过10ms时，话者会
感知到自己的回声；当回声的总延时超过数十毫秒时，话者会觉得回声干扰了自己说话。而实际通信中，信号延时往往会大于100ms以上。在语音/
音频信号处理领域，已经有了很多经典而成熟的回声消除方法。其中LMS、NLMS、GSPAPA等
算法已经广泛应用于话机、
移动电话，保证了优质的通话音质，使用户免受回声的干扰。
[0003] 然而，对于免提通话以及视频会议系统等开放式通信场景中，特别是对于麦克风与喇叭安置在同一机器上的免提式通信设备，麦克风拾取的声信号中，回声
能量经常会多于话者能量及SER(Signal-to-Acoustic Echoes-Ratio)小于0。而SER越小，回声消除算法需要抑制的回声等级就越大。为了提升回声的抑制等级，回声消除算法需要提高自适应
滤波器的精确度。而过度提高精确度导致自适应速率变慢。考虑到各方面因数，传统回声消除算法在上述环境中无法干净的消除拾取信号中的回声信号。
[0004] 在不得已的情况下，多数免提设备限制了扬声器输出保证的最低限度的SER。还有一些设备在扬声器播放大功率声音时动态降低输入灵敏度，保留了扬声器音量，保证了单工通信
质量。然而面对多人会议环境中，不免会出现远端与近端话者同时讲话的情况。此时，远端话者会接收到忽高忽低的拾音信号，严重影了全双工通信的自然通话质量。
[0005] 总而言之，如何能够在免提环境下提供更有效的回声消除方法，保证自然地全双工通信成为了本领域技术人员迫切解决的问题。
发明内容
[0006] 本发明的目的在于提供一种基于麦克风阵列的回声消除方法与系统，其可以大幅提高回声抑制等级，有效消除回声，保证了全双工通信的自然通话质量。
[0007] 为实现上述目的，本发明采用的技术方案是：
[0008] 一种基于麦克风阵列的回声消除方法，其具体包括以下步骤：
[0009] 步骤1、对多路输入的时域音频信号xM进行分
帧，短时时频变换得到频域
输入信号Xn；对参考信号ref进行分帧，短时时频变换得到频域参考信号REFn；
[0010] 步骤2、根据频域输入信号Xn，通过
波速形成
迭代算法更新自适应波速形成增益WABF；
[0011] 步骤3、通过步骤1得到的频域输入信号Xn与步骤2得到的波速形成增益WABF，计算波速形成后的频域信号YABF,
[0012] YABF[k]＝WABF[k]*Xn[k] (1)
[0013] 其中，k＝0,…,N/2-1；N＝2R，R为帧长；
[0014] 步骤4、根据步骤3得到的波速形成后的频域信号YABF，与步骤1得到的频域参考信号REFn，通过自适应回声消除得到频域回声消除残留信号EAEC；
[0015] 步骤5、根据步骤4得到的频域回声消除残留信号EAEC，与步骤1得到的频域参考信号REFn，通过幅度自适应回声消除得到回声消除
输出信号YAEC；
[0016] 步骤6、将回声消除输出信号YAEC进行时频逆变换的得到时域输出信号。
[0017] 所述步骤1中的短时时频变换为重
叠加窗短时傅里叶变换，窗函数选择汉宁窗；
[0018] 所述频域输入信号Xn为：
[0019]
[0020] 所述频域参考信号REFn为：
[0021]
[0022] 其中，k＝0,…,N/2-1；N＝2R，R为帧长；m＝0,…,M-1，M为麦克风阵列的麦克风总数；n为帧序号。
[0023] 所述步骤2具体如下：
[0024] 步骤2.1、根据步骤1得到的频域输入信号xM得到频域输入矩阵Xm，[0025]
[0026] 其中，BNABF为波速形成的分
块数目；k＝0,…,N/2-1；N＝2R，R为帧长；m＝0,…,M-1，M为麦克风阵列的麦克风总数；n为帧序号；
[0027] 步骤2.2、根据步骤2.1得到的频域输入矩阵Xm与自适应波速形成增益WABF，计算[0028]
[0029] 其中， 为上一帧回声消除后的残留信号，PA为波速形成的投影矩阵，μABF为波速形成自适应步长；
[0030] 步骤2.3、计算更新后的自适应波速形成增益WABF，
[0031]
[0032] 其中， 为自适应增益， 为约束增益。
[0033] 所述步骤2中的波速形成的投影矩阵PA与约束增益 根据麦克风摆设及拾音方位确定，具体如下：
[0034] 步骤2.2.1、根据麦克风摆设
位置 与拾音方位 计算约束矩阵计算距离dm，[0035]
[0036] 其中，m＝1，2，…，M；x，y，z是空间坐标；
[0037] 步骤2.2.2、根据频点
频率f计算各频率
波长λ，
[0038] λ[k]＝c/f[k] (25)
[0039] 其中，c为环境中声音传播速度；
[0040] 步骤2.2.3、根据麦克风距离dm与波长λ计算各个频率的转向矩阵a，[0041]
[0042] 其中，d0为参考位置的坐标；
[0043] 步骤2.2.4、根据转向矩阵a计算约束矩阵C，
[0044] C[k]＝{a[k]，sp[k]} (27)
[0045] 其中，sp为喇叭方位的转向矩阵；
[0046] 步骤2.2.5、根据约束矩阵C计算各个频率的约束增益
[0047]
[0048] 其中，
[0049] 步骤2.2.6、根据约束矩阵C计算各个频率的投影矩阵PA，
[0050] PA[k]＝I-C[k]/(C*[k]C[k])C*[k] (29)
[0051] 其中，I为单位矩阵。
[0052] 所述步骤4具体为：
[0053] 步骤4.1、根据频域参考信号REFn计算回声消除的参考矩阵XFm，
[0054]
[0055] BNAEC为回声消除的分块数目，k＝0，…，N/2-1；N＝2R，R为帧长；n为帧序号；
[0056] 步骤4.2、根据回声消除的参考矩阵XFm与上一帧迭代得到的回声消除自适应增益计算AAEC，
[0057]
[0058] 其中，BNABF为波速形成的分块数目；
[0059] 步骤4.3、根据步骤3得到的波速形成后的频域信号YABF与步骤4.2得到的AAEC，计算回声消除后的残留信号EAEC，
[0060] EAEC[k]＝YABF[k]-AAEC[k] (9)
[0061] 步骤4.4、根据步骤3得到的波速形成后的频域信号YABF计算回声消除的谱功率PAEC，
[0062]
[0063] 其中，α为平滑系数，α∈[0,1]；
[0064] 步骤4.5、根据步骤4.4得到的谱功率PAEC与步骤4.3得到的回声消除后的残留信号EAEC，计算回声消除的归一化残留信号
[0065]
[0066] 步骤4.6、根据步骤4.5得到的回声消除的归一化残留信号 计算回声消除的迭代因子PP，
[0067]
[0068] 其中，n＝0,…,BNAEC-1,thAEC为归一化
阈值；
[0069] 步骤4.7、根据回声消除的迭代因子PP更新回声消除自适应增益WAEC，[0070]
[0071] 其中，μAEC为回声消除的迭代步长。
[0072] 所述步骤5中，幅频自适应回声消除具体如下：
[0073] 步骤5.1、根据步骤1得到的频域参考信号REFn计算幅度回声消除的参考矩阵XFmAEC_Mag，
[0074]
[0075] 步骤5.2、根据幅度回声消除的参考矩阵XFmAEC_Mag与上一帧迭代得到的幅度回声消除的自适应增益 计算AAEC_Mag，
[0076]
[0077] 步骤5.3、根据AAEC_Mag与步骤5得到的回声消除后的残留信号EAEC得到幅度回声消除的后的残留信号EAEC_Mag，
[0078] EAEC_Mag[k]＝|EAEC[k]|-AAEC_Mag[k] (16)
[0079] 步骤5.4、根据回声消除后的残留信号EAEC计算幅度回声消除的谱功率PAEC_Mag，[0080]
[0081] 其中，α为平滑系数，α∈[0，1]；
[0082] 步骤5.5、根据谱功率PAEC_Mag与幅度回声消除后的残留信号EAEC_Mag计算归一化残留信号
[0083]
[0084] 步骤5.6、根据上述归一化残留信号 计算幅度回声消除的迭代因子PPMag，[0085]
[0086] 其中，n＝0，…，BNAEC-1；
[0087] 步骤5.7、根据幅度回声消除的迭代因子PPMag更新幅度回声消除自适应增益WAEC_Mag，
[0088]
[0089] 其中，μAEC_Mag为迭代步长；
[0090] 步骤5.8、根据步骤5.4得到的幅度回声消除后的残留信号EAEC_Mag与步骤4得到的幅度回声消除后的残留信号EAEC，计算幅度回声消除非线性增益
[0091]
[0092] 其中，MINGAIN为非线性最小增益；
[0093] 步骤5.9、根据幅度回声消除非线性增益 与步骤4得到的回声消除后的残留信号EAEC得到回声消除输出信号YAEC，
[0094]
[0095] 其中， 为幅度回声消除非线性增益，EAEC为回声消除后的残留信号。
[0096] 所述步骤6中的幅频逆变换为短时傅里叶逆变换，得到的时域输出信号为：
[0097]
[0098] 其中，r＝0,…,N-1，N＝2R，n为帧序号，时域输出信号y初始值为0。
[0099] 相应地，本发明还提供一种基于麦克风阵列的回声消除系统，其包括短时时频变换模块、话者方位检测模块、波速形成重置模块、波速形成自适应模块、波速形成模块、线性回声消除模块、非线性回声消除模块和短时时频逆变换模块；
[0100] 所述短时时频变换模块，用于将输入的时域音频信号进行分帧并变换到频域，并进行缓冲处理；
[0101] 所述话者方位检测模块，用于输出实时检测的话者方位坐标；
[0102] 所述波速形成重置模块，用户设定的指向或检测到的用户方位变更时重新计算约束增益与投影矩阵；
[0103] 所述波速形成自适应模块，根据回声与干扰声最小化算法，自适应变化波速形成增益；
[0104] 所述波速形成模块，用于对频域输入信号应用波速形成增益，消除回声的直达部分；
[0105] 所述线性回声消除模块，根据参考信号自适应学习，对波速形成后的残留回声进行抑制，消除回声的反射部分；
[0106] 所述非线性回声消除模块，利用参考信号与残留回声信号的长时幅度相关度进行增益抑制，进一步抑制回声能量；
[0107] 所述短时时频逆变换模块，用于将回声消除信号进行时频逆变换，以获取其频域信号。
[0108] 本发明通过波速形成、线性回声消除与非线性回声消除三重结构抑制了拾取信号中的回声能量，其中，波束形成模块稳定抑制喇叭方位的直达声信号，不受环境变化影响；线性回声消除模块对波速形成后的残留回声进行抑制，消除回声的反射部分；而非线性回声消除模块利用参考信号与残留回声信号的长时幅度相关度进行增益抑制，进一步抑制回声能量。本发明在大幅提升对回声的抑制效果的同时保证了拾取信号的清晰、自然，与
现有技术相比有更强的性能以及更好的
稳定性。
[0109] 本发明在回声消除之前首先进行波速形成，确保了波速形成算法不受回声消除自适应变化带来的干扰，同时也避免了采用多路回声消除时产生的庞大计算量。
附图说明
[0112] 图3本发明具体实施例预设拾音方向示意图；
[0113] 图4为本发明具体实施例样本信号时幅图；
[0114] 图5为本发明具体实施例样本信号时频图；
[0115] 图6为传统回声消除方式的处理后的时幅图；
[0116] 图7为传统回声消除方式的处理后的时频图；
[0117] 图8为本发明具体实施例处理后的时幅图；
[0118] 图9为本发明具体实施例处理后的时频图。
具体实施方式
[0119] 如图1所示，本发明揭示了一种基于麦克风阵列的回声消除方法，其具体包括以下步骤：
[0120] 步骤1、对多路输入的时域音频信号xM进行分帧，短时时频域变换得到频域输入信号Xn；对参考信号ref进行分帧，短时时频变换得到频域参考信号REFn。
[0121] 其中，短时时频变换为重叠
加窗(overlap-add method)短时傅里叶变换(STFT)，窗函数选择汉宁窗，得到的频域输入信号Xn如式(2)所示：
[0122]
[0123] 其中，k＝0,…,N/2-1；N＝2R，R为帧长；m＝0,…,M-1，M为麦克风阵列的麦克风总数；n为帧序号。后续步骤中k均表示子带序号，m均表示麦克风序号，不再进行重复说明。
[0124] 而进行短时时频变换得到的频域参考信号REFn如式(3)所示：
[0125]
[0126] 步骤2、根据频域输入信号Xn，通过波速形成迭代算法更新自适应波速形成增益WABF；
[0127] 而波速形成迭代算法包括：
[0128] 步骤2.1、根据步骤1得到的频域输入信号xM得到频域输入矩阵Xm，如式(4)所示：
[0129]
[0130] 其中，BNABF为波速形成的分块数目。
[0131] 步骤2.2、根据步骤2.1得到的频域输入矩阵Xm，与自适应波速形成自增益WABF得到如式(5)所示：
[0132]
[0133] 其中， 为上一帧回声消除后的残留信号，PA为波速形成的投影矩阵，μABF为波速形成自适应步长。
[0134] 步骤2.3、根据步骤2.2得到的 与约束增益 相加得到更新后的自适应波速形成增益WABF，如式(6)所示：
[0135]
[0136] 步骤3、通过步骤1得到的频域输入信号Xn与步骤2得到的波速形成增益WABF，计算波速形成后的频域信号YABF,如式(1)所示；
[0137] YABF[k]＝WABF[k]*Xn[k] (1)
[0138] 步骤4、根据步骤3得到的波速形成后的频域信号YABF，与步骤1得到的频域参考信号REFn，通过自适应回声消除得到频域回声消除残留信号EAEC。
[0139] 步骤4.1、根据频域参考信号REFn计算回声消除的参考矩阵XFm，如式(7)所示：
[0140]
[0141] BNAEC为回声消除的分块数目。
[0142] 步骤4.2、根据回声消除的参考矩阵XFm与上一帧迭代得到的回声消除自适应增益计算回声信号AAEC，如式(8)所示：
[0143]
[0144] 步骤4.3、根据步骤3得到的波速形成后的频域信号YABF与步骤4.2得到的AAEC，计算回声消除后的残留信号EAEC，如式(9)所示：
[0145] EAEC[k]＝YABF[k]-AAEC[k] (9)
[0146] 步骤4.4、根据步骤3得到的波速形成后的频域信号YABF计算回声消除的谱功率PAEC，如式(10)所示：
[0147]
[0148] 其中，α为平滑系数，α∈[0,1]。
[0149] 步骤4.5、根据步骤4.4得到的谱功率PAEC与步骤4.3得到的回声消除后的残留信号EAEC，计算回声消除的归一化残留信号 如式(11)所示：
[0150]
[0151] 步骤4.6、根据步骤4.5得到的回声消除的归一化残留信号 计算回声消除的迭代因子PP，如式(12)所示：
[0152]
[0153] 其中，n＝0,…,BNAEC-1,thAEC为归一化阈值。
[0154] 步骤4.7、根据回声消除的迭代因子PP更新回声消除自适应增益WAEC，如式(13)所示：
[0155] WAEC[k][n]＝WAEC[k][n]+μAECPP[k][n] (13)
[0156] 其中，μAEC为回声消除的迭代步长。
[0157] 步骤5、根据步骤4得到的频域回声消除残留信号EAEC，与步骤1得到的频域参考信号REFn，通过幅度自适应回声消除得到回声消除输出信号YAEC。
[0158] 其中，幅频自适应回声消除具体如下：
[0159] 步骤5.1、根据步骤1得到的频域参考信号REFn计算幅度回声消除的参考矩阵XFmAEC_Mag，如式(14)所示：
[0160]
[0161] 步骤5.2、根据幅度回声消除的参考矩阵XFmAEC_Mag与上一帧迭代得到的幅度回声消除的自适应增益 计算AAEC_Mag，如式(15)所示：
[0162]
[0163] 步骤5.3、根据AAEC_Mag与步骤5得到的回声消除后的残留信号EAEC得到幅度回声消除的后的残留信号EAEC_Mag，如式(16)所示：
[0164] EAEC_Mag[k]＝|EAEC[k]|-AAEC_Mag[k] (16)
[0165] 步骤5.4、根据回声消除后的残留信号EAEC计算幅度回声消除的谱功率PAEC_Mag，如式(17)所示：
[0166]
[0167] 其中，α为平滑系数，α∈[0，1]。
[0168] 步骤5.5、根据谱功率PAEC_Mag与幅度回声消除后的残留信号EAEC_Mag计算归一化残留信号 如式(18)所示：
[0169]
[0170] 步骤5.6、根据上述归一化残留信号 计算幅度回声消除的迭代因子PPMag，如式(19)所示：
[0171]
[0172] 其中，n＝0，…，BNAEC-1。
[0173] 步骤5.7、根据幅度回声消除的迭代因子PPMag更新幅度回声消除自适应增益WAEC_Mag，如式(20)所示：
[0174]
[0175] 其中，μAEC_Mag为迭代步长。
[0176] 步骤5.8、根据步骤5.4得到的幅度回声消除后的残留信号EAEC_Mag与步骤4得到的幅度回声消除后的残留信号EAEC，计算幅度回声消除非线性增益 如式(21)所示：
[0177]
[0178] 其中，MINGAIN为非线性最小增益。
[0179] 步骤5.9、根据幅度回声消除非线性增益 与步骤4得到的回声消除后的残留信号EAEC得到回声消除输出信号YAEC，如式(22)所示：
[0180]
[0181] 步骤6、根据上述回声消除输出信号YAEC进行时频逆变换的得到时域输出信号。
[0182] 其中，幅频逆变换为短时傅里叶逆变换，如式(23)所示：
[0183]
[0184] 其中，r＝0,…,N-1，N＝2R，n为帧序号，时域输出信号y初始值为0。
[0185] 上述步骤2中的波速形成的投影矩阵PA与约束增益 可根据麦克风摆设及拾音方位确定，具体如下：
[0186] 步骤2.2.1、根据麦克风摆设位置 与拾音方位 计算约束矩阵计算距离dm，如式(24)所示：
[0187]
[0188] 其中，m＝1,2,…,M；x,y,z是空间坐标。
[0189] 步骤2.2.2、根据频点频率f计算各频率波长λ，如式(25)所示：
[0190] λ[k]＝c/f[k] (25)
[0191] 其中，c为环境中声音传播速度。
[0192] 步骤2.2.3、根据声延时τ计算各个频率的转向矩阵a，如式(26)所示：
[0193]
[0194] 其中，d0为参考位置的坐标，通常选择麦克风阵列的中心位置。
[0195] 步骤3.2.4、根据转向矩阵a计算约束矩阵C，如式(27)所示：
[0196] C[k]＝{a[k]，sp[k]} (27)
[0197] 其中，sp为喇叭方位的转向矩阵，同样可通过式(26)得出。
[0198] 步骤2.2.5、根据约束矩阵C计算各个频率的约束增益 如式(28)所示：
[0199]
[0200] 其中，
[0201] 步骤2.2.6、根据约束矩阵C计算各个频率的投影矩阵PA，如式(29)所示：
[0202] PA[k]＝I-C[k]/(C*[k]C[k])C*[k] (29)
[0203] 其中，I为单位矩阵。
[0204] 基于上述回声消除方法，本发明还揭示了一种基于麦克风阵列的回声消除系统，其包括短时时频变换模块、话者方位检测模块、波速形成重置模块、波速形成自适应模块、波速形成模块、线性回声消除模块、非线性回声消除模块和短时时频逆变换模块；
[0205] 其中，短时时频变换模块，用于将输入的时域音频信号进行分帧并变换到频域，并进行缓冲处理；
[0206] 话者方位检测模块，用于输出实时检测的话者方位坐标；
[0207] 波速形成重置模块，用户设定的指向或检测到的用户方位变更时重新计算约束增益与投影矩阵；
[0208] 波速形成自适应模块，根据回声与干扰声最小化算法，自适应变化波速形成增益；
[0209] 波速形成模块，用于对频域输入信号应用波速形成增益，消除回声的直达部分；
[0210] 线性回声消除模块，根据参考信号自适应学习，对波速形成后的残留回声进行抑制，消除回声的反射部分；
[0211] 非线性回声消除模块，波速形成与回声消除的迭代误差导致残留部分仍存在回声信号，非线性回声消除模块利用参考信号与残留回声信号的长时幅度相关度进行增益抑制，进一步抑制回声能量。
[0212] 所述短时时频逆变换模块，用于将回声消除信号进行时频逆变换，以获取其频域信号。
[0213] 在上述基于麦克风阵列的回声消除方法与系统的
基础上，如图2所示，本发明的具体实施例采用了6个麦克风1的环形摆设，即M＝6，麦克1风间距d为5.6CM；喇叭2摆设于设备中心位置；
采样率fs＝16KHz，帧长R＝256。在该实施例中，波速形成的分块数目BNABF＝10；波速形成自适应步长μABF＝0.05。
[0214] 如图3所示，该实施例预先规定了12个方向的拾音区，话者方位检测模块会根据这12个区域进行检测，返回话者方位，当话者所属区域发生变化时，触发波速形成重置模块，根据步骤2.2.1至步骤2.2.6重新计算投影矩阵PA与约束增益
[0215] 根据麦克风间距为d，步骤2.2.1中的麦克风摆设位置 如式(30)所示：
[0216]
[0217] 参考位置规定为拾音区中心
角度上距离1M，高度0.2M。因此，步骤2.2.1所述的拾音方位 根据12个拾音区分别为：
[0218]
[0219] 其中，l为0至11的整数。
[0220] 进一步的，步骤2.2.2中，频点频率f，可通
过采样率fs与帧长R得出，如式(32)所示：
[0221]
[0222] 本实施例中喇叭安置在设备中心位置与参考位置d0一致。为了得出步骤2.2.4的喇叭方位的转向矩阵sp，把0代入式(26)：
[0223] sp[m][k]＝e-2πj(0)/λ[k] (33)
[0224] 其中，指数部分有一项始终为0，因此，sp始终为1。
[0225] 本实施例中的参数设置如下：步骤4.2中的回声消除的分块数目BNAEC＝5；步骤4.4与步骤5.4的平滑系数α＝0.1；步骤4.6与步骤5.6所述的归一化阈值thAEC＝1.5e-6；步骤4.7的回声消除的迭代步长μAEC＝0.05；步骤5.7的幅度回声消除的迭代步长μAEC_Mag＝0.2；步骤6.8的为非线性最小增益MINGAIN＝0.01。
[0226] 现结合样本信号进行进一步说明，如图4所示，样本信号采样率为16KHz，长度1分60秒。在该样本信号中，前33秒为单独喇叭发生区域；33秒至1分03秒为单独话者发声区域；
1分03秒之后为喇叭与话者同时发声区域。可以看到，因麦克风非常靠近喇叭，麦克风拾取的回声信号明显大于话者信号；样本信号的时频图如图6所示。
[0227] 图8为样本采用本发明回声消除方法后的信号时幅图，图9为对应时频图。作为对比，图6为通过传统回声消除NLMS方法得到的输出信号时幅图，图7为对应时频图。
[0228] 通过对比可知，传统NLSM回声消除方法有明显的回声消除效果，然而，当回声信号大于话者信号时仍有明显的残留回声。而采用本发明回声消除方法后，明显降低了残留回声能量，大幅改善了拾音音质。
[0229] 再者，通过对比图7与图9的前10秒区间可知，在滤波器快速收敛阶段，本发明的处理效果明显优于传统的回声消除方法。这是因为本发明中波速形成算法优先对回声直达声进行约束，降低了系统对滤波器系数变化的敏感度，提高了整体稳定性。
[0230] 总之，本发明提出的基于麦克风阵列的回声消除方法对比传统的回声消除方法有效降低了残留回声能量，保证了免提设备的全双工通话质量的同时没有产生信号失真或谐波过大等不良效应。
[0231] 以上所述，仅是本发明实施例而已，并非对本发明的技术范围作任何限制，故凡是依据本发明的技术实质对以上实施例所作的任何细微
修改、等同变化与修饰，均仍属于本发明技术方案的范围内。
