fft.c

对FFT算法的基2时析型实现

#include <ansi_c.h>
#include <cvirte.h>		/* Needed if linking in external compiler; harmless otherwise */
#include <userint.h>
#include <math.h>
#include "20041001.h"
#define PI 3.141592654

static int panelHandle;
int n,L,N,z,i=0,j=0,k,p,h,z,x,ab,Y,K1,K2,K3,JJ;
static int n1,n2;
int dxk,dxt,dk;
double *c,*cj,*e,*ej,*f,*fj,*g,*gj;
static double *a,*aj,*b,*bj;
double d,jh,xs1,xs2,xs3;
double w,wj,zh,zhj,zz,zzj;
int main (int argc, char *argv[])
{
	if (InitCVIRTE (0, argv, 0) == 0)	/* Needed if linking in external compiler; harmless otherwise */
		return -1;	/* out of memory */
	if ((panelHandle = LoadPanel (0, "20041001.uir", PANEL)) < 0)
		return -1;
	DisplayPanel (panelHandle);
	RunUserInterface ();
	return 0;
}

int CVICALLBACK Canshu (int panel, int control, int event,
		void *callbackData, int eventData1, int eventData2)
{
	switch (event)
		{
		case EVENT_COMMIT:
			DeleteGraphPlot (panelHandle, PANEL_GRAPH, -1, VAL_IMMEDIATE_DRAW);
		    GetCtrlVal(panelHandle,PANEL_SRXZ,&ab);
		    GetCtrlVal(panelHandle,PANEL_DIANSHU,&n);
		    if(n<=0)
		    {
		    n=16;
		    SetCtrlVal(panelHandle,PANEL_DIANSHU,n);
		    }
		    if(ab==0)
		    {
		    n1=n;
		    }
		    else
		    {
		    n2=n;
		    }
		    d=n;
			h=0;
			while(d>1)
			{
			 d=d/2;
			 h++;
			}
			if(d<1.0)
			{
			z=1;
			for(k=0;k<h;k++)
			 z=z*2;
			 n=z;
			}
			SetCtrlVal(panelHandle,PANEL_DIANSHU,n);          //使n为2的整数倍
			if(ab==0)
			{
			a=malloc(n*sizeof(double));
            aj=malloc(n*sizeof(double));
            }
            else
            {
            b=malloc(n*sizeof(double));
            bj=malloc(n*sizeof(double));
            }
            if(ab==0)
            {
             for(k=0;k<n;k++)
             {
             a[k]=0;aj[k]=0;
             } 
            }
            else
            {
             for(k=0;k<n;k++)
             {
             b[k]=0;bj[k]=0;
             } 
            }
            i=0;j=0; 
			break;
		}
	return 0;
}

int CVICALLBACK INPUT (int panel, int control, int event,
		void *callbackData, int eventData1, int eventData2)
{
	switch (event)
		{
		case EVENT_COMMIT:
			
			if(ab==0)
				{   
				GetCtrlVal(panelHandle,PANEL_SHIBU,&a[i]);
				GetCtrlVal(panelHandle,PANEL_XUBU,&aj[i]);
				i++;
			    SetCtrlVal(panelHandle,PANEL_YISHUDIANSHU,i);
				}
			else
				{
				GetCtrlVal(panelHandle,PANEL_SHIBU,&b[j]);
				GetCtrlVal(panelHandle,PANEL_XUBU,&bj[j]);
				j++;
				SetCtrlVal(panelHandle,PANEL_YISHUDIANSHU,j); 
				}
			break;
			 
		}
	return 0;
}   

