实验一 基于DFT的信号识别系统
实验时间：5月8日 实验地点：电子楼316
实验人员：12微电子1228402027 许一飞、
1228402045 祁响、 1228402046 刘秦华
一、实验目的
（1） 通过实验巩固离散傅立叶变换DFT的认识和理解。 （2） 熟练掌握应用DFT进行频谱分析的方法。 （3） 理解DFT离散频谱分析的应用价值。
二、实验内容
在语音识别、雷达信号处理、生物医学信号检测与识别等应用领域广泛使用基于离散傅立叶变换的谱分析技术。一个典型的信号识别系统如图：
x(n) 短时窗 DFT 峰值检测 分类器 y(n)
设系统的输入信号x(n)是具有单一频谱峰值的正弦信号，短时矩形窗将信号截短为有限长，经过DFT变换得到频谱，频率检测器检测频谱最大峰值的位置，即对应的频率，然后由分类器识别信号的类别。分类器的分类判决规则为：
第一类：最大峰值频率分布范围（Hz）为0≤f≤200。 第二类：最大峰值频率分布范围（Hz）为200≤f≤500。 第三类：最大峰值频率分布范围（Hz）为500≤f≤1000。 第四类：最大峰值频率分布范围（Hz）为f≥1000。
设采样频率fs=10000Hz，短时矩形窗宽度为N=1000，短时加窗信号经过DFT可以得到连续频谱在0≤ω＜2范围内的1000个取样点。 （1）编程实现该系统 （2）输入信号x（n）=1.2sin(0.08πn)，理论计算并画出0≤f≤fs范围的幅度谱，标出峰值频率，观察系统的实际识别结果，分析其正确性。 （3）输入信号x（n）=1.5+3cos(0.5πn)，理论计算并画出0≤f≤fs范围的幅度谱，标出峰值频率，观察系统的实际识别结果，分析其正确性。 （4）输入信号x（n）=0.7sin(0.14πn)，理论计算并画出0≤f≤fs范围的幅度谱，
标出峰值频率，观察系统的实际识别结果，分析其正确性。
（5）输入信号x（n）=1.2cos(0.5πn)+ 9.5sin(0.02πn)，理论计算并画出0≤f≤fs范围的幅度谱，标出峰值频率，观察系统的实际识别结果，分析其正确性。 （6）输入信号x（n）=cos(0.102πn)，理论计算并画出0≤f≤fs范围的幅度谱，标出峰值频率，观察系统的实际识别结果，分析其正确性。
三、实验分析
设x(n)是长度为N的有限长信号（注意这个前提），即信号仅仅分布在[0,N-1]区间，其余时间均为0，那么，该信号的离散傅立叶变换定义如下：
X(k)x(n)en0N12jknN ， k=0~N-1、
又因为ω=2πf/fs，所以有2πf=ωfs。以此为依据对输入信号作如下理论分析： (1) 输入信号为x(n)=1.2sin(0.08πn)时，峰值出现在f=ωfs/2/π =400Hz
(2) 输入信号为x(n)=1.5+3cos(0.5πn)时，峰值出现在f=ωfs/2/π=0Hz 或f=2500Hz，而因为满足A2>2C，所以，应该出现在0Hz处。
(3) 输入信号为x(n)=0.7sin(0.14πn)时，峰值出现在f=ωfs/2/π=700Hz
(4) 输入信号为x(n)=1.2cos(0.5πn)+ 9.5sin(0.02πn)时，峰值出现在f=ωfs/2/π=2500Hz 或f=100Hz，因为A1>A2，所以应该出现在100Hz处。 (5) 输入信号为x(n)=cos(0.102πn)时，峰值出现在f=ωfs/2/π=510Hz
四、实验结果
（1）MATLAB 程序如下：
function result=DFTexp(A,a,B,b,C) fs=10000; N=1000; n=0:(N-1); X=zeros(1,N);
x=A*cos(a*pi*n)+B*sin(b*pi*n)+C; for k=0:(N-1)
X(k+1)=sum(x.*exp(-j*2*pi*k*n/N)); end k=0:(N-1); f=fs*k/N;
X_mag=abs(X)/(N/2); X_mag(1)=X_mag(1)/2;
[peak_mag,peak_k]=max(X_mag(1:N/2)); peak_k=peak_k-1; peak_f=fs*peak_k/N; if (peak_f<200) type=1;
else if (peak_f<500)
type=2;
else if (peak_f<1000) type=3; else type=4; end end end
result={peak_f,peak_mag,type}; subplot(2,1,1);
plot(n(1:N/2),x(1:N/2)); title('signal'); xlabel('n'); ylabel('x(n)'); subplot(2,1,2);
plot(f(1:N),X_mag(1:N)); title('spectrum'); xlabel('frequency/Hz'); ylabel('magnitude');
axis([0,fs,0,1.2*peak_mag]);
peak_str=sprintf('type %d\\n Fp=%.2fHz\\n Mp=%.2f',type,peak_f,peak_mag); text(peak_f+fs/50,peak_mag,peak_str);
(2)实验仿真结果
① X(n)=1.2sin(0.08πn)
频率：400 属于第2类
幅值1.2
② X(n)=3cos(0.5πn)+1.5
频率：2500 幅值3 属于第4类
③ X(n)=0.7sin(0.14πn)
频率：700 幅值0.7 属于第3类
④ X(n)=1.2cos(0.5πn)+9.5sin(0.02πn)
频率：100 幅值9.5 属于第1类
⑤ X(n)=cos(0.102πn)
频率：510 属于第3类
幅值1
五、思考题
1.当矩形窗长度比1000小，例如32，以上实验内容（6）可能出现什么情况？
答：N取值过小，频率分辨率过低，出现失真现象。如下图所示：
2.当输入信号x(n)=cos(0.199πn)+0.9sin(0.204πn)时，系统能够得到正确识别结果吗？为什么？ 答：由信号频率量值可以看出两个频率分量很靠近，因此在采样点不够的情况下很难分辨出两个峰值，下面给出在采样点不同的时候仿真结果：
N=1000
N=2000
N=5000
N=10000
从图中还可以看出，采样点并不是越多越好，不同的采样点个数对应的对峰值的识别也不会达到很满意的效果。
3. 如果输入信号x(n)中含有叠加性宽带噪声e(n)会影响识别结果么？为什么？
答：可能会影响识别结果。
设混叠后的信号为s(n)，则s(n)=x(n)+e(n)；所以S(K)=X(K)+E(K)，因此S(K)^2=S(K)S(K)*=X(K)^2+E(K)^2。两边求算术平方根得S(K)=X(K)^2E(K)^2
所以如果输入信号中含有叠加性宽带噪声，且噪声较大，则有可能会影响识别结果。
在噪声比较小的时候，仿真结果如下图所示：
4. 如果系统中的DFT需要更新为FFT，并且短时窗不变，则FFT计算时应该做哪儿些考虑？对识别结果会产生什么影响？
当DFT更新为FFT时，要短时窗不变，在计算的时候应该考虑以下几个方面： ① 首先应该将信号序列分为偶序列X0(K)和奇序列X1(K)两个长度相等的序列； ② 奇序列X1(K)和偶序列X0(K)的短时窗长度都将变为N/2=500点，而不再是N=1000；
③ 奇序列X1(K)和偶序列X0(K)的K值取值范围为K=0~N/2-1；，当N/2≤K＜N时，应该利用X0(K)、X1(K)的周期性特征，即X0(K)= X0(K-N/2)，X1(K)=X1(K-N/2)。
对识别结果不会产生什么影响。
※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※※ 参考资料:
◎：基于DFT的信号识别系统
◎《数字信号处理》第二版 俞一彪 孙兵/著
