Differences between revisions 3 and 14 (spanning 11 versions)
Revision 3 as of 2020-03-21 13:11:40
Size: 1026
Comment:
Revision 14 as of 2021-08-11 15:28:39
Size: 6556
Comment:
Deletions are marked like this. Additions are marked like this.
Line 3: Line 3:

=== 複数の時系列を縦に並べてplotする ===
[[https://gwpy.github.io/docs/stable/timeseries/plot.html#plotting-multiple-timeseries-together|ここ]]に書いてある、gwpy.plotのPlot()がオススメ。<<BR>>
GUIで拡大等をしたとき、他のaxesも一緒に動く。<<BR>>
注意点は、TimeSeriesのlistやTimeSeriesDictを引数にしてまとめて渡した場合は ''separate=True'' をしても分割されないという謎仕様。面倒だが1つ1つ手で並べる必要がある。

{{{#!python
from gwpy.plot import Plot
plot = Plot(l1hoft, h1hoft, separate=True, sharex=True)
plot.show()
}}}
{{https://gwpy.github.io/docs/stable/_images/plot-7.png}}

戻り値 (ここでは ''plot'')は、[[https://matplotlib.org/3.2.1/api/_as_gen/matplotlib.pyplot.figure.html|matplotlibのfigure]]だと思ってしまってよい。<<BR>>
なので、それぞれのy軸ラベルや全体のタイトルを付けるには以下のようにする。
{{{#!python
import matplotlib.pyplot as plt
axes = plot.get_axes()
ylabel = ['y1','y2']
for i, ax in enumerate(axes): ax.set_ylabel(ylabel[i])
plt.tight_layout()
plt.subplots_adjust(top=0.95)
plot.suptitle('title')
}}}
Line 10: Line 34:
もしくは[[https://gwpy.github.io/docs/stable/plot/gps.html|ここ]]参照
Line 12: Line 37:
''.add_segments_bar()'' を使う([[https://gwpy.github.io/docs/latest/examples/miscellaneous/open-data-spectrogram.html|ここ]]を参照)。 ''.add_segments_bar()'' を使う([[https://gwpy.github.io/docs/latest/examples/miscellaneous/open-data-spectrogram.html|ここ]]を参照)。<<BR>>
Line 14: Line 39:
Line 16: Line 42:

=== スペクトルのパーセンタイル/σバンドを付ける ===
[[https://gwpy.github.io/docs/stable/examples/frequencyseries/percentiles.html|ここ]]を参照。手順としては、
 1. まずTimeSeriesから[[https://gwpy.github.io/docs/stable/api/gwpy.timeseries.TimeSeries.html#gwpy.timeseries.TimeSeries.spectrogram2|spectrogram2()]]を使ってnon-averaged power Spectrogramを作る。
 2. .percentile() でパーセンタイルのスペクトルを3つ作る。
 3. plot_mmm() で描画するとバンド付きでplotされる。

{{{#!python
from gwpy.timeseries import TimeSeries
from gwpy.plot import Plot

hoft = TimeSeries. ...
sg = hoft.spectrogram2(fftlength=4, overlap=2, window='hanning') ** (1/2.)
median = sg.percentile(50)
low = sg.percentile(5)
high = sg.percentile(95)

plot = Plot()
ax = plot.gca(xscale='log', xlim=(10, 1500), xlabel='Frequency [Hz]', yscale='log', ylim=(3e-24, 2e-20), ylabel=r'Strain noise [1/$\sqrt{\mathrm{Hz}}$]')
ax.plot_mmm(median, low, high, color='gwpy:ligo-hanford')
ax.set_title('LIGO-Hanford strain noise variation around GW170817', fontsize=16)
plot.show()
}}}
{{https://gwpy.github.io/docs/stable/_images/percentiles-4.png}}

パーセンタイルではなくmeanとσで表したい場合は、Spectrogramの2次元配列に対して時間軸(0軸)方向にmean, stdを取ればよい。
{{{#!python
mean = FrequencySeries(data=sg.mean(0), frequencies=sg.frequencies)
sigma = FrequencySeries(data=sg.std(0), frequencies=sg.frequencies)
}}}
Line 17: Line 74:
[[https://qiita.com/simonritchie/items/da54ff0879ad8155f441|プロットの複雑なレイアウトはGridSpecが便利かも、という話]]  * [[https://qiita.com/simonritchie/items/da54ff0879ad8155f441|プロットの複雑なレイアウトはGridSpecが便利かも、という話]]
 * [[https://bunseki-train.com/axvspan-and-axhspan/|matplotlibで一定区間に背景色をつける方法]]
 * [[https://qiita.com/nkay/items/d1eb91e33b9d6469ef51|matplotlibのめっちゃまとめ]]
Line 19: Line 78:
[[https://bunseki-train.com/axvspan-and-axhspan/|matplotlibで一定区間に背景色をつける方法]]
=== imshow を使って2Dヒストグラムを書く時の軸の取り方 ===
以下のような2変数f_in, f_outの間の相関データを、imshowを使って2Dヒストグラムとして描画したいとしよう。
{{{#!python
matrix = [[0, 1, 0, 1],
          [1, 0, 1, 0],
          [0, 1, 0, 1],
          [1, 0, 2, 0],
          [0, 1, 0, 2],
          [1, 0, 2, 0]]
f_in = [1,2,3,4]
f_out = [1,2,3,4,5,6]
}}}
素直に以下のように書くと、縦軸が反転して(上から下になる)しまい、bin幅もごみになる。
{{{#!python
fig = Plot(figsize=(5, 5))
ax = fig.gca()
ax.imshow(matrix, extent=[f_in[0], f_in[-1], f_out[0], f_out[-1]], aspect='auto', interpolation='nearest')
ax.set_xlabel("Injected Noise Frequency $f'$ [Hz]")
ax.set_ylabel('Obsurved Frequency $f$ [Hz]')
fig.show()
}}}
{{attachment:imshow_bad.png}}

縦軸の向きは origin='lower' を入れると下から上になる。<<BR>>
bin幅については、上下左右に半分ずつ追加することで変数の値(今の場合はf_in, f_out)がbin中心になるようにできる。
{{{#!python
df_in = f_in[1] - f_in[0]
df_out = f_out[1] - f_out[0]

f_L = f_in[0] - df_in*0.5
f_R = f_in[-1]+ df_in*0.5
f_B = f_out[0] - df_out*0.5
f_T = f_out[-1]+ df_out*0.5

fig = Plot(figsize=(5, 5))
ax = fig.gca()
ax.imshow(matrix, extent=[f_L, f_R, f_B, f_T], aspect='auto', origin='lower', interpolation='nearest')
ax.set_xlabel("Injected Noise Frequency $f'$ [Hz]")
ax.set_ylabel('Obsurved Frequency $f$ [Hz]')
fig.show()
}}}

{{attachment:imshow_good.png}}



=== TitleやLabelなどにアンダーバー(_)が含まれてると怒られる ===
環境によっては起きない問題だが、例えば神岡のK1sum*とかで発生する。<<BR>>
pythonやmatplotlibのバージョンの問題化と思いきや、virtualenvで最新環境にしても生じるので原因はよくわからない。

とりあえず、頭に以下のおまじないを付けると解決する。
{{{#!python
import matplotlib.pyplot as plt
plt.rc('text.latex', preamble=r'\usepackage{underscore}')
}}}


=== 正負両方の対数プロット ===
[[https://helve-blog.com/posts/python/matplotlib-object-oriented-logarithmic-graph/#%E6%AD%A3%E8%B2%A0%E4%B8%A1%E6%96%B9%E3%81%AE%E5%AF%BE%E6%95%B0%E3%83%97%E3%83%AD%E3%83%83%E3%83%88symlog]]
{{{#!python
x = np.arange(-10, 10)
y = np.arange(-10, 10)

fig, ax = plt.subplots()
ax.plot(x, y)
ax.set_yscale("symlog", linthreshy=1)
ax.grid()
plt.show()
}}}
{{attachment:xscale_symlog.png}}

Plotの描画関連

複数の時系列を縦に並べてplotする

ここに書いてある、gwpy.plotのPlot()がオススメ。
GUIで拡大等をしたとき、他のaxesも一緒に動く。
注意点は、TimeSeriesのlistやTimeSeriesDictを引数にしてまとめて渡した場合は separate=True をしても分割されないという謎仕様。面倒だが1つ1つ手で並べる必要がある。

   1 from gwpy.plot import Plot
   2 plot = Plot(l1hoft, h1hoft, separate=True, sharex=True)
   3 plot.show()

https://gwpy.github.io/docs/stable/_images/plot-7.png

戻り値 (ここでは plot)は、matplotlibのfigureだと思ってしまってよい。
なので、それぞれのy軸ラベルや全体のタイトルを付けるには以下のようにする。

   1 import matplotlib.pyplot as plt
   2 axes = plot.get_axes()
   3 ylabel = ['y1','y2']
   4 for i, ax in enumerate(axes): ax.set_ylabel(ylabel[i])
   5 plt.tight_layout()
   6 plt.subplots_adjust(top=0.95)
   7 plot.suptitle('title')

時刻の基準を手で設定する

ここを参照

   1 ax = plot.gca()
   2 ax.set_epoch(1126259462) # GPS timeを入れる

もしくはここ参照

DQ Flagを添える

.add_segments_bar() を使う(ここを参照)。
他にも .add_state_segments().add_dataqualityflag() といった似たものがある。詳細はPlotting APIを参照。

https://gwpy.github.io/docs/latest/_images/open-data-spectrogram-5.png

スペクトルのパーセンタイル/σバンドを付ける

ここを参照。手順としては、

  1. まずTimeSeriesからspectrogram2()を使ってnon-averaged power Spectrogramを作る。

  2. .percentile() でパーセンタイルのスペクトルを3つ作る。
  3. plot_mmm() で描画するとバンド付きでplotされる。

   1 from gwpy.timeseries import TimeSeries
   2 from gwpy.plot import Plot
   3 
   4 hoft = TimeSeries. ...
   5 sg = hoft.spectrogram2(fftlength=4, overlap=2, window='hanning') ** (1/2.)
   6 median = sg.percentile(50)
   7 low = sg.percentile(5)
   8 high = sg.percentile(95)
   9 
  10 plot = Plot()
  11 ax = plot.gca(xscale='log', xlim=(10, 1500), xlabel='Frequency [Hz]',  yscale='log', ylim=(3e-24, 2e-20),  ylabel=r'Strain noise [1/$\sqrt{\mathrm{Hz}}$]')
  12 ax.plot_mmm(median, low, high, color='gwpy:ligo-hanford')
  13 ax.set_title('LIGO-Hanford strain noise variation around GW170817', fontsize=16)
  14 plot.show()

https://gwpy.github.io/docs/stable/_images/percentiles-4.png

パーセンタイルではなくmeanとσで表したい場合は、Spectrogramの2次元配列に対して時間軸(0軸)方向にmean, stdを取ればよい。

   1 mean  = FrequencySeries(data=sg.mean(0), frequencies=sg.frequencies)
   2 sigma = FrequencySeries(data=sg.std(0),  frequencies=sg.frequencies)

matplotlib色々

imshow を使って2Dヒストグラムを書く時の軸の取り方

以下のような2変数f_in, f_outの間の相関データを、imshowを使って2Dヒストグラムとして描画したいとしよう。

   1 matrix = [[0, 1, 0, 1],
   2           [1, 0, 1, 0],
   3           [0, 1, 0, 1],
   4           [1, 0, 2, 0],
   5           [0, 1, 0, 2],
   6           [1, 0, 2, 0]]
   7 f_in  = [1,2,3,4]
   8 f_out = [1,2,3,4,5,6]

素直に以下のように書くと、縦軸が反転して(上から下になる)しまい、bin幅もごみになる。

   1 fig = Plot(figsize=(5, 5))
   2 ax = fig.gca()
   3 ax.imshow(matrix, extent=[f_in[0], f_in[-1], f_out[0], f_out[-1]], aspect='auto', interpolation='nearest')
   4 ax.set_xlabel("Injected Noise Frequency $f'$ [Hz]")
   5 ax.set_ylabel('Obsurved Frequency $f$ [Hz]')
   6 fig.show()

imshow_bad.png

縦軸の向きは origin='lower' を入れると下から上になる。
bin幅については、上下左右に半分ずつ追加することで変数の値(今の場合はf_in, f_out)がbin中心になるようにできる。

   1 df_in  = f_in[1]  - f_in[0]
   2 df_out = f_out[1] - f_out[0]
   3 
   4 f_L = f_in[0] - df_in*0.5
   5 f_R = f_in[-1]+ df_in*0.5
   6 f_B = f_out[0] - df_out*0.5
   7 f_T = f_out[-1]+ df_out*0.5
   8 
   9 fig = Plot(figsize=(5, 5))
  10 ax = fig.gca()
  11 ax.imshow(matrix, extent=[f_L, f_R, f_B, f_T], aspect='auto', origin='lower', interpolation='nearest')
  12 ax.set_xlabel("Injected Noise Frequency $f'$ [Hz]")
  13 ax.set_ylabel('Obsurved Frequency $f$ [Hz]')
  14 fig.show()

imshow_good.png

TitleやLabelなどにアンダーバー(_)が含まれてると怒られる

環境によっては起きない問題だが、例えば神岡のK1sum*とかで発生する。
pythonやmatplotlibのバージョンの問題化と思いきや、virtualenvで最新環境にしても生じるので原因はよくわからない。

とりあえず、頭に以下のおまじないを付けると解決する。

   1 import matplotlib.pyplot as plt
   2 plt.rc('text.latex', preamble=r'\usepackage{underscore}')

正負両方の対数プロット

https://helve-blog.com/posts/python/matplotlib-object-oriented-logarithmic-graph/#%E6%AD%A3%E8%B2%A0%E4%B8%A1%E6%96%B9%E3%81%AE%E5%AF%BE%E6%95%B0%E3%83%97%E3%83%AD%E3%83%83%E3%83%88symlog

   1 x = np.arange(-10, 10)
   2 y = np.arange(-10, 10)
   3 
   4 fig, ax = plt.subplots()
   5 ax.plot(x, y)
   6 ax.set_yscale("symlog", linthreshy=1)
   7 ax.grid()
   8 plt.show()

xscale_symlog.png

KAGRA/Subgroups/PEM/PythonMemoJP/plot (last edited 2021-08-11 15:28:39 by tatsuki.washimi)