工程师&程序员的自我修养 No.3 啸叫信号与啸叫控制方法

 

在信号与系统Ⅱ的学习过程中，有约一个月的研讨汇报时间，供同学们在信号处理方面的知识进行拓展、研究，并总结自己的学习成果在课堂上进行汇报。此处我将总结一些自己学到的啸叫信号的产生原理及控制方法。并大致描述我的研讨内容。

---

一、什么是啸叫信号：

在平时跟别人用电脑语音聊天，声音传给你→再传给对方→再传给你，循环往复且不消失，这个就是啸叫了。专业一些的说法就是：

在声学场景中，当形成声反馈闭合回路时容易出现啸叫现象

 其原因是：话筒采集的声音信号包含了经过扬声器放大的声音信号，该信号在声反馈回路中不断叠加放大，正反馈产生震荡，产生啸叫。

 

二、啸叫形成原理：

其中，v是近端语音信号，u是扬声器聩给信号，y是传声器采集信号，x是经过传递函数作用产生的反馈信号。G是放音系统传递函数，C主要是A/D、D/A产生的时延，K是扬声器增益。

因而可以建立系统函数：

这里我们可以得到： 

要想产生自激振荡，可以通过自激振荡的产生条件对啸叫产生系统进行运用，得到幅频条件和相频条件如下所示：

三、啸叫信号的三种控制方式：

1、移频移相法：

利用频率也即相位上的特性，设计反馈系统抑制啸叫信号，这种方法会很大程度上破坏相位特性，引发失真[1]。

2、陷波抑制法：

 （1）陷波抑制法就是通过窄带滤波器/自适应滤波器进行特定频率的滤波，前提是找到这个频率，这就需要先进行啸叫检测：

（2）检测原则可以通过峰值/均值比等参数准则得出。检测出成分之后，利用陷波滤波器[3]：

3、自适应反馈抵消法：

因为扬声器的信号是已知的，这就是一个先验知识（也就是desired signal），从而可以利用adaptive filter，该类方法复杂度高[2]：

四、啸叫抑制的具体实现：

此处感谢wufankun在codeforge上分享的部分代码，给我了很多启示和思考。

1、移频移相法：

clear all
close all
clc

[x,fs] = audioread('man.wav');
x = x(1:fs*2);
g = load('path.txt');

