加窗插值FFT算法的matlab 512点仿真

修改程序能够实现加窗插值FFT算法的matlab仿真，已经得到实验验证，注解清晰详细，代码可读性强，以512点为例进行仿真计算，对与其他点数，如256或是1024可以通过简单修改即可

fs=6400;N=512;   %采样频率和数据点数
f0=50;%原信号频率
n=0:N-1;
t=n/fs;   %时间序列
s=10*sin（100*pi*t+pi/3）+0.3*sin（200*pi*t+pi/4）+...
+0.6*sin（300*pi*t+pi/6）+0.1*sin（400*pi*t+pi/3）+...
+0.5*sin（500*pi*t+pi/4）+0.1*sin（600*pi*t+pi/5）+...
0.2*sin（700*pi*t+pi/10）; %信号
w = 0.5-0.5*cos（2*pi*n./N）;%汉宁窗函数
x = s.*w;%加窗
figure（1）;
subplot（211）;
plot（tx）;
xlabel（‘t‘）;
ylabel（‘x‘）;
title（‘信号时域波形‘）;
grid;
y=fft（xN）;    %对信号进行快速Fourier变换
%wd=atan（imag（y（5））/real（y（5）））;
vz=abs（y）/256*2;%未经插值的谐波系数%在FFT后为了求幅值一般 进行计算abs（y）*2/N，N是FFT的点数。2是hanning窗的恢复系数。
f=（0:length（y）-1）‘*fs/length（y）;
figure（1）;
subplot（212）;
plot（fvz）;
axis（[01000010]）;
xlabel（‘频率（Hz）‘）;
ylabel（‘幅值‘）;
title（‘信号幅频谱图‘）;
grid;
u=abs（y）;%准备插值
A=zeros（17）;%修正后幅值存储
Ph1=zeros（17）;%有改动Ph
Ph2=zeros（17）;%有改动Ph
Ph=zeros（17）;%有改动Ph
%插值运算
for I=0:6
    y1=u（floor（4*I+5））;
    y2=u（floor（4*I+6））;
    yy2=y2;
    y3=u（floor（4*I+4））;
    yy3=y3;
    max=y2;
    if y3>y2
        max=y3;
    end
    if max==y3;
        t=y1;
        y1=max;
        y2=t; 
    end  
    b=（y2-y1）/（y2+y1）;
    a=1.5*b;
    A（I+1）=（y1+y2）*（2.35619403+1.15543682*a^2+0.32607873*a^4+0.07891461*a^6）/512;%谐波波幅值
        if yy2>=yy3
   Ph1（I+1）=atan（imag（y（4*I+6））/real（y（4*I+6）））+pi/2-pi*（a+0.5）;%谐波相角值 左减右加      
        end
        if yy2   Ph2（I+1）=atan（imag（y（4*I+4））/real（y（4*I+4）））+pi/2-pi*（a-0.5）;%谐波相角值 左减右加
        end
   Ph=Ph1+Ph2;
end

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       1573  2014-03-07 09:44  jiachuangchazhi512.m

----------- ---------  ---------- -----  ----

                 1573                    1