void GFFT(double *c,double *cj,int N,int Y)
{
            //h=log(N)/log(2);
            h=0;
            d=N;
            while(d>1)
            {
            d=d/2;
            h++;
            }
            k=1;
		    for(i=1;i<N;i++)
		    {
		     z=0;
		     for(j=1;j<=h;j++)
		     {
		      p=pow(2,j);
		      z=z+N/p;
		      if(k<=z)
		      {
		       k=k+2*N/p-z;
		       if(i<k)
		        {
		         jh=c[i],c[i]=c[k-1],c[k-1]=jh;
		         jh=cj[i],cj[i]=cj[k-1],cj[k-1]=jh;
		        }
		        break;
		      }
		     }
			}											//排序正确
		    for(i=1;i<=h;i++)
		    {
		      dxk=N/pow(2,i);
		      dxt=pow(2,i-1);
		      dk=pow(2,i);
		      x=0;
		      for(j=0;j<dxk;j++) 
		      {
		        for(k=0;k<dxt;k++)
		        {
		         w=cos(2*PI*k/dk);
		         if(Y==0)
		         wj=(-1)*sin(2*PI*k/dk);
		         else
		         wj=sin(2*PI*k/dk);
		         zh=c[x+dxt]*w-cj[x+dxt]*wj;
		         zhj=c[x+dxt]*wj+cj[x+dxt]*w;
		         zz=c[x];
		         zzj=cj[x];
		         c[x]=zz+zh;
		         cj[x]=zzj+zhj;
		         c[x+dxt]=zz-zh;
		         cj[x+dxt]=zzj-zhj;
		          if((k+1)<dxt)
		          {
		          x++;
		          }
		        }
		        x=x+dxt+1;
		      } 
		    }
		     if(Y==1)
             for(z=0;z<N;z++)
		     {
		     c[z]=c[z]/N;
		     cj[z]=cj[z]/N;
		     } 
}
int CVICALLBACK MFFT (int panel, int control, int event,
		void *callbackData, int eventData1, int eventData2)
{
	switch (event)
		{
		case EVENT_COMMIT:
			 JJ=0;
             Y=0;
             if(ab==0)
             {
			 GFFT(a,aj,n,Y);
		     }
		     else
		     {
			 GFFT(b,bj,n,Y);
		     }
			break;
		}
	return 0;
}

int CVICALLBACK MIFFT (int panel, int control, int event,
		void *callbackData, int eventData1, int eventData2)
{
	switch (event)
		{
		case EVENT_COMMIT:
			 JJ=0;
			 Y=1;
             if(ab==0)
             {
			 GFFT(a,aj,n,Y);
		     }
		     else
		     {
			 GFFT(b,bj,n,Y);
		     }
            
			break;
		}
	return 0;
}

int CVICALLBACK ValueSHOW (int panel, int control, int event,
		void *callbackData, int eventData1, int eventData2)
{
	switch (event)
		{
		case EVENT_COMMIT:
		GetCtrlVal(panelHandle,PANEL_NUMERIC1,&K1);
		GetCtrlVal(panelHandle,PANEL_NUMERIC2,&K2);
		GetCtrlVal(panelHandle,PANEL_NUMERIC3,&K3); 
	    if(JJ==1)
		{
		if(K1<0||K1>=L)
		  {
		  K1=0;
		  SetCtrlVal(panelHandle,PANEL_NUMERIC1,K1);
		  }
		if(K2<0||K2>=L)
		  {
		  K2=0;
		  SetCtrlVal(panelHandle,PANEL_NUMERIC2,K2);
		  }
		if(K1<0||K1>=L)
		  {
		  K3=0;
		  SetCtrlVal(panelHandle,PANEL_NUMERIC3,K3);
		  }
		 xs1=g[K1]*g[K1]+gj[K1]*gj[K1];
		 xs2=g[K2]*g[K2]+gj[K2]*gj[K2];
		 xs3=g[K3]*g[K3]+gj[K3]*gj[K3];
		}
		else
		{
		if(K1<0||K1>=n)
		  {
		  K1=0;
		  SetCtrlVal(panelHandle,PANEL_NUMERIC1,K1);
		  }
		if(K2<0||K2>=n)
		  {
		  K2=0;
		  SetCtrlVal(panelHandle,PANEL_NUMERIC2,K2);
		  }
		if(K1<0||K1>=n)
		  {
		  K3=0;
		  SetCtrlVal(panelHandle,PANEL_NUMERIC3,K3);
		  }
		if(ab==0)
		 {
		 xs1=a[K1]*a[K1]+aj[K1]*aj[K1];
		 xs2=a[K2]*a[K2]+aj[K2]*aj[K2];
		 xs3=a[K3]*a[K3]+aj[K3]*aj[K3];
		 }
		else
		 {
		 xs1=b[K1]*b[K1]+bj[K1]*bj[K1];
		 xs2=b[K2]*b[K2]+bj[K2]*bj[K2];
		 xs3=b[K3]*b[K3]+bj[K3]*bj[K3];
		 }
		}
		 SetCtrlVal(panelHandle,PANEL_NUMERIC4,xs1);
		 SetCtrlVal(panelHandle,PANEL_NUMERIC5,xs2);
		 SetCtrlVal(panelHandle,PANEL_NUMERIC6,xs3);
			break;
		}
	return 0;
}

