信号与系统部分实验matlab代码供参考解析.docx
《信号与系统部分实验matlab代码供参考解析.docx》由会员分享，可在线阅读，更多相关《信号与系统部分实验matlab代码供参考解析.docx（31页珍藏版）》请在冰豆网上搜索。
信号与系统部分实验matlab代码供参考解析
信号与系统第三版112学期部分实验代码
实验一常用信号的实现和时域变换部分实验答案
（1）
……………………………………………………………………………………………
t=0:
0.1:
30;
x=exp（-0.1*t）.*sin（2/3*t）;
plot（t,x）;
axis（[030-11]）;
grid;
xlabel（'Time（sec）'）;
ylabel（'f（t）'）;
……………………………………………………………………………………………
对上面的f（t）进行抽样f（n）=f（t）|t=nT=x（nT）就可以得到抽样信号，将上述命令里的时间增量改为1s，plot（t,x）命令改为可以得到下图。
4．已知离散序列x（n）如图所示，试用Matlab编程绘出以下的离散序列波形。
（1）x（-n+3）
clc;
closeall;
clearall;
n=-3:
2
x=[333321];
[x1,n1]=xlpy（x,n,-3）;
[x2,n2]=xlfz（x1,n1）;
subplot（311）
stem（n,x,'filled'）;
axis（[min（n）-1,max（n）+1,min（x）-0.5,max（x）+0.5]）
subplot（312）
stem（n1,x1,'filled'）;
axis（[min（n1）-1,max（n1）+1,min（x1）-0.5,max（x1）+0.5]）
subplot（313）
stem（n2,x2,'filled'）;
axis（[min（n2）-1,max（n2）+1,min（x2）-0.5,max（x2）+0.5]）
（2）x（n-3）x（n+2）
clc;
closeall;
clearall;
n=-3:
2
x=[333321];
[x1,n1]=xlpy（x,n,3）;
[x2,n2]=xlpy（x,n,-2）;
nn=-5:
5;
x3=[00000x1];
x4=[x2zeros（1,5）];
x=x3.*x4;
subplot（311）
stem（nn,x3,'filled'）;
subplot（312）
stem（nn,x4,'filled'）;
subplot（313）
stem（nn,x,'filled'）;
实验二LTI时间系统的时域分析
6.对于因果和稳定的LTI系统，对于下列二阶微分方程确定其单位冲激响应是否是欠阻尼、过阻尼或临界阻尼，画出系统的h（t）和频率响应模的波特图。
（1）
（2）
（3）
（4）
……………………………………………………………………………………………
clc;
closeall;
clearall;
num=[1];
den=[144];
H=tf（num,den）;
bode（H）;
holdon;
num=[7];
den=[545];
H=tf（num,den）;
bode（H）;
holdon;
num=[1/37];
den=[545];
H=tf（num,den）;
bode（H）;
……………………………………………………………………………………………
实验四连续信号与系统的频域分析
1.利用那个fourier函数求下列信号的傅里叶变换F（jw），并用ezplot函数绘出其幅度谱和相位谱.
（1）
closeall;
clearall;
symstvwphaseimre;%定义变量t,v,w,phase,imre
f=exp（-3*abs（t））*sin（2*t）;%
Fw=fourier（f）;
subplot（311）;
ezplot（f）;%画-2*pi到2*pi内函数
axis（[-0.01200.5]）;
subplot（312）;
ezplot（abs（Fw））;
im=imag（Fw）;
re=real（Fw）;
phase=atan（im/re）;
subplot（313）;
ezplot（phase）;
axis（[-66-0.50.5]）;
（2）
closeall;
clearall;
symstvwphaseimre;%定义变量t,v,w,phase,imre
f=t*exp（-2*t）*sin（4*t）*sym（'Heaviside（t）'）;%
Fw=fourier（f）;
subplot（311）;
ezplot（f）;
subplot（312）;
ezplot（abs（Fw））;
title（'幅度谱'）;
im=imag（Fw）;
re=real（Fw）;
phase=atan（im/re）;
subplot（313）;
ezplot（phase）;
title（'相位谱'）;
（3）
closeall;
clearall;
symstvwphaseimre;%定义变量t,v,w,phase,imre
f=sin（pi*t）/（pi*t）*（sin（2*pi*（t-1））/（2*pi*（t-1）））;%
Fw=fourier（f）;
subplot（311）;
ezplot（f）;
axis（[-55-0.20.2]）;
subplot（312）;
ezplot（abs（Fw））;
title（'幅度谱'）;
im=imag（Fw）;
re=real（Fw）;
phase=atan（im/re）;
subplot（313）;
ezplot（phase）;
title（'相位谱'）;
2.用ifourier函数求下列傅里叶变换的逆变换，并画出其时域波形。
（1）
symstw
F=0.5*（sin（0.25*w）/（0.25*w））^2;
f=ifourier（F,t）
ezplot（f）
（2）
closeall
clearall;
clc;
symstw
F=cos（4*w+pi/3）;
f=ifourier（F,t）
ezplot（f）
结果：
f=
1/4*Dirac（t+4）+1/4*i*3^（1/2）*Dirac（t+4）+1/4*Dirac（t-4）-1/4*i*3^（1/2）*Dirac（t-4）
无波形图
定义如下Dirac冲激函数：
……………………………………………………………………………………………
functionf=Dirac（t）
f=10*（t==0）
……………………………………………………………………………………………
closeall
clearall;
clc;
symstw
F=cos（4*w+pi/3）;
f=ifourier（F,t）
ezplot（f）
figure;
t=-6:
0.001:
6;
f=1/4*Dirac（t+4）+1/4*i*3^（1/2）*Dirac（t+4）+1/4*Dirac（t-4）-1/4*i*3^（1/2）*Dirac（t-4）
plot（t,f）;
axis（[-66-1111]）;
……………………………………………………………………………………………
……………………………………………………………………………………………
（3）
closeall
clearall;
clc;
symstw
F=sym（'sin（3*（w-2*pi））/（w-2*pi）'）
f=ifourier（F,t）
ezplot（abs（f）,[-66]）;
结果：
F=sin（3*w）/（w-2*pi）
f=
1/4*exp（2*i*pi*（t+3））*Heaviside（t+3）-1/4*exp（2*i*pi*（t+3））*Heaviside（-t-3）-1/4*exp（2*i*pi*（t-3））*Heaviside（t-3）+1/4*exp（2*i*pi*（t-3））*Heaviside（-t+3）
3.已知信号f1（t）的波形如图所示，其傅里叶变换为F1（jw）。
现有信号f（t）=f1（t）*f1（t），试用Matlab求f（t）的傅里叶变换F（jw），分别绘出f1（t）和f（t）的时域波形及F1（w）和F（w）的频谱曲线，验证时域卷积定理。
……………………………………………………………………………………………
closeall;
clearall;
symstvwphaseimre;%定义变量t,v,w,phase,imre
f=sym（'Heaviside（t+1）-Heaviside（t-1）'）;%
Fw=fourier（f）;
subplot（221）;
ezplot（f）;gridon;
title（'f（t）'）;
holdon;
plot（[-1-1],[01]）;
plot（[11],[01]）;
axis（[-22-0.42]）;
holdoff;
subplot（222）;
FFw=maple（'convert',Fw,'piecewise'）;%将傅里叶变换的符号表达式转换为样条曲线
ezplot（abs（FFw））;gridon;
title（'f（t）幅度谱'）;
subplot（223）;
ff=ifourier（FFw*FFw）;
ezplot（ff）;gridon;
title（'f（t）*f（t）'）;
subplot（224）;
Fw1=FFw*FFw
FFww=maple（'convert',Fw1,'piecewise'）;%将傅里叶变换的符号表达式转换为样条曲线
ezplot（abs（FFww））;
title（'f（t）*f（t）幅度谱'）;
gridon;
……………………………………………………………………………………………
4.要求产生一个时间从0到250ms的含有噪声的频率为50Hz和120Hz的时域信号，噪声的标准差为2，离散信号的时间间隔为2ms。
绘出信号的时域波形和频谱进行分析。
（建议用fft函数进行频谱分析）。
如何画正弦信号的频谱图
……………………………………………………………………………………………
clc;
closeall;
clearall;
ti=0.002;
t=0:
ti:
0.25;%时间范围从0到250ms
fs=500;%时间间隔2ms，对应采样频率为500Hz
y=sin（2*pi*50*t）;%正弦信号频率50Hz
subplot（2,1,1）;
plot（t,y）;
xlabel（'t'）;
ylabel（'y=sin（2*pi*50*t'）;
subplot（212）;
N=512;%FFT变换的点数，一般取2的n次方
Y=fft（y,N）;%对离散序列
f=-0.5*fs:
fs*1/N:
0.5*fs-1/N;%
plot（f,fftshift（abs（Y）））;%
ylabel（'频谱'）;
xlabel（'f'）;
……………………………………………………………………………………………
画法一：
……………………………………………………………………………………………
clearall
closeall
t=0:
0.001:
0.25;
x=sin（2*pi*50*t）+sin（2*pi*120*t）;%正弦信号加了噪声
y=x+2*randn（size（t））;
figure
（1）
plot（y,'g'）%画图
fs=1000;
f=-0.5*fs:
fs*1/256:
0.5*fs-fs/256;
ffty=fft（y,256）;%快速傅里叶变换
figure
（2）
plot（f,fftshift（abs（ffty）.^2））;%画正弦信号频谱
……………………………………………………………………………………………
画法二：
……………………………………………………………………………………………
closeall;
clearall;
t=0:
0.001:
0.25;
x=sin（2*pi*50*t）+sin（2*pi*120*t）;
y=x+2*randn（size（t））;
figure;
plot（y（1:
50））;
title（'NoisyTimedomainsignal'）;
figure;
Y=fft（y,256）;
pyy=Y.*conj（Y）/256;
f=1000/256*（0:
127）;
plot（f,pyy（1:
128））;
title（'powerspectraldensity'）
……………………………………………………………………………………………
5.图所示电路为二阶低通滤波器。
设
L=0.4H,C=0.05F,R=2Ω，试用Matlab编程绘出该系统频率响应
的幅频响应和相频响应曲线，并求出H（jw）的截止频率。
UR（t）
系统的过阻尼、欠阻尼和临界阻尼状态
%计算连续系统输出响应y''（t）+0.5y'（t）+2y（t）=x（t）,x（t）=u（t）实部为负共轭根—欠阻尼
%计算连续系统输出响应
a=[10.52];%微分方程左端系数
b=[1];%微分方程右端系数
%计算脉冲响应
impulse（b,a）
……………………………………………………………………………………………
%计算连续系统输出响应y''（t）+2y'（t）+4y（t）=x（t）,x（t）=u（t）一对重实根—临界阻尼
%计算连续系统输出响应
a=[124];%微分方程左端系数
b=[1];%微分方程右端系数
%计算脉冲响应
impulse（b,a）
……………………………………………………………………………………………
%计算连续系统输出响应y''（t）+4y'（t）+3y（t）=x（t）,x（t）=u（t）一对互异实根—过阻尼
%计算连续系统输出响应
a=[143];%微分方程左端系数
b=[1];%微分方程右端系数
%计算脉冲响应
impulse（b,a）
……………………………………………………………………………………………
%计算连续系统输出响应y''（t）+3y（t）=x（t）,x（t）=u（t）一对纯虚根—零阻尼
%计算连续系统输出响应
a=[103];%微分方程左端系数
b=[1];%微分方程右端系数
%计算脉冲响应
impulse（b,a）
%计算连续系统输出响应y''（t）-0.5y'（t）+2y（t）=x（t）,x（t）=u（t）实部为正共轭根—负阻尼
%计算连续系统输出响应
a=[1-0.52];%微分方程左端系数
b=[1];%微分方程右端系数
%计算脉冲响应
impulse（b,a）
求e-tu（t）信号的频谱的方法。
1.……已知其傅立叶变换求模和相位……………………………………………………………
closeall;
clearall;
w=-50:
0.2:
50;
Y=1./（j*w+1）;
plot（w,abs（Y））;
figure;
plot（w,angle（Y））;
……………………………………………………………………………………………
2.………………利用fourier函数求其傅立叶变换…………………………………………
closeall;
clearall;
symst;
f=exp（-t）*sym（'Heaviside（t）'）;
F=fourier（f）;
ezplot（abs（F））;
figure;
%ezplot（angle（F））;thisiserror
re=real（F）;
im=imag（F）;
phase=atan（im/re）;
ezplot（phase）;
……………………………………………………………………………………………
3.…………利用傅立叶变换的定义式求其傅立叶变换………………………………………
closeall;
clearall;
symstwfft;%定义符号变量
f=exp（-t）*sym（'Heaviside（t）'）;%单边指数信号
%F=fourier（f）;%调用fourier函数计算傅立叶变换
ft=f*exp（-j*w*t）;%计算被积函数
F=int（ft,t,-Inf,Inf）;%计算傅立叶变换
ezplot（abs（F））;
figure;
%ezplot（angle（F））;thisiserror
re=real（F）;
im=imag（F）;
phase=atan（im/re）;
ezplot（phase）;
……………………………………………………………………………………………
4.…………利用FFT求其傅立叶变换…………………………………………………………
closeall;
clearall;
clc;
ti=0.1
t=0:
ti:
20;
f=exp（-t）;
N=256;
F=fft（f,N）;
fs=1/ti;
f=-0.5*fs:
fs/N:
0.5*fs-fs/N;
plot（f,fftshift（abs（F）））;
figure;
plot（f,（angle（F）））;%?
?
……………………………………………………………………………
实验七连续时间信号的采样与重构
二、实验任务
连续信号f（t）=cos（8*pi*t）+2*sin（40*pi*t）+cos（24*pi*t），
（1）计算该信号的奈奎斯特频率；
（2）以不同的采样频率对该信号进行采样，画出采样前和采样后信号的频谱，对比观察和分析信号临界采样、过采样和欠采样情况下，信号频谱有何变化；（3）经过理想抽样后得到抽样信号fs（t），通过理想低通滤波器后重构信号f（t）。
三、实验方法
1．确定f（t）的最高频率fm。
对于无限带宽信号，确定最高频率fm的方法：
设其频谱的模降到10-5左右时的频率为fm。
2．确定Nyquist抽样间隔TN。
选定两个抽样时间：
TSTN。
3．MATLAB的理想抽样为
n=-200:
200;nTs=n*Ts;或nTs=-0.04:
Ts:
0.04
4．抽样信号通过理想低通滤波器的响应
理想低通滤波器的冲激响应为
系统响应为
由于
所以
（1）
MATLAB计算为
ft=fs*Ts*wc/pi*sinc（（wc/pi）*（ones（length（nTs）,1）*t-nTs'*ones（1,length（t））））;
式
（1）即用fs（nTs）恢复y（t）的表达式，它表明，连续信号y（t）可以展开成正交采样函数（Sa函数）的无穷级数，该级数的系数等于采样值f（nTs）。
也就是说，若在采样信号fs（t）的每个样点处，画一个最大峰值为f（nTs）的Sa函数波形，那么其合成的波形就是原信号y（t）.
四、实验要求（画出6幅图）
1.当TS（1）在一幅图中画原连续信号f（t）和抽样信号fS（t）。
f（t）是包络线，fS（t）是离散信号。
（2）画出重构的信号y（t）。
（3）画出误差图，即error=abs（f（t）-y（t））的波形。
2.当TS>TN时同样可画出3幅图
……………………………………………………………………………………………
%时域采样定理
display（'奈奎斯特周期1/24s，Ts<1/24s过采样频谱不混叠；Ts>1/24s欠采样频谱混叠'）;
display（'pleaseinputthevalueofsampleperiod'）;
Ts=input（'Ts='）;
%绘制有限长余弦信号
t=1:
0.01:
2;
y=cos（8*pi*t）+1*sin（40*pi*t）+cos（24*pi*t）;
subplot（221）;
plot（t,y）;
%axis（[06-1.11.1]）;
xlabel（'t单位:
s','Fontsize',8）;
title（'f（t）'）;
%line（[06],[00],'color',[000]）;
%数值求解余弦信号的频谱
N=300;
W=2*pi*50;%50=1/dt
k=-N:
N;
w=k*W/N;
Y=0.01*y*exp（-j*t'*w）;%求f（t）的傅里叶变换F（jw）
Y=abs（Y）;
subplot（222）
plot（w/pi,Y）;
%axis（[-220pi*7+0.2]）;
title（'F（j\omega）'）;
xlabel（'\omega单位:
pi'）;
%采样后的余弦信号
subplot（223）;
plot（t,y,'b:
'）;%绘制包络
holdon;
t2=1:
Ts:
2;
y2=cos（8*pi*t2）+2*sin（40*pi*t2）+cos（24*pi*t2）;
stem（t2,y2）;
%axis（[06-1.11.1]）;
xlabel（'t单位:
s','Fontsize',8）;
title（'fs（t）'）;
holdoff;
%采样后余弦新的频谱
Y2=Ts*y2*exp（-j*t2'*w）;
Y2=abs（Y2）;
subplot（224）;
plot（w/pi,Y,'b'）%蓝色绘制原始信号频谱
xlabel（'\omega单位:
pi'）;
title（'Fs（j\omega）'）;
holdon;
plot（w/pi,Y2,'r'）;%红色绘制采样信号频谱
%axis（[-220pi*10]）;
holdoff;
……………………………………………………………………………………………
Ts=0.02过采样
Ts=0.05欠采样
……………………………………………………………………………………………
%时域采样定理
clc;
closeall;
clearall;
display（'奈奎斯特周期1/24s，Ts<1/24s过采样频谱不混叠；Ts>1/24s欠采样频谱混叠'）;
display（'pleaseinputthevalueofsampleperiod'）;
Ts=input（'Ts='）;
%绘制有限长余弦信号
t=1:
0.01:
2;
y=cos（8*pi*t）+1*sin（40*pi*t）+cos（24*pi*t）;
subplot（221）;
plot（t,y）;
%axis（[06-1.11.1]）;
xlabel（'t单位:
s','Fontsize',8）;
title（'f（t）'）;
升级会员 