WEBRTC AEC3算法原理

延迟估计算法
原理
AEC3的延迟估计算法与AEC的非线性处理的延迟估计算法思想一致，因为回声能量是呈指数衰减，所以计算滤波器能量最大块作为延迟估计值，但是比AEC的延迟估计算法复杂的多
AEC3延迟估计模块由步长为0.7的5个时域NLMS自适应滤波器组成，每个NLMS滤波器默认32块，每块16个sample共 512点，5个滤波器之间互相重叠8块，这里的重叠指的是输入的信号在时间上重叠
滤波器的输入信号是经过分频后的0~16kHz低频段信号然后经过4倍下采样相当于采样率为4000Hz的信号，5个nlms滤波器可以估计32165 - 8164 =2048个samples，也就是说该延迟估计算法最多能估计2048/4000 = 512ms的延迟

流程
1、滤波器更新
滤波器更新算法：NLMS更新公式
滤波器更新条件：

近端语音帧没有超过32000或-32000
远端能量大于一定阈值
2、计算每个滤波器的延时值
延时值：将每个滤波器系数能量最大块作为延时值
延时值可信度：

延时值落在滤波器>2且<10范围内
误差能量小于近端能量的0.2倍
滤波器是否更新
3、从5个滤波器中挑出最佳延时值估计
挑选方法：滤波器更新且可信度为真，且accuracy最大作为延时值估计
accuracy：近端能量-误差能量

4、延时值过滤
使用大小为250的历史窗口统计历史延时值，统计历史窗口中出现次数最多的延时值，当出现次数大于20次认为这是一个质量较高的估计值，反之认为是较差的估计值

5、时钟漂移检测
质量好的延时值会用来做时钟漂移检测，时钟漂移检测使用大小为3的历史窗口判断时钟漂移是否产生，窗口内的值时前3个历史时刻延时估计值：
判断标准：d1，d2，d3分别时历史窗口历史值与当前时刻延时值之差
增量漂移：d1==-1且d2==-2 或d1==-2且d2==-1 d3==-3确认产生反之可能产生

const bool probable_drift_up =
      (d1 == -1 && d2 == -2) || (d1 == -2 && d2 == -1);
const bool drift_up = probable_drift_up && d3 == -3;

减量漂移：d11且d22 或 d12且d21 d3==3确认产生反之可能产生

const bool probable_drift_down = (d1 == 1 && d2 == 2) || (d1 == 2 && d2 == 1);
const bool drift_down = probable_drift_down && d3 == 3;

目前aec3只做了时钟漂移检测，没有做相应处理

6、对齐远端数据
质量好的延时值和质量差的延时值都会用来进行远端对齐

另外，当上一次延时值和当前延时值的质量都是高的情况下会有一个迟滞量，即延迟增量小于某一阈值时认为延时相较于上一次没有改变

7、处理回声路径变化(延迟改变，增益改变，时钟漂移)
主辅滤波器处理声路径变化
将非线性增益设置为初始状态
自适应回声消除
原理
线性回声消除使用主辅滤波器结构，总体来说，主滤波器可以看作是一个精简版的频域分块的kalman滤波器，只不过精简了系数误差协方差的计算，计算的是滤波器抽头系数误差协方差的平均，而不是计算每块滤波器抽头系数的误差协方差，回声路径跟踪能力和双工效果比PBFDAF效果要好一些，但比起kalman滤波器又差了点，相当于在计算量上和效果上做了权衡。辅滤波器是PBFDAF，会一直更新，能快速对回声路径进行追踪，当辅滤波器发散时主滤波器系数拷贝至辅滤波器。因此这样的双滤波器结构既能够有效消除线性回声，有着不错的双工效果，又能够快速的跟踪回声路径的变化。
一般来说，滤波器的抽头长度应该设计的足够长，以完全逼近回声路径，但是为了快速的收敛以及降低计算量，AEC3自适应滤波器抽头长度设计的较短，只有12块，所以AEC3的自适应滤波器主要是消除早期的回声，后期的混响回声是通过混响模型进行估计然后用NLP进行消除。

流程
1、分析远端信号中的窄带信号
窄带信号检测算法：

远端信号中某个频率点的幅度大于左右频点幅度的3倍，则认为是弱窄带频率分量
远端信号中某个频率点的幅度大于左右频点的100倍，则认为该频率点是强窄带分量
自适应滤波器对窄带信号处理：

