logo

ジュリアで逆フーリエ変換を直接計算する方法 📂Julia

ジュリアで逆フーリエ変換を直接計算する方法

概要

Julia高速フーリエ変換を使うにはFFTW.jlパッケージを使えばよいが、既に実装されたフーリエ逆変換 ifft の代わりに三角関数を直接扱うならフーリエ級数を数式的に理解し係数を扱える必要がある。

コード

遅いフーリエ逆変換

当然、速度は最適化された ifft よりずっと遅いが、我々が理論的に知っているフーリエ変換方式そのまま計算される分だけ直感的で介入の余地が多い。

using DataFrames, FFTW, Plots

data = CSV.read("https://freshrimpsushi.github.io/ko/posts/2765/example.csv", DataFrame).x

Fs = 8
T = 120/8
fften = fft(data[1:round(Int64, Fs*T):end])

L = length(fften)
t_ = 1:L
a_ = 2real(fften[2:div(L,2)]) / L
b_ = -2imag(fften[2:div(L,2)]) / L

y = fill(real(fften[1]) / L, round(Int64, L))
for k in eachindex(a_)
    _sincos = sincospi.((2k/L)*t_)
    _sin = first.(_sincos)
    _cos = last.(_sincos)
    y .+= (a_[k] .* _cos) + (b_[k] .* _sin)
end
plot(data)
plot!(t_*(Fs*T), y, lw = 2)

サンプリング間隔

次はサンプリング間隔 $T$ を変えながらフーリエ逆変換を求めた例だ。

$T = 120/8$

120

$T = 12/8$

12

$T = 1/8$

1