geopulse.earth

Earth conductivity models and the surface-impedance abstraction.

Every EarthModel produces an Impedance object. The efield module calls impedance.apply(Bx_f, By_f) and never branches on model dimensionality — polymorphism is the whole point.

class geopulse.earth.CoastalCorrection2D(land_model, ocean_model, distance_from_coast_m, tm_decay_ratio=0.7)[source]

Bases: EarthModel

Semi-analytic coast-effect 2-D Earth model.

Interpolates between two 1-D reference models (land and ocean) via a tanh-shaped coupling function of distance from the coast, producing a TensorImpedance at each frequency.

Parameters:
  • land_model (Layered1D) – Layered 1-D model for the land side (well inland). Provides the far-field d +∞ asymptote.

  • ocean_model (Layered1D) – Layered 1-D model for the ocean side (well offshore). Typically a thin high-conductivity seawater layer over the same deep conductivity structure. Provides the d −∞ asymptote.

  • distance_from_coast_m (float) – Distance of the observation site from the coastline, in metres. Positive inland, negative seaward. 0 sits exactly at the coast.

  • tm_decay_ratio (float) – Ratio of TM to TE decay length L_TM / L_TE. Smaller → TM transitions more sharply than TE. Default 0.7 matches the Boteler-Pirjola (1998) observation that TM is more perturbed than TE near the interface.

Notes

tier == 2. Coordinate frame: coast along y (east); x (north) points inland. The tensor is anti-diagonal.

Examples

>>> import numpy as np
>>> from geopulse.earth.base import ConductivityLayer
>>> from geopulse.earth.layered_1d import Layered1D
>>> land = Layered1D([ConductivityLayer(np.inf, 0.001)])
>>> ocean = Layered1D([
...     ConductivityLayer(1_000.0, 3.3),
...     ConductivityLayer(np.inf, 0.001),
... ])
>>> coast = CoastalCorrection2D(land, ocean, distance_from_coast_m=0.0)
>>> imp = coast.compute_impedance(np.array([1e-3, 1e-2, 1e-1]))
>>> imp.rank
2
property tier: int

Return 2.

compute_impedance(freqs_Hz)[source]

Compute the surface impedance tensor per the coast-effect model.

Parameters:

freqs_Hz (ndarray) – Frequencies in Hz. Shape (n_freqs,). Must be non-negative; DC returns a zero tensor (no plane-wave response at ω = 0).

Return type:

Impedance

Returns:

TensorImpedance – Anti-diagonal 2×2 impedance tensor per frequency.

class geopulse.earth.ConductivityLayer(thickness_m, conductivity_Sm)[source]

Bases: object

A single layer in a 1-D conductivity model.

thickness_m

Layer thickness in meters. Use numpy.inf for the terminating half-space (bottom layer).

Type:

float

conductivity_Sm

Electrical conductivity in S/m.

Type:

float

Notes

Resistivity ρ = 1/σ is NOT stored — compute it when needed to avoid redundancy and potential inconsistency.

Raises:

ValueError – If conductivity_Sm is non-positive, or thickness_m is non-positive and not infinite.

Parameters:
thickness_m: float
conductivity_Sm: float
class geopulse.earth.EarthModel[source]

Bases: ABC

Abstract base class for Earth conductivity models.

Subclasses represent different spatial complexity:

  • Layered1D — horizontally layered half-space (Wait recursion).

  • Structured2D — 2-D cross-section (finite-difference).

  • Unstructured3D — full volumetric (ModEM / SimPEG wrapper).

All produce an Impedance object via compute_impedance().

abstractmethod compute_impedance(freqs_Hz)[source]

Compute the surface impedance at the given frequencies.

Parameters:

freqs_Hz (ndarray) – Frequency array in Hz. Shape: (n_freqs,). Must be non-negative. DC (0 Hz) is allowed but may return NaN because the plane-wave impedance is undefined at ω = 0.

Return type:

Impedance

