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.