# フーリエ変換 ## 高速フーリエ変換 高速フーリエ変換(Fast Fourier Transformation)を利用することができます。\ 例えば、データとしてmydataがあるとき、 ``` fftresult = np.fft.fft(mydata) ``` で、mydataの離散フーリエ変換結果が手に入ります。\ fftの計算は ![](https://1559620454-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-MEGDe0jc4Gum8QxA8RN%2F-MehADo0wjwMtYDpKaTa%2F-MehAlR3QYcJYmPI1vfx%2Ffft.svg?alt=media\&token=559af1fc-3afd-4a50-856a-51c37a01b1f4) なので（参照： arrayです。プロットする際には、np.real(), np.imag(), np.abs(), np.angle()などを使って、実部, 虚部, 絶対値, 位相などを取り出すと良いでしょう。 結果を解釈するには（例えばスペクトルをプロットしたいとき）各データがどのような周波数に対応するのかという情報が必要です。これは、mydataがそもそも何秒間隔のデータなのかということに依存します。もし0.01秒間隔であれば ``` freqarray = np.fft.fftfreq(xdata.size, d=0.01) ``` とすることで、freqarrayという名前の、周波数軸に対応するarrayを作ることができます。 ## Power Spectral Density ここではPower Spectral Densityの詳しい説明は省きますが、大雑把には、「フーリエ変換の絶対値の自乗を適切に規格化したもの」といえます。「適切に」の意味合いは以下のようになります。対象とするデータを $$g(t)$$ とするとき $$(g(t))^2$$ のことを、慣習的に「パワー」と呼びます（これは電気回路の解析に影響された用語です。実際に電気回路で電圧 $$V(t)$$ を測定している場合は、その回路で消費されるパワーは $$(V(t))^2$$ で書くことができます）。この瞬間的なパワーの値の統計的な平均値を、周波数ごとに分割したものがPSDです。前述のようにフーリエ変換の絶対値の自乗に比例するような量で $$ \lim\_{T \to \infty} \dfrac{1}{T}\int\_{-\infty}^{\infty}g^2(t) dt = \int\_{-\infty}^{\infty}{S\_g}(f) df= \dfrac{1}{2\pi}\int\_{-\infty}^{\infty}{S\_g}(\omega) d\omega $$ という性質を満たします。この $$S\_g(f)$$ をtwo-sided power spectral densityと呼びます。 測定値 $$g(t)$$ とす$$S\_g(-f)=S\_g(f)$$ という性質を満たすので、新しく $$S'\_g(f)=2S\_g(f)$$ というone-sided power spectral densityを定義すると $$ \lim\_{T \to \infty} \dfrac{1}{T}\int\_{-\infty}^{\infty}g^2(t) dt = \int\_{0}^{\infty}{S'*g}(f) df= \dfrac{1}{2\pi}\int*{0}^{\infty}{S'\_g}(\omega) d\omega $$ のように正の周波数領域だけを気にすれば良くなります。 pythonではscipyモジュールを使うことでpower spectral densityを計算できます。mydataというnumpy arrayのone-sided power spectral densityを計算するには ``` import scipy.signal as signal fps = 100 segmentlength=1000 mypsd=signal.welch(mydata,fps,'boxcar',segmentlength) ``` とします。mypsd\[0]に周波数のarrayが、mypsd\[1]にPSDのarrayが格納されます。\ 実際にやっていることは以下の図のようになっています。 ![](https://1559620454-files.gitbook.io/~/files/v0/b/gitbook-legacy-files/o/assets%2F-MEGDe0jc4Gum8QxA8RN%2F-MehADo0wjwMtYDpKaTa%2F-MehE_0BO3JEfq6wX3gz%2Fpsd.png?alt=media\&token=02ae42b2-d3bc-4e92-bf2c-3d3536365e52) segmentlength/fpsがフーリエ変換の時間範囲に対応し、周波数分解能を決めます。したがって、segmentlengthは十分な周波数分解能が得られるような値を選ぶ必要があります。窓関数やオーバーラップの具合などもオプションで設定可能です。詳細は\ \ を参照してください