Returns:

Impedance – An Impedance subclass matching this model’s dimensionality.

abstract property tier: int

1 for 1-D, 2 for 2-D, 3 for 3-D.

Type:

Dimensionality tier

classmethod from_library(name)[source]

Load a named model from the built-in library.

Parameters:

name (str) – Model name. See geopulse.earth.library.list_models() for the registry.

Return type:

EarthModel

Returns:

EarthModel – The requested model.

Raises:

KeyError – If name is not registered.

class geopulse.earth.Impedance(freqs_Hz)[source]

Bases: ABC

Abstract base class for surface impedance.

Parameters:

freqs_Hz (ndarray) – Frequency array in Hz. Cast internally to float64. Shape: (n_freqs,).

freqs_Hz

The stored frequency array.

Type:

numpy.ndarray

abstractmethod apply(Bx_f, By_f)[source]

Apply the impedance to frequency-domain B-field, returning E-field.

Parameters:
  • Bx_f (ndarray) – FFT of Bx in Tesla. Shape: (n_freqs,).

  • By_f (ndarray) – FFT of By in Tesla. Shape: (n_freqs,).

Return type:

tuple[ndarray, ndarray]

Returns:

  • Ex_f (numpy.ndarray) – FFT of Ex in V/m. Shape matches subclass semantics.

  • Ey_f (numpy.ndarray) – FFT of Ey in V/m. Shape matches subclass semantics.

abstract property rank: int

0 for scalar (1-D), 2 for tensor (2-D), 4 for kernel (3-D).

Type:

Rank marker

abstractmethod to_hdf5(group)[source]

Serialise impedance data to an HDF5 group.

Return type:

None

Parameters:

group (Group)

abstractmethod classmethod from_hdf5(group)[source]

Deserialise impedance data from an HDF5 group.

Return type:

Impedance

Parameters:

group (Group)

class geopulse.earth.ScalarImpedance(freqs_Hz, Z_values)[source]

Bases: Impedance

1-D scalar impedance Z(ω).

For a horizontally layered Earth under plane-wave excitation, the surface impedance is a single complex number per frequency. This is the standard magnetotelluric relation for the 1-D case.

Parameters:
  • freqs_Hz (ndarray) – Frequency array in Hz. Shape: (n_freqs,).

  • Z_values (ndarray) – Complex impedance values in Ohms. Shape: (n_freqs,). Cast to complex128 on construction.

Raises:

ShapeMismatchError – If Z_values.shape does not equal freqs_Hz.shape.

Notes

The applied relation is [1]_:

E_x(ω) =  (Z(ω) / μ₀) · B_y(ω)
E_y(ω) = -(Z(ω) / μ₀) · B_x(ω)

Examples

>>> import numpy as np
>>> freqs = np.array([1e-3, 1e-2, 1e-1])
>>> Z = np.array([1+1j, 3+3j, 10+10j])
>>> imp = ScalarImpedance(freqs, Z)
>>> Ex, Ey = imp.apply(np.ones(3), np.ones(3))
>>> Ex.shape
(3,)
apply(Bx_f, By_f)[source]

Apply the 1-D plane-wave MT relation E = (Z/μ₀) · B_perp.

Parameters:
  • Bx_f (ndarray) – FFT of the horizontal B-field components in Tesla. Shape: (n_freqs,).

  • By_f (ndarray) – FFT of the horizontal B-field components in Tesla. Shape: (n_freqs,).

Return type:

tuple[ndarray, ndarray]

Returns:

Ex_f, Ey_f (numpy.ndarray) – FFT of the horizontal E-field components in V/m. Shape: (n_freqs,).

Raises:

ShapeMismatchError – If Bx_f or By_f do not have the same length as freqs_Hz.

property rank: int

Return 0 — scalar per frequency.

to_hdf5(group)[source]

Serialise to an HDF5 group with a schema-version attribute.

Layout:

