回响去除装置、回响去除方法、回响去除程序、以及记录介质

摘要

在声源模型存储单元存储将从声源发出的音频信号作为概率密度函数来表示的声源模型。接收音频信号而得到的观测信号被变换为与多个频域的各自对应的按频率观测信号。然后，基于用于表示音频信号、观测信号以及回响去除滤波器的关系的回响模型和声源模型，利用各按频率观测信号来估计与各频带对应的回响去除滤波器。对各按频率观测信号应用各回响去除滤波器，从而求出与各频带对应的按频率目标信号，并将其合并。

说明书

技术领域

本发明涉及从所观测的观测信号去除回响信号的回响去除装置、回响去除方法、回响去除程序、以及记录介质。

背景技术

在以下说明中，将从声源发出的信号设为音频(audio)信号，有回响的房间发出音频信号，并通过多个接收部件(例如麦克风)接收而得到的信号设为观测信号。观测信号是对音频信号叠加了回响信号的信号。因此，从观测信号提取原来的音频信号的性质变得困难，且声音的清晰度也降低。相对于此，回响去除处理从观测信号除去所叠加的回响信号，从而容易提取音频信号原来的性质，并能够恢复声音的清晰度。这是有关通过作为其他各种音频信号处理系统的要素技术来使用，从而提高该系统整体的性能的技术。将回响去除处理作为要素技术能够有助于提高性能的音频信号处理系统例如有以下系统。

(1)将回响信号去除作为前处理来使用的声音识别系统

(2)通过回响信号去除提高声音的清晰度的电视会议系统等通信系统

(3)通过去除在讲演的录音中所包含的回响信号，从而提高录音的声音的清晰度的再现系统

(4)通过去除回响信号，从而提高听取容易度的助听器

(5)响应于人发出的声音，从而对机器传递命令的机器控制接口、以及机械和人之间的对话装置

(6)接收音频内容时包括回响信号，改善所接收的音频信号的音质的后期制作(postproduction)系统

(7)去除或附加音乐内容的回响信号，从而进行音乐内容的音频控制的音频效果器(effecter)。

图1表示以往的回响去除装置100的功能结构例(以下，称为“现有技术1”)。回响去除装置100包括估计单元104、去除单元106、声源模型存储单元108。在声源模型存储单元108，通过有限状态机使不包括回响信号的音频信号的短时间区间的波形模型化，并存储通过信号的自相关函数表现各状态的波形的特征的声源模型。此外，基于在时域对观测信号应用了回响去除滤波器的运算和上述声源模型，定义表现了从观测信号去除回响的信号(理想的目标信号)的似然性(likelihood)的最佳化函数。该最佳化函数涉及为将回响去除滤波器系数和声源模型的状态时系列作为参数，通过提供更加合适的滤波器系数和状态时系列，取得更大的值的函数。

以下的说明中，将输入的时域的观测信号设为x_t⁽¹⁾、…、x_t^(q)、…、x_t^(Q)。其中，x的右下方的下标“t”表示离散时刻的指示符，右上的q(q＝1、…、Q)表示接收部件(例如“麦克风”)的指示符。以下，将指示符为q的麦克风表示为第q信道的麦克风。这在下面也相同。

若输入观测信号x_t^(q)，则估计单元104使用该观测信号x_t^(q)和上述最佳化函数，估计回响去除滤波器。具体地说，估计单元104通过求出使上述最佳化函数的值最大的参数，从而估计回响去除滤波器。去除单元106通过对估计的回响去除滤波器卷积观测信号，从而输出从观测信号去除了回响信号的信号。该信号被称为目标信号。

图2表示以往的回响去除装置200的功能结构例(以下，称为“现有技术2”)。回响去除装置200包括用于将观测信号分割为U个频带的分割单元202、每个该频带的存储单元204_u(u＝0、…、U-1)、每个该频带的去除单元206_u、合并单元208。

分割单元202通过对观测信号进行子带分割，从而求出对U个频带的每一个分割的子带信号。该被分割的子带信号是时域的信号。此外，在子带分割时，有时进行下采样(downsample)(模型的间隔剔除)。在以下的说明中，将子带信号设为x’_n，u^(q)。其中，将n设为下采样后的采样的指示符，将u设为频带的指示符(u＝0、…、U-1)。以下，说明通过第q信道的麦克风接收的观测信号x_t^(q)的第u个频带的子带信号x’_n，u^(q)。

如上所述，在U个频带的每一个，设置去除单元206_u(u＝0、…、U-1)、存储单元204_u。在存储单元204_u存储回响去除滤波器。回响去除滤波器基于最小二乘法，预先决定回响去除滤波器的系数，使得利用从预先测定的声源至各麦克风的室内传达函数，在依次应用了该室内传达函数、分割单元202的子带分割处理、去除单元206_u的回响去除处理、合并单元208的合并处理的各处理时得到的系统整体的输入输出函数尽量成为单位脉冲(impulse)函数。

去除单元206_u通过对子带信号x’_n，u^(q)叠加回响去除滤波器，从子带信号去除回响信号。将从子带信号去除了回响信号的每个频带的子带信号设为按频率目标信号s~_n，u。然后，合并单元208合并按频率目标信号s_n，u~(u＝0、…、U-1)，从而求出目标信号s_t~。

回响去除装置100、200的细节记载在非专利文献1、2、3中。

(非专利文献1)T.Nakatani，B.H.Juang，T.Yoshioka，K.Kinoshita，M.Delcroix，and M.Miyoshi，“Study on speech dereverberation with autocorrelation codebook，”Proc.IEEE International Conference on Acoustics，Speech，and Signal Processing

(ICASSP-2007)，vol.I，pp.193-196，April 2007.

(非专利文献2)T.Nakatani，B.H.Juang，T.Yoshioka，K.Kinoshita，M.Miyoshi，“Importance of energy and spectral features in Gaussian source model for speech dereverberation，”WASPAA-2007，2007

(非专利文献3)N.D.Gaubitch，M.R.P.Thomas，P.A.Naylor，“Subband Method for Multichannel Least Squares Equalization of Room Transfer Functions，”Proc，IEEE Workshop on Applications of Signal Processing to Audio and Acoustics(WASPAA-2007)，pp.14-17，2007.

发明内容

在上述现有技术1的回响去除装置100中，为了通过最佳的方法来利用音频信号的时变性特征，需要计算非常大尺寸的协方差矩阵，以进行上述最佳化函数的值的最大化。因此，最佳化函数的值的最大化需要庞大的计算时间。说明协方差矩阵的尺寸变大的理由。有关在现有技术1中处理的观测信号的协方差矩阵H(r)通过以下的式(1)来表示。

$H\left(r\right)=\underset{t}{\Sigma }{X}_{t-1}^T{r}_t^{-1}{X}_{t-1}---\left(1\right)$

另外，在以下的说明中，协方差矩阵H(r)表示有关在现有技术1中处理的观测信号的协方差矩阵。其中，若设通过2个麦克风接收1个音频信号，则X_t-1＝[x^-_t-1⁽¹⁾，…，x^-_t-K⁽¹⁾，x^-_t-1⁽²⁾，…，x^-_t-K⁽²⁾]，x^-_t⁽¹⁾是由x_t⁽¹⁾的长度N的短时间帧构成的列向量x^-_t-1⁽¹⁾＝[x_t⁽¹⁾，x_t+1⁽¹⁾，…，x_t+N-1⁽¹⁾]^T，且设x_t⁽¹⁾、x_t⁽²⁾分别设为通过第1信道、第2信道的麦克风接收的观测信号。T表示矩阵、向量的转置。K是预测滤波器(估计的回响去除滤波器)的长度。此外，r_t是有关由音频信号的短时间帧构成的列向量s^-_t＝[s_t，s_t+1，…，s_t+N-1]^T的协方差矩阵r_t＝E{s^-_ts^-_t^T}，E{·}表示期待值函数。一般由于r_t不是已知的，因此基于在声源模型存储单元108中存储的声源模型，由估计单元104求出的估计值来代替。

