一到N阶全通系统的传输函数和相位特点

理解全通系统的传输函数和相位特点。

（1）设计并实现一个稳定的一阶全通IIR滤波器（零极点自行安排），再级联成N阶（N较大），计算此N阶全通滤波器的相位时延和群时延；讨论N的大小对时延的影响。

（2）产生一段音频（注意设置频率等参数），通过这个N阶滤波器，输出结果，分析相位延迟的影响。

设计：

一阶全通系统的频率响应：

 

其中，a为极点。

N阶全通系统由N个一阶全通系统级联而成，所以其频率响应为：

我们使用filter函数对音频信号进行滤波操作，对与N阶滤波器的系统响应函数的分子分母系数，我们采用二项式系数的计算方法来得到。

即：

结果分析：

实验中，我们设置a=A=0.5，即极点为0.5；

取一段14s的音频，采样频率为48kHz,我们截取8s的音频信号进行处理。

由

可知，N阶全通滤波器的相移函数为一阶全通滤波器的相移函数的N倍。群时延也是一阶全通滤波器的N倍。接下来，我们依次设置N的大小，来观察滤波器相时延和群时延函数的变化。

N=3；

N=10；

N=20；

频响函数没有变化，相频函数在随着N的增大成倍增加。

N=29;

在滤波器阶数达到29时，滤波器的相时延和群时延就开始有异常，与一阶全通滤波器的相位时延倍数不再N。

N=36;

此时，频谱函数没有改变，群时延和相时延异常，此时声音信号还没有完全失真。

N=37；

N为37 的时候，高阶全通滤波器相位和群时延都严重失真，此时，因为计算机的有限字长效应，频谱信号被全部滤除，音频信号也完全失真。

代码：

% ——————输入信号——————

clear all;

clc

N=37;%高阶滤波器阶数为N

A=0.5;%零极点，极点为A

filename = 'E:\音频.wav';

[x,Fs] = audioread(filename);

sample=[1,8.*Fs];

clear x Fs

[x,Fs]=audioread(filename,sample);

x(:,2)=[];

n1=length(x);

ww=[0:n1-1].*2.*pi./n1;

%——————filter滤波——————

k = 1:N;

%由于为n个一阶全通系统级联，所以分子分母涉及到二项式系数，保存到index里面；

index = [1 cumprod((N-k+1)./(k))];

for i=1:N+1

    a(i)=index(i)*((-A)^(N-i+1));%传递函数分子系数

    b(i)=index(i)*((-A).^(i-1));%传递函数分母系数

end

y1=filter([-A 1],[1 -A],x);

yn=filter(a,b,x);

Y1=fftshift(fft(y1));

X=fftshift(fft(x));

YN=fftshift(fft(yn));

%——————计算相位时延和群时延———————

[tp1,wp1]=phasedelay([-A,1],[1,-A],length(x));

[tpn,wpn]=phasedelay(a,b,length(x));

[tg1,wg1]=grpdelay([-A,1],[1,-A],length(x));

[tgn,wgn]=grpdelay(a,b,length(x));

%——————保存滤波音频——————

audiowrite('E:\god\七七的敢死队友\dsp\homework\第三次discuss\音频y1.wav',y1,Fs);

audiowrite('E:\god\七七的敢死队友\dsp\homework\第三次discuss\音频yn.wav',yn,Fs);

%————————画图————————

figure,subplot(221),plot(ww,abs(X)),title('原信号频谱图');

subplot(222),plot(ww,abs(Y1),'r'),title('一阶全通滤波信号频谱图');

subplot(223),plot(ww,abs(YN),'g'),title('N阶全通滤波信号频谱图');

figure,subplot(221),plot(ww,abs(x)),title('原信号时域图');

subplot(222),plot(ww,abs(y1),'r'),title('一阶全通滤波信号时域图');

subplot(223),plot(ww,abs(yn),'g'),title('N阶全通滤波信号时域图');

figure,plot(wp1,tp1,'r'),title('一阶全通系统相延迟','r');hold on;

plot(wpn,tpn,'g'),title('1（红）和N（绿）阶全通系统相延迟');

figure,plot(wg1,tg1,'r'),title('一阶全通系统群延迟');hold on;

plot(wgn,tgn,'g'),title('1（红）和N（绿）阶全通系统群延迟');