Ruby - （離散）フーリエ変換！

前回、「離散フーリエ変換」の C++ での実装に関する記事を紹介しました。

• C++ - （離散）フーリエ変換

今回は、同じアルゴリズムを Ruby で実装してみました。
実際、ほとんど同じです。

0. 前提条件

• Linux Mint 14 Nadia (64bit) での作業を想定。
• Ruby 2.0.0-p195

1. （離散）フーリエ変換について

前回の記事を参照ください。

• C++ - （離散）フーリエ変換

2. Ruby スクリプト作成

• 変換元の周期関数は、 \(f(t) = 2 * \sin(4t) + 3 \times \cos(2t)\) とした。
• t の範囲は \(0 \sim 2\pi\) に限定している。
• 分割数は N の値を変更して対応する。（以下の例では N=100 としている）
• 「元データ作成」、「元データを離散フーリエ変換」、「離散フーリエ変換されたデータを逆離散フーリエ変換」としている。

File: discrete_fourier_transform.cpp

#! /usr/local/bin/ruby
#*********************************************
# 離散フーリエ変換
#   f(t) = 2 * sin(4 * t) + 3 * cos(2 * t)
#          ( 0 <= t < 2 * pi )
#*********************************************/
#
# 計算クラス
class DiscreteFourierTransform
  N = 100                # 分割数
  CSV_DFT  = "DFT.csv"   # 出力ファイル (DFT)
  CSV_IDFT = "IDFT.csv"  # 出力ファイル (IDFT)

  def initialize
    @src_re, @src_im = Array.new, Array.new
    @dft_re, @dft_im = Array.new, Array.new
  end

  # 元データ作成
  def make_source_data
    N.times do |i|
      @src_re << 2 * Math.sin(4 * (2 * Math::PI / N) * i) \
               + 3 * Math.cos(2 * (2 * Math::PI / N) * i)
      @src_im << 0.0
    end
  rescue => e
    raise
  end

  # 離散フーリエ変換
  def execute_DFT
    # 出力ファイルOPEN
    open(CSV_DFT, "w") do |f|
      # ヘッダ出力 ( k, 角周波数, 元データ(実部), 元データ(虚部), DFT(実部), DFT(虚部) )
      f.puts "k,f,x_re,x_im,X_re,X_im"
      # 計算・結果出力
      N.times do |k|
        dft_re = 0.0
        dft_im = 0.0
        N.times do |n|
          dft_re += @src_re[n] * ( Math.cos((2 * Math::PI / N) * k * n))
                  + @src_im[n] * ( Math.sin((2 * Math::PI / N) * k * n))
          dft_im += @src_re[n] * (-Math.sin((2 * Math::PI / N) * k * n))
                  + @src_im[n] * ( Math.cos((2 * Math::PI / N) * k * n))
        end
        @dft_re << dft_re
        @dft_im << dft_im
        f.puts "%d,%.6f,%.6f,%.6f,%.6f,%.6f" %
          [k, (2 * Math::PI / N) * k, @src_re[k], @src_im[k], dft_re, dft_im]
      end
    end
  rescue => e
    raise
  end

  # 逆離散フーリエ変換
  def execute_IDFT
    # 出力ファイルOPEN
    open(CSV_IDFT, "w") do |f|
      # ヘッダ出力 ( k, 角周波数, DFT(実部), DFT(虚部), IDFT(実部), IDFT(虚部) )
      f.puts "k,f,X_re,X_im,x_re,x_im"
      # 計算・結果出力
      N.times do |n|
        idft_re, idft_im  = 0.0, 0.0
        N.times do |k|
          idft_re += @dft_re[k] * (Math.cos((2 * Math::PI / N) * k * n)) \
                   - @dft_im[k] * (Math.sin((2 * Math::PI / N) * k * n))
          idft_im += @dft_re[k] * (Math.sin((2 * Math::PI / N) * k * n)) \
                   + @dft_im[k] * (Math.cos((2 * Math::PI / N) * k * n))
        end
        idft_re /= N
        idft_im /= N
        f.puts "%d,%.6f,%.6f,%.6f,%.6f,%.6f" %
          [n, (2 * Math::PI / N) * n, @dft_re[n], @dft_im[n], idft_re, idft_im]
      end
    end
  rescue => e
    raise
  end
end