一般，预测滤波器长度K理论上至少必须与室内脉冲响应长度相同。从而，协方差矩阵H(r)的尺寸变得非常大。另一方面，若假设音频信号为恒定信号，则由于上述的协方差矩阵可近似为自相关矩阵，因此能够使用快速傅立叶变换等快速的计算方法，但若对声音信号等时变信号使用该假设，则存在回响去除的计算精度变差的问题。这样，为了在回响去除装置100中进行高精度的回响去除需要非常庞大的计算时间，此外为了快速地进行回响去除，存在音频信号是时变信号的情况下回响去除的精度变差的问题。

此外，在上述现有技术2的回响去除装置200中，必须预先估计回响去除滤波器(室内传达函数的反滤波器)，需要预先求出室内传达函数。而且，利用室内传达函数的反滤波器进行回响去除的处理方法，若对室内传达函数的误差非常敏感而且在某种程度以上的误差包含在室内传达函数内的情况下，则存在通过回响去除处理反而音频信号的失真增大的问题。而且，室内传达函数对声源的位置和室温的变化敏感，在声源的位置和室温事先无法正确确定的情况下，不能决定正确的室内传达函数。这样，在回响去除装置200中，需要预先准备高精度的室内传达函数，而且，在一个条件下求出的室内传达函数只能在非常有限的条件下用于回响去除。

因此，本发明如下进行回响去除。在存储单元存储用于将音频信号作为概率密度函数来表现的声源模型。接收音频信号而获得的观测信号被变换为与多个频带的各自对应的按频率观测信号。然后，基于用于表示各频带中的音频信号、观测信号、回响去除滤波器之间关系的回响模型和声源模型，利用各按频率观测信号估计与各频带对应的回响去除滤波器。对各按频率观测信号应用各回响去除滤波器，从而求出与各频带对应的按频率目标信号，并合并各按频率目标信号。

附图说明

图1是表示了现有技术1的回响去除装置的功能结构例的方框图。

图2是表示了现有技术2的回响去除或者的功能结构例的方框图。

图3是表示了实施例1的回响去除装置的功能结构例的方框图。

图4是表示了实施例1的回响去除装置的主要处理的流程图。

图5是表示了实施例2的回响去除装置的功能结构例的方框图。

图6是表示了实施例2的回响去除装置的主要处理的流程图。

图7是表示了实施例3的回响去除装置的功能结构例的方框图。

图8是表示了实施例4的回响去除装置的功能结构例的方框图。

图9是表示实验结果的图。

图10A是表示利用一个麦克风基于实施例4证实了回响去除的效果的实验中的观测信号的声谱的图。

图10B是表示利用一个麦克风基于实施例4证实了回响去除的效果的实验结果的声谱的图。

具体实施方式

以下，表示用于实施发明的最佳方式。另外，对具有相同功能的结构单元和进行相同处理的过程附加相同标号，并省略重复说明。

实施例1

图3表示实施例1的回响去除装置300的方框图，图4表示回响去除装置300的主要处理的流程。如图3所示，实施例1的回响去除装置300包括对U个频带的每一个分割观测信号的分割单元302、声源模型存储单元304、每个频带的估计单元306_u(u＝0、…、U-1)、每个频带的去除单元308_u、合并单元310。

分割单元302一边对每个频带分割观测信号，一边减少观测信号的采样数，从而输出按频率观测信号。该实施例1的分割单元302对观测信号时间偏移并应用短时间分析窗，并变换为频域从而分割给每个频带。

在声源模型存储单元304存储对每个频带表现按频率观测信号的特征的声源模型。

对每个频带设置估计单元306_u，估计单元306_u基于与声源模型相关联定义的观测信号的最佳化函数，从按频率观测信号估计回响去除滤波器。

对每个频带设置去除单元308_u，利用按频率观测信号和回响去除滤波器，对每个频带求出按频率目标信号。该实施例1的去除单元308_u通过将回响去除滤波器卷积在所述按频率观测信号，从而求出按频率目标信号。

合并单元310通过合并按频率目标信号，从而输出后述的目标信号。该实施例1的合并单元310合并按频率目标信号，并变换为将全部频带汇总为一个的时域信号，从而输出后述的目标信号。

首先，说明音频信号s_t和观测信号x_t^(q)的关系。假设从声源至各麦克风的室内传递函数不具有共通零点，将最接近声源的麦克风设为q＝1(第1信道的麦克风)。音频信号和观测信号的关系可以如下式(11)那样表示。此外，细节记载在“M.Miyoshi，“Estimating AR parameter-sets for linear-recurrent signals in convolutive mixtures，”Proc.ICA-2003，pp.585-589，2003.”。

$x_t^{\left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =1}{\overset{K}{\Sigma }}{c}_{\tau }^{\left(q\right)}{x}_{t-\tau }^{\left(q\right)}+h_0^1{s}_t---\left(11\right)$

h_o⁽¹⁾是从声源至q＝1的麦克风的室内脉冲响应的第一抽头(tap)值，c_t^(q)称为预测系数，是在估计单元306_u估计的回响去除滤波器的系数，τ是离散时间的指示符，如上述那样，K是预测滤波器长度(在现有技术1中估计的回响去除滤波器的尺寸)。

其中，若忽略音频信号的增益，则右边的第2项h₀⁽¹⁾s_t是音频信号s_t的固定倍数的信号，因此能够将该信号看做应估计的音频信号s_t。由此，式(11)可以改写成以下式(12)。

$x_t^{\left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =1}{\overset{K}{\Sigma }}{c}_{\tau }^{\left(q\right)}{x}_{t-\tau }^{\left(q\right)}+s_t---\left(12\right)$

在式(12)中，当前的观测信号x_t^(q)从过去的观测信号的时间序列x_t-τ^(q)预测，音频信号s_t看做预测的残差信号。另外，作为式(12)的前提，假设第1信道的麦克风(q＝1)最接近声源，但在不满足该条件的情况下，也同样可以利用式(12)，表示观测信号和音频信号的关系。即，通过对第1信道的麦克风(q＝1)以外的麦克风的观测信号导入充分的延迟，从而能够假设来自声源的声音最初到达的麦克风设为第1信道的麦克风(q＝1)，能够将第1信道的麦克风作为最接近声源的麦克风来处理。因此，例如，若将导入麦克风q的延迟时间设为d^(q)抽头，则保留式(12)的表示，考虑为对q＝1以外的预测系数{c₁^(q)，c₂^(q)，…，c_K^(q)}的第一个开始对d^(q)抽头代入固定值0，从而能够将观测信号与音频信号之间的关系表示为与上述(12)相同。

若对分割单元302输入观测信号x_t^(q)，则一边对每个频带分割观测信号，一边减少观测信号的采样数，从而输出按频率观测信号(步骤S2)。然后，实施例1的分割单元302对观测信号时间偏移并应用短时间分析窗，并变换为频域，从而对每个频带分割。例如，分割单元302进行短时间傅立叶变换。以下，设为分割单元302进行短时间内傅立叶变换，从而具体说明。

