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.

Looking for a simple core FFT algorithm code

Status
Not open for further replies.

fala

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.
 

cydi

Full Member level 3
Joined
Dec 12, 2006
Messages
165
Helped
14
Reputation
28
Reaction score
3
Trophy points
1,298
Activity points
2,004
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.
 

    fala

    Points: 2
    Helpful Answer Positive Rating

amit_8561

Member level 2
Joined
Apr 18, 2006
Messages
43
Helped
11
Reputation
22
Reaction score
3
Trophy points
1,288
Activity points
1,467
simple FFT

you can search for some help on dspguide.com
 

echo47

Advanced Member level 5
Joined
Apr 7, 2002
Messages
3,942
Helped
637
Reputation
1,272
Reaction score
88
Trophy points
1,328
Location
USA
Activity points
33,178
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;
  }
}
 

    fala

    Points: 2
    Helpful Answer Positive Rating

cedance

Advanced Member level 2
Joined
Oct 24, 2003
Messages
551
Helped
30
Reputation
60
Reaction score
7
Trophy points
1,298
Location
Germany
Activity points
4,622
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.
 

    fala

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

Similar threads

Part and Inventory Search

Welcome to EDABoard.com

Sponsor

Top