连续5个block远端信号都含有弱窄带频率分量后，滤波器在该频率点的左右2个频率点的系数都不会更新
连续10个block远端信号都含有弱窄带频率分量后，滤波器系数全部都不更新
如果远端信号含有强窄带频率分量，则非线性处理的高频部分抑制系数限制为0.001
那为什么要对这种窄带信号做处理呢？
因为窄带信号会导致自适应滤波器朝着错误的方向收敛，最终导致滤波器发散和回声泄漏。
假设远端信号为正弦信号，自适应滤波器不会朝着真实的回声路径去收敛，而是朝着单一频率方向收敛，最终导致自适应滤波器的效果像一个陷波器，基于LMS算法的滤波器如NLMS，FDAF等都有这个问题，如下图所示。

2、主滤波器更新
更新方式1：通过计算主滤波器系数的失调增益调整滤波器系数

    mis = mis + 0.1*( e2/y2 - mis)
    scale = 2 / sqrt(mis)
    H(t)=H(t) * scale

通过公式可以看出：如果误差信号一直大于某个阈值会导致失调量持续增加，当mis大于10时，计算失调增益scale来调整滤波器系数
总体来说该算法会一直限制主滤波器误差输出在与近端能量有关的一定范围内，即滤波器不会发散

更新方式2：通过系数更新公式更新滤波器系数

    Erl(t) = sum{H(n) * H(n)}
    mu(t) = H_error(t) / (0.5 * H_error(t) * X2(t) + N * E2(t)) —— (X2(t) > noise_gate)
    mu(t) = 0  —— (X2(t) <= noise_gate)      
    H_error(t) = H_error(t) - 0.5 * mu(t) * X2(t) * H_error(t)
    G(t) = mu(t) * E(t)
    factor = 0.00005f —— (E2_shadow>=E2_main)
    factor = 0.05f —— (E2_shadow<E2_main)
    H_error(t+1) = H_error(t) + factor * Erl(t)
    H(t+1)=H(t)+G(t) * conj(X(t))

注意：这个公式是我在代码中的注释，并不规范，只是记录了aec3滤波器的思想
通过与kalman滤波器仔细对比可以发现，其实这个滤波器是kalman滤波器的精简版，其更新公式与kalman滤波器是一样的，只不过精简了系数误差协方差的计算，计算的是滤波器抽头系数误差协方差的平均，而不是计算每块滤波器抽头系数的误差协方差，把mu(t)带入到 G(t) = mu(t)E(t)中，会发现G(t)就是kalman滤波器增益K的计算公式，H_error(t)的计算也与kalman滤波器的系数协方差更新公式一样，NE2(t)是测量噪声，factor*Erl(t)是过程噪声

另外还可以看出：

当远端信号能量小于预设的噪音门限值时，mu=0，代表不更新滤波器
当主滤波器误差小于辅滤波器误差时，主滤波器收敛较好，factor=0.00005，设置了比较小的过程噪声，当主滤波器误差大于辅滤波器误差时，主滤波器收敛较差，factor=0.05，设置了比较大的过程噪声，这也是保证滤波器不会发散的一个手段
3、辅滤波器更新
当主滤波器误差连续超过5个block小于辅滤波器误差时，认为辅滤波器发散，将主滤波器系数拷贝至辅滤波器，更新辅滤波器
按PBFDAF公式更新辅滤波器
4、选择线性回声消除结果
选择线性回声消除结果使用主滤波器误差还是使用辅滤波器误差，然后进行平滑过渡
假如配置了使用辅助滤波器的输出，则满足以下情况使用辅滤波器误差输出：

当近端信号大于一定阈值且辅滤波器误差输出能量小于0.9倍主滤波器误差输出能量
主滤波器输出信号能量大于一定阈值或者辅滤波器输出能量大于一定阈值，即辅助滤波器检测到的回声大于一定量
假如没有配置使用辅助滤波器的输出，则满足以下情况使用辅助滤波器误差输出

辅滤波器误差能量小于主滤波器误差能量且近端信号能量小于主滤波器误差能量
非线性回声消除
原理
AEC3的非线性回声消除原理主要思想是通过估计残余回声以及噪声，设计一个滤波器，以达到去除残余回声的目的，这实际上与降噪算法非常类似。

流程
1、AEC状态更新和参数计算
1、根据主辅滤波器的误差与近端信号分析滤波器收敛与发散状态

