Progressing with my attempt to record and process ultrasonic sounds, I have spent this evening with a trial version of a fast fourier transform (FFT) on the KL25z ARM board. My hope is that if I can get it running fast enough and well enough, I can offload that part of the work to this "co-processor", leaving the Raspberry Pi to deal with all the rest of the UI, communications, pattern-matching and so on.

There are several FFT implementations in the mbed code repository. Unfortunately, all of them bar one failed to compile on my system. I can only assume that this is because they were originally written for the early mbed boards, and have not been brought into line with the new KL25z board.

Luckily, one version did compile. Unlike most of the others this was not a stand-alone application, but a library, so I had to write my own "wrapper" around it. After a few hours of fiddling, I am now at a position where I can at least see it basically working.

#include "mbed.h"
#include "math.h"
#include "stdio.h"
#include "FFT.h"

#define POINTS 256
#define BKSP 8
#define SCALE 4

const double PI = 3.14159265358979323846;
const double CIRCLE = 2 * PI;

float data[POINTS];

float height(int i) {
    return abs(data[i]);
}

void run(float freq) {
    float step = freq * CIRCLE / POINTS;
    for (int i = 0; i < POINTS; ++i) {
        data[i] = (sin(i * step) * SCALE * log(freq)) + SCALE;
    }

    Timer t;
    t.start();
    vRealFFT(data, POINTS);
    t.stop();
    printf("FFT with %d points took %f seconds\r\n", POINTS, t.read());

    float max = 0;
    for (int i = 0; i < POINTS/2; ++i) {
        if (height(i) > max) max = data[i];
    }

    for (float level = max; level > 0; level -= max/10) {
        for (int i = 0; i < POINTS/2; i += 2) {
            putchar( (height(i) > level) ? '|' : ' ');
        }
        printf("\r\n");
    }
    for (int i = 0; i < POINTS/4; ++i) {
        putchar('-');
    }
    printf("\r\n");
}

int main() {
    for (;;) {
        printf("Enter frequency:\r\n");
        int number = 0;
        for (;;) {
            int c = getchar();
            putchar(c);
            if (c == BKSP) {
                number /= 10;
                continue;
            }
            if (c < '0' || c > '9') {
                printf("\r\n");
                break;
            }
            number *= 10;
            number += c - '0';
        }
       
        if (number > 0) run(1.0 * number);
    }
}

The code is pretty simple at the moment. It asks the user for a "frequency", then generates data representing a sine wave which fits that specified number of cycles in the array. This array is then passed to the FFT routine, which replaces the data with the frequency-domain values.

By experimentation I have found that only the odd-numbered frequency levels seem to make any sense, so I ignore all the even ones. I have also found that the frequency levels seem to drop off the higher the input frequency, so this version of the code applies a logarithmic "fiddle-factor" which sort of corrects it. There's still work to do on that, as I don't yet understand why it happens.

The code times how long the FFT took, then dumps out the frequency values in a vaguely pretty "ascii art" form:

It's clear from the picture that the FFT is essentially doing its job. Changing the input frequency moves the line on the chart.

This version is not as fast as I had hoped it might be, though. Even for only 256 samples per transform it can't even do 50 transforms per second. It's a shame that I could not get any of the other versions working, as the ones which didn't compile might potentially have been faster or more efficient. One was coded in assembler, which ought to be able to squeeze maximum performance out of the CPU. I might look deeper into why they won't compile another time.

This test is very simplistic. I need to check that more complex signals than a simple sine wave don't confuse it too much, and I also need to experiment with actually reading from an ADC pin to prepare the data. The documentation suggests that this can be done using DMA. With any luck I can get the board to carry on recording the next chunk of data while it is transforming the previous chunk.