数字信号处理a双语dspa实验报告.docx
《数字信号处理a双语dspa实验报告.docx》由会员分享，可在线阅读，更多相关《数字信号处理a双语dspa实验报告.docx（73页珍藏版）》请在冰豆网上搜索。
数字信号处理a双语dspa实验报告
实验报告
实验名称：
1、离散时间信号与系统的时、频域表示
2、离散傅立叶变换和z变换
3、数字滤波器的频域分析和实现
4、数字滤波器的设计
课程名称数字信号处理A（双语）
班级学号
姓名
指导老师
实验教室
开课时间2016-2017年学年第一学期
实验一：
离散时间信号与系统的时、频域表示
一、实验目的和任务：
●主要内容：
序列的生成、序列DTFT谱计算、简单滤波器的仿真，对应实验指导书Page1－34。
●实验要求：
●
（1）按照要求产生基本序列和复杂序列，对序列进行运算，求序列的DTFT幅度谱和相位谱。
●
（2）计算简单的离散时间系统的输出，包括简单滤波器的仿真、计算单个系统和互连系统的输出。
二、实验内容：
●具体包括
●基本序列产生和运算：
Q1.1－Q1.3，Q1.23，Q1.30－Q1.33
●离散时间系统仿真：
Q2.1－Q2.3
●LTI系统：
Q2.19、Q2.21，Q2.28
●DTFT：
Q3.1，Q3.2，Q3.4
三、实验过程与结果分析：
Q1.1运行P1_1产生单位样本序列u[n]的程序与显示的波形如下:
clf;
n=-10:
20;
u=[zeros（1,10）1zeros（1,20）];
stem（n,u）;
xlabel（'时间序列n'）;ylabel（'振幅'）;
title（'单位样本序列'）;
axis（[-102001.2]）;
结果：
Q1.2
clf清除波形axis设置坐标轴范围，可读比例等title给图形加标题xlabel给x加标注
Ylabel给y加标注
Q1.3
clf;
n=-10:
20;
u=[zeros（1,10）1zeros（1,20）];
stem（n+11,u）;
xlabel（'时间序列n'）;ylabel（'振幅'）;
title（'单位样本序列'）;
axis（[03201.2]）;
结果：
Q1.23
n=0:
50;
f=0.08;
phase=pi/2;
A=2.5;
arg=2*pi*f*n-phase;
x=A*cos（arg）;
clf;
stem（n,x）;
axis（[050-33]）;
grid;
title（'正弦序列'）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
axis;
结果：
Q1.30
s[n]是线性增加伴随着实指数缓慢衰减的图像加性噪声d[n]是均匀分布在-0.4和+0.4之间的随机序列
Q1.31
不能，因为d是列向量，s是行向量
Q1.32
x1是x的延时，x2和x相等，x3超前于x
Q1.33
产生图例说明
Q2.1
clf;
n=0:
100;
s1=cos（2*pi*0.05*n）;
s2=cos（2*pi*0.47*n）;
x=s1+s2;
M=input（'滤波器所需的长度='）;
num=ones（1,M）;
y=filter（num,1,x）/M;
subplot（2,2,1）;
plot（n,s1）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'信号#1'）;
subplot（2,2,2）;
plot（n,s2）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'信号#2'）;
subplot（2,2,3）;
plot（n,x）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'输入信号'）;
subplot（2,2,4）;
plot（n,y）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'输出信号'）;
axis;
结果：
S[n]被离散时间系统抑制
Q2.2
n=0:
100;
s1=cos（2*pi*0.05*n）;
s2=cos（2*pi*0.47*n）;
x=s1+s2;
M=input（'滤波器所需长度='）;
num=（-1）.^[0:
M-1];
y=filter（num,1,x）/M;
clf;
subplot（2,2,1）;
plot（n,s1）;
axis（[0,100,-2,2]）;
xlabel（'时间序号n'）;ylabel（'振幅'）;
title（'信号#1'）;
subplot（2,2,2）;
plot（n,s2）;
axis（[0,100,-2,2]）;
xlabel（'时间序号n'）;ylabel（'振幅'）;
title（'信号#2'）;
subplot（2,2,3）;
plot（n,x）;
axis（[0,100,-2,2]）;
xlabel（'时间序号n'）;ylabel（'振幅'）;
title（'输入信号'）;
subplot（2,2,4）;
plot（n,y）;
axis（[0,100,-2,2]）;
xlabel（'时间序号n'）;ylabel（'振幅'）;
title（'输出信号'）;
axis;
显示的波形如下:
改变LTI系统对输入的影响是,系统现在是一个高通滤波器。
它通过高频输入组件s2来替代低频输入组件s1.
Q2.3
clf;
n=0:
100;
s1=cos（2*pi*0.04*n）;
s2=cos（2*pi*0.6*n）;
x=s1+s2;
M=input（'滤波器所需的长度='）;
num=ones（1,M）;
y=filter（num,1,x）/M;
subplot（2,2,1）;
plot（n,s1）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'信号#1'）;
subplot（2,2,2）;
plot（n,s2）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'信号#2'）;
subplot（2,2,3）;
plot（n,x）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'输入信号'）;
subplot（2,2,4）;
plot（n,y）;
axis（[0,100,-2,2]）;
xlabel（'时间序列n'）;
ylabel（'振幅'）;
title（'输出信号'）;
axis;
结果：
Q2.19
clf;
N=40;
num=[2.24032.49082.2403];
den=[1-0.40.75];
y=impz（num,den,N）;
stem（y）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
title（'冲激响应'）;
grid;
结果：
Q2.21
clf;
N=40;
num=[0.9-0.450.350.002];
den=[1.00.71-0.46-0.62];
x=[1zeros（1,N-1）];
y=filter（num,den,x）;
stem（y）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
title（'冲激响应'）;
grid;
结果：
Q2.28
clf;
h=[321-210-403];
x=[1-23-4321];
y=conv（h,x）;
n=0:
14;
subplot（2,1,1）;
stem（n,y）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
title（'用卷积得到的输出'）;
grid;
x1=[xzeros（1,8）];
y1=filter（h,1,x1）;
subplot（2,1,2）;
stem（n,y1）;
xlabel（'时间序号n'）
ylabel（'振幅'）;
title（'由滤波生成的输出'）;
grid;
结果：
Q3.1
原始序列是：
Pause:
暂停命令
Q3.2
clf;
w=-4*pi:
8*pi/511:
4*pi;
num=[21];
den=[1-0.6];
h=freqz（num,den,w）;
subplot（2,1,1）
plot（w/pi,real（h））;
grid
title（'H（e^{j\omega}的实部'）
xlabel（'\omega/\pi'）;
ylabel（'振幅'）;
subplot（2,1,2）
plot（w/pi,imag（h））;
grid
title（'H（e^{j\omega}的虚部'）
xlabel（'\omega/\pi'）;
ylabel（'振幅'）;
pause
subplot（2,1,1）
plot（w/pi,abs（h））;
grid
title（'|H（e^{j\omega}|幅度谱'）
xlabel（'\omega/\pi'）;
ylabel（'振幅'）;
subplot（2,1,2）
plot（w/pi,angle（h））;
grid
title（'相位谱arg[H（e^{j\omega}]'）
xlabel（'\omega/\pi'）;
ylabel（'以弧度为单位的相位'）;
结果：
是w的周期。
周期2π
实部是2π为周期偶对称的
虚部是2π为周期奇对称的
幅度是2π为周期偶对称的
相位是2π为周期奇对称的
Q3.4
clf;
w=-4*pi:
8*pi/511:
4*pi;
num=[1357911131517];
den=1;
h=freqz（num,den,w）;
subplot（2,1,1）
plot（w/pi,real（h））;
grid
title（'H（e^{j\omega}的实部'）
xlabel（'\omega/\pi'）;
ylabel（'振幅'）;
subplot（2,1,2）
plot（w/pi,imag（h））;
grid
title（'H（e^{j\omega}的虚部'）
xlabel（'\omega/\pi'）;
ylabel（'振幅'）;
pause
subplot（2,1,1）
plot（w/pi,abs（h））;
grid
title（'|H（e^{j\omega}|幅度谱'）
xlabel（'\omega/\pi'）;
ylabel（'振幅'）;
subplot（2,1,2）
plot（w/pi,angle（h））;
grid
title（'相位谱arg[H（e^{j\omega}]'）
xlabel（'\omega/\pi'）;
ylabel（'以弧度为单位的相位'）;
结果：
实验二：
离散傅里叶变换和z变换
一、实验目的和任务：
●实验主要内容和要求（对应实验指导书3.5和3.6节）：
●
（1）DFT计算、常见性质的验证、用DFT实现快速卷积（应用于简单的含噪信号滤波）；
●
（2）实现z变换和逆z变换，（用于：
根据差分方程求系统的传递函数，进一步求冲激响应、零极点图）。
二、实验内容：
●具体包括——
●DFT：
Q3.23
●DFT性质：
●圆周移位子函数（Q3.26、Q3.27）
●圆周卷积子函数（Q3.28、Q3.29）
●圆周移位（Q3.30－35）
●圆周卷积（Q3.36－40）
●Z变换：
Q3.46、Q3.47、Q3.48
●逆Z变换：
Q3.50、Q3.51
三、实验过程与结果分析：
Q3.23
clf;
N=200;
L=256;
nn=[0:
N-1];
kk=[0:
L-1];
xR=[0.1*（1:
100）zeros（1,N-100）];
xl=[zeros（1,N）];
x=xR+i*xl;
XF=fft（x,L）;
subplot（3,2,1）;
grid;
plot（nn,xR）;
grid;
title（'实\{x[n]\}'）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
subplot（3,2,2）;
plot（nn,xl）;
grid;
title（'虚\{x[n]\}'）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
subplot（3,2,3）;
plot（kk,real（XF））;
grid;
title（'实\{x[k]\}'）;
xlabel（'频率指数k'）;
ylabel（'振幅'）;
subplot（3,2,4）;
plot（kk,imag（XF））;
grid;
title（'虚\{x[k]\}'）;
xlabel（'频率指数k'）;
ylabel（'振幅'）;
xx=ifft（XF,L）;
subplot（3,2,5）;
plot（kk,real（xx））;
grid;
title（'IDFT\{x[k]\}实部'）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
subplot（3,2,6）;
plot（kk,imag（xx））;
grid;
title（'IDFT\{x[k]\}虚部'）;
xlabel（'时间序号n'）;
ylabel（'振幅'）;
结果：
Q3.26在函数circshift中，命令rem的作用是什么？
R=rem（X,Y）求余数函数
Q3.27解释函数circshift怎样实现圆周移位运算。
输入序列x是循环左移M位。
如果M>0,那么circshift删除左边的元素向量x，并且附加他们到剩下的元素右边来获得循环转移序列。
如果如果M<0,然后circshift首先补充的x的长度,最右边的长度（x）-m样品从x中移走并且附加在剩下的M样本右边来得到循环转移序列。
Q3.28在函数circshift中，运算符~=的作用是什么？
如果A和B不相等返回值1
如果A和B相等返回值0
Q3.29解释函数circonv怎样实现圆周卷积运算。
函数circonv操作如下:
输入的是两个相等长度为L的两个向量x1和x2.，为了理解circonv是如何工作的,从x2的周期延拓角度来考虑很有用。
让x2p作为x2的无限长的周期延拓。
从概念上讲,常规时间反转x2p并且让x2tr通过x2p的时间反转等于元素1。
输出向量y元素1到L是通过x1和一个长度L的通过循环右移一个时间反转序列x2tr得到的序列sh之间的内积来获得的。
对于输出样例y[n],1≤n≤L、正确的循环移位是n-1点。
Q3.30
clf;
M=6;
a=[0123456789];
b=circshift（a,M）;
L=length（a）-1;
n=0:
L;
subplot（2,1,1）;
stem（n,a）;
axis（[0,L,min（a）,max（a）]）;
title（'原序列'）;
xlabel（'时间序号n'）;
ylabel（'a[n]'）;
subplot（2,1,2）;
stem（n,b）;
axis（[0,L,min（a）,max（a）]）;
title（'通过循环移位得到的序列',num2str（M）,'样本'）;
xlabel（'时间序号n'）;
ylabel（'b[n]'）;
决定时移的数量的部分是M
如果时移的数量大于序列长度，实际实现的循环时移是rem（M,length（a））点左移，相当于循环移动的M点（不止一次）,也相当于通过M点周期延拓的左移。
Q3.31
结果：
序列的长度是10，并且M=12。
这可能被解释为一个12点的循环左移（不止一次）,作为一个2点循环左移,或者作为一个2的线性左移，或者序列的12点周期延拓。
Q3.32&Q3.33.
程序：
clf;
x=[0246810121416];
N=length（x）-1;n=0:
N;
y=circshift（x,5）;
XF=fft（x）;
YF=fft（y）;
subplot（2,2,1）;
stem（n,abs（XF））;grid;
title（'原序列的DFT的幅度'）;
xlabel（'频率序号k'）;
ylabel（'|X[k]|'）;
subplot（2,2,2）;
stem（n,abs（YF））;grid;
title（'圆周位移后序列的DFT幅度'）;
xlabel（'频率序号k'）;
ylabel（'|Y[k]|'）;
subplot（2,2,3）;
stem（n,angle（XF））;grid;
title（'原序列的DFT的幅度'）;
xlabel（'频率序号k'）;
ylabel（'arg（X[k]）'）;
subplot（2,2,4）;
stem（n,angle（YF））;grid;
title（'圆周位移后序列的DFT相位'）;
xlabel（'频率序号k'）;
ylabel（'arg（Y[k]）'）;
结果：
序列的长度N=8并且时移是五个样品提前转移到左边。
相位是
。
Q3.34修改M=2
Q3.35序列为长度14，结果图如下：
序列长度为16，上述程序结果图如下：
Q3.36
g1=[123456];g2=[1-233-21];
ycir=cconv（g1,g2）;
disp（'循环卷积图像='）;disp（ycir）
G1=fft（g1）;G2=fft（g2）;
yc=real（ifft（G1.*G2））;
disp（'DFT变换乘积的IDFT变换的图像='）;disp（yc）
结果：
循环卷积图像=
Columns1through10
1.000002.00007.000010.000014.000011.000028.000012.0000-7.0000
Column11
6.0000
DFT变换乘积的IDFT变换的图像=
12281401614
Q3.37
g1=[1234];g2=[1-233];
ycir=cconv（g1,g2）;
disp（'循环卷积图像='）;disp（ycir）
G1=fft（g1）;G2=fft（g2）;
yc=real（ifft（G1.*G2））;
disp（'DFT变换乘积的IDFT变换的图像='）;disp（yc）
循环卷积图像=
1.00000.00002.00007.00007.000021.000012.0000
DFT变换乘积的IDFT变换的图像=
821147
g1=[12];g2=[1-2];
ycir=cconv（g1,g2）;
disp（'循环卷积图像='）;disp（ycir）
G1=fft（g1）;G2=fft（g2）;
yc=real（ifft（G1.*G2））;
disp（'DFT变换乘积的IDFT变换的图像='）;disp（yc）
循环卷积图像=
1.0000-0.0000-4.0000
DFT变换乘积的IDFT变换的图像=
-30
Q3.38
g1=[12345];g2=[22011];
g1e=[g1zeros（1,length（g2）-1）];
g2e=[g2zeros（1,length（g1）-1）];
ylin=cconv（g1e,g2e）;
disp（'通过圆周卷积的线性卷积='）;disp（ylin）;
y=conv（g1,g2）;
disp（'直接线性卷积='）;disp（y）
结果：
通过圆周卷积的线性卷积=
Columns1through10
2.00006.000010.000015.000021.000015.00007.00009.00005.00000.0000
Columns11through17
0.00000.00000.0000000.0000-0.0000
直接线性卷积=
2610152115795
观察可得:
零填充适当的长度确实可以实现用循环卷积实现线性卷积。
Q3.39
g1=[123];g2=[220];
g1e=[g1zeros（1,length（g2）-1）];
g2e=[g2zeros（1,length（g1）-1）];
ylin=cconv（g1e,g2e）;
disp（'通过圆周卷积的线性卷积='）;disp（ylin）;
y=conv（g1,g2）;
disp（'直接线性卷积='）;disp（y）
通过圆周卷积的线性卷积=
2.00006.000010.00006.000000.00000.00000.00000.0000
直接线性卷积=
261060
g1=[12];g2=[22];
g1e=[g1zeros（1,length（g2）-1）];
g2e=[g2zeros（1,length（g1）-1）];
ylin=cconv（g1e,g2e）;
disp（'通过圆周卷积的线性卷积='）;disp（ylin）;
y=conv（g1,g2）;
disp（'直接线性卷积='）;disp（y）
通过圆周卷积的线性卷积=
2.00006.00004.000000.0000
直接线性卷积=
264
Q3.40
g1=[12345];
g2=[22011];
g1e=[g1zeros（1,length（g2）-1）];
g2e=[g2zeros（1,length（g1）-1）];
G1EF=fft（g1e）;
G2EF=fft（g2e）;
ylin=real（ifft（G1EF.*G2EF））;
disp（'通过DFT的线性卷积='）;disp（ylin）;
结果：
通过DFT的线性卷积=
2.00006.000010.000015.000021.000015.00007.00009.00005.0000
Q3.46
clf;
w=0:
pi/51:
pi;
num=[259
升级会员 