信息工程概论

1、信息工程概论作业窗口傅里叶变换姓名：白子轩学号： 2130602008班级：信计31一、传统的傅里叶变换我们都知道， 信号分析中最重要的两个参数是时间和频率，而我们一般所得到的信号表示形式都是 f t 的形式，而我们可以通过传统的傅里叶变换f?( )f (t )e j t dt ，可以把信号变为频域表示f（。）但是，传统的傅里叶变换只对平稳的信号有用。 对于非平稳的信号需要用时间和品率的联合函数来表示信号。因此，我们需要短时傅里叶变换，也就是窗口傅里叶变换。二、窗口傅里叶变换对于信号的频率是随时间变化的信号。为了获得它的随时间变化的频谱，最采用的处理办法是加窗技术对信号截取，然后对截取的局部信

2、号作 Fourier 变换。然后不断地移动窗口函数中心的位置，就可以得到信号的局部区域的瞬时频率，因此，对于连续的信号，它的窗口傅里叶变换为：Sf (u, )f (t )g(tu)e i t dt窗口傅里叶逆变换为：f (t )1Sf (u, )g(t u)ei t d du2而对于离散的信号，它的窗口傅里叶变换为：N1i 2 lnSf m,lf n g nm expn0N窗口傅里叶逆变换为：f n1 N 1N 1Sf m, l g nm exp i 2 lnN m 0l 0N三、窗口函数要进行窗口傅里叶变换，首先要要选择窗口函数，窗口函数有很多，例如高斯窗、hamming 窗和 Hannin

3、g 窗等等。其中高斯窗函数被设计为了分析瞬态信号，Hamming 和 Hann窗函数被设计为了分析窄带信号，Kaiser-Bessel 窗函数可用于更好地分离两个频率成分非常接近但振幅完全不同的信号。在 matlab 中我们也可以直接调用一些窗口函数，调用的方法如下：四、实验实验 1：题目：在这里我先做了一下书上的例子，对线性调频信号ftsin 22t 2 进行频谱分析。题目分析：这个例子有两种做法，第一种方法是直接调用matlab 中的 spectrogram 函数，第二种方法是按照定义选取窗口函数，然后对每一小段的做快速傅里叶变换就可以了。方法一：源程序：cleart=0:0.001:10

4、;t1=t;f1=sin(2*pi*2*power(t,2);subplot(2,2,1);plot(t,f1);subplot(2,2,2);g=1/6*exp(-0.5*power(t,2).*(t>=-3 & t<=3)+0.*(t>3 | t<-3);t=-4:0.01:4;g1=1/6*exp(-0.5*power(t,2).*(t>=-3 & t<=3)+0.*(t>3 | t<-3);plot(t,g1);subplot(2,2,3);S,F,T,P = spectrogram(f1,gausswin(600),58

5、0,600,1E3);surf(T,F,10*log10(abs(P),'edgecolor', 'none');axis(0 10 0 50)view(0,90);xlabel('Time (Seconds)');ylabel('Hz');subplot(2,2,4);surf(T,F,10*log10(abs(P)+80,'edgecolor', 'none');axis(0 10 0 50 0 200)%axis tightxlabel('Time (Seconds)');yl

6、abel('Hz');zlabel('enargy')得到结果：结果分析：我们用的是高斯窗口，得到了一个很好的结果，无论是2D 图还是 3D 图，都与书上的图十分相似。但有一个十分大的缺陷，就是无法重构原来的信号，因为我们是直接调用的 spectrogram 函数，并不太知道里面的具体程序是什么样的，所以无法还原原信号，也无法计算误差。因此我们就需要第二种方法。方法二：源程序：cleart=0:0.001:10;t1=t;f1=sin(2*pi*2*power(t,2);subplot(3,2,1);plot(t,f1);subplot(3,2,2);g=1/6

7、*exp(-0.5*power(t,2).*(t>=-3 & t<=3)+0.*(t>3 | t<-3);t=-4:0.01:4;g1=1/6*exp(-0.5*power(t,2).*(t>=-3 & t<=3)+0.*(t>3 | t<-3);plot(t,g1);N=length(f1);Nw=20;%窗函数长window lengthL=19;%窗函数每次移动的样点数, 重叠宽度Ts=round(N-Nw)/L)+1;%计算把数据x 共分成多少段nfft=2ceil(log2(Nw) ;%FFT的长度TF=zeros(Ts

8、,nfft);%将存放三维谱图，先清零%for i=1:Tsi=0;flag=0;whileflag=0;i=i+1;if(i-1)*Nw+Nw<N%y(i-1)*L+1:i*L+L)%hamming(Nw)xw=f1(i-1)*L+1:i*L+L)*hamming(2*L);%取一段数据temp=fft(xw,nfft);%FFT变换temp=fftshift(abs(temp);%频谱以 0 频为中心%TF(i,:);TF(i,:)=temp;%把谱图存放在TF 中elseflag=1;endendsubplot(3,2,3);%mesh(abs(TF);%view(0,90);%a

9、xis tight%imagesc(TF);contour(abs(TF);xlabel(' 时间 ' );ylabel(' 频率 ' )subplot(3,2,4);mesh(abs(TF);%三维绘图axistighttitle('STFT' );xlabel(' 时间 ' );ylabel(' 频率 ' );%短时傅里叶变换X=fft(f1);X=fftshift(X);%重构x1=ifftshift(X);x1=ifft(x1);subplot(3,2,5);plot(t1,x1);xlabel('

10、 时间 ' );%x轴ylabel(' 振幅 ' );%y轴%误差e=x1-f1;subplot(3,2,6);plot(t1,abs(e);xlabel(' 时间 ' );%x轴ylabel(' 振幅 ' );%y轴得到结果：结果分析：我们发现得到的结果非常不理想， 经过多次试验， 我并没有找到是哪里出的问题， 因此在接下来的实验中将要放弃这种方法。实验 2：题目：录一段自己的声音，对这段声音进行频谱分析。实验分析：这个题目和实验1 类似，只不过一个用的是连续的窗口傅里叶变换，一个用的是离散的窗口傅里叶变换，所以在这里我们仍然运用mat

11、lab 中的 spectrogram函数即可。在这里我们需要注意一下， mp3是双音轨的，所以数字化后是一个是一个二维数组，而我们只需要分析其中的一组即可。于是我就录了自己的一段笑声，并抽取其中的一组数组，进行分析。源程序：clearx,fs=audioread(小波分析笑 .mp3' );y=x(:,1)'subplot(3,2,1);plot(y);title(' 信号波形图 ' );%图名xlabel(' 时间 ' );%x轴ylabel(' 振幅 ' );%y轴holdonN=length(y);Nw=600;%窗函数长w

12、indow lengthL=Nw/2;%窗函数每次移动的样点数, 重叠宽度%窗口函数subplot(3,2,2);g=hamming(Nw);j=-Nw/2;fori=1:Nwt(i)=j;j=j+1;endplot(t,g,'.');subplot(3,2,3);S,F,T,P = spectrogram(y,hamming(600),580,600,1E3);surf(T,F,10*log10(P),'edgecolor', 'none' );axistightview(0,90);xlabel('Time (Seconds)');ylabel('Hz' );subplot(3,2,4);surf(T,F,10*log10(P),'edgecolor' , 'none' );axis tightxlabel('Time (Seconds)');ylabel('Hz' );zlabel('enargy')Y=fft(y);X=fftshift(Y);%重构x1=ifftshift(X);x1=ifft(x1);subplot(3,2,5);plot(x