if __FILE__ == $0
  begin
    obj = DiscreteFourierTransform.new  # 計算クラスインスタンス化
    obj.make_source_data                # 元データ作成
    obj.execute_DFT                     # 離散フーリエ変換
    obj.execute_IDFT                    # 逆離散フーリエ変換
  rescue => e
    $stderr.puts "[#{e.class}] #{e.message}"
    e.backtrace.each{ |tr| $stderr.puts "\t#{tr}" }
  end
end

• Gist - Ruby script to compute discrete Fourier transform.

離散フーリエ変換と逆離散フーリエ変換の処理は、実質符号を変えるだけなので、１つにまとめても良いかもしれない。

3. 実行

まず、実行権限を付与。

$ chmod +x discrete_fourier_transform.rb

そして、実行。

$ ./discrete_fourier_transform

コンソールには特に何も表示しない。
アプリと同じディレクトリに “DFT.csv”, “IDFT.csv” という CSV ファイルが作成される。
内容は以下のようになっているはず。

File: DFT.csv

k,f,x_re,x_im,X_re,X_im
0,0.000000,3.000000,0.000000,0.000000,0.000000
1,0.062832,3.473724,0.000000,0.000000,-0.000000
2,0.125664,3.869257,0.000000,150.000000,-0.000000
3,0.188496,4.158424,0.000000,0.000000,-0.000000
4,0.251327,4.317576,0.000000,0.000000,-100.000000
5,0.314159,4.329164,0.000000,-0.000000,0.000000
6,0.376991,4.182959,0.000000,0.000000,0.000000
7,0.439823,3.876846,0.000000,0.000000,0.000000
8,0.502655,3.417134,0.000000,0.000000,0.000000
9,0.565487,2.818364,0.000000,0.000000,-0.000000
         :
====< 途中省略 >====
         :
90,5.654867,-0.248520,0.000000,-0.000000,-0.000000
91,5.717699,-0.263689,0.000000,-0.000000,0.000000
92,5.780530,-0.202174,0.000000,-0.000000,-0.000000
93,5.843362,-0.052303,0.000000,-0.000000,0.000000
94,5.906194,0.190852,0.000000,0.000000,0.000000
95,5.969026,0.524938,0.000000,-0.000000,-0.000000
96,6.031858,0.940264,0.000000,0.000000,100.000000
97,6.094690,1.420235,0.000000,-0.000000,-0.000000
98,6.157522,1.942242,0.000000,150.000000,-0.000000
99,6.220353,2.478964,0.000000,0.000000,-0.000000

File: IDFT.csv

k,f,X_re,X_im,x_re,x_im
0,0.000000,0.000000,0.000000,3.000000,-0.000000
1,0.062832,0.000000,-0.000000,3.473724,-0.000000
2,0.125664,150.000000,-0.000000,3.869257,-0.000000
3,0.188496,0.000000,-0.000000,4.158424,-0.000000
4,0.251327,0.000000,-100.000000,4.317576,-0.000000
5,0.314159,-0.000000,0.000000,4.329164,-0.000000
6,0.376991,0.000000,0.000000,4.182959,-0.000000
7,0.439823,0.000000,0.000000,3.876846,-0.000000
8,0.502655,0.000000,0.000000,3.417134,-0.000000
9,0.565487,0.000000,-0.000000,2.818364,-0.000000
         :
====< 途中省略 >====
         :
90,5.654867,-0.000000,-0.000000,-0.248520,0.000000
91,5.717699,-0.000000,0.000000,-0.263689,-0.000000
92,5.780530,-0.000000,-0.000000,-0.202174,-0.000000
93,5.843362,-0.000000,0.000000,-0.052303,0.000000
94,5.906194,0.000000,0.000000,0.190852,-0.000000
95,5.969026,-0.000000,-0.000000,0.524938,-0.000000
96,6.031858,0.000000,100.000000,0.940264,0.000000
97,6.094690,-0.000000,-0.000000,1.420235,0.000000
98,6.157522,150.000000,-0.000000,1.942242,0.000000
99,6.220353,0.000000,-0.000000,2.478964,-0.000000

離散フーリエ変換後のデータを逆離散フーリエ変換して、元のデータに戻っていることが分かる。
また、C++ 版と同じ結果になっていることも分かる。（結果が同じなのでグラフ出力は省略）

---

電気工学、音響学、振動、光学等でよく使用する重要な概念です。応用範囲は広いので他にも利用できるかと思います。

以上。