Created
May 10, 2024 09:44
-
-
Save aadm/8ebe73b5db97a5a711ede51e35dbd8ec to your computer and use it in GitHub Desktop.
Simple interactive viewer for seismic data built with python+streamlit+segyio
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
# seismic_viewer_app_demo | |
#-------------------------------------------------------------------- | |
# | |
# testing streamlit for a seismic viewer app | |
# using Near stack subcube from Per Avseth's QSI dataset | |
# | |
# to run: | |
# $ streamlit run seismic_viewer_app_demo.py | |
#-------------------------------------------------------------------- | |
import numpy as np | |
import matplotlib.pyplot as plt | |
import streamlit as st | |
import os | |
import requests | |
from tempfile import NamedTemporaryFile | |
import segyio | |
@st.cache_data | |
def st_load_seismic(): | |
remote_file = 'https://github.com/aadm/geophysical_notes/raw/master/3d_nearstack.sgy' | |
response = requests.get(remote_file) | |
# raise an HTTPError if the download fails | |
if response.raise_for_status() is False: | |
st.print('error while downloading {}'.format(remote_file)) | |
return None, None, None, None | |
# create a temporary file to save the downloaded SEGY file | |
with NamedTemporaryFile(delete=False) as tmp: | |
tmp.write(response.content) | |
local_file = tmp.name | |
with segyio.open(local_file, iline=41, xline=21) as f: | |
ntraces = f.tracecount | |
n_inl = f.ilines.size | |
n_crl = f.xlines.size | |
sr = segyio.tools.dt(f)/1e3 | |
nsamples = f.samples.size | |
twt = f.samples | |
data = segyio.tools.cube(f) | |
os.remove(local_file) # remove local files | |
seis = data.T | |
z, crl, inl = seis.shape | |
return seis, nsamples, n_crl, n_inl | |
@st.cache_data | |
def st_get_section(data, line_type, line): | |
if line_type == 'inline': | |
section = data[:,:,line] | |
else: | |
section = data[:,line,:] | |
return section | |
@st.cache_data | |
def st_plot_section(data, line_type, line, z, ztop=None, zbot=None, | |
cmap='bwr_r', clip=None, interpolation='spline16', | |
add_colorbar=True, square=True): | |
''' modified from squit.plot.model_seis() ''' | |
n_traces = data.shape[1] | |
if ztop is None: ztop = z.min() | |
if zbot is None: zbot = z.max() | |
# calculate upper/lower clips value for color palette | |
if clip is None: | |
clip = np.nanpercentile(data, 100) | |
c1, c2 = -clip, +clip | |
elif isinstance(clip, list) or isinstance(clip, tuple): | |
c1, c2 = clip | |
else: | |
c1, c2 = -clip, +clip | |
# make plot | |
fig, ax = plt.subplots(ncols=1, constrained_layout=True) | |
extents = [0, n_traces, z[-1], z[0]] | |
opt = dict(extent=extents, cmap=cmap, | |
vmin=c1, vmax=c2, | |
aspect='auto', interpolation=interpolation) | |
im = ax.imshow(data, **opt) | |
if add_colorbar: | |
cbar = plt.colorbar(im, ax=ax, orientation='horizontal') | |
cbar.ax.tick_params(labelsize='x-small') | |
title = '{} {}'.format(line_type, line) | |
ax.set_title(title, fontsize='large') | |
ax.tick_params(axis='both', labelsize='x-small') | |
ax.set_ylim(ztop, zbot) | |
ax.invert_yaxis() | |
if square: | |
ax.set_aspect('equal', 'box') | |
opt = dict(ztop=ztop, zbot=zbot, cmap=colormap, clip=clip) | |
st.pyplot(fig=fig, use_container_width=True) | |
st.set_page_config(page_title='Simple Seismic Viewer', layout="centered") | |
colormap_list = [ | |
'RdBu', 'RdGy', 'bwr_r', 'coolwarm_r', 'seismic_r', | |
'gray_r' or 'gist_gray_r' | |
] | |
#--- BEGIN APP | |
st.header('Simple Seismic Viewer') | |
seis, z, crl, inl = st_load_seismic() | |
widg_load = st.columns(2) | |
with widg_load[0]: | |
line_type = st.selectbox('inline or crossline', ('inline', 'crossline')) | |
with widg_load[1]: | |
min_line = 0 | |
if line_type == 'inline': | |
max_line = inl | |
else: | |
max_line = crl | |
opt_line = dict(min_value=min_line, max_value=max_line, format='%i') | |
line = st.slider('{} #'.format(line_type), **opt_line) | |
section = st_get_section(seis, line_type, line) | |
z_vect = np.arange(z) | |
clip000, clip100 = np.percentile(np.abs(section), (1, 100)) | |
clip099 = np.percentile(np.abs(section), 99) | |
widg_a = st.columns(3) | |
with widg_a[0]: | |
colormap = st.selectbox('colormap', colormap_list, label_visibility="collapsed",placeholder='colormap') | |
with widg_a[1]: | |
opt_clip = dict(min_value=clip000, max_value=clip100, value=clip099) | |
clip = st.number_input('clip', **opt_clip, label_visibility="collapsed") | |
with widg_a[2]: | |
square = st.checkbox('1:1 aspect ratio') | |
widg_b = st.columns(2) | |
opt_z = dict(min_value=z_vect.min(), max_value=z_vect.max(), step=10, format='%i') | |
with widg_b[0]: | |
ztop = st.slider('top window', value=z_vect.min(), **opt_z) | |
with widg_b[1]: | |
zbot = st.slider('bot window', value=z_vect.max(), **opt_z) | |
if zbot <= ztop: | |
zbot = ztop+100 | |
opt = dict(ztop=ztop, zbot=zbot, cmap=colormap, clip=clip, square=square) | |
st_plot_section(section, line_type, line, z_vect, **opt) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment