DSP第一次实验报告

更新时间:2023-08-24 13:24:01 阅读量: 教育文库 文档下载

说明:文章内容仅供预览,部分内容可能不全。下载后的文档,内容与下面显示的完全一致。下载之前请确认下面内容是否您想要的,是否完整无缺。

数字信号处理课程第一次实验报告

数字信号处理实验报告

实验名称:快速傅立叶变换(FFT)及其应用

一.实验目的:

1.在理论学习的基础上,通过本实验,加深对FFT的理解,熟悉MATLAB中的有关函数。 2.应用FFT对典型信号进行频谱分析。

3.了解应用FFT进行信号频谱分析过程中可能出现的问题,以便在实际中正确应用FFT。 4.应用FFT实现序列的线性卷积和相关。

二.实验原理:

快速傅氏变换(FFT),是离散傅里叶变换的快速算法,它是根据离散傅氏变换的奇,偶, 虚,实等特性,对离散傅里叶变换的算法进行改进获得的,根据不同的情况又分为按时间抽 取的FFT 和频率抽取的FFT,同时还包括N 为任意复合数的算法以及Chirp-z 变化算法。 下面对本实验要用到的几点进行简单的说明:

(1) 混叠:采样序列的频谱是被采样信号频谱的周期延拓,当采样频率不满足奈奎 斯特采样定理的时候,就会发生混叠,使得刺痒后的序列信号的频谱不能真实 的反映原采样信号的频谱。

(2) 泄露:根据理论分析,一个时间的信号其频带宽度为无限,一个时间无限的信 号其频带宽度则为有限。因此对一个时间有限的信号,应用DFT 进行分析,频 谱混叠难以避免。对一个时间无限的信号虽然频带有限,但在实际运算中,时 间总是取有限值,在将信号截断的过程中,出现了分散的扩展谱线的现象,称 之为频谱泄露或功率泄露。

(3) 栅栏效应:DFT 是对单位圆上Z 变换的均匀采样,所以它不可能将频谱视为一 个连续函数,就在一定意义上看,用DFT 来观察频谱就好象通过一个栅栏来观 看一个景象一样,只能在离散点上看到真实的频谱,这样就有可能发生一些频 谱的峰点和谷点被“尖桩的栅栏”所挡住,不能被我们观察到。

(4) 圆周卷积:把序列X(N)分布在N 等份的圆周上,而序列Y(N)经反摺后也 分布在另一个具有N 等份的同心圆的圆周上。两圆上对应的数两量两相乘求和, 就得到全部卷积序列。这个卷积过程称做圆周卷积。

(5) 互相关函数反映了两个序列X(N)和Y(N) 的相似程度,用FFT 可以很快 的计算互相关函数。

三、 实验内容:

实验中用到的信号序列: a) 高斯序列

e

(n p)2

q

0≤n≤15

xa(n)

0 b) 衰减正弦序列

其它

e ansin(2 fn)

0≤n≤15

数字信号处理课程第一次实验报告

xb(n)

0 c) 三角波序列

n

其它 0≤n≤3 4≤n≤7 其它 0≤n≤3 4≤n≤7

xc(n)

8-n

0 d) 反三角波序列

4-n

xd(n)

n-4

0 其它

上机实验内容:

1.观察高斯序列的时域和幅频特性,固定信号xa(n)中参数p=8,改变q的值,使q分别等于2、4、8,观察它们的时域和幅频特性,了解当q取不同值时,对信号序列的时域和幅频特性的影响;固定q=8,改变p,使p分别等于8、13、14,观察参数p变化对信号序列的时域及幅频特性影响,注意p等于多少时,会发生明显的泄漏现象,混叠是否也随之出现?记录实验中观察到的现象,绘出相应的时域序列和幅频特性曲线。 程序:(p=8,q=2时的程序) clc; clear; close all; p=8; q=2; n=[0:15];

xa=exp(-((n-p).^2)/q); subplot(2,1,1); stem(n,xa);

xlabel('n');ylabel('x(n)');

title('p=8,q=2 时域特性曲线'); xk=abs(fft(xa)); subplot(2,1,2); stem(n,xk);

xlabel('k');ylabel('X(k)');

title('p=8,q=2 幅频特性曲线'); 实验结果:

固定信号xa(n)中参数p=8: q=2时的时域及幅频特性:

数字信号处理课程第一次实验报告

q=4时的时域及幅频特性:

数字信号处理课程第一次实验报告

q=8时的时域及幅频特性:

固定信号xa(n)中参数q=8: p=8时的时域及幅频特性:

数字信号处理课程第一次实验报告

p=13时的时域及幅频特性:

p=14时的时域及幅频特性:

