We need to know the spatial and temporal behaviors of earthquake’s rupture in order to study about the physics of earthquake sources.
This is one of the fundamental methods proposed by Professor M. Kikuchi and Professor Hiroo Kanamori [Link] to calculate teleseismic body-wave inversion for modelling the earthquake source.
Nowadays, there are so many methods and approaches to get more information about the detailed rupture process. A very good paper review has been done by Professor Satoshi Ide [Link 1, Link 2] entitled “Slip Inversion” as a book chapter in Treatise on Geophysics [Link]. This paper gives us more detail explanation and development history of the methods/theories about slip inversion of earthquake.
There is also a very helpful book written by Prof. Agustin Udias, Prof. Raul Madariaga, and Prof. Elisa Buforn, entitled “Source Mechanisms of Earthquakes, Theory and Practices” [Link]. This book also gives us the more detail theories, methodology and some technical issues to provide the source mechanism of the earthquake.
You can also see some finite fault model database that I have compiled before. [Click Here].
Please prepare this programs first:
- Install Rdseed [Link]
- Install Seismic Analysis Code (SAC) [Link]
- Install Numerical Recipes [Link]. Numerical Recipes is available by open access to users in developing countries, including Indonesia. In this case, we need Numerical Recipes in Fortran.
- Kikuchi’s and Kanamori program from website maintained by ERI, Univ. of Tokyo [Link].
Step by step:
- Download package of Numerical Recipes
- Extract (easily using tar -xzvf)
- Go to directory of Numerical Recipes
- Enter directory “f”
- Delete a file named libnrf.a
- Open file named jacobi.f and ludcmp.f
- In those files, change value in line 4 (NMAX = 500) to be NMAX = 10000 (nmax must be set at a larger value than 1000)
- Close those files
- Open file Makefile
- In line 86 change $(FLIB) to only flib
- Close those files
- Open terminal, type: make flib
- A new libnrf.a is available now
- Copy this file (libnrf.a)
After we download the inversion program from the website. You can also download the sample files. Extract both of Program and Sample to directories.
- Go to directory Program
- Make a new directory named Recipes (easily using command: mkdir Recipes)
- Paste the libnrf.a but rename it to be recipes.a
- Next, go to directory Subroutine
- Open a file named Readme
- This is how to update our library file
- Copy all of the commands in the Readme file
- Open terminal, paste, and it will be run
Next go to directory Sample, cut the folder Input and Output and then paste in directory Program. Now, directory Program have 5 sub-directories i.e. Main, Recipes, Subroutine, Input and Output.
Now, we will to compile all of our fortran (main programs). Easy, there is a Makefile in directory Program. But, we need to change the Makefile to be executable. Open terminal (Ctrl+Alt+T), type sudo su, enter your Linux password, chmod 777 Makefile), close terminal.
Then you can easily compile the fortran program that you need. For example type: make green to get a executable program to have calculate green function using given velocity model.
In this stage, you can also modify the .bashrc file to write a path to this directory, thus you can use the program anywhere in your linux environment.
Preparing Our Seismic Data
- Select our case study (earthquake).
- We can download seismic data from some institutions, e.g. IRIS, GFZ, Japan’s network, BMKG network, etc. I have compiled some of data resources in this website [Link].
- In this tutorial, we use IRIS’s Wilber 3 facility.
- Easy, select your earthquake.
- Select stations with epicentral distance 30-90 degree.
- Why? Because for data from stations with distance less than 30 degree, the waves suffer strong perturbations by the structure of the upper mantle. For data from stations larger than 90 degree, the waves are diffracted because they reach the nucleus of the Earth.
- Download the data with some length for example 5 minute before P arrival and 60 minutes after P arrival or 3 minutes after S arrival. Select waveform with sampling each 0.05 seconds (20 Hz) or usually named as BH* channel. Seismic data usually configured to 5 levels of samples: HH a 100Hz; SH a 50Hz; BH a 20Hz; LH a 1Hz & VH a 0.1Hz [more information, click here: Link].
- Note: waveforms from IRIS wilber 3 have event informations inside its header (evla, evlo and evdp). If you use waveforms from other agencies (e.g. GFZ or BMKG), please check whether the event information have been put in the data header or not. If not yet, easy, you can use command: ch (change header) to fill the event informations. IRIS’s waveforms also have arrival time data on its header.
- Download the seismic data in format of seed (because we also need the RESP file) [Link]
- Make a new directory named Data, put the seed file inside it.
- Extract this rdseed data using command: rdseed -f filename.seed -R -d -o 1 -p [see this tutorial: Link]
- Then the folder contains a lot of files including SAC files, RESP files and pole-zeros files.
- Select waveforms with same (uniform) data, for example only select 00.BH* data. Remove all other data except data with 00.BH* (00.BHE, 00.BHN and 00.BHZ)
- Make a script named “mk_conv.farm.pl” according to the manual from the website.
- Change the value of high_pass_co, low_pass_co, Duration, pre_ev and sample (depend on your purposes).
- Don’t forget to provide files named “hypo” and “sacmacro“
- In line 61, change the coda of the SAC file, according to your data, for example, if we download the data from IRIS now, the coda will be “M.SAC” instead of “Q.SAC“.
- In the end of “sacmacro” file, add command “quit”
- Or you can also add some additional commands such as rmean, rtr and taper in sacmacro file (read the function of this command in manual of SAC).
- Run the perl program (make sure “perl” have been installed in your Linux)
- Use command: perl mk_conv.farm.pl
- Now two new files are available i.e. response and i_conv.farm
- Make sure there is no “error” in every step!
- Go to our main directory
- Copy some input data e.g. jb.ptime, jb.stime, jb.table to directory Data.
- In directory Data, run ./conv.sac.farm < i_conv.farm
- This will result in a new file named fort.22
- Then we will convert fort.22 to fort.1 (main input of inversion program)
- Use command: ./rotSH to produce fort.1
- Again, make sure there is no error in every step!
This is the main purpose of this tutorial, but if you have prepare your input data very well, you can continue easily to this step by reading the manual. There are there different sets of inversion program i.e.:
Set 1: Inversion allowing mechanism changes.
Set 2: Fine tuning of the source time function and slip distribution.
Set 3: Inversion with a fixed fault mechanism.
For each step, you can see the input, stdin and output in the manual which is available in the website. Use command: ./program name < stdin [general format]. [more information click here: Link]
To run the program:
- ./green < i_green
- ./inversion < i_inversion
- ./graphics < i_graphics
The parameters that you need to set in i_green, i_inversion and i_graphics depend on your cases and situations. You can read more descriptions in the manual.
For more description of theory and methodology, read these papers:
Kikuchi,M., & Kanamori,H., Bull. Seism. Soc. Am., 72, 491-506, 1982. [Link]
Kikuchi,M., & Kanamori,H., Phys. of the Earth and Plan. Int., 43, 205-222, 1986. [Link]
Kikuchi,M., & Kanamori,H., Bull. Seism. Soc. Am., 81, 2335-2350, 1991. [Link]
Kikuchi,M., Kanamori,H. & Satake,K., J. Geophys. Res., 98, 15797-15808, 1993. [Link]
If you have further questions, don’t hesitate to send me a message (email@example.com).
Special thanks to Fadiah Khairina and Renhard Sipayung who helped me to handle Numerical Recipes!