int CVICALLBACK SHOW (int panel, int control, int event,
		void *callbackData, int eventData1, int eventData2)
{
	switch (event)
	{
		case EVENT_COMMIT:
		DeleteGraphPlot (panelHandle, PANEL_GRAPH, -1, VAL_IMMEDIATE_DRAW);
		if(ab==0)
        {
		PlotY (panelHandle, PANEL_GRAPH, a, n, VAL_DOUBLE,VAL_THIN_LINE, VAL_EMPTY_SQUARE, VAL_SOLID, 1, VAL_RED);
		PlotY (panelHandle, PANEL_GRAPH, aj, n, VAL_DOUBLE,VAL_THIN_LINE, VAL_EMPTY_SQUARE, VAL_SOLID, 1, VAL_YELLOW);
		}
		else
		{
		PlotY (panelHandle, PANEL_GRAPH, b, n, VAL_DOUBLE,VAL_THIN_LINE, VAL_EMPTY_SQUARE, VAL_SOLID, 1, VAL_RED);
		PlotY (panelHandle, PANEL_GRAPH, bj, n, VAL_DOUBLE,VAL_THIN_LINE, VAL_EMPTY_SQUARE, VAL_SOLID, 1, VAL_YELLOW);
		}
			break;
	}
	return 0;
}

int CVICALLBACK JUANJI (int panel, int control, int event,void *callbackData, int eventData1, int eventData2)
{
	switch (event)
		{
		case EVENT_COMMIT:
			JJ=1;
		    L=n1+n2-1;
			d=L;
			h=0;
			while(d>1)
			{
			 d=d/2;
			 h++;
			}
			if(d<1.0)
			{
			z=1;
			for(k=0;k<h;k++)
			 z=z*2;
			 L=z;
			}
			e=malloc(L*sizeof(double));
            ej=malloc(L*sizeof(double));
            f=malloc(L*sizeof(double));
            fj=malloc(L*sizeof(double));
            g=malloc(L*sizeof(double));
            gj=malloc(L*sizeof(double));
            for(i=0;i<L;i++)
			{
			e[i]=0; ej[i]=0; f[i]=0; fj[i]=0;		 //得到的n1,n2值正确
			}
			for(i=0;i<n1;i++)
			{
			e[i]=a[i];
			ej[i]=aj[i];
			}
			for(i=0;i<n2;i++)
			{
			f[i]=b[i];
			fj[i]=bj[i];
			}									   //得到的e,ej,f,fj正确
			Y=0;
			GFFT(e,ej,L,Y);
			GFFT(f,fj,L,Y);
			for(i=0;i<L;i++)
			{
			 g[i]=e[i]*f[i]-ej[i]*fj[i];
			 gj[i]=e[i]*fj[i]+ej[i]*f[i];
			}
			Y=1;
			GFFT(g,gj,L,Y);
			DeleteGraphPlot (panelHandle, PANEL_GRAPH, -1, VAL_IMMEDIATE_DRAW);
			PlotY (panelHandle, PANEL_GRAPH, g, L, VAL_DOUBLE,VAL_THIN_LINE, VAL_EMPTY_SQUARE, VAL_SOLID, 1, VAL_RED);
		    PlotY (panelHandle, PANEL_GRAPH, gj, L, VAL_DOUBLE,VAL_THIN_LINE, VAL_EMPTY_SQUARE, VAL_SOLID, 1, VAL_YELLOW);
			break;
		}
	return 0;
}

int CVICALLBACK QUIT (int panel, int control, int event,
		void *callbackData, int eventData1, int eventData2)
{
	switch (event)
		{
		case EVENT_COMMIT:
			QuitUserInterface (0);
			break;
		}
	return 0;
}