JTransforms FFT in Android from PCM data

Question

I've been playing with this now for sometime, I cant work out what I am meant to be doing here.

I am reading in PCM audio data into an audioData array:

 recorder.read(audioData,0,bufferSize);     //read the PCM audio data into the audioData array

I want to use Piotr Wendykier's JTransform library in order to preform an FFT on my PCM data in order to obtain the frequency.

import edu.emory.mathcs.jtransforms.fft.DoubleFFT_1D;

At the moment I have this:

       DoubleFFT_1D fft = new DoubleFFT_1D(1024); // 1024 is size of array

for (int i = 0; i < 1023; i++) {
           a[i]= audioData[i];               
           if (audioData[i] != 0)
           Log.v(TAG, "audiodata=" + audioData[i] + " fft= " + a[i]);
       }
       fft.complexForward(a);

I cant make sense of how to work this, can somebody give me some pointers? Will i have to perform any calculations after this?

I'm sure I'm way off, anything would be greatly appreciated!

Ben

Shawnee

Sorry I should have said what complexForward does: Computes 1D forward DFT of complex data leaving the result in a. Complex number is stored as two double values in sequence: the real and imaginary part, i.e. the size of the input array must be greater or equal 2*n. The physical layout of the input data has to be as follows: a[2*k] = Re[k], a[2*k+1] = Im[k], 0<=k<n – Shawnee

give me any example... – Impose

Answer 1

If you're just looking for the frequency of a single sinusoidal tone in the input waveform then you need to find the FFT peak with the largest magnitude, where:

Magnitude = sqrt(re*re + im*im)

The index i of this largest magnitude peak will tell you the approximate frequency of your sinusoid:

Frequency = Fs * i / N

where:

Fs = sample rate (Hz)
i = index of peak
N = number of points in FFT (1024 in this case)

Formularize

That sounds very much like what I am trying to do. However I am unsure how to get the real and imaginary values out of the transform. Also I'm not sure as to whether I should be using a DoubleFFT_1D class? as opposed to the DoubleFFT_2D class, or the FloatFFT_1D class etc – Shawnee

You only need a 1D FFT for audio (2D FFTs are usually for 2D data such as images etc). From what you said in your question above it sounds like the real and imag components are interleaved, so bin i has real component at index 2*i and imag component at index 2*i + 1. – Formularize

Thanks Paul, thats great, a couple of things I dont understand... do I pass the audioData array like this fft.complexForward(audioData);? and if I do, does it manipulate the array? it has no return type. Thanks – Shawnee

I'm not familiar with this specific FFT but yes it sounds like it's what is known as an "in place" implementation, so you pass time domain data in and it gets overwritten with the frequency domain result. For the input you set the real parts (2*i) to your sample data values and the imaginary parts (2*i+1) to 0. – Formularize

So should i double each index in the audio data array before passing it to the FFT? ie (2*i)? – Shawnee

Yes, when copying your audio data to the FFT then you would do e.g. fft[2*i] = pcm[i]; fft[2*i+1] = 0; – Formularize

Answer 2

Since I've spent some hours on getting this to work here's a complete implementation in Java:

import org.jtransforms.fft.DoubleFFT_1D;

public class FrequencyScanner {
    private double[] window;

    public FrequencyScanner() {
        window = null;
    }

    /** extract the dominant frequency from 16bit PCM data.
     * @param sampleData an array containing the raw 16bit PCM data.
     * @param sampleRate the sample rate (in HZ) of sampleData
     * @return an approximation of the dominant frequency in sampleData
     */
    public double extractFrequency(short[] sampleData, int sampleRate) {
        /* sampleData + zero padding */
        DoubleFFT_1D fft = new DoubleFFT_1D(sampleData.length + 24 * sampleData.length);
        double[] a = new double[(sampleData.length + 24 * sampleData.length) * 2];

        System.arraycopy(applyWindow(sampleData), 0, a, 0, sampleData.length);
        fft.realForward(a);

        /* find the peak magnitude and it's index */
        double maxMag = Double.NEGATIVE_INFINITY;
        int maxInd = -1;

        for(int i = 0; i < a.length / 2; ++i) {
            double re  = a[2*i];
            double im  = a[2*i+1];
            double mag = Math.sqrt(re * re + im * im);

            if(mag > maxMag) {
                maxMag = mag;
                maxInd = i;
            }
        }

        /* calculate the frequency */
        return (double)sampleRate * maxInd / (a.length / 2);
    }

    /** build a Hamming window filter for samples of a given size
     * See http://www.labbookpages.co.uk/audio/firWindowing.html#windows
     * @param size the sample size for which the filter will be created
     */
    private void buildHammWindow(int size) {
        if(window != null && window.length == size) {
            return;
        }
        window = new double[size];
        for(int i = 0; i < size; ++i) {
            window[i] = .54 - .46 * Math.cos(2 * Math.PI * i / (size - 1.0));
        }
    }

    /** apply a Hamming window filter to raw input data
     * @param input an array containing unfiltered input data
     * @return a double array containing the filtered data
     */
    private double[] applyWindow(short[] input) {
        double[] res = new double[input.length];

        buildHammWindow(input.length);
        for(int i = 0; i < input.length; ++i) {
            res[i] = (double)input[i] * window[i];
        }
        return res;
    }
}

FrequencyScanner will return an approximation of the dominant frequency in the presented sample data. It applies a Hamming window to it's input to allow passing in arbitrary samples from an audio stream. Precision is achieved by internally zero padding the sample data before doing the FFT transform. (I know there are better - and far more complex - ways to do this but the padding approach is sufficient for my personal needs).

I testet it against raw 16bit PCM samples created from reference sounds for 220hz and 440hz and the results match.

Stockholm

I think you should call realforwardfull to calculate the fft results. Realforward has a different interpretation of the results. if you check the doc of the api. – Stevenstevena

Answer 3

Yes you need to use realForward function instead of complexForward, because you pass it a real array and not a complex array from doc.

EDIT:

Or you can get the real part and perform complex to complex fft like this :

double[] in = new double[N];
read ...
double[] fft = new double[N * 2];

for(int i = 0; i < ffsize; ++i)
{
  fft[2*i] = mic[i];
  fft[2*i+1] = 0.0;
}
fft1d.complexForward(fft);

I try and I compare results with matlab, and I don't get same results... (magnitude)

Erme

Answer 4

If you're looking for the FFT of an Audio input ( 1D, real data ), should'nt you be using the 1D REAL Fft?

Magic

its something i will experiment with – Shawnee