Click here to Skip to main content
15,063,556 members
Please Sign up or sign in to vote.
0.00/5 (No votes)
Hello guys,
I'm having a problem with my coding. Currently I am working on a evaluation tool. This tool has to do some frequency analysis. Therefore I bought a book and read quite some articles about the FFT in C# and finally got some code to calculate it.

For a small amount of samples, e.g. when I calculate a sin and put it in the algorithm, everything works fine. BUT if I put real data in it (ca. 93k samples @ 500 Hz) the output is weird...

Here is the result of my FFT vs. the result of the FFT of NI DIAdem:
<img src="http://imgur.com/9w8i4au">

Here the Code:

C#
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace SignalAnalysis
{
    public class FFT
    {
        // Declaration of Variables
        private Complex[] data;
        const double PI = Math.PI;

        // new FFT Operator
        public FFT(ref Complex[] Input_Data)
        {
            data = new Complex[(ulong)Input_Data.LongCount()];
            int i = 0;
            foreach (Complex clx in Input_Data)
            {
                data[i] = new Complex(clx.re, clx.im);
                i++;
            }
            //Array.Copy(Input_Data, data, Input_Data.LongCount());
            if (IsPowerOfTwo((ulong)data.LongCount()) == false)
                data = MakePowerOfTwo(data);
        }
        
        // Check if number is power of 2
        bool IsPowerOfTwo(ulong x)
        {
            return (x & (x - 1)) == 0;
        }

        // Make power of 2
        Complex[] MakePowerOfTwo(Complex[] input)
        {
            double root = Math.Log(input.LongCount(),2.00);
            int pow = (int)root + 1;
            ulong new_cnt = (ulong)Math.Pow(2, pow);

            Complex[] output = new Complex[new_cnt];
            Complex clx = new Complex(0, 0);
            output = Enumerable.Repeat(clx, output.Count()).ToArray();
            Array.Copy(input, output, input.LongCount());

            //for (long i = input.LongCount(); i < output.LongCount(); i++)
            //    output[i] = new Complex(0, 0);

            return output;
        }

        // Calculate the FFT
        public bool calculate()
        {
            int Abtastwerte = data.Count();
            double tempr = 0;
            double tempi = 0;
            double wreal = 0;
            double wimag = 0;
            double real1 = 0;
            double imag1 = 0;
            double real2 = 0;
            double imag2 = 0;
            int i = 0;
            int j = 0;
            int k = 0;
            int stufen = 0;
            int sprung = 0;
            int schritt = 0;
            int element = 0;

            for (j = 0; j < Abtastwerte - 1; j++)
            {
                if (j < i)
                {
                    tempr = data[j].re;
                    tempi = data[j].im;
                    data[j].re = data[i].re;
                    data[j].im = data[i].im;
                    data[i].re = tempr;
                    data[i].im = tempi;
                }
                k = Abtastwerte / 2;
                while (k <= i)
                {
                    i -= k;
                    k /= 2;
                }
                i += k;
            }
            
            stufen = (int)(Math.Log10((double)Abtastwerte) / Math.Log10((double)2));
            sprung = 2;

            for (i = 0; i < stufen; i++)
            {
                element = 0;
                for (j = Abtastwerte / sprung; j >= 1; j--)
                {
                    schritt = sprung / 2;
                    for (k = 0; k < schritt; k++)
                    {
                        wreal = Math.Cos(k * 2.0 * PI / sprung);
                        wimag = Math.Sin(k * 2.0 * PI / sprung);

                        real1 = data[element + k].re +
                                (wreal * data[element + k + schritt].re -
                                wimag * data[element + k + schritt].im);

                        imag1 = data[element + k].im +
                                (wreal * data[element + k + schritt].im +
                                wimag * data[element + k + schritt].re);

                        wreal *= -1.0;
                        wimag *= -1.0;

                        real2 = data[element + k].re +
                                (wreal * data[element + k + schritt].re -
                                wimag * data[element + k + schritt].im);

                        imag2 = data[element + k].im +
                                (wreal * data[element + k + schritt].im +
                                wimag * data[element + k + schritt].re);

                        data[element + k].re = real1;
                        data[element + k].im = imag1;
                        data[element + k + schritt].re = real2;
                        data[element + k + schritt].im = imag2;
                    }
                    element += sprung;
                }
                sprung *= 2;
            }
            return true;
        }

        // Retrieve the Result
        public Complex[] getResult()
        {
            return data;
        }

        // Retrieve the frequency channel
        public double[] getFrequencyCh(double SamplingRate)
        {
            double[] frequency = new double[data.Count()];
            double coefficient = (SamplingRate / 2.00) / ((double)data.Count() / 2.00);
            for (int i = 0; i < data.Count(); i++)
                frequency[i] = i * coefficient;
            return frequency;
        }

    }
}
Posted
Comments
Jochen Arndt 22-Jan-16 9:29am
   
Did you implement the code yourself?

Than it is a quite time consuming task to find where the codes differs from the used algorithm. So you might show that too and somebody here has the time to check it. But you should not rely on that and try to find the error yourself.

Another option is using existing code (library or source). Just search for "c# fft [code|library]".

This content, along with any associated source code and files, is licensed under The Code Project Open License (CPOL)

  Print Answers RSS
Top Experts
Last 24hrsThis month



CodeProject, 20 Bay Street, 11th Floor Toronto, Ontario, Canada M5J 2N8 +1 (416) 849-8900