1.本发明属于语音信号处理技术领域，具体涉及一种回声消除方法。
背景技术：
2.在语音交互系统中，回声消除质量直接影响用户语音识别交互体验，因此回声消除功能起到重要作用，如何有效增强回声抑制能力，是提升语音交互质量的关键性问题。
3.目前常见的回声消除方案是采用基于递推最小二乘法的自适应滤波器来估计回声路径传递函数，从而消除线性回声，并结合深度学习方法实现非线性回声抑制功能。这类方法的回声抑制能力容易受到复杂环境的影响，例如语音和回声的信噪比较低时，这类方法的回声抑制能力下降，语音失真加重。
技术实现要素：
4.为克服现有技术存在的缺陷，本发明公开了一种回声消除方法。
5.本发明所述回声消除方法，包括以下步骤:s1.将模拟参考信号和模拟麦克风信号转换成数字参考信号和数字麦克风信号；s2.通过短时傅里叶变换将数字麦克风信号x(n)和数字参考信号y(n)的时域信号分解成频域麦克风信号x(l,k)和频域参考信号y(l,k)；频域麦克风信号x(l,k)表示为；式中k为频率刻度，l为时域帧索引，d(l,k)为第l帧第k个频点的目标语音信号，为扬声器播放参考信号音频后麦克风拾取到的第l帧第k个频点的线性回声信号，v(l,k)为第l帧第k个频点的周围环境噪声信号；频域线性回声信号；上式中第l帧第k个频点的频域参考信号向量yv(l,k)=[y(l,k),y(1-1,k)...y(l-l 1,k)]y(l,k)为当前帧频域参考信号，y(l-1,k)
…
y(l-l 1,k)分别为当前帧之前的l-1帧，即第l-1，
…
，l-l 1帧的频域参考信号；wv(l,k)表示参考信号从扬声器传播到麦克风的理想频域回声路径传递函数，即wv(l,k)=[w
l
(0,k),w
l
(1,k)...w
l
(l-1,k)]，l为频域参考信号向量和理想频域回声路径传递函数的长度；s3. 将当前帧和过去l-1帧的频域参考信号组成的频域参考信号向量yv(l,k)=[y(l,k),y(1-1,k)...y(l-l 1,k)]存储到缓存器中；s4.对频域参考信号向量进行滤波，输出回声估计信号
；式中上标h表示共轭转置，为第l帧第k个频点的回声估计信号，yv(l,k)为第l帧第k个频点的频域参考信号向量,gv(l,k)为第l帧第k个频点的估计回声路径传递函数；对第一帧估计回声路径传递函数，取任意不为零的初值，对后续帧的估计回声路径传递函数gv(l,k)，更新公式如下gv(l,k)=gv(l-1,k) kv(l,k)e1(l,k)式中kv(l,k)为第l帧第k个频点的卡尔曼增益，由所述频域麦克风信号减去回声估计信号得到第一残余信号e1(l,k)，；s5.利用第二滤波器对第一残余信号e1(l,k)滤波，输出第二残余信号；首先根据测量噪声协方差矩阵r(l,k)和回声估计信号确定后验信噪比γ(l,k)=r(l,k)/t(l,k)其中t(l,k)为第l帧第k个频点的经过平滑处理后的平滑回声估计能量；r(l,k)为第l帧第k个频点的测量噪声协方差矩阵；；为s4步骤中求得的第l帧第k个频点的回声估计信号；β为第二平滑参数，为第l帧第k个频点的原始回声估计能量；e{}表示数学期望 ；接着由以下方程组，根据平滑回声估计能量t(l,k)和后验信噪比γ（l,k），进一步循环计算得到第二残余信号能量r
e2,e2
和先验信噪比；先验信噪比；第二残余信号e2(l-1,k)=g(l-1,k)e1(l,k)第二残余信号能量；第二残余信号平滑后能量；式中max{}表示取最大值，g(l,k)为第二滤波器参数，s(l,k)为第l帧第k个频点的经过平滑处理后的第二残余信号平滑后能量，e2(l-1,k)为第l-1帧第k个频点的第二残余信号，λ为权重因子，c为经验常数，η为第三平滑参数，初始状态下，l=1,第二滤波器参数和第二残余信号平滑后能量的初始值取设定的初值，计算先验信噪比，然后根据先验信噪比更新第二滤波器参数；
；式中g(l,k)为第l帧第k个频点的第二滤波器参数，更新后计算第二残余信号能量r
e2,e2
和第二残余信号平滑后能量s(l,k)，再重复循环计算；最后输出第二残余信号e2(l,k)=g(l,k)e1(l,k)；将第二残余信号转换成时域信号，完成回声消除。
[0006]
优选的，所述s2步骤中，通过短时傅里叶变换将数字麦克风信号x(n)和数字参考信号y(n)的时域信号分解成频域麦克风信号和频域参考信号。
[0007]
优选的，所述s4步骤中，卡尔曼增益的计算方式如下；式中上标t表示转置，p-(l,k)为第l帧第k个频点的先验误差协方差矩阵，p(l,k)为第l帧第k个频点的后验误差协方差矩阵，f为状态转换矩阵，q为过程噪声协方差矩阵，r(l,k)为第l帧第k个频点的测量噪声协方差矩阵，i为单位矩阵。
[0008]
优选的，所述测量噪声协方差矩阵r(l,k)计算方式如下:；式中e{
·
}表示求数学期望，α为第一平滑参数，r
e,e
(l,k)为第l帧第k个频点的第一残余信号的自相关矩阵。
[0009]
本发明利用优化的回声消除算法，可以有效提高回声抑制能力，提升语音交互质量。
附图说明
[0010]
图1是本发明所述回声消除方法的一个具体流程示意图；
[0011]
图2上方为本发明一个具体实施例中的输入信号时域波形示意图，图2下方为采用现有技术对输入信号进行回声消除的输出时域波形示意图；
[0012]
图3上方为本发明一个具体实施例中的输入信号时域波形示意图，图3下方为采用本发明所述方法对输入信号进行回声消除的输出时域波形示意图。
具体实施方式
[0013]
下面结合附图对本发明的具体实施方式进行详细说明，但并不用来限制本发明的保护范围。
[0014]
本发明的一个具体实施方式流程如图1所示。
[0015]
s1. 通过模拟数字转换器将模拟参考信号和模拟麦克风信号转换成数字参考信号和数字麦克风信号。
[0016]
可通过模拟数字转换器adc将时间连续、幅值连续的模拟信号转换为时间离散、幅值离散的数字信号,模拟信号转换为数字信号之后，进行后续的数字信号处理。
[0017]
s2. 通过短时傅里叶变换将数字麦克风信号x(n)和数字参考信号y(n)的时域信号分解成频域麦克风信号x(l,k)和频域参考信号y(l,k)。
[0018]
x(l,k)表示为；式中k为频率刻度，l为时域帧索引，d(l,k)为第l帧第k个频点的目标语音信号，为扬声器播放参考信号音频后麦克风拾取到的第l帧第k个频点的线性回声信号，v(l,k)为第l帧第k个频点的周围环境噪声信号；频域线性回声信号；
[0019]
上式中第l帧第k个频点的频域参考信号向量yv(l,k)=[y(l,k),y(1-1,k)...y(l-l 1,k)]y(l,k)为当前帧频域参考信号，y(l-1,k)
…
y(l-l 1,k)分别为当前帧之前的l-1帧，即第l-1，
…
，l-l 1帧的频域参考信号；wv(l,k)表示参考信号从扬声器传播到麦克风的理想频域回声路径传递函数，即wv(l,k)=[w
l
(0,k),w
l
(1,k)...w
l
(l-1,k)]，w
l
(0,k)，w
l
(1,k)，
…
，w
l
(l-1,k)分别为作用于第l帧第k个频点参考信号向量y(l,k)，y(l-1,k)，
…
，y(l-l 1,k)的传递函数系数。
[0020]
l为频域参考信号向量和理想频域回声路径传递函数的长度。
[0021]
s3. 将当前帧和过去l-1帧的频域参考信号组成的频域参考信号向量yv(l,k)=[y(l,k),y(1-1,k)...y(l-l 1,k)]存储到缓存器中。
[0022]
根据频域麦克风信号和频域参考信号向量，可对传播通过回声路径的参考信号进行线性预测，向量长度l取任意正整数值，通常取值l∈[1,4]。
[0023]
l值越大，回声消除效果越好，但内存和算力也会随之增大，所需要的缓存器也越大，将频域参考信号向量存储到缓存器中，便于进行后续s4-s5步骤的数据调用和计算，l值根据硬件条件和回声消除效果要求综合考虑决定。
[0024]
所述的频域参考信号向量是由当前帧频域参考信号y(l,k)和过去l-1帧的频域参考信号y(l-1,k)
…
y(l-l 1,k)所组成的。
[0025]
s4.采用线性预测的方法，根据频域参考信号向量和频域麦克风信号来估计扬声器到麦克风的回声路径传递函数。
[0026]
l=1时，第一帧参考信号向量中过去帧的频域参考信号取设定的初值，例如均为零，可取yv(l,k)=[y(l,k),0,...0]。
[0027]
利用频域卡尔曼滤波器作为第一滤波器对频域参考信号向量滤波，输出回声估计信号
[0028]
；式中上标h表示共轭转置，为第l帧第k个频点的回声估计信号，yv(l,k)为第l帧第k个频点的频域参考信号向量,gv(l,k)为第l帧第k个频点的估计回声路径传递函数；
[0029]
l=1时，第一帧估计回声路径传递函数可取任意不为零的初值，例如可取gv(l,k)=[0.1,0.1,...0.1]。
[0030]
由所述频域麦克风信号减去回声估计信号得到第一残余信号e1(l,k)，
[0031][0032]
估计回声路径传递函数gv(l,k)的更新公式如下
[0033]gv
(l,k)=gv(l-1,k) kv(l,k)e1(l,k)式中kv(l,k)为第l帧第k个频点的卡尔曼增益，
[0034]
其中，卡尔曼增益的计算方式如下
[0035][0036]
式中上标t表示转置，p-(l,k)为第l帧第k个频点的先验误差协方差矩阵，p(l,k)为第l帧第k个频点的后验误差协方差矩阵，f为状态转换矩阵，q为过程噪声协方差矩阵，r(l,k)为第l帧第k个频点的测量噪声协方差矩阵，i为单位矩阵。
[0037]
状态转移矩阵f和过程噪声协方差矩阵q为可以取得优良效果的经验值，例如可取f=i，q=i，初始后验误差协方差矩阵可取任意不为零的对角矩阵。
[0038]
一个具体实施例中，可采用优化的测量噪声协方差矩阵，其计算方式如下
[0039][0040]
式中e{
·
}表示求数学期望，α为第一平滑参数，本具体实施方式中取α=0.89，r
e,e
(l,k)为第l帧第k个频点的第一残余信号的自相关矩阵。
[0041]
若测量噪声协方差矩阵的值过大，滤波响应速度会变慢，若测量噪声协方差矩阵值太小，可能引起滤波发散。采用上述优化的测量噪声协方差矩阵由当前帧第一残余信号自相关矩阵和上一帧测量噪声协方差矩阵组成，进一步平衡了回声抑制能力和滤波响应速
度。
[0042]
s5. 利用第二滤波器对第一残余信号滤波，输出第二残余信号。
[0043]
首先根据步骤s4中所述的测量噪声协方差矩阵和回声估计信号确定后验信噪比
[0044]
γ(l,k)=r(l,k)/t(l,k)
[0045]
为第l帧第k个频点的回声估计信号；
[0046]
式中t(l,k)为第l帧第k个频点的经过平滑处理后的第一回声估计能量，β为第二平滑参数，本具体实施方式中取β=0.99，为第l帧第k个频点的回声估计能量，其中初始平滑回声估计能量取设定的初值，例如可取t(l,k)=0.0001。
[0047]
接着由以下方程组，根据平滑回声估计能量t(l,k)和后验信噪比γ（l,k），进一步循环计算得到第二残余信号能量r
e2,e2
和先验信噪比；
[0048]
先验信噪比；第二残余信号e2(l-1,k)=g(l-1,k)e1(l,k)第二残余信号能量；第二残余信号平滑后能量；
[0049]
式中max{}表示取最大值，g(l,k)为第二滤波器参数，s(l,k)为第l帧第k个频点的经过平滑处理后的第二残余信号平滑后能量，e2(l-1,k)为第l-1帧第k个频点的第二残余信号，λ为权重因子，c为经验常数，η为第三平滑参数，具体实施方式中取λ=0.7，c=0.08，η=0.85
[0050]
初始状态下，l=1,第二滤波器参数取不为零的初值，例如可取g(l,k)=0.01，第二残余信号平滑后能量的初始值取设定的初值，例如可取s(l,k)=0.0001。计算先验信噪比，然后根据先验信噪比更新第二滤波器参数；
[0051][0052]
式中g(l,k)为第l帧第k个频点的第二滤波器参数，更新后计算第二残余信号能量r
e2,e2
和第二残余信号平滑后能量s(l,k)，再重复循环计算；
[0053]
利用第二滤波器对第一残余信号滤波，输出第二残余信号e2(l,k)=g(l,k)e1(l,k)。
[0054]
最后通过短时傅里叶逆变换将经过回声消除处理后的频域信号转换成时域信号。
[0055]
本发明相对现有技术，采用基于卡尔曼滤波的回声消除方法，结合优化设计的时变测量噪声协方差矩阵，由当前帧第一残余信号自相关矩阵和上一帧测量噪声协方差矩阵
组成，第一残余信号自相关矩阵通过引入当前残余回声信号来反应对测量值的信任度，使建模更加精准，提高了回声抑制能力；采用优选的测量噪声协方差矩阵可避免更新的测量噪声协方差矩阵的值变化过大，导致滤波发散，进一步平衡了回声抑制能力和滤波响应速度。此外，还联合第二滤波器，对残留回声做出进一步的抑制。
[0056]
如图2所示，上方为本发明具体实施例中的输入信号时域波形示意图，下方为采用现有技术对输入信号进行回声消除后输出的时域波形示意图；如图3所示给出具体实施例中，上方为输入信号时域波形示意图，下方为采用本发明所述方法对输入信号进行回声消除后输出的时域波形示意图。从图2和图3可见，经回声消除处理后，时域波形中存在的多个电压幅值较大部分为保留的语音信号，电压幅值较小部分为残留回声信号，显然，经本发明所述方法处理输出的语音信号与残留信号之间幅值差别更大，表明本发明所述方法的回声抑制能力更强。
[0057]
前文所述的为本发明的各个优选实施例，各个优选实施例中的优选实施方式如果不是明显自相矛盾或以某一优选实施方式为前提，各个优选实施方式都可以任意叠加组合使用，所述实施例以及实施例中的具体参数仅是为了清楚表述发明人的发明验证过程，并非用以限制本发明的专利保护范围，本发明的专利保护范围仍然以其权利要求书为准，凡是运用本发明的说明书及附图内容所作的等同结构变化，同理均应包含在本发明的保护范围内。