Get Started in 10 Minutes¶
Altay Sansal
Sep 22, 2026
8 min read
In this page we will be showing basic capabilities of MDIO.
For demonstration purposes, we will ingest a public raw shot record from the Soda Lake 2010 3D/3C seismic survey (Geothermal Data Repository / AWS Open Data, CC-BY-4.0).
We use one ~7.5 MiB field record: F7733R1.SGY. No AWS credentials are required.
The file is a 3C common-source gather: shot 7733, 958 traces (1 aux + 3 × 319
receivers), 2000 samples at 2 ms.
Warning
For plotting and remote ingestion the notebook requires Matplotlib and aiohttp
as a dependency. Please install it before executing using pip install matplotlib aiohttp or
conda install matplotlib aiohttp.
Defining the SEG-Y Dataset¶
Since MDIO 0.8 we can directly ingest remote SEG-Y files. To make the header scan performant we can also set up an environment variable for MDIO. See here for when to use this: Buffered Reads in Ingestion.
The gather is a single live shot, so the live mask is all True. The first trace is an
auxiliary channel at (0, 0) and is much louder than the seismic. Plots below drop it
and keep one of the three components.
The following environment variables are essential here:
MDIO__IMPORT__CLOUD_NATIVEtells MDIO to do buffered reads for headers due to remote file.MDIO__IMPORT__SAVE_SEGY_FILE_HEADERmakes MDIO save the SEG-Y specific file headers (text, binary) which is not strictly necessary for consumption and is disabled by default.
import os
os.environ["MDIO__IMPORT__CLOUD_NATIVE"] = "true"
os.environ["MDIO__IMPORT__SAVE_SEGY_FILE_HEADER"] = "true"
input_url = "https://gdr-data-lake.s3.us-west-2.amazonaws.com/soda_lake/raw_seismic/2010/v1.0.0/F7733R1.SGY"
Ingesting to MDIO¶
To do this, we can use the convenient SEG-Y to MDIO converter.
The shot identifier is the field record number at byte 9, and the channel number is at byte 13.
MDIO uses TGSAI/segy to parse the SEG-Y; the field names conform to its canonical keys defined in SEGY Binary Header Keys and SEGY Trace Header Keys. Since MDIO v1 we also introduced templates for common seismic data types. For this shot gather we use the StreamerShotGathers2D template here, which expects shot_point and channel.
We will also specify the units for the time domain. The spatial units will be automatically parsed from SEG-Y binary header. However, there may be a case where it is corrupt in the file, for that see the Fixing X/Y Units Issues section.
In summary, we will use the byte locations as defined for ingestion.
import matplotlib.pyplot as plt
from segy.schema import HeaderField
from segy.standards import get_segy_standard
from mdio import segy_to_mdio
from mdio.builder.schemas.v1.units import TimeUnitModel
from mdio.builder.template_registry import get_template
shot_trace_headers = [
HeaderField(name="shot_point", byte=9, format="int32"),
HeaderField(name="channel", byte=13, format="int32"),
]
rev1_segy_spec = get_segy_standard(1.0)
soda_lake_segy_spec = rev1_segy_spec.customize(trace_header_fields=shot_trace_headers)
mdio_template = get_template("StreamerShotGathers2D")
unit_ms = TimeUnitModel(time="ms")
mdio_template.add_units({"time": unit_ms})
segy_to_mdio(
input_path=input_url,
output_path="soda_lake.mdio",
segy_spec=soda_lake_segy_spec,
mdio_template=mdio_template,
overwrite=True,
)
It only took a few seconds to ingest, since this is a very small file.
However, MDIO scales up to TB (that’s ~1,000 GB) sized volumes!
Opening the Ingested MDIO File¶
Let’s open the MDIO file with the open_mdio function. This will return a pretty xarray representation with our standardized format.
from mdio import open_mdio
dataset = open_mdio("soda_lake.mdio")
dataset
<xarray.Dataset> Size: 8MB
Dimensions: (shot_point: 1, channel: 958, time: 2000)
Coordinates:
* shot_point (shot_point) int32 4B 7733
source_coord_x (shot_point) float64 8B ...
source_coord_y (shot_point) float64 8B ...
* channel (channel) int32 4kB 1 2 3 4 5 6 ... 954 955 956 957 958
group_coord_y (shot_point, channel) float64 8kB ...
group_coord_x (shot_point, channel) float64 8kB ...
* time (time) int32 8kB 0 2 4 6 8 10 ... 3990 3992 3994 3996 3998
Data variables:
segy_file_header <U1 4B ...
trace_mask (shot_point, channel) bool 958B ...
headers (shot_point, channel) [('trace_seq_num_line', '<i4'), ('trace_seq_num_reel', '<i4'), ('energy_source_point_num', '<i4'), ('ensemble_num', '<i4'), ('trace_num_ensemble', '<i4'), ('trace_id_code', '<i2'), ('vertically_summed_traces', '<i2'), ('horizontally_stacked_traces', '<i2'), ('data_use', '<i2'), ('source_to_receiver_distance', '<i4'), ('receiver_group_elevation', '<i4'), ('source_surface_elevation', '<i4'), ('source_depth_below_surface', '<i4'), ('receiver_datum_elevation', '<i4'), ('source_datum_elevation', '<i4'), ('source_water_depth', '<i4'), ('receiver_water_depth', '<i4'), ('elevation_depth_scalar', '<i2'), ('coordinate_scalar', '<i2'), ('source_coord_x', '<i4'), ('source_coord_y', '<i4'), ('group_coord_x', '<i4'), ('group_coord_y', '<i4'), ('coordinate_unit', '<i2'), ('weathering_velocity', '<i2'), ('subweathering_velocity', '<i2'), ('source_uphole_time', '<i2'), ('group_uphole_time', '<i2'), ('source_static_correction', '<i2'), ('receiver_static_correction', '<i2'), ('total_static_applied', '<i2'), ('lag_time_a', '<i2'), ('lag_time_b', '<i2'), ('delay_recording_time', '<i2'), ('mute_time_start', '<i2'), ('mute_time_end', '<i2'), ('samples_per_trace', '<i2'), ('sample_interval', '<i2'), ('instrument_gain_type', '<i2'), ('instrument_gain_const', '<i2'), ('instrument_gain_initial', '<i2'), ('correlated_data', '<i2'), ('sweep_freq_start', '<i2'), ('sweep_freq_end', '<i2'), ('sweep_length', '<i2'), ('sweep_type', '<i2'), ('sweep_taper_start', '<i2'), ('sweep_taper_end', '<i2'), ('taper_type', '<i2'), ('alias_filter_freq', '<i2'), ('alias_filter_slope', '<i2'), ('notch_filter_freq', '<i2'), ('notch_filter_slope', '<i2'), ('low_cut_freq', '<i2'), ('high_cut_freq', '<i2'), ('low_cut_slope', '<i2'), ('high_cut_slope', '<i2'), ('year_recorded', '<i2'), ('day_of_year', '<i2'), ('hour_of_day', '<i2'), ('minute_of_hour', '<i2'), ('second_of_minute', '<i2'), ('time_basis_code', '<i2'), ('trace_weighting_factor', '<i2'), ('group_num_roll_switch', '<i2'), ('group_num_first_trace', '<i2'), ('group_num_last_trace', '<i2'), ('gap_size', '<i2'), ('taper_overtravel', '<i2'), ('cdp_x', '<i4'), ('cdp_y', '<i4'), ('inline', '<i4'), ('crossline', '<i4'), ('shot_point', '<i4'), ('shot_point_scalar', '<i2'), ('trace_value_unit', '<i2'), ('transduction_const_mantissa', '<i4'), ('transduction_const_exponent', '<i2'), ('transduction_unit', '<i2'), ('device_trace_id', '<i2'), ('times_scalar', '<i2'), ('source_type_orientation', '<i2'), ('source_energy_dir_mantissa', '<i4'), ('source_energy_dir_exponent', '<i2'), ('source_measurement_mantissa', '<i4'), ('source_measurement_exponent', '<i2'), ('source_measurement_unit', '<i2'), ('channel', '<i4')] 218kB ...
amplitude (shot_point, channel, time) float32 8MB ...
Attributes:
apiVersion: 1.2.3
createdOn: 2026-09-22 21:03:51.528919+00:00
name: StreamerShotGathers2D
attributes: {'surveyType': '2D', 'gatherType': 'common_source', 'default...Querying Metadata¶
Now let’s look at the file text header saved in the segy_file_header metadata variable.
print(dataset["segy_file_header"].attrs["textHeader"])
C 1 CLIENT COMPANY CREW NO
C 2 LINE AREA MAP ID
C 3 REEL NO DAY-START OF REEL YEAR OBSERVER
C 4 INSTRUMENT: MFG MODEL SERIAL NO
C 5 DATA TRACES/RECORD AUXILIARY TRACES/RECORD CDP FOLD
C 6 SAMPLE INTERVAL SAMPLES/TRACE BITS/IN BYTES/SAMPLE
C 7 RECORDING FORMAT FORMAT THIS REEL MEASUREMENT SYSTEM
C 8 SAMPLE CODE: FLOATING PT FIXED PT FIXED PT-GAIN CORRELATED
C 9 GAIN TYPE: FIXED BINARY FLOATING POINT OTHER
C10 FILTERS: ALIAS HZ NOTCH HZ BAND - HZ SLOPE - DB/OCT
C11 SOURCE: TYPE NUMBER/POINT POINT INTERVAL
C12 PATTERN: LENGTH WIDTH
C13 SWEEP: START HZ END HZ LENGTH MS CHANNEL NO TYPE
C14 TAPER: START LENGTH MS END LENGTH MS TYPE
C15 SPREAD: OFFSET MAX DISTANCE GROUP INTERVAL
C16 GEOPHONES: PER GROUP SPACING FREQUENCY MFG MODEL
C17 PATTERN: LENGTH WIDTH
C18 TRACES SORTED BY: RECORD CDP OTHER
C19 AMPLITUDE RECOVERY: NONE SPHERICAL DIV AGC OTHER
C20 MAP PROJECTION ZONE ID COORDINATE UNITS
C21 PROCESSING:
C22 PROCESSING:
C23
C24
C25
C26
C27
C28
C29
C30
C31
C32
C33
C34
C35
C36
C37
C38
C39
C40 END EBCDIC
Since we saved the binary header, we can look at that as well.
dataset["segy_file_header"].attrs["binaryHeader"]
{'job_id': 18909,
'line_num': 1,
'reel_num': 1,
'data_traces_per_ensemble': 957,
'aux_traces_per_ensemble': 1,
'sample_interval': 2000,
'orig_sample_interval': 2000,
'samples_per_trace': 2000,
'orig_samples_per_trace': 2000,
'data_sample_format': 1,
'ensemble_fold': 0,
'trace_sorting_code': 1,
'vertical_sum_code': 1,
'sweep_freq_start': 0,
'sweep_freq_end': 0,
'sweep_length': 0,
'sweep_type_code': 4,
'sweep_trace_num': 0,
'sweep_taper_start': 0,
'sweep_taper_end': 0,
'taper_type_code': 3,
'correlated_data_code': 1,
'binary_gain_code': 2,
'amp_recovery_code': 1,
'measurement_system_code': 2,
'impulse_polarity_code': 0,
'vibratory_polarity_code': 0,
'fixed_length_trace_flag': 1,
'num_extended_text_headers': 0,
'segy_revision_major': 0,
'segy_revision_minor': 0}
Fetching Data and Plotting¶
We keep one component (transduction_const_mantissa == 2), clip by the 98th percentile,
and plot against source-receiver offset. Channel order follows the 2D receiver patch,
not offset.
Since MDIO v1 we are using Xarray under the hood, so we can use its convenient indexing. It also handles the plotting with proper dimension coordinate labels.
import numpy as np
COMPONENT_A = 2
shot = dataset.sel(shot_point=7733)
mantissa = shot.headers.to_numpy()["transduction_const_mantissa"]
comp = shot.isel(channel=mantissa == COMPONENT_A)
offset = np.hypot(
comp.group_coord_x.to_numpy() - float(comp.source_coord_x),
comp.group_coord_y.to_numpy() - float(comp.source_coord_y),
)
comp = comp.assign_coords(offset=("channel", offset)).sortby("offset")
amp = comp["amplitude"]
clip = float(np.percentile(np.abs(amp.to_numpy()), 98))
amp.plot(x="offset", y="time", vmin=-clip, vmax=clip, cmap="gray_r", yincrease=False);
Receiver coordinates in map view, then the same points colored by amplitude at 1,000 ms.
fig, ax = plt.subplots(figsize=(6, 5))
ax.scatter(comp.group_coord_x, comp.group_coord_y, s=16, c="0.3", label="receivers")
ax.scatter(
float(comp.source_coord_x),
float(comp.source_coord_y),
marker="*",
s=180,
c="tab:red",
label="source",
zorder=3,
)
ax.set_xlabel("group_coord_x")
ax.set_ylabel("group_coord_y")
ax.set_aspect("equal", adjustable="box")
ax.legend()
ax.set_title("component A receiver patch");
Color the same map by amplitude at 1,000 ms.
twt = amp.sel(time=1000)
fig, ax = plt.subplots(figsize=(6, 5))
sc = ax.scatter(
comp.group_coord_x,
comp.group_coord_y,
c=twt,
s=16,
cmap="gray_r",
vmin=-clip,
vmax=clip,
)
ax.scatter(
float(comp.source_coord_x),
float(comp.source_coord_y),
marker="*",
s=180,
c="tab:red",
label="source",
zorder=3,
)
ax.set_xlabel("group_coord_x")
ax.set_ylabel("group_coord_y")
ax.set_aspect("equal", adjustable="box")
ax.legend()
plt.colorbar(sc, ax=ax, label="amplitude @ 1000 ms");
Query Headers¶
We can query headers for the whole dataset very quickly because they are separated from the seismic wavefield.
Let’s get all the headers for source and receiver coordinates in this dataset.
Note that the header maps will still share the geometry/grid of the dataset!
The compute property fetches the lazily opened MDIO values.
dataset.headers["group_coord_x"].compute()
<xarray.DataArray 'group_coord_x' (shot_point: 1, channel: 958)> Size: 8kB
array([[ 0., 7734489., 7735044., 7735667., 7736158., 7736714.,
7737268., 7737825., 7738382., 7738938., 7739494., 7741720.,
7742273., 7742831., 7743326., 7743945., 7744504., 7745053.,
7745611., 7746169., 7746725., 7735425., 7735981., 7736537.,
7737093., 7737650., 7738204., 7738762., 7739319., 7739874.,
7740432., 7740987., 7741544., 7742098., 7742657., 7743212.,
7743767., 7744323., 7744879., 7745435., 7745947., 7746549.,
7747660., 7736361., 7736909., 7737472., 7738028., 7738587.,
7739142., 7739699., 7740255., 7740810., 7741366., 7741925.,
7742480., 7743037., 7743591., 7744147., 7744704., 7745258.,
7745818., 7746372., 7746875., 7737298., 7737856., 7738411.,
7738969., 7739523., 7740086., 7740636., 7741192., 7741747.,
7742321., 7742858., 7743416., 7743973., 7744529., 7745083.,
7745638., 7746197., 7746717., 7748979., 7749532., 7738232.,
7738789., 7739347., 7739904., 7740457., 7741018., 7741565.,
7742129., 7742683., 7743236., 7743790., 7744350., 7744908.,
7745463., 7746023., 7746577., 7747132., 7747690., 7748244.,
7748800., 7749360., 7749914., 7750469., 7739172., 7739720.,
7740280., 7740842., 7741404., 7741955., 7742505., 7743068.,
7743616., 7744177., 7744732., 7745294., 7745847., 7746406.,
...
7745143., 7745699., 7746254., 7746808., 7747367., 7747923.,
7748480., 7749035., 7749591., 7750145., 7750704., 7751260.,
7751817., 7752373., 7752928., 7753486., 7754042., 7754596.,
7755154., 7743853., 7744409., 7744965., 7745521., 7746078.,
7746633., 7747192., 7747747., 7748303., 7748859., 7749416.,
7749972., 7750529., 7751083., 7751639., 7752199., 7752752.,
7753307., 7753864., 7754421., 7754977., 7755534., 7756089.,
7744789., 7745346., 7745903., 7746459., 7747015., 7747571.,
7748128., 7748684., 7749237., 7749796., 7750351., 7750908.,
7751464., 7752021., 7752577., 7753135., 7753688., 7754246.,
7754802., 7755357., 7755914., 7756470., 7757026., 7745728.,
7746283., 7746840., 7747394., 7747951., 7748506., 7750175.,
7750731., 7751289., 7751845., 7752402., 7752958., 7753514.,
7754070., 7754625., 7755170., 7755736., 7756295., 7756852.,
7757405., 7757961., 7747221., 7747775., 7748332., 7748888.,
7752779., 7753338., 7753893., 7754450., 7755005., 7755564.,
7756119., 7756676., 7757229., 7757787., 7758343., 7758899.,
7748711., 7749269., 7749826., 7753161., 7753717., 7754272.,
7754831., 7755385., 7755940., 7756501., 7757053., 7757613.,
7758169., 7758724., 7759277., 7759836.]])
Coordinates:
* shot_point (shot_point) int32 4B 7733
source_coord_x (shot_point) float64 8B 7.735e+06
source_coord_y (shot_point) float64 8B 4.534e+07
* channel (channel) int32 4kB 1 2 3 4 5 6 ... 953 954 955 956 957 958
group_coord_y (shot_point, channel) float64 8kB 0.0 ... 4.535e+07
group_coord_x (shot_point, channel) float64 8kB 0.0 7.734e+06 ... 7.76e+06
Attributes:
unitsV1: {'length': 'ft'}dataset.headers["group_coord_y"].compute()
<xarray.DataArray 'group_coord_y' (shot_point: 1, channel: 958)> Size: 8kB
array([[ 0., 45338731., 45338358., 45337937., 45337606., 45337232.,
45336858., 45336483., 45336108., 45335735., 45335358., 45333861.,
45333488., 45333112., 45332779., 45332367., 45332013., 45331615.,
45331238., 45330866., 45330488., 45340121., 45339748., 45339376.,
45338997., 45338622., 45338248., 45337872., 45337501., 45337122.,
45336752., 45336440., 45336001., 45335627., 45335252., 45334878.,
45334503., 45334128., 45333754., 45333380., 45332941., 45332630.,
45331881., 45341515., 45341140., 45340762., 45340388., 45340016.,
45339638., 45339266., 45338892., 45338516., 45338142., 45337767.,
45337393., 45337019., 45336643., 45336269., 45335893., 45335518.,
45335144., 45334769., 45334431., 45342905., 45342531., 45342145.,
45341778., 45341405., 45341034., 45340657., 45340279., 45339906.,
45339520., 45339156., 45338783., 45338406., 45338034., 45337660.,
45337284., 45336907., 45336559., 45335038., 45334661., 45344294.,
45343919., 45343547., 45343172., 45342797., 45342423., 45342048.,
45341668., 45341298., 45340919., 45340557., 45340165., 45339796.,
45339422., 45339048., 45338672., 45338298., 45337925., 45337549.,
45337174., 45336802., 45336426., 45336052., 45345684., 45345307.,
45344935., 45344564., 45344186., 45343811., 45343437., 45343065.,
45342690., 45342315., 45341943., 45341561., 45341191., 45340814.,
...
45349747., 45349373., 45348997., 45348623., 45348247., 45347873.,
45347499., 45347124., 45346750., 45346375., 45346001., 45345628.,
45345251., 45344878., 45344501., 45344127., 45343754., 45343379.,
45343004., 45352638., 45352261., 45351886., 45351512., 45351138.,
45350761., 45350389., 45350015., 45349640., 45349265., 45348890.,
45348515., 45348142., 45347766., 45347392., 45347016., 45346642.,
45346268., 45345892., 45345519., 45345145., 45344769., 45344393.,
45354027., 45353652., 45353276., 45352902., 45352527., 45352152.,
45351779., 45351403., 45351030., 45350656., 45350280., 45349907.,
45349531., 45349158., 45348781., 45348409., 45348032., 45347657.,
45347284., 45346909., 45346533., 45346158., 45345785., 45355416.,
45355043., 45354668., 45354294., 45353918., 45353543., 45352417.,
45352045., 45351670., 45351295., 45350923., 45350547., 45350171.,
45349797., 45349424., 45349051., 45348675., 45348301., 45347925.,
45347550., 45347174., 45356434., 45356059., 45355685., 45355310.,
45352687., 45352312., 45351937., 45351564., 45351189., 45350814.,
45350440., 45350064., 45349691., 45349315., 45348941., 45348565.,
45357449., 45357074., 45356700., 45354450., 45354078., 45353704.,
45353328., 45352954., 45352578., 45352202., 45351828., 45351458.,
45351080., 45350706., 45350331., 45349955.]])
Coordinates:
* shot_point (shot_point) int32 4B 7733
source_coord_x (shot_point) float64 8B 7.735e+06
source_coord_y (shot_point) float64 8B 4.534e+07
* channel (channel) int32 4kB 1 2 3 4 5 6 ... 953 954 955 956 957 958
group_coord_y (shot_point, channel) float64 8kB 0.0 ... 4.535e+07
group_coord_x (shot_point, channel) float64 8kB 0.0 7.734e+06 ... 7.76e+06
Attributes:
unitsV1: {'length': 'ft'}As we mentioned before, we can also get specific dataset slices of headers while fetching a slice.
Let’s overlay receiver coordinates on the offset-sorted gather.
group_x_header = comp["group_coord_x"]
group_y_header = comp["group_coord_y"]
image = comp["amplitude"]
gs_kw = {"height_ratios": (1, 8)}
fig, (hdr_ax, img_ax) = plt.subplots(2, 1, figsize=(7, 6), gridspec_kw=gs_kw, sharex="all")
hdr_ax2 = hdr_ax.twinx()
group_x_header.plot(ax=hdr_ax, x="offset", c="red")
group_y_header.plot(ax=hdr_ax2, x="offset", c="black")
image.plot(
ax=img_ax,
x="offset",
y="time",
vmin=-clip,
vmax=clip,
cmap="gray_r",
yincrease=False,
add_colorbar=False,
)
img_ax.set_title("")
hdr_ax.set_xlabel("")
plt.tight_layout()
MDIO to SEG-Y Conversion¶
Finally, let’s demonstrate going back to SEG-Y.
We will use the convenient mdio_to_segy function and write it out as a round-trip file.
The output spec can be modified if we want to write things to different byte locations, etc, but we will use the same one as before.
from mdio import mdio_to_segy
mdio_to_segy(
input_path="soda_lake.mdio",
output_path="soda_lake_roundtrip.sgy",
segy_spec=soda_lake_segy_spec,
)
Validate Round-Trip SEG-Y File¶
We can validate if the round-trip SEG-Y file is matching the original using TGSAI/segy.
Step by step:
Open original file
Open round-trip file
Compare text headers
Compare binary headers
Compare 100 random headers and traces
import numpy as np
from segy import SegyFile
original_segy = SegyFile(input_url)
roundtrip_segy = SegyFile("soda_lake_roundtrip.sgy")
# Compare text header
assert original_segy.text_header == roundtrip_segy.text_header
# Compare bin header
assert original_segy.binary_header == roundtrip_segy.binary_header
# Compare 100 random trace headers and samples.
# IBM float32 cannot preserve denormals through IEEE decode/encode.
rng = np.random.default_rng()
rand_indices = rng.integers(low=0, high=original_segy.num_traces, size=100)
for idx in rand_indices:
np.testing.assert_equal(original_segy.header[idx], roundtrip_segy.header[idx])
np.testing.assert_allclose(original_segy.sample[idx], roundtrip_segy.sample[idx], rtol=0, atol=1e-37)
print("Files identical!")
Files identical!