#include <math.h>
#include <complex.h>
#include <stdint.h>
#include <stdlib.h>

complex float* fft(complex float* p, uint32_t N) {
    complex float w = cexpf(-2 * I * M_PI / N);

    if (N == 1)
        return p;

    complex float* pe = (complex float*)malloc(sizeof(complex float) * (N / 2));
    complex float* po = (complex float*)malloc(sizeof(complex float) * (N / 2));
    for (int i = 0; i < N - 1; i += 2) {
       pe[i / 2] = p[i];
       po[i / 2] = p[i + 1];
    }

    complex float* xe = fft(pe, N / 2);
    complex float* xo = fft(po, N / 2);

    complex float* ps = (complex float*)malloc(sizeof(complex float) * N);

    for (int i = 0; i < N / 2; i++) {
        ps[i] = xe[i] + cpowf(w, (complex float)i) * xo[i];
        ps[i + N / 2] = xe[i] - cpowf(w, (complex float)i) * xo[i];
    }
    
    return ps;
}


complex float* ifft(complex float* p, uint32_t N) {
    complex float w = ((complex float)1 / N) * cexpf(2 * I * M_PI / N);

    if (N == 1)
        return p;

    complex float* pe = (complex float*)malloc(sizeof(complex float) * (N / 2));
    complex float* po = (complex float*)malloc(sizeof(complex float) * (N / 2));
    for (int i = 0; i < N - 1; i += 2) {
       pe[i / 2] = p[i];
       po[i / 2] = p[i + 1];
    }

    complex float* xe = fft(pe, N / 2);
    complex float* xo = fft(po, N / 2);

    complex float* ps = (complex float*)malloc(sizeof(complex float) * N);

    for (int i = 0; i < N / 2; i++) {
        ps[i] = xe[i] + cpowf(w, (complex float)i) * xo[i];
        ps[i + N / 2] = xe[i] - cpowf(w, (complex float)i) * xo[i];
    }
    
    return ps;
}