原创 MATLAB信号处理仿真-短时傅里叶变换,时、频分辨力的跷跷板

2012-10-9 22:48 19977 19 31 分类: 处理器与DSP 文集: MATLAB信号处理仿真

本次,我们将展示另外一个信号分析的有力工具,短时傅里叶变换。

对于探讨经典信号处理的问题,我的感觉相对舒畅,因为从信号与系统到离散时间信号处理,建立在无穷时域前提下的各种结论,完美的仿佛传说中的理想国与乌托邦。尽管DFT与FFT的时域有限长采样造成了频谱泄露,但是我们似乎仍然心存希望,至少可以通过增加采样长度来降低泄露,至少还有这一途径可以让我们感到些许温暖,对于频率成分固定不变的信号,这种增加采样长度的方法可以让我们逼近理想主义的完美。但是,对于频率成分随时间发生改变的信号,以上的方法陷入一个两难境地,如果采集很长的一段数据进行DFT分析,可以很精确的得知出现了那些频率成分,但是却无法知道这个成分是在何时出现的,或者,采集很短的数据做DFT分析,可以相对准确的知道在什么时候出现了某个频率成分,然而对于这个成分的具体频率值可能就探测的不那么准确,于是,对于频谱成分随时间变化的信号,我们的傅里叶分析真的变成了一个跷跷板问题,一头好了,另一头就糟糕。不管怎样,有总比没有强,这就是本次的主题,短时傅里叶变换,在奥本海姆的书中称为“依时傅里叶变换”。

具体的数学推导此处就不再赘述了,并且没有数学排版工具也不方便展示公式,维基百科英文版的短时傅里叶变换写的也比较生动,此处我们定性的讲一讲,就是如果有一个很长的时间信号序列,我们每次从中取出一个长度为L的片段,加窗做DFT,看看频谱,这就是短时傅里叶变换,由于多引入了一个时间的维度,所以存在一个问题,就是两次做DFT的L点数据帧,它们在原始信号中的相对位置是怎样的,是否允许有些重叠呢?我们思考一下,如果两帧信号完全不重叠,那么由于加窗的原因(想想,窗函数的边界处通常都是接近0的数据)交界处的信号突变会被窗函数抹杀掉,所以通常来说两帧数据是要有重叠的,而且重叠的部分越多,对于捕获突发的信号越有帮助,这样做的代价呢,就是计算量的增加。

那么可能您会有疑问,在做短时傅里叶变换的时候,增加DFT的帧长度L会不会是个一本万利的事情呢?很遗憾,不是,因为对于L点的DFT的分析中找到的谱峰,其有效意义对应于整个的L点数据,所以无法分辨出是在L点数据中的那个位置出现的,也许你会想,可以把两帧数据的重叠点数增加啊,这样是否利用观测某个谱峰在那个时间出现呢?的确,提高重叠点数是有帮助,但是,得到的分析结果仍然是一种不清晰的结果。口说无凭,下面我们通过实验现象来观察。

首先我们观察一个DTMF信号离散双音多频信号,说白了就是电话拨号音,这个是我手机号码的电话拨号音,如果你能够译码出我的手机号,欢迎发短信告诉我,请注明是在信号处理博客上看到的,我的博客ID是“duweitao”,
请不要把我的手机号公布在网站上,免得给大家添麻烦

首先看这个拨号音的时域信号,是若干的时域脉冲的组合,脉冲之间是短暂的静默信号。下图第一行是全局波形,第二行是放大的脉冲局部,每一个拨号数字是用2个频率分量来表示的,不同数字的两个频率分量的组合不同。

 

20121008225612512.jpg


然后我们把全部的时域信号拿来做个DFT,我们看到似乎是有6根谱线,但是你能说出那根谱线在何时出现么,或者某个谱线总共出现了几次呢?

20121008230021733.jpg


接下来,采用256点的短时DFT进行分析,用光谱图spectrogram来观察结果,可以清楚的看到在时间轴上各个断续的音调,但是频率轴上音调的中心位置就不那么清晰了。

20121009210528202.jpg

然后,我们采用1024点的短时DFT进行分析,这一次可以清楚的看到频率的中心位置,但是,你发现了吗,有些音调之间本来存在的空隙,在这个光谱图上看不见了,可是,事实上它们是存在的。

20121009210608216.jpg


所以说时域和频域的解析能力是一个跷跷板的两头,你要优化这一头,那么另一头就变得糟糕。如果你对“spectrogram”函数不了解,请自行在matlab中doc 之,或者,请看下面这张图,光谱图是以下这个3D谱图的俯视图,即从上往下看,能量越强,离观察者越近的谱峰,颜色越热。

20121008230156515.jpg

把同样的方法,应用于附件的“白杨礼赞”普通话播音文件,我们可以看到如下的图片
首先是语谱图

 

20121009214432814.jpg

然后是3D视图

20121009214505462.jpg

从此中,我们可以看到,人类的语音中,充满了时变的频谱分量,人类语音信号中有一个重要参数叫做基音周期,亦即基波的频率的倒数,这个基频的频率和其各次倍频的能量关系,以及它们的时变规律,决定了一个人的声音特征。

语言,多么普通,我们天天在使用,用来传递含义或是表达感情。语音分析与合成是最早使用数字信号处理作为研究工具的领域之一,尽管人类的语音识别、合成的算法日益成熟乃至可以商用,然而,在如此平常的信号里却蕴藏了我们永远也无法完全探究的奥秘,没有两个人的语谱图特征是完全一样的,这东西就像指纹一样,在刑侦学科上,被称作“声纹”,说实话,每次看这种语谱图我的心中都充满敬畏,为什么会是这样? 我不知道,没人知道,也许只有上帝知道。

 

 

