MATLAB-巴特沃兹滤波器.doc

相位是信号的一个重要参数，所以一般在滤波中希望实现零相位偏移或具有线性相位偏移。现在大致有两类方法实现零相位偏移：1、将滤波器与一个全通滤波器级联，用后者的相位响应矫正前者。需要附加全通滤波器的设计运算量以及额外的矫正算法运算量；2、使用两个完全相同的滤波器，将通过第一个滤波器后的信号反转再通过第二个滤波器，再将第二个滤波器的输出序列反转作为整个系统的输出。综上我们采用第二种方法，用MATLAB实现：clearclcclose all%生成一个信号Fs=1000;t=(0:1/Fs:1);y=sin(2*pi*5*t+pi/4)+sin(2*pi*50*t+pi/5)+sin(2*pi*60*t+pi/12); yy=sin(2*pi*5*t+pi/4);figure(1)subplot(211)plot(t,y),grid ontitle(输入信号波形)xlabel(t)ylabel(magnitude)subplot(212)plot(t,angle(y)title(输入信号相位)xlabel(t)ylabel(angle)%设计一个低通滤波器fp=10/500;%通带截止频率fs=30/500;%阻带截止频率wp=3;%通带最大衰减ws=30;%阻带最大衰减n,fc=buttord(fp,fs,wp,ws);b,a=butter(n,fc);figure(4)freqz(b,a);%将生成的信号通过低通滤波器z=filter(b,a,y);figure(2)subplot(221)plot(t,z)title(通过第一个滤波器后的波形)xlabel(t)ylabel(output)axis(0 1 -1 1)subplot(222)plot(t,angle(z)title(通过第一个滤波器后的相位)xlabel(t)ylabel(angle)subplot(223)plot(t,yy)title(理论上过滤后的波形)xlabel(t)ylabel(output)axis(0 1 -1 1)subplot(224)plot(t,angle(yy);title(原相位)xlabel(t);ylabel(angle)zf=fliplr(z); %对经过第一个滤波器后的信号序列反转zfz=filter(b,a,zf); %对经过第一个滤波器后的信号序列反转再通过滤波器zz=fliplr(zfz); %对经过两次滤波器后的输出反转作为最后输出figure(3)subplot(221)plot(t,zz)title(通过两个滤波器后的波形)xlabel(t)ylabel(output)axis(0 1 -1 1)subplot(222)plot(t,angle(zz)title(通过两个滤波器后的相位)xlabel(t)ylabel(angle)subplot(223)plot(t,yy)title(理论上过滤后的波形)xlabel(t)ylabel(output)axis(0 1 -1 1)subplot(224)plot(t,angle(yy);title(原相位)xlabel(t);ylabel(angle)%波形的前端和后端有失真，所以对信号进行延拓for i=1:1001 ytleft(i)=y(1002-i); ytright(i)=y(1002-i);endyt=ytleft y ytright;tt=(0:3/(length(yt)-1):3);figure(5)subplot(211)plot(tt,yt)title(通过镜像延拓后的输入信号波形)xlabel(t)ylabel(input)subplot(212)plot(tt,angle(yt)title(通过镜像延拓后的输入信号的相位)xlabel(t)ylabel(angle)%将生成的信号通过低通滤波器z=filter(b,a,yt);figure(6)subplot(221)plot(tt,z)title(通过第一个滤波器后的波形)xlabel(t)ylabel(output)axis(1 2 -1 1)subplot(222)plot(tt,angle(z)title(通过第一个滤波器后的相位)xlabel(t)ylabel(angle)axis(1 2 0 4)subplot(223)plot(t,yy)title(理论上过滤后的波形)xlabel(t)ylabel(output)axis(0 1 -1 1)subplot(224)plot(t,angle(yy);title(原相位)xlabel(t);ylabel(angle)axis(0 1 0 4)zf=fliplr(z); zfz=filter(b,a,zf); %对经过第一个滤波器后的信号序列翻转再通过滤波器zz=fliplr(zfz);figure(7)subplot(221)for i=1:1001output(i)=zz(1001+i);endplot(t,output)title(通过两个滤波器后的波形)xlabel(t)ylabel(output)axis(0 1 -1 1)subplot(222)plot(t,angle(output)title(通过两个滤波器后的相位)xlabel(t)ylabel(angle)axis(0 1 0 4)subplot(223)plot(t,yy)title(理论上过滤后的波