CN109151237A - 基于零吸引的仿射投影自适应回声消除方法 - Google Patents

Abstract

本发明公开了一种基于零吸引的仿射投影自适应回声消除方法，其步骤如下：A、远端信号处理，将远端传来的远端信号采样得到为当前时刻n的远端信号离散值x(n)，进而构建当前时刻的自适应滤波器输入向量X(n)和仿射投影矩阵U(n)；B、回声估计，将滤波器输入向量X(n)通过自适应滤波器得到输出值y(n)、也即回声的估计值y(n)；C、回声抵消，近端麦克风拾取出当前时刻n的带回声的近端信号d(n)，将当前时刻n的带回声的近端信号d(n)与当前时刻n的自适应滤波器的输出值y(n)相减得到残差信号e(n)，D、成比例仿射投影矩阵的构建；E、零吸引因子计算；F、权系数向量的更新；G、重复。该方法的辨识能力强，收敛速度快、稳态误差低，回声消除效果好。

Description

基于零吸引的仿射投影自适应回声消除方法

技术领域

本发明属于通信的自适应回声对消技术领域。

背景技术

在以语音为主的通信系统(比如免提电话、电视电话会议系统等)中，语音质量通常会受到以声学回声为主的回声影响。以电视电话会议为例，因为扬声器和麦克风被置于同一空间，本地扬声器发出的远端语音被本地近端麦克风接收并传回远端，导致远端说话者听到自己的声音。由于回声对通信质量的严重影响，研究人员提出了多种方法来消除回声。其中，自适应回声消除技术成本低效果好，得到研究学者们的一致认可，是目前国际上公认的最有前景的回声消除技术，也是回声消除采用的主要技术。

自适应回声消除技术是用自适应滤波器来估计回声，并在近端信号中减去回声的估计值以消除回声。自适应回声消除技术的核心是自适应回声消除算法，因此该研究领域最主要工作是找到性能卓越的自适应回声消除算法。

在目前的自适应回声消除应用中，性能较好的自适应回声消除算法有文献1的方法：“An Efficient Proportionate Affine Projection Algorithm for EchoCancellation”(Constantin Paleologu,Member,IEEE,Silviu Ciochina,and JacobBenesty,IEEE Signal Processing Letters,Vol.17,No.2,pp.165-168,2010.)。该方法通过记忆成比例的方式，即：权系数的更新不仅与当前时刻的权系数成正相关，还与前p-1时刻的权系数成正相关，产生多个时刻的记忆效应，更好的消除了背景噪声的干扰，以较低的计算复杂度，提高了自适应滤波器的收敛速度。然而，对于稀疏度较高、零值的权系数较多的回声消除系统，它对小的权系数没有加速变小的运算，导致其收敛速度仍然较慢。

发明内容

本发明的发明目的就是提供一种基于零吸引的新型记忆成比例仿射投影自适应回声消除方法，该方法对通信系统的声学回声的消除效果好，收敛速度快，稳态误差小。

本发明实现其发明目的所采用的技术方案是，基于零吸引的仿射投影自适应回声消除方法，其步骤如下：

A、远端信号处理

将远端传来的远端信号采样得到当前时刻n的远端信号离散值x(n)，将当前时刻n与之前的L-1时刻的远端信号离散值x(n)，x(n-1),...,x(n-L+1)，构成当前时刻n的自适应滤波器输入向量X(n)，X(n)＝[x(n)，x(n-1),...,x(n-L+1)]^T，其中L＝512，是滤波器抽头数；T代表转置运算；

将当前时刻n与之前的p-1时刻的自适应滤波器输入向量X(n)，X(n-1),...,X(n-p+1))，构成当前时刻n的仿射投影输入矩阵U(n)，U(n)＝[X(n),X(n-1),...,X(n-p+1)]，其中p为仿射投影阶数，p＝2～8；

B、回声估计

