Make your own free website on Tripod.com

To see the Dsp1.h Header File
To see the dft.h Header File
To see the makefile
To see the input
To see the output
To go back to index
To go back to homepage
/*------------------------------------------------------------------------*/
/*									  */
/*			Program to compute DFT coefficients directly	  */
/*			3 other functions are used			  */
/*									  */
/*									  */
/*			Program design by: E C Ifeachor, July, 1992 	  */
/*									  */
/*			Program written  and modified 			  */
/*			by: Ronald Bazillion, 7/22/96	     		  */
/*									  */
/*------------------------------------------------------------------------*/

#include "dsp1.h"
#include "dft.h"

main()
{
	extern long npt;
	extern int inv;
	printf("select type of transform\n");
	printf("\n");
	printf("0 for forward DFT\n");
	printf("1 for inverse DFT\n");
	scanf("%d",&inv);
	read_data();
	dft();
	save_data();
	exit();
}
/*--------------------------------------------------------------*/
/* 	Function to read the data from a file			*/
/*								*/
/*  	Program designed by E C Ifeachor July, 1992		*/
/*								*/
/*	Written by Ronald Bazillion 7/22/96			*/
/*								*/
/*--------------------------------------------------------------*/

void read_data()
{
	extern long npt;
	int n;
	extern complex x[size];

	for (n=0;n< size; ++n)
	{
		x[n].real=0;
		x[n].imag=0;
	}
	if ((in=fopen("coeff.dat","r"))==NULL)
	{
		printf("cannot open file coeff.dat\n");
		exit(1);
	}
	fscanf(in,"%ld",&npt);
	for(n=1; n<=npt; n++)
	{
		fscanf(in,"%lf %lf",&x[n].real,&x[n].imag);
	}
	fclose(in);
}

/*--------------------------------------------------------------*/
/* 	Function to save the data to a file 			*/
/*								*/
/*  	Program designed by E C Ifeachor July, 1992		*/
/*								*/
/*	Written by Ronald Bazillion 7/22/96			*/
/*								*/
/*--------------------------------------------------------------*/

void save_data()
{
	long k;
	int k1;
	extern  long npt;
	extern complex x[size];
	if((out=fopen("dftout.dat","w"))==NULL)
	{
		printf("Cannot open file dftout.dat\n");
		exit(1);
	}
	fprintf(out,"k \t\t XR(k) \t\t XI(k) \n");
	fprintf(out,"\n");
	for (k=1; k<=npt; ++k)
	{
		k1=k-1;
		fprintf(out, "%d \t %f \t %f \n", k1, x[k].real, x[k].imag);
	}
	fclose(out);
}
/*------------------------------------------------------------------*/
/*								    */
/*	Function to compute the DFT of a discrete-time              */
/*	sequence directly					    */
/*								    */
/*								    */
/*      Program designed by: E C Ifeachor 31.10.91		    */
/*								    */
/*	Program written by: Ronald Bazillion 7/22/96		    */
/*								    */
/*------------------------------------------------------------------*/
void dft()
{
	extern	int inv;
	extern  long npt;
	long k,n;
	double WN, wk, c, s, XR[size], XI[size];
	extern complex x[size];

	WN=2*pi/npt;
	if(inv==1)
		WN=-WN;
	for (k=0;k<npt;++k)
	{
		XR[k]=0.0;
		XI[k]=0.0;
		wk=k*WN;
		for (n=0; n<npt;++n)
		{
		c=cos(n*wk);
		s=sin(n*wk);
		XR[k]=XR[k]+x[n+1].real*c+x[n+1].imag*s;
		XI[k]=XI[k]-x[n+1].real*s+x[n+1].imag*c;
	      	}
		if (inv==1)
		{	/* divide by N for IDFT */
			
			XR[k]=XR[k]/npt;
			XI[k]=XI[k]/npt;
		}
	}
	for (k=1;k<=npt; ++k)
	{
		/* store transformed data in x */
	
		x[k].real=XR[k-1];
		x[k].imag=XI[k-1];
	}
}