数字信号处理 实验报告 打印.docx
《数字信号处理 实验报告 打印.docx》由会员分享,可在线阅读,更多相关《数字信号处理 实验报告 打印.docx(18页珍藏版)》请在冰豆网上搜索。
数字信号处理实验报告打印
例一:
单位脉冲序列的产生
n=[-3:
3];xn=[(n-0)==0];%序列的起点为-3,终点为3,在n=0处生成一个单位脉冲;
subplot(3,2,1);stem(n,xn,'.');%stem:
离散序列绘图指令;
line([-3,3],[0,0])%创建直线;
axis([-3,3,0,1]);%控制坐标轴比例和外观;
xlabel('n');ylabel('x(n)')%X轴和Y轴标注;
例二:
单位阶跃序列的产生
n=[-10:
10];xn=[(n-0)>=0];%序列的起点为-10,终点为10,在n=0处生成一个单位阶跃;
subplot(3,2,1);stem(n,xn,'.');line([-10,10],[0,0])
axis([-10,10,0,1]);
xlabel('n');ylabel('x(n)')
例三:
矩形序列的产生
clear,closeall%清除所有程序;
n=[-10:
10];xn1=[(n-0)>=0];
xn2=[(n-5)>=0];%定义两个阶跃序列;
xn=xn1-xn2;%两个阶跃序列之差得到矩形序列;
subplot(3,2,1);stem(n,xn,'.');line([-10,10],[0,0])
axis([-10,10,0,1]);
xlabel('n');ylabel('x(n)')
例四:
幅度为1,脉冲宽度占空比duty=50,%周期T=5的周期方波脉冲信号的产生。
T=5,t=-2*T:
0.01:
2*T;duty=50;%duty—信号为正的区域在一个周期内所占的百分比;
x=square(t,duty);
subplot(2,2,1)
plot(t,x)
title('连续时间函数-周期方波脉冲')%title—标出图名;
axis([-10,10,-1.2,1.2])
line([-10,10],[0,0])
例五:
幅度为±1,周期T=4的三角波脉冲信号的产生
T=4;t=-2.5*T:
0.01:
2.5*T;
x=sawtooth(t-2*T/3,0.5);%sawtooth(x,width):
三角波函数,width为该位置的横坐标和周期的比值(0-1的标量)。
width=1时产生锯齿波,width=0.5时产生标准对称三角波;
subplot(3,3,1);
plot(t,x);
title('周期三角波脉冲,T=4')
axis([-10,10,-1.2,1.2]);
line([-10,10],[0,0])
实验二
clear;closeall;
N=64;n=[0:
N-1];
A=[1,-0.9];B=[0.05,0.05];%系统差分方程系数向量;
x=[n==0];%产生deta(n);
x1n=[11111111zeros(1,56)];%产生
x2n=ones(1,64);%产生
hn=filter(B,A,x);%filter求系统对序列x的响应;
subplot(2,2,1);
stem(n,hn,'.');
line([0,64],[0,0]);
axis([0,64,0,0.1]);
title('(a)系统单位脉冲响应h(n)')
y1n=filter(B,A,x1n)%求系统对序列x1n的响应;
subplot(2,2,2);
stem(n,y1n,'.');line([0,64],[0,0]);
axis([0,64,0,0.6]);
title('(b)系统对R8(n)的响应y1(n)')
y2n=filter(B,A,x2n);%求系统对序列x2n的响应;
subplot(2,2,4);
stem(n,y2n,'.');line([0,64],[0,0]);
axis([0,64,0,1]);
title('(c)系统对u(n)的响应y2(n)')
closeall;clearall
x1n=[11111111];N=8;M1=10;M2=4;
h1n=ones(1,10);
h2n=[12.52.51];
y1n=conv(x1n,h1n);%conv求卷积;
y2n=conv(x1n,h2n);
L1=N+M1-1;L2=N+M2-1;
n1=0:
L1-1;n2=0:
L2-1;
subplot(2,2,1);
stem(n1,y1n,'.');
line([0,20],[0,0]);
axis([0,20,0,10]);
title('(d)h1(n)与R8(n)的卷积');
subplot(2,2,2);
stem(n2,y2n,'.');
line([0,15],[0,0]);
axis([0,15,0,10]);
title('(e)h2(n)与R8(n)的卷积')
un=ones(1,256);%产生信号un
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);%谐振器对un的响应y31n
y32n=filter(B,A,xsin); %谐振器对正弦信号的响应y32n
figure(3)
subplot(2,1,1);stem(n,y31n,’.’)
title(′(f)谐振器对u(n)的响应y31(n)′)
subplot(2,1,2);stem(n,y32n,’.’)
title(′(g)谐振器对正弦信号的响应y32(n)′)
实验三
clearall;closeall;
x1n=ones(1,4);
N1=8;N2=16;k1=0:
N1-1;k2=0:
N2-1;w1=k1*2/N1;w2=k2*2/N2;%k1和k2是频域取样点的取值区间,w1和w2是将k值转换为归一化频率;
M=8;xa=1:
(M/2);xb=(M/2):
-1:
1;x2n=[xa,xb];%构造x2n
x3n=[xb,xa];
X1k8=fft(x1n,8);
X1k16=fft(x1n,16);
X2k8=fft(x2n,8);
X2k16=fft(x2n,16);
X3k8=fft(x3n,8);
X3k16=fft(x3n,16);
subplot(2,2,1);stem(w1,abs(X1k8),'.');
title('(1a)8点DFT[x1(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X1k8))]);
subplot(2,2,2);stem(w2,abs(X1k16),'.');
title('(1b)16点DFT[x1(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X1k16))]);
figure
(2)%创建图形窗
subplot(2,2,1);stem(w1,abs(X2k8),'.');
title('(1c)8点DFT[x2(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X2k8))]);
subplot(2,2,2);stem(w2,abs(X2k16),'.');
title('(1d)16点DFT[x2(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X2k16))]);
subplot(2,2,3);stem(w1,abs(X3k8),'.');
title('(1e)8点DFT[x3(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X3k8))]);
subplot(2,2,4);stem(w2,abs(X3k16),'.');
title('(1f)16点DFT[x3(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X3k16))]);
clearall;closeall;
N1=8;N2=16;k1=0:
N1-1;k2=0:
N2-1;w1=k1*2/N1;w2=k2*2/N2;n=0:
15;
x1n=cos(pi*n/4);
x2n=cos(pi*n/4)+cos(pi*n/8);
X1k8=fft(x1n,8);
X1k16=fft(x1n,16);
X2k8=fft(x2n,8);
X2k16=fft(x2n,16);
subplot(2,2,1);stem(w1,abs(X1k8),'.');
title('(2a)8点DFT[x1(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X1k8))]);
subplot(2,2,2);stem(w2,abs(X1k16),'.');
title('(2b)16点DFT[x1(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X1k16))]);
subplot(2,2,3);stem(w1,abs(X2k8),'.');
title('(2c)8点DFT[x2(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X2k8))]);
subplot(2,2,4);stem(w2,abs(X2k16),'.');
title('(2d)16点DFT[x2(n)]');
xlabel('w/pi');ylabel('幅度');
axis([0,2,0,1.2*max(abs(X2k16))]);
clearall;closeall;
Fs=64;T=1/Fs;
N=16;n=0:
N-1;
xnT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x(t)取样;
Xk16=fft(xnT,16);
Xk16=fftshift(Xk16);%将0频率移到频谱中心;
Tp=N*T;F=1/Tp;%求频率分辨率F
k=-N/2:
N/2-1;fk=k*F;
subplot(3,1,1);
stem(fk,abs(Xk16),'.');
title('(3a)16点|DFT[x(nT)]|');
xlabel('f(Hz)');ylabel('幅度');
axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(Xk16))]);
N=32;n=0:
N-1;
xnT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);
Xk32=fft(xnT,32);
Xk32=fftshift(Xk32);
Tp=N*T;F=1/Tp;
k=-N/2:
N/2-1;fk=k*F;
subplot(3,1,2);
stem(fk,abs(Xk32),'.');
title('(3b)32点|DFT[x(nT)]|');
xlabel('f(Hz)');ylabel('幅度');
axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(Xk32))]);
N=64;n=0:
N-1;
xnT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);
Xk64=fft(xnT,64);
Xk64=fftshift(Xk64);
Tp=N*T;F=1/Tp;
k=-N/2:
N/2-1;fk=k*F;
subplot(3,1,3);
stem(fk,abs(Xk64),'.');
title('(3c)64点|DFT[x(nT)]|');
xlabel('f(Hz)');ylabel('幅度');
axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(Xk64))]);
实验四
functionst=mstg%编写程序mstg,供后面调用;
N=1600;Fs=10000;T=1/Fs;Tp=N*T;%N为信号s(t)的长度;采样频率Fs=10kHz;
t=0:
T:
(N-1)*T;
k=0:
N-1;
f=k/Tp;
fc1=Fs/10;%第一路调幅信号的载波频率fc1=1kHz;
fm1=fc1/10;%第一路调幅信号的调制信号频率fm1=100Hz;
fc2=Fs/20;%第二路调幅信号的载波频率fc2=500Hz;
fm2=fc2/10%第二路调幅信号的调制信号频率fm2=50Hz;
fc3=Fs/40;%第三路调幅信号的载波频率fc3=250Hz;
fm3=fc3/10;%第三路调幅信号的调制信号频率fm3=25Hz;
xt1=cos(2*pi*fm1*t).*cos(2*pi*fc1*t);%产生第一路调幅信号;
xt2=cos(2*pi*fm2*t).*cos(2*pi*fc2*t);%产生第二路调幅信号;
xt3=cos(2*pi*fm3*t).*cos(2*pi*fc3*t);%产生第三路调幅信号;
st=xt1+xt2+xt3;%三路信号叠加形成时域混叠信号;
fxt=fft(st,N);%计算信号st的频谱
subplot(3,1,1);
plot(t,st);grid;%plot线性x-y坐标绘图;图上加坐标网格;
xlabel('t/s');ylabel('s(t)');
axis([0,Tp/2,min(st),max(st)]);
title('(a)s(t)的波形');
subplot(3,1,3);
stem(f,abs(fxt)/max(abs(fxt)),'.');grid;%stem离散序列绘图;
title('(b)s(t)的频谱');
axis([0,Fs/5,0,1.2]);
xlabel('f/Hz');ylabel('幅度');
clearall;closeall;
Fs=10000;T=1/Fs;M=1600;t=0:
T:
(M-1)*T;
st=mstg;%调用子程序mstg;
%==============以下为低通滤波器的设计与实现======================
fp=280;fs=450;
wp=2*fp/Fs;ws=2*fs/Fs;rp=0.1;rs=60;%获取数字滤波器的指标;
[N,wp]=ellipord(wp,ws,rp,rs);%计算椭圆DF的阶数N;
[B,A]=ellip(N,rp,rs,wp);%计算椭圆低通DF系统函数系数向量B和A;
[h1w,w]=freqz(B,A);%计算DF的频率响应
y1t=filter(B,A,st);%对信号st进行滤波
figure
(2);subplot(3,1,1);%figure创建图形窗;
plot(w/pi,20*log10(abs(h1w)));grid%grid图上加坐标网络;
axis([0,1,-100,10]);
xlabel('w/pi');ylabel('幅度(dB)');
subplot(3,1,3);
plot(t,y1t);grid
axis([0,0.1,-1.5,1.5]);
xlabel('t/s');ylabel('y1(t)');
%===============以下为带通滤波器的设计与实现=========================
fpl=440;fpu=560;fsl=275;fsu=900;rp=0.1;rs=60;
wp=[2*fpl/Fs,2*fpu/Fs];ws=[2*fsl/Fs,2*fsu/Fs];
[N,wp]=ellipord(wp,ws,rp,rs);
[B,A]=ellip(N,rp,rs,wp);
[h2w,w]=freqz(B,A);
y2t=filter(B,A,st);
figure(3);subplot(3,1,1);
plot(w/pi,20*log10(abs(h2w)));grid
axis([0,1,-100,10]);
xlabel('w/pi');ylabel('幅度(dB)');
subplot(3,1,3);
plot(t,y2t);grid
axis([0,0.1,-1.5,1.5]);
xlabel('t/s');ylabel('y2(t)');
%===============以下为高通滤波器的设计与实现====================
fp=890;fs=600;
wp=2*fp/Fs;ws=2*fs/Fs;rp=0.1;rs=60;
[N,wp]=ellipord(wp,ws,rp,rs);
[B,A]=ellip(N,rp,rs,wp,'high');
[h3w,w]=freqz(B,A);
y3t=filter(B,A,st);
figure(4);subplot(3,1,1);
plot(w/pi,20*log10(abs(h3w)));grid
axis([0,1,-100,10]);
xlabel('w/pi');ylabel('幅度(dB)');
subplot(3,1,3);
plot(t,y3t);grid
axis([0,0.1,-1.5,1.5]);
xlabel('t/s');ylabel('y3(t)');
实验五
functionxt=xtg(N)%信号x(t)产生函数,并显示信号的幅频特性曲线
%xt=xtg产生一个长度为N,有加性高频噪声的单频调幅信号xt,采样频率Fs=1khz
%载波频率fc=Fs/10=100Hz,调制正弦波频率f0=Fc/10=10Hz
N=2000;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;rs=70;%滤波器指标
fb=[fp,fs];m=[0,1];%计算remenzord函数所需参数f,m,dev
dev=[10^(-rs/20),(10^(rp/20)-1)/(10^(rp/20)+1)];
[n,fo,mo,w]=remezord(fb,m,dev,Fs);%确定remenz函数所需参数
hn=remez(n,fo,mo,w);%调用remenz函数进行设计,用于滤除噪声nt中的低频成分
yt=filter(hn,1,10*nt);%滤除随机噪声中低频成分,生成高频噪声yt
%===================================================================
xt=xt+yt;%噪声加信号
fst=fft(xt,N);k=0:
N-1;f=k/Tp;
subplot(3,2,1);plot(t,xt);grid;
xlabel('t/s');ylabel('x(t)');
axis([0,0.5,min(xt),max(xt)]);
title('(a)信号加噪声波形')
subplot(3,2,2);plot(f,abs(fst)/max(abs(fst)));grid;title('(b)信号加噪声的频谱')
axis([0,Fs/2,0,1.2]);xlabel('f/Hz');ylabel('幅度')
clearall;closeall;
%=====调用xtg产生信号xt,xt长度N=1000,并显示xt及其频谱===========
N=1000;xt=xtg(N);
fp=120;fs=150;rp=0.2;rs=60;%输入给定指标
Fs=1000;T=1/Fs;Tp=N*T;k=0:
N-1;f1=k/Tp;
%
(1)用窗函数法设计滤波器
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));
ywt=fftfilt(hn,xt,N);f=[f1,zeros(1,length(ywt)-1)];%调用函数fftfilt对xt滤波
t=0:
T:
[length(ywt)-1]*T;
subplot(3,2,5);plot(t,ywt);grid;title('(c)滤除噪声后的信号波形')
axis([0,0.5,-1,1]);xlabel('t/s');ylabel('yw(t)')
%=======
(2)用等波纹逼近法设计滤波器========================
fb=[fp,fs];m=[1,0];%确定remenzord函数所需参数f,m,dev
dev=[(10^(rp/20)-1)/(10^(rp/20)+1),10^(-rs/20)];
[Ne,fo,mo,W]=remezord(fb,m,dev,Fs);%确定remenz函数所需参数
hn=remez(Ne,fo,mo,W);%调用remenz函数进行设计
yet=fftfilt(hn,xt,N);%调用函数fftfilt对xt滤波
f=[f1,zeros(1,length(yet)-1)];
subplot(3,2,6);plot(t,yet);grid;title('(d)滤除噪声后的信号波形')
axis([0,0.5,-1,1]);xlabel('t/s');ylabel('ye(t)')
实验二
实验三
实验四
实验五