Note

このページは Jupyter Notebook から生成されました。 ノートブックをダウンロード (.ipynb)

gwexpy Plotting Introduction

Open In Colab

このノートブックでは、gwexpy の強化されたプロット機能を紹介します。 TimeSeriesMatrixFrequencySeriesMatrix などの行列データ構造、および SpectrogramList などのコレクションに対する、「いい感じ」のデフォルト設定(figsize自動調整、対数軸、レイアウト、軸ラベル)についてデモを行います。

[1]:
import matplotlib.pyplot as plt
import numpy as np
from astropy import units as u

from gwexpy.frequencyseries import (
    FrequencySeries,
    FrequencySeriesList,
    FrequencySeriesMatrix,
)
from gwexpy.spectrogram import Spectrogram, SpectrogramList
from gwexpy.timeseries import TimeSeries, TimeSeriesMatrix

#
np.random.seed(42)

1. TimeSeriesMatrix Plotting

TimeSeriesMatrix は (Row, Col, Time) の3次元構造を持つ時系列データを管理します。 .plot() を呼び出すと、以下の最適化が行われます:

  • 軸ラベル: determine_xlabel により “Time [s]” や “Time [s] from …” が自動設定されます。

  • Auto-GPS: 時間軸の場合、デフォルトで auto-gps スケールが適用され、GPS時刻(大きなオフセット)が適切に処理されます。

  • レイアウト: 行・列の構造に合わせて figsize が自動調整されます。

[2]:
# 2x2 TimeSeriesMatrix Create
# GPS ( : 1356998418 = 2023-01-01 00:00:00 UTC) from 10seconds
t0 = 1356998418.0
times = t0 + np.linspace(0, 10, 1000)
data = np.random.randn(2, 2, 1000)
# : Sensor A, B / : Axis X, Y
tsm = TimeSeriesMatrix(
    data,
    times=times,
    unit="m",
    name="Velocity",
    rows=["Sensor A", "Sensor B"],
    cols=["X-axis", "Y-axis"],
)

# Plot
plot = tsm.plot()
plot.show();
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_3_0.png

2. FrequencySeriesMatrix Plotting

FrequencySeriesMatrix は (Row, Col, Frequency) の3次元構造を持つ周波数データを管理します。 一般的なスペクトル解析の結果などを保持します。 .plot() の特徴:

  • 対数軸: データサイズ等の条件に応じて、X軸(周波数)やY軸(振幅)が自動的に対数スケールになります。

  • 軸ラベル: “Frequency [Hz]” や “[unit]” が自動設定されます。

[3]:
# FrequencySeriesMatrix Create
# 10Hz from 1kHz
f = np.logspace(1, 3, 500)
# 1/f Noise Data Generate
data = np.zeros((2, 2, 500))
for i in range(2):
    for j in range(2):
        data[i, j, :] = 10 / f * (1 + 0.1 * np.random.randn(500)) + (i + j)

# unit='m' (Amplitude Spectral Density do m/sqrt(Hz) 、 )
fsm = FrequencySeriesMatrix(
    data,
    frequencies=f,
    unit="m",
    name="Displacement",
    rows=["Sensor A", "Sensor B"],
    cols=["X-axis", "Y-axis"],
)

# Plot
# log-log Plot Verify
plot = fsm.plot()
plot.show();
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_5_0.png

3. FrequencySeriesList Plotting (参考)

FrequencySeriesList (または Dict) は、グリッド構造を持たない単純なリストです。 これも同様にプロット可能です。

[4]:
# FrequencySeriesList Create
fs1 = FrequencySeries(1 / f, frequencies=f, unit=u.m, name="Series 1")
fs2 = FrequencySeries(2 / f, frequencies=f, unit=u.m, name="Series 2")

fs_list = FrequencySeriesList([fs1, fs2])

fs_list.plot();
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_7_0.png
[5]:
del times, data, tsm, fsm, fs1, fs2, fs_list

import gc

gc.collect()
[5]:
18018

4. SpectrogramList Plotting

SpectrogramList のプロットは、デフォルトで geometry=(N, 1) の縦並びレイアウト、対数カラーバーなどが適用されます。 スペクトログラムも時間軸を持つため、GPS時刻を与えると auto-gps が機能します。

[6]:
# Spectrogram Create
# TimeSeriesMatrix GPS
t = t0 + np.linspace(0, 10, 100)
# 10Hz - 1000Hz
f = np.logspace(1, 3, 50)