数字信号处理课程第一次实验报告

分析:

(1)分析当q取不同值时,对信号序列的时域和幅频特性的影响

答:p固定,q增大的过程中,从时域波形来看,图像接近峰值处越平缓。从频域波形来看,高频分量减少,凹陷部分变宽,接近零的点个数增多,频谱越窄,越不容易产生混叠。 (2)察参数p变化对信号序列的时域及幅频特性影响,注意p等于多少时,会发生明显的泄漏现象,混叠是否也随之出现?

答:q固定,p增大的过程中,时域波形右移,在规定的窗口内有效值被截断的越多。因为窗口截断会造成窗口泄露,所以我们可以在幅频特性图中看到,随着p值的变大,高频分量会增加。易出现泄露,当p=13时,特别是p=14时,产生了明显的泄露与混叠。

通过这个实验,我们知道了固定信号的参数改变会对信号的时域和频域波形造成影响,还有可能出现泄漏与混叠。

2.观察衰减正弦序列xb(n)的时域和幅频特性,a=0.1,f=0.0625,检查谱峰出现位置是否正确,注意频谱的形状,绘出幅频特性曲线,改变f,使f分别等于0.4375和0.5625,观察这两种情况下频谱的形状和谱峰出现位置,有无混叠和泄漏现象?说明产生现象的原因。 程序:(a=0.1,f=0.0625时的程序) clc; clear; close all; a=0.1; f=0.0625; n=[0:15];

xb=exp(-a.*n).*sin(2*pi.*f.*n);

数字信号处理课程第一次实验报告

subplot(2,1,1); stem(n,xb);

xlabel('n');ylabel('x(n)');

title('a=0.1,f=0.0625 时域特性曲线'); xk=abs fft(xb)); subplot(2,1,2); stem(n,xk);

xlabel('k');ylabel('X(k)');

title('a=0.1,f=0.0625 幅频特性曲线');

实验结果:

a=0.1,f=0.0625的时域和幅频特性:

a=0.1,f=0.4375的时域和幅频特性:

数字信号处理课程第一次实验报告

a=0.1,f=0.5625的时域和幅频特性:

分析:

数字信号处理课程第一次实验报告

观察衰减正弦序列xb(n)的时域和幅频特性,a=0.1,f=0.0625,检查谱峰出现位置是否正确,注意频谱的形状,绘出幅频特性曲线,改变f,使f分别等于0.4375和0.5625,观察这两种情况下频谱的形状和谱峰出现位置,有无混叠和泄漏现象?说明产生现象的原因。 答:该实验中f=F/fs,F为固有频率,而fs为采样频率,统一做归一化处理使fs=1。而本题的频率分辨率为fs/N,当N=16时也为0.0625.

a=0.1,f=0.0625,sin(2πfn)为一个周期的波形,频域为单一谱线。而fs/2=0.5,故频谱关于 fs/2成镜像对称。

而a=0.1,f=0.4375时,sin(2πfn)为6.5个周期波形。a=0.1,f=0.5625时,sin(2πfn)为8.5个周期波形,频域中出现了泄漏。因为时域波形中出现了半个周期,频域中得不到单一谱线。这时由于加窗截断造成的。

而当a=0.1,f=0.5625时,由于此时不满足奈奎斯特采样定理,导致有混叠现象的产生, 并且可以看出,图中后两个序列的时域图:因为0.4375+0.5625=1,满足如下等式(此情况只适用于正弦序列),Xb(n)|f=0.4375=-Xb(n)|0.5625,即sin(2pi*fn)=-sin[2pi(1-f)n],其幅频特性是完全相同的。

通过这个实验,我们知道了频率分辨率的实际意义,以及信号频率不满足奎斯特采样定理,会有混叠现象的产生,而加窗截断会导致泄漏现象的出现。

3.观察三角波和反三角波序列的时域和幅频特性,用N=8点FFT分析信号序列xc(n) 和xd(n)的幅频特性,观察两者的序列形状和频谱曲线有什么异同?绘出两序列及其幅频特性曲线。在xc(n) 和xd(n)的末尾补零,用N=32点FFT分析这两个信号的幅频特性,观察幅频特性发生了什么变化?两种情况下的FFT频谱还有相同之处吗?这些变化说明了什么?

程序1: clc; clear; close all; n=0:7; n1=0:3; n2=4:7; x1=n1; x2=8-n2; xc=[x1,x2];

yc=abs(fft(xc,8)); subplot(2,2,1); stem(n,xc);

xlabel('t/T');ylabel('xc(n)'); title('三角波序列 时域特性'); subplot(2,2,2); stem(n,yc);

xlabel('k');ylabel('Xc(k)'); title('三角波序列 幅频特性'); x3=4-n1; x4=n2-4;