当主滤波器误差小于近端信号的一半且近端信号低于收敛阈值时，认为主滤波器收敛
当辅滤波器误差小于近端信号的0.05倍且近端信号低于收敛阈值时，认为辅滤波器收敛
当两个滤波器中较小误差大于1.5倍的近端信号且近端信号大于某阈值时，认为两个滤波器都发散
2、根据主滤波器的时域响应分析主滤波器的属性，主要分析如下参数：

寻找滤波器响应峰值
计算峰值的延迟filter_delays_blocks_
计算峰值的增益，即回声路径最大增益值
计算第二峰值分析滤波器响应是否随时间变化一致
any_filter_consistent： 主滤波器时域响应是否随时间变化一致
max_echo_path_gain： 最大回声路径增益值
3、计算相对于滤波器开始的直接路径延迟
4、更新远近端语音活动统计
5、根据回声混响模型，计算远端信号加上混响后的功率
6、检测回声是否过饱和：通过最大回声增益值和远端信号最大值来判断

7、更新erle和erl：根据erle的估计计算线性滤波器的好坏，以及通过erle计算残余回声

ERL : echo return loss = Y2 / X2 近端功率/远端功率 高表示回声强，低表示回声弱
ERLE : echo return loss enhancement = Y2 / E2 近端功率/误差功率 衡量AEC的性能，计算残余回声
RERL : Residual Echo Return Loss = ERL + ERLE
8、更新AEC工作状态标志，在AEC由初始状态转为正常工作状态erle估计器会重置

通过远端语音是否活动和近端语音是否过饱和判断AEC的状态，初始状态或正常工作状态
从初始状态进入正常工作状态会将线性滤波和非线性处理模块的状态都从初始状态转换为正常工作状态
9、transparent_mode检测：检测滤波器是否工作正常

总的来说就是滤波器工作情况良好的情况下 transparent_mode = false
在滤波器工作不合理且近端语音未过饱和超过6秒 transparent_mode = true
在transparent_mode = true 时，滤波器性能很差，在非线性处理时会使用不同函数计算残余回声
10、计算滤波器工作性能，来决定是否使用线性滤波的结果作为非线性处理的输入，以及使用不同的非线性抑制函数
在滤波器工作良好时，使用线性滤波结果计算残余回声，工作良好判定需同时满足以下条件：

初始化后滤波器更新次数超过0.4秒
观测到任一滤波器收敛 或 延时估计模块估计出了延时
非transparent mode
在滤波器工作较差时，非线性处理的输入为近端信号，并且使用不同的函数计算残余回声
11、估计后期回声混响模型

2、计算频谱
加窗、FFT，计算频谱
计算近端信号、误差信号、线性回声信号频谱

3、估计残余回声
滤波器工作良好情况下：

如果回声过饱和，则R2=Y2，即认为误差输入中都是残余回声
如果回声未过饱和，则R2 = S2_linear / ERLE + Reverb
滤波器工作较差情况下：

如果回声过饱和，则R2=Y2，即认为误差输入中都是残余回声
如果回声未过饱和，则R2 = X2 * 0.0001 + Reverb
这里Reverb就是根据估计出的混响模型计算出来的自适应滤波器未消除的后期混响回声

4、估计舒适噪声
更新公式：

N2 = Y2_smoothed < N2 ? (0.9*Y2_smoothed + 0.1*N2)* 1.0002f : N2*1.0002f

缓慢更新N2噪声能量，当近端语音能量小于噪音能量时，噪音能量快速降低至近端能量水平，当近端语音能量大于噪音能量时，噪音能量N2缓慢增加

5、近端语音活动检测
如果近端(误差)能量显著大于残余回声能量和舒适噪声能量，且连续超过12个block时则认为产生了近端语音，最短保持时间设置为50个block

6、选择计算非线性增益的输入参数
根据近端功率谱、线性回声估计功率谱、残余回声估计功率谱、舒适噪声功率谱计算非线性增益值

如果使用线性滤波的结果作为非线性处理的输入，则非线性处理的近端功率谱使用线性滤波的误差输出，回声功率谱使用线性滤波估计的线性回声功率谱
如果不使用线性滤波的结果作为非线性处理的输入，则非线性处理的近端功率谱使用近端功率谱，回声功率谱使用估计的残余回声功率谱
7、计算非线性增益值
低频段非线性增益值核心计算函数GainToNoAudibleEcho代码如下:

