Fast two-point seismic ray tracing in layered media.
LayTracer is an open-source Python package for computing ray paths, travel times, and amplitude attributes in horizontally layered (1D) velocity models with constant layer velocities. It also computes sparse fixed-topology derivatives of traveltime and ray parameter with respect to layer velocities, interface depths, and source/receiver coordinates for sensitivity analysis and inversion. The ray tracing core algorithm is based on the dimensionless ray parameter method of Fang & Chen (2019), achieving rapid convergence.
Documentation: https://danikiev.github.io/LayTracer
| Category | Capability |
|---|---|
| Ray tracing | Second-order dimensionless-q iteration with a checked bracketed fallback and independently verified endpoint residual |
| Phase itineraries | Explicit prescribed direct, reflected, transmitted, converted, and multiple P–SV/SH paths |
| Diagnostics | Opt-in solve route, iterations, residual, conditioning, and criticality margin for every ray |
| Sensitivities | Sparse fixed-topology derivatives of traveltime and physical ray parameter with respect to layer velocities, interface depths, and endpoints |
| Travel time | Layer-by-layer travel time summation from the solved ray parameter |
| Attenuation | Intrinsic absorption operator t* from quality factors Q |
| Relative amplitude attributes | Geometrical spreading, attenuation, and angle-dependent P-SV/SH coefficient products; these are not complete synthetic amplitudes |
| Brewster angles | Automatic detection of Brewster-like zeros in R/T coefficient curves |
| Parallel execution | Multi-ray tracing with joblib / loky backend for large surveys |
| Visualisation | 2-D ray path plots (matplotlib) and interactive 3-D viewer (Plotly) |
| Documentation | Comprehensive Sphinx docs with extensive theory, gallery examples, and API reference |
- Python 3.9 or newer and
pip - Conda package manager (recommended for full reproducible setup via miniforge)
# 1. Create the environment base (Python and pip) from environment.yml
conda env create -f environment.yml
# 2. Activate it
conda activate laytracer
# 3. Install LayTracer in editable mode with the development dependencies defined in pyproject.toml
python -m pip install -e ".[dev]"For a minimal runtime-only installation, use python -m pip install -e . instead.
python -m pip install --upgrade pip
pip install laytracerpython -c "import laytracer; print(laytracer.__version__)"Alternative:
python -m pip show laytracerUse this mode when you want a published stable release. For development or latest unreleased changes, install from the repository with pip install -e ..
install.batchmod +x install.sh
./install.sh# Linux / macOS
python -m venv .venv
source .venv/bin/activate
python -m pip install --upgrade pip
pip install -e .
# Windows (PowerShell)
python -m venv .venv
.\.venv\Scripts\Activate.ps1
python -m pip install --upgrade pip
pip install -e .pip-only installation is suitable for running LayTracer. In either workflow, add the optional extras defined in pyproject.toml for development (".[dev]"), documentation (".[docs]") or both (".[dev,docs]").
- Python ≥ 3.9
- NumPy (≥ 2), SciPy, Pandas
- Matplotlib, Plotly, cmcrameri
- psutil, joblib (≥ 1.4)
When NumPy is imported, its OpenBLAS library may reserve a working buffer for every thread it can use. On Windows this memory is committed at once. With conda-forge's OpenBLAS it is about 135 MB per thread, so on a machine with 22 logical CPUs a Python process commits about 3 GB before doing any work. LayTracer's parallel workers already use one BLAS thread each, because ray tracing makes no BLAS calls. For your own Python process, limit the threads, for example to four:
conda env config vars set OPENBLAS_NUM_THREADS=4 -n laytracer
conda deactivate
conda activate laytracerThe setting applies whenever the environment is activated. Remove it with conda env config vars unset OPENBLAS_NUM_THREADS -n laytracer. Without conda, set the environment variable OPENBLAS_NUM_THREADS before Python starts.
import laytracer
import numpy as np
import pandas as pd
vel_df = pd.DataFrame({
"Depth": [0.0, 1000.0, 2000.0, 3500.0],
"Vp": [3000.0, 4500.0, 5500.0, 6500.0],
"Vs": [1500.0, 2250.0, 2750.0, 3250.0],
"Rho": [2200.0, 2500.0, 2700.0, 2900.0],
"Qp": [200.0, 400.0, 600.0, 800.0],
"Qs": [100.0, 200.0, 300.0, 400.0],
})source = np.array([0.0, 0.0, 3000.0])
receiver = np.array([5000.0, 1000.0, 0.0])
result = laytracer.trace_rays(
source,
receiver,
vel_df,
requested={"travel_times", "rays", "ray_parameters"},
)
print(f"Travel time: {result.travel_times[0]:.4f} s")
print(f"Ray parameter: {result.ray_parameters[0]:.6e} s/m")itinerary = laytracer.RayItinerary(
source_phase="P",
interactions=[
laytracer.Interaction(3500.0, "reflect", "SV"),
],
)
converted = laytracer.trace_rays(
source,
receiver,
vel_df,
itinerary=itinerary,
requested={
"travel_times",
"rays",
"ray_parameters",
"diagnostics",
"sensitivities",
},
)
certificate = converted.diagnostics[0]
kernel = converted.sensitivities[0]
print(certificate.method, certificate.absolute_offset_residual)
print(kernel.vp_layer_indices, kernel.dtravel_time_dvp)
print(kernel.vs_layer_indices, kernel.dtravel_time_dvs)The itinerary is fixed by the caller. LayTracer does not search for phase
branches or guarantee a global first arrival. The legacy reflection= and
refraction= tuple arguments remain available; in that API refraction
means a prescribed transmission or mode conversion, not a head wave.
approximator = laytracer.TravelTimeApproximator.fit(
sources,
receivers,
vel_df,
source_max_distance=100.0,
receiver_max_distance=None, # keep every receiver as an exact anchor
source_phase="P",
)
prediction = approximator.predict(order=2)
times = prediction.travel_time_matrix
valid = prediction.valid_matrix
print(approximator.exact_ray_count, times.shape)Anchor placement accepts arbitrary 3-D point clouds. Predictions use only
nearby anchors with the same ordered layer, phase, and vertical-direction
topology. Unsupported targets are returned as NaN with a validity mask and
reason code; LayTracer does not silently retrace them. order=1 is the
default linear endpoint predictor. order=2 adds analytic fixed-topology
endpoint curvature and uses the exact straight-path expression when both
endpoints lie in one homogeneous layer.
src = np.array([0.0, 0.0, 3000.0])
rcvs = np.column_stack([
np.arange(500, 15001, 500),
np.zeros(29),
np.zeros(29),
])
result = laytracer.trace_rays(
sources=src,
receivers=rcvs,
velocity_df=vel_df,
source_phase="P", # also accepts "SV", "SH", "S", or a list such as ["P", "SH", "SV"]
requested={"travel_times", "rays", "ray_parameters", "tstar", "spreading", "trans_product"},
transcoef_method="standard", # standard Zoeppritz coefficients without normalization
)
# Access results
print(result.travel_times) # travel times (s)
print(result.tstar) # attenuation operator t*
print(result.spreading) # geometrical spreading
print(result.trans_product) # nonnegative coefficient-magnitude product# 2-D ray paths over velocity cross-section
laytracer.plot.rays_2d(vel_df, rays=[r.ray_path for r in ...])
# 1-D velocity profile
laytracer.plot.velocity_profile(vel_df, param="Vp")
# Interactive 3-D viewer
fig = laytracer.plot.rays_3d(vel_df, rays=result.rays, sources=src, receivers=rcvs)
fig.show()| Symbol | Description |
|---|---|
LayerStack |
Data class holding layer thicknesses, velocities (Vp, Vs), densities, and Q-factors |
ModelArrays |
Pre-extracted NumPy arrays from a velocity DataFrame for efficient repeated tracing |
build_layer_stack(vel_model, z_src, z_rcv) |
Extract the traversed layer stack between source and receiver depths (accepts DataFrame or ModelArrays) |
| Symbol | Description |
|---|---|
solve(h, v, segments, interactions, ...) |
Low-level two-point solve for a compiled path |
RayResult |
Result container: travel time, ray path, ray parameter, optional attributes, and diagnostics |
SolveDiagnostics |
Solve route, iteration counts, transformed and physical parameters, endpoint residual, conditioning, and criticality margin |
offset(q, h, lmd) |
Total horizontal offset X(q) |
offset_dq(q, h, lmd) |
First derivative dX/dq |
offset_dq2(q, h, lmd) |
Second derivative d²X/dq² |
q_from_p(p, vmax) / p_from_q(q, vmax) |
Convert between slowness p and dimensionless q |
initial_q(X_target, h, lmd) |
Asymptotic initial estimate for Newton iteration |
newton_step(q_i, X_target, h, lmd) |
One quadratic Newton step |
| Symbol | Description |
|---|---|
psv_rt_coefficients(p, vp1, vs1, rho1, vp2, vs2, rho2) |
All 8 P-SV reflection/transmission coefficients (Zoeppritz) |
sh_rt_coefficients(p, vs1, rho1, vs2, rho2) |
Decoupled SH-SH reflection/transmission coefficients |
normalize_rt_coefficient(coeff, p, v_in, rho_in, v_out, rho_out) |
Energy-flux normalization of R/T coefficients (Červený, 2001) |
find_brewster_angles(rt_coefficients, angles, ...) |
Detect Brewster-like zeros in R/T curves |
| Symbol | Description |
|---|---|
trace_rays(sources, receivers, velocity_df, ...) |
Trace all source–receiver pairs with optional parallelism; multi-phase input returns a phase-keyed result dictionary |
Interaction / RayItinerary |
Ordered, fixed-topology reflection and transmission/mode-conversion description |
TraceResult |
Container: source phase, travel times, ray paths, ray parameters, optional attributes, diagnostics, and sensitivities |
RaySensitivity |
Sparse fixed-topology derivatives with respect to layer velocities, interface depths, and endpoint coordinates |
| Symbol | Description |
|---|---|
linearized_ray_change(...) |
Apply one sparse sensitivity record to model and endpoint perturbations |
select_anchors(points, max_distance, groups=None) |
Deterministic grouped point-cloud anchor coverage |
TravelTimeApproximator |
Fit exact topology-compatible anchors and predict nearby traveltimes at first or second endpoint-Taylor order |
TravelTimePrediction |
Predicted times, validity/reason arrays, and anchor-assignment diagnostics |
| Function | Description |
|---|---|
velocity_profile(vel_df, ...) |
1-D velocity–depth step profile (matplotlib) |
rays_2d(vel_df, rays, ...) |
2-D ray paths over layered velocity cross-section (matplotlib) |
rays_3d(vel_df, rays, ...) |
Interactive 3-D ray visualisation (Plotly) |
Full documentation is built with Sphinx and includes:
- Getting Started — installation and environment setup
- Methodology — complete mathematical derivations (dimensionless parameter, Newton iteration, travel time, t*, geometrical spreading, Zoeppritz coefficients, critical & Brewster angles, 3-D extension)
- Examples Gallery — runnable scripts rendered with Sphinx-Gallery
- API Reference — auto-generated from docstrings with numpydoc
Install the optional documentation dependencies defined in pyproject.toml first:
conda activate laytracer
python -m pip install -e ".[docs]"Then build:
conda activate laytracer
# Windows
build-docs.bat
# Linux / macOS
chmod +x build-docs.sh
./build-docs.shBuild docs with PDF output:
conda activate laytracer
# Windows
build-docs.bat -pdf
# Linux / macOS
chmod +x build-docs.sh
./build-docs.sh -pdfYou can do also using make commands:
conda activate laytracer
cd docs
# Build HTML
make html
# Build PDF
make latexpdf
# Run a local server to view HTML
cd build/html
python -m http.serverThis repository includes a GitHub Actions workflow at .github/workflows/docs.yml that runs on every push to main and:
- builds Sphinx HTML docs,
- builds the PDF (
laytracer.pdf), - copies the PDF into the published site (
_static/laytracer.pdf), - deploys HTML docs to GitHub Pages,
- uploads the PDF as a workflow artifact.
One-time GitHub setup:
- Open Settings → Pages in your GitHub repository.
- Set Source to GitHub Actions.
- Push to
main.
Published docs URL:
https://danikiev.github.io/LayTracer
LayTracer implements the method of Fang & Chen (2019) for two-point ray tracing in horizontally layered media:
-
Dimensionless ray parameter q = p · v_max / √(1 − p² · v²_max) maps the full range of take-off angles to [0, ∞) without singularities.
-
Offset equation X(q) is a smooth, monotonically increasing function — ideal for Newton iteration.
-
Safeguarded transformed-parameter solver attempts rapid Newton iteration, falls back to a checked bracketed solve when needed, and accepts a solution only after independently recomputing the endpoint offset.
-
Amplitude attributes are computed inline:
- Travel time from vertical slowness summation
- Attenuation t* from per-layer Q-factors
- Geometrical spreading from analytic ∂X/∂p
- Full P-SV R/T coefficients and decoupled SH-SH coefficients (Lay & Wallace (1995) formulation)
- Fang, X. & Chen, X. (2019). A fast and robust two-point ray tracing method in layered media with constant or linearly varying layer velocity. Geophysical Prospecting, 67(7), 1811–1824. doi:10.1111/1365-2478.12799
- Aki, K. & Richards, P.G. (2002). Quantitative Seismology. 2nd ed., University Science Books.
- Lay, T. & Wallace, T.C. (1995). Modern Global Seismology. Academic Press.
- Červený, V. (2001). Seismic Ray Theory. Cambridge University Press. doi:10.1017/CBO9780511529399
LayTracer includes a comprehensive test suite covering the solver, amplitude calculations, API, and physical symmetries.
conda activate laytracer
pytestTest modules:
test_solver.py— Newton convergence, Snell's law, travel time accuracytest_amplitude.py— Zoeppritz coefficients, energy-flux normalization, Brewster detectiontest_api.py— multi-ray tracing interfacetest_sensitivity.py— analytic traveltime and ray-parameter derivativestest_generalized.py— generalized layered-media validation casestest_homogeneous_equivalence.py— homogeneous-medium equivalence checkstest_plot.py— plotting helpers and visualisation optionstest_symmetry.py— reciprocity and physical consistency checks
LayTracer/
├── laytracer/ # Main package
│ ├── __init__.py # Public API exports
│ ├── model.py # LayerStack, ModelArrays, build_layer_stack
│ ├── solver.py # Core ray tracing solver (q-parameter + Newton)
│ ├── sensitivity.py # Sparse fixed-topology analytic sensitivities
│ ├── amplitude.py # Transmission coefficients, Zoeppritz, Brewster
│ ├── api.py # High-level multi-ray interface (trace_rays)
│ └── plot.py # Visualisation (2-D, 3-D, velocity profiles)
├── examples/ # Sphinx-Gallery example scripts
│ ├── 01_basic_raytracing.py
│ ├── 02_paper_examples.py
│ ├── 03_reflection_transmission.py
│ ├── 04_amplitude_analysis.py
│ ├── 05_homogeneous_equivalence.py
│ ├── 06_visualizing_sensitivities.py
│ └── README.txt
├── branding/ # Branding assets and logo source
│ └── logo/ # Logo generator, vendored Poppins font, and usage docs
├── pytests/ # Test suite
│ ├── test_solver.py
│ ├── test_amplitude.py
│ ├── test_api.py
│ ├── test_sensitivity.py
│ ├── test_generalized.py
│ ├── test_homogeneous_equivalence.py
│ ├── test_plot.py
│ └── test_symmetry.py
├── docs/ # Sphinx documentation
│ └── source/
| ├── api/ # Generated API reference pages
│ | └── index.rst # API reference landing page
| ├── changelog.rst # Changelog page with a placeholder for changelog content generated from CHANGELOG.md
│ ├── citing.rst # Citation information and BibTeX entries
│ ├── conf.py # Sphinx configuration
│ ├── credits.rst # Author and contributor acknowledgements and library credits
│ ├── getting_started.rst # Installation and quick start
│ ├── index.rst # Main index page with links to all sections
│ ├── methodology.rst # Full mathematical derivations
│ └── references.bib # Bibliography for the docs
├── pyproject.toml # Build configuration (setuptools + setuptools-scm)
├── environment.yml # Conda environment base (Python and pip)
├── pytest.ini # Pytest configuration
├── CHANGELOG.md # Changelog
├── LICENSE # MIT License file
├── RELEASE.md # Release guidelines and GitHub Actions workflow documentation
└── README.md # This README file
LayTracer is released under the MIT License.
Denis Anikiev — danikiev@gmail.com
If you use a LayTracer in your research, please cite:
@software{Anikiev2026LayTracer,
author = {Anikiev, Denis},
title = {{LayTracer}: {F}ast two-point seismic ray tracing in layered media},
year = {2026},
publisher = {Zenodo},
doi = {10.5281/zenodo.18850919},
url = {https://github.com/danikiev/LayTracer}
}The DOI in the above citation points to the latest version of LayTracer. If you use a specific version, please cite that version's DOI (available on the Zenodo release page), e.g.:
@software{Anikiev2026LayTracerVersion,
author = {Anikiev, Denis},
title = {{LayTracer}: {F}ast two-point seismic ray tracing in layered media},
year = {2026},
publisher = {Zenodo},
version = {0.4.0},
license = {MIT},
doi = {10.5281/zenodo.19983599},
url = {https://github.com/danikiev/LayTracer}
}