FFT

#define _CRT_SECURE_NO_WARNINGS
#define _USE_MATH_DEFINES
#include <math.h>
#include <stdio.h>
#define N 4096
#define FFT 1
#define IFFT -1
void data_in(int, float *, float*);
void fft_sub(int, int, float *, float *);
float f[N];
void main(void)
{
  int i, n, flag;
  float fr[N], fi[N];
  time_t t1, t2;
  printf("データ個数を入力 n=");
  scanf("%d", &n);
  printf("\n");
  data_in(n, fr, fi);
  time(&t1);
  flag = FFT;
  fft_sub(n, flag, fr, fi);
  time(&t2);
  printf("time = %f\n", difftime(t2, t1));
  printf("データ個数 n=%d \n", n);
  printf("k   f   fr   fi\n");
  for (i = 0; i < n; i++)
  {
    printf("%3d %8.8f %8.8f %8.8f\n", i, f[i], fr[i], fi[i]);
  }
  printf("%d\n", calCount);
}
void data_in(int n, float fr[], float fi[])
{
  //方形波
  int i;
  /*  
  for (i = 0; i < n; i++)
  {
    fr[i] = 0.0;
    fi[i] = 0.0;
    if (i < n / 4)
      fr[i] = 1.0;
    if (i == n / 4)
      fr[i] = 0.5;
    if (i == 3 * n / 4)
      fr[i] = 0.5;
    if (i > 3 * n / 4)
      fr[i] = 1.0;
    if (i > n / 4 && i < 3 * n / 4)
      fr[i] = 0.0;
    f[i] = fr[i];
  }
  */  
  //4_2の波形
}
void fft_sub(int n, int flag, float *fr, float *fi)
{
  int i, m, iw, j, k, k0, k2, l, lp, lp2, n2, kk;
  double arg, darg, c, s, wr, wi, xa, ya;
  //Bit reversal
  j = 0;
  for (i = 0; i < n - 1; i++)
  {
    if (i < j)
    {
      xa = fr[i];
      fr[i] = fr[j];
      fr[j] = xa;
      ya = fi[i];
      fi[i] = fi[j];
      fi[j] = ya;
    }
    n2 = n / 2;
    while (j >= n2)
    {
      j = j - n2;
      n2 = n2 / 2;
    }
    j += n2;
  }
  //fft
  m = 0;
  n2 = n;
  while (n2 != 1)
  {
    m++;
    n2 = n2 / 2;
  }
  for (l = 1; l <= m; l++)//段
  {
    lp = pow(2.0, (double)l);
    lp2 = lp / 2;
    darg = -(double)flag*M_PI / (double)lp2;
    arg = 0;
    for (j = 0; j < lp2; j++)
    {
      c = cos(arg);
      s = sin(arg);
      arg += darg;
      for (i = j; i < n; i += lp)
      {
        iw = i + lp2;
        /*バタフライ演算部分を自作する*/
      }
    }
  }
  if (flag == FFT)
  {
    for (i = 0; i < n; i++)
    {
      fr[i] /= (double)n;
      fi[i] /= (double)n;
    }
  }
}

ソースコード

ソースコードファイル

課題

パソコン実践5.1を演習4.2の問題の波形に対して行った結果を報告せよ.DFTとFFTの値を比較すること.提出レポートには,ソースコード,結果を示すこと.

提出方法

締切日時:6月19日授業時
提出物:課題の問題とプログラムコード,結果を含めたレポートを印刷したもの.