Q&A

地震時動水圧の観測データとフィルタ処理について
水利基盤チーム

はじめにラインにおいて、地震時動水圧1）の観測を実施しています。図－1および図－2に示しているのは、地震が発2017ᖺ2᭶27᪥10:53:24

0.010ືỈᅽ

水利基盤チームでは、供用されている農業用パイプ

[MPa]

0.02

生した際に観測されたパイプラインの水圧の時刻歴デ

-0.01
-0.02
0ータです。どちらのデータも観測地点近傍の震度観測

50

100

150

200

250300᫬㛫 t [sec.]

点で震度2の地震動が記録された際に取得されまし

図－1

地震時動水圧の波形た。地震が発生するとパイプライン内には0.1～2.0秒周期の圧力の増減が発生します。

認できます。地震時動水圧に着目したい場合、その他の余計な圧力波形を除外する必要があります。このよ
うに不要な周波数（周期）の成分を除去し、目的の成分を抽出する処理をフィルタ処理といいます。
本Q&Aでは、図－1および図－2に示す地震時動水圧観測データについて解説するとともに、図－2のデータに対してフィルタ処理を行う一手法を紹介します。

[MPa]が30～40秒程度の水圧の増減が重畳していることが確0.02

2015ᖺ8᭶14᪥13:09:27

0.010ືỈᅽ

さて、図－2を見ると地震時動水圧に加えて、周期-0.01

࿘ᮇ30~40sࡢỈᧁᅽ-0.02
050

100

150

200

250300᫬㛫 t [sec.]

図－2 地震時動水圧＋水撃圧の波形（標本点数N＝30,000、標本点間隔Δt＝0.01s）

また、本稿に記載しているコードは、以下の環境で動作することを確認しています。
Python 3.8.0、

NumPy 1.18.0

水撃圧、地震時動水圧はどちらも水撃現象により発生する水圧の変化と考えられます。
水撃圧はパイプラインの水を使用するためにバルブ

Ｑ１：図－2で確認できる周期30～40秒の圧力振動はを開閉する際、あるいはポンプを動作させるなどで流速が変化した際に発生します。例えば、バルブを閉め何ですか？ていくと流速の減少に伴い圧力が増加していきます。
Ａ１：この圧力振動は水撃圧だと考えられます。水撃そして全閉状態となると圧力振動が発生します。圧力圧は、例えば水を利用しようとする際にバルブを開閉

振動周期T0は、することで発生します。そのため、営農作業等が行わ（1）れる期間の日中には水撃圧が頻繁に観測されます。図－1のような冬期に取得したデータには、水撃圧はほととなります2）。ここで、Lは管路延長、aは圧力伝播速んどみられません。

度です。圧力伝播速度は管の物性や形状などから計算されます。このように水撃圧の振動周期はパイプライ

Ｑ２：水撃圧と地震時動水圧にはどのような違いがあンの物性や配置、規模により決まります。
一方で、地震時動水圧は地震動によりパイプラインりますか？

が振動し、内部の水が管壁により強制的に変位させら
Ａ２：管路内の水の流速が急に変化したことにより圧れることで発生します。そのため、地震時動水圧の振力が変化する現象を水撃現象といいます。すなわち、動周期は地震の特性に依存することになります。

寒地土木研究所月報№820

2021年７月 35

実際に供用されるパイプラインはLが大きいため、水撃圧の周期は地震時動水圧より長くなります。
Ｑ３：図－2のデータからどのように水撃圧を除外できますか？
Ａ３：水撃圧は地震時動水圧より周期が長い（周波数が低い）ため、図－2の水圧データ（以下、p(t)とします）に含まれる低周波数成分を除去することで水撃圧を削除することができます。このような処理は、ある

図－3周波数より高い周波数成分を残すことからハイパスフ

フィルタリング手順

ィルタ（ローカットフィルタ）といいます。周波数に基づいてデータを操作したいため、p(t)を周波数領域に変換する作業が必要です。
計算機を用いて時系列データを周波数成分に変換する手法として、離散フーリエ変換が挙げられます。標import numpy as np
def high_pass(P, dt, fcut):
N = len(P)freq = np.fft.rfftfreq(N, dt)

本点数N、標本点間隔Δtの実数からなる時系列データ

FP = np.fft.rfft(P)#ձ 㞳ᩓࣇ࣮࢚ࣜኚ᥮-> P(ω)ࢆィ⟬x(t)＝｛x0, x1,…, xN－1｝を離散フーリエ変換することによ

FQ = FP.copy()
FQ[freq <= fcut] = 0.0#ղ fcut௨ୗࡢ࿘Ἴᩘᡂศࢆ0࡟ࡍࡿ

