numpy, scipyなどのメモ

工事中

TimeSeriesの計算小技・注意事項

四則演算、指数

gwpyのTimeSeriesや FrequencySeriesは、四則演算や指数の演算がそのまま使える。横軸の配列(timesやfrequencies)が同一なら、TimeSeries同士の演算も可能。

   1 t = np.linspace(start=-10, stop=10, num=201)
   2 A = TimeSeries(np.sin(t), t0=1261872018)                                                                                                                                 
   3 plot = Plot(A, A+1, 2*A, A**2)
   4 ax = plot.gca().legend(['A', 'A+1', '2*A', 'A**2'], fontsize=18)
   5 plot.show()

FFT length の注意

TimeSeriesからスペクトルを計算する際、fftlengthの値が小さいと下図左のように謎の冪スペクトルに埋もれてしまうことがある。(ドリフト+窓関数の影響？)
十分なfftlengthをとればいいのだが、平均回数を稼ぎたいときや短時間のみを見たい場合にはそうもいかない。
そんなときは先に適当なhighpassフィルターをかけることで、下図右のように回避できる。

FFT の窓関数について

窓関数について参考になるページ
Heinzel-Rudiger-Schilling
SciPy で使用可能な窓関数
小野精器計測コラム

window	周波数分解能	主な用途
hanning	△	連続的な、一般的な信号 (とりあえずコレ)
flattop	×	線スペクトルのピーク値に読み取り
boxcar	〇	過渡信号、ハンマリング試験

gwpyのFFTでは、窓関数による減衰は補正されているっぽい。(Note参照)
試しに同じ時系列波形に対してhanningとflattopでASDを作ってband limited RMSを計算してみたら、同じ値になった。(補正なしだと2倍ほど変わるはず)

KAGRA/Subgroups/PEM/PythonMemoJP/calculation (last edited 2021-07-29 15:34:38 by tatsuki.washimi)

-  ⇤ ← Revision 5 as of 2020-06-14 17:04:58 → 
  Size: 1659
  Editor: tatsuki.washimi
  Comment:
+   ← Revision 6 as of 2020-09-10 10:16:37 → ⇥
  Size: 2822
  Editor: tatsuki.washimi
  Comment:
-Deletions are marked like this.
+Additions are marked like this.
 Line 9:
-= TimeSeriesなどの計算小技・注意事項 =
+= TimeSeriesの計算小技・注意事項 =
 Line 32:
+== FFT の窓関数について ==
窓関数について参考になるページ<<BR>>
[[https://holometer.fnal.gov/GH_FFT.pdf|Heinzel-Rudiger-Schilling]] <<BR>>
[[https://org-technology.com/posts/scipy-window-function.html|SciPy で使用可能な窓関数]] <<BR>>
[[https://www.onosokki.co.jp/HP-WK/eMM_back/emm150.pdf|小野精器 計測コラム]]

||'''window''' || '''周波数分解能''' || '''主な用途'''                              ||
|| hanning     || △                 || 連続的な、一般的な信号 (とりあえずコレ) ||
|| flattop     || ×                 || 線スペクトルのピーク値に読み取り            ||
|| boxcar      || 〇                 || 過渡信号、ハンマリング試験                  ||


gwpyのFFTでは、窓関数による減衰は補正されているっぽい。([[https://gwpy.github.io/docs/stable/api/gwpy.timeseries.TimeSeries.html#gwpy.timeseries.TimeSeries.fft|Note参照]])<<BR>>
試しに同じ時系列波形に対してhanningとflattopでASDを作ってband limited RMSを計算してみたら、同じ値になった。(補正なしだと2倍ほど変わるはず)