ymaeda_opentoolsの目的別ガイド
火山性長周期・超長周期地震の力源解析
(Guide of ymaeda_opentools by purposes:
Source analyses of volcanic long-period and very-long-period events)
私がこれまでに作成したコンピュータツールの中でも特に利用のリクエストを受ける頻度が高いのが
火山性長周期・超長周期地震の力源解析のための波形逆解析のコードです。
以下では波形逆解析やその下準備のためのプログラムを紹介します。
The waveform inversion code for the source analysis of
volcanic long-period (LP) and very-long-period (VLP) events
is one of the most frequently requested computer tools.
This code and several additional tools for preparation of the analysis
are introduced below.
◆準備
(Preparation)
波形逆解析のコードに入力として与える地震波形データは
Seismic Analysis Code (SAC)形式で与える必要があります。
その用意については
地震データの基本的な処理
をご覧ください。
The input seismic waveform data for waveform inversion
must be given in Seismic Analysis Code (SAC) format.
See Basic treatments of seismic data
for the preparation of the SAC data.
ymaeda_opentoolsの波形逆解析のアルゴリズムでは直交座標系を使用します。
そのため緯度・経度で与えられている観測点座標を直交座標に換算する必要があります。
ymaeda_opentoolsの
latlon2xyコマンド
を用いると緯度・経度から直交座標への換算を行うことができます。
ここで重要な点は換算の基準点(直交座標の原点)を自由に選べることです。
座標系の変換だけであれば
国土地理院のWEBサイト上のツール
などを用いて行うこともできますが、
変換の基準点
はいくつかの都道府県をまとめた地域毎に固定値で与えられています。
基準点から離れるほど球面である地表を平面近似することによる誤差が大きくなり、
また得られる座標値も直感的に理解しにくい値になります。
ymaeda_opentoolsのlatlon2xyコマンドでは国土地理院と同じ換算式を使用していますが
座標変換の基準点をユーザが指定できます。
火山地域の解析であれば例えば火山の山頂や噴火口の中心点を基準点にすることができます。
これによって平面近似の誤差が小さくなるとともに、
得られる座標も例えば
「山頂から東∗∗km、北∗∗km」
のように理解しやすい値となります。
Because the waveform inversion algorithm of ymaeda_opentools
is based on a cartesian coordinate system,
the coordinates of stations given by the latitudes and longitudes
must be converted to cartesian coordinates.
This conversion is implemented by
latlon2xy command
of ymaeda_opentools.
An advantage of this tool over existing conversion tools, for example
a WEB tool of the Geospatial Information Authority of Japan (GSI),
is that the reference point (the origin of the cartesian coordinate)
can be chosen arbitrarily.
In the GSI, the reference point is given as
a fixed value for each region composed of several prefectures.
The larger the distance from the reference point,
the larger the approximation error of the spherical ground surface
by a flat plane.
In addition, the values of cartesian coordinates are not intuitive.
The conversion of coordinates by
latlon2xy command of ymaeda_opentools
uses the same equation as GSI but with an arbitrary reference point
that can be specified by the user.
For example, the user can use the summit or the center of an active crater
of a volcano as the reference point.
By this, the approximation error by the flat plane is minimized
and the derived coordinates are intuitive, for example
“∗∗ km east and
∗∗ km north from the summit”.
ymaeda_opentoolsの波形逆解析アルゴリズムの特徴の1つは3次元地形を計算に入れられることです。
そのために3次元地形データ(数値標高モデル)を用意する必要があります。
日本国内の解析であれば数値標高モデルは
国土地理院のWEBサイト
からダウンロードできます。
そのときに必要なデータの範囲をメッシュ番号で指定する必要があります。
ymaeda_opentoolsの
latlon2JapanMeshCodeコマンド
を用いると、解析対象地域の代表的な地点の緯度経度を与えて
その地点が含まれるメッシュ番号を調べることができます。
A feature of the waveform inversion algorithm of ymaeda_opentools
is that 3-D topography can be taken into account.
To realize this, a 3-D topography data (a digital elevation model; DEM)
must be prepared.
The DEM data in Japan can be downloaded from
the WEB site of GSI.
When downloading the data,
the user must specify the spatial range of the data by mesh numbers.
The mesh number that includes a given latitude and longitude
can easily be surveyed by
latlon2JapanMeshCode command of ymaeda_opentools.
国土地理院のWEBサイトからダウンロードしたオリジナルの数値標高モデルはXML形式になっています。
これを緯度・経度・標高の組の一覧表形式のテキストデータに変換することで
解析での利用が可能になり、プロットも容易になります。
この変換にはymaeda_opentoolsの
GSIdem2latlonDataコマンド
が利用できます。
The original DEM data downloaded from the GSI's website is in XML format.
It has to be converted to a table text format
composed of the combination of a latitude, a longitude,
and the elevation of the ground surface at that point.
For this conversion, use
GSIdem2latlonData command of ymaeda_opentools.
◆順問題
(Forward problem)
逆問題を解くためには前提として順問題が解けている必要があります。
波形逆解析の場合であれば順問題に相当するのは
地震波動ソースを与えて波動場を求める理論波形計算です。
波形逆解析のプロの間で通称「グリーン関数」と呼ばれているものですが、
実際にはグリーン関数(インパルス応答)ではなく、与えた震源時間関数のもとでの理論波形です。
ymaeda_opentoolsでは想定される全てのソース位置の候補と観測点のペアについて
前もって理論波形(通称「グリーン関数」)を計算しておきます。
そのために
waterPMLコマンド
を使用します。
このコマンドでは3次元地形を使用し、
Perfectly Matched Layer (PML)を用いて理論波形を計算できます。
並進運動(速度)だけでなく傾斜変動も計算できる点が特色です。
また水域も計算に含めることができます。
アルゴリズムとして
Maeda et al. (2011), Maeda and Kumagai (2013)
で提唱したものを用いております。
Solving an inverse problem requires
the corresponding forward problem to have been solved.
In waveform inversion, the corresponding forward problem
is the computation of synthetic waveforms from a given source.
The synthetic waveforms for this purpose
are called “Green's functions”
among the professionals of waveform inversion,
but they are actually not Green's functions (responses to an impulsive source)
but synthetic waveforms caused by a given source time function.
In ymaeda_opentools, the synthetic waveforms
(so-called “Green's functions”)
for all pairs of potential source locations and stations
have to be computed before performing inversion.
The waterPML command
is available for this computation.
It computes synthetic waveforms using 3-D topography
and Perfectly Matched layers (PMLs)
and can incorporate water-filled regions.
The algorithm is based on Maeda et al. (2011) and Maeda and Kumagai (2013).
waterPMLコマンドはパラメータも多く、一見難しそうに見えるかもしれません。
同コマンドのマニュアル
内にプログラムを動かすためのbashスクリプトを具体的に示した例題を用意しておりますので、
それを見ながらプログラムを動かしてみると理解が早いかと思われます。
Beginners may feel waterPML command to be difficult to use
because it consists of many parameters.
Therefore, the
documentation of waterPML command
includes several examples with explicit bash scripts
to execute the program.
Exercise following these examples is probably the best way
to learn the usage of the program.
理論波形計算の結果を検証する1つのアプローチは
「正解」が分かっている単純な媒質においてwaterPMLコマンドにより理論波形を計算し、
結果をその「正解」と比較することです。
例えば無限等方均質媒質であればAki and Richards (2002)で解析解が与えられていますので
それと比較することができます。
ymaeda_opentoolsではAki and Richards (2002)の解析解を
指定した地震波動ソースと観測点の組合せについて計算するプログラムとして
WIHMコマンド
を用意しております。
waterPMLコマンドとWIHMコマンドの計算は独立です。
したがって無限等方均質媒質中の理論波形をwaterPMLコマンドとWIHMコマンドで計算・比較すれば
waterPMLコマンドによる計算結果の1つの検証になります。
なおwaterPMLコマンドでは速度波形、WIHMコマンドでは変位波形を出力する点に留意してください。
速度波形から変位波形への変換(積分)には
sacfile_integralコマンド
や
sequencefile_integralコマンド
を利用できます。また変位波形から速度波形への変換(微分)には
sacfile_differentiateコマンド
や
sequencefile_differentiateコマンド
を利用できます。
An approach to examine the computed synthetic waveforms
is to calculate them with the waterPML command
assuming a simple medium for which the “correct answer” is known
and compare them.
For example, an analytical solution for the synthetic wave field
in an infinite, isotropic, and homogeneous medium
is given by Aki and Richards (2002).
This analytical solution,
for a specified pair of a seismic wave source and a station,
can be computed with
WIHM command
of ymaeda_opentools.
The computations by waterPML and WIHM commands are independent.
Therefore, a comparison of synthetic waveforms
in an infinite, isotropic, and homogeneous medium
calculated with waterPML and WIHM commands
is an approach to examine the outputs from waterPML command.
Note that waterPML command outputs the velocity waveforms,
whereas WIHM command outputs the displacement waveforms.
A velocity waveform can be converted (integrated) to a displacement waveform
using
sacfile_integral
or
sequencefile_integral
command, and
a displacement waveform can be converted (differentiated)
to a velocity waveform using
sacfile_differentiate
or
sequencefile_differentiate
command.
◆逆問題
(Inverse problem)
以上の準備のもとでようやく逆問題の解析に取り掛かれます。
ymaeda_opentoolsの波形逆解析のプログラムは
winvコマンド
になります。
これも多数のパラメータを含んでおり初学者には難しく感じられるかもしれませんが、
waterPMLコマンドと同様にプログラムを動かすための具体的なbashスクリプトを含めて
例題を載せておりますのでご参照ください。
The inverse problem can be solved after these preparations.
A computer program to conduct waveform inversion in ymaeda_opentools
is winv command.
This command may also look difficult to beginners
because of many parameters,
but users can learn the usage of the program through examples,
available in the documentation,
in which explicit bash scripts to execute the program are given.