Continue to Site

Welcome to EDAboard.com

Welcome to our site! EDAboard.com is an international Electronics Discussion Forum focused on EDA software, circuits, schematics, books, theory, papers, asic, pld, 8051, DSP, Network, RF, Analog Design, PCB, Service Manuals... and a whole lot more! To participate you need to register. Registration is free. Click here to register now.

Register Log in

Looking for a simple core FFT algorithm code

Status
Not open for further replies.
F

fala

Full Member level 5
Full Member level 5
Joined
Sep 18, 2005
Messages
249
Helped
19
Reputation
38
Reaction score
4
Trophy points
1,298
Activity points
3,569
Hello, I'm developing an FFT algorithm for impedance calculations. I searched the net and found very long and sophisticated codes and libraries for that like https://www.fftw.org I don't want to use pre-written libraries. I just want to know the very basic idea so I can find out what differences with DFT it has so hopefully I can write an algorithm in C with a few dozen lines.
Can someone give me a hint what keyword or what site I should look for the core idea of FFT?

thanks alot.
 

simple fft algorithm

HI

there are many algorithms, look at sites such as this **broken link removed**

as long as I know, the fastest is https://www.fftw.org/

Sal
 

    F

    fala

    Points: 2
    Helpful Answer Positive Rating
simple fft

Check out DSP.com, and get the e-book Scientists and engineers guide to DSP. Yuo will get pseudo-code in it. Its pretty simple.
 

    F

    fala

    Points: 2
    Helpful Answer Positive Rating
simple FFT

you can search for some help on dspguide.com
 

simple FFT

Here's a simple (low performance) 16-bit complex FFT. I forget where I found it, and I may have rewritten some of it. You can change NFFT to some other power of two. It assumes 32-bit int and 16-bit short. The output array overwrites the input array. The array's real and imaginary values are interleaved i0, j0, i1, j1, i2, j2, ...

Code: 
#include <math.h>
#define PI 3.141592653589793238462643
#define NFFT 256

void fft(short *x)
{
  static short once=1, sint[3 * NFFT / 4];
  int n1, n2, ie, i, j, k, xt, yt;
  short *xi, *xl, *sinus;

  if (once)
  {
    for (i=0; i < 3*NFFT/4; i++)    /* sine table */
      sint[i] = 0x7FFF * sin(2 * PI * i / NFFT);
    once = 0;
  }

  for (n2=NFFT, ie=1; ie<NFFT; ie=ie<<1)
  {
    n1 = n2;
    n2 = n2 >> 1;
    sinus = sint;
    for (j=0; j < n2; j++)
    {
      for (i=j<<1; i < 2*NFFT; i+=n1<<1)
      {
        xi = x + i;
        xl = xi + n1;
        xt = (xi[0] - xl[0]) >> 1;
        yt = (xi[1] - xl[1]) >> 1;
        *(int*)xi = (((xi[0] + xl[0]) << 15) & 0xFFFF0000) | (unsigned short)((xi[1] + xl[1]) >> 1);
        *(int*)xl = (((xt * sinus[NFFT/4] + yt * sinus[0]) << 1) & 0xFFFF0000) | (unsigned short)((yt * sinus[NFFT/4] - xt * sinus[0]) >> 15);
      }
      sinus = sinus + ie;
    }
  }

  for (j=i=0; i < 2*NFFT; i+=2)   /* bit reverse */
  {
    if (i < j)
    {
      ie = *(int*)(x+i);
      *(int*)(x+i) = *(int*)(x+j);
      *(int*)(x+j) = ie;
    }
    k = NFFT;
    while (j >= k)
    {
      j -= k;
      k = k >> 1;
    }
    j += k;
  }
}
 

    F

    fala

    Points: 2
    Helpful Answer Positive Rating
Re: simple FFT

here is another 1D fft, i had to write for my image processing project. as for excellent explanations of algorithms (many), both in DIT and DIF, check out Inside fft blackbox book from CRC press.

