时域采样与频域采样定理的验证实验Word下载.docx
《时域采样与频域采样定理的验证实验Word下载.docx》由会员分享，可在线阅读，更多相关《时域采样与频域采样定理的验证实验Word下载.docx（33页珍藏版）》请在冰豆网上搜索。
为使用DFT，首先用下面的公式产生时域离散信号，对三种采样频率，采样序列按顺序用x1（n）、x2（n）、x3（n）表示。
因为采样频率不同，得到的x1（n）、x2（n）、x3（n）的长度不同，长度（点数）用公式N=Tp×
Fs计算。
选FFT的变换点数为M=64，序列长度不够64的尾部加零。
X（k）=FFT[x（n）]，k=0,1,2,3,…,M－1
式中,k代表的频率为
要求：
编写实验程序，计算x1（n）、x2（n）和x3（n）的幅度特性，并绘图显示。
观察分析频谱混叠失真。
（2）频域采样理论的验证。
给定信号如下：
编写程序分别对频谱函数X（ejω）=FT［x（n）］在区间［0,2π］上等间隔采样32点和16点，得到X32（k）和X16（k）：
再分别对X32（k）和X16（k）进行32点和16点IFFT，得到x32（n）和x16（n）：
分别画出、X32（k）和Ｘ16（k）的幅度谱，并绘图显示x（n）、x32（n）和x16（n）的波形，进行对比和分析，验证总结频域采样理论。
时域：
A=444.128;
a=64*sqrt
（2）*pi;
w0=64*sqrt
（2）*pi;
Tp=64/1000;
F1=1000;
F2=300;
F3=200;
T1=1/F1;
T2=1/F2;
T3=1/F3;
n1=0:
Tp*F1-1;
n2=0:
Tp*F2-1;
n3=0:
Tp*F3-1;
x1=A*exp（-a*n1*T1）.*sin（w0*n1*T1）;
x2=A*exp（-a*n2*T2）.*sin（w0*n2*T2）;
x3=A*exp（-a*n3*T3）.*sin（w0*n3*T3）;
f1=fft（x1,length（n1））;
f2=fft（x2,length（n2））;
f3=fft（x3,length（n3））;
k1=0:
length（f1）-1;
fk1=k1/Tp;
k2=0:
length（f2）-1;
fk2=k2/Tp;
k3=0:
length（f3）-1;
fk3=k3/Tp;
subplot（3,2,1）
stem（n1,x1,'
.'
）
title（'
（a）Fs=1000HZ'
）;
xlabel（'
n'
ylabel（'
x1（n）'
subplot（3,2,3）
stem（n2,x2,'
（b）Fs=300HZ'
x2（n）'
subplot（3,2,5）
stem（n3,x3,'
（c）Fs=200HZ'
x3（n）'
subplot（3,2,2）
>
plot（fk1,abs（f1））;
title（'
（a）FT[xa（nT）],Fs=1000Hz'
xlabel（'
f（Hz）'
ylabel（'
幅度'
subplot（3,2,4）;
plot（fk2,abs（f2））;
（a）FT[xa（nT）],Fs=300Hz'
plot（fk3,abs（f3））;
（a）FT[xa（nT）],Fs=200Hz'
频域：
M=27;
N=32;
n=0:
M;
xa=0:
floor（M/2）;
xb=ceil（M/2）-1:
-1:
0;
xn=[xa,xb];
Xk=fft（xn,1024）;
X32k=fft（xn,32）;
X32n=ifft（X32k）;
X16k=X32k（1:
2:
N）;
X16n=ifft（X16k,N/2）;
subplot（3,2,2）;
stem（n,xn,'
boxon
（b）三角形序列x（n）'
x（n）'
axis（[0,32,0,20]）
k=0:
1023;
wk=2*k/1024;
subplot（3,2,1）;
plot（wk,abs（Xk））;
（a）FT[x（n）]'
\omega^pi'
|X（e^j^\omega）|'
axis（[0,1,0,200]）
N/2-1;
subplot（3,2,3）;
stem（k,abs（X16k）,'
boxon;
（c）16点频率采样'
k'
|X_1_6（k）|'
axis（[0,8,0,200]）;
subplot（3,2,4）;
stem（n1,abs（X16n）,'
（d）16IDFT|X_1_6（k）|'
|X_1_6（n）|'
axis（[0,32,0,20]）;
N-1;
subplot（3,2,5）;
stem（k,abs（X32k）,'
（e）32点频率采样'
|X_3_2（k）|'
axis（[0,16,0,200]）;
subplot（3,2,6）;
stem（n1,abs（X32n）,'
（f）32点IDFT|X_3_2（k）|'
|X_3_2（n）|'
axis（[0,16,0,20]）
4．思考题
如果序列x（n）的长度为M，希望得到其频谱在［0,2π］上的N点等间隔采样，当N<
M时，如何用一次最少点数的DFT得到该频谱采样？
答：
从实验的结果中，可以看出对于求频域采样点数N小于原时域序列长度M的N点离散频谱是，可以先对元序列x（n）以N为周期进行周期延拓后去其主值区序列
再计算N点DFT则得到N点频域采样：
其所求的N点离散频谱对应的时域离散序列是元序列x（n）以N为周期进行周期延拓后取主值区序列，而不是原序列x（n）.
5、实验小结
通过这次实验，我对时域采样和频域采样的理论、定理的理解更加透彻，以前只是课堂上老师说的，现在通过自己亲手做了，也知道了其内在的原理是怎么样的了。
在这次实验中，无论是在时域还是频域，对信号采样必须仔细考虑采样的参数：
采样频谱、采样周期、采样点数，因为我在实验过程中，把stem（）函数的参数设置错了，其直接导致结果出错，我对代码检查了好几遍才改过来的，所以，我们要做，就是选取好采样的参数，避免另一个域周期延拓时发生混叠，否则，我们采样所得的数据肯定丢失一部分原信号的信息，我们便无法对原信号对原信号进行恢复和正确分析。
此次实验所遇到的问题：
主要是时域和频域的理论通过自己亲自设计的时候，对一些小知识点理解的不够彻底，使得在用MATLAB命令时，有点难下手。
所以，理论应用于实践，实践源于理论，理论很重要的。
设计二正余弦信号的谱分析
1．实验目的
学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便正确应用FFT。
2.实验原理
用FFT对信号作频谱分析是学习数字信号处理的重要内容。
经常需要进行谱分析的信号是模拟信号和时域离散信号。
对信号进行谱分析的重要问题是频谱分辨率F和分析误差。
频谱分辨率直接和FFT的变换区间N有关，因为FFT能够实现的频率分辨率是2π/N，因此要求2π/N≤F。
可以根据此式选择FFT的变换区间N。
误差主要来自于用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时,离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。
周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。
如果不知道信号周期，可以尽量选择信号的观察时间长一些。
对模拟信号进行谱分析时，首先要按照采样定理将其变成时域离散信号。
如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。
3．实验步骤及内容
（1）对一个频率为10Hz，采样频率为64Hz的32点余弦序列进行谱分析,画出其频谱图；
若将频率改为11Hz，其他参数不变，重新画出该序列的频谱图，观察频率泄漏现象，分析原因。
M文件：
31;
x1=cos（pi*20*n/64）;
%
x2=cos（pi*22*n/64）;
subplot（2,2,1）;
stem（n,x1）%绘制序列的幅频特性曲线
10HZ余弦序列'
subplot（2,2,2）;
stem（n,x2）
11HZ余弦序列'
X1=abs（fft（x1,32））;
%求x1余弦序列
subplot（2,2,3）
stem（k,X1）%绘制序列的幅频特性
X1（k）'
10HZ时32点FFT幅频曲线'
X2=abs（fft（x2,32））;
%求x2余弦序列
subplot（2,2,4）
stem（k,X2）
X2（k）'
11HZ时32点FFT幅频曲线'
泄漏的原因:
通过图可看出：
频率为10Hz的余弦曲线DFT只有两个点不等于零，位于k=5和
k=27处，k=5对应于频率10Hz，k=27对应于频率54Hz。
这样DFT确实正确的分辨了余弦信号的频率。
将频率改为11Hz，采样频率和窗长度依然为32点，计算图像可看出：
频谱图上k=5和k=27处都有较大的峰值，而其它的点上幅度不再为零。
这两个峰值对应的频率为10Hz和12Hz所以信号的峰值位于两者之间本来是单一的11Hz频率的能量会分不到许多DFT频率上，这种现象叫频率泄露，来源于截断效应。
（2）考察DFT的长度对双频率信号频谱分析的影响。
设待分析的信号为
令两个长度为16的正余弦序列的数字频率为
及
。
取N为四个不同值16，32，64，128。
画出四个DFT幅频图。
分析DFT长度对频谱分辨率的影响。
程序：
16;
x=0.5*sin（pi*2*0.22*n）+sin（pi*2*0.34*n）;
X1=abs（fft（x,16））;
%求余弦序列的16点FFT
subplot（2,2,1）
15;
stem（k,X1）%绘制序列的幅频特性曲线
16点FFT幅频曲线'
）
X2=abs（fft（x,32））;
%求余弦序列的32点FFT
subplot（2,2,2）
stem（k,X2）%绘制序列的幅频特性曲线
32点FFT幅频曲线'
X3=abs（fft（x,64））;
%求余弦序列的64点FFT
63;
stem（k,X3）%绘制序列的幅频特性曲线
X3（k）'
64点FFT幅频曲线'
X4=abs（fft（x,128））;
%求余弦序列的128点FFT
127;
stem（k,X4）%绘制序列的幅频特性曲线
X4（k）'
128点FFT幅频曲线'
DFT长度对频谱分辨率的影响:
DFT样本值就是其DTFT在相应位置的采样。
在图中很难看
出两个峰值，因此要提高它的分辨率，故把R增大，逐渐可以看出它有两个峰值，将k换算
成数字频率f=w/2*pi=k/R.这样可确定峰值的位置大体在f=0.21和0.35之附近与信号
的给定频率有一定的误差，这也是截断和泄露带来的问题，在这图上还可以看到一些较小
的峰，这是很难判断是输入信号固有的还是由泄露引起的。
这说明了增加DFT长度R减小
了相邻样本间的频率间距，提高频谱的视在分辨率，因而可以提高样本未知的测定精度。
（3）在上题中若把两个正弦波的频率取得较近，令
，
，试问怎样选择FFT参数才能在频谱分析中分辨出这两个分量？
要能分清两个频率，分辨率至少应达到f=0.03。
因为此处的数字频率是对采样频率Fs
进行归一化，因此总的样本数至少要达到1/0.03=33。
加窗以后可以使频谱函数更加光滑便于分辨峰值位置和准确的数值，为了提高实际分辨率，应该尽量增加信号的长度N及DF长度R,当受到条件限制不能提高N，则单独提高R可以提高视在分辨率。
（1）对于周期序列，如果周期不知道，如何用FFT进行谱分析？
如果X（n）的周期预先不知道，可先截取M点进行DFT，即
0<
=K<
=M-1
再将截取长度扩大一倍，截取
=2M-1
比较
和
，如果二者的主谱差别满足分析误差要求，则以
或
近似表示
的频谱，否则继续将截取长度加倍，直至前后再次分析所得主谱频率差别满足误差要求，设最后截取长度为iM，则
表示
点的谱线强度
（2）如何选择FFT的变换区间？
（包括非周期信号和周期信号）
对于非周期信号：
有频谱分辨率F，而频谱分辨率直接和FFT的变换区间有关，因为FFT能够实现的频率分辨率是2π/N...因此有最小的N>
2π/F。
就可以根据此式选择FFT
的变换区间。
对于周期信号，周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。
设计三语音信号滤波处理
1.实验目的
（1）了解语音信号的产生、采集，能绘制语音信号的频率响应曲线及频谱图；
（2）学会用MATLAB对语音信号进行分析和处理；
（3）掌握用滤波器去除语音信号噪声的方法，观察去噪前后的语音信号。
2.实验原理
（1）语音信号的采集
在MATLAB软件平台下，利用函数wavread（）对语音信号采集，并记录采样频率和采样点数。
将语音信号转换成计算机能够运算的有限长序列。
wavread函数的调用格式如下：
y=wavread（file）
读取file所规定的wav文件，返回采样值放在向量y中。
[y,fs]=wavread（file）
采样值放在向量y中，fs表示采样频率（Hz）。
（2）用FFT作谱分析
FFT即快速傅立叶变换，它是从DFT运算中发展起来的，利用系数
的对称性和周期性减少运算量。
长度为N的序列直接计算DFT需要
次复乘和
次复加，而用FFT进行运算一般需要
次复加，从而使DFT的运算大大简化。
用FFT对连续信号进行谱分析的步骤如下：
图3.1连续信号谱分析过程
引入前置低通滤波器LPF是为了消除或减少时域连续信号转换成序列时可能出现的频谱混叠现象。
表示时域有限的窗函数。
（3）设计滤波器去除语音信号的噪声
通过wavread（）函数将语音信号读入，通过频率采样及fft（）产生信号，并对之加噪，通过窗函数法设计滤波器滤掉该语音信号的噪声，并对比滤波前后的语音波形和频谱。
3.实验内容
（1）利用Windows下的录音机录制一段自己的话音，时间在1s内。
然后在Matlab软件平台下，利用wavread函数对语音信号进行采样，记住采样频率和采样点数；
x=wavread（'
D:
\test.wav'
）;
[x,fs]=wavread（'
%把语音信号进行加载入Matlab仿真软件平台中
sound（x,fs）;
%对加载的语音信号进行回放
stem（x）;
语音信号的时域波形'
（2）画出语音信号的时域波形，对采样后的语音信号进行快速傅立叶变换，得到信号的频谱特性，画出采样信号的时域波形和频谱图；
Fs=8000;
%给出抽样频率
[x,fs,nbits]=wavread（'
sound（x,fs,nbits）;
n=length（x）;
%求出语音信号的长度
X=fft（x,n）;
%傅里叶变换
plot（x）;
原始信号时域波形'
plot（abs（X））;
原始域波形'
[x,fs,nbits]=wavread（'
noise=0.04*sin（10000*pi*x）;
%sin函数产生正弦噪声
s=x+noise;
%语音信号加入正弦噪声
sound（s）;
plot（s）;
加正弦噪语音信号的时域波形'
S=fft（s）;
plot（abs（S））;
加正弦噪语音信号的频域波形'
n=length（x）;
noise=rand（size（x/5））;
%产生白噪声
%语音信号加入白噪声
加白噪声语音信号的时域波形'
加白噪声语音信号的频域波形'
（3）根据对语音信号谱分析的结果，确定滤除噪声所需滤波器的技术指标，设计合适的数字滤波器，并画出滤波器的频域响应；
根据对加噪语音信号谱的分析得滤除噪声所需的滤波器应为低通滤波器。
技术指标如下：
Fp=1000HZ;
Fs=1200HZ;
rp=0.1;
rs=60
Fp=1000;
Fs=1200;
wp=2*pi*Fp;
ws=2*pi*Fs;
rp=0.1;
rs=60
[n,wn]=ellipord（wp,ws,rp,rs,'
s'
[b,a]=ellip（n,rp,rs,wn,'
low'
'
fk=linspace（1,3000,100）
w=2*pi*fk
h=freqs（b,a,w）
magh=abs（h）
subplot（2,1,1）;
plot（fk,magh）;
低通幅频响应'
phah=unwrap（angle（h））
subplot（2,1,2）;
plot（fk,phah）;
低通相频响应'
（4）用所设计的滤波器对采集的信号进行滤波，在同一个窗口画出滤波前后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化；
答：
程序:
Ft=8000;
wp=2*pi*Fp/Ft;
ws=2*pi*Fs/Ft;
fp=2*Ft*tan（wp/2）;
fs=2*Fs*tan（wp/2）;
[n11,wn11]=buttord（wp,ws,1,50,'
%求低通滤波器的阶数和截止频率
[b11,a11]=butter（n11,wn11,'
%求S域的频率响应的参数
[num11,den11]=bilinear（b11,a11,0.5）;
%利用双线性变换实现频率响应S域到Z域的变换
[x,fs,nbits]
升级会员 