この例では、長さ 8 までのバイト・スライスしか完全に復元できません。
FFT を理解しょうと思うも、手を動かしてみないとガマンできない検証猿のため、手を出したものの、何かを掴んだ気がするも引き出せない、何がなんだか分からずパニくったのでクールダウン中の備忘録として。勉強して適宜アップデートしていきます。
Go言語(golang)で FFT してみたかった
高速フーリエ変換(FFT, Fast Fourier Transform)と逆高速フーリエ変換(IFFT, Inverse Fast Fourier Transform)のシンプルなサンプルです。
とりあえず、バイトのスライス(可変長配列)をフーリエ変換して、それを元に逆フーリエ変換するだけです。
ポイント:
「複雑な値の変化を、複数の『シンプルな波』の集合体とみなした場合、FFT は各々の波に分解する(周波数成分ごとに分ける)処理である」と解釈した場合、 sin
や cos
の円形運動(三角関数)だけではダメで、いったん複素数(complex number)への変換が必要となる。
package main
import (
"fmt"
"math"
"math/cmplx"
)
func main() {
byteData := []byte{1, 2, 3, 4, 5, 6, 7, 8}
initialStep := 1
// FFT(入力データを周波数成分に変換)
complexData := bytesToComplex128(byteData)
transformed := fft(complexData, len(complexData), initialStep)
fmt.Println("FFT結果(実部=cos成分、虚部=sin成分、強度=振幅):")
for i, v := range transformed {
if i != 0 {
fmt.Printf(" #%d: 周波数成分(実部: %f, 虚部: % 2.6f), 強度: %f\n",
i+1, real(v), imag(v), cmplx.Abs(v))
continue
}
fmt.Printf(" #%d: 直流成分(入力データの平均値)(実部: %f, 虚部: % 2.6f), 強度: %f\n",
i+1, real(v), imag(v), cmplx.Abs(v))
}
// IFFT(周波数成分を元のデータに戻す)
restoredComplex := ifft(transformed, len(transformed), initialStep)
restoredBytes := complex128ToBytes(restoredComplex)
fmt.Println("Original bytes:", byteData)
fmt.Println("Restored bytes:", restoredBytes)
//
// Output:
// FFT結果(実部=cos成分、虚部=sin成分、強度=振幅):
// #1: 直流成分(入力データの平均値)(実部: 36.000000, 虚部: 0.000000), 強度: 36.000000
// #2: 周波数成分(実部: -4.000000, 虚部: 9.656854), 強度: 10.452504
// #3: 周波数成分(実部: -4.000000, 虚部: 4.000000), 強度: 5.656854
// #4: 周波数成分(実部: -4.000000, 虚部: 1.656854), 強度: 4.329569
// #5: 周波数成分(実部: -4.000000, 虚部: 0.000000), 強度: 4.000000
// #6: 周波数成分(実部: -4.000000, 虚部: -1.656854), 強度: 4.329569
// #7: 周波数成分(実部: -4.000000, 虚部: -4.000000), 強度: 5.656854
// #8: 周波数成分(実部: -4.000000, 虚部: -9.656854), 強度: 10.452504
// Original bytes: [1 2 3 4 5 6 7 8]
// Restored bytes: [1 2 3 4 5 6 7 8]
}
// bytesToComplex128 converts a byte slice to a complex128 slice.
func bytesToComplex128(data []byte) []complex128 {
complexData := make([]complex128, len(data))
for i, v := range data {
complexData[i] = complex(float64(v), 0)
}
return complexData
}
// complex128ToBytes converts a complex128 slice to a byte slice.
func complex128ToBytes(data []complex128) []byte {
byteData := make([]byte, len(data))
for i, v := range data {
byteData[i] = byte(math.Round(real(v)))
}
return byteData
}
// fft returns a slice of complex128 values that represent the FFT of the input
// data.
func fft(data []complex128, size, step int) []complex128 {
if size == 1 {
return []complex128{data[0]}
}
angle := -2 * math.Pi / float64(size) // θ
exp := cmplx.Rect(1, angle) // 回転因子(単位円上の角度から計算)
wn := complex(1, 0) // 回転因子の初期化
a0 := fft(data, size/2, 2*step) // 偶数番目の要素を処理
a1 := fft(data[step:], size/2, 2*step) // 奇数番目の要素を処理
y := make([]complex128, size)
for k := 0; k < size/2; k++ {
y[k] = a0[k] + wn*a1[k]
y[k+size/2] = a0[k] - wn*a1[k]
wn *= exp // 回転因子の更新
}
return y
}
// ifft returns a slice of complex128 values that represent the IFFT of the input
// data.
func ifft(data []complex128, size, step int) []complex128 {
if size == 1 {
return []complex128{data[0]}
}
angle := 2 * math.Pi / float64(size) // θ
exp := cmplx.Rect(1, angle) // 回転因子(単位円上の角度から計算)
wn := complex(1, 0) // 回転因子の初期化
a0 := ifft(data, size/2, 2*step) // 偶数番目の要素を処理
a1 := ifft(data[step:], size/2, 2*step) // 奇数番目の要素を処理
y := make([]complex128, size)
for k := 0; k < size/2; k++ {
y[k] = (a0[k] + wn*a1[k]) / 2
y[k+size/2] = (a0[k] - wn*a1[k]) / 2
wn *= exp // 回転因子の更新
}
return y
}
- オンラインで動作を見る @ GoPlayground
みんな、よくこんなん計算式だけで理解できんな。すげぇよ。動いたとはいえ、意味わかんねぇよ。