加窗傅里叶变换,加hamming窗,序列是离散的,存在数值A中可能有1000左右的数值在A中,怎么做它的加窗傅里叶变换。做傅里叶变换时不能利用fft快速傅里叶变换,因为不是二的整次幂。我写了个小程序,可是有问题,时频图也画不出来
clc;clear;
fid=fopen('C:\Documents and Settings\Administrator\桌面\Matlab数据程序\必需基因\M65145.txt','r');%里面是1072个碱基序列。
b=fscanf(fid,'%s');
%l=length(a);
%b=a(1:1638);
N=length(b);
fclose(fid);
A=zeros(1,N);
for i=1:1:N
if b(i)=='a';
A(i)=0.1260;
else
if b(i)=='t';
A(i)=0.1335;
else
if b(i)=='g';
A(i)=0.0806;
else
if b(i)=='c';
A(i)=0.1340;
end
end
end
end
end
M=109;
%Lap=M/2;
Lap=90;
Lshift=M-Lap;
w=hamming(M);
Nshift=floor((length(A)-Lap)/Lshift);% 计算窗函数移动次数
tf=zeros(M/2,Nshift);% 生成时频矩阵
for i=0:Nshift-1;
xi=A((i*Lshift+1) i*Lshift+M));% 截取数据段
xw=xi.*w'; %数据相乘
n=0 M-1);k=0 M-1);
temp=xw*exp(-j*2*pi/M).^(n'*k);%离散傅里叶变换
tf(:,i+1)=temp(1:M/2);% 将频率为0-fs/2范围内的频率点数傅里叶值存入时频矩阵
end
tf_t=(0:Nshift-1)*Lshift;% 计算时间
tf_f=(0:M/2-1)/M;% 计算频率
[T,F]=meshgrid(tf_t,tf_f);
contour(T,F,abs(tf));
![]()
未命名.jpg |