将当前时刻n的滤波器输入向量X(n)通过自适应滤波器得到当前时刻n的输出值y(n)，也即当前时刻n的回声的估计值y(n)，y(n)＝W^T(n)X(n)；其中，W(n)为当前时刻n的自适应滤波器的抽头权向量，W(n)＝[w₁(n),w₂(n),...,w_l(n),...,w_L(n)]^T；其长度为L，其初始值为零向量即W(0)＝0，w_l(n)为当前时刻n的自适应滤波器的第l个抽头权系数；

C、回声抵消

近端麦克风拾取出当前时刻n的带回声的近端信号d(n)；将当前时刻n的带回声的近端信号d(n)与当前时刻n的自适应滤波器的输出值y(n)相减，得到当前时刻n的残差信号e(n)，e(n)＝d(n)-y(n)；再将当前时刻n的残差信号e(n)回送给远端；

D、成比例仿射投影矩阵的构建

D1、归一化平滑变化率的计算

将当前时刻n的平滑时间窗的近半段的第l个抽头权系数w_l(n)的绝对值加和，得到当前时刻n的第l个抽头权系数的近半段记忆值Q_l(n)，Q_l(n)＝|w_l(n)|+|w_l(n-1)|+...+|w_l(n-J/2+1)|；

其中J为平滑时间窗的长度，其取值为10-20间的偶数；

将当前时刻n的平滑时间窗远半段的第l个抽头权系数w_l(n)的绝对值加和，得到当前时刻n的第l个抽头权系数的远半段记忆值

进而得到当前时刻n的第l个抽头权系数w_l(n)的归一化平滑变化率V_l(n)

其中，min[·]为求最小值运算；

D2、成比例因子参数的计算

由下式计算出当前时刻n的的第l个抽头权系数w_l(n)的成比例因子参数γ_l(n)

其中，α是成比例控制因子，其取值范围为-1～1，max(·)表示取最大值运算，ξ是正则化参数，其取值为0.001～0.01，||·||₁表示欧几里得一范数运算，为归一化平滑变化率阈值，取值为0.01～0.1；

D3、成比例增益系数的计算

用下式计算出当前时刻n的的第l个抽头权系数w_l(n)的成比例增益系数g_l(n)，

将当前时刻n的所有抽头权系数w_l(n)的成比例增益系数g_l(n)，组成当前时刻n的成比例增益向量G(n)，G(n)＝[g₁(n),g₂(n),...,g_L(n)]^T；

D4、成比例仿射投影矩阵的构建

再由下式得出当前时刻n的成比例仿射投影矩阵P(n)：

P(n)＝[G(n)⊙X(n),P_-1(n-1)]

其中，⊙表示两个向量的对应元素相乘的运算；P_-1(n-1)为前一时刻n-1的成比例仿射投影矩阵P(n-1)的前p-1列；成比例仿射投影矩阵的初始值为零矩阵，即P(0)＝0_L×p

E、零吸引因子计算

由下式计算出当前时刻n的零吸引因子Z(n)，

其中ρ为零吸引参数，取值为10^-8～10^-5，ε为零吸引控制因子，取值为1～20，sgn(·)为取符号运算；

F、滤波器抽头权向量更新

根据当前时刻n的：仿射投影输入矩阵U(n)、自适应滤波器的抽头权向量W(n)、零吸引因子Z(n)和成比例仿射投影矩阵P(n)，计算出下一时刻n+1的自适应滤波器的抽头权向量W(n+1)：

W(n+1)＝W(n)+μP(n)(U^T(n)P(n)+δI_P)^-1e(n)-Z(n)

其中，μ为自适应滤波器的步长参数，其取值范围为0＜μ＜2，δ为防止矩阵求逆计算困难的常数，其取值为0.001～0.01，I_p为p×p的单位矩阵；

G、令n＝n+1,重复A～F的步骤,直至通话结束。

与现有技术相比，本发明的有益效果是：