接着，将上述式(12)一般化，考虑以下的式(12’)。

$x_t^{\left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =d}{\overset{K}{\Sigma }}{c}_{\tau }^{\left(q\right)}{x}_{t-\tau }^{\left(q\right)}+{\tilde{s}}_t---\left({12}^{\prime }\right)$

这里，d是对用于预测当前的观测信号的过去的观测信号导入延迟的常数。在设d＝1时，与式(12)相同。另一方面，在设d＞1时，上述式(12’)不能严格表示观测信号和音频信号的关系。这是因为，从当前的时刻t追溯到过去，由来于d抽头的音频信号的信号不包含在上述式(12’)的右边的过去的信号序列中，因此由来于该时间区间的音频信号而包含在当前的观测信号中的回响信号不能通过过去的观测信号的线性结合来表示。“由来于该时间区间的音频信号而包含在当前的观测信号中的回响信号”对应于与室内脉冲响应的最初的d抽头相当的初始反射音。从而，假设在上述式(12’)中，在残差信号中除了音频信号以外包含该初始反射音。为了明确该事实将残差信号记述为s_t~。另外，在该说明书中，记号A_α~表示在记号A的正上方附加了记号～的组合字符。

<对于频率信号的卷积运算>

接着，说明对频域信号进行相当于上述式(12’)的右边第1项中包含的时域中的卷积的运算的计算方法。首先，将在时域对某音频信号x_t卷积了滤波器长度为K的某回响去除滤波器c_t而得到的信号设为y_t。在时刻t0开始的y_t通过窗函数以时间窗切出了短时间帧的信号，能够通过z变换区域可以表示为式(13)那样。

W_N(y(z)z^t0)＝W_N(c(z)·x(z)z^t0)    (13)

其中，设y(z)＝c(z)·x(z)，“·”表示卷积，W()设为相当于时域中的长度N的窗函数的函数。W_N(c(z))从c(z)中的-N+1次获取0次的项，并与窗的形状成比例地变更各系数，并去除窗外的项。z^t0是将在时刻t0开始的短时间帧移动至窗函数中的时间偏移操作。

而且，将从时刻t中的滤波器系数c_t切出长度M的帧表示为设W_M^R()表示长度M的短时间分析窗(方形窗)。则，明显成为c(z)＝∑_τc_τM，M(z)。上述式(13)能够改写为如下。

$W_N\left(y_{t0,N}\left(z\right)\right)=W_N\left(\underset{\tau =0}{\overset{K_R}{\Sigma }}{c}_{\mathrm{\tau M},M}\left(z\right)z^{-\mathrm{\tau M}}x\left(z\right)z^{t0}\right)---\left(14\right)$

$=\underset{\tau =0}{\overset{K_R}{\Sigma }}{W}_N{c}_{\mathrm{\tau M},M}\left(z\right)x\left(z\right)z^{t0-\mathrm{\tau M}}---\left(15\right)$

$=\underset{\tau =0}{\overset{K_R}{\Sigma }}{W}_N\left(c_{\mathrm{\tau M},M}\left(z\right)x_{t0-M+1-\mathrm{\tau M},M+N-1}\left(z\right)z^{M-1}\right)---\left(16\right)$

其中，式(14)中的∑_τc_τM，M(z)z^-τM对应于c(z)(参照式(13))，式(16)中的x_{t0-M+1-τM，M+N-1}(z)对应于x(z)(参照式(13))。

此外，设K_R＝<K/M>，<K/M>表示K/M以上的最小的整数。K_R成为在估计单元306_u估计的回响去除滤波器的滤波器长度(抽头数)。在式(15)中，在窗函数的自变量中包含的项中，通过去除窗外的项，导出式(16)。

式(16)中的“c_τM，M(z)x_{t0-M+1-τM，M+N-1}(z)”表示在时域中的从滤波器系数c_τ的第τM抽头切出了长度M的帧、在时域的观测信号x_t的时刻t0-M+1-τM时间点开始切出了长度M的帧、在z域将它们相乘。z域中的乘法相当于卷积运算，因此表示上述观测信号x_t和滤波器系数c_t的各帧的时域中的卷积运算。此外，c_τM，M(z)的帧长度是M，x_{t0-M+1-τM，M+N-1}(z)的帧长度是M+N-1，因此在短时间傅立叶变换的点数(频带的数)U为U≥2M+N-2时，时域的卷积通过短时间傅立叶变换域的积严格表示，其中，利用在音频信号处理中屡次使用的近似。即，短时间分析窗中包含的信号和滤波器的卷积在该滤波器的长度M与短时间分析窗N相比充分短时，在短时间傅立叶变换域，能够通过该信号和滤波器的积来近似。若利用该近似，则式(16)可以在z域的单位圆上(相当于短时间傅立叶变换域)改写为以下的式(17)。

$W_N\left(y_{t0,N}\left(z\right)\right)\approx \underset{\tau =0}{\overset{K_R}{\Sigma }}{W}_N^R\left(c_{\mathrm{\tau M},M}\left(z\right)\right)W_N\left(x_{t0-\mathrm{\tau M},N}\left(z\right)\right)---\left(17\right)$

若使用短时间傅立叶变换表示，则式(17)成为以下的式(18)那样。

$Y_n\approx \underset{\tau =0}{\overset{K_R}{\Sigma }}\mathrm{diag}\left(X_{n-\tau }\right)C_{\tau }---\left(18\right)$

其中，n和τ是短时间帧的指示符，Y_n、C_n、X_n分别表示作为元素而具有从对应于y(z)、c(z)、x(z)的时域信号通过时间窗切出的信号的短时间傅立叶变换后的各频带的值的矢量，diag(X)是作为对焦分量而具有矢量X的分量的对角矩阵。另外，在本说明书中，如下表示短时间傅立叶变换。其中，t_τ表示帧τ的最初的采样的离散时间指示符。

$X_{\tau ,u}=\underset{\tau =0}{\overset{U-1}{\Sigma }}{x}_{t_{\tau }+t}\exp (-j2\mathrm{\pi ut}/U)---\left(19\right)$

X_r＝[X_r，0 X_r，1…X_r，U-1]^T                   (20)

通过式(18)，时域的卷积运算能够作为按频率观测信号的每个频带的卷积运算来计算。在式(17)中，M是相当于帧偏移的值，因此在该近似计算中，帧偏移M需要与窗函数W_N()的窗长度N相比充分小的值。

以上<补充说明：对于频率信号的卷积运算>结束。

在对上述式(12’)的两边使用式(16)，例如实施短时间傅立叶变换，从而得到以下的式(22)。

$X_n^{\left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =D}{\overset{K_R}{\Sigma }}\mathrm{diag}\left(X_{n-\tau }^{\left(q\right)}\right)C_{\tau }^{\left(q\right)}+{\tilde{s}}_n---\left(22\right)$

式(22)等价于式(22a)。

$X_{n,u}^{\left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =D}{\overset{K_R}{\Sigma }}{X}_{n-\tau ,u}^{\left(q\right)}{C}_{\tau ,u}^{\left(q\right)}+{\tilde{S}}_{n,u}---\left(22a\right)$