h = fir2(200,[0,0.48,0.5,1],[1,1,0,0]); h = h(:); h = h.*exp(2*pi*1i*(1:length(h))'/4);
                                        % 希尔伯特变换所用滤波器，方法：先得到低通，然后移频
h_dummy = zeros(size(h)); h_dummy((end+1)/2) = 1;

K = 0.2;                                % 增益

g = g(:);                               % 反馈声学路径g
c = [0,0,0,0,1]';                       % 扩音系统内部传递路径c

xs1 = zeros(size(c));
xs2 = zeros(size(g));
xs3 = zeros(size(h_dummy));

y1 = zeros(size(x));                    % 先分配y1和y2空间，避免运行中临时分配空间占用大量的运算量
y2 = zeros(size(x));
temp = 0;

for i = 1:length(x)                     % 卷积形成反馈回路
    xs1 = [x(i)+temp; xs1(1:end-1)];    % 等待与c卷积的信号缓存
    y1(i) = K*(xs1'*c);                 % 馈给扬声器的信号
    
    xs3 = [y1(i); xs3(1:end-1)];        % 通过一个只有时延的滤波器，为了将y1与y2的群时延条件相同
    y1(i) = xs3' * h_dummy;
    
    y1(i) = min(1,y1(i));               % 幅度约束，啸叫则出现截止
    y1(i) = max(-1,y1(i));
    xs2 = [y1(i); xs2(1:end-1)];        % 等待与g卷积的信号缓存
    temp = xs2'*g;                      % temp作为单样点缓存，待下一采样点处理时与输入信号混合
end


xs1 = zeros(size(c));
xs2 = zeros(size(g));
xs3 = zeros(size(h));
temp = 0;
f_shift = 3;                            % 移频频率为3Hz

for i = 1:length(x)
    xs1 = [x(i)+temp; xs1(1:end-1)];
    y2(i) = K*(xs1'*c);
    
    xs3 = [y2(i); xs3(1:end-1)];
    y2(i) = xs3' * h;                   % 通过滤波器得到信号频谱的正半轴部分
    y2(i) = y2(i)*exp(2*pi*1i*i/fs*f_shift);    % 频移f_shift
    y2(i) = real(y2(i));                % 取实部，恢复出频谱在负半轴部分的信号

    y2(i) = min(1,y2(i));
    y2(i) = max(-1,y2(i));
    xs2 = [y2(i); xs2(1:end-1)];
    temp = xs2'*g;
end

figure
subplot(211),plot(y1)           % 未处理信号
subplot(212),plot(y2)           % 啸叫抑制信号

 可以明显看到通过一个简单的H移相器的作用，把让原混叠的功率尖峰大大缩减，最终得到一个啸叫情况极小的音频信号。这里有很大程度的失真，不难理解这不适用于实时啸叫控制系统的设计。

2、陷波抑制法：

clear all
close all
clc

[x,fs] = audioread('man.wav');
% x = x(1:fs);
g = load('path.txt');

K = 0.2;                                % 增益

g = g(:);                               % 反馈声学路径g
c = [0,0,0,0,1]';                       % 扩音系统内部传递路径c

xs1 = zeros(size(c));
xs2 = zeros(size(g));

y1 = zeros(size(x));                    % 先分配y1和y2空间，避免运行中临时分配空间占用大量的运算量
y2 = zeros(size(x));
temp = 0;

for i = 1:length(x)                     % 卷积形成反馈回路
    xs1 = [x(i)+temp; xs1(1:end-1)];    % 等待与c卷积的信号缓存
    y1(i) = K*(xs1'*c);                 % 馈给扬声器的信号
    y1(i) = min(1,y1(i));               % 幅度约束，啸叫则出现截止
    y1(i) = max(-1,y1(i));
    xs2 = [y1(i); xs2(1:end-1)];        % 等待与g卷积的信号缓存
    temp = xs2'*g;                      % temp作为单样点缓存，待下一采样点处理时与输入信号混合
end

d  = fdesign.notch('N,F0,Q,Ap',2,922/(fs/2),1,1);   % 设计陷波器，中心频率为922Hz
Hd = design(d);
iir_coef1 = Hd.sosMatrix;

d  = fdesign.notch('N,F0,Q,Ap',2,4534/(fs/2),1,1);  % 设计陷波器，中心频率为4534Hz
Hd = design(d);
iir_coef2 = Hd.sosMatrix;

iir_coef = [iir_coef1; iir_coef2];  % 两个陷波器级联成陷波器组


iir_buffer = zeros(size(iir_coef,1),5);

xs1 = zeros(size(c));
xs2 = zeros(size(g));
temp = 0;

for i = 1:length(x)
    xs1 = [x(i)+temp; xs1(1:end-1)];
    y2(i) = K*(xs1'*c);
    y2(i) = min(1,y2(i));
    y2(i) = max(-1,y2(i));
    d = y2(i);                      % iir滤波输入
    for k = 1:size(iir_coef,1)      % iir滤波过程
        iir_buffer(k,1:3) = [d, iir_buffer(k,1:2)];
        d = iir_buffer(k,1:3) * iir_coef(k,1:3)' - iir_buffer(k,4:5) * iir_coef(k,5:6)';
        iir_buffer(k,4:5) = [d, iir_buffer(k,4)];
    end
    y2(i) = d;
    xs2 = [y2(i); xs2(1:end-1)];
    temp = xs2'*g;
end

[b,a] = sos2tf(iir_coef,ones(size(iir_coef,1)+1,1));
figure,freqz(b,a,1:fs/2,fs)     % 观察陷波器组的频响

figure
subplot(211),plot(y1)           % 未处理信号
subplot(212),plot(y2)           % 啸叫抑制信号

陷波抑制法通过一个滤波器的作用让啸叫严重的频段大大衰减，从而得到我们想要的啸叫情况极弱的音频信号。这种方法失真是低于上面的移频移相法的。

---

参考文献：

 [1] 徐国伟, 王凤刚, 彭启琮. 频率采样法设计陷波器及其在啸叫抑制中的应用[J]. 西部广播电视, 2005, 000(010):20-21.

[2] 周璐. 影响自适应反馈抵消啸叫抑制算法性能的声学因素分析[D]. 2012.

[3] 张金亮, 王彦. 啸叫检测和抑制方法及其装置, 2016.