一、本发明的抽头权向量更新公式是：W(n+1)＝W(n)+μP(n)(U^T(n)P(n)+δI_P)^-1e(n)-Z(n)。即在权向量更新时增加了零吸引因子作为减项。该零吸引因子减项，对小于ε的抽头权系数w_l(n)，该抽头权系数w_l(n)越小，其变化越快，得出的零吸引因子越迅速的接近于公式中的被减项抽头权系数w_l(n)本身，使得该抽头权系数加速向零靠近；对于稀疏度较高、零值的权系数较多的回声消除系统，这种对小的抽头权系数的加速变小运算，能够有效提高其收敛速度。

二、本发明抽头权系数w_l(n)的归一化平滑变化率V_l(n)小于等于阈值时，其对应的成比例增益系数g_l(n)为成比例增益向量G(n)的平均值，不再与该抽头权系数w_l(n)本身有关，从而能够抑制背景噪声的干扰，有效提高收敛速度，降低稳态误差。相反，抽头权系数w_l(n)的归一化平滑变化率V_l(n)大于阈值时，其对应的成比例增益系数g_l(n)与抽头权系数w_l(n)本身有关，从而能够有效利用大抽头权系数含有的有效信息，加快收敛速度，降低稳态误差。

下面结合附图和具体实施方式对本发明进行详细说明。

附图说明

图1是本发明仿真实验的稀疏信道图。

图2是本发明仿真实验中的近端信号和远端信号图。

图3是文献1和本发明仿真实验的归一化稳态失调曲线。

下面结合附图和具体实施方式对本发明进行详细说明。

具体实施方式

实施例

本发明的一种具体实施方式是，基于零吸引的仿射投影自适应回声消除方法，其步骤如下：

A、远端信号处理

将远端传来的远端信号采样得到当前时刻n的远端信号离散值x(n)，将当前时刻n与之前的L-1时刻的远端信号离散值x(n)，x(n-1),...,x(n-L+1)，构成当前时刻n的自适应滤波器输入向量X(n)，X(n)＝[x(n)，x(n-1),...,x(n-L+1)]^T，其中L＝512，是滤波器抽头数；T代表转置运算；

将当前时刻n与之前的p-1时刻的自适应滤波器输入向量X(n)，X(n-1),...,X(n-p+1))，构成当前时刻n的仿射投影输入矩阵U(n)，U(n)＝[X(n),X(n-1),...,X(n-p+1)]，其中p为仿射投影阶数，p＝2～8；

B、回声估计

将当前时刻n的滤波器输入向量X(n)通过自适应滤波器得到当前时刻n的输出值y(n)，也即当前时刻n的回声的估计值y(n)，y(n)＝W^T(n)X(n)；其中，W(n)为当前时刻n的自适应滤波器的抽头权向量，W(n)＝[w₁(n),w₂(n),...,w_l(n),...,w_L(n)]^T；其长度为L，其初始值为零向量即W(0)＝0，w_l(n)为当前时刻n的自适应滤波器的第l个抽头权系数；

C、回声抵消

近端麦克风拾取出当前时刻n的带回声的近端信号d(n)；将当前时刻n的带回声的近端信号d(n)与当前时刻n的自适应滤波器的输出值y(n)相减，得到当前时刻n的残差信号e(n)，e(n)＝d(n)-y(n)；再将当前时刻n的残差信号e(n)回送给远端；

D、成比例仿射投影矩阵的构建

D1、归一化平滑变化率的计算

将当前时刻n的平滑时间窗的近半段的第l个抽头权系数w_l(n)的绝对值加和，得到当前时刻n的第l个抽头权系数的近半段记忆值Q_l(n)，Q_l(n)＝|w_l(n)|+|w_l(n-1)|+...+|w_l(n-J/2+1)|；

其中J为平滑时间窗的长度，其取值为10-20间的偶数；

将当前时刻n的平滑时间窗远半段的第l个抽头权系数w_l(n)的绝对值加和，得到当前时刻n的第l个抽头权系数的远半段记忆值

进而得到当前时刻n的第l个抽头权系数w_l(n)的归一化平滑变化率V_l(n)

其中，min[·]为求最小值运算；

D2、成比例因子参数的计算

