使用数字全通滤波器对IIR滤波器进行相位补偿

IIR滤波器在通带内不同频率的相移既不是线性的，也不是一个常值，使IIR滤波器的输出信号不同频率成分的延迟不可控制和调整。为了解决这个问题，可以在IIR滤波器后级联一个全通滤波器来调整输出信号的相移。

IIR全通滤波器函数如下：
函数：iirgrpdelay
功能：设计一个在已知频率范围内接近规定群延迟的全相位数字滤波器
调用格式：
[b,a]= iirgrpdelay(N,F,Edges,Gd)
说明：N是全通滤波器的阶数（N必须是偶数）；F和Gd是指在某一频率区间内群延迟的指标，F和Gd都是矢量，Edges是通带的边缘值；b和a是数字滤波器的系数。

案例、设计一个巴特沃斯数字低通滤波器,采样频率为4000Hz,滤波器通带400Hz，阻带600Hz，通带波纹和阻带衰减为Rp=3和Rs=20。在该滤波器之后级联一个全极点滤波器，使在通带内的时延接近于一个常数值。程序如下：

clear all; close all; clc;

fc=400;fb=600;                   % 设置通带和阻带频率
Rp=3;Rs=20;                      % 设置通带波纹和阻带衰减
Fs=4000; Fs2=Fs/2;               % 采样频率和奈奎斯特频率
Wp=fc/Fs2; Ws=fb/Fs2;            % 通带和阻带归一化频率
[N,Wn]=buttord(Wp,Ws,Rp,Rs);     % 设计巴特沃斯原型滤波器
[bn,an] = butter(N,Wn);          % 求出滤波器系数bn,an
[H1,w]=freqz(bn,an);             % 计算响应曲线
Hgd=grpdelay(bn,an);             % 计算群延迟曲线

F = 0:0.001:Wp;                  % 通带区间
g = grpdelay(bn,an,F,2);         % 求出通带群延迟
Gd = max(g)-g;                   % 给出一个反向群延迟值
% 设计一个IIR全通滤波器
[num,den,tau]=iirgrpdelay(4, F, [0 0.2], Gd);  

B=conv(num,bn);                  % 两滤波器级联后系数
A=conv(den,an);
[Ho,wo]=freqz(B,A);              % 计算级联滤波器响应曲线
[Hogd,wgd]=grpdelay(B,A);        % 计算级联滤波器群延迟曲线
% 作图
subplot 221; plot(w*Fs/2/pi,20*log10(abs(H1)),'k'); 
xlabel('频率/Hz'); ylabel('幅值/dB');
title('(a)巴特沃斯滤波器幅频响应'); axis([0 2000 -100 10]);

subplot 222; plot(w*Fs/2/pi,Hgd,'k'); xlim([0 2000]);
xlabel('频率/Hz率'); ylabel('延迟量/样点数');
title('(b)巴特沃斯滤波器群延迟')
subplot 223; plot(wo/pi*Fs2,20*log10(abs(Ho)),'k');
xlabel('频率/Hz'); ylabel('幅值/dB');
title('(c)级联滤波器幅频响应'); axis([0 2000 -100 10]);
subplot 224; plot(wgd/pi*Fs2,Hogd,'k');
xlabel('频率/Hz'); ylabel('延迟量/样点数');
title('(d)级联滤波器群延迟')
set(gcf,'color','w'); 

运行结果如下：

在程序中先设计了IIR的巴特沃斯低通滤波器，它的幅频响应曲线和群延迟曲线分别如图（a）和（b）所示。再进一步设计了IIR全通滤波器，它与低通滤波器级联。级联后滤波器幅频响应曲线和群延迟曲线如图（c）和（d）所示。从图中可看出，只有低通滤波器时群延迟在通带内为一条曲线（见图（b）），而级联了全通滤波器后，群延迟在通带内接近于一条直线（见图（d））。

在低通滤波器与全通滤波器级后群延迟只是接近于一个常数值，但它并不是零延迟，所以不能把这种级联方法作为零延迟来处理。

用全通滤波器和IIR滤波器级联后，相位可以得到补偿，群延迟也接近于一个常数值。但要记住，只是接近于常数值，有时还会有起伏，可以适当调整全通滤波器的阶数以进一步观察。

参考文献：MATLAB数字信号处理85个实用案例精讲——入门到进阶；宋知用（编著）