CN103840794B - 简化子带结构无次级路径有源控制方法 - Google Patents

Abstract

本发明公布了一种简化子带结构无次级路径有源控制方法。该方法对时延修正参考信号，正交时延修正参考信号以及误差信号做子带分解，全带自适应滤波器不做子带划分。无需次级路径建模，依次判断每个子带合适参考信号类型及该子带相应更新方向，直接更新全带自适应滤波器系数。其显著优点：无需次级路径建模，直接更新全带自适应滤波器系数，系统结构设计简单。

Description

简化子带结构无次级路径有源控制方法

一.技术领域

本发明涉及无次级路径建模的有源控制方法，针对宽带信号提出了一种简化子带结构无次级路径的有源控制方法。

二.背景技术

大多有源控制算法需要次级声源和误差传声器间的传递函数，即次级路径模型。前馈有源控制算法常采用滤波最小均方(FxLMS)算法，当真实次级路径传递函数和次级路径模型间差别较大时，FxLMS算法可能不稳定。Zhou等为避免对次级路径建模，提出一种无次级路径建模算法(D.Zhou and V. DeBrunner,“A new active noise controlalgorithm that requires no secondary path identification based on SPRproperty,”IEEE Trans.Signal Processing,vo1.55，no.5，pp.1719-1729，May2007)。这种算法无需对次级路径建模，在更新自适应滤波器系数的过程中，通过判断残余噪声能量降低与否，在0°和180°两个方向选择更新过程中步长前的符号，以降低算法残差。但当信号频率对应的次级路径传递函数的相位与自适应滤波器系数更新方向间的相位差接近+90°时，算法收敛较慢。Wu等对之做了改进，将算法改在频域实施，且引入±90°两个更新方向(M.Wu，G.Chen，and X.Qiu，“An improved active noise control algorithm withoutsecondary path identification based on the frequency-domain subbandarchitecture,”IEEE Trans.on Speech andAudio Processing，vol.16，pp.1409-1419，Nov2008)，提高了这种情况下Zhou等算法的收敛速度。

当参考信号为宽带信号时，Zhou等和Wu等的算法均指出可采用子带结构将参考信号、误差信号以及自适应滤波器均做子带分解，并通过子带参考信号和子带误差信号调整对应子带自适应滤波器权系数，之后通过权变换的方法将各子带自适应滤波器系数堆砌为全频带自适应滤波器系数。

Morgan等提出了无延迟子带自适应滤波算法(D.R.Morgan and J.C.Thi，“Adelayless subband adaptive filter architecture,”IEEE Trans.Signal Processing，vol.43，pp.1819-1830，1995)，在他们的算法中，控制声源的输入信号为全通带自适应滤波器的输出信号而非采用的子带自适应滤波器输出信号，从而避免信号通路中引入的延迟，并同时保证子带自适应滤波的优点。

Courville等(M.Courville and P.Duhamel，“Adaptive flitering in subbandsusing a weighted criterion,”IEEE Trans.Signal Processing，vol.56，pp.2359-2371，1998.)曾提出过一种子带结构，对全带参考信号、期望信号以及误差信号做子带分解，而全带自适应滤波器系数不做子带分解，并采用子带信号更新全带自适应滤波器系数，避免了堆砌过程。

DeBrunner等(V.DeBrunner，L.DeBrunner,and L.Wang，“Subband adaptivefilering with delay compensation for active noise control,”IEEE Trans.SignalProcessing，vo1.52，pp.2932-2937，2004)将该子带结构应用到有源控制系统，通过子带滤波参考信号和子带误差信号更新全带自适应滤波器系数，分析子带误差信号组成，从而补偿由次级路径和分析滤波器组对误差信号造成的时延。

中国专利CN101393736B公开一种无次级路径建模的有源控制方法，该方法在频域实施，在0°，180°，90°，-90°四个方向搜索合适的更新方向，对单频和窄带信号均有很好的效果，系统结构简单，易于实现。然而当参考信号为宽带噪声时，该方法可能无法搜索到合适的更新方向，且对如何采用无次级路径建模有源控制方法实现宽带有源降噪并未提及。

