import cmath
import math

def bit_reversal_permutation(x):
    """
    配列xの要素をビット反転順に並べ替える関数
    例: N=8 なら [0, 1, 2, 3, 4, 5, 6, 7] -> [0, 4, 2, 6, 1, 5, 3, 7]
    """
    N = len(x)
    # Nが2の何乗かを求める（ビット長）
    num_bits = N.bit_length() - 1 
    
    result = [0] * N
    for i in range(N):
        # i (例: 3 = 011) をビット反転させる (110 = 6)
        # '{:0{width}b}' で2進数文字列にし、[::-1] で反転し、int(..., 2) で数値に戻す
        reversed_i = int('{:0{width}b}'.format(i, width=num_bits)[::-1], 2)
        result[reversed_i] = x[i]
        
    return result

def fft(x):
    """
    非再帰型 Cooley-Tukey FFT (Radix-2)
    入力 x の長さ N は 2の累乗でなければならない。
    """
    N = len(x)
    if N <= 1:
        return x
        
    # Nが2の累乗かどうかのチェック (ビット演算を利用した簡単な判定)
    if (N & (N - 1)) != 0:
        raise ValueError("データ長 N は 2の累乗である必要があります。")

    # 1. 限界まで分解した状態（ビット反転順）からスタート
    X = bit_reversal_permutation(x)

    # 2. 3重ループによるバタフライ演算（合体プロセス）
    
    # 外側のループ: ブロックサイズ L を 2, 4, 8, ... N と倍増させていく（トーナメントの階層）
    L = 2
    while L <= N:
        # この階層における「基本となる回転因子」
        # 1ブロックを 1周 (2pi) に対応させるので、1歩あたりの角度は -2pi / L
        w_L = cmath.exp(-2j * math.pi / L)
        
        # 内側のループ(1): 配列全体をブロックサイズ L ごとにジャンプして進む
        for start in range(0, N, L):
            # ブロック内で使う回転因子の初期値 (W^0 = 1)
            w = 1
            
            # 内側のループ(2): ブロックの前半と後半をペアにしてバタフライ演算（回転因子の更新）
            # ペアの距離 (half_L) は常にブロックサイズの半分
            half_L = L // 2
            
            for k in range(half_L):
                # 上のハコ A のインデックス (これが前半の k 番目)
                index_A = start + k
                # 下のハコ B のインデックス (これが後半の k 番目、上のハコから half_L だけ離れた位置)
                index_B = start + k + half_L
                
                # 現在のハコの中身を取り出す
                A = X[index_A]
                B = X[index_B]
                
                # バタフライ演算の適用
                # ルール1: 回転させるのは「下のハコ (後半の k 番目の奇数 DFT)」だけ
                W_B = w * B
                
                # ルール2: クロスして足し引きする
                # 前半 k 番目の DFT = k 番目の偶数 DFT (A) + W^k * k 番目の奇数 DFT (W_B)
                X[index_A] = A + W_B  
                # 後半 k 番目の DFT = k 番目の偶数 DFT (A) - W^k * k 番目の奇数 DFT (W_B)
                X[index_B] = A - W_B  
                
                # 次のペアのために、回転因子を1歩進める (w = w * w_L)
                w = w * w_L
                
        # 階層が1つ終わったら、ブロックサイズを倍にする
        L *= 2

    return X

# N=8 の簡単なデータ
data = [3.2,3.5,4.5,3.0,8.0,5.6,7.8,4.0]

print("Input:")
print([round(val, 2) for val in data])

result = fft(data)

print("\nFFT Output:")
for i, val in enumerate(result):
    # 複素数を見やすく丸めて表示
    real = round(val.real, 2)
    imag = round(val.imag, 2)
    print(f"X[{i}] = {real} + {imag}j")