# -*- coding: utf-8 -*-
"""
Created on Wed Oct 26 01:09:48 2016

@author:
"""

import numpy as np
#from time import clock

def randomnumbers(N): #creates an array of N random complex numbers
    a = np.zeros(N)
    b = np.zeros(N)
    for i in range(0,N):
        a[i] = 2 * np.random.random() - 1
        b[i] = 2 * np.random.random() - 1
    return a + 1j * b

def ft(f, multi): #fourier transformation - optional multiplier for periodic extension
    N = f.size
    ft  = np.zeros(multi*N, dtype=complex)
#    ft_ = np.zeros(N,       dtype=complex) #only for slower version
    for n in range(0, multi*N):
        ft_ = np.exp((-2j*np.pi*(n - N//2)*np.arange(-N//2, N//2))/N) * f
#        Slower version of line above:
#        for i in range(0, N):
#            ft_[i] = np.exp((-2j*np.pi*(n - N//2)*(i - N//2))/N) * f[i]
        ft[n] = np.sum(ft_)/N
    return ft

def invft(ft, multi): #inverse fourier transformation - optional multiplier for periodic extension
    N = ft.size
    f  = np.zeros(multi*N, dtype=complex)
#    f_ = np.zeros(N,       dtype=complex) #only for slower version
    for i in range(0, multi*N):
        f_ = np.exp((2j*np.pi*(i - N//2)*np.arange(-N//2, N//2))/N) * ft
#        Slower version of line above:
#        for n in range(0, N):
#            f_[n] = np.exp((2j*np.pi*(n - N//2)*(i - N//2))/N) * ft[n]
        f[i] = np.sum(f_)
    return f

#f = randomnumbers(5000)
#time = clock()
#invft(ft(f,1),1)
#print(clock()-time)

randomnumbers(10).tofile("numbers")

f = np.fromfile("numbers", dtype=complex)
N = f.size
print("Original", N, "numbers:\n", f)

ft = ft(f, 1)
ft.tofile("ftnumbers")
print("\nFourier transformed numbers:\n", ft)
ft_ = np.fft.fft(f)
print("\nComparison with FFT:\n", ft_)

f2 = invft(ft, 1)
print("\nInverse transformed numbers:\n", f2)
f2_ = np.fft.ifft(ft_)
print("\nComparison with IFFT:\n", f2_)