其中，D对应于式(22)的延迟d，是通过帧数表示对频率信号中的过去的观测信号导入的延迟的量。相邻的帧的频率信号在时域相互具有重叠。因此，帧n的观测信号(式(22)的左边X_n⁽¹⁾)中包含的音频信号的一部分还包含在与其正前的帧对应的观测信号中。从而，在式(22)中，若利用包括正前的帧的过去的观测信号预测X_n⁽¹⁾，则也能够预测音频信号的一部分。观测信号的可预测的部分不包含在残差信号中，因此音频信号的一部分通过回响去除而被去除。为了防止上述情况，在利用频率信号的本发明中，如式(22)所示，在预测当前的观测信号时不应用正前的帧的观测信号，而只使用某程度的延迟D以上分离的过去的观测信号。另外在满足d＝DM时，上述式(12’)和式(22)一致。在以下的说明中，将式(22)作为用于表示观测信号和音频信号的关系的式来说明本实施方式。在式(22)中，X_n^(q)相当于与在第q信道的麦克风接收的时域信号有关的短时间傅立叶变换。短时间傅立叶变换遵照式(19)、(20)。其中，n表示帧号。此外，频带u(u＝0、…、U-1)中的按频率观测信号表示为X_n，u^(q)。为了求出该X_n，u^(q)，分割单元302每M采样，时间偏移并应用短时间分析窗，并变换为频域。由此，得到对每个频带分割的按频率观测信号X_n，u^(q)。

此外，将在后面叙述细节的估计单元306_u估计用于从按频率观测信号X_n，u^(q)去除回响的回响去除滤波器。若得到作为该回响去除滤波器的系数的预测系数C_τ^(q)，则能够如下估计目标信号(包括初始反射音的音频信号)S~_n。

${\tilde{S}}_n{X}_n^{\left(1\right)}-\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =D}{\overset{K_R}{\Sigma }}\mathrm{diag}\left(X_{n-\tau }^{\left(q\right)}\right)C_{\tau }^{\left(q\right)}---\left(23\right)$

式(23)按照频带对于S_n~＝[S_n，0~，S_n，1~，…，S_n，U-1~]的各元素，还能够表示为以下的式(24)那样。

${\tilde{S}}_{n,u}=X_{n,u}^{\left(1\right)}-\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =D}{\overset{K_R}{\Sigma }}{X}_{n-\tau ,u}^{\left(q\right)}{C}_{\tau ,u}^{\left(q\right)}---\left(24\right)$

这里若利用式(25)-(28)，则式(24)还能够表示为式(29)那样。

$C_u=[C_u^{\left(1\right)},C_u^{\left(2\right)}···C_u^{\left(Q\right)}]---\left(25\right)$

$C_u^{\left(q\right)}=[C_{D,u}^{\left(q\right)},C_{D+1,u}^{\left(q\right)}···C_{K_R,u}^{\left(q\right)}]---\left(26\right)$

$B_{n-D,u}=[B_{n-D,u}^{\left(1\right)},B_{n-D,u}^{\left(2\right)}···B_{n-D,u}^{\left(Q\right)}]---\left(27\right)$

$B_{n-D,u}^{\left(q\right)}=[X_{n-D,u}^{\left(q\right)},X_{n-D-1,u}^{\left(q\right)}···X_{n-K_R,u}^{\left(q\right)}]---\left(28\right)$

${\tilde{S}}_{n,u}=X_{n,u}^{\left(1\right)}-B_{n-D,u}{C}_u^T---\left(29\right)$

其中T表示矢量和矩阵的转置。在本实施例中，将C_u称为第u个频带的回响去除滤波器。另外，式(29)中的B_n-D，uC_u^T相当于将在每个信道卷积了B_n-D，u和C_u^T而得到的信号，关于全部q进行了加法运算的信号。在估计单元306_u估计回响去除滤波器C_u，去除单元308_u基于式(29)去除回响信号。

此外，若将0_D-1设为全部元素为0的D-1维的行矢量，则回响去除滤波器W_u还能够定义为以下那样。

W_u＝[1，0_D-1，C_u⁽¹⁾，0，0_D-1，C_u⁽²⁾，…，0，0_D-1，C_u^(Q)]

此时，通过去除单元308_u基于以下式去除回响信号。

${\tilde{S}}_{n,u}={\zeta }_{n,u}{W}_u^T---\left(30\right)$

${\zeta }_{n,u}=\left(\begin{array}{cccc}{\zeta }_{n,u}^{\left(1\right)}& {\zeta }_{n,u}^{\left(2\right)}& ···& {\zeta }_{n,u}^{\left(Q\right)}\end{array}\right)$

${\zeta }_{n,u}^{\left(q\right)}=\left(\begin{array}{cccc}{X}_{n,u}^{\left(q\right)}& X_{n-1,u}^{\left(q\right)}& ···& X_{n-K_R,u}^{\left(q\right)}\end{array}\right)$

如上那样，若估计单元306_u能够估计回响去除滤波器C_u或W_u，则去除单元308_u基于式(29)或式(30)，能够去除回响信号。接着，在说明回响去除滤波器的估计之前，说明声源模型。

在声源模型存储单元304存储对每个频带表示了按频率观测信号的特征的声源模型。

该实施例的声源模型将音频信号能够取的值的倾向通过概率分布来表示。此外，基于该概率分布定义最佳化函数。关于有效的声源模型例如时变正规分布，且如下定义所求的按频率信号S_n~的概率密度函数。

p(S_n~)＝N(S_n~；0，Ψ_n)        (31)

Ψ_n∈Ω_Ψ        (32)

其中，N(S_n~；0，Ψ_n)是平均为0且声源模型的协方差矩阵Ψ_n＝E(S_n~(S_n~)^*T)的多维复数正态分布，Ψ_n按每个短时间帧n取不同的值，或是相同的值。在以下的说明中，Ψ_n称为模型协方差矩阵，假设模型协方差矩阵Ψ_n是每个短时间帧n取不同的值的对角矩阵。此外，“*”表示复数共轭。Ω_Ψ表示包括Ψ_n可取的全部值的集合(即，Ψ_n的参数空间)。若设ψ_n，u²＝E(S_n，u~S_n，u~^*T)表示Ψ_n的第u个对角元素，则由于Ψ_n是对角矩阵，因此概率分布函数能够在各频带独立而成为

p(S_n，u~)＝N(S_n，u~；0，ψ_n，u²)        (33)

                                              。

每个频带的估计单元306_u基于与声源模型关联定义的观测信号的最佳化函数，从按频率观测信号估计回响去除滤波器(步骤S4)。具体说明回响去除滤波器的估计的细节。

如上述式(25)所示，回响去除滤波器C_u通过由与全部的麦克风有关的观测信号的预测系数C_u^(q)构成的矢量来表示。预测系数C_u^(q)是频域预测系数。ψ_u²表示模型协方差矩阵的第u对角元素的时间序列，并设ψ_u²＝{ψ_n，U²}。此外，设θ_u＝{C_u、ψ_u²}表示估计参数的集合。而且，将全部频带的估计参数全体的集合表示为θ＝{θ₀，θ₁，…，θ_U-1}。此外，如下那样，作为每个频带的最佳化函数定义对数似然函数L_u(θ_u)，作为全频带的最佳化函数定义对数似然函数L(θ)。

$L\left(\theta \right)=\underset{u}{\Sigma }{L}_u\left({\theta }_u\right)---\left(35\right)$

式(34)基于式(29)和(33)，能够表示为以下的式(36)。

通过估计使式(35)的左边最大化的参数，从而能够求出回响去除滤波器的预测系数C_u^(q)。式(35)的最大化能够通过以下的最佳算法来实现。

1.关于全部频带u，将初始值例如以下的式(37)那样定义。

$C_{n,u}^{\left(q\right)}=0---\left(37\right)$

2.重复以下的2个式直到收敛。