%///////////////////////////////////////////////////////////
% short time DFT analyse wav file
%///////////////////////////////////////////////////////////
close all;
clear;
clc;
% set target file at latter line
wav_fname              = 'baiyang_test.wav' ;
wav_fname              = 'dial_tone.wav'    ;

[wav_data, fs, q_bits] = wavread(wav_fname) ;
data_len               = length(wav_data)   ;
kaiser_beta            =  8                 ;

stdft_len_win          = 256               ;
stdft_len_overlap      = stdft_len_win-16  ;
stdft_len_dft          = stdft_len_win     ;

figure;
subplot(2,1,1);
plot(wav_data); title('Total signal');
subplot(2,1,2);
sect_len = 800;
plot(wav_data(1:sect_len));title('First part');

kaiser_win_spectrum_plot(fs, wav_data, kaiser_beta);
y_min = -160; y_max = 10;
ylim([y_min, y_max])    ;   % set y-axis range

title('Input Wave File, Normalized Spectrum', 'fontsize', 14);

dscp_str   = sprintf('fs %.2f KHz, L %d, q %d bits, beta %d',...
fs/1E3, data_len, q_bits, kaiser_beta);

text(0.5*(-fs/2),0.8*y_min,dscp_str, ...
     'FontSize',14, 'color', 'g','BackgroundColor','k') ;

len_w = stdft_len_win      ;   
len_o = stdft_len_overlap  ;  
len_d = stdft_len_dft      ;
figure  ;
% use tic toc get running time
tic     ;
spectrogram(wav_data, len_w, len_o, len_d, fs);
title('Spectrogram', 'fontsize', 14);
toc     ;
text_x = 0.7*(fs/2);
text_y = 0.5* (data_len/fs);
dscp_str = sprintf(' data len %d \n win len %d\n overlap %d\n dft len %d',...
           data_len, len_w, len_o, len_d);
text(text_x, text_y,dscp_str, ...
     'FontSize',14, 'color', 'g','BackgroundColor','k') ;
colormap('jet');
% plot 3d-mesh
S = spectrogram(wav_data, len_w, len_o, len_d, fs);
figure; mesh(log10(1E-4+abs(S))*10);
xlabel('Time', 'FontSize',14);
ylabel('Frequency', 'FontSize',14);

 

文章评论12条评论)

登录后参与讨论

用户1871089 2016-2-17 18:04

谢谢

用户377235 2016-1-4 15:54

太好了,通俗易懂

用户377235 2015-8-10 11:20

太感动了写得真好!醍醐灌顶!

用户377235 2015-6-1 11:18

你使用的什么窗函数来构造图谱的?

用户1836868 2015-5-4 09:22

以前没学过信号方面的,正好下载来看看,多谢了!

用户617985 2014-11-12 15:40

多谢楼主,学习了

用户377728 2013-11-8 00:41

谢谢分享

用户452017 2013-11-4 17:36

太棒了 正好需要这个呢 哈哈 谢谢分享了

用户433433 2013-3-18 19:30

楼主辛苦了,学习了

用户1469633 2012-11-12 13:33

楼主帅气不解释!
相关推荐阅读
用户424825 2013-12-04 21:13
嵌入式环境下的算法开发之学习建议
         说到在嵌入式Linux的平台上设计算法,目的无非是2个,一是Linux,这东西无孔不入,可以塞到各种板子上。二是算法代码执行的效率,嵌入式的平台不像是PC各种资源还是比较紧张...
用户424825 2013-10-27 15:56
博客测试记录
根据EDN的网页温馨提示,添加了Host DNS的信息,可以从学校访问了,哈利路亚,党的光芒暖人心 可以从传媒大学校内访问,2013年4月22日,于实验室 无法从传媒大学校内访问,2013...
用户424825 2013-10-06 18:34
博客页面测试
博客测试 插图测试     不管怎么说,EDN的博客是越做越好,而且编辑器的功能也日益强大。 edn的确是用心在做这个事情了,打从两年前网站重新改版了一次之后,折腾的有点乱,现在似乎是又变得好起...
用户424825 2013-08-20 09:56
高性能高频收发机设计:业余无线电爱好者的福音
之前在EDN上发表过的文章《21世纪的业余无线电》引起了人们诸多的兴趣,人们对这个业余爱好的各个方面也提出了许多问题。本文撰写的动机也是源于这方面的兴趣,但不只限于业余无线电。本文讨论了工程师在设计各...
用户424825 2013-08-15 09:04
评论:@用芯创造未来 博客中提到的“【博客大赛】心得分享”
此文有深意,收藏之 ---------------- 这段时间经历太多踏入职场的第一次,第一次引荐别人,也第一次经历被引荐的人离开,第一 次与上层深入探讨行业发展并交互多重信息,第一次考博...
用户424825 2013-08-14 06:55
EDN 好文记录帖
  “高性能高频收发机设计:业余无线电爱好者的福音 ”。 介绍接收机体系结构的文章,其中末尾部分对于80MHz,14比特带通直采 与 传统混频方案的噪声对比分析比较精彩。值得收藏   ...
我要评论
12
19
关闭 站长推荐上一条 /2 下一条