The resistivity model

Everything in pycsamt.interp starts from a ResistivityModel — a 2-D section of \(\log_{10}\rho\) on a (depth, distance) grid, plus station positions and metadata. It usually comes from an inversion (ResistivityModel.from_occam2d), but you can build one directly from arrays, which is what this whole gallery does so every deliverable traces back to a known model.

This first example builds the shared synthetic section and reads it as 1-D resistivity-depth profiles.

Build a model from arrays

ResistivityModel.from_array takes a (n_z, n_x) array of \(\log_{10}\rho\) with cell-centre depth and distance axes. Here the section is a classic hydrogeological sequence: a resistive dry overburden, a conductive sand aquifer, a very conductive clay aquitard, and a resistive basement.

import matplotlib.pyplot as plt
import numpy as np
from _interp_data import demo_model

# Use the resistivity section (2nd figure) as the card thumbnail.

rm = demo_model()
print(rm)
print(
    "grid (n_z, n_x):",
    rm.rho_2d.shape,
    f" depth 0-{rm.z_centers.max():.0f} m,  distance 0-{rm.x_centers.max():.0f} m",
)

ResistivityModel(x_centers=ndarray(shape=(44,), dtype=float64), z_centers=ndarray(shape=(34,), dtype=float64), rho_2d=ndarray(shape=(34, 44), dtype=float64), station_x=ndarray(shape=(11,), dtype=float64), station_names=list(len=11, sample=['S00', 'S01', 'S02', ...]), method='synthetic', rms=nan, metadata=dict(len=0, keys=[]))
grid (n_z, n_x): (34, 44)  depth 0-300 m,  distance 0-2000 m

The section itself

A quick look at the model with an EM-style colour scale: cool colours are conductive (water, clay), warm colours resistive (dry ground, basement).

fig, ax = plt.subplots(figsize=(10, 4.2), constrained_layout=True)
im = ax.pcolormesh(
    rm.x_centers, rm.z_centers, rm.rho_2d, cmap="Spectral", shading="auto"
)
ax.plot(
    rm.station_x,
    np.full_like(rm.station_x, rm.z_centers.min()),
    "kv",
    ms=6,
    clip_on=False,
)
ax.invert_yaxis()
ax.set_xlabel("distance (m)")
ax.set_ylabel("depth (m)")
ax.set_title("Synthetic resistivity section (station markers on top)")
fig.colorbar(im, ax=ax, label=r"$\log_{10}\rho$  ($\Omega\cdot$m)")

<matplotlib.colorbar.Colorbar object at 0x7f2e07d1db50>

1-D resistivity-depth profiles

PlotResistivityDepthProfile extracts the sounding beneath a station and plots resistivity against depth — the most direct way to pick layer boundaries. Passing a ResistivityModel shows the raw log; later examples pass an EMHydroResult to add zone shading.

from pycsamt.interp.plot import PlotResistivityDepthProfile

PlotResistivityDepthProfile(rm).plot()

<Figure size 450x900 with 1 Axes>

Reading it. The profile drops sharply into the aquifer (~40 Ohm-m near 20 m), bottoms out in the clay (~12 Ohm-m), then rises steeply into the resistive basement (>1000 Ohm-m). Those four segments are exactly what the lithology and hydro examples turn into named units and aquifer zones.