IGU-BD3C-5型地震仪SAC数据一致性检验

最新推荐文章于 2024-07-23 17:15:35 发布

李百柏

最新推荐文章于 2024-07-23 17:15:35 发布

阅读量133

点赞数

分类专栏：地震数据简单处理文章标签： python3.11

本文链接：https://blog.csdn.net/aoaoooo/article/details/132692337

版权

地震数据简单处理专栏收录该内容

2 篇文章 0 订阅

订阅专栏

IGU-BD3C-5型号地震仪是教学工作中常用的地震仪。在地震仪布线埋设之前，常常需要进行仪器响应的一致性检验。

现在给出以下Shell+SAC+Python3联合编写的代码，实现利用pz1文件结合SAC软件数据去除仪器响应、（手动选取）截取时间窗口计算相关系数、将仪器一致性数据导出成txt文件、绘制相关系数柱状图和平均相关系数柱状图：

#! /bin/bash
#! /usr/bin/python3
#Create by ljw in 2023/7/22, China University of Geoscience(Wuhan)
rm -rf *.func
rm -rf sac_temp
rm -rf sac/*E.sac
rm -rf sac/*N.sac
rm -rf result.txt
rm -rf average.txt
cp -R sac sac_temp
for file in sac/*.Z.sac
do
echo $file >> number.txt
for asc in sac_temp/*.Z.sac
do
sac << EOF
cut 2000 2600
read ${asc} ${file}
cor norm
cor master 2
w append .func
q
EOF
mv sac/*.func .
mv sac_temp/*.func .
done
saclst depmax depmax depmin f *.func | awk '{system("echo "$2" >> result.txt")}'
saclst depmax depmax depmin f *.func | awk '{system("echo "$2" >> average.txt")}'
rm -rf *.func
python3 << EOF
import matplotlib.pyplot as plt
import re
import os
from pylab import mpl
path = r'./sac'
filename= os.listdir(path)
x_temp = []
x = []
for name in filename:
        x_temp.append(int(name[5:9]))
x_temp.sort()
for x_value in x_temp:
        x.append(str(x_value))
mpl.rcParams['font.sans-serif'] = ['simhei']
mpl.rcParams['axes.unicode_minus'] = False
in_file = open("result.txt",mode='r')
txt = in_file.read()
txt_list = re.split('[\t \n]',txt.strip())
txt_out = map(float,txt_list)
y = list(txt_out)
with open("number.txt","r",encoding='utf-8') as f:
        z_temp = f.read()
        z = z_temp[9:13]
fig, ax = plt.subplots(figsize=(10, 7))
ax.bar(x=x, height=y)
for a,b in zip(x,y):
        ax.text(a,b,'%.3f'%b,ha = 'center',va = 'bottom',fontsize = 12)
        ax.set_title('与'+z+'号仪器进行一致性检验', fontsize=25)
plt.xlabel("仪器编号",fontsize=15)
plt.ylabel("相关系数",fontsize=15)
plt.ylim((0,1.1))
plt.savefig(z+'.png')
EOF
rm -rf number.txt
rm -rf result.txt
done
python3 << EOF
import matplotlib.pyplot as plt
import re
import os
import math
from pylab import mpl
path = r'./sac'
filename= os.listdir(path)
x_temp = []
x = []
for name in filename:

        x_temp.append(int(name[5:9]))
x_temp.sort()
for x_value in x_temp:

        x.append(str(x_value))
mpl.rcParams['font.sans-serif'] = ['simhei']
mpl.rcParams['axes.unicode_minus'] = False
in_file = open("average.txt",mode='r')
txt = in_file.read()
txt_list = re.split('[\t \n]',txt.strip())
txt_out = map(float,txt_list)
y_temp = list(txt_out)
length_all = len(y_temp)
length = int(math.sqrt(length_all))
y = [0 for num in range(0,length)]
i = 0
for value_1 in y_temp:

        j = int(i % length)

        y[j] = y[j] + value_1
        i = i + 1
k = 0
for value_2 in y:
        y[k] = value_2/length
        k = k + 1
fig, ax = plt.subplots(figsize=(10, 7))
ax.bar(x=x, height=y)
for a,b in zip(x,y):
        ax.text(a,b,'%.3f'%b,ha = 'center',va = 'bottom',fontsize = 12)
        ax.set_title("平均相关系数", fontsize=25)
plt.xlabel("仪器编号",fontsize=15)
plt.ylabel("相关系数",fontsize=15)
plt.ylim((0,1.1))
plt.savefig("平均相关系数.png")
EOF
rm -rf average.txt
rm -rf sac_temp