const auto& p = dominant_nearend_detector_->IsNearendState() ? nearend_params_
                                                               : normal_params_;
for (size_t k = 0; k < gain->size(); ++k) {
    float enr = echo[k] / (nearend[k] + 1.f);  // Echo-to-nearend ratio.
    float emr = echo[k] / (masker[k] + 1.f);   // Echo-to-masker (noise) ratio.
    float g = 1.0f;
    if (enr > p.enr_transparent_[k] && emr > p.emr_transparent_[k]) {
        g = (p.enr_suppress_[k] - enr) /
             (p.enr_suppress_[k] - p.enr_transparent_[k]);
        g = std::max(g, p.emr_transparent_[k] / emr);
    }
    (*gain)[k] = g;
}

enr为残余回声功率和误差功率比值，emr为残余回声功率和噪声功率比值
可以看出来：只有在enr和emr大于一定阈值的情况下才会计算增益，低于这个阈值增益为1，相当于不做处理，这是为了减小抑制效果，提升近端语音“透明度”，即双工效果。

8、对误差信号进行非线性抑制
使用非线性增益与误差信号频谱相乘
加上舒适噪声
IFFT
重叠相加获得去窗后的时域平滑的信号
与AEC对比
AEC3的延迟估计算法与AEC非线性处理延迟估计算法类似，但AEC就简单的多，虽然AEC也有延迟估计算法，但实际测试还是有不少机器会出现回声消除效果不好甚至无法消除的情况，但用AEC3就不会出现完全无法消除的情况，我们使用了大量机型测试，AEC3都能够保证算法能正常工作，回声能正常消除，我认为这得益于AEC3的延迟估计算法。
AEC的自适应滤波为了快速对回声路径跟踪，采用了大的固定步长的滤波器，在双工的情况下滤波器会发散，发散后只是简单的重置滤波器，而AEC3的使用双滤波器结构，而且对滤波器发散进行了处理，对比AEC的单个滤波器要强不少
AEC没有初始状态控制，在AEC刚开始工作时，算法还在学习阶段，这时候回声消除处理的不好，容易产生残留回声，AEC3在回声路径变化或者初始化时，会增大非线性抑制，以确保不会产生回声
总结
暂时只是对webrtc aec3的算法粗略的了解了一下，可以发现aec3算法对工程化做了大量的处理，例如延迟估计算法、处理回声路径变化、初始状态设置、线性回声消除以及非线性回声消除，其中的每一部分内还有更多复杂的细节及逻辑，因此不难看出Google推出aec3就是为了针对总类繁多的webrtc设备终端一站式解决适配问题。
有任何疑问，欢迎加微信交流：xu19315519614
最后附上我单独提取测试的AEC3代码： https://github.com/ewan-xu/AEC3

问题记录
AEC3的线性滤波器只有在远端信号达到一定阈值才会更新滤波器，这样可以防止远端的噪音信号干扰滤波器导致滤波器错误的跟踪回声路径，但也带来一个缺点，那就是在初始化时，如果远端的语音信号非常小，一直低于阈值，那么AEC3会认为是噪声，导致滤波器一直不收敛，因而回声无法消除的情况
AEC3算法舒适噪声估计和AEC一样，通过近端信号能量产生舒适噪声，但是都会产生一个问题，那就是当近端持续一直很大时（比如播放音乐），这时候产生的噪声就会非常大，在实际的通话中，会造成正向激励的环路，导致即使两边都没有通话，但是噪音会持续增大的情况
AEC3算法对于音乐的回声消除效果较差，因为对于音乐，有许多的乐器会产生这种类似于窄带信号的频率分量，在这种情况下AEC3的自适应滤波器基本上是不更新的，这也导致了自适应滤波器没法快速跟踪回声路径，导致回声消除效果较差
AEC3延时估计的参考信号是经过下采样至频率只有0~2000Hz的信号，对于人声来说没问题，但对于音乐或者某些信号2000Hz以下没有频率分量，这时候AEC3的延时估计算法就失效了
AEC3的非线性抑制算法为了提升双工效果，即近端语音“透明度”，在对每个子频带计算增益值时，enr和emr小于一定值就不进行抑制，即抑制系数设置为1，这样导致了AEC3在双工时会有残留的回声，而且残留的回声听感不平滑