由下式计算出当前时刻n的的第l个抽头权系数w_l(n)的成比例因子参数γ_l(n)

其中，α是成比例控制因子，其取值范围为-1～1，max(·)表示取最大值运算，ξ是正则化参数，其取值为0.001～0.01，||·||₁表示欧几里得一范数运算，为归一化平滑变化率阈值，取值为0.01～0.1；

D3、成比例增益系数的计算

用下式计算出当前时刻n的的第l个抽头权系数w_l(n)的成比例增益系数g_l(n)，

将当前时刻n的所有抽头权系数w_l(n)的成比例增益系数g_l(n)，组成当前时刻n的成比例增益向量G(n)，G(n)＝[g₁(n),g₂(n),...,g_L(n)]^T；

D4、成比例仿射投影矩阵的构建

再由下式得出当前时刻n的成比例仿射投影矩阵P(n)：

P(n)＝[G(n)⊙X(n),P_-1(n-1)]

其中，⊙表示两个向量的对应元素相乘的运算；P_-1(n-1)为前一时刻n-1的成比例仿射投影矩阵P(n-1)的前p-1列；成比例仿射投影矩阵的初始值为零矩阵，即P(0)＝0_L×p

E、零吸引因子计算

由下式计算出当前时刻n的零吸引因子Z(n)，

其中ρ为零吸引参数，取值为10^-8～10^-5，ε为零吸引控制因子，取值为1～20，sgn(·)为取符号运算；

F、滤波器抽头权向量更新

根据当前时刻n的：仿射投影输入矩阵U(n)、自适应滤波器的抽头权向量W(n)、零吸引因子Z(n)和成比例仿射投影矩阵P(n)，计算出下一时刻n+1的自适应滤波器的抽头权向量W(n+1)：

W(n+1)＝W(n)+μP(n)(U^T(n)P(n)+δI_P)^-1e(n)-Z(n)

其中，μ为自适应滤波器的步长参数，其取值范围为0＜μ＜2，δ为防止矩阵求逆计算困难的常数，其取值为0.001～0.01，I_p为p×p的单位矩阵；

G、令n＝n+1,重复A～F的步骤,直至通话结束。

仿真实验：

为了验证本发明的有效性，进行了仿真实验，并与现有的文献1算法进行了对比。

仿真实验的远端信号u(n)为语音信号，采样频率为8000Hz，采样点个数40000。回声信道脉冲响应在高2.5m，宽3.75m，长6.25m，温度20℃，湿度50％的安静密闭房间内获得，脉冲响应长度即滤波器抽头数L＝512。实验的背景噪声为高斯白噪声，信噪比为30dB。并且在麦克风接收到的近端信号中加入采集得到的冲击干扰。

按照以上实验条件，用本发明方法与现有的两种方法进行回声消除实验。各种方法的参数具体取值如表1。

表1实验各算法的最优参数近似取值

图1是实验用的安静密闭房间构成的通信系统的稀疏信道图。图2是实验中近端信号图。图3是文献1和本发明的实验的归一化稳态失调曲线。

从图3可以看出：本发明在约6000个采样时刻收敛，稳态误差约在-38.5dB；而文献1则在约9000个采样时刻收敛，稳态误差约在-35dB；本发明比文献1稳态误差减小了-3dB；同时，本发明比文献1收敛速度增加了近一倍。

Claims (1)

1.基于零吸引的仿射投影自适应回声消除方法，其步骤如下：

A、远端信号处理

将远端传来的远端信号采样得到当前时刻n的远端信号离散值x(n)，将当前时刻n与之前的L-1时刻的远端信号离散值x(n)，x(n-1),...,x(n-L+1)，构成当前时刻n的自适应滤波器输入向量X(n)，X(n)＝[x(n)，x(n-1),...,x(n-L+1)]^T，其中L＝512，是滤波器抽头数；T代表转置运算；

将当前时刻n与之前的p-1时刻的自适应滤波器输入向量X(n)，X(n-1),...,X(n-p+1))，构成当前时刻n的仿射投影输入矩阵U(n)，U(n)＝[X(n),X(n-1),...,X(n-p+1)]，其中p为仿射投影阶数，p＝2～8；

