数字信号处理实验第一次报告实验三快速傅立叶变换及其应用Word文档下载推荐.docx
《数字信号处理实验第一次报告实验三快速傅立叶变换及其应用Word文档下载推荐.docx》由会员分享,可在线阅读,更多相关《数字信号处理实验第一次报告实验三快速傅立叶变换及其应用Word文档下载推荐.docx(19页珍藏版)》请在冰豆网上搜索。
subplot(2,1,1);
plot(n,xa);
ya=fft(xa);
ya=abs(ya);
subplot(2,1,2);
stem(n,ya);
p=8,q=2(注:
上面是时域,下面是频域)
P=8,q=4
P=8,q=8
P=13,q=8
P=14,q=8
结论:
X(n)中的参数p为高斯序列的峰值位置,q则表示高斯序列峰的尖锐度,(即峰值边沿的陡峭度)。
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;
xb=exp(-a*n).*sin(2*pi*f*n);
plot(n,xb);
yb=fft(xb);
yb=abs(yb);
stem(n,yb);
f=0.0625;
f=0.4375;
f=0.5625
结论:
该实验中f=F/fs(F—固有频率fs—采样频率,统一做归一化处理fs=1)
图中的幅频特性图:
当f=0.0625时,没有产生明显的混叠和泄露;
当f=0.4375和f=0.5625时,产生了混叠,是因为不满足奈奎斯特采样定理的缘故;
图中后两个序列的时域图:
因为0.4375+0.5625=1,满足如下等式(此情况只适用于正弦序列),Xb(n)|f=0.4375=-Xb(n)|f=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频谱还有相同之处吗?
这些变化说明了什么?
n1=0:
3;
xc1=n1;
xd1=4-n1;
n2=4:
7;
xc2=8-n2;
xd2=n2-4;
xc=[xc1,xc2];
xd=[xd1,xd2];
subplot(2,2,1);
31;
n2=0:
plot(n2,xc);
yc=fft(xc,n);
yc=abs(yc);
subplot(2,2,2);
stem(n1,yc);
subplot(2,2,3);
plot(n2,xd);
subplot(2,2,4);
yd=fft(xd,n);
yd=abs(yd);
stem(n1,yd);
n=8;
(左边是时域,右边是频域,下同)
n=32;
反三角波的边沿比较陡峭,因此它的幅频特性曲线中高频分量比较多。
由图知:
当N=8时,正反三角波的幅频特性相同,因为两者的时域只差一个相位;
当N=16时,正,反三角波的幅频特性不同。
这是因为栅栏效应,当N=8时,一些谱线被挡住。
通过在原序列的末端补零,N=16,即增加采样的点数和改变采样的位置,使这些被挡住的谱线显露出来,弱化了栅栏效应。
3.一个连续信号含两个频率分量,经采样得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-1;
x=sin(2*pi*0.125*n)+cos(2*pi*(0.125+f)*n);
plot(n,x);
y=fft(x);
y=abs(y);
stem(n,y);
N=16,f=1/16
N=16,f=1/64
N=128,f=1/16
N=128,f=1/64
当N=16,f=1/16,N=128,f=1/16以及N=128,f=1/64时,均反应了真实的频谱;
只有当N=16,f=1/64时,频谱发生了严重的栅栏效应。
这是由于分辨率等于1/N,当f>
=1/N时,能分辨,不会发生栅栏效应;
当f=<
1/N时,不能分辨,会发生栅栏效应。
4.用FFT分别计算xa(n)(p=8,q=2)和xb(n)(a=0.1,f=0.0625)的16点循环卷积和线形卷积。
程序:
p=8;
q=2;
a=0.1;
yb=abs(yb);
y1=ya.*yb;
stem(n,y1);
yaa=fft(xa,32);
yaa=abs(yaa);
ybb=fft(xb,32);
ybb=abs(ybb);
y2=yaa.*ybb;
stem(n,y2);
(上图是循环卷积,下图是线性卷积)
比较图中线性卷积与圆周卷积序列:
Xa(n)(序列长度为N1)与Xb(n)(序列长度为N2)的N点圆周卷积序列(当N<
N1+N2-1),即为将Xa(n)与Xb(n)线性卷积序列中序号从N到N1+N2-1的序列叠加到原序列序号从0到N-1的地方。
5.产生一512点的随即序列xe(n)并用xc(n)和xe(n)做线形卷积,观察卷积前后xe(n)频谱的变化。
要求将xe(n)分成8段,分别采用重叠相加法和重叠保留法。
用重叠保留法和重叠相加法实现线形卷积的过程为:
Xc(n)序列长度为8,Xe(n)序列长度为512,分Xe(n)序列为8段,每段长度为64,则每段序列与Xc(n)序列卷积后的长度为72,总长度为520。
(凑成2的整数倍)
(重叠相加法)
e=rand(1,512);
yc=fft(xc,72);
//将短序列补零后做72点的FFT
xe1=xe(1:
64);
ye1=fft(xe1,72);
//对长序列第一段做72点的FFT
y1=ye1.*yc;
//将上述两个FFT相乘
y1=[y1,zeros(1,448)];
//补上448个零,以便相加,以下7段重复上述过程
xe2=xe(65:
128);
ye2=fft(xe2,72);
y2=ye2.*yc;
y2=[zeros(1,64),y2,zeros(1,384)];
xe3=xe(129:
192);
ye3=fft(xe3,72);
y3=ye3.*yc;
y3=[zeros(1,128),y3,zeros(1,320)];
xe4=xe(193:
256);
ye4=fft(xe4,72);
y4=ye4.*yc;
y4=[zeros(1,192),y4,zeros(1,256)];
xe5=xe(257:
320);
ye5=fft(xe5,72);
y5=ye5.*yc;
y5=[zeros(1,256),y5,zeros(1,192)];
xe6=xe(321:
384);
ye6=fft(xe6,72);
y6=ye6.*yc;
y6=[zeros(1,320),y6,zeros(1,128)];
xe7=xe(385:
448);
ye7=fft(xe7,72);
y7=ye7.*yc;
y7=[zeros(1,384),y7,zeros(1,64)];
xe8=xe(449:
512);
ye8=fft(xe8,72);
y8=ye8.*yc;
y8=[zeros(1,448),y8];
y=y1+y2+y3+y4+y5+y6+y7+y8;
//将这8个序列相加,便可得到最终的结果。
n=1:
520;
plot(n,y)
(重叠保留法)
xe=rand(1,512);
xe=[zeros(1,8),xe,zeros(1,56)]//长序列前添8个零,后添56个零,构成576点的序列
xc1=n1;
对短序列做72点的FFT
72);
//将所得序列分成8段,每段序列长度为72
对长序列的第一段做72点的FFT
将上述两段序列相乘
y1=y1(9:
取第一段所得结果的后64点,以下七段同上述布骤。
xe2=xe(73:
144);
y2=y2(9:
xe3=xe(145:
216);
y3=y3(9:
xe4=xe(216:
288);
y4=y4(9:
xe5=xe(289:
360);
y5=y5(9:
xe6=xe(361:
432);
y6=y6(9:
xe7=xe(433:
504);
y7=y7(9:
xe8=xe(505:
576);
y=[y1,y2,y3,y4,y5,y6,y7,y8]//将上述8段合并,变可得到最终结果
比较图中序列的线形卷积频谱:
原序列的频谱曲线较线性卷积序列的频谱曲线陡峭,即一个长序列与一个短序列作线性卷积,短序列就相当于一个低通滤波器,滤除长序列的一部分高频分量;
6.用FFT分别计算xa(n)(p=8,q=2)和xb(n)(a=0.1,f=0.0625)的16点循环相关和线形相关,问一共有多少种结果,他们之间有何异同点。
functiony=t27
N=16;
m=(-N+1):
(N-1);
xa=exp(-(n-8).^2/2);
xb=exp(-0.1*n).*sin(2*pi*0.0625*n);
Xa1=abs(fft(xa,N));
Xb1=abs(fft(xb,N));
Xa2=fft(xa,2*N);
Xb2=fft(xb,2*N);
rm1=real(ifft(conj(Xa1).*Xb1));
rm2=real(ifft(conj(Xa2).*Xb2));
rm2=[rm2(N+2:
2*N)rm2(1:
N)];
subplot(2,1,1)
stem(n,rm1)
subplot(2,1,2)
stem(m,rm2)
(上面是循环相关,下面是线性相关)
由上图可以看到,16点的循环相关由于高斯噪声的干扰,衰减发生了微小的变化,时间位置不对了。
并且线形相关32点,循环相关只有16点。
7.用FFT分别计算xa(n)(p=8,q=2)和xb(n)(a=0.1,f=0.0625)的自相关函数。
n=0:
p=8;
q=2;
k1=length(xa);
k2=length(xb);
xak=fft(xa,2*k1);
xbk=fft(xb,2*k2);
rma=real(ifft(conj(xak).*xak));
rma=[rma(k1+2:
2*k1)rma(1:
k1)];
rmb=real(ifft(conj(xbk).*xbk));
rmb=[rmb(k2+2:
2*k2)rmb(1:
k2)];
m1=(-k1+1):
(k1-1);
m2=(-k2+1):
(k2-1);
stem(m1,rma);
stem(m2,rmb);
(上面是Xa下面是Xb)
由图可以看出最大值出现在0点,这是因为自相关的两个序列是完全一样的,之间不存在延迟。
而且,两个序列的互相关与自相关不相同。
五.思考题:
(1)实验中的信号序列Xc(n)和Xd(n),在单位圆上的变化频谱|Xc(exp(jw)|和|Xd(exp(jw))|会相同吗?
如果不同,说出哪一个低频分量更多一些,为什么?
答:
他们的单位圆上的Z变换频谱不同,Xc(n)的时域波形比较平缓,顾其低频分量会多一些。
(2)对于一个有限长度序列进行DFT等价于将该序列周期延拓后进行DFT展开,因为DFS也知识取其中一个周期来运算,所以FFT在一定条件下也可以用以分析周期信号序列。
如果实正弦信号sin(2pi*fn)f=0.1用16点FFT来做DFS运算,得到的频谱是信号本身的真实谱吗?
为什么?
不是信号本身的真实谱。
因为原信号的周期是10,因而进行16点FFT时,其16点的周期延拓后的波形已经不再是原来的波形了。