2
$\begingroup$

While working with FFT, I have a strange case with a current experimental setup. I am working with beat frequencies (intermediate frequency output) and using a standard FFT algorithm (complex to complex). My sampled data are real-valued so I created a double size array with alternating real and complex values (with complex values set to zero) as in the figure below. The particular implementation of FFT that I used (LomontFFT) requires that to work around the lack of complex-value-type handling.

enter image description here

Here the ADC samples 400 data points and then 112 bits are zero padded after Hanning to make 512 data points.

Here is the code I have written (C#):

public void fftElavualtion (double[] adcValues)
{
        int N = 512;
        double Fs = 195312.5;
        double[] data_sample = Enumerable.Repeat(0.0, N).ToArray();
        double[] data_sampleHanned = new double[adcValues.Length];
        double[] fftData = new double[N * 2];
        double[] magnitude = new double[N];
        double max_magnitude = double.NegativeInfinity;
        int max_index = -1;
        double fundamentalFrequency;

        data_sampleHanned = HanningWindow(adcValues);

        Array.Copy(data_sampleHanned, data_sample, data_sampleHanned.Length);

        for (int i = 0; i <= (N - 1); i++)
        {
            fftData[2 * i] = data_sample[i];
            fftData[(2 * i) + 1] = 0;
        }

        var fftMethod = new Lomont.LomontFFT();
        fftMethod.FFT(fftData, true);

        for (int y = 0; y <= (N - 1); y++)
        {
            double re = fftData[2 * y];
            double im = fftData[2 * y + 1];
            magnitude[y] = Math.Sqrt(re * re + im * im);
        }

        for (int k = 1; k <= (N - 1); k++)
        {
            if (magnitude[k] > max_magnitude)
            {
                max_magnitude = magnitude[k];
                max_index = k;
            }
        }

        fundamentalFrequency = max_index * Fs / N;

        textBoxFundaFreq.Text = "Freq = " + fundamentalFrequency.ToString();
        for (int x = 1; x < N; x++)
        {
            this.chart2.Series[0].Points.AddXY(x, magnitude[x]);
        }
    }

The problems are:

  1. Each frequency bin represents a bandwidth of (Sample Rate)/(total number of samples) so 200,000/512. i.e., approx 390 Hz. Is this too high for a bin? If so how do I manage as I can have only 400 samples and also can't change sampling rate (system constraint). Does zero padding improve?

  2. To find the fundamental frequency I am locating the index of highest magnitude (max_index) and multiplying it with badwidth of each bin. This yields me the same frequency for all trials as second place of magnitude is always high irespective to what signal is fed to the processing system.

What is a possible reason for always the second bin having the highest energy or highest magnitude?

Note: I am ignoring the case of a problem in FFT, as LomontFFT is trusted. Info: (Example double valued sampled with ADC - 400 data points in CSV)

 2047.5,2063.88,2082.99,2099.37,2115.75,2140.32,2149.875,2156.7,2167.62,2178.54,2194.92,2205.84,2214.03,2222.22,2231.775,2238.6,2184,2115.75,2061.15,2013.375,1951.95,1897.35,1835.925,1774.5,1842.75,1911,2033.85,2115.75,2184,2293.2,2375.1,2457,2375.1,2252.25,2129.4,1992.9,1842.75,1706.25,1610.7,1501.5,1638,1777.23,1919.19,2059.785,2199.015,2336.88,2476.11,2613.975,2457,2297.295,2170.35,1904.175,1636.635,1449.63,1365,1236.69,1569.75,1835.925,1979.25,2080.26,2306.85,2457,2605.785,2852.85,2730,2661.75,2530.71,2429.7,2306.85,2182.635,2115.75,2047.5,1911,1774.5,1638,1501.5,1365,1228.5,1092,955.5,1092,1228.5,1365,1501.5,1638,1774.5,1911,2047.5,2197.65,2347.8,2497.95,2648.1,2798.25,2948.4,3091.725,3228.225,3084.9,2934.75,2784.6,2634.45,2484.3,2334.15,2190.825,2047.5,1774.5,1528.8,1228.5,1092,962.325,832.65,709.8,655.2,791.7,941.85,1090.635,1283.1,1528.8,1760.85,1911,2047.5,2293.2,2455.635,2593.5,2852.85,3016.65,3289.65,3412.5,3549,3412.5,3262.35,3001.635,2728.635,2455.635,2319.135,2115.75,2047.5,1794.975,1542.45,1365,1078.35,941.85,764.4,477.75,273,526.4805,682.5,832.65,1090.635,1337.7,1501.5,1909.635,2047.5,2184,2443.35,2661.75,2989.35,3274.635,3562.65,3701.88,3842.475,3685.5,3535.35,3207.75,2934.75,2730,2525.25,2306.85,2047.5,1842.75,1490.58,1119.3,982.8,709.8,436.8,163.8,40.95,204.75,436.8,764.4,1078.35,1300.845,1610.7,1883.7,2047.5,2320.5,2613.975,2907.45,3194.1,3480.75,3753.75,3890.25,4093.635,3549,2866.5,2184,1842.75,1706.25,1228.5,682.5,1.092,271.635,464.1,750.75,955.5,1228.5,1569.75,1842.75,2047.5,2320.5,2593.5,2866.5,3139.5,3412.5,3685.5,3822,3999.45,3822,3685.5,3412.5,3139.5,2866.5,2593.5,2320.5,2047.5,1774.5,1501.5,1228.5,955.5,682.5,518.7,313.95,135.135,313.95,518.7,682.5,955.5,1228.5,1501.5,1774.5,2047.5,2252.25,2525.25,2757.3,3057.6,3344.25,3494.4,3671.85,3842.475,3671.85,3494.4,3344.25,3057.6,2757.3,2525.25,2252.25,2047.5,1774.5,1501.5,1228.5,955.5,819,696.15,477.75,375.375,477.75,696.15,819,955.5,1228.5,1501.5,1774.5,2047.5,2320.5,2593.5,2866.5,3016.65,3139.5,3276,3426.15,3617.25,3426.15,3276,3139.5,3016.65,2866.5,2593.5,2320.5,2047.5,1829.1,1556.1,1255.8,1146.6,955.5,873.6,764.4,679.77,846.3,955.5,1201.2,1351.35,1501.5,1706.25,1842.75,2047.5,2184,2347.8,2525.25,2702.7,2866.5,3003,3139.5,3288.285,3139.5,3003,2866.5,2702.7,2525.25,2347.8,2184,2047.5,1911,1774.5,1638,1501.5,1365,1228.5,1092,1044.225,1092,1228.5,1365,1501.5,1638,1774.5,1911,2046.135,2184,2317.77,2454.27,2525.25,2593.5,2689.05,2852.85,2921.1,2852.85,2689.05,2593.5,2525.25,2454.27,2317.77,2184,2047.5,1965.6,1815.45,1733.55,1636.635,1474.2,1228.5,1351.35,1481.025,1562.925,1644.825,1726.725,1808.625,1890.525,1945.125,2027.025,2047.5,2095.275,2145.78,2197.65,2245.425,2293.2,2327.325,2361.45,2388.75,2361.45,2327.325,2293.2,2245.425,2197.65,2145.78,2095.275,2047.5,2020.2,1992.9,1965.6,1938.3,1911,1883.7,1829.1,1774.5,1842.75,1911,1979.25,2047.5,2115.75,2184,2115.75,2047.5
Improve this question
$\endgroup$
8
  • $\begingroup$ Your illustration says you're doing a 512-point FFT, but it looks like it should be 1024 based on how you're doing the processing. Is that a typo? Also, if you're doing this in MATLAB or other similar environment, it might help to provide the source code for your simulation to see if you have any bugs. $\endgroup$
    – Jason R
    Jul 14, 2015 at 15:18
  • $\begingroup$ @JasonR I am using a 512 point FFT, but since I cange the real values to real complex pairs, I just double the FFT length. 2 consequite value (Real and Comples) makes a single value or data point. So the FFT calculation limits to 512 points. $\endgroup$
    – dDebug
    Jul 14, 2015 at 15:24
  • $\begingroup$ Can you reword your title as a question and rework your explanation a bit: I don't understand what the problem is. I also don't understand if it is you how does the "pseudo complex" transformation. If you do, DON'T: that's done by your programming language. What you describe is a possible memory layout for complex numbers, but the algorithm need to know what is the real and what is the imaginary part. $\endgroup$
    – user13706
    Jul 14, 2015 at 15:42
  • 1
    $\begingroup$ Also, in your new, reworked, explanation, put some code or describe what you are doing as opposed to what you suppose the computer is doing. You could also provide a couple of plot to show the problem before you explain why it is a problem. $\endgroup$
    – user13706
    Jul 14, 2015 at 15:43
  • $\begingroup$ @dDebug: I understand now; I was confused by the way you presented the problem above. I wouldn't say you have 1024 points at any given time. Instead, you have 512 points, either real or complex. $\endgroup$
    – Jason R
    Jul 14, 2015 at 17:05

1 Answer 1

Reset to default
0
$\begingroup$

Ok, so I tried to perform the same calculation (but forgive me, I did it in Python, but it should be readable).

import numpy as np                    # Linear algebra module
from matplotlib import pyplot as plt  # Plot module
from scipy import signal              # One common DSP module

# Spectral density calculation
data = np.array([2047.5,2063.88,2082.99,2099.37,2115.75,2140.32,2149.875,2156.7,2167.62,2178.54,2194.92,2205.84,2214.03,2222.22,2231.775,2238.6,2184,2115.75,2061.15,2013.375,1951.95,1897.35,1835.925,1774.5,1842.75,1911,2033.85,2115.75,2184,2293.2,2375.1,2457,2375.1,2252.25,2129.4,1992.9,1842.75,1706.25,1610.7,1501.5,1638,1777.23,1919.19,2059.785,2199.015,2336.88,2476.11,2613.975,2457,2297.295,2170.35,1904.175,1636.635,1449.63,1365,1236.69,1569.75,1835.925,1979.25,2080.26,2306.85,2457,2605.785,2852.85,2730,2661.75,2530.71,2429.7,2306.85,2182.635,2115.75,2047.5,1911,1774.5,1638,1501.5,1365,1228.5,1092,955.5,1092,1228.5,1365,1501.5,1638,1774.5,1911,2047.5,2197.65,2347.8,2497.95,2648.1,2798.25,2948.4,3091.725,3228.225,3084.9,2934.75,2784.6,2634.45,2484.3,2334.15,2190.825,2047.5,1774.5,1528.8,1228.5,1092,962.325,832.65,709.8,655.2,791.7,941.85,1090.635,1283.1,1528.8,1760.85,1911,2047.5,2293.2,2455.635,2593.5,2852.85,3016.65,3289.65,3412.5,3549,3412.5,3262.35,3001.635,2728.635,2455.635,2319.135,2115.75,2047.5,1794.975,1542.45,1365,1078.35,941.85,764.4,477.75,273,526.4805,682.5,832.65,1090.635,1337.7,1501.5,1909.635,2047.5,2184,2443.35,2661.75,2989.35,3274.635,3562.65,3701.88,3842.475,3685.5,3535.35,3207.75,2934.75,2730,2525.25,2306.85,2047.5,1842.75,1490.58,1119.3,982.8,709.8,436.8,163.8,40.95,204.75,436.8,764.4,1078.35,1300.845,1610.7,1883.7,2047.5,2320.5,2613.975,2907.45,3194.1,3480.75,3753.75,3890.25,4093.635,3549,2866.5,2184,1842.75,1706.25,1228.5,682.5,1.092,271.635,464.1,750.75,955.5,1228.5,1569.75,1842.75,2047.5,2320.5,2593.5,2866.5,3139.5,3412.5,3685.5,3822,3999.45,3822,3685.5,3412.5,3139.5,2866.5,2593.5,2320.5,2047.5,1774.5,1501.5,1228.5,955.5,682.5,518.7,313.95,135.135,313.95,518.7,682.5,955.5,1228.5,1501.5,1774.5,2047.5,2252.25,2525.25,2757.3,3057.6,3344.25,3494.4,3671.85,3842.475,3671.85,3494.4,3344.25,3057.6,2757.3,2525.25,2252.25,2047.5,1774.5,1501.5,1228.5,955.5,819,696.15,477.75,375.375,477.75,696.15,819,955.5,1228.5,1501.5,1774.5,2047.5,2320.5,2593.5,2866.5,3016.65,3139.5,3276,3426.15,3617.25,3426.15,3276,3139.5,3016.65,2866.5,2593.5,2320.5,2047.5,1829.1,1556.1,1255.8,1146.6,955.5,873.6,764.4,679.77,846.3,955.5,1201.2,1351.35,1501.5,1706.25,1842.75,2047.5,2184,2347.8,2525.25,2702.7,2866.5,3003,3139.5,3288.285,3139.5,3003,2866.5,2702.7,2525.25,2347.8,2184,2047.5,1911,1774.5,1638,1501.5,1365,1228.5,1092,1044.225,1092,1228.5,1365,1501.5,1638,1774.5,1911,2046.135,2184,2317.77,2454.27,2525.25,2593.5,2689.05,2852.85,2921.1,2852.85,2689.05,2593.5,2525.25,2454.27,2317.77,2184,2047.5,1965.6,1815.45,1733.55,1636.635,1474.2,1228.5,1351.35,1481.025,1562.925,1644.825,1726.725,1808.625,1890.525,1945.125,2027.025,2047.5,2095.275,2145.78,2197.65,2245.425,2293.2,2327.325,2361.45,2388.75,2361.45,2327.325,2293.2,2245.425,2197.65,2145.78,2095.275,2047.5,2020.2,1992.9,1965.6,1938.3,1911,1883.7,1829.1,1774.5,1842.75,1911,1979.25,2047.5,2115.75,2184,2115.75,2047.5])
freq, spectrum = signal.periodogram(data, window="hann", fs=195312.5, nfft=1024) # using Hann window and padding with 1024 - 400 zeros

# Plotting
plt.figure()
plt.xlabel("Frequency")
plt.ylabel("Magnitude")
plt.title("Spectral density")
plt.plot(freq, spectrum)

## zoom and plot points
plt.figure()
plt.xlim([0, 10000])
plt.xlabel("Frequency")
plt.ylabel("Magnitude")
plt.title("Spectral density")
plt.plot(freq, spectrum)
plt.plot(freq, spectrum, "ko")

plt.show()

I got this:

enter image description here

enter image description here

First thing then: you should use the RealFFT since your input is real-valued. From N samples, you get N//2+1 samples. This is faster and less error prone. You also need to get rid of redundancy everywhere: this is a mess.

As far as your problem is concerned, you forgot to remove the DC component (aka detrend-ing) which is obviously leaking there. Removing the mean value of your data should be enough to get rid of it.

Improve this answer
$\endgroup$
2
  • $\begingroup$ First, I didn't understand: your program seemed correct... But the problem was not a bug in what you've written but in what you forgot to write. That would have been identified in mere seconds if you had just plotted your data. Please do that in the future. $\endgroup$
    – user13706
    Jul 15, 2015 at 18:41
  • $\begingroup$ I used your Python code to test other samples and got acceptable results. But using the RealFFT from LomotFFT also didn't work. I have asolution with the Python code but the true reason for this behaviour is not clear to me. $\endgroup$
    – dDebug
    Jul 16, 2015 at 15:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.