离散信号时域分析的MALAB实现.docx
《离散信号时域分析的MALAB实现.docx》由会员分享，可在线阅读，更多相关《离散信号时域分析的MALAB实现.docx（10页珍藏版）》请在冰豆网上搜索。
离散信号时域分析的MALAB实现
实验十一离散信号时域分析的MATLAB实现
一、实验目的
1.熟悉MATLAB编程方法、常用语句和可视化绘图技术；
2.掌握序列时域运算的MATLAB编程方法。
二、实验原理
在用MATLAB表示离散信号并将其可视化时，由于矩阵元素个数是有限的，因此无法表示无限长序列；另外，离散信号无法进行符号运算。
在MATLAB中，绘制离散序列波形图的专用命令为stem（）。
其格式有：
（1）stem（k,f）
在图形窗口中，绘制出样值顶部为空心圆的序列f（k）波形图。
（2）stem（k,f,’fill’）
在图形窗口中，绘制出样值顶部为实心圆的序列f（k）波形图。
下面介绍离散序列的MATLAB表示、基本运算（相加、相乘、平移、反转、尺度变换）、卷积和的实现及其图形显示方法。
1.单位序列（k）
单位序列的定义：
下面为绘制（k-k0）波形图的子程序：
functionimpseq（k1,k2,k0）%单位序列（k-k0）,k0为时移量
图11-1
k=k1:
k2;%k1,k2为序列的起止序列号
n=length（k）;
x=zeros（1,n）;
x（1,k0-k1+1）=1;%在k0时刻信号赋值为1
stem（k,x,'fill'）
axis（[k1,k2,0,1.1]）
title（'单位序列d（k-k0）'）
输入如下命令，则可获得单位序列（k-3）的波形图，如图11-1所示。
impseq（-1,5,3）
2.单位阶跃序列（k）
单位序列的定义：
下面为绘制（k-k0）波形图的MATLAB子程序。
functionstepseq（k1,k2,k0）%单位阶跃序列,k0为时移量
k=k1:
k0-1;%k1,k2为序列的起止序列号
kk=length（k）;
图11-2
x=zeros（1,kk）;%k0前信号赋值为0
stem（k,x,'fill'）%绘出k1~k0-1的波形（0值）
holdon
n=k0:
k2;
nn=length（n）;
x=ones（1,nn）;%k0后信号赋值为1
stem（n,x,'fill'）%绘出k0~k2的波形（1值）
holdoff
axis（[k1,k2,0,1.1]）
title（'单位阶跃序列'）
运行如下命令，则可获得单位序列（k-3）的波形图，如图11-2所示。
stepseq（-1,10,3）
3.序列的相加（减）、相乘运算
对序列向量f1（k）、f2（k）相加或相乘，可以通过补零的方式使f1（k）、f2（k）成为具有相同维数的序列向量s1（k）、s2（k），然后对s1（k）、s2（k）相加或相乘。
因此，序列向量f1（k）、f2（k）的维数可以不同。
以下函数可实现序列向量f1（k）、f2（k）的相加或相乘运算。
function[f,k]=sigadd（f1,k1,f2,k2）
%实现序列f1,f2的相加,相减,相乘,可据实际需要作选择
%f1,k1;f2,k2是参加运算的序列向量及其时间向量
%f,k作为返回的和（差,积）序列及其时间向量
%将f1,f2转换成等长序列s1,s2
k=min（min（k1）,min（k2））:
max（max（k1）,max（k2））;
s1=zeros（1,length（k））;s2=s1;%初始化序列
s1（find（（k>=min（k1））&（k<=max（k1））==1））=f1;
s2（find（（k>=min（k2））&（k<=max（k2））==1））=f2;
f=s1+s2;%序列相加
%f=s1-s2;%序列相减
%f=s1.*s2;%序列相乘
stem（k,f,'fill'）
axis（[（min（min（k1）,min（k2））-1）,（max（max（k1）,max（k2））+1）,（min（f）-0.5）,（max（f）+0.5）]）
例11-1.已知序列
。
编写M文件求
。
解：
运行如下M文件，可实现
，结果如图11-1所示。
图11-1
k1=-1:
3;
f1=[-32312];
k2=-1:
1;
kk=length（k2）;
f2=ones（1,kk）;
subplot（2,2,1）;stem（k1,f1,'fill'）;title（'f1（k）'）;
subplot（2,2,2）;stem（k2,f2,'fill'）;title（'f2（k）'）;
subplot（2,2,3）;[f,k]=sigadd（f1,k1,f2,k2）;title（'f1（k）+f2（k）'）
若要实现序列f1（k）、f2（k）的相乘或相减运算，只需将xlyunxuan（f1,k1,f2,k2）子程序中的相乘或相减语句设置为有效即可。
4.序列的平移、反转
（1）序列的平移
序列的平移可以看作是将序列的时间序号向量平移，而对应原时间序号的序列样值不变。
要将序列左移k0个单位时，则将时间序号向量都减小k0个单位；若要右移k0个单位时，则将时间序号向量都增大k0个单位。
实现序列平移的子函数如下：
function[x,n]=sigshift（f,k,k0）%实现序列平移:
x（k）=f（k-k0）
n=k+2;x=f;
subplot（1,2,1）;stem（k,f,'fill'）;title（'f（k）'）;xlabel（'k'）;
subplot（1,2,2）;stem（n,x,'fill'）;title（'f（k-k0）'）;xlabel（'k'）;
例11-2.已知指数序列
，绘出f（k-2）的波形图。
解：
运行如下M文件，可得如图11-2所示的结果。
图11-2
k=0:
5;
f=（0.5）.^k;%定义序列f（k）
[x,k]=sigshift（f,k,2）%调用平移子函数
（2）序列的反转
序列的反转可用MATLAB中的fliplr函数实现。
以下是实现序列反转及其结果可视化的函数。
function[x,n]=sigfold（f,k）%实现序列反转:
x（k）=f（-k）
x=fliplr（f）;
n=-fliplr（k）;
subplot（1,2,1）;stem（k,f,'fill'）;title（'f（k）'）;xlabel（'k'）;
subplot（1,2,2）;stem（n,x,'fill'）;title（'f（-k）'）;xlabel（'k'）;
5.序列的卷积运算
序列f1（k）、f2（k）的卷积和运算f（k）=f1（k）*f2（k），可由MATLAB的conv（）函数实现，调用格式为：
f=conv（f1,f2）
如：
已知序列：
，运行如下M文件可求其卷积和：
k1=-2:
2;
f1=ones（1,length（k1））;
k2=0:
3;
f2=2.^k2;
f=conv（f1,f2）
结果为：
f=
137151514128
可见，conv（）函数不需要给定f1（k）、f2（k）的非零样值的时间序号，也不返回卷积和序列f（k）的时间序号；此外，conv（）假定f1（k）、f2（k）都是从k=0开始，这就限制了它的应用范围。
因此，要对从任意k值开始的序列进行卷积和运算，同时正确标识出函数conv（）的计算结果各量f，还须构造序列f1（k）、f2（k）和f（k）的对应序号向量。
下面是求序列f1（k）、f2（k）卷积和的实用函数dconv（），它可实现序号向量的返回。
function[f,k]=dconv（f1,k1,f2,k2）%求卷积和:
f（k）=f1（k）*f2（k）
f=conv（f1,f2）
k0=k1
（1）+k2
（1）;%计算序列f非零样值的起点位置k0
k3=length（k1）+length（k2）-2;%计算序列f非零样值的宽度
k=k0:
k0+k3;%确定序列f非零样值的序号向量
subplot（2,2,1）;stem（k1,f1,'fill'）;title（'f1（k）'）;xlabel（'k'）;
subplot（2,2,2）;stem（k2,f2,'fill'）;title（'f2（k）'）;xlabel（'k'）;
subplot（2,2,3）;stem（k,f,'fill'）;title（'f（k）=f1（k）*f2（k）'）;xlabel（'k'）;
h=get（gca,'position'）;
h（3）=2.5*h（3）;
set（gca,'position',h）
例11-3.已知序列：
，求其卷积和。
图11-3
解：
运行如下M文件
k1=-2:
2;
f1=0.5.*k1;%定义序列f1（k）
k2=0:
3;
f2=2.^k2;%定义序列f2（k）
[f,k]=dconv（f1,k1,f2,k2）;%求卷积和
图形结果则如图11-3所示，文本结果如下：
f=
Columns1through6
-1.0000-2.5000-5.0000-9.5000-2.00004.0000
Columns7through8
8.00008.0000
6.离散信号相关函数的计算
在数字信号处理中，广泛用到信号的相关运算。
两个序列的互相关用来度量这两个序列的相似程度。
给定两个长度相同、能量受限的序列x（k）和y（k），它们的互相关函数是另一个序列，定义为：
其中，n称为滞后参数。
当x（k）=y（k）时，称为x（k）的自相关函数，它提供了序列位置不同基准情况下自相似程度的度量。
由于信号的自相关函数与信号的功率谱密度是一对傅里叶变换对，这为信号的功率谱计算提供了另一重要途径。
由互相关的定义式可知，x（k）和y（k）的互相关函数可以由序列x（n）和y（-n）的卷积和求得，实现互相关函数correlation（x,k1,y,k2）运算的程序如下：
function[Rxy,n]=correlation（x,k1,y,k2）
[y0,k0]=sigfold（y,k2）;%y（k）反转,生成y（-k）
[Rxy,n]=dconv（x,k1,y0,k0）;
holdoff
subplot（3,1,1）;stem（k1,x,'fill'）;xlabel（'k'）;ylabel（'x（k）'）;
subplot（3,1,2）;stem（k2,y,'fill'）;xlabel（'k'）;ylabel（'y（k）'）;
subplot（3,1,3）;stem（n,Rxy,'fill'）;xlabel（'n'）;ylabel（'Rxy（n）'）;
例11-4.设
为原型序列，y（k）为x（k）加入噪声干扰并移位后的序列：
y（k）=x（k-1）+w（k）
其中，w（k）为具有零均值和单位方差的高斯序列。
计算序列y（k）和x（k）的互相关。
图11-4
解：
运行如下M文件，可得如图11-4的结果。
k1=-3:
3;
x=[2351-131];%生成x（k）
[y2,k2]=sigshift（x,k1,1）;%生成x（k-1）
w=randn（1,length（y2））;%模拟噪声信号w（k）
kw=k2;
[y,k]=sigadd（y2,k2,w,kw）;%生成y（k）=x（k-1）+w（k）
[Rxy,n]=correlation（x,k1,y,k）%互相关运算
三.实验内容与步骤
1.编写M文件，绘制复指数序列
在
范围内x（k）的实部、虚部、模值和相角的波形图。
2.已知序列f1（k）如图11-5（A）所示，编写M文件，绘出如下序列的波形：
（1）f1（-k+2），
（2）f1（k-2）（k-1）。
3.已知序列f1（k）、f2（k）如图11-5所示，编写M文件，求离散卷积和f1（k）*f2（k），并绘出其波形。
4.编写M文件，求如图11-5所示序列f1（k）、f2（k）的互相关函数。
升级会员 