Code: 
void hiwi::fft1(Vec_IO_DP &inpR, Vec_IO_DP &inpI, int Inv)

{
	
	int NumberOfProblems = 1;
	int ProblemSize = inpR.size();
	int N = inpR.size();
	int Distance = 1;
	bool NotSwitchInput = 1;
	int i;
	int JFirst = 0;
	int J, K, JTwiddle;
	double WR, WI;
	Vec_DP wR(N), wI(N), aR(N), aI(N), bR(N), bI(N);
	
	for(i=0; i<N; i++) {
		wR[i] = 0;
		wI[i] = 0;
		aR[i] = inpR[i];
		aI[i] = inpI[i];
		bR[i] = 0;
		bI[i] = 0;
	}
	for(i=0; i<N; i++) {
		wR[i] = cos(2* pi * i / N);
		if(Inv == 1)
			wI[i] = -sin(2* pi * i / N);
		else
			wI[i] = +sin(2* pi * i / N);
	}
	while(ProblemSize > 1) {
		if(NotSwitchInput) {
			for(JFirst = 0; JFirst < NumberOfProblems; JFirst++) {
				J = JFirst;
				JTwiddle = 0;
				K = JFirst;
				while( J < N - 1 ) {
					WR = wR[JTwiddle];
					WI = wI[JTwiddle];
					bR[J] = aR[K] + aR[K + N/2];
					bI[J] = aI[K] + aI[K + N/2];
					bR[J + Distance] = (aR[K] - aR[K + N/2]) * WR - (aI[K] - aI[K + N/2]) * WI;
					bI[J + Distance] = (aR[K] - aR[K + N/2]) * WI + (aI[K] - aI[K + N/2]) * WR;
					JTwiddle += NumberOfProblems;
					J += 2 * NumberOfProblems;
					K += NumberOfProblems;
				}
			}
			NotSwitchInput = 0;
		}
		else {
			for(JFirst = 0; JFirst < NumberOfProblems;JFirst++) {
				J = JFirst;
				JTwiddle = 0;
				K = JFirst;
				while(J < N - 1) {
					WR = wR[JTwiddle];
					WI = wI[JTwiddle];
					aR[J] = bR[K] + bR[K + N/2];
					aI[J] = bI[K] + bI[K + N/2];
					aR[J + Distance] = (bR[K] - bR[K+N/2]) * WR - (bI[K] - bI[K + N/2]) * WI;
					aI[J + Distance] = (bR[K] - bR[K+N/2]) * WI + (bI[K] - bI[K + N/2]) * WR;
					JTwiddle += NumberOfProblems;
					J += 2 * NumberOfProblems;
					K += NumberOfProblems;
				}
			}
			NotSwitchInput = 1;
		}
		NumberOfProblems *= 2;
		ProblemSize /= 2;
		Distance *= 2;
	}
	double factor = log((double) N) / log((double) 2);
	if ( fmod(factor,(double) 2) == 0) {
		for(i = 0; i<N; i++) {
			inpR[i] = aR[i];
			inpI[i] = aI[i];
		}
	}
	else {
		for(i = 0; i<N; i++) {
			inpR[i] = bR[i];
			inpI[i] = bI[i];
		}
	}
	
}

the first and the second parameters are the input vector real and imaginary part and the third parameter tells whether to calculate fft or ifft. a 1 calculates fft and 0 or any other number calculates ifft.

2D fft should be straight forward from this as u can take the 1D fft of the rows and from the resulting matrix, take 1D fft of columns. hope this might be of use to some one. however, the function input specified here uses matrix templates which i wrote. so, specify the i/p and vector/matrix declaration yourself, as to how u usually do it.

regards,
cedance.
 

    F

    fala

    Points: 2
    Helpful Answer Positive Rating
Status
Not open for further replies.

Similar threads

Part and Inventory Search

Welcome to EDABoard.com

Sponsor

Back
Top