Skip to content

WP8: SEG-Y and LAS adapters, the iwt command, and the critical assessment - #9

Merged
ingehap merged 1 commit into
mainfrom
claude/blissful-gauss-3y1u5w
Oct 15, 2026
Merged

ingehap merged 1 commit into
mainfrom
claude/blissful-gauss-3y1u5w

Conversation

@ingehap

@ingehap ingehap commented Oct 15, 2026

Copy link
Copy Markdown
Owner

Summary

Work package 8 of docs/IMPLEMENTATION_PLAN.md, the last one. It adds three things:

  • file adapters, so recorded data can replace the synthetic set;
  • the iwt command, which runs the ZO-VSP → model → tie chain on files;
  • the README's consolidated assessment of the article.

New code

  • iwt/io.py:
    • SEG-Y in numpy, revisions 1 and 2 with fixed-length traces: IBM, IEEE and integer samples, both byte orders, EBCDIC or ASCII text, and scaled elevation and coordinate words. The core needs no new dependency;
    • ZO-VSP and CMP-gather files, with a configurable receiver-depth word, sign and datum shift;
    • LAS through lasio: common mnemonics; us/ft, us/m, g/cm³ and kg/m³; depth from a TVDSS curve or MD − KB; missing samples dropped; optional resampling;
    • CSV logs, layer models (with Q and labels) and time-depth tables.
  • iwt/vsp/processing.py: the WP3 ZO-VSP chain in one call (picks, separation, 30 Hz times, interval velocities, t*, interval Q, corridor stack). The stratigraphic synthetic is driven by the recorded first arrival.
  • iwt/cli.py: iwt demo | vsp | model | tie, plus python -m iwt and a console entry point.
  • tests/test_io_cli.py (14 tests):
    • IBM floats against the textbook value −118.625 = 0xC276A000;
    • the byte layout against fixed byte positions;
    • segyio in both directions (segyio added to the dev extras only; wheels exist for cp39 and cp313);
    • LAS unit and depth handling;
    • the processing chain on a model with known Q;
    • the command line end to end.
  • README: a command-line section, WP8 results and a real-data guide. It also adds the article's claims tested on the synthetic well, the consolidated critical assessment and the errata (one new: the article does not give the inputs of its sonic-to-seismic correction).
  • Plan: WP8 marked done, with as-built notes. Version 0.8.0.

Findings

  • Intrinsic Q depends on the wavelet driving the stratigraphic synthetic:

    • with the recorded first arrival: 45.5/59.0/120/72 (true 45.5/58.8/80/70);
    • with a Ricker at the first arrivals' peak frequency (47 Hz, raised by the ghost): the middle block is 20 % high.

    The recorded arrival beats even the true source wavelet.

  • iwt model: the model's times are within 0.83 ms RMS of the VSP. What remains is the bias of the reference pick plus the ray-theory against effective-medium scale effect.

  • iwt tie:

    • corridor stack 0.944;
    • seismic trace 0.950 with Q in the synthetic, and 0.909 with a 16° phase rotation without it;
    • an estimated wavelet scores 0.986, but that tie is in-sample and is documented as a diagnostic only.

Not built (documented)

DLIS, deviation surveys, 3C rotation, walkaway picking, and command-line access to walkaway inversion and CMP migration (library only). Whether the Hanssen data are released was not checked.

Test plan

  • pytest -q -W error::RuntimeWarning -W error::DeprecationWarning (Python 3.13): 168 passed
  • pytest -q (Python 3.9, without lasio and segyio): 161 passed, 4 skipped
  • pyflakes iwt tests examples
  • pytest -W error tests/test_io_cli.py: clean
  • README quick start runs; pip install -e . gives a working iwt command (--version, every --help)
  • Demo chain iwt demo → vsp → model → tie reproduces the README numbers

🤖 Generated with Claude Code

https://claude.ai/code/session_019vaUDUPicoFeQaUhfhnoYq


Generated by Claude Code

…ment

- iwt.io: numpy SEG-Y reader/writer (rev 1/2, fixed-length traces, IBM/IEEE/int
  samples, both byte orders, scaled header words), ZO-VSP and gather files,
  LAS via lasio (mnemonics, us/ft-us/m, g/cc-kg/m3, TVDSS or MD - KB), CSV logs,
  layer models and time-depth tables
- iwt.vsp.processing: the WP3 ZO-VSP chain in one call; the stratigraphic t* is
  driven by the recorded first arrival (a peak-frequency Ricker biased Q by 20 %)
- iwt.cli: iwt demo / vsp / model / tie, python -m iwt, console entry point
- tests/test_io_cli.py (14 tests): IBM floats, byte layout, segyio cross-check
  both ways, LAS units and depths, processing chain, CLI end to end