数字信号处理课程第一次实验报告

xd=[x3,x4];

yd=abs(fft(xd,8)); subplot(2,2,3); stem(n,xd);

xlabel('t/T');ylabel('xd(n)');

title('反三角波序列 时域特性'); subplot(2,2,4); stem(n,yd);

xlabel('k');ylabel('Xd(k)');

title('反三角波序列 幅频特性');

实验结果:

分析:

(1)观察两者的序列形状和频谱曲线有什么异同?

答:三角波序列和反三角波序列在时域上波形是不相同的,两者相差一个相位。事实上此两种波形为圆周移位关系。但是频域上,当N=8时,因为两者的时域只差一个相位,所以正反三角波的幅频特性相同。

程序2:

clc; clear; close all; n=0:31;

数字信号处理课程第一次实验报告

n1=0:3; n2=4:7; x1=n1; x2=8-n2;

x3=[0,0,0,0, 0,0,0,0, 0,0,0,0, 0,0,0,0, 0,0,0,0, 0,0,0,0] xc2=[x1,x2,x3];

yc2= abs(fft(xc2,32)); subplot(2,2,1); stem(n,xc2);

xlabel('t/T');ylabel('xd(n)'); title('三角波序列 时域特性'); subplot(2,2,2); stem(n,yc2);

xlabel('k');ylabel('Xd(k)'); title('三角波序列 幅频特性'); x4=4-n1; x5=n2-4;

xd2=[x4,x5,x3];

yd2=abs(fft(xd2,32)); subplot(2,2,3); stem(n,xd2);

xlabel('t/T');ylabel('xd(n)');

title('反三角波序列 时域特性'); subplot(2,2,4); stem(n,yd2);

xlabel('k');ylabel('Xd(k)');

title('反三角波序列 幅频特性');

实验结果:

数字信号处理课程第一次实验报告

分析:

(2)观察幅频特性发生了什么变化?两种情况下的FFT频谱还有相同之处吗?这些变化说明了什么?

答:改用32点FFT分析后,三角波序列和反三角波序列的幅频特性不再相同。由于栅栏效应,即由于DFT是对单位圆上Z变换的均匀采样,当N=8时,采样点较少,一些谱线被“栅栏”挡住,而通过在原序列的末端补零,N=32,即增加采样的点数和改变采样的位置,使这些被挡住的谱线显露出来,弱化了栅栏效应。于是我们知道,减小栅栏效应的一个方法是在原序列的末端填补一些零值,从而变动DFT的点数。

通过这个实验,我们充分理解了栅栏效应。同时也说明了,时域上的补零等价于频域上的线性内插。

4.一个连续信号含两个频率分量,经采样得

x(n) sin 2 0.125n cos 2 (0.125 f)n n=0,1, ,N-1

已知N=16,△f分别为1/16和1/64,观察其频谱;当N=128时,△f不变,其结果有何不同,为什么?

