一、speex回声消除aec小析：
频域自适应算法采用了分块处理的思想，以减少高阶自适应滤波器的计算复杂度，多延迟自适应滤波器（MDF)则更一般可以分M块来处理，其中块的长度比自适应滤波器的阶数更小。
后置滤波器系数一直在更新，单讲时前置滤波器用后置滤波器系数更新，双讲时前置滤波器系数不更新。
如果近端远端的数据线性比较好的话，用speex处理回声效果还不错。若是非线性严重效果就差了。另外memath认为相比其它的回声消除算法，speex对在滤波器更新时做了远端检测，回声估计比较好些。
以前写的文章，分享下，不当之处，请指正。
二、源码解读：
标记：N=2*frame_size 其中frame_size为帧长
M= filter_length/frame_size，其中filter_length为回声消除拖尾
1、按照公式（1）把近端信号in[frame_size]直流去除掉，存入缓存input[frame_size]。
2、按照公式（2）对远端与近端信号分别进行预加重处理，远端信号存入x[frame_size]，并对近端延迟的信号进行frame_size的交叠，以下均按预加重处理后信号处理。
3、频域X[M*frame_size]后面M-1个频域段的历史数据后移，然后对远端信号x使用FFT变换到频域X[M*frame_size]前面部分。
4、求远端信号x时域的能量Sxx和以及其功率谱Xf[frame_size+1]。
5、前置滤波器foreground与远端信号x的输出再与近端信号相减得到误差信号存入e[N]前半部分，并求其能量方差Sff。
6、求取后置滤波器W各段的系数幅值
并进行归一化处理
，
其中初始化max_sum=1。
更新系数w(i)=w(i)+PHI(i)，PHI由远端信号频域X(i)与误差E(i)复数相乘并乘以power_1*prop(i)计算得到，其中power_1为自适应步长因子。
7、为防止循环卷积，对于第1段以及当前帧数对M-1求余为j-1的j(j=0...M-1)段作处理。（把W[j*N]段做反傅里叶变换后的后一半部分全部置0，再做傅里叶变换回来）
8、计算远端信号与后置滤波器的输出并反傅里叶变换到时域y，前置、后置滤波器输出之间的误差与它们误差的方差Dbf；计算近端信号减去后置滤波器输出差值的方差See。
9、分别计算平均能量差别与平均方差差别Davg1= 0.6*Davg1 + 0.4*(Sff-See); Davg2 = 0.85*Davg2 + 0.15*(Sff-See);Dvar1 =0.36*Dvar1 + 0.16*Sff*Dbf；Dvar2 = 0.7225*Dvar2 + 0.0225*Sff*Dbf；以便判断是否要更新前置后置滤波器系数。
如果
;其中在公式(3)中d1=1，d2=0.5，d3=0.25任一个条件成立，就要更新前置滤波器系数。若是要更新，则把后置滤波器的系数复制给前置滤波器，并平滑前置滤波器的输出误差与后置滤波器的输出误差得到总的误差。Davg1，Davg2，Dvar1，Dvar2重置为0。
10、若
在公式(4)中d1=4，d2=4，d3=4任一个条件成立，就要更新后置滤波器系数。则把前置滤波器的系数复制给后置滤波器，并把前置滤波器的输出误差复制给后置滤波器的输出误差，再更新后置滤波器输出前部分。把See=Sff,Davg1，Davg2，Dvar1，Dvar2重置为0。
11、自适应滤波器的输出out[frame_size]由近端信号减去前置滤波器的输出，并且去预加重处理公式（4）；保存前置滤波器生成人造回声信号并去预加重处理公式（5）。
12、误差缓存区前部分拷贝给后部分，并置前部分为0。
13、计算前置滤波器输出与前置滤波器的输出的互相关Sey，后置滤波器输出的自相关Syy，近端信号的自相关Sdd；计算分别前置、后置滤波器输出误差的频域数据E与Y，及其功率谱Rf与Yf。
14、若!(Syy>=0 且Sxx>=0且 See >= 0)|| !(Sff < N*1e9且 Syy < N*1e9且 Sxx< N*1e9),则screwed_up=screwed_up+50;否则Sff>Sdd+10000N时，screwed_up=screwed_up+1;否则screwed_up=0；当screwed_up>=50时，重置全部参数。
15、计算远端信号的能量Sxx及其功率谱Xf;平滑远端信号的功率
16、计算相关系数及标准差
Eht=Rf(i)-Eh(i);
Yht=Yf(i)-Yh(i);
Peyt=Peyt+Eht*Yht;
Pyyt=Pyyt+Yht*Yht;
更新：Eh(i)=(1- N/Fs)*Eh(i)+N/Fs*Rf(i)
Yh(i)=(1- N/Fs)*Yh(i)+N/Fs*Yf(i)
最后Pyyt=sqrt(Pyyt), Peyt=Peyt/Pyyt;
17、计算相关系数更新率alpha=min(Syy*2*frame_size/Fs,See*0.5*frame_size/Fs)/See; alpha_1=1-alpha。
18、更新相关系数Pey=alpha_1*Pey+ alpha*Peyt; Pyy=alpha_1*Pyy+ alpha*Pyyt;
19、Pyy=max(Pyy,1);Pey=max(Pey,0.005*Pyy);Pey=min(Pey,Pyy)
20、泄漏估计leak_estimate=Pey/Pyy,若leak_estimate>16383,则令leak_estimate=32767。
21、计算错误比率RER = (0.0001*Sxx + 3*leak_estimate*Syy)/ See;对RER作上下限限制。
RER=max(RER, Sey*Sey/(1+See*Syy);RER=min(RER,0.5);
22、若sum_adapt>M并且leak_estimate>0.03,则更新自适应滤波器步长因子:先求后置滤波器的回声消除功率 r=min(leak_estimate*Yf(i),0.5(Rf(i)+1))*0.7+RER*(Rf(i)+1)*0.3，
最优步长为power_1=r/((Rf(i)+1)*(power(i)+10))；否则，adapt_rate=0,若Sxx>1000N时，
adapt_rate=min(Sxx, See)*0.25/See;更新power_1= adapt_rate/(power+10)，
总自适应次数sum_adapt=sum_adapt+adapt_rate。
三、speex回声消除流程图：
memath回声消除系列文章：
speex与webrtc回声消除小结
QQ、YY与webRTC回声消除效果对比分析与展望
webrtc aec3效果对比aec与aecm（webrtc M64 20180115版本）