B、回声估计

将当前时刻n的滤波器输入向量X(n)通过自适应滤波器得到当前时刻n的输出值y(n)，也即当前时刻n的回声的估计值y(n)，y(n)＝W^T(n)X(n)；其中，W(n)为当前时刻n的自适应滤波器的抽头权向量，W(n)＝[w₁(n),w₂(n),...,w_l(n),...,w_L(n)]^T；其长度为L，其初始值为零向量即W(0)＝0，w_l(n)为当前时刻n的自适应滤波器的第l个抽头权系数；

C、回声抵消

近端麦克风拾取出当前时刻n的带回声的近端信号d(n)；将当前时刻n的带回声的近端信号d(n)与当前时刻n的自适应滤波器的输出值y(n)相减，得到当前时刻n的残差信号e(n)，e(n)＝d(n)-y(n)；再将当前时刻n的残差信号e(n)回送给远端；

D、成比例仿射投影矩阵的构建

D1、归一化平滑变化率的计算

将当前时刻n的平滑时间窗的近半段的第l个抽头权系数w_l(n)的绝对值加和，得到当前时刻n的第l个抽头权系数的近半段记忆值Q_l(n),Q_l(n)＝|w_l(n)|+|w_l(n-1)|+...+|w_l(n-J/2+1)|；

其中J为平滑时间窗的长度，其取值为10-20间的偶数；

将当前时刻n的平滑时间窗远半段的第l个抽头权系数w_l(n)的绝对值加和，得到当前时刻n的第l个抽头权系数的远半段记忆值

进而得到当前时刻n的第l个抽头权系数w_l(n)的归一化平滑变化率V_l(n)

其中，min[·]为求最小值运算；

D2、成比例因子参数的计算

由下式计算出当前时刻n的的第l个抽头权系数w_l(n)的成比例因子参数γ_l(n)

其中，α是成比例控制因子，其取值范围为-1～1，max(·)表示取最大值运算，ξ是正则化参数，其取值为0.001～0.01，||·||₁表示欧几里得一范数运算，为归一化平滑变化率阈值，取值为0.01～0.1；

D3、成比例增益系数的计算

用下式计算出当前时刻n的的第l个抽头权系数w_l(n)的成比例增益系数g_l(n)，

将当前时刻n的所有抽头权系数w_l(n)的成比例增益系数g_l(n)，组成当前时刻n的成比例增益向量G(n)，G(n)＝[g₁(n),g₂(n),...,g_L(n)]^T；

D4、成比例仿射投影矩阵的构建

再由下式得出当前时刻n的成比例仿射投影矩阵P(n)：

P(n)＝[G(n)⊙X(n),P_-1(n-1)]

其中，⊙表示两个向量的对应元素相乘的运算；P_-1(n-1)为前一时刻n-1的成比例仿射投影矩阵P(n-1)的前p-1列；成比例仿射投影矩阵的初始值为零矩阵，即P(0)＝0_L×p

E、零吸引因子计算

由下式计算出当前时刻n的零吸引因子Z(n)，

其中ρ为零吸引参数，取值为10^-8～10^-5，ε为零吸引控制因子，取值为1～20，sgn(·)为取符号运算；

F、滤波器抽头权向量更新

根据当前时刻n的：仿射投影输入矩阵U(n)、自适应滤波器的抽头权向量W(n)、零吸引因子Z(n)和成比例仿射投影矩阵P(n)，计算出下一时刻n+1的自适应滤波器的抽头权向量W(n+1)：

W(n+1)＝W(n)+μP(n)(U^T(n)P(n)+δI_P)^-1e(n)-Z(n)

其中，μ为自适应滤波器的步长参数，其取值范围为0＜μ＜2，δ为防止矩阵求逆计算困难的常数，其取值为0.001～0.01，I_p为p×p的单位矩阵；

G、令n＝n+1,重复A～F的步骤,直至通话结束。