#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <complex.h>
#include "matrixlib.h"

void fft(int N,int s,complex x[N],complex b[N]){
	if(N==1) {
		b[0]=x[0];
		return;
	}
	if(N%2){
		printf("N not divisible by 2!\n");
		exit(1);
	}
	int K=N/2;
	fft(K,2*s,&x[0],&b[0]);
	fft(K,2*s,&x[s],&b[K]);
	for(int k=0;k<K;k++){
		complex even=b[k],odd=b[k+K];
		complex w=cexp(-I*2*M_PI*k/N);
		b[k]=even+w*odd;
		b[k+K]=even-w*odd;
	}
}

#define FFTSIZE 16
int main(){
	complex X[FFTSIZE],B[FFTSIZE];
	for(int i=0;i<FFTSIZE;i++) X[i]=i;
	fft(FFTSIZE,1,X,B);
	cvecprint(FFTSIZE,B);
	return 0;
}