CN109102794A - 基于凸组合的m估计成比例类仿射投影的回声消除方法 - Google Patents
基于凸组合的m估计成比例类仿射投影的回声消除方法 Download PDFInfo
- Publication number
- CN109102794A CN109102794A CN201810816058.1A CN201810816058A CN109102794A CN 109102794 A CN109102794 A CN 109102794A CN 201810816058 A CN201810816058 A CN 201810816058A CN 109102794 A CN109102794 A CN 109102794A
- Authority
- CN
- China
- Prior art keywords
- length
- current time
- filter
- value
- vector
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Withdrawn
Links
- 238000000034 method Methods 0.000 title claims abstract description 25
- PXFBZOLANLWPMH-UHFFFAOYSA-N 16-Epiaffinine Natural products C1C(C2=CC=CC=C2N2)=C2C(=O)CC2C(=CC)CN(C)C1C2CO PXFBZOLANLWPMH-UHFFFAOYSA-N 0.000 title claims abstract description 23
- 230000011664 signaling Effects 0.000 claims abstract description 5
- 230000003044 adaptive effect Effects 0.000 claims description 15
- 239000011159 matrix material Substances 0.000 claims description 9
- 230000000694 effects Effects 0.000 claims description 7
- 238000001914 filtration Methods 0.000 claims description 5
- 238000004364 calculation method Methods 0.000 claims description 3
- 230000017105 transposition Effects 0.000 claims description 3
- 238000010276 construction Methods 0.000 claims 1
- 238000004422 calculation algorithm Methods 0.000 description 33
- 238000002474 experimental method Methods 0.000 description 15
- 238000004891 communication Methods 0.000 description 8
- 230000001276 controlling effect Effects 0.000 description 6
- 238000005516 engineering process Methods 0.000 description 6
- 230000002596 correlated effect Effects 0.000 description 2
- 238000013507 mapping Methods 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 230000003068 static effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G10—MUSICAL INSTRUMENTS; ACOUSTICS
- G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
- G10L21/00—Speech or voice signal processing techniques to produce another audible or non-audible signal, e.g. visual or tactile, in order to modify its quality or its intelligibility
- G10L21/02—Speech enhancement, e.g. noise reduction or echo cancellation
- G10L21/0208—Noise filtering
-
- G—PHYSICS
- G10—MUSICAL INSTRUMENTS; ACOUSTICS
- G10K—SOUND-PRODUCING DEVICES; METHODS OR DEVICES FOR PROTECTING AGAINST, OR FOR DAMPING, NOISE OR OTHER ACOUSTIC WAVES IN GENERAL; ACOUSTICS NOT OTHERWISE PROVIDED FOR
- G10K11/00—Methods or devices for transmitting, conducting or directing sound in general; Methods or devices for protecting against, or for damping, noise or other acoustic waves in general
- G10K11/16—Methods or devices for protecting against, or for damping, noise or other acoustic waves in general
- G10K11/175—Methods or devices for protecting against, or for damping, noise or other acoustic waves in general using interference effects; Masking sound
- G10K11/178—Methods or devices for protecting against, or for damping, noise or other acoustic waves in general using interference effects; Masking sound by electro-acoustically regenerating the original acoustic waves in anti-phase
- G10K11/1785—Methods, e.g. algorithms; Devices
- G10K11/17853—Methods, e.g. algorithms; Devices of the filter
- G10K11/17854—Methods, e.g. algorithms; Devices of the filter the filter being an adaptive filter
-
- G—PHYSICS
- G10—MUSICAL INSTRUMENTS; ACOUSTICS
- G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
- G10L19/00—Speech or audio signals analysis-synthesis techniques for redundancy reduction, e.g. in vocoders; Coding or decoding of speech or audio signals, using source filter models or psychoacoustic analysis
- G10L19/04—Speech or audio signals analysis-synthesis techniques for redundancy reduction, e.g. in vocoders; Coding or decoding of speech or audio signals, using source filter models or psychoacoustic analysis using predictive techniques
- G10L19/26—Pre-filtering or post-filtering
-
- G—PHYSICS
- G10—MUSICAL INSTRUMENTS; ACOUSTICS
- G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
- G10L21/00—Speech or voice signal processing techniques to produce another audible or non-audible signal, e.g. visual or tactile, in order to modify its quality or its intelligibility
- G10L21/02—Speech enhancement, e.g. noise reduction or echo cancellation
- G10L21/0208—Noise filtering
- G10L2021/02082—Noise filtering the noise being echo, reverberation of the speech
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Acoustics & Sound (AREA)
- Multimedia (AREA)
- Computational Linguistics (AREA)
- Signal Processing (AREA)
- Health & Medical Sciences (AREA)
- Audiology, Speech & Language Pathology (AREA)
- Human Computer Interaction (AREA)
- Quality & Reliability (AREA)
- Filters That Use Time-Delay Elements (AREA)
- Cable Transmission Systems, Equalization Of Radio And Reduction Of Echo (AREA)
Abstract
一种基于凸组合的M估计成比例类仿射投影回声消除方法，其步骤如下：A、远端信号采样；B、凸组合；将当前时刻n的大步长滤波器值y1(n)和小步长滤波器值y2(n)通过当前时刻n的大步长滤波器权重λ(n)进行凸组合得到当前时刻n的组合滤波器值y(n)。C、回声信号消除，将近端麦克风拾取的带回声的近端信号d(n)与自适应滤波器的输出值y(n)相减后再回送给远端，回送信号为残差信号e(n)，e(n)＝d(n)‑y(n)；D、滤波器抽头权向量更新，使用基于凸组合的M估计成比例类仿射投影的方法，计算自适应滤波器抽头权向量的更新；E、滤波器权重的更新；F、滤波器权重的限定；G、令n＝n+1,重复A、B、C、D、E、F、G的步骤,直至迭代结束。
Description
技术领域
本发明涉及一种用于通信系统的自适应回声消除方法。
背景技术
在以语音为主的通信系统(比如免提电话、电视电话会议系统等)中，语音质量通常会受到以声学回声为主的回声影响，严重影响通话质量。回声，即声音或信号经过延时或形变被反射回信号源的一种现象，在语音通信、数据通信、卫星通信、免提电话、电话会议系统等通信系统中，都不同程度的存在回声现象。以电视电话会议为例，因为扬声器和麦克风被置于同一空间，本地扬声器发出的远端语音被本地近端麦克风接收并传回远端，导致远端说话者听到自己的声音。因此，必须采取有效的措施来消除回声信号，减轻其影响，提高语音通话质量。目前，在众多回声消除方法中，自适应回声消除技术具有逐步调节性能，应用成本低，收敛速度快，回声残差小，是目前国际上公认的最有前景的回声消除技术，也是回声消除目前采用的主流技术。
所谓自适应回声消除技术，是利用自适应滤波器的一种回声消除方法。自适应滤波器是能够根据输入信号自主调整性能进行数字信号处理的数字滤波器，作为对比，传统的一些非自适应滤波器(如FIR滤波器，IIR滤波器)的参数是静态的预先设定不可调整的。对于一些未知的环境来说，自适应滤波器能够在工作过程中逐步学习出所需的统计特性，并以此为依据自动调整滤波器系数，以达到最佳滤波的效果；一旦输入信号的统计特性发生变化，它又可以跟踪这种变化，自动调整滤波器参数，使滤波性能重新达到最佳。自适应回声消除技术的本质是通过自适应滤波器来估计回声，并在近端信号中减去回声的估计值以消除回声。自适应回声消除技术的核心是自适应回声消除算法。
目前的自适应回声消除方法有参考文献1“类仿射投影算法的改进研究”(周彤，33-41)提出的IPAPLM算法。该方法使输入向量乘以跟输入向量正相关的成比例矩阵，使得输入信号分配到不同的自适应滤波器权值更新参数与输入向量正相关，从而加快了算法辨识稀疏系统时的收敛速度；通过M估计的思想，设置一个阀值参数，当误差小于阀值参数时，权值向量正常更新，当误差大于阀值参数时，阀值参数代替误差进行权值向量更新，使算法具有良好的抗冲激噪声能力；但该方法只考虑了滤波器权值更新参数的改进，但其自适应滤波器的更新步长是固定的，采用的步长大，收敛速度快，稳态误差也大；采用的步长小，则稳态误差小，收敛速度慢；无法兼顾收敛速度和稳态误差的问题。
发明内容
本发明的发明目的就是提供一种基于凸组合的M估计成比例类仿射投影自适应回声消除方法，该方法对通信系统的收敛速度快，稳态误差小，回声消除效果好。
本发明实现其发明目的所采用的技术方案是，一种基于凸组合的M估计成比例类仿射投影自适应回声消除方法，其步骤如下：
A、远端信号采样
A1、将当前时刻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代表转置运算；
A2、将当前时刻n的滤波器输入向量X(n)通过凸组合自适应回声消除滤波器滤波后分别得到当前时刻n的大步长滤波器值y1(n)，y1(n)＝W1 T(n)X(n)和当前时刻n的小步长滤波器值y2(n)，y2(n)＝W2 T(n)X(n)；其中，W1(n)＝[w1,0(n),w1,1(n),...w1,l(n)...,w1,L-1(n)]T，为当前时刻n的凸组合自适应回声消除滤波器中的大步长滤波器的抽头权向量，其初始值为零向量；w1,l(n)为大步长滤波器的第l个抽头权系数，也是大步长滤波器向量中的第l个元素；W2(n)＝[w2,0(n),w2,1(n),...w2,l(n)...,w2,L-1(n)]T为当前时刻n的小步长滤波器的抽头权向量，其初始值为零向量；w2,l(n)为小步长滤波器的第l个抽头权系数，也是小步长滤波器向量中的第l个元素；
B、凸组合
将当前时刻n的大步长滤波器值y1(n)和小步长滤波器值y2(n)通过当前时刻n的大步长滤波器权重λ(n)进行凸组合得到当前时刻n的组合滤波器值y(n)，y(n)＝λ(n)y1(n)+(1-λ(n))y2(n)。
其中，大步长滤波器权重λ(n)的计算公式为取值范围为0～1；α(n)为当前时刻n的混合参数，其初始值为0；
C、回声信号消除
将近端麦克风拾取的当前时刻n的带回声的近端信号d(n)，与自适应滤波器的当前时刻n的输出值y(n)相减后得到当前时刻n的残差信号e(n)，e(n)＝d(n)-y(n)；再将残差信号e(n)回送给远端；
同时，将当前时刻n的近端信号d(n)与当前时刻n的大步长滤波器值y1(n)和小步长滤波器值y2(n)相减，分别得到当前时刻n的大步长残差信号e1(n)和小步长残差信号e2(n)：
e1(n)＝d(n)-y1(n)，e2(n)＝d(n)-y2(n)；
D、滤波器抽头权向量更新
D1、M估计函数计算
将当前时刻n到时刻n-Nw+1之间的大步长残差信号平方e1 2(n)，e1 2(n-1)，…,e1 2(n-Nw+1)，构成当前时刻n估计窗内的大步长残差信号平方序列Ae1(n)：
Ae1(n)＝[e1 2(n),e1 2(n-1),...,e1 2(n-Nw+1)]
其中，Nw为估计窗的长度，其取值范围为5～15；
再由当前时刻n估计窗内的大步长残差信号平方序列Ae1(n)，计算出当前时刻n的大步长残差信号的方差
其中，λ为遗忘因子，其取值范围为0.800～0.999,C为常数，C＝1.483(1+5/(Nw-1))，med(·)表示取中间值的运算；
然后由当前时刻n的大步长残差信号的方差得出当前时刻n的大步长M估计阀值参数ξ1(n)，
由下式计算出滤波器当前时刻n的大步长抽头权向量更新的M估计函数值
其中，的初始值即为0；
由当前时刻n到时刻n-Nw+1之间的小步长残差信号平方值e2 2(n),e2 2(n-1),…,e2 2(n-Nw+1)，构成当前时刻n估计窗内的小步长残差信号平方序列Ae2(n)：
Ae2(n)＝[e2 2(n),e2 2(n-1),...,e2 2(n-Nw+1)]
再由当前时刻n估计窗内的小步长残差信号平方序列Ae2(n)，计算出当前时刻n的小步长残差信号的方差
然后由当前时刻n的小步长残差信号的方差得出当前时刻n的小步长M估计阀值参数ξ2(n)，
由下式计算出滤波器当前时刻n的小步长抽头权向量更新的M估计函数值
其中，的初始值即为0；
D2、类仿射投影
将当前时刻n到时刻n-P+1的自适应滤波器输入向量X(n),X(n-1),..,X(n-P+1)，构成当前时刻n的类仿射投影输入矩阵
将当前时刻n到n-P+1时刻的大步长抽头权向量更新的M估计函数值构成大步长抽头权向量更新的类仿射投影M估计向量ψ1(n)，
将当前时刻n到n-P+1时刻的小步长抽头权向量更新的M估计函数值构成小步长抽头权向量更新的类仿射投影M估计向量ψ2(n)，
其中，P为投影阶数，其取值为2、4、8；
D3、成比例控制因子计算
由当前时刻n大步长滤波器的第l个抽头权系数w1,l(n)，得出当前时刻n第l个大步长成比例控制因子g1,l(n)，
其中，κ为成比例调整参数，其取值范围为：-1≤κ＜1；ε为成比例限制参数，取值为0.01-0.001，作用是防止公式分母为0；
再计算出，当前时刻n的大步长成比例矩阵G1(n)，
G1(n)＝diag[g1,0(n),g1,1(n),...,g1,l(n),...,g1,L-1(n)]，其中diag[·]表示构造对角矩阵；
由当前时刻n小步长滤波器的第l个抽头权系数w2,l(n)，得出当前时刻n第l个小步长成比例控制因子g2,l(n)，
再计算出，当前时刻n的小步长成比例矩阵G2(n)，
G2(n)＝diag[g2,0(n),g2,1(n),...,g2,l(n),...,g2,L-1(n)]；
D3、组合滤波器抽头权向量的更新
如果当前时刻n的混合参数α(n)大于等于阈值σ，所述的阈值σ的取值为4；
则下一时刻n+1的大步长滤波器的抽头权向量W1(n+1)和小步长滤波器的抽头权向量W2(n+1)分别为：
如果当前时刻n的混合参数α(n)小于阈值σ，则n+1时刻的大步长滤波器的抽头权向量W1(n+1)和小步长滤波器的抽头权向量W2(n+1)分别为：
其中，μ1为大步长滤波器的步长，其取值为0.20～0.85；μ2为小步长滤波器的步长，其取值为0.020～0.085；δ为滤波器抽头权向量的正则化参数，其取值为0.001～0.01；
E、滤波器的权重更新
更新得到下一时刻n+1的混合参数值α(n+1)：
其中，μα是一个常数，取值为10；
再更新得到下一时刻n+1的大步长滤波器权重λ(n+1)：
F、滤波器的权重限定
如果下一时刻n+1的混合参数α(n+1)小于阈值σ的负值，则令α(n+1)＝-σ，λ(n+1)＝0；如果下一时刻n+1的混合参数α(n+1)大于等于阈值σ，则令α(n+1)＝σ，λ(n+1)＝1；
G、重复
令n＝n+1，重复步骤A-F的操作，直至通话结束。
与现有技术相比，本发明的有益效果是：
一、本发明将一个大步长滤波器和一个小步长滤波器结合形成凸组合滤波器，能同时保证滤波器具有较快的收敛速度和较低的稳态误差(残差)，解决固定步长滤波器的不能同时兼顾收敛速度和稳态误差之间的矛盾；
二、两个滤波器的权重更新，基于对应滤波器的误差及输出值，即滤波器的权重与其误差及输出值成正相关，从而进一步降低了算法的稳态误差，调高了回声消除效果。
下面结合附图和具体实施方式对本发明进行详细的说明。
附图说明
图1是本发明仿真实验的声学回声信道图。
图2是本发明仿真实验的网络回声信道图。
图3是本发明仿真实验中的语音信号。
图4是大步长类仿射投影算法的改进方法(IPAPLMμ＝1)、小步长类仿射投影算法的改进方法(IPAPLMμ＝0.2)和本发明(CIPAPLMWT)方法在声学回声信道中的仿真实验归一化稳态失调(MSD)曲线。
图5是大步长类仿射投影算法的改进方法(IPAPLMμ＝1)、小步长类仿射投影算法的改进方法(IPAPLMμ＝0.2)和本发明(CIPAPLMWT)方法在网络回声信道中的仿真实验归一化稳态失调(MSD)曲线。
具体实施方式
实施例
本发明的一种具体实施方式是，一种基于凸组合的M估计成比例类仿射投影回声消除方法，其步骤如下：
A、远端信号采样
A1、将当前时刻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代表转置运算；
A2、将当前时刻n的滤波器输入向量X(n)通过凸组合自适应回声消除滤波器滤波后分别得到当前时刻n的大步长滤波器值y1(n)，y1(n)＝W1 T(n)X(n)和当前时刻n的小步长滤波器值y2(n)，y2(n)＝W2 T(n)X(n)；其中，W1(n)＝[w1,0(n),w1,1(n),...w1,l(n)...,w1,L-1(n)]T，为当前时刻n的凸组合自适应回声消除滤波器中的大步长滤波器的抽头权向量，其初始值为零向量；w1,l(n)为大步长滤波器的第l个抽头权系数，也是大步长滤波器向量中的第l个元素；W2(n)＝[w2,0(n),w2,1(n),...w2,l(n)...,w2,L-1(n)]T为当前时刻n的小步长滤波器的抽头权向量，其初始值为零向量；w2,l(n)为小步长滤波器的第l个抽头权系数，也是小步长滤波器向量中的第l个元素；
B、凸组合
将当前时刻n的大步长滤波器值y1(n)和小步长滤波器值y2(n)通过当前时刻n的大步长滤波器权重λ(n)进行凸组合得到当前时刻n的组合滤波器值y(n)，y(n)＝λ(n)y1(n)+(1-λ(n))y2(n)。
其中，大步长滤波器权重λ(n)的计算公式为取值范围为0～1；α(n)为当前时刻n的混合参数，其初始值为0；
C、回声信号消除
将近端麦克风拾取的当前时刻n的带回声的近端信号d(n)，与自适应滤波器的当前时刻n的输出值y(n)相减后得到当前时刻n的残差信号e(n)，e(n)＝d(n)-y(n)；再将残差信号e(n)回送给远端；
同时，将当前时刻n的近端信号d(n)与当前时刻n的大步长滤波器值y1(n)和小步长滤波器值y2(n)相减，分别得到当前时刻n的大步长残差信号e1(n)和小步长残差信号e2(n)：
e1(n)＝d(n)-y1(n)，e2(n)＝d(n)-y2(n)；
D、滤波器抽头权向量更新
D1、M估计函数计算
将当前时刻n到时刻n-Nw+1之间的大步长残差信号平方e1 2(n)，e1 2(n-1)，…,e1 2(n-Nw+1)，构成当前时刻n估计窗内的大步长残差信号平方序列Ae1(n)：
Ae1(n)＝[e1 2(n),e1 2(n-1),...,e1 2(n-Nw+1)]
其中，Nw为估计窗的长度，其取值范围为5～15；
再由当前时刻n估计窗内的大步长残差信号平方序列Ae1(n)，计算出当前时刻n的大步长残差信号的方差
其中，λ为遗忘因子，其取值范围为0.800～0.999,C为常数，C＝1.483(1+5/(Nw-1))，med(·)表示取中间值的运算；
然后由当前时刻n的大步长残差信号的方差得出当前时刻n的大步长M估计阀值参数ξ1(n)，
由下式计算出滤波器当前时刻n的大步长抽头权向量更新的M估计函数值
其中，的初始值即为0；
由当前时刻n到时刻n-Nw+1之间的小步长残差信号平方值e2 2(n),e2 2(n-1),…,e2 2(n-Nw+1)，构成当前时刻n估计窗内的小步长残差信号平方序列Ae2(n)：
Ae2(n)＝[e2 2(n),e2 2(n-1),...,e2 2(n-Nw+1)]
再由当前时刻n估计窗内的小步长残差信号平方序列Ae2(n)，计算出当前时刻n的小步长残差信号的方差
然后由当前时刻n的小步长残差信号的方差得出当前时刻n的小步长M估计阀值参数ξ2(n)，
由下式计算出滤波器当前时刻n的小步长抽头权向量更新的M估计函数值
其中，的初始值即为0；
D2、类仿射投影
将当前时刻n到时刻n-P+1的自适应滤波器输入向量X(n),X(n-1),..,X(n-P+1)，构成当前时刻n的类仿射投影输入矩阵
将当前时刻n到n-P+1时刻的大步长抽头权向量更新的M估计函数值构成大步长抽头权向量更新的类仿射投影M估计向量ψ1(n)，
将当前时刻n到n-P+1时刻的小步长抽头权向量更新的M估计函数值构成小步长抽头权向量更新的类仿射投影M估计向量ψ2(n)，
其中，P为投影阶数，其取值为2、4、8；
D3、成比例控制因子计算
由当前时刻n大步长滤波器的第l个抽头权系数w1,l(n)，得出当前时刻n第l个大步长成比例控制因子g1,l(n)，
其中，κ为成比例调整参数，其取值范围为：-1≤κ＜1；ε为成比例限制参数，取值为0.01-0.001，作用是防止公式分母为0；
再计算出，当前时刻n的大步长成比例矩阵G1(n)，
G1(n)＝diag[g1,0(n),g1,1(n),...,g1,l(n),...,g1,L-1(n)]，其中diag[·]表示构造对角矩阵；
由当前时刻n小步长滤波器的第l个抽头权系数w2,l(n)，得出当前时刻n第l个小步长成比例控制因子g2,l(n)，
再计算出，当前时刻n的小步长成比例矩阵G2(n)，
G2(n)＝diag[g2,0(n),g2,1(n),...,g2,l(n),...,g2,L-1(n)]；
D3、组合滤波器抽头权向量的更新
如果当前时刻n的混合参数α(n)大于等于阈值σ，所述的阈值σ的取值为4；
则下一时刻n+1的大步长滤波器的抽头权向量W1(n+1)和小步长滤波器的抽头权向量W2(n+1)分别为：
如果当前时刻n的混合参数α(n)小于阈值σ，则n+1时刻的大步长滤波器的抽头权向量W1(n+1)和小步长滤波器的抽头权向量W2(n+1)分别为：
其中，μ1为大步长滤波器的步长，其取值为0.20～0.85；μ2为小步长滤波器的步长，其取值为0.020～0.085；δ为滤波器抽头权向量的正则化参数，其取值为0.001～0.01；
E、滤波器的权重更新
更新得到下一时刻n+1的混合参数值α(n+1)：
其中，μα是一个常数，取值为10；
再更新得到下一时刻n+1的大步长滤波器权重λ(n+1)：
F、滤波器的权重限定
如果下一时刻n+1的混合参数α(n+1)小于阈值σ的负值，则令α(n+1)＝-σ，λ(n+1)＝0；如果下一时刻n+1的混合参数α(n+1)大于等于阈值σ，则令α(n+1)＝σ，λ(n+1)＝1；
G、重复
令n＝n+1，重复步骤A-F的操作，直至通话结束。
仿真实验
为了验证本发明的有效性，进行了仿真实验，并与文献1中大步长IPAPLM算法和小步长IPAPLM算法进行对比。
仿真实验的远端信号x(n)为图3的语音信号，采样点个数800000。脉冲响应长度即滤波器抽头数L＝512。实验的背景噪声为高斯白噪声，信噪比为20dB。并且在麦克风接收到的近端信号中加入采集得到的冲击干扰。
图1是仿真实验的声学回声信道图，它是安静密闭房间构成的通信系统的稀疏信道。图2是仿真实验的网络回声信道图。
仿真实验时，三种方法的参数具体取值如表1。
表1各算法仿真实验的参数
| 大步长IPAPLM | μ＝1,δ＝0.001,ρ＝0.01,δ<sub>p</sub>＝0.001 |
| 小步长IPAPLM | μ＝0.2,δ＝0.001,ρ＝0.01,δ<sub>p</sub>＝0.001 |
| 本发明(CIPAPLMWT) | μ<sub>1</sub>＝1,μ<sub>2</sub>＝0.2,δ＝0.001,ρ＝0.01,δ<sub>p</sub>＝0.001 |
仿真结果通过独立运行50次平均得到。
图4是大步长IPAPLM算法和小步长IPAPLM算法和本发明(CIPAPLMWT)方法在声学回声信道中的仿真实验归一化稳态失调(MSD)曲线。由图4可以看出，本发明对语音具有很好的稳定性，与大步长IPAPLM算法和小步长IPAPLM算法相比。本发明在2×105迭代处收敛开始达到稳态，大步长IPAPLM算法在2×105迭代处收敛开始达到稳态，小步长IPAPLM算法在4×105迭代处收敛开始达到稳态；本发明的收敛达到稳态后，如6×105迭代时，大步长IPAPLM算法MSD的稳态误差值为-72dB，小步长IPAPLM算法MSD的稳态误差值为-60dB，而本发明的MSD的稳态误差值仅为-78。
图5是大步长IPAPLM算法和小步长IPAPLM算法和本发明(CIPAPLMWT)方法在网络回声信道中的仿真实验归一化稳态失调(MSD)曲线。由图5可以看出，当信道为网络回声信道时，本发明对语音也具有很好的稳定性，本发明对语音具有很好的稳定性，与大步长IPAPLM算法和小步长IPAPLM算法相比。本发明在4×105迭代处收敛开始达到稳态，大步长IPAPLM算法在4×105迭代处收敛开始达到稳态，小步长IPAPLM算法在6×105迭代处收敛开始达到稳态；本发明的收敛达到稳态后，如6×105迭代时，大步长IPAPLM算法MSD的稳态误差值为-58dB，小步长IPAPLM算法MSD的稳态误差值为-48dB，而本发明的MSD的稳态误差值仅为-64。
由图4和图5分析可知，本发明与大步长IPAPLM算法的收敛速度基本相同，比小步长IPAPLM算法的收敛速度快；本发明的稳态误差低于比大步长IPAPLM算法及小步长IPAPLM算法的稳态误差；其回声消除效果好。
Claims (1)
1.一种基于凸组合的M估计成比例类仿射投影的回声消除方法，其步骤如下：
A、远端信号采样
A1、将当前时刻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代表转置运算；
A2、将当前时刻n的滤波器输入向量X(n)通过凸组合自适应回声消除滤波器滤波后分别得到当前时刻n的大步长滤波器值y1(n)，y1(n)＝W1 T(n)X(n)和当前时刻n的小步长滤波器值y2(n)，y2(n)＝W2 T(n)X(n)；其中，W1(n)＝[w1,0(n),w1,1(n),...w1,l(n)...,w1,L-1(n)]T，为当前时刻n的凸组合自适应回声消除滤波器中的大步长滤波器的抽头权向量，其初始值为零向量；w1,l(n)为大步长滤波器的第l个抽头权系数，也是大步长滤波器向量中的第l个元素；W2(n)＝[w2,0(n),w2,1(n),...w2,l(n)...,w2,L-1(n)]T为当前时刻n的小步长滤波器的抽头权向量，其初始值为零向量；w2,l(n)为小步长滤波器的第l个抽头权系数，也是小步长滤波器向量中的第l个元素；
B、凸组合
将当前时刻n的大步长滤波器值y1(n)和小步长滤波器值y2(n)，通过当前时刻n的大步长滤波器权重λ(n)进行凸组合得到当前时刻n的组合滤波器值y(n)，y(n)＝λ(n)y1(n)+(1-λ(n))y2(n)；
其中，大步长滤波器权重λ(n)的计算公式为取值范围为0～1；α(n)为当前时刻n的混合参数，其初始值为0；
C、回声信号消除
将近端麦克风拾取的当前时刻n的带回声的近端信号d(n)，与自适应滤波器的当前时刻n的输出值y(n)相减后得到当前时刻n的残差信号e(n)，e(n)＝d(n)-y(n)；再将残差信号e(n)回送给远端；
同时，将当前时刻n的近端信号d(n)与当前时刻n的大步长滤波器值y1(n)和小步长滤波器值y2(n)相减，分别得到当前时刻n的大步长残差信号e1(n)和小步长残差信号e2(n)：
e1(n)＝d(n)-y1(n)，e2(n)＝d(n)-y2(n)；
D、滤波器抽头权向量更新
D1、M估计函数计算
将当前时刻n到时刻n-Nw+1之间的大步长残差信号平方e1 2(n)，e1 2(n-1)，…,e1 2(n-Nw+1)，构成当前时刻n估计窗内的大步长残差信号平方序列Ae1(n)：
Ae1(n)＝[e1 2(n),e1 2(n-1),...,e1 2(n-Nw+1)]
其中，Nw为估计窗的长度，其取值范围为5～15；
再由当前时刻n估计窗内的大步长残差信号平方序列Ae1(n)，计算出当前时刻n的大步长残差信号的方差
其中，λ为遗忘因子，其取值范围为0.800～0.999,C为常数，C＝1.483(1+5/(Nw-1))，med(·)表示取中间值的运算；
然后由当前时刻n的大步长残差信号的方差得出当前时刻n的大步长M估计阀值参数ξ1(n)，
由下式计算出滤波器当前时刻n的大步长抽头权向量更新的M估计函数值
其中，的初始值即为0；
由当前时刻n到时刻n-Nw+1之间的小步长残差信号平方值e2 2(n)，e2 2(n-1)，…，e2 2(n-Nw+1)，构成当前时刻n估计窗内的小步长残差信号平方序列Ae2(n)：
Ae2(n)＝[e2 2(n),e2 2(n-1),...,e2 2(n-Nw+1)]
再由当前时刻n估计窗内的小步长残差信号平方序列Ae2(n)，计算出当前时刻n的小步长残差信号的方差
然后由当前时刻n的小步长残差信号的方差得出当前时刻n的小步长M估计阀值参数ξ2(n)，
由下式计算出滤波器当前时刻n的小步长抽头权向量更新的M估计函数值
其中，的初始值即为0；
D2、类仿射投影
将当前时刻n到时刻n-P+1的自适应滤波器输入向量X(n),X(n-1),..,X(n-P+1)，构成当前时刻n的类仿射投影输入矩阵
将当前时刻n到n-P+1时刻的大步长抽头权向量更新的M估计函数值构成大步长抽头权向量更新的类仿射投影M估计向量ψ1(n)，
将当前时刻n到n-P+1时刻的小步长抽头权向量更新的M估计函数值构成小步长抽头权向量更新的类仿射投影M估计向量ψ2(n)，
其中，P为投影阶数，其取值为2、4、8；
D3、成比例控制因子计算
由当前时刻n大步长滤波器的第l个抽头权系数w1,l(n)，得出当前时刻n第l个大步长成比例控制因子g1,l(n)，
其中，κ为成比例调整参数，其取值范围为：-1≤κ＜1；ε为成比例限制参数，取值为0.01-0.001，作用是防止公式分母为0；
再计算出，当前时刻n的大步长成比例矩阵G1(n)，
G1(n)＝diag[g1,0(n),g1,1(n),...,g1,l(n),...,g1,L-1(n)]，其中diag[·]表示构造对角矩阵；
由当前时刻n小步长滤波器的第l个抽头权系数w2,l(n)，得出当前时刻n第l个小步长成比例控制因子g2,l(n)，
再计算出，当前时刻n的小步长成比例矩阵G2(n)，
G2(n)＝diag[g2,0(n),g2,1(n),...,g2,l(n),...,g2,L-1(n)]；
D3、组合滤波器抽头权向量的更新
如果当前时刻n的混合参数α(n)大于等于阈值σ，所述的阈值σ的取值为4；
则下一时刻n+1的大步长滤波器的抽头权向量W1(n+1)和小步长滤波器的抽头权向量W2(n+1)分别为：
如果当前时刻n的混合参数α(n)小于阈值σ，则n+1时刻的大步长滤波器的抽头权向量W1(n+1)和小步长滤波器的抽头权向量W2(n+1)分别为：
其中，μ1为大步长滤波器的步长，其取值为0.20～0.85；μ2为小步长滤波器的步长，其取值为0.020～0.085；δ为滤波器抽头权向量的正则化参数，其取值为0.001～0.01；
E、滤波器的权重更新
更新得到下一时刻n+1的混合参数值α(n+1)：
其中，μα是一个常数，取值为10；
再更新得到下一时刻n+1的大步长滤波器权重λ(n+1)：
F、滤波器的权重限定
如果下一时刻n+1的混合参数α(n+1)小于阈值σ的负值，则令α(n+1)＝-σ，λ(n+1)＝0；如果下一时刻n+1的混合参数α(n+1)大于等于阈值σ，则令α(n+1)＝σ，λ(n+1)＝1；
G、重复
令n＝n+1，重复步骤A-F的操作，直至通话结束。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810816058.1A CN109102794A (zh) | 2018-07-24 | 2018-07-24 | 基于凸组合的m估计成比例类仿射投影的回声消除方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810816058.1A CN109102794A (zh) | 2018-07-24 | 2018-07-24 | 基于凸组合的m估计成比例类仿射投影的回声消除方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN109102794A true CN109102794A (zh) | 2018-12-28 |
Family
ID=64847371
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810816058.1A Withdrawn CN109102794A (zh) | 2018-07-24 | 2018-07-24 | 基于凸组合的m估计成比例类仿射投影的回声消除方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109102794A (zh) |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110138990A (zh) * | 2019-05-14 | 2019-08-16 | 浙江工业大学 | 一种消除移动设备VoIP电话回声的方法 |
| CN110191245A (zh) * | 2019-07-10 | 2019-08-30 | 西南交通大学 | 一种基于时变参数的自适应回声消除方法 |
| CN110572525A (zh) * | 2019-10-30 | 2019-12-13 | 西南交通大学 | 一种用于语音通信的自适应通信回声消除方法 |
| CN111243568A (zh) * | 2020-01-15 | 2020-06-05 | 西南交通大学 | 一种凸约束自适应回声消除方法 |
| CN113113038A (zh) * | 2021-05-19 | 2021-07-13 | 北京字节跳动网络技术有限公司 | 回声消除方法、装置及电子设备 |
| CN117789691A (zh) * | 2024-02-27 | 2024-03-29 | 苏州大学 | 一种用于电话通信的回声消除方法及系统 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104683614A (zh) * | 2015-03-24 | 2015-06-03 | 西南交通大学 | 基于m估计的记忆成比例仿射投影自适应回声消除方法 |
| CN106060295A (zh) * | 2016-05-17 | 2016-10-26 | 西南交通大学 | 一种凸组合系数差的成比例仿射投影回声消除方法 |
| CN109040497A (zh) * | 2018-07-24 | 2018-12-18 | 西南交通大学 | 一种基于m估计的成比例类仿射投影自适应回声消除方法 |
-
2018
- 2018-07-24 CN CN201810816058.1A patent/CN109102794A/zh not_active Withdrawn
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104683614A (zh) * | 2015-03-24 | 2015-06-03 | 西南交通大学 | 基于m估计的记忆成比例仿射投影自适应回声消除方法 |
| CN106060295A (zh) * | 2016-05-17 | 2016-10-26 | 西南交通大学 | 一种凸组合系数差的成比例仿射投影回声消除方法 |
| CN109040497A (zh) * | 2018-07-24 | 2018-12-18 | 西南交通大学 | 一种基于m估计的成比例类仿射投影自适应回声消除方法 |
Non-Patent Citations (4)
| Title |
|---|
| P. SONG 等: "Affine-Projection-Like M-Estimate Adaptive Filter for Robust Filtering in Impulse Noise", 《IEEE TRANSACTIONS ON CIRCUITS AND SYSTEMS II: EXPRESS BRIEFS》 * |
| ZHAO HAIQUAN 等: "A novel sign adaptation scheme for convex combination of two adaptive filters", 《AEU - INTERNATIONAL JOURNAL OF ELECTRONICS AND COMMUNICATIONS》 * |
| 周彤: "类仿射投影算法的改进研究", 《中国优秀硕士学位论文全文数据库（信息科技辑）》 * |
| 夏可智: "基于M估计的鲁棒自适应滤波算法在回声消除中的应用", 《中国优秀硕士学位论文全文数据库（信息科技辑）》 * |
Cited By (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110138990A (zh) * | 2019-05-14 | 2019-08-16 | 浙江工业大学 | 一种消除移动设备VoIP电话回声的方法 |
| CN110191245A (zh) * | 2019-07-10 | 2019-08-30 | 西南交通大学 | 一种基于时变参数的自适应回声消除方法 |
| CN110191245B (zh) * | 2019-07-10 | 2021-06-22 | 西南交通大学 | 一种基于时变参数的自适应回声消除方法 |
| CN110572525A (zh) * | 2019-10-30 | 2019-12-13 | 西南交通大学 | 一种用于语音通信的自适应通信回声消除方法 |
| CN111243568A (zh) * | 2020-01-15 | 2020-06-05 | 西南交通大学 | 一种凸约束自适应回声消除方法 |
| CN111243568B (zh) * | 2020-01-15 | 2022-04-26 | 西南交通大学 | 一种凸约束自适应回声消除方法 |
| CN113113038A (zh) * | 2021-05-19 | 2021-07-13 | 北京字节跳动网络技术有限公司 | 回声消除方法、装置及电子设备 |
| CN117789691A (zh) * | 2024-02-27 | 2024-03-29 | 苏州大学 | 一种用于电话通信的回声消除方法及系统 |
| CN117789691B (zh) * | 2024-02-27 | 2024-05-03 | 苏州大学 | 一种用于电话通信的回声消除方法及系统 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109102794A (zh) | 基于凸组合的m估计成比例类仿射投影的回声消除方法 | |
| CN109040499B (zh) | 一种抗冲击干扰的自适应回声消除方法 | |
| CN104683614B (zh) | 基于m估计的记忆成比例仿射投影自适应回声消除方法 | |
| CN104506746B (zh) | 一种改进的凸组合解相关成比例自适应回声消除方法 | |
| CN110191245B (zh) | 一种基于时变参数的自适应回声消除方法 | |
| CN110572525B (zh) | 一种用于语音通信的自适应通信回声消除方法 | |
| CN109040497B (zh) | 一种基于m估计的成比例类仿射投影自适应回声消除方法 | |
| CN103369162A (zh) | 一种低复杂度的电话回声自适应消除方法 | |
| Huang et al. | Practically efficient nonlinear acoustic echo cancellers using cascaded block RLS and FLMS adaptive filters | |
| CN104144269B (zh) | 一种基于解相关的成比例自适应电话回声消除方法 | |
| CN107871510A (zh) | 一种零吸引变核宽度的最大熵回声消除方法 | |
| CN115457928A (zh) | 一种基于神经网络双讲检测的回声消除方法及系统 | |
| CN109089004B (zh) | 一种基于相关熵诱导的集员自适应回声消除方法 | |
| CN109151237A (zh) | 基于零吸引的仿射投影自适应回声消除方法 | |
| CN106170049B (zh) | 一种有偏补偿的归一化子带自适应回声消除方法 | |
| Shen et al. | An improved LMS adaptive filtering algorithm and its analysis | |
| CN106157965B (zh) | 一种基于权向量重用的零范数集员仿射投影自适应回声消除方法 | |
| CN113873090B (zh) | 一种稳健估计仿射投影样条自适应回声消除方法 | |
| CN107527622B (zh) | 一种基于最大熵的变步长符号子带回声消除方法 | |
| Kim et al. | Comparison of two proposed methods in adaptive noise canceling | |
| CN106060295B (zh) | 一种凸组合系数差的成比例仿射投影回声消除方法 | |
| CN113409806B (zh) | 一种基于反正切函数的零吸引回声消除方法 | |
| CN114400016A (zh) | 基于自适应解相关和变步长成比例m估计的回声消除方法 | |
| CN110767245B (zh) | 基于s型函数的语音通信自适应回声消除方法 | |
| CN107578782B (zh) | 一种变核宽度的最大熵子带回声消除方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| WW01 | Invention patent application withdrawn after publication | ||
| WW01 | Invention patent application withdrawn after publication |
Application publication date: 20181228 |