2-1.关于全部频带u，固定C_n，u^(q)，从而更新模型协方差矩阵Ψ_n，使得最佳化函数L(θ)最大化。

${\hat{\Psi }}_n=\arg \underset{\Psi \in {\Omega }_{\Psi }}{\max }L\left(\theta \right)\to {\Psi }_n---\left(38\right)$

2-2.固定Ψ_n，从而关于全部频带u，更新回响去除滤波器C_u，使得最佳化函数L_u(θ_u)最大化。

${\hat{C}}_u={\left(\underset{n}{\Sigma }\frac{B_{n-D,u}^{*T}{B}_{n-D,u}}{{\psi }_{n,u}^2}\right)}^{+}\underset{n}{\Sigma }\frac{B_{n-D,u}^{*T}{X}_{n,u}^{\left(1\right)}}{{\psi }_{n,u}^2}\to C_u---\left(39\right)$

其中，在上述算法的表示中，将参数A的值更新为B的操作表示为“A→B”。此外，“+”表示摩尔-彭罗斯(Moore-Penrose)的虚拟反转逆矩阵。另外，在上述算法中需要计算的有关观测信号的协方差矩阵H’(ψ_n，u²)成为如以下的式(40)那样。

$H^{\prime }\left({\psi }_{n,u}^2\right)=\underset{n}{\Sigma }\frac{B_{n-D,u}^{*T}{B}_{n-D,u}}{{\psi }_{n,u}^2}---\left(40\right)$

基于该最佳化算法，根据最终得到的C_u构成回响去除滤波器。去除单元308_u基于式(29)或式(30)，将该回响去除滤波器C_u或W_u卷积在按频率观测信号X_n，u^(q)，从而从X_n，u^(q)去除回响信号，求出按频率目标信号S_n，u~(步骤S12)。

此外，合并单元310合并每个频带的按频率信号S_n，u~，并变换为时域，从而输出目标信号S_t~(步骤S14)。具体地说，能够使用用于将短时间傅立叶变换的帧的时间序列变换为时域信号的一般方法。即，对各帧n，对S_n~＝[S_n，0~，S_n，1~，…，S_n，U-1~]适用短时间傅立叶反变换，得到各帧的时间信号，并对各帧的信号进行重叠加法运算(overlap-added)，从而得到目标信号S_t~。帧τ的短时间傅立叶反变换通过式(40a)来表示。重叠加法运算通过对应用短时间傅立叶反变换而得到的各帧的时间信号应用任意的时间窗，并以在分割单元利用的相同的帧偏移宽度M对信号进行加法运算而实现。具体的计算式通过式(40b)来表示。其中，w_t^I表示长度为N的时间窗，floor(a)表示a以下的最大的整数。

$x_{\tau ,1}=\frac{1}{U}\underset{u=0}{\overset{U-1}{\Sigma }}{X}_{r,u}\exp (j2\mathrm{\pi ut}/U)---\left(40a\right)$

$x_{\tau }=\underset{\tau =\mathrm{floor}((1-N)/M)+1}{\overset{\mathrm{floor}(t/M)}{\Sigma }}{w}_{t-\mathrm{\tau M}}^I{x}_{\tau ,t-\mathrm{\tau M}}---\left(40b\right)$

说明该实施例1的回响去除装置300的效果。根据该回响去除装置300的观测信号x_t^(q)(q＝1，…，Q)，将回响去除处理作为各频带的运算而能够近似计算。每M采样时间偏移并应用长度N的短时间分析窗，从而进行对频域信号的变换，从而能够缩短各频带的回响去除滤波器的长度。此外，能够缩小回响去除滤波器的估计所需的协方差矩阵的尺寸。若说明该理由，一般，回响去除滤波器的尺寸与用于求出该回响去除滤波器的协方差矩阵的尺寸相等。此外，由于每M采样，一边时间偏移一边(适用长度N的短时间分析窗)切取N采样量，从而进行频域变换处理，因此与现有技术1比较，卷积的回响去除滤波器的尺寸也变小。从而协方差矩阵的尺寸也变小。这从式(1)、式(40)也能够明白，即，若比较式(1)所示的协方差矩阵H(r)的尺寸和式(40)所示的协方差矩阵H’(ψ_n，u²)的尺寸，则现有技术1的协方差矩阵H(r)的尺寸依赖于预测滤波器长度(室内脉冲响应长度)K。但是，在本实施例1中使用的协方差矩阵H’(ψ_n，u²)依赖于K_R(即，<K/M>)。这是因为，如式(35)所示，构成协方差矩阵H’(ψ_n，u²)的B_n-D，u^(q)的元素的数量(抽头数)是K_R-D个。从而，与现有技术1相比，能够理解在本实施例1中使用的协方差矩阵的尺寸能够变小。在回响去除滤波器的估计中，除了需要协方差矩阵的计算，还需要其逆矩阵的计算。有关上述计算的计算成本占据回响去除处理整体的计算成本的大部分。而且，该两者的计算成本能够通过缩小协方差矩阵的尺寸而缩小。如上所述，在本实施例中，能够大幅削减回响去除处理整体的计算成本。

实施例2

在实施例1中，在各频带，通过对观测信号卷积回响去除滤波器，从而实现回响去除。一方面已知，估计回响信号并求出观测信号的能量和回响信号的能量之差即差信号方法比实施例1的回响去除方法能够实现不容易受到回响去除滤波器的估计误差的影响的回响去除。例如，记载在“k.Kinoshita，T.Nakatani，and M.Miyoshi，“Spectral subtraction steered by multi-step forward linear prediction for single channel speech dereverberation，”Proc.ICASSP-2206，vol.1，pp.817-820，May，2006.”。在本实施例2中应用该想法。

说明实施例2的回响去除装置400。图5表示回响去除装置400的功能结构例，图6表示主要处理的流程。回响去除装置400与回响去除装置300相比，不同之处在于去除单元308_u被去除单元407_u代替。去除单元407_u由每个频带的回响信号生成部件408_u、每个频带的按回响信号频率的功率生成部件410_u、每个频带的按观测信号频率的功率生成部件412_u、每个频带的减法运算部件414_u构成。

若通过分割单元302，观测信号被分割为每个频带(步骤S22)，通过估计单元306_u，估计每个频带的回响去除滤波器(步骤S4)，则回响信号生成部件408_u利用回响去除滤波器和按频率观测信号X_n，u^(q)，生成按频率回响信号R_n，u。(步骤S22)。具体地说，例如通过以下的式(41)求出按频率回响信号R_n，u。

$R_{n,u}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =D}{\overset{K_R}{\Sigma }}\mathrm{diag}\left(X_{n-\tau ,u}^{\left(q\right)}\right)C_{\tau ,u}^{\left(q\right)}---\left(41\right)$

按回响信号频率的功率生成部件410_u求出按频率回响信号R_n，u的按频率功率|R_n，u|²(步骤S24)。另一方面，按观测信号频率的功率生成部件412_u例如求出通过第1信道的麦克风接收的按频率观测信号的按频率功率|X⁽¹⁾_n，u|²(步骤S26)。然后，减法运算部件414_u通过计算按频率回响信号的按频率功率和按频率观测信号的按频率功率之差，从而求出差信号|X⁽¹⁾_n，u|²-|R_n，u|²，并基于用于计算该差信号的按频率观测信号X⁽¹⁾_n，u和该差信号，求出按频率目标信号(步骤S28)。例如基于以下的式来求出按频率目标信号S_n，u~。

$S_{n,u}^{~}=G_{n,u}{X}_{n,u}^{\left(1\right)}$

