matlab-基于matlab的IIR数字陷波器设计和仿真-源码

clc; clear; close all; warning off; addpath(genpath(pwd)); %**************************************************************************** %更多关于matlab和fpga的搜索“fpga和matlab”的CSDN博客: %matlab/FPGA项目开发合作 %https://blog.csdn.net/ccsss22?type=blog %**************************************************************************** %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %IIR数字陷波器设计 i=[1 2 3 4 5 6]; %wn=[0.025*pi 0.075*pi 0.125*pi]; wn=[5/320*pi 15/320*pi 25/320*pi]; w=zeros(1,6); for j=1:6 k=floor((j+1)/2); w(j)=wn(k)-(1/2)*(1-(-1).^mod(j,2))*0.002*pi; end %w=[0.023*pi 0.027*pi 0.073*pi 0.077*pi 0.0123*pi 0.0127*pi]; %o1=[-0.5 -1.5 -2.5 -3.5 -4.5 -5.5]; o1=-(2*floor((i+1)/2)-1)+(1/2)*(1-(-1).^mod(i,2))*0.5; o=o1*pi; B=(o+4*w)/2; P=tan(B); %P=P'; Q=zeros(6,6); for ii=1:6 for kk=1:6 Q(kk,ii)=sin(kk*w(ii))-tan(B(ii))*cos(kk*w(ii)) end end Q1=pinv(Q) s=P*Q1; a=[2 2*s] b=[s(6)+1 s(5)+s(1) s(4)+s(2) s(3)+s(3) s(2)+s(4) s(1)+s(5) s(6)+1]; freqz(b,a) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% load c1.mat y=c; Fs=6400; figure(2) yy=filter(b,a,y); subplot(211) plot(y(1000+800:1511+800)) xlabel('(a)原波形') subplot(212) plot(yy(1000+800:1511+800)) xlabel('(b)陷波后') figure(8) M1=fft(y(1000+800:1255+800)) M=20*log10(abs(M1)); f1=Fs*(0:127)/256; plot(f1,M(1:128)) M2=fft(yy(1000+800:1255+800)) MM=20*log10(abs(M2)); hold on plot(f1,MM(1:128),'r') xlabel('陷波前后信号的频谱图') axis([0 1000 -50 40]) figure(3) N=128; A1=fft(yy(3688:3815),N); %A1=fft(yy(684:811),N); f=Fs*(-N/2:(N/2-1))/N; AA1=abs(fftshift(A1)); subplot(211) bar(f,AA1,'g'); xlabel('没有调制信号陷波后') %AA=A1.*conj(A1) %db1=10*log10(AA); %plot(f,db1,'r') A2=fft(yy(3815:3942),N); %A2=fft(yy(812:939),N); AA2=abs(fftshift(A2)); subplot(212) bar(f,AA2,'g'); xlabel('有调制信号陷波后') figure(4) %z1=yy(3815:3942)-yy(3688:3815) %z1=y(3432:3559)-y(3688:3815) z1=yy(940:1067)-yy(812:939) subplot(211) plot(z1) Z=fft(z1,N); Z2=abs(fftshift(Z)); subplot(212) bar(f,Z2,'g'); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% load c1.mat y1=c; Fs=6400; figure(5) yy1=filter(b,a,y1); subplot(211) plot(y1) subplot(212) plot(yy1) figure(6) z11=y1(3817:3944)-y1(3689:3816) %z1=y(3432:3559)-y(3688:3815) subplot(211) plot(z11) Z1=fft(z11,N); Z22=abs(fftshift(Z1)); subplot(212) bar(f,Z22,'g');