sub fft {
  my $PI = 3.141593;
  my ($F, $LN) = @_;
  my @f = @{$F};

  my $tmp;
  my @t = (0.0, 0.0);
  my @u = (0.0, 0.0);
  my @w = (0.0, 0.0);

  my ($N, $NV2, $NM1, $i, $j, $k, $l, $ip, $le, $le1);
  $N = 2 ** $LN;
  $NV2 = $N / 2;
  $NM1 = $N - 1;
  $j = 1;

  for ($i=1; $i<=$NM1; ++$i) {
    if ($i < $j) {
      $t[0] = ${$f[$j-1]}[0];
      $t[1] = ${$f[$j-1]}[1];
      ${$f[$j-1]}[0] = ${$f[$i-1]}[0];
      ${$f[$j-1]}[1] = ${$f[$i-1]}[1];
      ${$f[$i-1]}[0] = $t[0];
      ${$f[$i-1]}[1] = $t[1];
    }
    $k = $NV2;
    while ($k < $j) {
      $j -= $k;
      $k /= 2;
    }
    $j += $k;
  }

  for ($l=1; $l<=$LN; ++$l) {
    $le = 2 ** $l;
    $le1 = $le / 2;
  
    $u[0] = 1.0;
    $u[1] = 0.0;
  
    $w[0] = cos($PI / $le1);
    $w[1] = - sin($PI / $le1);
print "angle for $le1\n";
    for ($j = 1; $j <= $le1; ++$j) {
      for ($i = $j; $le < 0 ? $i >= $N : $i <= $N; $i += $le) {
        $ip = $i + $le1;
        $t[0] = ${$f[$ip-1]}[0] * $u[0] -
                      ${$f[$ip-1]}[1] * $u[1];
        $t[1] = ${$f[$ip-1]}[0] * $u[1] +
                      ${$f[$ip-1]}[1] * $u[0];
        ${$f[$ip-1]}[0] = ${$f[$i-1]}[0] - $t[0];
        ${$f[$ip-1]}[1] = ${$f[$i-1]}[1] - $t[1];
        ${$f[$i-1]}[0] += $t[0];
        ${$f[$i-1]}[1] += $t[1];
      }
    $tmp = $u[0];
    $u[0] = $u[0] * $w[0] - $u[1] * $w[1];
    $u[1] = $tmp * $w[1] + $u[1] * $w[0];
    }
    $tmp = $u[0];
    $u[0] = $u[0] * $w[0] - $u[1] * $w[1];
    $u[1] = $tmp * $w[1] + $u[1] * $w[0];
  }
  
  for ($i=1; $i<=$N; $i++) {
    ${$f[$i-1]}[0] /= $N;
    ${$f[$i-1]}[1] /= $N;
  }
}


my @F = ([2.0, 0.0], [3.0, 0.0], [4.0, 0.0], [4.0, 0.0]);

foreach $e (@F) {
  print "${$e}[0] + i ${$e}[1]\n";
}
fft(\@F, 2);
foreach $e (@F) {
  print "${$e}[0] + i ${$e}[1]\n";
}