$G_{n,u}=\max \{\frac{{\left|X_{n,u}^{\left(1\right)}\right|}^2-{\left|R_{n,u}\right|}^2}{{\left|X_{n,u}^{\left(1\right)}\right|}^2},G_0\}$

其中，设max{A，B}是选择A、B中的大的函数，G₀是G₀＞0，且决定用于通过功率减法运算来抑制信号的能量的下限的地面系数(flooring coefficient)。此外，合并单元416将该按频率目标信号变换为时域，从而求出目标信号s_t~(步骤S30)。

该回响去除装置400即使通过实施例1的回响去除装置300而在回响去除滤波器中包含估计误差也能够进行音质的劣化少的回响去除。

此外，现有技术的回响去除处理只能在时域工作。但是，在实施例1、2中说明的回响去除装置300、400能够在频域工作，因此能够与盲声源分离(blind source separation)和维纳滤波器(Wiener filter)等在频域工作的其他多个有用的声音增强技术进行组合。

实施例3

图7表示实施例3的回响去除装置500的功能结构例。若说明与实施例1的回响去除装置300的主要不同点则如下：(1)相对于回响去除装置300的分割单元302通过一边时间偏移时域的观测信号一般变换为频域从而分割给频带，回响去除装置500的分割单元502通过子带分割而分割给频带；(2)回响去除装置300的去除单元、合并单元的处理在频域进行，相对于此，本实施例的回响去除装置500的去除单元、合并单元的处理在时域进行。

将子带分割的信号设为子带信号，将子带的数量设为V，将子带的指示符设为v(v＝0、…、V-1)。估计单元506_v对各子带信号估计回响去除滤波器，去除单元508_v对各子带信号去除回响。通过合并单元510进行合并从而求出目标信号s_t~。分割单元502的子带分割处理、合并单元510的合并处理记载在“M.R.Portnoff，“Implementation of the digital phase vocoder using the fast Fourier transform，”IEEE Trans.ASSP，vol.24，No.3，pp.243-248，1976.(以下称为“非专利文献A”。)”和“J.P.Reilly，M.Wilbur，M.Seibert，and N.Ahmadvand，“The complex subband decomposition and its application to the decimation of large adaptive filtering problems，”IEEE Trans.Signal Processing，vol.50，no.11，pp.2730-2743，Nov.2002”等。以下的说明中，利用非专利文献A的技术进行说明。在该非专利文献A中，记载了后述的式(50)。此外，主要的处理的流程与图4相同，因此省略。

首先，说明音频信号与观测信号的关系。分割单元502对观测信号进行子带分割，分割给V个频带的每一个(子带)。若该分割根据非专利文献A的定义通过式来表示则成为以下的式(50)那样。

$x_{t,v}^{\left(q\right)}=\underset{\tau =-N_h}{\overset{N_h}{\Sigma }}{x}_t^{\left(q\right)}{h}_{t-\tau }{e}^{-j2\mathrm{\pi v\tau }/V}---\left(50\right)$

其中，在各子带中，将应用观测信号的频率偏移和低通滤波器而得到的信号的采样指示符设为t(与被子带处理的之前的观测信号的离散时刻)相同，将与通过第q信道的麦克风接收的观测信号有关的第v(v＝0、…、V-1)个子带的第t个采样设为x_t，v^(q)。e^-j2πvτ/V是与第v个子带对应的频率偏移运算符，h_τ是长度为2N_h+1的低通滤波器的系数。此外，若将式(50)应用于上述式(12’)的两边则得到以下的式。

$x_{t,v}^{\left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =d}{\overset{K}{\Sigma }}{c}_{\tau }^{\left(q\right)}{x}_{t-\tau ,v}^{\left(q\right)}+{\tilde{s}}_{t,v}---\left(51\right)$

其中式(51)的右边的s_t，v~是对包含初始反射音的音频信号应用分割处理而得到的信号。在本实施例中将s_t，v~作为应求出的目标信号来处理。此外，分割单元502进行子带分割并对各子带信号进行下采样。例如将对通过第1信道的麦克风接收的观测信号x_t，v(1)和音频信号s_t，v的各时间序列以γ个采样间隔进行了下采样(采样的间隔剔除(thinning out))的信号的采样的指示符表示为b，将下采样后得到的子带信号表示为x_b，v^’(q)或s_b，v~’。将与下采样的信号的采样指示符b对应的、下采样之前的信号的采样指示符设为t_b。则，能够表示为以下的式(52)那样。

