1， 2， 3， 4， 5， 6， 7， 8，
系统响应及系统稳定性 信号的时域采样与频域采样 用FFT对信号作频谱分析 FIR数字滤波器的特性
频率采样法设计FIR数字滤波器 FIR数字滤波器设计综合实验 双线性交换法设计IIR数字滤波器 IIR数字滤波器设计综合实验
实验八
%产生信号序列向量st，并显示st的时域波形和频谱
%st=mstg 返回三路调幅信号相加形成的混合信号，长度N=800 function st=mstg
N=800;%N为信号st的长度
Fs=10000;T=1/Fs;Tp=N*T;%采样频率Fs=10kHz，Tp=为采样时间 t=0:T:(N-1)*T;k=0:N-1;f=k/Tp;
fc1=Fs/10;%第1路调幅信号的载波信号fc1=1000Hz fm1=fc1/10;%第1路调幅信号的载波信号fm1=100Hz fc2=Fs/20;%第2路调幅信号的载波信号fc2=500Hz fm2=fc2/10;%第2路调幅信号的载波信号fm2=50Hz fc3=Fs/40;%第3路调幅信号的载波信号fc3=250Hz fm3=fc3/10;%第3路调幅信号的载波信号fm3=25Hz xt1=cos(2*pi*fm1*t).*cos(2*pi*fc1*t);%产生第1路调幅信号 xt2=cos(2*pi*fm2*t).*cos(2*pi*fc2*t);%产生第2路调幅信号 xt3=cos(2*pi*fm3*t).*cos(2*pi*fc3*t);%产生第3路调幅信号 st=xt1+xt2+xt3; %三路调幅信号相加 fxt=fft(st,N); %计算条幅信号st的频谱 %以下为绘图部分，绘制st的时域波形和幅频特性曲线 subplot(211),plot(t,st);grid;xlabel('t/s');ylabel('s(t)'); axis([0,Tp/8,min(st),max(st)]);title('(a)s(t)的波形');
subplot(212),stem(f,abs(fxt)/max(abs(fxt)),'.');grid;title('(b)s(t)的频谱') axis([0 Fs/5 0 1.2]); xlabel('f/Hz');ylabel('幅度');
clear all;close all;clf;clc; Fs=10000;T=1/Fs; %采样频率
%调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st st=mstg;N=800;k=(N-1)/Fs;t=(0:N-1)/Fs;%低通滤波器设计与实现
fp=280;fs=450;wp=2*fp/Fs;ws=2*fs/Fs; %计算滤波器边界频率关于pi归一化 rp=0.1;rs=60; [N,wp]=ellipord(wp,ws,rp,rs);
%调用ellipord计算椭圆DF阶数N和通带截止频率wp
[B,A]=ellip(N,rp,rs,wp); %调用ellip计算椭圆带通DF系统函数系数向量B和A y1t=filter(B,A,st); %数字滤波器软件实现w=0:T*pi:pi;[h,w]=freqz(B,A,w); % 低通滤波器设计与实现绘图部分
figure(2);subplot(2,1,1);plot(w/pi,20*log10(abs(h)),'k'); %调用绘图函数plot绘制损耗函数曲线 axis([0,1,-80,10]);xlabel('\\omega/\\pi');ylabel('幅度/dB');title('低通滤波器损耗函数曲线'); subplot(2,1,2);
plot(t,y1t); %调用绘图函数plot绘制滤波器输出波形
axis([0 k -1.2 1.2]);xlabel('t/s');ylabel('y_1(t)');title('低通滤波器滤除后的调幅信号y_1(t)'); %带通滤波器设计与实现fpl=440;fpu=560;fsl=275;fsu=900;
wp=[2*fpl/Fs,2*fpu/Fs];ws=[2*fsl/Fs,2*fsu/Fs];%滤波器边界频率关于pi归一化
rp=0.1;rs=60; [N,wp]=ellipord(wp,ws,rp,rs); %调用ellipord计算椭圆DF阶数N和通带截止频率wp [B,A]=ellip(N,rp,rs,wp); %调用ellip计算椭圆带通DF系统函数系数向量B和A y2t=filter(B,A,st); %滤波器软件实现 % 带通滤波器设计与实现绘图部分
w=0:T*pi:pi;[h,w]=freqz(B,A,w);figure(3);subplot(2,1,1);
plot(w/pi,20*log10(abs(h)),'k'); %调用绘图函数plot绘制损耗函数曲线
axis([0,1,-80,10]);xlabel('\\omega/\\pi');ylabel('幅度/dB');title('带通滤波器损耗函数曲线'); subplot(2,1,2);
plot(t,y2t); %调用绘图函数plot绘制滤波器输出波形
axis([0 k -1.2 1.2])xlabel('t/s');ylabel('y_2(t)');title('带通滤波器滤除后的调幅信号y_2(t)'); %高通滤波器设计与实现
fp=890;fs=600;wp=2*fp/Fs;ws=2*fs/Fs;rp=0.1;rs=60; %DF指标（低通滤波器的通、阻带边界频） [N,wp]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N,rp,rs,wp,'high');y3t=filter(B,A,st);w=0:T*pi:pi;
[h,w]=freqz(B,A,w);figure(4);subplot(2,1,1);plot(w/pi,20*log10(abs(h)),'k'); %调用绘图函数plot绘制损耗函数曲线 axis([0,1,-80,10]);
xlabel('\\omega/\\pi');ylabel('幅度/dB'); title('高通滤波器损耗函数曲线'); subplot(2,1,2);
plot(t,y3t); %调用绘图函数plot绘制滤波器输出波形 axis([0 k -1.2 1.2]); xlabel('t/s');ylabel('y_3(t)');
title('高通滤波器滤除后的调幅信号y_3(t)');
实验7 clear all;clc;clf; T=1;Fs=1/T; wpz=0.2;wsz=0.3; rp=1;rs=15;
Wp=2/T*tan(wpz*pi/2); Ws=2/T*tan(wsz*pi/2); [N,Wc]=buttord(Wp,Ws,rp,rs,'s'); [B,A]=butter(N,Wc,'s'); [Bz,Az]=bilinear(B,A,Fs); w=0:0.001*pi:pi; [h,w]=freqz(Bz,Az,w);
subplot(2,1,1);plot(w/pi,20*log10(abs(h)),'k'); grid on;axis([0,1,-200,10]);xlabel('\\omegal/\\pi'); ylabel('幅度/dB');title('低通滤波器（T=1s）'); [Nd,Wdc]=buttord(wpz,wsz,rp,rs); [Bdz,Adz]=butter(Nd,Wdc); [h,w]=freqz(Bdz,Adz,w);
subplot(212);plot(w/pi,20*log10(abs(h)),'k'); grid on;axis([0,1,-200,10]);
xlabel('\\omega/\\pi');ylabel('幅度/dB'); title('直接设计法（T=1s）');
clear all;clc;clf; T=2;Fs=1/T; wpz=0.2;wsz=0.3; rp=1;rs=15;
Wp=2/T*tan(wpz*pi/2); Ws=2/T*tan(wsz*pi/2); [N,Wc]=buttord(Wp,Ws,rp,rs,'s'); [B,A]=butter(N,Wc,'high','s'); [Bz,Az]=bilinear(B,A,Fs); w=0:0.001*pi:pi; [h,w]=freqz(Bz,Az,w);
subplot(2,1,1);plot(w/pi,20*log10(abs(h)),'k'); grid on;axis([0,1,-200,10]);xlabel('\\omegal/\\pi'); ylabel('幅度/dB');title('高通滤波器（T=2s）'); [Nd,Wdc]=buttord(wpz,wsz,rp,rs); [Bdz,Adz]=butter(Nd,Wdc,'high'); [h,w]=freqz(Bdz,Adz,w);
subplot(212);plot(w/pi,20*log10(abs(h)),'k'); grid on;axis([0,1,-200,10]);
xlabel('\\omega/\\pi');ylabel('幅度/dB'); title('直接设计法（T=2s）');
clc;clf;clear;clear all;
fsl=1500;fsu=2700;fpl=2025;fpu=2225; Fs=8000;T=1/Fs;rp=1;rs=40;
wsz=[2*fsl/Fs,2*fsu/Fs];wpz=[2*fpl/Fs,2*fpu/Fs]; wp=2/T*tan(wpz*pi/2);ws=2/T*tan(wsz*pi/2); [N,wpo]=ellipord(wp,ws,rp,rs,'s'); [B,A]=ellip(N,rp,rs,wpo,'s'); [Bz,Az]=bilinear(B,A,Fs); w=0:T*pi:pi; [h,w]=freqz(Bz,Az,w);
subplot(111),plot(w/pi,20*log10(abs(h)),'k'); grid on;
axis([0.3,0.8,-60,10]);xlabel('\\omega/\\pi'); set(gca,'xtick',[0.4 0.5 0.6 0.7]); ylabel('幅度/dB');title('带通滤波器');
ws=[2*fsl/Fs,2*fsu/Fs];wp=[2*fpl/Fs,2*fpu/Fs]; wp=2/T*tan(wpz*pi/2);ws=2/T*tan(wsz*pi/2); w=0:T*pi:pi;
[Nd,wpod]=ellipord(wp,ws,rp,rs,'s'); [Bz,Az]=ellip(Nd,rp,rs,wpod,'stop','s'); [Bdz,Adz]=bilinear(Bz,Az,Fs); [h1,w]=freqz(Bdz,Adz,w); figure(2)
subplot(111),plot(w/pi,20*log10(abs(h1)),'k'); grid on;axis([0.3,0.8,-60,10]); set(gca,'xtick',[0.3 0.4 0.5 0.6 0.7]); xlabel('\\omega/\\pi');ylabel('幅度/dB'); title('带阻滤波器');
实验6
function xt=xtg(N)
%实验五信号x(t)产生,并显示信号的幅频特性曲线
%xt=xtg(N) 产生一个长度为N,有加性高频噪声的单频调幅信号xt,采样频率Fs=1000Hz %载波频率fc=Fs/10=100Hz,调制正弦波频率f0=fc/10=10Hz. Fs=1000;T=1/Fs;Tp=N*T;t=0:T:(N-1)*T;
fc=Fs/10;f0=fc/10; %载波频率fc=Fs/10，单频调制信号频率为f0=Fc/10; mt=cos(2*pi*f0*t); %产生单频正弦波调制信号mt，频率为f0 ct=cos(2*pi*fc*t); %产生载波正弦波信号ct，频率为fc xt=mt.*ct; %相乘产生单频调制信号xt nt=2*rand(1,N)-1; %产生随机噪声nt
%=======设计高通滤波器hn,用于滤除噪声nt中的低频成分,生成高通噪声======= fp=150; fs=200;Rp=0.1;As=70;% 滤波器指标
fb=[fp,fs];m=[0,1]; %计算remezord函数所需参数f,m,dev dev=[10^(-As/20),(10^(Rp/20)-1)/(10^(Rp/20)+1)]; [n,fo,mo,W]=remezord(fb,m,dev,Fs);
% 确定remez函数所需参数
hn=remez(n,fo,mo,W);%调用remez函数进行设计,用于滤除噪声nt中的低频成分 yt=filter(hn,1,10*nt); %滤除随机噪声中低频成分，生成高通噪声yt
%================================================================ xt=xt+yt; %噪声加信号 fst=fft(xt,N);k=0:N-1;f=k/Tp;
subplot(4,1,1);plot(t,xt);grid;xlabel('t/s');ylabel('x(t)'); axis([0,Tp/5,min(xt),max(xt)]);title('(a) 信号加噪声波形')
subplot(4,1,2);plot(f,abs(fst)/max(abs(fst)));grid;title('(b) 信号加噪声的频谱'); axis([0,Fs/2,0,1.2]);xlabel('f/Hz');ylabel('幅度');
clear all;close all;clc;
%==调用xtg产生信号xt, xt长度N=1000,并显示xt及其频谱,========= N=1000;xt=xtg(N);
fp=120; fs=150;Rp=0.2;As=60;Fs=1000; % 用窗函数法设计滤波器
wc=(fp+fs)/Fs; %理想低通滤波器截止频率(关于pi归一化） B=2*pi*(fs-fp)/Fs; %过渡带宽度指标 Nb=ceil(11*pi/B); %blackman窗的长度N hn=fir1(Nb-1,wc,blackman(Nb)); Hw=abs(fft(hn,1024));
% 求设计的滤波器频率特性
ywt=fftfilt(hn,xt,N); %调用函数fftfilt对xt滤波 subplot(4,1,3);plot(Hw);axis([0 1000 -0.2 1.2]); subplot(4,1,4);plot(ywt);
% 输入给定指标
实验3
clear all;close all
%实验内容(1)=================================================== x1n=[ones(1,4)];%产生序列向量x1(n)=R4(n)
M=8;xa=1:(M/2); xb=(M/2):-1:1; x2n=[xa,xb];%产生长度为8的三角波序列x2(n) X1k8=fft(x1n,8); %计算x1n的8点DFT X1k16=fft(x1n,16); %计算x1n的16点DFT X2k8=fft(x2n,8); %计算x1n的8点DFT X2k16=fft(x2n,16); %计算x1n的16点DFT %以下绘制幅频特性曲线 clc;clf;clear;
k=0:1023;w=k/1024;subplot(3,1,1);plot(w,abs(fft(x1n,1024)));axis([0,1,0,1.2*max(abs(X1k8))]); set(gca,'Xtick',[0,0.25,0.5,0.75,1]);title('x1n的DFT');xlabel('wk');ylabel('幅度'); subplot(3,1,2);stem(abs(X1k8),'.'); %绘制8点DFT的幅频特性图 title('(1a) 8点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度');
axis([0,8,0,1.2*max(abs(X1k8))]);subplot(3,1,3); stem(abs(X1k16),'.'); %绘制16点DFT的幅频特性 title('(1b)16点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度');axis([0,16,0,1.2*max(abs(X1k16))]); figure(2)
k=0:1023;w=k/1024;subplot(3,1,1);plot(w,abs(fft(x2n,1024)));axis([0,1,0,1.2*max(abs(X1k8))]); set(gca,'Xtick',[0,0.25,0.5,0.75,1]);title('x2n的DFT');xlabel('w');ylabel('幅度'); subplot(3,1,2);stem(abs(X2k8),'.'); %绘制8点DFT的幅频特性图 title('(2a) 8点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度');
axis([0,8,0,1.2*max(abs(X2k8))]);subplot(3,1,3);stem(abs(X2k16),'.'); %绘制16点DFT的幅频特性 title('(2b)16点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度'); axis([0,16,0,1.2*max(abs(X2k16))]); N=8;n=0:N-1; %FFT的变换区间N=8 x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8); X4k8=fft(x4n); %计算x4n的8点DFT X5k8=fft(x5n); %计算x5n的8点DFT N=16;n=0:N-1; %FFT的变换区间N=16 x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8); X4k16=fft(x4n); %计算x4n的16点DFT X5k16=fft(x5n); %计算x5n的16点DFT k1=0:7;k2=0:15;figure(3)
subplot(2,2,1);stem(1/8*k1,abs(X4k8),'.'); %绘制8点DFT的幅频特性图 title('(4a) 8点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度');
axis([0,1,0,1.2*max(abs(X4k8))])subplot(2,2,3);stem(1/16*k2,abs(X4k16),'.'); %绘制16点DFT title('(4b)16点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度');axis([0,1,0,1.2*max(abs(X4k16))]) subplot(2,2,2);stem(1/8*k1,abs(X5k8),'.'); %绘制8点DFT的幅频特性图 title('(5a) 8点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度');
axis([0,1,0,1.2*max(abs(X5k8))])subplot(2,2,4);stem(1/16*k2,abs(X5k16),'.'); %绘制16点DFT的 title('(5b)16点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度'); axis([0,1,0,1.2*max(abs(X5k16))])
实验5
function[Hr,w,b,L]=Hr_type2(h) M=length(h);L=M/2; b=2*h(L:-1:1);n=[1:1:L]; n=n-0.5;w=[0:1:511]'*6*pi/511; Hr=cos(w*n)*b';
clear all;close all;clc;clf;clear; h=[-4,1,-1,-2,5,6,5,-2,-1,1,-4]; M=length(h);n=0:M-1; [Hr,w,a,L]=Hr_Type2(h); amax=max(a)+1;amin=min(a)-1; subplot(2,2,1);stem(n,h,'.'); axis([-1,2*L+1,amin,amax]); xlabel('n');ylabel('h(n)'); title('冲激响应');
subplot(2,2,3);stem(0:L,a'.'); axis([-1,2*L+1,amin,amax])
xlabel('n'),ylabel('a(n)');title('a(n) 系数'); subplot(2,2,2);plot(w/pi,Hr);grid;
xlabel('w(单位:π)');ylabel('Hr');title('1型振幅相应'); subplot(2,2,4); zplane(h,L'.'); clear all;close all;clc;clf; M=20;tao=(M-1)/2; k=0:M-1;w1=(2*pi/M)*k; Hrs=[1,1,1,zeros(1,15),1,1]; k1=0:floor((M-1)/2); k2=floor((M-1)/2)+1:M-1;
angH=[-tao*(2*pi)/M*k1,tao*(2*pi)/M*(M-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H,M));
[db,mag,pha,grd,w]=freqz_m(h,1); [Hr,ww,a,L]=Hr_type2(h);
Hdr=[1,1,0,0,1,1];wdL=[0,0.25,0.25,1.75,1.75,2]; subplot(2,2,1);plot(w1/pi,Hrs,'o',wdL,Hdr);
axis([0,2,-0.1,1.1]);grid on;set(gca,'xtick',[0,0.2,0.3,1.8,1.9,2]); title('M=20');
subplot(2,2,2);stem(k,h,'.');axis([-1,20,-0.1,0.3]); title('脉冲响应');
subplot(2,2,3);plot(ww/pi,Hr,w1(1:11)/pi,Hrs(1:11),'o'); axis([0,1,-0.2,1.2]);set(gca,'xtick',[0,0.2,0.3,1]); title('振幅响）');
subplot(2,2,4);plot(w/pi,db);
axis([0,1,-60,10]);set(gca,'xtick',[0,0.2,0.3,1]); title('幅度响应（db）');
实验1
close all;clear all
%======内容1：调用filter解差分方程，由系统对u(n)的响应判断稳定性====== A=[1,-0.9];B=[0.05,0.05]; %系统差分方程系数向量B和A x1n=[1 1 1 1 1 1 1 1 zeros(1,50)]; %产生信号x1(n)=R8(n) x2n=ones(1,128); %产生信号x2(n)=u(n) hn=impz(B,A,58); %求系统单位脉冲响应h(n)
subplot(2,2,1);y='h(n)';stem(hn, '.','b'); %调用函数stem绘图 title('(a) 系统单位脉冲响应h(n)');
y1n=filter(B,A,x1n); %求系统对x1(n)的响应y1(n) subplot(2,2,2);y='y1(n)';stem(y1n, '.','b'); title('(b) 系统对R8(n)的响应y1(n)');
y2n=filter(B,A,x2n); %求系统对x2(n)的响应y2(n) subplot(2,1,2);y='y2(n)';stem(y2n,'.', 'b'); title('(c) 系统对u(n)的响应y2(n)');
%===内容2：调用conv函数计算卷积============================ x1n=[1 1 1 1 1 1 1 1 ]; %产生信号x1(n)=R8(n) h1n=[ones(1,10) zeros(1,10)]; h2n=[1 2.5 2.5 1 zeros(1,10)]; y21n=conv(h1n,x1n); y22n=conv(h2n,x1n); figure(2)
subplot(2,2,1);y='h1(n)';stem(h1n, '.', 'b'); %调用函数tstem绘图 title('(d) 系统单位脉冲响应h1(n)'); subplot(2,2,2);y='y21(n)'; stem(y21n,'.', 'b'); title('(e) h1(n)与R8(n)的卷积y21(n)');
subplot(2,2,3);y='h2(n)'; stem(h2n, '.', 'b'); %调用函数tstem绘图 title('(f) 系统单位脉冲响应h2(n)'); subplot(2,2,4);y='y22(n)';stem(y22n,'.', 'b'); title('(g) h2(n)与R8(n)的卷积y22(n)');
%=========内容3：谐振器分析======================== un=ones(1,256); %产生信号u(n) n=0:255;
xsin=sin(0.014*n)+sin(0.4*n); %产生正弦信号
A=[1,-1.8237,0.9801];B=[1/100.49,0,-1/100.49]; %系统差分方程系数向量B和A y31n=filter(B,A,un); %谐振器对u(n)的响应y31(n) y32n=filter(B,A,xsin); %谐振器对u(n)的响应y31(n) figure(3)
subplot(2,1,1);y='y31(n)';stem(y31n, '.', 'b'); title('(h) 谐振器对u(n)的响应y31(n)'); subplot(2,1,2);y='y32(n)';stem(y32n, '.', 'b'); title('(i) 谐振器对正弦信号的响应y32(n)');
实验2
% 时域采样理论验证程序exp2a.m
Tp=64/1000; %观察时间Tp=64微秒 %产生M长采样序列x(n) % Fs=1000;T=1/Fs; Fs=1000;T=1/Fs; M=Tp*Fs;n=0:M-1;
A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5; xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); %M点FFT[xnt)] subplot(2,1,1);
stem(n,xnt,'.'); box on;title('(a) Fs=1000Hz'); k=0:M-1;fk=k/Tp;
subplot(2,,2);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=1000Hz'); xlabel('f(Hz)');ylabel('幅度');axis([0,Fs,0,1.2*max(abs(Xk))])
M=27;N=32;n=0:M; %产生M长三角波序列x(n)
xa=0:floor(M/2); xb= ceil(M/2)-1:-1:0; xn=[xa,xb];
Xk=fft(xn,1024); %1024点FFT[x(n)], 用于近似序列x(n)的TF X32k=fft(xn,32) ;%32点FFT[x(n)]
x32n=ifft(X32k); %32点IFFT[X32(k)]得到x32(n) X16k=X32k(1:2:N); %隔点抽取X32k得到X16(K) x16n=ifft(X16k,N/2); %16点IFFT[X16(k)]得到x16(n) subplot(3,2,2);stem(n,xn,'.');box on
title('(b) 三角波序列x(n)');xlabel('n');ylabel('x(n)');axis([0,32,0,20]) k=0:1023;wk=1*k/1024; %
subplot(3,2,1);plot(wk,abs(Xk));title('(a)FT[x(n)]'); xlabel('omega/pi');ylabel('|X(e^j^omega)|');axis([0,1,0,200]) k=0:N/2-1;
subplot(3,2,3);stem(k,abs(X16k),'.');box on
title('(c) 16点频域采样');xlabel('k');ylabel('|X_1_6(k)|');axis([0,16,0,200]) n1=0:N/2-1;
subplot(3,2,4);stem(n1,x16n,'.');box on
title('(d) 16点IDFT[X_1_6(k)]');xlabel('n');ylabel('x_1_6(n)');axis([0,32,0,20]) k=0:N-1;
subplot(3,2,5);stem(k,abs(X32k),'.');box on
title('(e) 32点频域采样');xlabel('k');ylabel('|X_3_2(k)|');axis([0,32,0,200]) n1=0:N-1;
subplot(3,2,6);stem(n1,x32n,'.');box on
title('(f) 32点IDFT[X_3_2(k)]');xlabel('n');ylabel('x_3_2(n)');axis([0,32,0,20])
实验4
function[Hr,w,a,L]=Hr_type1(h) M=length(h); L=(M-1)/2;
a=[h(L+1),2*h(L:-1:1)]; n=[0:1:L];
w=[0:1:511]'*pi/512; Hr=cos(w*n)*a';
clear all;close all;clc;clf; h=[-4,1,-1,-2,5,6,5,-2,-1,1,-4]; M=length(h); n=0:M-1;
[Hr,w,a,L]=Hr_type1(h); a,L
amax=max(a)+1;amin=min(a)-1; subplot(2,2,1);stem(n,h,'.'); axis([-1,2*L+1,amin,amax]); xlabel('n');ylabel('h(n)'); title('冲击响应');
subplot(2,2,3);stem(0:L,a,'.'); axis([-1,2*L+1,amin,amax]);
xlabel('n');ylabel('a(n)');title('a(n)系数') subplot(2,2,2);plot(w/pi,Hr);grid;
xlabel('w(单位：π)');ylabel('Hr');title('I型振幅响应'); subplot(2,2,4);zplane(h,1);
function[Hr,w,d,L]=Hr_type4(h) M=length(h); L=M/2; d=2*[h(L:-1:1)]; n=[1:1:L]; n=n-0.5;
w=[0:1:511]'*6*pi/511; Hr=sin(w*n)*d';
clear all;close all;clc;clf;
h=[-4,1,-1,-2,5,6,6,5,-2,-1,1,-4]; M=length(h); n=0:M-1;
[Hr,w,d,L]=Hr_type4(h); d,L
dmax=max(d)+1;dmin=min(d)-1; subplot(2,2,1);stem(n,h,'.'); axis([-1,2*L+1,dmin,dmax]); xlabel('n');ylabel('h(n)'); title('冲击响应');
subplot(2,2,3);stem(1:L,d,'.'); axis([-1,2*L+1,dmin,dmax]);
xlabel('n');ylabel('d(n)');title('d(n)系数') subplot(2,2,2);plot(w/pi,Hr);grid;
xlabel('w(单位：π)');ylabel('Hr');title('I型振幅响应'); subplot(2,2,4);zplane(h,1);
function[Hr,w,c,L]=Hr_type3(h) M=length(h); L=(M-1)/2;
c=[2*h(L+1:-1:1)]; n=[0:1:L];
w=[0:1:511]'*6*pi/511; Hr=sin(w*n)*c';
clear all;close all;clc;clf;
h=[-4,1,-1,-2,5,6,6,5,-2,-1,1,-4]; M=length(h); n=0:M-1;
[Hr,w,c,L]=Hr_type3(h); c,L
cmax=max(c)+1;cmin=min(c)-1; subplot(2,2,1);stem(n,h,'.'); axis([-1,2*L+1,cmin,cmax]); xlabel('n');ylabel('h(n)'); title('冲击响应');
subplot(2,2,3);stem(0:L,c,'.'); axis([-1,2*L+1,cmin,cmax]);
xlabel('n');ylabel('b(n)');title('b(n)系数') subplot(2,2,2);plot(w/pi,Hr);grid;
xlabel('w(单位：π)');ylabel('Hr');title('I型振幅响应'); subplot(2,2,4);zplane(h,1);
function[Hr,w,b,L]=Hr_type2(h) M=length(h); L=M/2;
b=2*[h(L:-1:1)]; n=[1:1:L]; n=n-0.5;
w=[0:1:511]'*6*pi/511; Hr=cos(w*n)*b';
clear all;close all;clc;clf;
h=[-4,1,-1,-2,5,6,6,5,-2,-1,1,-4]; M=length(h); n=0:M-1;
[Hr,w,b,L]=Hr_type2(h); b,L
bmax=max(b)+1;bmin=min(b)-1; subplot(2,2,1);stem(n,h,'.'); axis([-1,2*L+1,bmin,bmax]); xlabel('n');ylabel('h(n)'); title('冲击响应');
subplot(2,2,3);stem(1:L,b,'.'); axis([-1,2*L+1,bmin,bmax]);
xlabel('n');ylabel('b(n)');title('b(n)系数') subplot(2,2,2);plot(w/pi,Hr);grid;
xlabel('w(单位：π)');ylabel('Hr');title('I型振幅响应'); subplot(2,2,4);zplane(h,1);
因篇幅问题不能全部显示，请点此查看更多更全内容