fft.cpp

基于VC编写正弦信号的傅立叶变换程序

#include "stdafx.h"
#include <math.h>
#include < stdio.h>
#include "fft.h"
#define pi 3.14159265358979323846

COMPX EE(COMPX b1,COMPX b2)
{
	COMPX b3;
	b3.real=b1.real*b2.real-b1.imag*b2.imag;
	b3.imag=b1.real*b2. imag +b1.imag*b2.real;
	return(b3);
}

void FFT(COMPX *xin,int N)
{
	int m,LH,nm,I,k,J,M,K;
	float p,ps;
	int B,N1;
	COMPX w,T;
	M=log10(N)/log10(2);
	LH=N/2;
	J=LH;
	N1=N-2;
	for(I=1;I<=N1;I++)
	{
		if(I<J)
		{
			T=xin[I];
			xin[I]= xin[J]; 
			xin[J]=T;
		}
		
		K=LH;
		while(J>=K)
		{
				J=J-K;
				K=K/2;
		}
		J=J+K;
	}
	
	for(m=1;m<=M;m++)
	{
		B=pow(2,m-1);
		for(J=0;J<=B-1;J++)
		{
			p=pow(2,M-m)*J;
			ps=2*pi/N*p;
			w.real=cos(ps);
			w.imag=-sin(ps);
			
			for(k=J;k<=N-1;k=k+pow(2,m))
			{
				T=EE(xin[k+B],w);
				xin[k+B].real=xin[k].real-T.real;
				xin[k+B].imag=xin[k]. imag -T.imag;
				xin[k].real=xin[k].real+T.real;
				xin[k].imag=xin[k]. imag +T.imag;
			}
		}
	}
}

void sine_generate(COMPX *xin,int N,float npp,float amplitude)
{
	for(int i=0;i<N;i++)
		xin[i].imag=0;
	if(npp-2.00<0.00001)
	{
		for(i=0;i<N;i++)
			xin[i].real=amplitude * pow(-1,i);
	}
	else
	{
	
		for(i=0;i<N;i++)
			xin[i].real=amplitude * sin(2*pi*i/npp);
	}

}