$x_{b,v}^{\prime \left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{\tau =d}{\overset{K}{\Sigma }}{c}_{\tau }^{\left(q\right)}{x}_{t_b-\tau ,v}^{\left(q\right)}+{\tilde{s}}_{b,v}^{\prime }---\left(52\right)$

另一方面，由于h_τ是低通滤波器因此在以该低通滤波器的截止频率的2倍以上的采样频率进行下采样的情况下，能够通过上采样高精度地恢复到下采样之前的信号。该上采样例如以以下的步骤进行。

步骤1.在下采样后的信号的各采样之间，插入γ-1个“0”。

步骤2.应用低通滤波器。

在步骤2.中一般使用有限长度脉冲响应滤波器。这表示通过上采样恢复的信号能够通过下采样后的信号的线性组合来表示。

若使用该关系则式(52)的右边的记载x_tb-τ，v^(q)能够如以下的式(53)那样表示。

β_τ，k相当于与上采样中的低通滤波器的系数对应而决定的系数，k₀相当于用于上采样的低通滤波器的滤波的延迟，k₀+k₁+1相当于用于上采样的低通滤波器的滤波器长度。若将式(53)代入式(52)进行整理，则得到以下的式(54)。

$x_{b,v}^{\prime \left(1\right)}=\underset{q=1}{\overset{Q}{\Sigma }}\underset{k=d^{\prime }}{\overset{K^{\prime }}{\Sigma }}{\alpha }_{k,v}^{\left(q\right)}{x}_{b-k,v}^{\prime \left(q\right)}+{\tilde{s}}_{b,v}^{\prime }---\left(54\right)$

其中，α_k，v^(q)表示在将式(53)代入式(52)而整理时，成为x’_b-k，v^(q)的项的系数。d’表示基于α_k，v^(q)的滤波的延迟，K’表示基于α_k，v^(q)的滤波的滤波器长度。基于式(52)、(53)以及间隔剔除间隔γ的关系，能够决定为d’≒d/γ-k₀、K’≒K/γ+k₁。在d’≥1时，式(54)表示在对个子带信号以α_k，v^(q)作为预测系数(在估计单元506_v估计的回响去除滤波器的系数)从过去的观测信号预测了当前的观测信号时，残差信号成为包括初始反射音的音频信号的关系。在以下的说明中，将式(54)作为表示各子带信号中的观测信号和音频信号的关系的式进行处理。

这里，定义式(55)-(58)。

${\alpha }_v=\left(\begin{array}{ccccc}{\alpha }_v^{\left(1\right)}& ···& {\alpha }_v^{\left(q\right)}& ···& {\alpha }_v^{\left(Q\right)}\end{array}\right)---\left(55\right)$

${\alpha }_v^{\left(q\right)}=[{\alpha }_{d^{\prime },v}^{\left(q\right)},{\alpha }_{d^{\prime }+1,v}^{\left(q\right)}···{\alpha }_{K^{\prime },v}^{\left(q\right)}]---\left(56\right)$

$F_{b-d^{\prime },v}=\left(\begin{array}{ccccc}{F}_{b-d^{\prime },v}^{\left(1\right)}& ···& F_{b-d^{\prime },v}^{\left(q\right)}& ···& F_{b-d^{\prime },v}^{\left(Q\right)}\end{array}\right)---\left(57\right)$

$F_{b-d^{\prime },v}^{\left(q\right)}=[x_{b-d^{\prime },v}^{\prime \left(q\right)},x_{b-d^{\prime }-1,v}^{\prime \left(q\right)}···x_{b-K^{\prime },v}^{\prime \left(q\right)}]---\left(58\right)$

此时，式(54)能够表示为式(59)那样。

${\tilde{s}}_{b,v}^{\prime }=x_{b,v}^{\prime \left(1\right)}-F_{b-d^{\prime },v}·{\alpha }_v^T---\left(59\right)$

在实施例3中，将α_v设为对于第v个子带信号的回响去除滤波器，去除单元508_v基于上述式(59)进行回响信号的去除。另外，若将0_d’-1设为全部元素为0的d’-1维的行向量，则回响去除滤波器w_v还能够表示为以下的式(60)那样。

$w_v=\left(\begin{array}{ccccccccccc}1& 0_{d^{\prime }-1}& {\alpha }_v^{\left(1\right)}& ···& 0& 0_{d^{\prime }-1}& {\alpha }_v^{\left(q\right)}& ···& 0& 0_{d^{\prime }-1}& {\alpha }_v^{\left(Q\right)}\end{array}\right)---\left(60\right)$

此时，去除单元508_v基于式(61)，进行回响信号的去除。

${\tilde{s}}_{b,v}^{\prime }={\xi }_{b,v}{w}_v^T---\left(61\right)$

${\xi }_{b,v}=\left(\begin{array}{ccccc}{\xi }_{b,v}^{\left(1\right)}& ···& {\xi }_{b,v}^{\left(q\right)}& ···& {\xi }_{b,v}^{\left(Q\right)}\end{array}\right)$

${\xi }_{b,v}^{\left(q\right)}=\left(\begin{array}{cccc}{x}_{b,v}^{\left(q\right)}& x_{b-1,v}^{\left(q\right)}& ···& x_{b-K^{\prime },v}^{\left(q\right)}\end{array}\right)$

接着，说明估计单元506_v的回响去除滤波器的估计方法。在该实施例的声源模型存储单元504中存储的声源模型与实施例1、2相同，通过概率分布来表示音频信号可取的倾向，基于此来定义最佳化函数。作为有效的声源模型，例如时变正规分布。在以下的说明中，作为最简单的声源模型，导入在各子带间信号独立的模型。此外，假设各子带信号是频谱平坦，且只有信号的能量随时间变化的时变白色正态过程。

与上述式(31)、(32)一样，定义参数空间，如下变更。此时，可以如下那样定义s_b~’＝[s_b，0~’，s_b，1~’，…，s_b，V-1~’]^T的概率密度函数。

p(s_b~’)＝N(s_b~’；0，Ψ_b’)         (31’)

Ψ_b’∈Ω_Ψ’           (32’)

其中，N(s_b~’；0，Ψ_b’)表示平均为0且声源模型的协方差矩阵Ψ_b’＝E(s_b~’(s_b~’)^*T)的多维复数正态分布，Ψ_b’可以每采样b取不同的值，或者取相同的值。在以下的说明中，将Ψ_b’称为模型协方差矩阵，假设模型协方差矩阵Ψ_b’是每个采样取不同的值的对角矩阵。Ω_Ψ’表示Ψ_b’可取的全部值的集合(即，Ψ_b’的参数空间)。ψ_b，v’²＝E(s_b，v~’(s_b，v~’)*)是Ψ_b’的第v个对角元素。由于Ψ_b’是对角矩阵，因此概率密度函数能够对各子带独立地成为p(s_b，v~’)＝N(s_b，v~’；0，ψ_bv’²)。ψ_v’²表示模型协方差矩阵的第v对角元素的时间序列，表示为ψ_v’²＝{ψ_b，v’²}。此外，设θ_v＝{α_v，ψ_v’²}表示有关子带v的估计参数的集合。而且，将全部子带的全部估计参数的集合表示为θ’＝{θ₀，θ₁，…，θ_V-1}。此外，如下那样，作为每个频带的最佳化函数定义对数似然函数L_v(θ_v)，作为全频带的最佳化函数定义对数似然函数L’(θ’)。

$L^{\prime }\left({\theta }^{\prime }\right)=\underset{v}{\Sigma }{L}_v\left({\theta }_v\right)---\left({35}^{\prime }\right)$

基于式(59)和式(31’)，可以将式(63)表示为式(64)那样。

通过估计使式(64)最大化的参数，能够得到回响去除滤波器的系数的估计值。式(64)的最大化可通过以下的最佳化算法来实现。

1.关于全部的子带v，如以下式(65)那样决定初始值。

${\alpha }_{b,v}^{\left(q\right)}=0---\left(65\right)$

2.重复以下的2个式直到收敛。

2-1.关于全部的子带v，固定α_b，v^(q)，从而更新模型协方差矩阵Ψ_b’，使得最佳化函数L’(θ’)最大化。

${{\hat{\Psi }}_b}^{\prime }=\arg \underset{{{\Psi }_v}^{\prime }\in {\Omega }_{{\Psi }^{\prime }}}{\max }{L}^{\prime }\left({\theta }^{\prime }\right)\to {{\Psi }_b}^{\prime }---\left(66\right)$

2-2.固定Ψ_b’，从而关于全部子带v，更新回响去除滤波器系数α_v，使得最佳化函数L_v(θ_v)最大化。

${\hat{\alpha }}_v={\left(\underset{b}{\Sigma }\frac{F_{b-d^{\prime },v}^{*T}{F}_{b-d^{\prime },v}}{{\psi }_{b,v}^2}\right)}^{+1}\underset{b}{\Sigma }\frac{F_{b-d^{\prime },v}^{*T}{x}_{b,v}^{\prime \left(1\right)}}{{\psi }_{b,v}^{\prime 2}}\to {\alpha }_v---\left(67\right)$

基于最终得到的α_v，估计单元506_v构成回响去除滤波器，去除单元508_v通过上述式该回响去除滤波器基于上述式(59)或(61)去除回响信号，从而求出按频率目标信号s_b，v~’。此外，合并单元510通过与将按频率目标信号s_b，v~’上采样处理一同合并各子带信号，从而求出目标信号s_t~。

如以上说明那样，在子带处理中，将观测信号分割为每个子带的时域信号后，以γ个间隔进行下采样，从而能够使各频带的时域信号的采样频率成为1/γ。

在本实施例中，对各频带的时域信号分别进行回响去除处理，并将它们进行合并，从而实现全频域的回响去除。对于时域信号，若比较进行下采样的情况和不进行的情况，则进行下采样的情况能够将用于回响去除滤波器的估计的协方差矩阵的尺寸变小。这是因为协方差矩阵的尺寸由回响去除滤波器的滤波器长度决定，与房间的脉冲响应的抽头数对应地决定回响去除滤波器的滤波器长度K，且物理上，在采样频率越小则相同时间长度的脉冲响应成为越少抽头数。换言之，通过以γ个间隔进行下采样，回响去除滤波器的滤波器长度K’(＝K/γ+k₁)，从而比现有技术的回响去除滤波器的滤波器长度K小。

若回响去除滤波器的滤波器长度变小，则如上所述，可以将在回响去除滤波器估计时使用的协方差矩阵的尺寸变小，因此能够削减回响去除滤波器的估计处理的计算成本。

此外，在该下采样以低通滤波器的截止频率的2倍以上的采样频率进行的情况下，通过与该下采样处理一并进行的子带分割处理求得的子带信号具有通过上采样能够高精度地恢复的性质。从而，即使在合并单元510的合并处理时进行上采样，目标信号也不会变差。

实施例4

图8表示实施例4的回响去除装置600的功能结构例。回响去除装置600与回响去除装置500相比，不同点在于去除单元508_v被去除单元607_v代替。通过该代替，与回响去除装置500相比，能够进行不容易受到回响去除滤波器的估计误差的影响的回响去除。理由如在实施例2中说明那样。去除单元607_v对应于在实施例2中说明的去除单元407_v。去除单元607_v包括每个频带的回响信号生成部件608_v、每个频带的按回响信号频率的功率生成部件610_v、每个频带的按观测信号频率的功率生成部件612_v、每个频带的减法运算部件614_v。

回响信号生成部件608_v利用回响去除滤波器α_v和观测信号x_t，v^(q)，求出按频率回响信号r_b，v。具体来说通过以下的式(70)求出。

r_b，v＝F_b-d’，v·α_v^T        (70)

此外，按回响信号频率的功率生成部件610_v求出按频率回响信号的按频率功率|r_b，v|²。此外，按观测信号频率的功率生成部件612_v，求出通过第1信道的麦克风接收的观测信号x_b，v⁽¹⁾的按频率功率|x_b，v⁽¹⁾|²。然后，减法运算部件614_v通过计算按频率回响信号的按频率功率和按频率观测信号的按频率功率之差，从而求出差信号|x_b，v⁽¹⁾|²-|r_b，v|²，并基于用于该差信号的计算的按频率观测信号x_b，v⁽¹⁾和该差信号，求出按频率目标信号(步骤S28)。例如基于以下的式，求出按频率目标信号s_b，v~’。例如，通过以下的式求出按频率目标信号s_b，v~’。

${\tilde{s}}_{b,v}^{\prime }=G_{b,v}{x}_{b,v}^{\prime \left(1\right)}---\left(71\right)$

$G_{b,v}=\max \{\frac{{\left|x_{b,v}^{\prime \left(1\right)}\right|}^2-{\left|{\tilde{r}}_{b,v}\right|}^2}{{\left|x_{b,v}^{\prime \left(1\right)}\right|}^2},G_0\}---\left(72\right)$

其中，设max{A，B}是选择A和B中较大一方的函数，G₀是G₀＞0，且决定用于通过功率减法运算来抑制信号的能量的下限的地面系数(flooring coefficient)。

此外，通过合并单元510合并每个按频率目标信号s_b，v’~(v＝0，…，V-1)，并作为目标信号s_t~来输出。

通过回响去除装置600那样的结构，与回响去除装置500相比，不大受回响去除滤波器的估计误差的影响就能够进行回响信号的去除。

实施例5

在实施例1～4中说明的回响去除装置300～600中，以事先得到全部的信号的成批处理为前提而构成。作为实施例5，还可以对通过麦克风接收的观测信号依次进行回响信号的去除。例如，通过估计单元估计的回响去除滤波器以预先决定的时间间隔被估计、更新。在该更新时，通过对该时刻之前得到的观测信号的全部或一部分应用上述最佳化算法，从而估计回响去除滤波器。与该估计一同，回响去除装置300的估计单元306_u(参照图3)、回响去除装置400的回响信号生成部件408_u(参照图5)、回响去除装置的估计单元506_v(参照图7)、回响去除装置600的回响信号生成部件608_v(参照图8)可以对应于对在各时刻依次得到的观测信号，将其之前得到的最新的回响去除滤波器适用于该时刻的观测信号的结构。通过该依次处理，能够进行更正确的回响信号的去除。

[声源模型的具体例]

以下，关于实施例1至实施例5的声源模型的具体例，表示集合Ω_Ψ、Ω_Ψ’的例来进行说明。主要说明实施例1、2、5。对于实施例3、4，由于通过对以下的说明中的各标号进行以下的替换就能够构成具体例，因此省略说明。

Ω_Ψ→Ω_Ψ’

Ψ_u→Ψ_v’

ψ_n，u→ψ_b，v’

X_n，u^(q)→x_b，v^(q)’

S_n，u~→s_b，u~’

B_n，u→F_b，v

D→d’

C_u→α_v

i_n→i_b

式(38)→式(66)

式(39)→式(67)

306_u→506_v

(1)作为第一个具体例，设集合Ω_Ψ是由任意的正定对角矩阵构成的集合。这表示ψ_n，u²可以取任意正值。此时在上述最佳化算法中，式(38)的更新式可以置换为在全部频带中单独计算的以下的更新式(80)。另外，关于式(39)的更新式没有变更。

${\hat{\psi }}_{n,u}^2={(X_{n,u}^{\left(1\right)}-B_{n-D,u}{C}_u^T)(X_{n,u}^{\left(1\right)}-B_{n-D,u}{C}_u^T)}^{*}---\left(80\right)$

(2)说明第2个具体例。与非专利文献1中记载的技术相同，说明通过有限状态机将音频信号的波形模型化的情况。此时，集合Ω_Ψ成为由有限个正定对角矩阵构成的集合。各矩阵成为与对应于观测信号的短时间信号的频域信号可取的各个有限个状态的对应的协方差矩阵。这些有限个矩阵可以基于对事先在不包括回响的环境中接收的音频信号的频域信号或其协方差矩阵进行聚类(clustering)等方法来构成。此外，将有限个矩阵的数量设为Z，将指示符设为i(i＝1，…，Z)，将与状态i对应的协方差矩阵设为Ψ(i)。

据此，在上述重复算法中应估计的参数代替协方差矩阵而成为指示符的值。以下，将时刻n的状态设为i_n，将与状态i_n对应的协方差矩阵设为Ψ(i_n)，协方差矩阵Ψ(i_n)的对角元素设为ψ_u²(i_n)。各时刻的声源模型的状态i_n并不是在每个频带决定的值，是对全部频带决定一个的值。因此，基于对数似然函数决定的最佳化函数，对全部频带，可以如下的式(81)那样定义。

其中，设估计参数θ＝{C，I}由i_n的时间序列I＝{i₁，i₂…}和各频带的预测系数C＝{C₀，C₁，…，C_U-1}构成。基于该最佳化函数，在所述最佳化算法中，式(38)的更新式可以置换为关于全部频带的以下的更新式(82)。另外，对式(39)的更新式不进行变更。

${\hat{i}}_n=\arg \underset{i_n}{\max }\underset{u}{\Sigma }\log N(X_{n,u}^{\left(1\right)};B_{n-D,u}{C}_n^T,{\psi }_u^2\left(i_n\right))\to i_n---\left(82\right)$

通过从式(38)至式(82)的置换，估计单元306_u能够更正确地进行回响去除滤波器的估计。

(3)说明第3个具体例。通过将在(2)中说明的状态i_n假设为概率变量，从而能够构成基于更精密的声源模型的最佳化函数。作为一例，说明状态i_n通过一次马尔科夫过程就能够模型化的情况。通过马尔科夫过程的假设，能够成为p(I)＝p(i)II_np(i_n|i_n-1)。声源模型的参数是对于任意的状态i、j的p(i)、p(i|j)以及各状态的协方差矩阵Ψ(i)，这些参数和在不包括回响的环境接收的音频信号一并可以事先准备。此时用于去除回响信号的最佳化函数成为以下的式那样。