根据Wu等提出的频域子带结构无次级路径有源控制算法(M.Wu，G.Chen，andX.Qiu，“An improved active noise control algorithm without secondary pathidentification based on the frequency-domain subband architecture,”IEEETrans.on Speech andAudio Processing,vol.16，pp.1409-1419，Nov2008)，在该算法中可将宽带噪声做子带分解划分为多个窄带信号，依次对每个子带在0°，180°，90°，-90°四个方向中搜索合适的更新方向，并调整每个子带滤波器系数，然后将子带滤波器系数经过权变换堆砌为全带自适应滤波器系数，然而该算法中需要子带个数较多，且需要经过权变换合成全带滤波器系数，增加计算量。据此，本发明的方法是在Wu等算法上的一种改进，采用简化子带结构，将时延修正参考信号，正交时延修正参考信号和误差信号做子带分解，而全带自适应滤波器不做子带划分，依次判断每个子带参考信号类型，并搜索该子带更新方向，然后直接更新全带自适应滤波器系数，实现宽带噪声有源控制，系统结构简单，能有效减少子带个数，且计算量小。

三.发明内容

1.发明目的：提出一种简化子带结构无次级路径有源控制方法。该方法无需次级路径建模，根据各个子带参考信号类型以及该子带更新方向一起直接更新全带自适应滤波器系数，系统结构简单，易于实现。

2.技术方案

本发明的目的通过以下技术方案实现一种简化子带结构无次级路径有源控制方法，其步骤为：

(1)计算实验中所用电路板电时延τ₁＝Q/f_s，其中f_s为采样频率，Q为参考信号从控制器参考信号输入端到控制信号输出端延迟的采样点数。

(2)计算实验中声时延τ₂，计算公式为τ₂＝l/c，其中，l为次级声源与误差传感器之间的距离，c为声速，一般环境中取340m/s。

(3)基于步骤(1)电时延τ₁和步骤(2)声时延τ₂，计算总的时延τ＝τ₁+τ₂，总时延对应直达声时延点数为Δ，计算公式为Δ＝τf_s，其中f_s为采样频率。

(4)设计希尔伯特滤波器