- README: command line, WP8 results, using real data, the article's claims
  tested, critical assessment, errata; plan marked done; version 0.8.0

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_019vaUDUPicoFeQaUhfhnoYq
@ingehap
ingehap marked this pull request as ready for review October 15, 2026 06:38
Copilot AI balanced review requested due to automatic review settings October 15, 2026 06:38
@ingehap
ingehap merged commit 5255fbd into main Oct 15, 2026
4 checks passed

Copilot AI left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

🟡 Changes recommended

Recording-delay, geometry-unit and Q/datum handling can corrupt results, while accepted input cases can cause command failures.

11 open findings
What changed in this PR

Completes WP8 by connecting file-based data to the existing well-tie workflow and documenting the synthetic results and limitations.

Changes:

  • Adds SEG-Y, LAS and CSV adapters.
  • Introduces the ZO-VSP processing chain and iwt commands.
  • Adds interoperability tests and consolidates the article assessment.
File Description
tests/​test_io_cli.py Adapter, processing and CLI tests.
README.md Usage, results and critical assessment.
pyproject.toml CLI entry point, test dependency and version.
iwt/​vsp/​processing.py Integrated ZO-VSP processing chain.
iwt/​vsp/​__init__.py Documents the processing module.
iwt/​io.py SEG-Y, LAS and CSV adapters.
iwt/​cli.py Demo, VSP, model and tie commands.
iwt/​__main__.py Enables python -m iwt.
iwt/​__init__.py Exports adapters and updates version.
docs/​IMPLEMENTATION_PLAN.md Records WP8 completion and limitations.

🧠 Review effort: Balanced


💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.

Comment thread iwt/io.py
Comment on lines +345 to +348
depths = depth_sign * sg.headers[depth_field] + depth_shift
order = np.argsort(depths, kind="stable")
zs = float(np.median(sg.headers["source_depth"])) if source_depth is None else source_depth
xs = float(np.median(np.abs(sg.headers["offset"]))) if source_offset is None else source_offset
Comment thread iwt/io.py
Comment on lines +379 to +382
off = np.abs(sg.headers["offset"]).astype(float)
order = np.argsort(off, kind="stable")
return Gather(sg.t, off[order], sg.data[order], float(np.median(sg.headers["source_depth"])),
float(np.median(-sg.headers["receiver_elevation"])))
Comment thread iwt/vsp/processing.py
Comment on lines +133 to +136
syn = Z.synthetic_zovsp(strat_model, vsp.depths, vsp.source_depth, dt, t.size,
source_wavelet(t, down, pk), f_max=f_max, attenuation=False,
source="downgoing", workers=workers)
cs_strat = _chain(syn.t, syn.data, f_seismic, band)[4]
Comment thread iwt/cli.py
if args.vsp:
vsp = json.loads(Path(args.vsp).read_text())
if qp is None:
q = vsp["q"]
Comment thread iwt/cli.py
qv = np.where(inv > 0, 1.0 / np.where(inv > 0, inv, 1.0), np.inf)
qp = _q_curve(z, q["top"], q["bottom"], qv)
note.append(f"{q['kind']} Q from {args.vsp}")
dz, tz = np.asarray(vsp["depths"], float), np.asarray(vsp["vertical_time"], float)
Comment thread iwt/cli.py
sg = IO.read_segy(args.seismic)
obs = sg.data[args.trace]
t, dt = sg.t, sg.dt
i0 = int(round(sg.t0 / dt))
Comment thread iwt/cli.py
else:
spike = (np.array([0.0, dt]), np.array([1.0, 0.0]))
refl = normal_incidence(model, dt, nt, spike, **kw)[1][i0:]
tw, w = estimate_wavelet(obs[win], refl[win], dt, args.half_length, args.damping)
Comment thread iwt/io.py
samples = np.ascontiguousarray(raw[:, 240:]).view(order + kind).reshape(n_traces, ns)
data = ibm_to_ieee(samples) if fmt == 1 else samples.astype(np.float64)
dt_us = binary["interval"] & 0xFFFF or int(headers["sample_interval"][0])
t0 = float(headers["delay"][0]) * 1e-3 if n_traces else 0.0
Comment thread iwt/io.py
write_segy(path, vsp.data, vsp.dt, dict(
receiver_elevation=-vsp.depths, source_depth=np.full(n, vsp.source_depth),
offset=np.full(n, vsp.source_offset), field_record=np.ones(n),
trace_number=np.arange(1, n + 1)), text=text or "ZERO-OFFSET VSP WRITTEN BY IWT", fmt=fmt)
Comment thread iwt/vsp/processing.py
q = A.interval_q(vsp.depths, tv, dts, q_blocks)
corridor = None
if corridor_peak is not None:
dec = Z.deconvolve(t, up, down, tv, ricker(corridor_peak, vsp.dt))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants