USGS Catalog in Combination with Global CMT

从 USGS 地震目录 (https://earthquake.usgs.gov/earthquakes/search/) 下载的目录和 Global CMT 上下载的目录 (https://www.globalcmt.org) 对于震级较大的地震一般不一致。USGS 提供地震的发震时刻 (Origin time)，震中位置等，而 GCMT 提供地震的质心时刻（Centroid time）和质心位置，以及震源机制类型和标量矩等信息。本篇博客主要分享两个 shell 脚本，用来批量地将 USGS 和 GCMT 目录相匹配，并且在 GCMT 目录中判断断层机制类型。

USGS & GCMT

分别在 USGS 和 GCMT 网站下载地震目录: query.csv 和 GCMT.txt (GMT psvelomeca input 输出再手动创建文本文件)。由于二者的地震目录在震级上存在差异，在初始选择地震目录时，建议将 GCMT 的震级范围减两级；比如要匹配 7 级及以上地震时，在 GCMT 上下载 6.8 级及以上地震。

# USGS - query.csv format
time,latitude,longitude,depth,mag, ...
2020-07-22T06:12:44.719Z,55.0683,-158.5543,28,7.8,...

# GCMT - GCMT.txt format
lon lat str1 dip1 rake1 str2 dip2 rake2 sc iexp name
-159.27 54.84 242 17 90 62 73 90 7.402 27 202007220612A

根据上述数据的共性，选择时间和位置作为判别条件：如果两个事件的时间在 5 分钟以内，且位置差异在 5 deg 以内，则认为是同一事件：

#!/bin/csh

set pi = 3.1415926
rm result
sort -n -k 1 query.csv > tmp; mv tmp query.csv
set num = `wc -l query.csv | awk '{print $1}'`

@ i = 1
@ j = 0
while ($i < $num)
  @ i++
  set time_1 = `awk -F ',|T|Z|-|:' '{if (NR == '$i') print $1$2$3$4$5}' query.csv`
  @ k = -1
  while ($k < 0)
    @ j++

### Distance    
    set lon1 = `awk -F ',' '{if (NR == '$i') print $3}' query.csv`
    set lat1 = `awk -F ',' '{if (NR == '$i') print $2}' query.csv`
    set lon2 = `awk '{if (NR == '$j') print $1}' GCMT.txt`
    set lat2 = `awk '{if (NR == '$j') print $2}' GCMT.txt`
    set dtmp1 = `echo 'scale=5;c ( ( 90 - '$lat1' ) * '$pi' / 180 ) * c ((90 - '$lat2') * '$pi' / 180 ) ' | bc -l`         
    set dtmp2 = `echo 'scale=5;s ((90 - '$lat1' ) * '$pi' / 180) * s ((90 - '$lat2' ) * '$pi' / 180 )' | bc -l`
    set dtmp3 = `echo 'scale=5; '$dtmp2' * c (('$lon2' - '$lon1') * '$pi' / 180 )' | bc -l`    
    set tmp = `echo 'scale=0; ('$dtmp1' + '$dtmp3' + 1) * 1000 ' | bc -l`  
    if ($tmp == 0 ) then
      set dist_deg = `echo 'scale=8; (180 / '$pi') * '$pi'' | bc -l`
    else         
      set dist_deg = `echo 'scale=8; (180 / '$pi') * 2 * a ( sqrt (1 - ('$dtmp1' + '$dtmp3')^2 ) / (1 + '$dtmp1' + '$dtmp3' ))' | bc -l`
    endif
    set ttmp = `echo $dist_deg | awk '{if ($1 < 5) print 1; else print -1}'`    
    
### Time      
    set time_2 = `awk '{if (NR == '$j') print $11}' GCMT.txt | awk -F '[A|B]' '{print $1}'`
    set tmp = `echo $time_1 - $time_2 | bc -l`

    if ($tmp < 0) then
      @ k = 0; @ j--
    else if (($tmp == 0) && ($ttmp > 0)) then  
      echo $time_1 $time_2 $dist_deg >>  result
    else if (($tmp < 5) && ($ttmp > 0)) then
      echo $time_1 $time_2 $dist_deg "///" >>  result
    endif
  end   
end

Fault types

由 GCMT 地震目录中的滑动角，判断震源机制类型：

We parameterize the focal mechanism type using the rakes, r1 and r2 (degrees), of the two nodal planes using the following algorithm (J. Hardebeck, personal communication, 2005):if (abs(r1) >90) r1 = (180‐abs(r1))(r1/abs(r1)) if (abs(r2) >90) r2 = (180‐abs(r2))(r2/abs(r2)) if (abs(r1) < abs(r2)) then r = r1 else r = r2 end if fptype = r/90. The parameter fptype will vary from −1 (normal) to 0 (strike‐slip) to 1 (reverse) and has the advantage of providing a single scalar value for characterizing the faulting type (Shearer et al., JGR, 2006).

#!/bin/csh

rm GCMT_fault_types.txt

set num = `wc -l GCMT.txt | awk '{print $1}'`

@ i = 0

while ($i < $num)
  @ i++
  set slip1=`awk '{if (NR == '$i') print $5}' GCMT.txt`
  set slip2=`awk '{if (NR == '$i') print $8}' GCMT.txt`

  set slip1_abs=`echo $slip1 | awk '{if ($1 >= 0) print $1;else print -1*$1}'`
  set slip2_abs=`echo $slip2 | awk '{if ($1 >= 0) print $1;else print -1*$1}'`

  if ($slip1_abs > 90) then
    set slip1 = `echo "(180- $slip1_abs )* $slip1 / $slip1_abs" | bc -l`
  endif
  
  if ($slip2_abs > 90) then
    set slip2 = `echo "(180- $slip2_abs )* $slip2 / $slip2_abs" | bc -l`
  endif
  
  set slip1_abs=`echo $slip1 | awk '{if ($1 >= 0) print $1;else print -1*$1}'`
  set slip2_abs=`echo $slip2 | awk '{if ($1 >= 0) print $1;else print -1*$1}'`

  set tmp = `echo $slip1_abs $slip2_abs | awk '{if ($1<$2) print 1; else print -1}'`
  if ($tmp > 0) then
    set slip = `echo $slip1/90 | bc -l`
  else
    set slip = `echo $slip2/90 | bc -l`
  endif

  awk '{if (NR == '$i') print $0,'$slip'}' GCMT.txt >> GCMT_fault_types.txt
end

个人统计了 2008/01 至 2020/10 7 级及以上地震的断层类型，横坐标数字表示含义为： −1 (normal), 0 (strike‐slip) 和 1 (reverse).

# GMT plotting   
awk '{print $12}' GCMT_fault_types.txt  | gmt histogram -R-1/1/0/100 -JX8c \
-W0.2 -Ggray -Bxa0.5 -Bya20+l"Counts" -BWSne -D+f10p,4+o8p -L1p -png histVert

Contents

USGS & GCMT

Fault types