Travel time curve calculation for spherical, isotropic, homogeneous Earth model

A travel time curve is a graph of the time it takes a seismic ray to reach the station from an earthquake versus the distance of the earthquake to the station (epicentral distance).

Here, we calculate the travel time curve for a spherically homogeneous Earth. This means the seismic velocity is constant for the whole earth. We have considered the three cases,

  1. When the seismic ray travel through the body from the surface event to the station (body waves).
  2. When the seismic ray for a deeper event (300 km depth) reaches the station.

Fortran code for calculation of the Model:

program travel_time_calc
    !Program to calculate seismic travel times in a spherical, homogeneous Earth.
    implicit none
    real, parameter :: pi = 2*asin(1.0), R = 6371    ! radius of Earth (km
    real :: V, xx, h1, h, dx
    integer, parameter :: n = 89
    integer :: i
    real(kind=16), allocatable, dimension(:) :: x, xrad, Di, Ti, Tim, Dhi, Thi, Thim

    V = 12;        ! assumed constant seismic wave velocity (km/s)
    dx=2 !increment
     allocate(x(0:n),xrad(0:n),Di(0:n),Ti(0:n),Tim(0:n),Dhi(0:n),Thi(0:n),Thim(0:n)) !allocating the size of the arrays
    write(*,*) "Distance (in degrees)", "  ","T-time surface (in mins)", "  ", "T-time interior (in mins)"
    open(1,file="Tcalc.dat")
    xx=0
    do i=0,n !do loop
        x(i) = xx + dx    !distance
        xrad(i) = x(i)*pi/180 ! convert degrees to radians
        !Path along the body
        Di(i) = 2*R*sin(0.5*xrad(i))    ! distance of the surface event along the body in km
        Ti(i) = Di(i)/V               ! travel time in seconds
        Tim(i) = Ti(i)/60             ! travel time in minutes
        
        !For deeper events (300 km)
        h1=300
        h=R-h1
        Dhi(i) = sqrt(R**2 + h**2 - 2*R*h*cos(xrad(i)))    !distance of the interior events in km
        Thi(i) = Dhi(i)/V                ! travel time in seconds
        Thim(i)= Thi(i)/60                ! travel time in minutes

         write(1,10) x(i), Tim(i), Thim(i) !printing output on the terminal
         write(*,10) x(i), Tim(i), Thim(i) !printing output on the terminal
        xx = x(i)
     end do
    10 format(f7.2, "  ",f8.4, "  ",f8.4)
    deallocate(x, xrad, Di, Ti, Tim, Dhi, Thi, Thim) ! deallocates the array
end program travel_time_calc

For the above program, we have considered the Earth to be isotropic and homogeneous with constant velocity of 12 km/s. Then we parameterize for 89 events each with 2 degrees interval from 2 to 180 degrees. The radius of the earth is assumed to be 6371 km.

For compiling and plotting (bash script):

#!/bin/bash

gfortran -ffree-form travel_time_calc.f -o travel_time_calc
./travel_time_calc

gnuplot << EOF
set terminal postscript color solid #terminal type
set title 'Travel Time Curve' #figure title
set xlabel 'Delta (in degrees)' #xlabel
set ylabel 'Travel Time (in mins)' #ylabel
set output 'travel_time.ps' #output fig name
set grid #grid
plot [:] [:] "Tcalc.dat" using 1:2 title "Calculated body waves travel time (surface event)" with lines lc rgb "blue", \
    "Tcalc.dat" using 1:3 title "Calculated body wave travel time (300 km deep event)" with lines lc rgb "red"
EOF

travel_time

– Utpal Kumar (IES, Academia Sinica)

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s