spec_list = SpectrogramList(
    [
        Spectrogram(
            np.random.rand(100, 50), times=t, frequencies=f, unit=u.count, name="Spec 1"
        ),
        Spectrogram(
            np.random.rand(100, 50) * 10,
            times=t,
            frequencies=f,
            unit=u.count,
            name="Spec 2",
        ),
    ]
)
spec_list.plot();
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_10_0.png

4. Adaptive Plotting Optimization (Decimation)

大量のデータを持つ TimeSeries をプロットする際、描画パフォーマンスを維持するために自動的に間引き(Min-Max Decimation)を適用します。

  • 自動適用: デフォルトではサンプル数が 50,000 を超える場合に発動します。

  • ピーク維持: Min-Max アルゴリズムにより、波形のエンベロープ(最大・最小)を保持したままデータ点を削減(デフォルト 約10,000点)します。

  • カスタマイズ: decimate_thresholddecimate_points 引数で制御可能です。

[7]:
# 100,000,000 Sample Data Generate
ts_large = TimeSeries(np.random.randn(100000000), sample_rate=1024, name="Large TS")

print(f"Original samples: {len(ts_large)}")

# PlotRun ( )
plot_opt = ts_large.plot()

# donepoints Verify
ax = plot_opt.gca()
line = ax.get_lines()[0]
print(f"Plotted points: {len(line.get_xdata())}")
plot_opt.show();
Original samples: 100000000
Plotted points: 10000
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_12_1.png

コレクション・行列への適用とパラメータ指定

TimeSeriesListTimeSeriesMatrix でも透過的に機能します。また、閾値を下げて強制的に適用させることも可能です。

[8]:
# Threshold 10,000 、Plotpoints 2,000 do
ts_med = TimeSeries(np.random.randn(20000), sample_rate=16384, name="Medium TS")

plot_custom = ts_med.plot(decimate_threshold=10000, decimate_points=2000)

ax = plot_custom.gca()
print(f"Original: {len(ts_med)}, Plotted: {len(ax.get_lines()[0].get_xdata())}")
plot_custom.show();
Original: 20000, Plotted: 2000
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_14_1.png

5. Spectrogram Summary Plotting

plot_summary は、複数の Spectrogram を垂直に並べ、それぞれの右側に FrequencySeries の統計量(10%, 50%, 90% パーセンタイル)をプロットする機能です。 スペクトログラムの時間変化と、その全体的な周波数特性を並べて確認する際に非常に便利です。

  • 自動同期: スペクトログラムのカラー軸と FrequencySeries の縦軸が自動的に同期します。

  • MMMプロット: plot_mmm を用いて、Min/Median/Max(または任意のパーセンタイル)を塗りつぶし付きで表示します。

[9]:
from gwexpy.spectrogram import SpectrogramDict
from gwexpy.timeseries import TimeSeries


# Data Generate
def make_sg(name):
    ts = TimeSeries(np.random.randn(2048 * 10), sample_rate=2048, name=name, unit="m")
    return ts.spectrogram(stride=1, fftlength=1, overlap=0.5)


sg_dict = SpectrogramDict({"Channel A": make_sg("A"), "Channel B": make_sg("B")})

# ASDgram Plot Run
fig, axes = sg_dict.plot_summary(fmin=10, fmax=500, title="ASDgram Demonstration")
plt.show();
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_16_0.png

6. PairPlot Example

以下は gwexpy.plot.PairPlot の基本的な使用例です。

  • 2 つの TimeSeries を異なるサンプリングレートで作成

  • PairPlot が自動で最小サンプルレートにリサンプリングし、共通区間にトリミングします

  • corner=True(デフォルト)により左下三角だけが描画されます

実行すると、対角にヒストグラム、非対角に 2D ヒストグラムが表示されます。

[10]:
import numpy as np

from gwexpy.plot import PairPlot
from gwexpy.timeseries import TimeSeries

# 2 TimeSeries Create
ts1 = TimeSeries(np.random.randn(5000), sample_rate=1024, name="TS1")
ts2 = TimeSeries(np.random.randn(4000), sample_rate=2048, name="TS2")

# PairPlot Create(corner=True )
pair = PairPlot([ts1, ts2], corner=True)

#
pair.show()
../../../../_images/web_ja_user_guide_tutorials_intro_plotting_18_0.png
[10]:
<gwexpy.plot.pairplot.PairPlot at 0x7ed1287c7080>