Repository navigation
WP6: ray tracing and Q-compensated depth imaging (the Fig. 3 experiment) - #7
Merged
Merged
Conversation
- iwt.raytrace: exact qP group velocity; second derivative of q(p); reflection_table with two-way PP time, ray parameter, reflection angle, t* = sum t/Q and the stationary-phase (true-amplitude) spreading factor for every depth x offset in layered VTI, by Newton refinement on exact sums. - iwt.imaging: first-order Kjartansson constant-Q operator in numpy's FFT convention; inverse Q with Wang's (2002) stabilization set by the maximum gain; 1D true-amplitude prestack depth migration of CMP gathers with a per-output-point inverse Q filter (ray-based, Cavalca et al., 2013), amplitude+phase or phase-only, common high cut. - iwt.avo: event tracking on migrated gathers; A + B sin^2 fits. - examples/wp6_imaging.py and docs/figures/wp6_*.png: Alkhalifah-Tsvankin moveout, reflectivity Born gathers with and without Q, eight migration variants scored on the Sto-top AVA, gain limit against noise. - Tests (15 new), README (sections reordered by WP) and plan; version 0.6.0. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_019vaUDUPicoFeQaUhfhnoYq
There was a problem hiding this comment.
🟡 Changes recommended
Unresolved correctness issues affect valid reflection rays, Q filtering and event tracking.
4 open findings
What changed in this PR
Implements WP6, extending the well-tie workflow with layered VTI reflection rays and Q-compensated depth imaging for the Fig. 3 experiment.
Changes:
- Adds reflection tables, stabilized inverse-Q filtering, migration and AVA extraction.
- Adds numerical validation and an imaging comparison example.
- Documents WP6 results and updates the package version.
| File | Description |
|---|---|
| tests/test_imaging.py | Tests kinematics, Q operators, migration and AVA. |
| README.md | Documents WP6 usage, results and assumptions. |
| pyproject.toml | Updates version to 0.6.0. |
| iwt/raytrace.py | Adds group velocity and reflection tables. |
| iwt/imaging.py | Adds Q filtering and depth migration. |
| iwt/avo.py | Adds event tracking and AVA fitting. |
| iwt/__init__.py | Exports new modules and updates version. |
| examples/wp6_imaging.py | Demonstrates migration variants and gain sensitivity. |
| docs/IMPLEMENTATION_PLAN.md | Records WP6 completion and scope. |
🧠 Review effort: Balanced
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
Comment on lines
+298
to
+302
| q, dq, d2q = qp_vertical_slowness(p[:, None], *a, derivative=2) | ||
| X = -(wk * dq).sum(axis=1) | ||
| dX = -(wk * d2q).sum(axis=1) | ||
| p = np.clip(p - (X - xi) / dX, 0.0, p_cap[k - 1] * (1.0 - 1e-12)) | ||
| q, dq, d2q = qp_vertical_slowness(p[:, None], *a, derivative=2) |
| order = list(range(start, n)) + list(range(start - 1, -1, -1)) | ||
| for j in order: | ||
| if j == start - 1: | ||
| centre = out.depth[start] |
Comment on lines
+69
to
+73
| f = np.abs(np.asarray(f, dtype=float)) | ||
| tstar = np.asarray(tstar, dtype=float) | ||
| with np.errstate(divide="ignore", invalid="ignore"): | ||
| lnf = np.where(f > 0, np.log(np.where(f > 0, f, 1.0) / f_ref), 0.0) | ||
| return np.exp(-np.pi * f * tstar) * np.exp(2j * f * tstar * lnf) |
Comment on lines
+78
to
+90
| from scipy.optimize import brentq | ||
|
|
||
| if not gain_limit > 1.0: | ||
| raise ValueError("gain_limit must exceed 1") | ||
| if not np.isfinite(gain_limit): | ||
| return 0.0 | ||
|
|
||
| def peak(log_s): | ||
| s = np.exp(log_s) | ||
| r = np.sqrt(s * s + s) # beta* + s at the maximum | ||
| return r / (2.0 * s * (s + 1.0 - r)) - gain_limit | ||
|
|
||
| return float(np.exp(brentq(peak, -60.0, 10.0))) |
This file contains hidden or 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
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
This PR implements work package 6 of
docs/IMPLEMENTATION_PLAN.md: layered VTI ray tracing and ray-based Q-compensated depth migration (Cavalca et al., 2013). With them it reproduces the article's Fig. 3 experiment.What is added
iwt.raytracereflection_table: for every depth × offset in layered VTI it returns the two-way PP time, ray parameter, reflection angle and t* = Σt/Q. It also returns the stationary-phase spreading factor that turns a point-source amplitude into the plane-wave R(θ).iwt.imagingiwt.avo: event tracking across migrated gathers, andA + B sin²θfits.examples/wp6_imaging.pyanddocs/figures/wp6_*.png:Main findings
The "observed" gather has Q; the reference is the Q-free gather migrated with the true VTI model. AVA is compared over 20–44°.
Validation
🤖 Generated with Claude Code
https://claude.ai/code/session_019vaUDUPicoFeQaUhfhnoYq
Generated by Claude Code