本文实例为大家分享了python地震数据可视化的具体代码,供大家参考,具体内容如下
参考源码:
准备工作:
在windows10下安装python3.7,下载参考源码到本地。
1. demo绘图测试
demo绘图指令
cmd> python seisplot.py --demo
问题1)缺少依赖包
file "d:/desktop/python可视化/seisplot/seisplot.py", line 16, in
import yaml
modulenotfounderror: no module named 'yaml'
>pip3 install yaml
could not find a version that satisfies the requirement yaml (from versions: )
no matching distribution found for yaml
查看environment.yml:
channels:
- defaults
- conda-forge
dependencies:
- matplotlib
- numpy
- obspy
- pillow
- pyyaml
安装相关依赖包,如 obspy pillow pyyaml
cmd> pip3 install obspy
出现proxy error, 添加选项 --proxy server:port
安装ok
问题2)配置文件的编码异常
d:\desktop\python可视化\seisplot>python seisplot.py --demo
[91m
welcome to
┌─┐┌─┐┬┌─┐┌─┐┬ ┌─┐┌┬┐
└─┐├┤ │└─┐├─┘│ │ │ │
└─┘└─┘┴└─┘┴ ┴─┘└─┘ ┴
good luck[0m
traceback (most recent call last):
file "seisplot.py", line 421, in
cfg = yaml.load(f)
file "d:\programfiles\python36\lib\site-packages\yaml\__init__.py", line 70, in load
loader = loader(stream)
file "d:\programfiles\python36\lib\site-packages\yaml\loader.py", line 34, in __init__
reader.__init__(self, stream)
file "d:\programfiles\python36\lib\site-packages\yaml\reader.py", line 85, in __init__
self.determine_encoding()
file "d:\programfiles\python36\lib\site-packages\yaml\reader.py", line 124, in determine_encoding
self.update_raw()
file "d:\programfiles\python36\lib\site-packages\yaml\reader.py", line 178, in update_raw
data = self.stream.read(size)
unicodedecodeerror: 'gbk' codec can't decode byte 0x99 in position 2003: illegal multibyte sequence
分析源码中参数传递情况:
def load(stream, loader=loader):
"""
parse the first yaml document in a stream
and produce the corresponding python object.
"""
print(stream) ## debug by huanying03
loader = loader(stream)
try:
return loader.get_single_data()
finally:
loader.dispose()
打印结果如下:
编码不对?
65001 :utf-8代码页
936 :默认的gbk
437 :是美国英语
仔细检查发现config.yml中有异常字符,更正config.yml中错误字符:ok
问题3)python调试打印语句
print ("value=%d"%value)
demo数据绘图结果
如下:
2. 本地数据绘图测试
应用指令
cmd > python seisplot.py data/my.sgy
问题1)数据道的采样值异常
filename data2/t10.sgy
seismic.py:104:dt=4000
seismic.py:104:ns=512
seismic.py:144:header=b'c 1 client'
n_traces 500
n_samples 512
dt 0.004
t_start 0
t_end 2.044
max_val inf
min_val -inf
clip_val 180462501764003194804057887685476352
read data in 0.4 s
++++++++++++++++++++++++++++++++++++++++++++++++++++++++
plotting
plot width 17.00 in
plot height 11.00 in
traceback (most recent call last):
file "seisplot.py", line 442, in
main(t, cfg)
file "seisplot.py", line 214, in main
cfg)
file "d:\desktop\python可视化\seisplot\plotter.py", line 123, in plot_histogram
y, x, _ = ax.hist(np.ravel(data), bins=int(100.0 / (clip_val / largest)),
overflowerror: cannot convert float infinity to integer
解决办法:将数据做归一化处理
问题2)数据文件的卷头异常
filename data2/marmousi_vel.segy.hdrs.sgy.cdp.sgy
seismic.py:104:dt=10000
seismic.py:104:ns=350
seismic.py:144:header=b'\xc3@\xf1@\xc3\xd3\xc9\xc5\xd5\xe3'
traceback (most recent call last):
file "seisplot.py", line 442, in
main(t, cfg)
file "seisplot.py", line 40, in main
s = seismic.from_segy(target, params={'ndim': cfg['ndim']})
file "d:\desktop\python可视化\seisplot\seismic.py", line 164, in from_segy
return cls.from_obspy(stream, params=params)
file "d:\desktop\python可视化\seisplot\seismic.py", line 147, in from_obspy
x =np.array(list(stream.textual_file_header.decode()))
unicodedecodeerror: 'utf-8' codec can't decode byte 0xc3 in position 0: invalid continuation byte
解决办法:重写3200字节的segy卷头
marmousi相关数据绘图
速度模型绘图结果如下:
marmousi模型正演的炮集数据绘图结果如下:
测试备忘
1
cygwin执行指令,可以看到反馈信息为彩色,分析信息很方便。
相比之下,cmd终端执行指令后的反馈信息可读性较差。
2
本地数据出现异常后,通过cwp中的segyread, sushw, sugain, segywrite等程序进行修改,满足绘图程序要求。
1)3200字节的卷头必须规范
2)道头字中的cdpt取值必须规范(xline)
3)道采样值最好归一化处理。
以上就是本文的全部内容,希望对大家的学习有所帮助,也希望大家多多支持萬仟网。
希望与广大网友互动??
点此进行留言吧!