group/
  attrs:
    impedance_type = "scalar_1d"
    schema_version = 1
  freqs_Hz         dataset
  Z_values         dataset (complex128)
Return type:

None

Parameters:

group (Group)

classmethod from_hdf5(group)[source]

Read a ScalarImpedance from an HDF5 group.

Return type:

ScalarImpedance

Parameters:

group (Group)

class geopulse.earth.TensorImpedance(freqs_Hz, Z_tensor)[source]

Bases: Impedance

2-D magnetotelluric impedance tensor Z̃(ω) — 2×2 per frequency.

Applied relation (per frequency ω_i):

Ex(ω_i) = (Z[i,0,0] · Bx(ω_i) + Z[i,0,1] · By(ω_i)) / μ₀
Ey(ω_i) = (Z[i,1,0] · Bx(ω_i) + Z[i,1,1] · By(ω_i)) / μ₀

Shape convention: Z_tensor[i, :, :] is the 2×2 tensor at freqs_Hz[i]. A 1-D isotropic Earth is the special case Z_tensor[i] = [[0, Z], [-Z, 0]].

Parameters:
  • freqs_Hz (ndarray) – Frequency array in Hz. Shape (n_freqs,).

  • Z_tensor (ndarray) – Complex impedance tensor. Shape (n_freqs, 2, 2).

Raises:

ShapeMismatchError – If Z_tensor.shape does not equal (n_freqs, 2, 2).

apply(Bx_f, By_f)[source]

Apply the 2×2 tensor relation E(ω) = Z̃(ω) · B(ω) / μ₀.

Parameters:
  • Bx_f (ndarray) – FFT of horizontal B-field components in Tesla. Shape (n_freqs,).

  • By_f (ndarray) – FFT of horizontal B-field components in Tesla. Shape (n_freqs,).

Return type:

tuple[ndarray, ndarray]

Returns:

Ex_f, Ey_f (numpy.ndarray) – FFT of horizontal E-field components in V/m. Shape (n_freqs,).

Raises:

ShapeMismatchError – If either input does not match freqs_Hz in shape.

property rank: int

Return 2 — 2×2 tensor per frequency.

to_hdf5(group)[source]

Serialise to an HDF5 group with a schema-version attribute.

Layout:

group/
  attrs:
    impedance_type = "tensor_2d"
    schema_version = 1
  freqs_Hz      dataset
  Z_tensor      dataset (complex128, shape (n_freq, 2, 2))
Return type:

None

Parameters:

group (Group)

classmethod from_hdf5(group)[source]

Read a TensorImpedance from an HDF5 group.

Return type:

TensorImpedance

Parameters:

group (Group)

class geopulse.earth.KernelImpedance(freqs_Hz, Z_kernel)[source]

Bases: Impedance

3-D volumetric-response kernel Z̃(r, r', ω) — STUB (WP3).

Shape convention: (n_freqs, n_r, n_r, 2, 2). In practice the kernel will be sparse or block-structured and use chunked HDF5 storage.

Parameters:
  • freqs_Hz (np.ndarray)

  • Z_kernel (np.ndarray)

apply(Bx_f, By_f)[source]

Not yet implemented (WP3).

Return type:

tuple[ndarray, ndarray]

Parameters:
property rank: int

Return 4 — spatial kernel per frequency.

to_hdf5(group)[source]

Not yet implemented (WP3).

Return type:

None

Parameters:

group (Group)

classmethod from_hdf5(group)[source]

Not yet implemented (WP3).

Return type:

KernelImpedance

Parameters:

group (Group)

Modules

base

Abstract base class for Earth conductivity models.

impedance

Surface-impedance abstraction — the polymorphism spine of GeoPulse.

layered_1d

1-D layered-Earth surface impedance via Wait (1954) recursion.

library

Built-in 1-D Earth-model registry.

structured_2d

2-D structured Earth models.

unstructured_3d

3-D unstructured Earth model (ModEM / SimPEG wrapper) — STUB (WP3).