用matlab中firls函数设计一个长度为L的希尔伯特滤波器即通过运行firls(L’，f，g，’hilbert')，其中L'＝L-1，为希尔伯特滤波器阶数，数组f是以π为单位的频带边缘频率，且0＜f＜l，数组g为在f所给频率上的幅度响应。f和g数组的长度必须相同且是偶数，可得到希尔伯特滤波器时域脉冲响应h’＝[h’₀，h’₁，h’₂，...，h'_L-1]。

(5)设计原型滤波器

用matlab中firl函数设计一个长度为K的原型低通滤波器H₀(z)，即通过运行fir1(K'，W_n)，其中K'＝K-1，为原型低通滤波器阶数，W_n是以π为单位的归一化截止频率，可得到原型低通滤波器时域脉冲响应h＝[h₀，h₁，h₂，...，h_K']。

(6)产生时延修正参考信号

时延修正参考信号r(n)为参考信号x(n)经由纯时延滤波器z^-Δ滤波得到，计算公式为r(n)＝x(n-Δ)，其中Δ为估计次级路径直达声时延采样点数。

(7)产生正交时延修正参考信号

正交时延修正参考信号由参考信号x(n)通过希尔伯特滤波器得到，在希尔伯特滤波器通带频率范围内，与r(n)各频率点对应相位之间差值为90°。计算公式为：

$\hat{r}\left(n\right)=\underset{i=0}{\overset{L-1}{\mathrm{\Σ}}}{h}_i^{\text{'}}x(n-i)$

式中，L为希尔伯特滤波器长度，h’_i为希尔伯特滤波器的系数，且设计的希尔伯特滤波器的延时与直达声时延Δ相同。

(8)产生子带参考信号

第m个子带时延修正参考信号r_m(k)由时延修正参考信号r(n)经多相快速傅里叶变换得到，采用同样的方法可得到第m个子带误差信号e_m(k)和相应的子带正交时延修正参考信号计算公式分别为：

$r_m\left(k\right)=\underset{i=0}{\overset{K_L-1}{\mathrm{\Σ}}}{h}_i{e}^{j2\mathrm{\π}\frac{\mathrm{mi}}{M}}r(\mathrm{kD}-i)$

$e_m\left(k\right)=\underset{i=0}{\overset{K_L-1}{\mathrm{\Σ}}}{a}_i{e}^{j2\mathrm{\π}\frac{\mathrm{mi}}{M}}r(\mathrm{kD}-i)$

${\hat{r}}_m\left(k\right)=\underset{i=0}{\overset{K_L-1}{\mathrm{\Σ}}}{a}_i{e}^{j2\mathrm{\π}\frac{\mathrm{mi}}{M}}\hat{r}(\mathrm{kD}-i)$

式中，k为子带时间指数，h_i是长度为K的原型FIR滤波器的第i个元素值，D为降采样率，为简化文字，子带时延修正参考信号和子带正交时延修正参考信号简称子带参考信号和子带正交参考信号。

(9)子带滤波参考信号和更新方向的选择

步骤1.不更新自适应滤波器系数，获得最大误差信号幅度，迭代计算P个点的子带误差信号功率和ξ₁、子带参考信号功率和χ₁，子带正交参考信号功率和计算公式为：

e_max=max(e_max，|e_m(k)|)

ξ₁＝ξ₁+|e_m(k)|²

χ₁＝χ₁+|r_m(k)|²

${\widehat{\mathrm{\χ}}}_1={\widehat{\mathrm{\χ}}}_1+{\left|{\hat{r}}_m\left(k\right)\right|}^2$

步骤2.对时延修正参考信号r_m(k)和子带正交时延修正参考信号两种子带参考信号类型，分别在μ_m＝μ₀或-μ₀两种更新方向，重复下面步骤，寻找合适更新方向。

(a)迭代计算P个点的子带误差信号功率和ξ₂、子带参考信号功率和χ₂以及子带正交参考信号功率和

(b)如果|e_m(k)|＜(1+δ)e_max，采用下式更新自适应滤波器系数；否则停止更新。

$w_l(k+1)=w_l\left(k\right)-2{\mathrm{\μ}}_m\mathrm{Re}\left\{e_m\left(k\right)R_m^{*}(k-l)\right\}$

式中，Re表示取实部运算。表示R_m(k)的共轭，为子带时延修正参考信号r_m(k)或子带正交时延修正参考信号

步骤3.如果满足|e_m(k)|＞(1+δ)e_max或ξ₂/χ₂＞δ₁ξ₁/χ₁、其中δ＜1，δ₁＜1，则跳转至步骤5；

步骤4.若μ_m＝μ0和-μ₀，都满足|e_m(k)|＞(1+δ)e_max或ξ₂/χ₂＞δ₁ξ₁/χ₁，则μ_m＝μ_mα，α＜1，并返回步骤(1)；否则，μ_m＝μ_m/α，并返回步骤(1)；

步骤5.比较步骤2中四种情况下所得子带误差信号功率ξ₂，ξ₂最小值对应合适的子带参考信号类型以及更新方向；

步骤6.重复执行步骤1至步骤5，直至m＝M-1；

(10)所有子带信号一起调整全带自适应滤波器系数，直至误差降低达到稳态，更新全带自适应滤波器W(z)系数采用如下迭代公式：

$w_l(k+1)=w_l\left(k\right)-2\underset{m=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_m\mathrm{Re}\left\{e_m\left(k\right)R_m^{*}(k-l)\right\}$

式中，μ_m为第m个子带步长，w_l为自适应滤波器W(z)系数，l为0，1，2……N，N为自适应滤波器阶数，m取0，1，…，M-1，M为子带个数，Re表示取实部运算。表示R_m(k)的共轭，为子带时延修正参考信号r_m(k)或子带正交时延修正参考信号

3.有益效果

本发明与现有技术相比，其显著优点：(1)将全带参考信号通过一希尔伯特滤波器得到全带正交时延修正参考信号，该信号与全带时延参考信号在频带范围内，相位差为90°。(2)通过系统总时延，粗略估计次级路径直达声时延，并将直达声时延作为次级路径传递函数的一部分先进行补偿，然后对时延修正参考信号、正交时延修正参考信号和误差信号做子带划分，以减少无次级路径有源控制方法应用中子带个数。(3)无需次级路径建模，每个子带在0°或180°方向切换搜索合适更新方向，并选择合适的子带参考信号类型，调整全带自适应滤波器系数。(4)避免了频域运算，整个算法都在时域内实施，且无需采用权变换将子带自适应滤波器系数堆砌为全带自适应滤波器系数，系统结构简单，易于实现，减少计算量。

四.附图说明

图1算法结构框图

图2 MSE曲线比较

图3仿真降噪前后频谱

图4实验降噪前后频谱

图5步长为0.002时三个子带更新方向搜索图

图6三个子带正确更新方向图

五.具体实施方式

下面以两端端口均封闭的管道声场中实现简化子带结构无次级路径有源控制方法为例，对该方法的具体实施方式做详细说明。

1.系统构成

实验装置是端口直径为17cm的方形直管道，管道两端端口封闭，误差传声器位于次级声源的下游管道中，初级声源置于次级声源的上游管道中，测得初级声源与误差传声器之间的距离为136cm，次级声源与误差传声器之间的距离l为40cm。控制器选用ADSP21161EZKIT数字信号处理评估板，采样频率为48000Hz。图1为简化子带结构无次级路径的有源控制方法框图，x(n)为DSP生成参考信号，经过功率放大器驱动初级声源，实现过程将误差传声器采集到的信号降采样到2400Hz，为100Hz～500Hz的带限白噪声，误差传感器用于获取误差信号e(n)。

2.实施过程

(1)将控制器的输入端口和输出端口短接，以白噪声作为控制器的输入信号，采用LMS算法，测得控制器从输入端到输出端延迟采样点数。如AD21161评估板，降采样后，采样频率为2400Hz，测得控制器从输入端到输出端共延迟15个采样点，可得控制器电时延τ₁为6.25ms。(2)计算实验中声时延τ₂，计算公式为τ₂＝l/c，其中l为次级声源与误差传感器之间的距离，c为声速，一般环境中取340m/s。例如次级声源与误差传感器之间的距离l为40cm，计算声时延τ₂近似为1.2ms。

(3)基于步骤(1)电路板电时延τ₁和步骤(2)声时延τ₂，计算总的时延τ＝τ₁+τ₂，直达声时延Δ＝τf_s，f_s为采样频率。例如控制器电时延τ₁为6.25ms，τ₂近似为1.2ms，计算总时延τ为7.45ms，f_s为2400Hz，对应的直达声时延Δ为17。

(4)设计希尔伯特滤波器

用matlab中firls函数设计一个长度为L的希尔伯特滤波器即通过运行firls(L’，f，g，’hilbert')，其中L'＝L-1，为希尔伯特滤波器阶数，数组f是以π为单位的频带边缘频率，且0＜f＜1，数组g为在f所给频率上的幅度响应。f和g数组的长度必须相同且是偶数，可得到希尔伯特滤波器时域脉冲响应h’＝[h'₀，h'₁，h'₂，...，h'_L-1]。例如，L可取34，a可取0.01，b可取0.95，其中c和d均取1。采样频率f_s可取2400Hz，通过运行firls(60，[0.010.95]，[11]，’hilbert')，可得到通带频率范围为12Hz到1140Hz的希尔伯特滤波器时域脉冲响应h’＝[h'₀，h’₁，h’₂，...，h'_L-1]，该希尔伯特滤波器时延为L/2，即为17个采样点。

(5)设计原型滤波器

用matlab中firl函数设计一个长度为K的原型低通滤波器H₀(z)，即通过运行firl(K'，W_n)，其中K'＝K-1，为原型低通滤波器阶数，W_n是以π为单位的归一化截止频率，可得到原型低通滤波器时域脉冲响应h＝[h₀，h₁，h₂，...，h_K']。对于两端封闭的管道系统，K可取128，归一化截止频率W_n可取为1/8，若采样频率f_s取1200Hz，通过运行firl(127，1/8)，则可得到截止频率为200Hz的原型低通滤波器时域脉冲响应h＝[h₀，h₁，h₂，...，h_K-1]。

(6)产生时延修正参考信号

时延修正参考信号r(n)为参考信号x(n)经由纯时延滤波器z^-Δ滤波得到，计算公式为r(n)＝x(n-Δ)，其中Δ为估计次级路径直达声时延采样点数。

(7)产生正交时延修正参考信号

正交时延修正参考信号由参考信号x(n)通过希尔伯特滤波器得到，在希尔伯特滤波器通带频率范围内，与r(n)各频率点对应相位之间差值为90°。计算公式为：

$\hat{r}\left(n\right)=\underset{i=0}{\overset{L-1}{\mathrm{\Σ}}}{h}_i^{\text{'}}x(n-i)$

式中，L为希尔伯特滤波器长度，h'_i为希尔伯特滤波器的系数，且设计的希尔伯特滤波器的延时与直达声时延Δ相同。

(8)产生子带参考信号

第m个子带时延修正参考信号r_m(k)由时延修正参考信号r(n)经多相快速傅里叶变换得到，采用同样的方法可得到第m个子带误差信号e_m(k)和相应的子带正交时延修正参考信号计算公式分别为：

$r_m\left(k\right)=\underset{i=0}{\overset{K_L-1}{\mathrm{\Σ}}}{h}_i{e}^{j2\mathrm{\π}\frac{\mathrm{mi}}{M}}r(\mathrm{kD}-i)$

$e_m\left(k\right)=\underset{i=0}{\overset{K_L-1}{\mathrm{\Σ}}}{a}_i{e}^{j2\mathrm{\π}\frac{\mathrm{mi}}{M}}r(\mathrm{kD}-i)$

${\hat{r}}_m\left(k\right)=\underset{i=0}{\overset{K_L-1}{\mathrm{\Σ}}}{a}_i{e}^{j2\mathrm{\π}\frac{\mathrm{mi}}{M}}\hat{r}(\mathrm{kD}-i)$

式中，k为子带时间指数，h_i是长度为K的原型FIR滤波器的第i个元素值，D为降采样率，子带时延修正参考信号和子带正交时延修正参考信号以下简称子带参考信号和子带正交参考信号。例如两端封闭管道系统中，原型滤波器长度K为128，实验中采用子带个数M为3，因此D取值不超过M/2，D取值为1。

(9)子带滤波参考信号和更新方向的选择

P取值为2000，初始步长μ₀为0.002，α取值为0.5，子带个数M为3，具体实施如下：

步骤1.不更新自适应滤波器系数，获得最大误差信号幅度，迭代计算P个点的子带误差信号功率和ξ₁、子带参考信号功率和χ₁和子带正交参考信号功率和计算公式为

e_max＝max(e_max，|e_m(k)|)

ξ₁＝ξ₁+|e_m(k)|²

χ₁＝χ₁+|r_m(k)|²

${\widehat{\mathrm{\χ}}}_1={\widehat{\mathrm{\χ}}}_1+{\left|{\hat{r}}_m\left(k\right)\right|}^2$

步骤2.对时延修正参考信号r_m(k)和子带正交时延修正参考信号两种子带参考信号类型，分别在μ_m＝μ₀或-μ₀两种更新方向情况下，重复下面步骤，寻找合适更新方向。

(a)迭代计算P个点的子带误差信号功率和ξ₂、子带参考信号功率和χ₂以及子带正交参考信号功率和

(b)如果|e_m(k)|＜(1+δ)e_max，采用下式更新自适应滤波器系数；否则停止更新。

$w_l(k+1)=w_l\left(k\right)-2{\mathrm{\μ}}_m\mathrm{Re}\left\{e_m\left(k\right)R_m^{*}(k-l)\right\}$

式中，Re表示取实部运算。表示R_m(k)的共轭，为子带时延修正参考信号r_m(k)或子带正交时延修正参考信号

步骤3.如果满足|e_m(k)|＞(1+δ)e_max或ξ₂/χ₂＞δ₁ξ₁/χ₁、其中δ＜1，δ₁＜1，则跳转至步骤5；

步骤4.若μ_m＝μ₀和-μ₀，都满足|e_m(k)[＞(1+δ)e_max或&/χ₂＞δ₁ξ₁/χ₁，则μ_m＝μ_mα，α＜1，并返回步骤(1)；否则，μ_m＝μ_m/α，并返回步骤(1)；

步骤5.比较步骤2中四种情况下所得子带误差信号功率ξ₂，ξ₂最小值对应合适的子带参考信号类型以及更新方向；

步骤6.重复执行步骤1至步骤5，直至m＝M-1。

(10)所有子带信号一起调整全带自适应滤波器系数，直至误差降低达到稳态，更新全带自适应滤波器W（z)系数采用如下迭代公式：

$w_l(k+1)=w_l\left(k\right)-2\underset{m=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_m\mathrm{Re}\left\{e_m\left(k\right)R_m^{*}(k-l)\right\}$

式中，μ_m为第m个子带步长，w_l为自适应滤波器W(z)系数，l为0，1，2……N，N为自适应滤波器阶数，m取0，1，…，M-1，M为子带个数，Re表示取实部运算。表示R_m(k)的共轭，为子带时延修正参考信号r_m(k)或子带正交时延修正参考信号

3.实验结果

本发明方法(简称SSANC方法)与FxLMS方法的时域收敛曲线，如图2所示，达到稳态时，所提SSANC方法与FxLMS方法最小均方误差(MSE)一致，然而SSANC算法需要依次搜索每个子带合适的更新方向并判断每个子带参考信号类型，然后进入控制过程，因而占用一定时间，因此与FxLMS算法相比较，收敛速度稍慢。通过仿真和具体实验验证，两种方法降噪前后信号频谱比较图，如图3和图4所示，说明实验和仿真结果相一致，且SSANC方法能取得很好的降噪效果。对实验中选取小的步长，如μ＝0.002，SSANC方法中三个子带需要重复搜索判断过程，需要较长的时间，搜索到正确更新方向和子带参考信号类型，如图5所示。当步长μ＝0.02时，SSANC方法能在较短时间内找到三个子带正确更新方向和子带参考信号类型，如图6所示。

Claims (3)

1.一种简化子带结构无次级路径有源控制方法，其特征为：

(1)时域实现简化子带结构无次级路径有源控制方法，通过估计次级路径直达声时延，并将该时延作为次级路径的一部分，产生时延修正参考信号，为估计次级路径90°和-90°相位信息，采用希尔伯特变换获得正交时延修正参考信号，并对时延修正参考信号、正交时延修正参考信号以及误差信号做子带分解，且全频带自适应滤波器不做子带划分，时延修正参考信号r(n)为：

r(n)＝x(n-Δ)

其中，x(n)为参考信号，Δ为估计次级路径直达声时延采样点数；

(2)无需次级路径建模，依次判断每个子带参考信号类型，并搜索该子带合适的更新方向；

(3)更新全带自适应滤波器系数：

$w_l(k+1)=w_l\left(k\right)-2\underset{m=0}{\overset{M-1}{\Σ}}{\mathrm{\μ}}_m\mathrm{Re}\left\{e_m\left(k\right)R_m^{*}(k-l)\right\}$

其中，μ_m为第m个子带步长，k为子带时间指数，w_i为自适应滤波器W(z)系数，l为0，1，2......N，

N为自适应滤波器阶数，M为子带个数，Re表示取实部运算；表示R_m(k)的共轭，R_m(k)可为子带时延修正参考信号r_m(k)或子带正交时延修正参考信号

2.如权利要求1所述的一种简化子带结构无次级路径有源控制方法，其特征在于通过希尔伯特变换获得正交时延修正参考信号，从而对时域无次级路径有源控制方法中引入±90°相位，正交时延修正参考信号为：

$\hat{r}\left(n\right)=\underset{i=0}{\overset{L-1}{\Σ}}{g}_{i}x(n-i)$

其中，x(n)为参考信号，L为希尔伯特滤波器阶数，g_i为希尔伯特滤波器系数。

3.如权利要求1所述的一种简化子带结构无次级路径有源控制方法，其特征在于通过子带正交时延修正参考信号，实现时域无次级路径有源控制方法中次级路径90°和-90°相位要求，并在时域无次级路径有源控制方法实施每个子带次级路径四个方向相位搜索，判断每个子带参考信号类型以及搜索子带合适更新方向的方法为：分别对子带时延修正参考信号r_m(k)和子带正交时延修正参考信号在0°和180°两种可能更新方向中搜索判断，从而选择合适的更新方向以及相应子带参考信号类型。