程序:(以'N=16△f=1/16 为例) clc; clear; close all; N=15; n=0:N;

x1=sin(2*pi*0.125.*n)+cos(2*pi*(0.125+1/16).*n); y1=abs(fft(x1));

数字信号处理课程第一次实验报告

subplot(2,1,1); stem(n,x1);

xlabel('n');ylabel('x1(n)');

title('N=16△f=1/16 时域特性'); subplot(2,1,2); stem(n,y1);

xlabel('k');ylabel('X1(k)');

title(' N=16△f为1/16 频域特性');

实验结果:

数字信号处理课程第一次实验报告

数字信号处理课程第一次实验报告

分析:

已知N=16,△f分别为1/16和1/64,观察其频谱;当N=128时,△f不变,其结果有何不同,为什么?

答:该实验中f=F/fs,F为固有频率,而fs为采样频率,统一做归一化处理使fs=1。而本题的频率分辨率为fs/N,当N=16时为1/16.而当N=64时为1/128。

当N=16,f=1/16,N=128,f=1/16以及N=128,f=1/64时,均反应了真实的频谱;因为此时的f都大于等于频率分辨率,可以分辨清楚。只有当N=16,f=1/64时,频谱发生了泄漏。这是由于窗口的截取长度选取的不合理,因为在截断过程中X(k)的长度N应为信号周期的整数倍。而且由于栅栏效应故频谱不可以分辨。当N增加至128时,提高了DFT的频率分辨率,所以△f为1/16和1/64时频谱都没有发生泄漏。

通过这个实验,我们知道增加信号的截取长度N可以提高DFT的分辨率,并且对于连续周期信号,要选择截取长度N为信号周期的整数倍,否则可能会发生泄漏。而且当信号频率差下降时,要增加N,以提高频率分辨率。

5.用FFT分别计算xa(n) (p=8,q=2) 和xb(n) (a=0.1,f=0.0625)的16点循环卷积和线性卷积。

程序: clc; clear; close all; n=0:15;

xa=exp(-(n-8).^2/2);

xb=exp(-0.1*n).*sin(2*pi*0.0625*n);

数字信号处理课程第一次实验报告

f1=conv(xa,xb); subplot(2,1,1); stem(0:30,f1);

title('线性卷积结果');

xa2=fft(xa,16);xb2=fft(xb,16); f2=real(ifft(xa2.*xb2)); subplot(2,1,2); stem(0:15,f2);

title('循环卷积结果');

实验结果:

分析:

计算xa(n)和xb(n)的16点循环卷积可以直接用MATLAB中的conv函数求得xa(n)和xb(n)的卷积和,但是计算xa(n) 和xb(n)的16点线性卷积时必须先分别对xa(n) 和xb(n)做FFT变换,得到xa(k) 和xb(k),然后把xa(k) 和xb(k)相乘,最后进行FFT的反变换。因为两个周期序列的DFS的乘积等于这两个序列的周期卷积的DFS,而循环卷积相当于周期卷积后再取主值。Xa(n)(序列长度为N1)与Xb(n)(序列长度为N2)的N点圆周卷积序列(当N<N1+N2-1),即为将Xa(n)与Xb(n)线性卷积序列中序号从N到N1+N2-1的序列叠加到原序列序号从0到N-1的地方。即周期卷积是线性卷积的周期延拓。

通过这个实验,我们知道周期卷积是线性卷积的周期延拓,两者并不相等,只有在N>N1+N2-1时,才不产生混叠现象,此时循环卷积等于线性卷积。

数字信号处理课程第一次实验报告

6.产生一512点的随机序列xe(n),并用xc(n) 和xe(n)做线性卷积,观察卷积前后xe(n)频谱的变化。要求将xe(n)分成8段,分别采用重叠相加法和重叠保留法。

程序: clear all; n=0:7; m=0:518;

xc(1:4)=n(1:4); xc(5:8)=8-n(5:8); xe=rand(1,512); hk=fft(xc,128); Xe=fft(xe,519); Xc=fft(xc,519); subplot(3,1,1);

plot(m,abs(Xc.*Xe));axis([0,518,0,250]); xlabel('k');ylabel('频域'); title('直接卷积幅频特性') for j=1:8 %重叠相加法 x(j,:)=xe(64*(j-1)+1:64*j); xk(j,:)=fft(x(j,:),128);

f(j,(64*(j-1)+1):(64*(j-1)+128))=ifft(xk(j,:).*hk); end;

y=zeros(1,576); for i=1:8 y=y+f(i,:); end

subplot(3,1,2)

plot(m,abs(fft(y(1:519))));axis([0,518,0,250]);

xlabel('k');ylabel('频域');title('重叠相加法幅频特性') x2(1,1:71)=[zeros(1,7) xe(1:64)]; %重叠保留法 for j=2:8

x2(j,:)=xe((64*(j-1)-6):(64*j)); end

x2(9,:)=[xe(506:512),zeros(1,64)]; for j=1:9

xk2(j,:)=fft(x2(j,:),128); f2(j,:)=ifft(xk2(j,:).*hk); end;

y2=zeros(1,576); for i=1:9

y2(1,64*(i-1)+1:64*i)=f2(i,8:71); end

subplot(3,1,3)

数字信号处理课程第一次实验报告

plot(m,abs(fft(y2(1:519))));axis([0,518,0,250]);

xlabel('k');ylabel('频域');title('重叠保留法幅频特性')

实验结果:

数字信号处理课程第一次实验报告

分析:

比较图中序列的线形卷积频谱可得,卷积后的结果保留了原短序列的大致趋势,但原序列的频谱曲线较线性卷积序列的频谱曲线平滑,即一个长序列与一个短序列作线性卷积,短序列就相当于一个低通滤波器,滤除长序列的一部分高频分量。

通过这个实验,我们了解了可以用FFT算法计算线性卷积,并且了解了重叠相加法与重叠保留法的区别和联系。

7. 用FFT分别计算xa(n)(p=8,q=2)和xb(n)(a=0.1,f=0.0625)的16点循环相关和线形相关,问一共有多少种结果,他们之间有何异同点。 程序: clear all; N=16; n=0:N-1; p=8;q=2;

x=exp(-(n-p).^2/q); a=0.1;f=0.0625;

y=exp(-a*n).*sin(2*pi*f*n); k=length(x);

本文来源:https://www.bwwdw.com/article/u5ni.html

Top