CN115588438B - 一种基于双线性分解的wls多通道语音去混响方法 - Google Patents

一种基于双线性分解的wls多通道语音去混响方法 Download PDF

Info

Publication number
CN115588438B
CN115588438B CN202211587593.7A CN202211587593A CN115588438B CN 115588438 B CN115588438 B CN 115588438B CN 202211587593 A CN202211587593 A CN 202211587593A CN 115588438 B CN115588438 B CN 115588438B
Authority
CN
China
Prior art keywords
sub
filter
formula
equation
signal
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.)
Active
Application number
CN202211587593.7A
Other languages
English (en)
Other versions
CN115588438A (zh
Inventor
曾小金
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Chipintelli Technology Co Ltd
Original Assignee
Chipintelli Technology Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Chipintelli Technology Co Ltd filed Critical Chipintelli Technology Co Ltd
Priority to CN202211587593.7A priority Critical patent/CN115588438B/zh
Publication of CN115588438A publication Critical patent/CN115588438A/zh
Application granted granted Critical
Publication of CN115588438B publication Critical patent/CN115588438B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

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
    • G10L21/0216—Noise filtering characterised by the method used for estimating noise
    • G10L21/0224—Processing in the time domain
    • 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
    • G10L19/265—Pre-filtering, e.g. high frequency emphasis prior to encoding
    • 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
    • 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
    • G10L21/0216—Noise filtering characterised by the method used for estimating noise
    • G10L2021/02161—Number of inputs available containing the signal or the noise to be suppressed
    • G10L2021/02166—Microphone arrays; Beamforming

Landscapes

  • Engineering & Computer Science (AREA)
  • Computational Linguistics (AREA)
  • Signal Processing (AREA)
  • Health & Medical Sciences (AREA)
  • Audiology, Speech & Language Pathology (AREA)
  • Human Computer Interaction (AREA)
  • Physics & Mathematics (AREA)
  • Acoustics & Sound (AREA)
  • Multimedia (AREA)
  • Quality & Reliability (AREA)
  • Circuit For Audible Band Transducer (AREA)

Abstract

本发明公开了一种基于双线性分解的WLS多通道语音去混响方法,包括如下步骤:S1.采集混响语音时域信号；S2.以混响语音时域信号为基础描述参考麦克风信号模；S3.计算期望信号；S4.进行双线性分解：利用克罗内克乘积将预测滤波器系数向量分解为两个子滤波器；S5.求解预测滤波器系数向量；S6.正则化滤波器系数向量；S7.计算出去除混响之后的期望信号。本发明利用加权最小二乘方法，将最大似然估计问题转换为最小二乘问题进行处理。同时，由于矩阵在求逆过程中的计算复杂度较高，为了进一步降低计算量，引入了克罗内克积分解，将一个长的滤波器分解为两个短的子滤波器进行处理，从而使得算法在具有良好去混响性能的同时还拥有较低的计算复杂度。

Description

一种基于双线性分解的WLS多通道语音去混响方法
技术领域
本发明属于智能语音技术领域，涉及语音技术，具体涉及一种基于双线性分解的WLS多通道语音去混响方法。
背景技术
依据各类标准，语音去混响技术可以被划分为多个类别。其中一种分类标准是按照语音接收端麦克风的数量，可以将去混响方法分为单通道去混响方法和多通道去混响方法。由于单通道去混响技术对硬件设备的要求低，只需要单个麦克风，拥有较低的成本，因此也赢得了广泛的关注。然而，单个麦克风在接收语音信号时，无法获得语音信号传播过程中的空间信息，相比于多个麦克风而言，能够利用的信息更少。在此基础上，多通道去混响技术应运而生。
加权预测误差(WPE)方法是迄今为止最有效的去混响算法之一，它是一种经典的基于多通道线性预测的语音去混响方法，对于提高语音质量非常有效。该方法基于时变高斯模型对语音信号进行建模，利用最大似然估计来估计线性预测滤波器，通过方差与滤波器系数两者交替迭代计算得出去混响之后的语音信号。加权最小二乘法（WLS）是对原模型进行加权，使之成为一个新的不存在异方差性的模型，WLS直接应用于语音去混响，算力需求较高，限制了在自适应滤波领域的应用。
发明内容
为将去混响问题引入自适应滤波领域，本发明公开了一种基于双线性分解的WLS多通道语音去混响方法。
本发明所述基于双线性分解的WLS多通道语音去混响方法,如图1所示，包括如下步骤：
S1.利用麦克风阵列采集语音信号并转化为数字信号形式的混响语音时域信号xm(n)，表示为：
xm(n)=hm(n)*s(n)+vm(n);
式中*表示线性卷积，s(n)表示纯净语音信号，hm(n)表示的是声源与第m个麦克风之间的房间脉冲响应且其长度为Lh，vm(n)表示背景噪声，m和n分别表示第m个麦克风和第n个离散时间；
S2.令麦克风阵列中第一个麦克风为参考麦克风，假设背景噪声为零，以混响语音时域信号xm(n)为基础，将参考麦克风信号模型以公式3描述：
Figure 170287DEST_PATH_IMAGE001
---公式3
公式3中上标T和*分别表示矩阵的转置和共轭运算，d1(k)表示参考麦克风在时频域的期望信号，g(k)为预测滤波器系数向量，其上标^、*分别表示待估计量和共轭运算,Xτ(k)是在第k个频率处，经过时间延迟τ帧后的输入信号构成的矩阵,下标τ表示输入信号经过了τ帧的线性延迟；
S3. 利用公式3 ，得到期望信号d1(k)等于参考麦克风接收到的信号减去估计出来的中晚期混响分量：
Figure 493559DEST_PATH_IMAGE002
--公式4；
S4. 进行双线性分解：利用克罗内克乘积将预测滤波器系数向量g(k)分解为两个子滤波器如下：
Figure 737459DEST_PATH_IMAGE003
---公式5
其中上标^表示待估计量，
Figure 255028DEST_PATH_IMAGE004
和
Figure 420430DEST_PATH_IMAGE005
表示的是长度分别为L1和L2的第一子滤波器和第二子滤波器，
Figure 550322DEST_PATH_IMAGE006
表示克罗内克乘积，
Figure 383149DEST_PATH_IMAGE004
和
Figure 71619DEST_PATH_IMAGE007
下标中的分量序号ps表示的是子滤波器系数向量中第ps个分量，分量序号ps=1,2,…,P；P是分解之后每个子滤波器所含滤波器系数的个数；
根据公式5，将公式4分解可得
Figure 989897DEST_PATH_IMAGE008
---公式7
其中当i=1时表示第一子滤波器构成的期望信号向量d1,1(k)，取2时则表示第二子滤波器构成的期望信号d1,2(k),x1(k)为参考麦克风信号模型，*表示共轭计算，子滤波器系数向量由公式8表示：
Figure 422015DEST_PATH_IMAGE009
--公式8
i取不同的值表示不同的子滤波器，i=1，即
Figure 404621DEST_PATH_IMAGE010
表示第一子滤波器系数向量，i=2即
Figure 263993DEST_PATH_IMAGE011
则是第二子滤波器系数向量，
Figure 669566DEST_PATH_IMAGE012
在i=1,2时分别表示由输入混响语音信号与第一子滤波器和第二子滤波器构成的双线性分解矩阵；
S5. 预测滤波器系数向量的求解，
分别求解公式7中的第一子滤波器和第二子滤波器，通过期望信号时频域代价函数的最小化来估计获得，
第一子滤波器系数向量
Figure 170955DEST_PATH_IMAGE013
可由下式推导得到
Figure 948680DEST_PATH_IMAGE014
--公式9
式中上标^表示待估计量，上标H表示共轭转置,表示当括号里的表达式F(f)取最小值时，F(f)中变量f的值，||2表示求解绝对值的平方，
Figure 447795DEST_PATH_IMAGE015
为输入混响语音信号与第一子滤波器构成的双线性分解矩阵的内部向量，d1,1(k)为第一子滤波器构成的期望信号向量；
加权矩阵
Figure 871823DEST_PATH_IMAGE016
，diag[]表示对角阵，其中的对角线元素a1(n,k)是一个N维列向量，N为总帧数，即a1(n,k)=[w1(1,k), w1(2,k),... w1(N,k)]T；式中w1(n,k),n=1,2...N是第一子滤波器在第n帧第k个频率处的加权系数；
令
Figure 849006DEST_PATH_IMAGE017
表示代价函数，上标H表示共轭转置，d1,1(k)表示第一子滤波器构成的期望信号向量；
为了求取代价函数的最小值，结合公式9对代价函数
Figure 776511DEST_PATH_IMAGE018
进行求导并令求导之后的值为零，解出第一子滤波器系数向量
Figure 7378DEST_PATH_IMAGE013
的表达式如下
Figure 918703DEST_PATH_IMAGE019
式中[]-1表示矩阵求逆运算，
令自相关矩阵
Figure 699577DEST_PATH_IMAGE020
，互相关向量
Figure 747167DEST_PATH_IMAGE021
，则第一子滤波器系数向量
Figure 322505DEST_PATH_IMAGE013
的表达式可以简化为
Figure 425853DEST_PATH_IMAGE022
---公式91；
同理，可以得到第二子滤波器系数向量
Figure 807155DEST_PATH_IMAGE023
---公式92；
S7.利用步骤S5中的公式91和公式92，求得两个子滤波器系数向量
Figure 178094DEST_PATH_IMAGE013
和
Figure 924333DEST_PATH_IMAGE011
；结合公式5和公式8, 计算得到滤波器系数向量
Figure 747932DEST_PATH_IMAGE024
；带入公式4便可以计算出去除混响之后的期望信号d1(k)。
优选的,S5步骤之后还包括S6步骤，具体为：
S6.对公式91和公式92进行正则化，长度为L1的第一子滤波器系数向量改写为
Figure 136189DEST_PATH_IMAGE025
---公式93
式中δ是正则化因子，同理，可以得到长度为L2的第二子滤波器系数向量
Figure 656906DEST_PATH_IMAGE011
的表达式如下：
Figure 839626DEST_PATH_IMAGE026
---公式94
公式94式自相关矩阵
Figure 150522DEST_PATH_IMAGE027
，
而互相关向量
Figure 76889DEST_PATH_IMAGE028
；
所述S7步骤中,以公式93和公式94分别替换公式91和公式92进行计算。
本发明利用加权最小二乘方法，将最大似然估计问题转换为最小二乘问题进行处理。同时，由于矩阵在求逆过程中的计算复杂度较高，为了进一步降低计算量，引入了Kronecker积分解，将一个长的滤波器分解为两个短的子滤波器进行处理，从而使得算法在具有良好去混响性能的同时还拥有较低的计算复杂度。
附图说明
图1为本发明所述去混响方法的一种具体实施方式示意图；
图2为本发明一个具体实施例中未处理前的原始混响语音对应的波形图和语谱图；
图3为对图2所示原始混响语音采用现有技术进行处理得到的波形图和语谱图；
图4为对图2所示原始混响语音采用本发明进行处理得到的波形图和语谱图；
图2至图4中的（a）部分为原始混响语音对应的波形图，（b）部分为原始混响语音对应的语谱图；
图5为现有技术和本发明的算力对比示意图。
具体实施方式
下面对本发明的具体实施方式作进一步的详细说明。
本发明所述去混响方法，包括如下步骤：
S1. 利用麦克风阵列采集声源发出的语音信号，经过数模转换之后处理为数字信号形式的混响语音时域信号xm(n),表示为：
xm(n)=hm(n)*s(n)+vm(n)；
式中*表示线性卷积，s(n)表示纯净语音信号，hm(n)表示的是声源与第m个麦克风之间的房间脉冲响应且长度为Lh，vm(n)表示背景噪声，m和n分别表示第m个麦克风和第n个离散时间；
S2.利用短时傅里叶变换，将混响语音分为不同的频带进行处理：
对麦克风阵列设置预测滤波器，将第一个麦克风当作参考麦克风，参考麦克风最终矩阵形式的信号模型：
Figure 687999DEST_PATH_IMAGE001
---公式3
公式3中上标T和*分别表示矩阵的转置和共轭运算，d1(k)表示参考麦克风在第k个频率处的期望信号，Xτ(k)表示麦克风阵列的全部M个麦克风的阵列输入矩阵，下标τ表示输入信号经过了τ帧的线性延迟，g(k)为预测滤波器系数向量；
S3. 利用公式3，可以得到参考麦克风在时频域的期望信号d1(k)等于参考麦克风接收到的信号x1(k)减去估计出来的中晚期混响分量，即：
Figure 103937DEST_PATH_IMAGE002
---公式4
S4.利用克罗内克乘积将预测滤波器系数向量g(k)分解为两个子滤波器如下：
子滤波器系数向量
Figure 731490DEST_PATH_IMAGE029
--公式8
i取不同的值1,2表示不同的子滤波器，i=1，即
Figure 258286DEST_PATH_IMAGE010
表示第一子滤波器系数向量，i=2即
Figure 192744DEST_PATH_IMAGE011
则是第二子滤波器系数向量，上标^表示待估计量，gi,1（k）表示第i子滤波器在第k个频率处的第一个分量，其余以此类推，上标T表示转置；
并将公式4改写成
Figure 717266DEST_PATH_IMAGE008
---公式7
其中当i=1时表示第一子滤波器构成的期望信号向量d1,1(k)，取2时则表示第二子滤波器构成的期望信号d1,2(k)，两个期望信号都是由N个元素构成的列向量，在第n帧和第k个频率处的元素表示为d1,i(n,k),n=1,2...N ，x1(k)为参考麦克风接收到的信号;
公式7中，由输入混响语音信号与子滤波器系数向量所构成的分解矩阵
Figure 298027DEST_PATH_IMAGE030
是一个由N个
Figure 566197DEST_PATH_IMAGE031
维的列向量构成，其具体形式如下
Figure 355162DEST_PATH_IMAGE032
对于
Figure 316164DEST_PATH_IMAGE033
中的每个元素向量
Figure 823369DEST_PATH_IMAGE034
又进一步由P个
Figure 691968DEST_PATH_IMAGE035
的列向量构成，n-τ表示信号经过了τ帧的线性延迟：
Figure 368062DEST_PATH_IMAGE036
上标T表示矩阵转置，而构成
Figure 968808DEST_PATH_IMAGE034
向量的内部每一个Li维的混响信号矩阵列向量
Figure 228888DEST_PATH_IMAGE037
，ps=1，2，…，P则是由以下等式计算得出：
第一列向量
Figure 838861DEST_PATH_IMAGE038
第二列向量
Figure 664734DEST_PATH_IMAGE039
式中上标H表示共轭转置运算符，上标^表示待估计量，g1,ps(k)和g2,ps(k)分别表示两个子滤波器系数向量的其中第ps个分量，
Figure 701960DEST_PATH_IMAGE040
和
Figure 183757DEST_PATH_IMAGE041
则表示的是维度分别为L1和L2单位矩阵，符号
Figure 119394DEST_PATH_IMAGE006
表示的是克罗内克乘积运算，x(n-τ,k)表示经过了τ帧线性延迟的输入的混响语音信号。
而x(n-τ,k)是由M个Lg维的列向量构成,M是麦克风总个数，Lg是原始预测滤波器总长度；
Figure 799774DEST_PATH_IMAGE042
对于向量x(n-τ,k)中的每个Lg维的列向量由Lg个元素构成，以第m（m=1,2,…,M）个麦克风处的信号向量xm(n-τ,k)为例，其形式如下所示：
xm(n-τ,k)=[ xm(n-τ,k), xm(n-τ-1,k),... xm(n-τ-Lg+1,k)]T；
S5.预测滤波器系数向量的求解
第一子滤波器系数向量
Figure 866956DEST_PATH_IMAGE013
可由下式推导得到
Figure 836049DEST_PATH_IMAGE044
--公式9
式中上标^表示待估计量，上标H表示共轭转置,表示当括号里的表达式F(f)取最小值时，F(f)中变量f的值，||2表示求解绝对值的平方，
Figure 23710DEST_PATH_IMAGE015
为输入混响语音信号与第一子滤波器构成的双线性分解矩阵的内部向量，d1,1(k)为第一子滤波器构成的期望信号向量；
加权矩阵
Figure 27438DEST_PATH_IMAGE016
，diag[]表示对角阵，其中的对角线元素a1(n,k)是一个N维列向量，N为总帧数，即a1(n,k)=[w1(1,k), w1(2,k),... w1(N,k)]T；式中w1(n,k),n=1,2...N是第一子滤波器在第n帧第k个频率处的加权系数；
令
Figure 203205DEST_PATH_IMAGE017
表示代价函数，上标H表示共轭转置，d1,1(k)表示第一子滤波器构成的期望信号向量；
为了求取代价函数的最小值，结合公式9对代价函数
Figure 190752DEST_PATH_IMAGE018
进行求导并令求导之后的值为零，解出第一子滤波器系数向量
Figure 946218DEST_PATH_IMAGE013
的表达式如下
Figure 538874DEST_PATH_IMAGE019
式中[]-1表示矩阵求逆运算，
令自相关矩阵
Figure 852918DEST_PATH_IMAGE020
，互相关向量
Figure 593341DEST_PATH_IMAGE021
，则第一子滤波器系数向量
Figure 886919DEST_PATH_IMAGE013
的表达式可以简化为
Figure 130819DEST_PATH_IMAGE022
---公式91；
同理，可以得到第二子滤波器系数向量
Figure 851650DEST_PATH_IMAGE023
---公式92；
S6.为了避免在矩阵求逆过程中出现不可逆，需要对自相关矩阵R1进行正则化，于是公式91改写为
Figure 282631DEST_PATH_IMAGE045
式中δ是正则化因子，通常是一个很小的正常数。
同理，公式92改写如下：
Figure 678103DEST_PATH_IMAGE026
式中自相关矩阵，而互相关向量。
S7. 利用步骤S5中的公式91和公式92，或改写的公式93和公式94 ，可以求得两个子滤波器系数向量
Figure 245351DEST_PATH_IMAGE013
和
Figure 137083DEST_PATH_IMAGE046
；结合公式5和公式8, 计算得到滤波器系数向量
Figure 55361DEST_PATH_IMAGE047
；带入公式4便可以计算出去除混响之后的期望信号d1(k)。
具体为，将两个子滤波器系数向量
Figure 487479DEST_PATH_IMAGE013
和
Figure 470085DEST_PATH_IMAGE046
带入公式8：
Figure 798298DEST_PATH_IMAGE048
可得
Figure 938293DEST_PATH_IMAGE049
；再带入公式5，
i=1时，
Figure 439681DEST_PATH_IMAGE050
分别等于
Figure 981521DEST_PATH_IMAGE051
，ps=1,2...P；
i=2时，
Figure 215056DEST_PATH_IMAGE052
分别等于
Figure 343811DEST_PATH_IMAGE005
，ps=1,2...P；
计算得到滤波器系数向量
Figure 117732DEST_PATH_IMAGE053
；带入公式4便可以计算出去除混响之后的期望信号d1(k)。
一个具体实施方式为：
S1. 利用麦克风阵列采集语音信号，通过模数转换可以得到数字化后的语音信号。
S2. 由于语音信号的短时平稳特性，采用短时傅里叶变换将语音信号从时域变换到时频域进行分析，利用分帧加窗的思想，将信号划分为不同的频带进行处理。
每一帧语音的长度设置为32 ms，帧移为8 ms，窗函数选择汉宁窗，预测滤波器的长度设置为Lg为24，采用512点的傅里叶变换，考虑到语音信号进行快速傅里叶变换之后的频谱具有共轭对称性，因此为了提高计算效率只选取一半的频点进行计算，即K=256。
S3. 设置预测延迟τ=2，对每一个频带分别进行去混响处理，最后再将所有处理后的数据进行合成，从而输出去除混响之后的语音信号。
由于造成语音质量下降的主要原因是中晚期混响成分，而早期混响成分不仅不会对语音质量带来显著的影响，相反其有助于提高语音的清晰度，因此在进行去混响处理时，往往将早期混响与直达语音成分整合在一起，于是需要去除的便是中晚期混响分量。因此在处理过程中需要进行一个延迟线性预测，保留早期混响。
S4.利用延迟预测滤波器对中晚期混响成分进行估计，然后将该成分从参考语音信号中减去，便可以得到最终去除混响之后的信号，为了方便，将第一个麦克风接收到的信号x1(k)作为参考信号。
S5.分别求解公式7中的第一子滤波器和第二子滤波器，通过期望信号时频域代价函数的最小化来估计获得，
第一子滤波器系数向量
Figure 779658DEST_PATH_IMAGE013
可由下式推导得到
Figure 184094DEST_PATH_IMAGE054
--公式9
式中上标^表示待估计量，上标H表示共轭转置,表示当括号里的表达式F(f)取最小值时，F(f)中变量f的值，||2表示求解绝对值的平方，
Figure 298681DEST_PATH_IMAGE015
为输入混响语音信号与第一子滤波器构成的双线性分解矩阵的内部向量，d1,1(k)为第一子滤波器构成的期望信号向量；
加权矩阵
Figure 876293DEST_PATH_IMAGE016
，diag[]表示对角阵，其中的对角线元素a1(n,k)是一个N维列向量，N为总帧数，即a1(n,k)=[w1(1,k), w1(2,k),... w1(N,k)]T；式中w1(1,k),n=1,2...N是第一子滤波器在第n帧第k个频率处的加权系数；
令
Figure 891260DEST_PATH_IMAGE017
表示代价函数，上标H表示共轭转置，d1,1(k)表示第一子滤波器构成的期望信号向量；
为了求取代价函数的最小值，结合公式9对代价函数
Figure 732177DEST_PATH_IMAGE018
进行求导并令求导之后的值为零，解出第一子滤波器系数向量
Figure 68480DEST_PATH_IMAGE013
的表达式如下
Figure 653045DEST_PATH_IMAGE019
式中[]-1表示矩阵求逆运算，
令自相关矩阵
Figure 883038DEST_PATH_IMAGE020
，互相关向量
Figure 130742DEST_PATH_IMAGE021
，则第一子滤波器系数向量
Figure 954342DEST_PATH_IMAGE013
的表达式可以简化为
Figure 139336DEST_PATH_IMAGE022
---公式91；
同理，可以得到第二子滤波器系数向量
Figure 364781DEST_PATH_IMAGE023
---公式92。
S6. 得到了两个子滤波器系数向量之后，便能计算出分解之前的滤波器系数
Figure 281921DEST_PATH_IMAGE055
，于是可以得到最终的期望信号：
Figure 858396DEST_PATH_IMAGE056
。
S7. 将得到的期望信号进行短时傅里叶逆变换，从而将每一个去混响之后的子带信号进行合成，转换成时域信号进行输出。
如图2至图5所示给出了本发明一个具体实施例。图2是一段原始混响语音对应的波形图和语谱图，而图3和图4则分别表示的是在混响时间为300 ms的环境中，麦克风个数M为8,滤波器长度Lg为24时，去混响前后的波形图和语谱图对比。对比图2至图4可以发现，混响被去除，且图3和图4显示现有技术和本发明对原始混响语音去混响效果相当。
图5则是双线性分解前后两种算法的计算复杂度对比情况。图中横坐标P表示的是分解之后每个子滤波器所含滤波器系数的个数，纵坐标则是计算复杂度。图5种的计算复杂度是指算法中所涉及的实数乘法次数。定义两个实数相乘涉及一次乘法，一个实数与一个复数之间做乘积运算涉及2次乘法，而两个复数之间做乘积则需要4次乘法运算。纵坐标所示数值为以指数函数对乘法次数进行表示，均为无量纲单位。图例中的WLS表示现有技术的算法，D-WLS表示的是经过本发明所述双线性分解之后的WLS算法，括号中的内容指的是两个子滤波器的长度。以M=8（32-6）为例，其含义为麦克风个数为8，两个子滤波器长度分别为32和6。
结合图3和图4分析可知，分解之后的WLS算法不仅拥有良好的去混响效果，同时其计算复杂度也有显著降低。
前文所述的为本发明的各个优选实施例，各个优选实施例中的优选实施方式如果不是明显自相矛盾或以某一优选实施方式为前提，各个优选实施方式都可以任意叠加组合使用，所述实施例以及实施例中的具体参数仅是为了清楚表述发明人的发明验证过程，并非用以限制本发明的专利保护范围，本发明的专利保护范围仍然以其权利要求书为准，凡是运用本发明的说明书及附图内容所作的等同结构变化，同理均应包含在本发明的保护范围内。

Claims (2)

1.基于双线性分解的WLS多通道语音去混响方法,其特征在于,包括如下步骤：
S1.利用麦克风阵列采集语音信号并转化为数字信号形式的混响语音时域信号xm(n)，表示为：
xm(n)=hm(n)*s(n)+vm(n)；
式中*表示线性卷积，s(n)表示纯净语音信号，hm(n)表示的是声源与第m个麦克风之间的房间脉冲响应且其长度为Lh，vm(n)表示背景噪声，m和n分别表示第m个麦克风和第n个离散时间；
S2.令麦克风阵列中第一个麦克风为参考麦克风，假设背景噪声为零，以混响语音时域信号xm(n)为基础，将参考麦克风信号模型以公式3描述：
Figure QLYQS_1
---公式3
公式3中上标T和*分别表示矩阵的转置和共轭运算，d1(k)表示参考麦克风在时频域的期望信号，g(k)为预测滤波器系数向量, 其上标^、*分别表示待估计量和共轭运算，Xτ(k)是在第k个频率处，经过时间延迟τ帧后的输入信号构成的矩阵,下标τ表示输入信号经过了τ帧的线性延迟；
S3. 利用公式3 ，得到期望信号d1(k)等于参考麦克风接收到的信号减去估计出来的中晚期混响分量：
Figure QLYQS_2
--公式4；
S4. 进行双线性分解：利用克罗内克乘积将预测滤波器系数向量g(k)分解为两个子滤波器如下：
Figure QLYQS_3
---公式5
其中上标^表示待估计量，
Figure QLYQS_4
和
Figure QLYQS_5
表示的是长度分别为L1和L2的第一子滤波器和第二子滤波器，
Figure QLYQS_6
表示克罗内克乘积，
Figure QLYQS_7
和
Figure QLYQS_8
下标中的分量序号ps表示的是子滤波器系数向量中第ps个分量，分量序号ps=1,2,…,P；P是分解之后每个子滤波器所含滤波器系数的个数；
根据公式5，将公式4分解可得
Figure QLYQS_9
---公式7
其中当i=1时表示第一子滤波器构成的期望信号向量d1,1(k)，取2时则表示第二子滤波器构成的期望信号d1,2(k),x1(k)为参考麦克风信号模型，*表示共轭计算，子滤波器系数向量由公式8表示：
Figure QLYQS_10
--公式8
i取不同的值表示不同的子滤波器，i=1，即
Figure QLYQS_11
表示第一子滤波器系数向量，i=2即
Figure QLYQS_12
则是第二子滤波器系数向量，
Figure QLYQS_13
在i=1,2时分别表示由输入混响语音信号与第一子滤波器和第二子滤波器构成的双线性分解矩阵；
S5. 预测滤波器系数向量的求解，
分别求解公式7中的第一子滤波器和第二子滤波器，通过期望信号时频域代价函数的最小化来估计获得，
第一子滤波器系数向量
Figure QLYQS_14
可由下式推导得到
Figure QLYQS_15
--公式9
式中上标^表示待估计量，上标H表示共轭转置,表示当括号里的表达式F(f)取最小值时，F(f)中变量f的值，||2表示求解绝对值的平方，
Figure QLYQS_16
为输入混响语音信号与第一子滤波器构成的双线性分解矩阵的内部向量，d1,1(k)为第一子滤波器构成的期望信号向量；
加权矩阵
Figure QLYQS_17
，diag[]表示对角阵，其中的对角线元素a1(n,k)是一个N维列向量，N为总帧数，即a1(n,k)=[w1(1,k), w1(2,k),... w1(N,k)]T；式中w1(n,k),n=1,2...N是第一子滤波器在第n帧第k个频率处的加权系数；
令
Figure QLYQS_18
表示代价函数，上标H表示共轭转置，d1,1(k)表示第一子滤波器构成的期望信号向量；
为了求取代价函数的最小值，结合公式9对代价函数
Figure QLYQS_19
进行求导并令求导之后的值为零，解出第一子滤波器系数向量
Figure QLYQS_20
的表达式如下
Figure QLYQS_21
式中[]-1表示矩阵求逆运算，
令自相关矩阵
Figure QLYQS_22
，互相关向量
Figure QLYQS_23
，则第一子滤波器系数向量
Figure QLYQS_24
的表达式简化为
Figure QLYQS_25
---公式91；
同理，得到第二子滤波器系数向量
Figure QLYQS_26
---公式92；
S7.利用步骤S5中的公式91和公式92，求得两个子滤波器系数向量
Figure QLYQS_27
和
Figure QLYQS_28
；结合公式5和公式8, 计算得到滤波器系数向量
Figure QLYQS_29
；带入公式4计算出去除混响之后的期望信号d1(k)。
2.如权利要求1所述基于双线性分解的WLS多通道语音去混响方法,其特征在于,S5步骤之后还包括S6步骤，具体为：
S6.对公式91和公式92进行正则化，长度为L1的第一子滤波器系数向量改写为
Figure QLYQS_30
---公式93
式中δ是正则化因子，同理，得到长度为L2的第二子滤波器系数向量
Figure QLYQS_31
的表达式如下：
Figure QLYQS_32
---公式94
公式94式自相关矩阵
Figure QLYQS_33
，
而互相关向量
Figure QLYQS_34
；
所述S7步骤中,以公式93和公式94分别替换公式91和公式92进行计算。
CN202211587593.7A 2022-12-12 2022-12-12 一种基于双线性分解的wls多通道语音去混响方法 Active CN115588438B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211587593.7A CN115588438B (zh) 2022-12-12 2022-12-12 一种基于双线性分解的wls多通道语音去混响方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211587593.7A CN115588438B (zh) 2022-12-12 2022-12-12 一种基于双线性分解的wls多通道语音去混响方法

Publications (2)

Publication Number Publication Date
CN115588438A CN115588438A (zh) 2023-01-10
CN115588438B true CN115588438B (zh) 2023-03-10

Family

ID=84783065

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211587593.7A Active CN115588438B (zh) 2022-12-12 2022-12-12 一种基于双线性分解的wls多通道语音去混响方法

Country Status (1)

Country Link
CN (1) CN115588438B (zh)

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP2437257A1 (en) * 2006-10-16 2012-04-04 Fraunhofer-Gesellschaft zur Förderung der angewandten Forschung e.V. Apparatus and method for multi-channel parameter transformation
CN108141691A (zh) * 2015-10-14 2018-06-08 华为技术有限公司 自适应混响消除***

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9692939B2 (en) * 2013-05-29 2017-06-27 Yeda Research And Development Co. Ltd. Device, system, and method of blind deblurring and blind super-resolution utilizing internal patch recurrence
WO2021229197A1 (en) * 2020-05-12 2021-11-18 Queen Mary University Of London Time-varying and nonlinear audio processing using deep neural networks
CN115019818A (zh) * 2022-05-26 2022-09-06 重庆金美通信有限责任公司 基于渐进式混响感知网络的语音去混响方法、装置及设备

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP2437257A1 (en) * 2006-10-16 2012-04-04 Fraunhofer-Gesellschaft zur Förderung der angewandten Forschung e.V. Apparatus and method for multi-channel parameter transformation
CN108141691A (zh) * 2015-10-14 2018-06-08 华为技术有限公司 自适应混响消除***

Also Published As

Publication number Publication date
CN115588438A (zh) 2023-01-10

Similar Documents

Publication Publication Date Title
CN108172231B (zh) 一种基于卡尔曼滤波的去混响方法及***
Gannot et al. Subspace methods for multimicrophone speech dereverberation
JP5124014B2 (ja) 信号強調装置、その方法、プログラム及び記録媒体
JP4567655B2 (ja) 音声信号中の暗騒音を抑制する方法及び装置、並びに反響消去を伴う対応する装置
US8467538B2 (en) Dereverberation apparatus, dereverberation method, dereverberation program, and recording medium
JP5435204B2 (ja) 雑音抑圧の方法、装置、及びプログラム
JP2008203879A (ja) 雑音抑圧の方法及び装置並びにコンピュータプログラム
JP2007526511A (ja) 周波数領域で多重経路多チャネル混合信号のブラインド分離のための方法及びその装置
JP6987075B2 (ja) オーディオ源分離
CN110111802B (zh) 基于卡尔曼滤波的自适应去混响方法
US9078077B2 (en) Estimation of synthetic audio prototypes with frequency-based input signal decomposition
JP7348812B2 (ja) 雑音抑制装置、雑音抑制方法及び音声入力機器
CN115588438B (zh) 一种基于双线性分解的wls多通道语音去混响方法
EP2630812B1 (en) Estimation of synthetic audio prototypes
EP2730026B1 (en) Low-delay filtering
Shamsa et al. Noise reduction using multi-channel FIR warped Wiener filter
CN109074811B (zh) 音频源分离
CN114220453B (zh) 基于频域卷积传递函数的多通道非负矩阵分解方法及***
KR100863184B1 (ko) 간섭 및 반향신호 제거를 위한 다단계 암묵 디콘볼루션방법
WO2022215199A1 (ja) 情報処理装置、出力方法、及び出力プログラム
JP4714892B2 (ja) 耐高残響ブラインド信号分離装置及び方法
Masuyama et al. Causal and Relaxed-Distortionless Response Beamforming for Online Target Source Extraction
CN117877503A (zh) 一种基于自适应多通道线性预测的低复杂度语音去混响方法
Gaubitch et al. Subband method for multichannel least squares equalization of room transfer functions
CN114550739A (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
GR01 Patent grant
GR01 Patent grant