散フーリエ変換、逆変換はそれぞれ式（2）、（3）のよ
うになります。
（2）

Q = np.fft.irfft(FQ)return Q
if __name__ == "__main__":P = [p0, p1, ..., pN-1]
dt = 0.01
fcut = 0.1

2+1種類の周波数を持つ波に分解することができま

0.01

p(t)

0.02

1010リエ係数Xkは、周波数kΔt/Nの波の強度や位相の情報-0.01を持ちます（ただし0 ≤ k ≤ N/2です）。逆に言うと、

-0.02



00.02のように時系列データを周波数成分に分解し、各成分0.01

したがってp(t)から水撃圧の成分を取り除く手順
は、図－3に示す通りです。まず、離散フーリエ変換によりp(t)のフーリエ係数P(ω)を求め（①）、所定のq(t)波を足し合わせることで再現することができます。こ

0.1
0.01

元々のデータx(t)は、振幅が｜Xk｜、周波数がkΔt/Nのいます。

# ࣇ࢕ࣝࢱฎ⌮ᚋࡢᅽຊࢹ࣮ࢱ
# ᚋࡣฟຊࡋࡓࡾࠊࢢࣛࣇᥥ⏬ࡍࡿ

図－4 プログラム例（データの入出力等は省略しています）

離散フーリエ変換により、時系列データx(t)をN/の強度や位相などを調べることをスペクトル解析とい# ᅽຊࢹ࣮ࢱ㸦ࣇ࢓࢖ࣝ࡞࡝࠿ࡽㄞࡳ㎸ࡳ㸧# ᶆᮏⅬ㛫㝸[s]# ࢝ࢵࢺ࢜ࣇ࿘Ἴᩘ[Hz]Q = high_pass(P, dt, fcut)（3）

す。すなわち、式（2）によって計算できる複素フー#ճ 㞳ᩓࣇ࣮࢚ࣜ㏫ኚ᥮-> q(t)ࢆィ⟬|P(ω)|

ることができます。虚数単位をj＝√－1とすると、離100200

᫬㛫 t[sec.]

0.0010.001

3000.010.1

ṇつ໬ゅ࿘Ἴᩘ ω1[rad]

100



-0.01|Q(ω)|

り複素フーリエ係数 X(ω)＝｛X0, X1,…, XN－1｝を計算す

10.1
0.01

-0.02
0

100200᫬㛫 t[sec.]

図－5

300

0.001
0.001

0.010.1ṇつ໬ゅ࿘Ἴᩘ ω1[rad]フィルタ処理適用例

周波数より小さい周波数に相当するフーリエ係数を0とし（②）、フーリエ逆変換により時系列信号を復元する（③）ことで水撃圧を取り除いた圧力波形を得ることができます。
36

寒地土木研究所月報№820

2021年７月

Ｑ４：フィルタ処理の計算はどのように実装できます

おわりにか？
本稿では、現地で観測した水圧データから水撃圧をＡ４：離散フーリエ変換は高速フーリエ変換（FFT）簡易的に取り除く手法について解説しました。地震時により効率的に計算することができます。ここではプ

動水圧は頻繁に観測できるデータではないので、過3）ログラミング言語の一つであるPython を用いて計算4）するプログラムを示します。PythonではNumPy と大・過少評価を防ぐためにも上記のようなフィルタ処理は重要だと考えられます。

呼ばれる数値計算を行うためのモジュールを使用する（文責：寺田健司）ことができます。NumPyで提供されている離散フーリエ変換に関する機能を用いて、フィルタ処理を行う

参考文献

プログラムの例を図－4に示します。カットオフ周波数は0.1Hzとしています。言い換えると周期10秒以上の1） 大久保天、中村和正、今泉祐治、寺田健司、川口

成分をカットし、10秒以下の成分を残すということで

清美：農業用管水路で生じる地震時動水圧、農業す。

農村工学会論文集、310（88-1）、pp.I_135-I_144、

図－5はフィルタ処理の各段階のデータを描画したものです。図－5左上は現地観測で取得した水圧データ

2020．
2） 農林水産省農村振興局整備部設計課：土地改良事（図－2と同一）、左下がフィルタ処理後の水圧デー

業計画設計基準及び運用・解説

タです。水撃圧の成分が概ね除去できていることが確イン」、農業農村工学会、pp.218-263、2009．

認できます。

設計「パイプラ

3） Python、https://www.python.org/（2021/05/18参照）
4） NumPy API reference（Discrete FourierTransform）、https://numpy.org/doc/stable/reference/routines.fft.html（2021/05/18 参照）

寒地土木研究所月報№820

2021年７月 37