search

paudetseis / PyRaysum

Teleseismic body wave modeling through stacks of (dipping/anisotropic) layers
https://paudetseis.github.io/PyRaysum/
MIT License
43 stars 14 forks source link
anisotropy fortran geophysics modeling plane-waves python receiver-functions seismic-imaging seismology

readme

Software for modeling ray-theoretical body-wave propagation

This program generates sets of ray-theoretical seismograms for incident plane waves (teleseismic approximation) for seismic velocity models consisting of a stack of layers with planar but non-parallel (dipping) interfaces, allowing the possibility of anisotropy in the layers. Incident P and S waves are supported.

PyRaysum is a Python wrapper around the Fortran software Raysum, originally developed by Andrew Frederiksen. A trimmed down version of the Fortran code is supplied with PyRaysum. You can find the original version here.

DOI build codecov GitHub Code style: black

Authors: Wasja Bloch, Pascal Audet (Developers and Maintainers of PyRaysum) & Andrew Frederiksen & (Developer of original Fortran version)

Installation

PyRaysum can be installed from PyPI or from source.

To avoid conflicts with other programs, it is recommended to install PyRaysum inside a designated conda environment (called here prs) alongside its dependecies.

conda create -n prs "python<3.12" "numpy<1.23" setuptools fortran-compiler obspy -c conda-forge
conda activate prs

If you are using an alternative Python interpreter (e.g., IPython or a Jupyter notebook), include it into your installation as well. This ensures that the interpreter uses the correct Python version.

# IPython
conda install ipython

# Jupyter notebooks
conda install jupyter
Installing from PyPI
pip install pyraysum
Installing from source

The source code of PyRaysum can also be downloaded from GitHub and installed via pip.

git clone https://github.com/paudetseis/PyRaysum.git
cd PyRaysum
pip install .

Getting Started

To compute receiver functions for a range of back-azimuths considering a simple 1-layer over a half-space subsurface seismic velocity model, in a Python or iPython console, execute:

from pyraysum import Model, Geometry, Control, run

# Build a 1-layer-over-half-space subsurface model
model = Model(
    thickn=[32000, 0],  # m; half-space thickness is irrelevant
    rho=[2800, 3600],  # kg/m^3
    vp=[6400, 8100],  # m/s
    vs=[3600, 4650],  # m/s
)
model.plot()

# Define back-azimuth range and single horizontal slowness value
geometry = Geometry(baz=range(0, 360, 30), slow=0.07)
geometry.plot()

# Set sampling interval, number of points, alignment and ray-coordinate rotation
control = Control(dt=1./20., npts=800, align="P", rot="PVH")

# Run the simulation
result = run(model, geometry, control, rf=True)

# Filter below 2s period
result.filter("rfs", "lowpass", freq=0.5)

# Plot the results
result.plot("rfs")

Documentation

The complete API documentation, scripts and tutorials are described at https://paudetseis.github.io/PyRaysum/

Citing

If you use PyRaysum in your work, please cite the following references:

Contributing

All constructive contributions are welcome, e.g. bug reports, discussions or suggestions for new features. You can either open an issue on GitHub or make a pull request with your proposed changes. Before making a pull request, check if there is a corresponding issue opened and reference it in the pull request. If there isn't one, it is recommended to open one with your rationale for the change. New functionality or significant changes to the code that alter its behavior should come with corresponding tests and documentation. If you are new to contributing, you can open a work-in-progress pull request and have it iteratively reviewed.

Other examples of contributions include notebooks that describe published examples of PyRaysum usage and processing. Suggestions for improvements (speed, accuracy, plotting, etc.) are also welcome.