obspy去polezero类型仪器仪响应

  如果有polezero文件,那可以用attach_paz来读取。然后去仪器响应的脚本如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
import numpy as np
import obspy
import matplotlib.pyplot as plt
from obspy import read
from obspy.io.sac import attach_paz
from obspy.signal.invsim import corn_freq_2_paz
from pathlib import Path
fpath="/home/junxie/work/Snet/data/2016_08_16_01_00"
p=Path(fpath)
i=0
for file in p.rglob('*VX*.SAC'): #循环读取文件
st=read(file, debug_headers=True)
tr=st[0].copy()
aa=str(file).split("/")
pz=str(file)+"_PZ" # 这就是paz文件
attach_paz(tr,pz) # 贴到tr中
paz_1hz = corn_freq_2_paz(1.0, damp=0.707) # 1Hz instrument
paz_1hz['sensitivity'] = 1.0
st.simulate(paz_simulate=paz_1hz) #去仪器响应

  如果没有的话找到零点、极点和放大系数,赋值到paz里面,例如:

1
2
3
4
5
6
7
8
9
10
11
from obspy import read
from obspy.signal.invsim import corn_freq_2_paz
st = read()
paz_sts2 = {'poles': [-0.037004+0.037016j, -0.037004-0.037016j,
-251.33+0j,
-131.04-467.29j, -131.04+467.29j],
'zeros': [0j, 0j],
'gain': 60077000.0,
'sensitivity': 2516778400.0}
paz_1hz = corn_freq_2_paz(1.0, damp=0.707)
st.simulate(paz_remove=paz_sts2, paz_simulate=paz_1hz)