PyHydroGeophysX package#

Subpackages#

Submodules#

PyHydroGeophysX.em1d module#

Compatibility facade for PyHydroGeophysX.workflows.em1d.

exception PyHydroGeophysX.em1d.BackendUnavailable[source]#

Bases: RuntimeError

Raised when an optional numerical backend cannot be used.

PyHydroGeophysX.em1d.backend_status(method: str | None = None) → Dict[str, Any][source]#

Report whether the requested EM forward/inversion backend is usable.

The check imports the same method-specific forward class used by inversion. This prevents the UI and AQUAH from announcing a background inversion that cannot start because SimPEG or one of its runtime dependencies is missing.

PyHydroGeophysX.em1d.build_em_config(method: str, model: Dict[str, Any], geom: Dict[str, Any], inv: Dict[str, Any]) → Dict[str, Any][source]#

PyHydroGeophysX.em1d.calibrate_to_reference(path: str, method: str, geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], ref_resistivity: float, *, max_probe: int = 6, log: ~typing.Callable[[str], None] = <function noop>) → float[source]#

Find the data_scale that makes the recovered near-surface resistivity match a known/expected value.

The amplitude scale and the absolute resistivity level are degenerate — the EM data alone cannot fix the level (any data_scale fits the data, with the resistivity shifting to compensate). This breaks the degeneracy with EXTERNAL information: the user supplies ref_resistivity (a known background, e.g. from a borehole or regional geology), and this inverts a few probe soundings at two trial scales, fits the near-surface resistivity’s log-linear response to the scale, and solves for the scale that yields ref_resistivity. Returns the current data_scale unchanged if calibration is not possible.

PyHydroGeophysX.em1d.estimate_data_scale(path: str, method: str, geom: ~typing.Dict[str, ~typing.Any], *, max_soundings: int = 8, log: ~typing.Callable[[str], None] = <function noop>) → float[source]#

Estimate the amplitude calibration (data_scale) for normalized data.

Normalized airborne responses (e.g. moment-normalized dB/dt) differ from the studio’s 1D forward by a near-constant amplitude factor. This fits each sounding’s decay SHAPE to a grid of half-space forward responses at the current geometry and takes the geometric-mean amplitude ratio forward/observed at the best-fitting resistivity. Returns 1.0 if it cannot be estimated (so the caller can fall back to no scaling).

PyHydroGeophysX.em1d.example_catalog() → Dict[str, Dict[str, Any]][source]#

Return the desktop EM examples and their documented settings.

Prefer the source checkout for development and fall back to the compact package-data copy installed from a wheel.

PyHydroGeophysX.em1d.fdem_forward(model: ~typing.Dict[str, ~typing.Any], geom: ~typing.Dict[str, ~typing.Any], log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#: 1D FDEM forward response (secondary field, real/imag per frequency).

PyHydroGeophysX.em1d.fdem_invert(data: ~typing.Dict[str, ~typing.Any], geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#: Invert an FDEM sounding for a layered resistivity model (Occam 1D).

PyHydroGeophysX.em1d.invert_line(path: str, method: str, geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], *, spacing: float = 50.0, positions: ~numpy.ndarray | None = None, heights: ~numpy.ndarray | None = None, max_soundings: int = 12, lines: ~typing.Sequence[int] | None = None, doi_blank: bool = True, doi_factor: float = 0.5, ref_resistivity: float = 0.0, out_dir: ~pathlib.Path | None = None, initial_models: ~numpy.ndarray | None = None, log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#

Invert a line on a shared fixed-layer grid.

inv["lci_mode"] selects how the soundings are coupled:

simultaneous (the default whenever lateral_smoothness is positive): Solves the whole line as one system, with the lateral constraint part of what is being minimized. See PyHydroGeophysX.inversion.em1d_lci.
sequential: The older block-coordinate passes: each station is re-inverted on its own against the distance-weighted model its neighbours had at the end of the previous pass. Kept because it needs no analytic Jacobian, so it still runs against a forward operator that cannot supply one.
off: Independent 1D inversion per sounding, no lateral coupling.

The models are laid side by side to form a resistivity(position, depth) section ready for Model3DView.show_model(): edges = (ex, ey, ez) (ez is elevation, increasing upward) and model3d of shape (n_pos, 1, n_depth).

positions gives the along-line distance of each sounding (the section x-axis); heights overrides the sensor height per sounding. When doi_blank is set, cells below a per-sounding depth of investigation (a diffusion-depth estimate scaled by doi_factor) are blanked (NaN) so the unconstrained deep part of an early-time sounding is not shown as railed.

inv["robust_errors"] retains all imported gates and iteratively inflates effective errors for large residuals. It overrides hard rejection. The main chi2 uses ORIGINAL errors; result["robust"] records effective errors and a separate effective chi2. Import-time flags and QC still apply.

inv["auto_lambda"] re-solves the line at other smoothness weights to reach target_chi2. inv["reject_outliers"] drops the gates the converged model cannot explain (beyond outlier_threshold sigma, over outlier_passes cycles, never below min_data_fraction of the gates) and solves again; what it removed is reported under result["outliers"]. They address different causes, so they can be used together: relaxing the smoothness helps when the model is too stiff for the data, rejection helps when a minority of gates are simply wrong.

lines restricts the run to the named survey lines, so a line whose data is thinner than the rest can be given its own settings instead of one set having to suit every line. Passing None runs from the first station, as before. The lateral constraint already groups by line, so a line inverted on its own is tied exactly as it would be inside a whole-survey run; what changes is which settings reach it, and that the other lines are not re-solved. max_soundings then counts within the selection.

Stations arrive ordered by line, so a selection is a contiguous block. A set of lines that is not contiguous is refused rather than quietly widened to the span that encloses it, which would invert the lines in between under settings chosen for their neighbours.

PyHydroGeophysX.em1d.is_temcompany_source(path: str) → bool[source]#

Return whether path looks like a TEMcompany/TEM2Go export.

Both complete project directories and the self-describing *.xyz exports written by TEMImage are accepted.

PyHydroGeophysX.em1d.is_ttem_source(path: str) → bool[source]#: Return whether path contains a TEMcompany tTEM raw acquisition.

PyHydroGeophysX.em1d.load_line_geometry(path: str) → Dict[str, Any][source]#: Load per-sounding line geometry: along-line positions (m), optional sensor heights (m), and optional map coordinates x/y (e.g. easting/northing) for plan-view depth slices. Recognizes header names (distance/position for the position; alt/height for the height; easting/northing for the map coordinates, which also derive the distance when no distance column is present). A header-less file is read by column order (1 column = position; 2+ = position, height). positions is shifted to start at 0.

PyHydroGeophysX.em1d.load_sounding(path: str, method: str, sounding: int = 0, *, moment: str = 'HM', use_flags: bool = True, max_relative_std: float | None = None, gate_rejection: str = 'truncate', reject_negative: bool = False, min_gates_per_moment: Mapping[str, int] | None = None, ttem_loop_area: float | None = None, ttem_gex_path: str | None = None, ttem_tfi_path: str | None = None) → Dict[str, Any][source]#

Load one sounding from a sounding file.

The first column is the abscissa (FDEM: frequency Hz; TDEM: time s). The remaining columns hold the response(s) — a single sounding, or several stacked side by side so one file can carry a whole survey line (common for airborne EM exports). sounding picks which one (0-based):

TDEM: each extra column is one sounding’s response → column 1 + sounding.
FDEM: response columns come in (real, imag) pairs → one sounding is the pair starting at 1 + 2*sounding; a lone trailing real column gives imag = 0.

The returned dict also reports n_soundings so the caller can offer a picker. use_flags applies only to TEMcompany project databases; see load_temcompany_sounding().

PyHydroGeophysX.em1d.load_temcompany_sounding(path: str, sounding: int = 0, moment: str = 'HM', *, use_flags: bool = True, max_relative_std: float | None = None, gate_rejection: str = 'truncate', reject_negative: bool = False, min_gates_per_moment: Mapping[str, int] | None = None) → Dict[str, Any][source]#

Load a TEMcompany/TEM2Go sounding from a project folder or XYZ export.

The defaults reproduce the gate selection the project itself records, which was measured rather than assumed. Over 1,503 station-moment datasets of one project, the gates the stored inversion used (InversionModel.Datasets) are exactly the gates whose InUseFlags are set and whose value is finite and not a dummy, with 100 percent agreement. There is no further sign test: 87 low- moment and 251 high-moment datasets keep a non-positive gate. There is no further error cut either, and none is needed, because the largest relative error among the kept gates is exactly 0.250, so TEMImage applied that cut upstream when it wrote the flags. The selection is not even contiguous, so nor truncation: only 36 percent of the high-moment selections are a single run of gates.

So max_relative_std=None and reject_negative=False are the defaults. Both arguments remain, because a survey whose flags were written by an older an older acquisition release, or one being deliberately treated more strictly, still needs them. Note that the sign test only runs when max_relative_std is set: it condemns a gate alongside a noisy one rather than on its own.

use_flags=False ignores the project’s in-use flags and returns every gate with a finite, non-dummy value. Only a project database records those flags, so it makes no difference to an XYZ export.

min_gates_per_moment applies only to a joint LM+HM read of a project database, which is the only path where dropping one moment still leaves a sounding; see _load_temcompany_joint_database().

PyHydroGeophysX.em1d.load_ttem_sounding(path: str, sounding: int = 0, moment: str = 'LM+HM', *, max_relative_std: float | None = None, stack_seconds: float = 2.0, loop_area: float | None = None, gex_path: str | None = None, tfi_path: str | None = None) → Dict[str, Any][source]#

Load and stack one sounding from a raw TEMcompany tTEM survey.

Data are normalized by measured transmitter current and the loop area from the selected GEX (or an explicit UI override). A selected TFI is applied to each sign-corrected transient sequence before stacking.

PyHydroGeophysX.em1d.model_arrays(model: Dict[str, Any]) → tuple[ndarray, ndarray, ndarray][source]#: Return (thicknesses, resistivity, conductivity) from a model dict.

PyHydroGeophysX.em1d.model_depth_profile(thicknesses: ndarray, resistivity: ndarray, pad: float = 20.0) → tuple[ndarray, ndarray][source]#: Step profile (depth, resistivity) for plotting a layered model.

PyHydroGeophysX.em1d.save_inversion(result: Dict[str, Any], out_dir: Path) → List[str][source]#: Save recovered model + data fit to npy/csv; return written paths.

PyHydroGeophysX.em1d.tdem_forward(model: ~typing.Dict[str, ~typing.Any], geom: ~typing.Dict[str, ~typing.Any], log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#: 1D TDEM forward response (dB/dt or H per time channel).

PyHydroGeophysX.em1d.tdem_invert(data: ~typing.Dict[str, ~typing.Any], geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#: Invert a TDEM sounding for a layered resistivity model (Occam 1D).

PyHydroGeophysX.em1d.tdem_joint_invert(data: ~typing.Dict[str, ~typing.Any], geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#: Invert all available LM/HM gates at one station for one shared 1D model.

PyHydroGeophysX.gravmag module#

Compatibility facade for PyHydroGeophysX.workflows.gravmag.

exception PyHydroGeophysX.gravmag.InversionBackendUnavailable[source]#

Bases: BackendUnavailable

SimPEG / discretize / a usable solver could not be imported.

PyHydroGeophysX.gravmag.backend_status() → Dict[str, Any][source]#: Report whether the SimPEG potential-field inversion stack is available.

PyHydroGeophysX.gravmag.build_gravmag_config(kind: str, settings: Dict[str, Any], bodies: List[Dict[str, Any]], field: Dict[str, Any] | None = None) → Dict[str, Any][source]#

PyHydroGeophysX.gravmag.extract_profile(grid: Dict[str, ndarray], p1: Sequence[float], p2: Sequence[float], n: int = 200) → Dict[str, ndarray][source]#: Sample a gridded field along the line p1 -> p2 (bilinear).

PyHydroGeophysX.gravmag.forward_bodies(xobs: ~numpy.ndarray, yobs: ~numpy.ndarray, kind: str, bodies: ~typing.List[~typing.Dict[str, ~typing.Any]], field: ~typing.Dict[str, ~typing.Any] | None = None, log: ~typing.Callable[[str], None] = <function noop>) → ndarray[source]#: Sum the anomaly of a list of bodies. kind = ‘gravity’ or ‘magnetics’.

PyHydroGeophysX.gravmag.gravity_prism(xobs: ndarray, yobs: ndarray, body: Dict[str, Any]) → ndarray[source]#: Vertical gravity (mGal) of a right rectangular prism (Nagy 1966). z down.

PyHydroGeophysX.gravmag.gravity_sphere(xobs: ndarray, yobs: ndarray, body: Dict[str, Any]) → ndarray[source]#: Vertical gravity (mGal) of a buried sphere. z positive down, obs at z=0.

PyHydroGeophysX.gravmag.grid_data(x: ndarray, y: ndarray, value: ndarray, nx: int = 120, ny: int = 120, method: str = 'linear') → Dict[str, ndarray][source]#: Grid scattered station values onto a regular map. Returns xx, yy, zz.

PyHydroGeophysX.gravmag.invert_gravmag(x, y, value, kind: str, *, z: ~numpy.ndarray | None = None, field: ~typing.Dict[str, ~typing.Any] | None = None, detrend: int = 0, n_xy: int = 22, n_z: int = 12, max_iterations: int = 20, beta0_ratio: float = 1.0, max_stations: int = 600, relative_error: float = 0.03, noise_floor: float | None = None, solver: str = 'simpeg', auto_beta: bool = True, target_chi2: float = 1.0, chi2_tolerance: float = 0.2, max_beta_trials: int = 6, sensitivity_power: float = 1.0, out_dir: str | None = None, random_seed: int | None = 42, log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#

Run a SimPEG 3D potential-field inversion under the survey.

gravity recovers a density-contrast model (g/cc); magnetics recovers a susceptibility model (SI) and needs field = {inclination, declination, strength_nT}. z is optional per-station elevation (m, positive upward); a missing value falls back to 1 m. detrend (0..3) removes a polynomial regional trend before inversion. The returned grid uses elevation increasing upward. random_seed makes SimPEG’s eigenvalue-based beta estimate reproducible. Raises InversionBackendUnavailable if SimPEG is missing.

PyHydroGeophysX.gravmag.magnetic_dipole(xobs: ndarray, yobs: ndarray, body: Dict[str, Any], field: Dict[str, Any]) → ndarray[source]#: Total-field magnetic anomaly (nT) of an induced/magnetized sphere (a dipole).

PyHydroGeophysX.gravmag.qc_products(x: ndarray, y: ndarray, value: ndarray, *, detrend: int = 1, nx: int = 120, ny: int = 120) → Dict[str, Any][source]#: Calculate observed, regional and residual products for map/profile QC.

PyHydroGeophysX.gravmag.regional_residual(x: ndarray, y: ndarray, value: ndarray, degree: int = 1) → Tuple[ndarray, ndarray][source]#: Fit a polynomial regional trend of degree (1..3); return (regional, residual).

PyHydroGeophysX.gravmag.save_grid(grid: ~typing.Dict[str, ~numpy.ndarray], out_dir: ~pathlib.Path, name: str = 'anomaly', log: ~typing.Callable[[str], None] = <function noop>) → List[str][source]#: Save a grid to npy + CSV + VTK (best-effort). Return written paths.

PyHydroGeophysX.gravmag.spatially_balanced_indices(x: ndarray, y: ndarray, max_stations: int) → ndarray[source]#

Return deterministic farthest-point indices for a spatially balanced subset.

The previous evenly spaced file-row selection could over-sample a survey segment when input rows were ordered by flight line or acquisition order. Farthest-point selection starts near the survey centroid and repeatedly adds the station furthest from the selected set, preserving map coverage without a random seed.

PyHydroGeophysX.gui module#

Local GUI launcher for PyHydroGeophysX.

PyHydroGeophysX.gui.main(argv: Sequence[str] | None = None) → int[source]#

Launch the Streamlit GUI.

Parameters:: argv (sequence of str, optional) – Extra arguments passed to streamlit run.
Returns:: Streamlit exit code.
Return type:: int

PyHydroGeophysX.gui_mesh3d module#

Launcher for the 3D Mesh Builder Streamlit app.

PyHydroGeophysX.gui_mesh3d.main(argv: Sequence[str] | None = None) → int[source]#

Launch the 3D Mesh Builder GUI.

Parameters:: argv – Extra arguments forwarded to streamlit run.
Returns:: Streamlit exit code.
Return type:: int

PyHydroGeophysX.joint_api module#

Compatibility shim for PyHydroGeophysX.inversion.joint_api.

class PyHydroGeophysX.joint_api.JointInversionRequest(method_a: str, method_b: str, strategy: str, data: ~typing.Dict[str, ~typing.Any], parameters: ~typing.Dict[str, ~typing.Any] = <factory>, output_dir: str | ~pathlib.Path = 'results/joint_inversion', run_baseline: bool = True)[source]#

Bases: object

Input contract for a registered joint inversion runner.

data: Dict[str, Any]#

method_a: str#

method_b: str#

output_dir: str | Path = 'results/joint_inversion'#

parameters: Dict[str, Any]#

run_baseline: bool = True#

strategy: str#

class PyHydroGeophysX.joint_api.JointInversionResult(methods: ~typing.Tuple[str, str], strategy: str, models: ~typing.Dict[str, ~typing.Any] = <factory>, predicted: ~typing.Dict[str, ~typing.Any] = <factory>, coverage: ~typing.Dict[str, ~typing.Any] = <factory>, chi2: ~typing.Dict[str, float] = <factory>, history: ~typing.List[~typing.Dict[str, ~typing.Any]] = <factory>, baseline: ~typing.Dict[str, ~typing.Any] = <factory>, artifacts: ~typing.Dict[str, str] = <factory>, warnings: ~typing.List[str] = <factory>, meta: ~typing.Dict[str, ~typing.Any] = <factory>, status: str = 'success')[source]#

Bases: object

Method-neutral result returned by all registered joint runners.

artifacts: Dict[str, str]#

baseline: Dict[str, Any]#

chi2: Dict[str, float]#

coverage: Dict[str, Any]#

history: List[Dict[str, Any]]#

meta: Dict[str, Any]#

methods: Tuple[str, str]#

models: Dict[str, Any]#

predicted: Dict[str, Any]#

status: str = 'success'#

strategy: str#

summary() → Dict[str, Any][source]#: Return a compact JSON-ready representation without model arrays.

warnings: List[str]#

class PyHydroGeophysX.joint_api.JointMethodAdapter(*args, **kwargs)[source]#

Bases: Protocol

Interface implemented by a capability-specific request runner.

capability: JointPairCapability#

run(request: JointInversionRequest) → JointInversionResult[source]#: Execute a validated request and return method-neutral results.

class PyHydroGeophysX.joint_api.JointPairCapability(methods: Tuple[str, str], strategies: Mapping[str, str], dimension: str, model_parameter: str, implemented: bool, dependencies: Tuple[str, ...] = (), description: str = '', runner: str | None = None, backends: Tuple[str, ...] = ())[source]#

Bases: object

Describe the strategies available for one normalized method pair.

backends: Tuple[str, ...] = ()#

dependencies: Tuple[str, ...] = ()#

description: str = ''#

dimension: str#

implemented: bool#

methods: Tuple[str, str]#

model_parameter: str#

runner: str | None = None#

strategies: Mapping[str, str]#

PyHydroGeophysX.joint_api.get_joint_capabilities(include_planned: bool = True) → List[JointPairCapability][source]#: List implemented capabilities and, optionally, planned method pairs.

PyHydroGeophysX.joint_api.get_joint_capability(method_a: str, method_b: str) → JointPairCapability[source]#: Return one capability, including a planned placeholder if unsupported.

PyHydroGeophysX.joint_api.normalize_joint_pair(method_a: str, method_b: str) → Tuple[str, str][source]#: Return a stable pair key and reject duplicate methods.

PyHydroGeophysX.joint_api.pair_joint_soundings(f_value: Any, t_value: Any, parameters: Mapping[str, Any]) → List[Tuple[int, int, str]][source]#: Pair FDEM/TDEM soundings by coordinate or, when safe, by index.

PyHydroGeophysX.joint_api.split_joint_soundings(value: Any) → Tuple[List[Mapping[str, Any]], ndarray | None][source]#: Normalize one EM sounding or a line-of-soundings container.

PyHydroGeophysX.joint_api.validate_profile_interface(interface: Any, sensor_x: Any, sensor_z: Any) → Tuple[ndarray, ndarray][source]#: Validate a 2-D structural interface against an ERT surface profile.

PyHydroGeophysX.table_io module#

Compatibility shim for PyHydroGeophysX.data_processing.table_io.

PyHydroGeophysX.table_io.ensure_dir(path: str | Path) → Path[source]#: Create path and its parents if needed, then return it.

PyHydroGeophysX.table_io.load_2d_array(path: str | Path) → ndarray[source]#: Load an array from NPY, NPZ, CSV, TXT, or DAT input.

PyHydroGeophysX.table_io.load_xyz_table(path: str | Path, min_cols: int = 2) → ndarray[source]#: Load a two-dimensional table with at least min_cols columns.

PyHydroGeophysX.table_io.read_json(path: str | Path) → dict[str, Any] | None[source]#: Read JSON, returning None for a missing or malformed document.

PyHydroGeophysX.table_io.write_csv(path: str | Path, rows: Sequence[Sequence[Any]], header: Iterable[str] | None = None) → Path[source]#: Write rows to a CSV file.

PyHydroGeophysX.table_io.write_json(path: str | Path, obj: Any) → Path[source]#: Atomically write a JSON document.

Module contents#

Public, import-light API for PyHydroGeophysX.

The package root intentionally contains no eager scientific imports. Public objects are resolved on first access so importing a lightweight submodule does not also initialize plotting libraries or optional geophysics engines.

class PyHydroGeophysX.AgentCoordinator(api_key: str | None = None, output_dir: str = 'results/agents', llm_provider: str = 'openai')[source]#

Bases: object

Coordinates multiple agents to execute a complete workflow.

The coordinator manages cross-modal geophysical workflows such as: “load geophysical data → process → invert → convert to hydrologic parameters → report” with support for multiple data types (ERT, seismic, etc.) and LLM providers (GPT, Gemini, Claude).

execute_workflow(config: Dict[str, Any], dry_run: bool = False, resume: bool = False) → Any[source]#

Execute the complete workflow with registered agents.

Parameters:

config – Configuration dictionary containing: - data_file: Path to ERT data file - instrument: Instrument type (E4D, Syscal, etc.) - inversion_params: Parameters for inversion - petrophysical_params: Parameters for water content conversion - use_seismic: Whether to include seismic processing (default: False) - seismic_data: Optional seismic data file - use_climate: Whether to include climate data (default: False) - climate_config: Climate data configuration (coords/geometry, dates, etc.) - ert_timestamps: Timestamps for ERT acquisitions (for climate alignment)
dry_run – If True, return a preview plan without executing.
resume – If True, load checkpointed intermediate results and skip already-completed steps.

Returns:

Dictionary containing workflow results, or AgentResult when dry_run=True.

get_workflow_summary() → Dict[str, Any][source]#: Get summary of workflow execution.

preview_workflow(config: Dict[str, Any]) → AgentResult[source]#

Resolve and validate a workflow plan without running processing.

Parameters:: config (dict) – Workflow configuration. May include user_request or request for deterministic preview parsing.
Returns:: Preview result with resolved config, validation warnings, plan, and approximate LLM cost.
Return type:: AgentResult
Raises:: None –

Examples

>>> coordinator = AgentCoordinator(api_key=None)
>>> result = coordinator.preview_workflow({"data_file": "missing.ohm"})
>>> result["status"]
'failed'

register_agent(agent_name: str, agent_instance)[source]#

Register an agent with the coordinator.

Parameters:

agent_name – Unique identifier for the agent
agent_instance – Agent instance to register

class PyHydroGeophysX.ArtifactRef(artifact_id: str, kind: str, path: str, format: str, checksum: str = '', metadata: Mapping[str, ~typing.Any]=<factory>)[source]#

Bases: object

A serializable reference to a workflow input or output artifact.

artifact_id: str#

checksum: str = ''#

format: str#

classmethod from_dict(value: Mapping[str, Any]) → ArtifactRef[source]#

classmethod from_path(path: Path, *, artifact_id: str, kind: str, format: str | None = None, base_dir: Path | None = None, metadata: Mapping[str, Any] | None = None, checksum: bool = True) → ArtifactRef[source]#

kind: str#

metadata: Mapping[str, Any]#

path: str#

resolve(base_dir: Path) → Path[source]#

to_dict() → Dict[str, Any][source]#

exception PyHydroGeophysX.BackendUnavailable[source]#

Bases: RuntimeError

Raised when an optional numerical backend cannot be used.

class PyHydroGeophysX.BaseVelocityModel[source]#

Bases: object

Base class for seismic velocity models.

calculate_velocity(**kwargs) → ndarray[source]#

Calculate seismic velocity from rock properties.

Parameters:: **kwargs – Rock properties specific to each model
Returns:: Seismic velocity values (Vp, Vs, or both)

class PyHydroGeophysX.BrieModel(exponent: float = 3.0)[source]#

Bases: object

Brie’s model for calculating the effective bulk modulus of a partially saturated medium.

calculate_fluid_modulus(saturation: float, water_modulus: float = 2.0, gas_modulus: float = 0.01) → float[source]#

Calculate effective fluid bulk modulus using Brie’s equation.

Parameters:

saturation – Water saturation (0 to 1)
water_modulus – Bulk modulus of water (GPa, default: 2.0)
gas_modulus – Bulk modulus of gas (GPa, default: 0.01)

Returns:

Effective fluid bulk modulus (GPa)

calculate_saturated_modulus(dry_modulus: float, mineral_modulus: float, porosity: float, saturation: float, water_modulus: float = 2.0, gas_modulus: float = 0.01) → float[source]#

Calculate the saturated bulk modulus based on Brie’s equation.

Parameters:

dry_modulus – Bulk modulus of the dry rock (GPa)
mineral_modulus – Bulk modulus of the mineral matrix (GPa)
porosity – Porosity of the rock
saturation – Water saturation (0 to 1)
water_modulus – Bulk modulus of water (GPa, default: 2.0)
gas_modulus – Bulk modulus of gas (GPa, default: 0.01)

Returns:

Saturated bulk modulus (GPa)

class PyHydroGeophysX.DEMModel[source]#

Bases: BaseVelocityModel

Differential Effective Medium (DEM) model for calculating elastic properties and seismic velocities of porous rocks.

calculate_velocity(porosity: ndarray, saturation: ndarray, bulk_modulus: float, shear_modulus: float, mineral_density: float, aspect_ratio: float = 0.1) → Tuple[ndarray, ndarray, ndarray][source]#

Calculate P-wave velocity using the DEM model.

Parameters:

porosity – Porosity values (array)
saturation – Saturation values (array)
bulk_modulus – Initial bulk modulus of the solid matrix (GPa)
shear_modulus – Initial shear modulus of the solid matrix (GPa)
mineral_density – Density of the solid matrix (kg/m³)
aspect_ratio – Aspect ratio of pores (default: 0.1)

Returns:

Effective bulk modulus (GPa), effective shear modulus (GPa), and P-wave velocity (m/s)

class PyHydroGeophysX.ERTForwardModeling(mesh: pygimli.Mesh, data: pygimli.DataContainer | None = None)[source]#

Bases: object

Class for forward modeling of Electrical Resistivity Tomography (ERT) data.

create_synthetic_data(xpos: ndarray, ypos: ndarray | None = None, mesh: pygimli.Mesh | None = None, res_models: ndarray | None = None, schemeName: str = 'wa', noise_level: float = 0.05, absolute_error: float = 0.0, relative_error: float = 0.05, save_path: str | None = None, show_data: bool = False, seed: int | None = None, xbound: float = 100, ybound: float = 100) → Tuple[pygimli.DataContainer, pygimli.Mesh][source]#

Create synthetic ERT data using forward modeling.

This method simulates an ERT survey by placing electrodes, creating a measurement scheme, performing forward modeling to generate synthetic data, and adding noise.

Parameters:

xpos – X-coordinates of electrodes
ypos – Y-coordinates of electrodes (if None, uses flat surface)
mesh – Mesh for forward modeling
res_models – Resistivity model values
schemeName – Name of measurement scheme (‘wa’, ‘dd’, etc.)
noise_level – Level of Gaussian noise to add
absolute_error – Absolute error for data estimation
relative_error – Relative error for data estimation
save_path – Path to save synthetic data (if None, does not save)
show_data – Whether to display data after creation
seed – Random seed for noise generation
xbound – X boundary extension for mesh
ybound – Y boundary extension for mesh

Returns:

Tuple of (synthetic ERT data container, simulation mesh)

forward(resistivity_model: ndarray, log_transform: bool = True) → ndarray[source]#

Compute forward response for a given resistivity model.

Parameters:

resistivity_model – Resistivity model values
log_transform – Whether resistivity_model is log-transformed

Returns:

Forward response (apparent resistivity)

forward_and_jacobian(resistivity_model: ndarray, log_transform: bool = True) → Tuple[ndarray, ndarray][source]#

Compute forward response and Jacobian matrix.

Parameters:

resistivity_model – Resistivity model values
log_transform – Whether resistivity_model is log-transformed

Returns:

Tuple of (forward response, Jacobian matrix)

get_coverage(resistivity_model: ndarray, log_transform: bool = True) → ndarray[source]#

Compute coverage (resolution) for a given resistivity model.

Parameters:

resistivity_model – Resistivity model values
log_transform – Whether resistivity_model is log-transformed

Returns:

Coverage values for each cell

set_data(data: pygimli.DataContainer) → None[source]#

Set ERT data for forward modeling.

Parameters:: data – ERT data container

set_mesh(mesh: pygimli.Mesh) → None[source]#

Set mesh for forward modeling.

Parameters:: mesh – PyGIMLI mesh

class PyHydroGeophysX.ERTInversionAgent(api_key: str | None = None, model: str | None = None, llm_provider: str = 'openai')[source]#

Bases: BaseAgent

Agent specialized in ERT inversion.

Uses PyHydroGeophysX inversion module to perform resistivity inversion with optional structural constraints from seismic data.

execute(input_data: Dict[str, Any]) → Dict[str, Any][source]#

Perform ERT inversion (standard or time-lapse).

Parameters:: input_data – Dictionary containing: - ert_data: Loaded ERT data (for standard inversion) - inversion_mode: ‘standard’ or ‘time-lapse’ - time_lapse_data: List of ERT datasets (for time-lapse) - time_lapse_method: ‘difference’, ‘ratio’, or ‘joint’ (for time-lapse) - temporal_regularization: Temporal smoothing weight (for time-lapse) - inversion_params: Inversion parameters (lambda, max_iter, etc.) - use_structure_constraint: Whether to use seismic structure (default: False) - seismic_structure: Optional seismic structure data - output_dir: Directory for saving results
Returns:: Dictionary containing inversion results

class PyHydroGeophysX.ERTLoaderAgent(api_key: str | None = None, model: str | None = None, llm_provider: str = 'openai')[source]#

Bases: BaseAgent

Agent specialized in loading ERT data from various instruments.

Uses PyHydroGeophysX data_processing module to load, validate, and prepare ERT data for inversion.

execute(input_data: Dict[str, Any]) → Dict[str, Any][source]#

Load and process ERT data.

Parameters:: input_data – Dictionary containing: - data_file: Path to ERT data file - instrument: Instrument type (E4D, Syscal, ABEM, etc.) - project_dir: Project directory - crs: Coordinate reference system (‘local’ or EPSG code) - quality_check: Whether to perform quality checks (default: True)
Returns:: Dictionary containing loaded ERT data and quality metrics

class PyHydroGeophysX.FDEMForwardModeling(thicknesses: ndarray, survey_config: FDEMSurveyConfig | None = None, survey: simpeg.electromagnetics.frequency_domain.Survey | None = None)[source]#

Bases: object

Forward modeling of Frequency-Domain EM data using SimPEG.

Supports 1D layered-earth conductivity models.

forward(conductivity: ndarray) → ndarray[source]#: Compute FDEM response for a given conductivity model.

forward_with_noise(conductivity: ndarray, noise_level: float = 0.05, seed: int | None = None) → Tuple[ndarray, ndarray, ndarray][source]#: Compute noisy and clean FDEM responses with data uncertainties.

static hydro_to_fdem(water_content: ndarray, porosity: ndarray, layer_thicknesses: ndarray, **petro_params)[source]#: Convert hydrological properties to FDEM response via petrophysics.

class PyHydroGeophysX.FDEMSurveyConfig(source_location: ndarray = None, source_radius: float = 10.0, receiver_location: ndarray = None, receiver_orientation: str = 'z', receiver_component: str = 'secondary', frequencies: ndarray = None, waveform_type: str = 'dipole')[source]#

Bases: object

Configuration for FDEM survey geometry.

frequencies: ndarray = None#

receiver_component: str = 'secondary'#

receiver_location: ndarray = None#

receiver_orientation: str = 'z'#

source_location: ndarray = None#

source_radius: float = 10.0#

waveform_type: str = 'dipole'#

class PyHydroGeophysX.FirstBreakPick(source_id: int, receiver_id: int, time_s: float, source_x: float, source_z: float, receiver_x: float, receiver_z: float, field_record: int, trace_number: int, trace_index: int, amplitude: float)[source]#

Bases: object

One first-break pick.

amplitude: float#

field_record: int#

receiver_id: int#

receiver_x: float#

receiver_z: float#

source_id: int#

source_x: float#

source_z: float#

time_s: float#

to_dict() → Dict[str, Any][source]#: Return a CSV/JSON friendly representation.

trace_index: int#

trace_number: int#

class PyHydroGeophysX.GeophysicalInversion(method: str, **kwargs)[source]#

Bases: object

Unified factory for multi-method geophysical inversion.

Dispatches to the correct inversion engine.

SUPPORTED = {'ert', 'fdem', 'joint', 'joint_ert_srt', 'srt', 'tdem'}#

property engine#

run(**kwargs)[source]#

setup(**kwargs)[source]#

class PyHydroGeophysX.GeophysicalInversionAgent(api_key: str | None = None, model: str | None = None, llm_provider: str = 'openai')[source]#

Bases: BaseAgent

Agent for multi-method inversion orchestration.

Supports SRT and FDEM natively and delegates ERT to ERTInversionAgent.

execute(input_data: Dict[str, Any]) → Dict[str, Any][source]#

Execute the agent’s primary task.

Parameters:: input_data – Input data dictionary
Returns:: Dictionary containing execution results

class PyHydroGeophysX.GravityMagneticsJointResult(density: ~numpy.ndarray, susceptibility: ~numpy.ndarray, predicted_gravity: ~numpy.ndarray, predicted_magnetics: ~numpy.ndarray, coverage_gravity: ~numpy.ndarray, coverage_magnetics: ~numpy.ndarray, chi2_gravity: float, chi2_magnetics: float, convergence: ~typing.List[~typing.Dict[str, ~typing.Any]], edges: ~typing.Tuple[~numpy.ndarray, ~numpy.ndarray, ~numpy.ndarray], model_shape: ~typing.Tuple[int, int, int], cross_gradient: ~numpy.ndarray, baseline: ~typing.Dict[str, ~typing.Any] = <factory>, meta: ~typing.Dict[str, ~typing.Any] = <factory>)[source]#

Bases: object

Result of a shared-mesh gravity–magnetics inversion.

baseline: Dict[str, Any]#

chi2_gravity: float#

chi2_magnetics: float#

convergence: List[Dict[str, Any]]#

coverage_gravity: ndarray#

coverage_magnetics: ndarray#

cross_gradient: ndarray#

density: ndarray#

edges: Tuple[ndarray, ndarray, ndarray]#

meta: Dict[str, Any]#

model_shape: Tuple[int, int, int]#

predicted_gravity: ndarray#

predicted_magnetics: ndarray#

susceptibility: ndarray#

class PyHydroGeophysX.HertzMindlinModel(critical_porosity: float = 0.4, coordination_number: float = 4.0)[source]#

Bases: BaseVelocityModel

Hertz-Mindlin model and Hashin-Shtrikman bounds for seismic velocity in porous rocks.

calculate_velocity(porosity: ndarray, saturation: ndarray, bulk_modulus: float, shear_modulus: float, mineral_density: float, depth: float = 1.0) → Tuple[ndarray, ndarray][source]#

Calculate P-wave velocity for porous rocks.

Parameters:

porosity – Porosity values (array)
saturation – Saturation values (array)
bulk_modulus – Bulk modulus of the solid matrix (GPa)
shear_modulus – Shear modulus of the solid matrix (GPa)
mineral_density – Density of the solid matrix (kg/m³)
depth – Depth for pressure estimation (m, default: 1.0)

Returns:

Tuple of high-bound P-wave velocity (m/s) and low-bound P-wave velocity (m/s)

class PyHydroGeophysX.HydroModelOutput(model_directory: str)[source]#

Bases: ABC

Base class for all hydrological model outputs.

calculate_saturation(water_content: ndarray, porosity: float | ndarray) → ndarray[source]#

Calculate saturation from water content and porosity.

Parameters:

water_content – Water content array
porosity – Porosity value(s)

Returns:

Saturation array

abstractmethod get_timestep_info() → List[Tuple][source]#

Get information about each timestep.

Returns:: List of timestep information tuples

abstractmethod load_time_range(start_idx: int = 0, end_idx: int | None = None, **kwargs) → ndarray[source]#

Load data for a range of timesteps.

Parameters:

start_idx – Starting timestep index
end_idx – Ending timestep index (exclusive)
**kwargs – Additional parameters specific to the model type

Returns:

Data array for the specified timestep range

abstractmethod load_timestep(timestep_idx: int, **kwargs) → ndarray[source]#

Load data for a specific timestep.

Parameters:

timestep_idx – Index of the timestep to load
**kwargs – Additional parameters specific to the model type

Returns:

Data array for the specified timestep

class PyHydroGeophysX.JointGravityMagneticsInversion(gravity_data: Mapping[str, Any], magnetics_data: Mapping[str, Any], *, field: Mapping[str, Any] | None = None, n_xy: int = 12, n_z: int = 8, max_iterations: int = 10, max_stations: int = 600, gravity_relative_error: float = 0.03, magnetics_relative_error: float = 0.03, gravity_noise_floor: float = 0.5, magnetics_noise_floor: float = 2.0, gravity_weight: float = 1.0, magnetics_weight: float = 1.0, cross_gradient_weight: float = 2000000000000.0, beta0_ratio: float = 1.0, gravity_detrend: int = 0, magnetics_detrend: int = 0, run_baseline: bool = True, baseline_max_iterations: int | None = None, random_seed: int | None = 42, output_dir: str | None = None, progress_callback: Callable[[Dict[str, Any]], None] | None = None)[source]#

Bases: object

Jointly recover density and susceptibility with SimPEG cross-gradient coupling.

Inputs are mappings with x, y, value and optional z arrays. Gravity values use mGal and magnetic total-field values use nT. Both surveys are placed on one tensor mesh; the stacked inversion model is split by a simpeg.maps.Wires map into density (g/cc) and susceptibility (SI).

run() → GravityMagneticsJointResult[source]#: Run the native SimPEG similarity-measure inversion.

class PyHydroGeophysX.JointInversionRequest(method_a: str, method_b: str, strategy: str, data: ~typing.Dict[str, ~typing.Any], parameters: ~typing.Dict[str, ~typing.Any] = <factory>, output_dir: str | ~pathlib.Path = 'results/joint_inversion', run_baseline: bool = True)[source]#

Bases: object

Input contract for a registered joint inversion runner.

data: Dict[str, Any]#

method_a: str#

method_b: str#

output_dir: str | Path = 'results/joint_inversion'#

parameters: Dict[str, Any]#

run_baseline: bool = True#

strategy: str#

class PyHydroGeophysX.JointInversionResult(methods: ~typing.Tuple[str, str], strategy: str, models: ~typing.Dict[str, ~typing.Any] = <factory>, predicted: ~typing.Dict[str, ~typing.Any] = <factory>, coverage: ~typing.Dict[str, ~typing.Any] = <factory>, chi2: ~typing.Dict[str, float] = <factory>, history: ~typing.List[~typing.Dict[str, ~typing.Any]] = <factory>, baseline: ~typing.Dict[str, ~typing.Any] = <factory>, artifacts: ~typing.Dict[str, str] = <factory>, warnings: ~typing.List[str] = <factory>, meta: ~typing.Dict[str, ~typing.Any] = <factory>, status: str = 'success')[source]#

Bases: object

Method-neutral result returned by all registered joint runners.

artifacts: Dict[str, str]#

baseline: Dict[str, Any]#

chi2: Dict[str, float]#

coverage: Dict[str, Any]#

history: List[Dict[str, Any]]#

meta: Dict[str, Any]#

methods: Tuple[str, str]#

models: Dict[str, Any]#

predicted: Dict[str, Any]#

status: str = 'success'#

strategy: str#

summary() → Dict[str, Any][source]#: Return a compact JSON-ready representation without model arrays.

warnings: List[str]#

class PyHydroGeophysX.JointPairCapability(methods: Tuple[str, str], strategies: Mapping[str, str], dimension: str, model_parameter: str, implemented: bool, dependencies: Tuple[str, ...] = (), description: str = '', runner: str | None = None, backends: Tuple[str, ...] = ())[source]#

Bases: object

Describe the strategies available for one normalized method pair.

backends: Tuple[str, ...] = ()#

dependencies: Tuple[str, ...] = ()#

description: str = ''#

dimension: str#

implemented: bool#

methods: Tuple[str, str]#

model_parameter: str#

runner: str | None = None#

strategies: Mapping[str, str]#

class PyHydroGeophysX.MODFLOWPorosity(model_directory: str, model_name: str)[source]#

Bases: HydroModelOutput

Class for processing porosity data from MODFLOW simulations.

get_timestep_info() → List[Tuple][source]#

Get information about each timestep in the model. Returns a minimal placeholder since porosity doesn’t vary with time.

Returns:: List with single dummy timestep info

load_porosity() → ndarray[source]#

Load porosity data from MODFLOW model (supports both MODFLOW 6 and earlier versions).

Returns:: 3D array of porosity values (nlay, nrow, ncol)

load_time_range(start_idx: int = 0, end_idx: int | None = None, **kwargs) → ndarray[source]#

Load porosity for a range of timesteps. Since porosity is typically constant, this returns a stack of identical arrays.

Parameters:

start_idx – Starting timestep index (unused)
end_idx – Ending timestep index (unused)

Returns:

4D array of porosity values (nt, nlay, nrow, ncol) where all timesteps are identical

load_timestep(timestep_idx: int, **kwargs) → ndarray[source]#

Load porosity for a specific timestep. Note: For MODFLOW, porosity is typically constant over time, so this returns the same array regardless of timestep.

Parameters:: timestep_idx – Index of the timestep (unused)
Returns:: 3D array of porosity values

class PyHydroGeophysX.MODFLOWWaterContent(model_directory: str, idomain: ndarray)[source]#

Bases: HydroModelOutput

Class for processing water content data from MODFLOW simulations.

get_timestep_info() → List[Tuple[int, int, float, float]][source]#

Get information about each timestep in the WaterContent file.

Returns:: List of tuples (kstp, kper, pertim, totim) for each timestep

load_time_range(start_idx: int = 0, end_idx: int | None = None, nlay: int = 3) → ndarray[source]#

Load water content for a range of timesteps.

Parameters:

start_idx – Starting timestep index (default: 0)
end_idx – Ending timestep index (exclusive, default: None loads all)
nlay – Number of layers in the model (default: 3)

Returns:

Water content array with shape (timesteps, nlay, nrows, ncols)

load_timestep(timestep_idx: int, nlay: int = 3) → ndarray[source]#

Load water content for a specific timestep.

Parameters:

timestep_idx – Index of the timestep to load
nlay – Number of layers in the model

Returns:

Water content array with shape (nlay, nrows, ncols)

class PyHydroGeophysX.MeshCreator(quality: float = 28, area: float = 40)[source]#

Bases: object

Class for creating and managing meshes for geophysical inversion.

create_from_ert_data(data, max_depth: float = 30.0, quality: float = 34)[source]#

Create a mesh suitable for ERT inversion from ERT data.

Parameters:

data – PyGIMLI ERT data object
max_depth – Maximum depth of the mesh
quality – Mesh quality parameter

Returns:

PyGIMLI mesh for ERT inversion

create_from_layers(surface: ndarray, layers: List[ndarray], bottom_depth: float = 30.0, markers: List[int] = None) → pygimli.Mesh[source]#

Create a mesh from surface and layer boundaries.

Parameters:

surface – Surface coordinates [[x,z],…]
layers – List of layer boundary coordinates
bottom_depth – Depth below surface minimum for mesh bottom
markers – List of markers for each layer (default: [2, 3, 2, …])

Returns:

PyGIMLI mesh

class PyHydroGeophysX.ParflowPorosity(model_directory: str, run_name: str)[source]#

Bases: ParflowOutput

Processes porosity data from ParFlow simulations. Porosity in ParFlow is typically static (time-invariant) and stored in a single PFB file (e.g., <run_name>.out.porosity.pfb or similar).

get_timestep_info() → List[Tuple[int, float]][source]#

Returns timestep information, typically based on other ParFlow outputs (like saturation) as porosity itself is static.

Returns:

A list of (timestep_number, time_value) tuples,: derived from self.available_timesteps.

Return type:

List[Tuple[int, float]]

load_mask() → ndarray[source]#

Load the domain mask data from a ParFlow model. The mask file (.out.mask.pfb) indicates active (1) and inactive (0) cells.

Returns:

A 3D NumPy array representing the domain mask (nz, ny, nx).: Values are typically 0 or 1.

Return type:

np.ndarray

Raises:

FileNotFoundError – If no standard mask PFB file can be found.
ValueError – If there’s an error reading or processing the PFB file.

load_porosity() → ndarray[source]#

Load the static porosity data from the ParFlow model.

It searches for common ParFlow porosity filename patterns within the model directory.

Returns:

A 3D NumPy array of porosity values (nz, ny, nx).

Return type:

np.ndarray

Raises:

FileNotFoundError – If no standard porosity PFB file can be found.
ValueError – If there’s an error reading or processing the PFB file.

load_time_range(start_idx: int = 0, end_idx: int | None = None, **kwargs: Any) → ndarray[source]#

Load porosity data for a conceptual range of timesteps. Since porosity is time-invariant, this method returns a 4D array where the static 3D porosity data is repeated along the time axis.

The number of repetitions along the time axis (nt) is determined by the length of self.available_timesteps (discovered from saturation/pressure files) if end_idx is None, or by min(end_idx - start_idx, len(available_timesteps)). A minimum of 1 repetition is ensured if any timesteps are notionally available.

Parameters:

start_idx (int, optional) – Starting timestep index (used to determine nt). Defaults to 0.
end_idx (Optional[int], optional) – Ending timestep index (exclusive, used for nt). Defaults to None (use all available timesteps).
**kwargs (Any) – Additional keyword arguments (not used).

Returns:

A 4D NumPy array of porosity values (nt, nz, ny, nx).: All slices along the time dimension are identical.

Return type:

np.ndarray

load_timestep(timestep_idx: int, **kwargs: Any) → ndarray[source]#

Load porosity data. For ParFlow, porosity is typically time-invariant. This method returns the static porosity array, ignoring timestep_idx.

Parameters:

timestep_idx (int) – Index of the timestep (ignored, as porosity is static).
**kwargs (Any) – Additional keyword arguments (not used).

Returns:

A 3D NumPy array of porosity values (nz, ny, nx).

Return type:

np.ndarray

class PyHydroGeophysX.ParflowSaturation(model_directory: str, run_name: str)[source]#

Bases: ParflowOutput

Processes saturation data from ParFlow simulations (.out.satur.*.pfb files).

get_timestep_info() → List[Tuple[int, float]][source]#

Provides information about available ParFlow timesteps.

For ParFlow, the timestep number from the filename often directly corresponds to the simulation time (e.g., if output is every 1 hour, timestep 24 is 24 hours). This method returns a list of tuples: (timestep_number, simulation_time). Currently, simulation_time is simply cast from timestep_number. More accurate time mapping would require parsing ParFlow timing files if complex.

Returns:

A list where each tuple is (timestep_number, time_value).: Time_value is float representation of timestep_number.

Return type:

List[Tuple[int, float]]

load_time_range(start_idx: int = 0, end_idx: int | None = None, **kwargs: Any) → ndarray[source]#

Load saturation data for a specified range of zero-based timestep indices.

Parameters:

start_idx (int, optional) – Starting zero-based timestep index. Defaults to 0.
end_idx (Optional[int], optional) – Ending zero-based timestep index (exclusive). If None, loads up to the last available timestep. Defaults to None.
**kwargs (Any) – Additional keyword arguments (not used).

Returns:

A 4D NumPy array of saturation values (num_timesteps, nz, ny, nx).: Returns an empty 4D array if the range is invalid or no data is found.

Return type:

np.ndarray

Raises:

ValueError – If no timesteps are available, or if the specified range is invalid (e.g., start_idx out of bounds, end_idx <= start_idx leading to empty range).

load_timestep(timestep_idx: int, **kwargs: Any) → ndarray[source]#

Load saturation data for a specific, zero-based timestep index.

Parameters:

timestep_idx (int) – The zero-based index of the timestep to load from the list of available timesteps discovered during initialization.
**kwargs (Any) – Additional keyword arguments (not used by this method).

Returns:

A 3D NumPy array of saturation values (nz, ny, nx).

Return type:

np.ndarray

Raises:

ValueError – If no timesteps are available or if timestep_idx is out of range.

class PyHydroGeophysX.ProfileInterpolator(point1: List[int], point2: List[int], surface_data: ndarray, origin_x: float = 0.0, origin_y: float = 0.0, pixel_width: float = 1.0, pixel_height: float = -1.0, num_points: int = 200)[source]#

Bases: object

Class for handling interpolation of data to/from profiles.

interpolate_3d_data(data: ndarray) → ndarray[source]#

Interpolate 3D data (n_layers, ny, nx) to profile.

Parameters:: data – 3D array of values
Returns:: Array of interpolated values (n_layers, n_profile_points)

interpolate_layer_data(layer_data: List[ndarray]) → ndarray[source]#

Interpolate multiple layer data to profile.

Parameters:: layer_data – List of 2D arrays for each layer
Returns:: Array of interpolated values (n_layers, n_profile_points)

interpolate_to_mesh(property_values: ndarray, depth_values: ndarray, mesh_x: ndarray, mesh_y: ndarray, mesh_markers: ndarray, ID: ndarray, layer_markers: list = [3, 0, 2]) → ndarray[source]#

Interpolate property values from profile to mesh with layer-specific handling.

Parameters:

property_values – Property values array (n_points or n_layers, n_points)
depth_values – Depth values array (n_layers, n_points)
mesh_x – Coordinates of mesh cells
mesh_y – Coordinates of mesh cells
mesh_markers – Markers indicating different layers in mesh
layer_markers – List of marker values for each layer

Returns:

Interpolated values for mesh cells

class PyHydroGeophysX.RunContext(project_root: ~pathlib.Path = <factory>, output_dir: ~pathlib.Path = <factory>, progress: ~typing.Callable[[str], None] = <function _noop_progress>, cancelled: ~typing.Callable[[], bool] = <function _never_cancelled>, object_cache: ~typing.MutableMapping[str, ~typing.Any] = <factory>)[source]#

Bases: object

Process-local execution settings that are deliberately absent from recipes.

static cache_key(ref: ArtifactRef) → str[source]#

cancelled() → bool#

load_object(ref: ArtifactRef, loader: Callable[[Path], Any]) → Any[source]#: Load and cache an artifact by artifact_id + checksum.

object_cache: MutableMapping[str, Any]#

output_dir: Path#

prepare() → None[source]#

progress() → None#

project_root: Path#

resolve_artifact(ref: ArtifactRef) → Path[source]#

class PyHydroGeophysX.SeismicAgent(api_key: str | None = None, model: str | None = None, llm_provider: str = 'openai')[source]#

Bases: BaseAgent

Agent specialized in seismic refraction tomography (SRT) processing.

Uses PyGIMLI and PyHydroGeophysX seismic processing modules to invert seismic travel time data and extract velocity interfaces for structural constraints.

Supports two modes: - ‘inversion’: Load seismic data file and run SRT inversion - ‘interface’: Extract velocity interfaces from existing velocity model

Example

>>> agent = SeismicAgent()
>>> result = agent.execute({
...     'seismic_file': 'seismic_data.dat',
...     'velocity_threshold': 1200,
...     'output_dir': 'results/seismic'
... })

execute(input_data: Dict[str, Any]) → Dict[str, Any][source]#

Process seismic data and extract velocity structure.

Parameters:: input_data – Dictionary containing seismic_file or pre-loaded seismic_data, optional velocity thresholds, inversion parameters, an output directory, and an extract_interfaces flag. Supported inversion parameters include lam, zWeight, vTop, vBottom, paraDepth, paraMaxCellSize, and limits.
Returns:: Dictionary containing velocity model, mesh, interfaces, and visualizations

class PyHydroGeophysX.VRHModel[source]#

Bases: BaseVelocityModel

Voigt-Reuss-Hill (VRH) mixing model for effective elastic properties of composites.

calculate_properties(fractions: List[float], bulk_moduli: List[float], shear_moduli: List[float], densities: List[float]) → Tuple[float, float, float][source]#

Calculate effective elastic properties using the VRH model.

Parameters:

fractions – Volume fractions of each mineral (must sum to 1)
bulk_moduli – Bulk moduli of each mineral (GPa)
shear_moduli – Shear moduli of each mineral (GPa)
densities – Densities of each mineral (kg/m³)

Returns:

Effective bulk modulus (GPa), effective shear modulus (GPa), and effective density (kg/m³)

calculate_velocity(fractions: List[float], bulk_moduli: List[float], shear_moduli: List[float], densities: List[float]) → Tuple[float, float][source]#

Calculate P-wave and S-wave velocities using the VRH model.

Parameters:

fractions – Volume fractions of each mineral (must sum to 1)
bulk_moduli – Bulk moduli of each mineral (GPa)
shear_moduli – Shear moduli of each mineral (GPa)
densities – Densities of each mineral (kg/m³)

Returns:

P-wave velocity (m/s) and S-wave velocity (m/s)

PyHydroGeophysX.VRH_model(f: Any = [0.35, 0.25, 0.2, 0.125, 0.075], K: Any = [55.4, 36.6, 75.6, 46.7, 50.4], G: Any = [28.1, 45, 25.6, 23.65, 27.4], rho: Any = [2560, 2650, 2630, 2540, 3050]) → Any[source]#

Implements the Voigt-Reuss-Hill (VRH) mixing model to estimate the effective bulk modulus (Km), shear modulus (Gm), and density (rho_b) of a composite material made from various minerals.

Parameters: f (list): Fraction of each mineral in the composite (must sum to 1). K (list): Bulk modulus of each mineral (GPa). G (list): Shear modulus of each mineral (GPa). rho (list): Density of each mineral (kg/m^3).

Returns: Km (float): Effective bulk modulus of the composite material (GPa). Gm (float): Effective shear modulus of the composite material (GPa). rho_b (float): Effective density of the composite material (kg/m^3).

class PyHydroGeophysX.WaterContentAgent(api_key: str | None = None, model: str | None = None, llm_provider: str = 'openai')[source]#

Bases: BaseAgent

Agent specialized in converting resistivity to water content.

Uses PyHydroGeophysX petrophysical models and Monte Carlo uncertainty quantification to estimate water content from resistivity.

execute(input_data: Dict[str, Any]) → Dict[str, Any][source]#

Convert resistivity to water content.

Parameters:: input_data – Dictionary containing: - inversion_results: ERT inversion results - petrophysical_params: Parameters for each layer (rhos, n, porosity, etc.) - uncertainty_analysis: Whether to run Monte Carlo (default: False) - n_realizations: Number of MC realizations (default: 100) - output_dir: Directory for saving results
Returns:: Dictionary containing water content estimates and uncertainties

class PyHydroGeophysX.WorkflowRunResult(status: str, summary: Mapping[str, ~typing.Any]=<factory>, metrics: Mapping[str, ~typing.Any]=<factory>, artifacts: List[ArtifactRef] = <factory>, warnings: List[str] = <factory>, provenance: Mapping[str, ~typing.Any]=<factory>, objects: MutableMapping[str, ~typing.Any]=<factory>)[source]#

Bases: object

Public serializable result plus an optional same-process object channel.

artifacts: List[ArtifactRef]#

classmethod from_dict(value: Mapping[str, Any]) → WorkflowRunResult[source]#

legacy_payload() → Dict[str, Any][source]#

Return a same-process mapping for existing viewers during migration.

Unlike to_dict(), this deliberately merges the live objects channel back into the summary. It must never be persisted or emitted by CLI/bridge code.

metrics: Mapping[str, Any]#

objects: MutableMapping[str, Any]#

provenance: Mapping[str, Any]#

status: str#

summary: Mapping[str, Any]#

to_dict() → Dict[str, Any][source]#: Return the process-safe view. objects is intentionally excluded.

warnings: List[str]#

class PyHydroGeophysX.WorkflowSpec(workflow_id: str, inputs: Mapping[str, ~typing.Any]=<factory>, parameters: Mapping[str, ~typing.Any]=<factory>, seed: int | None = None, dependencies: List[str] = <factory>, metadata: Mapping[str, ~typing.Any]=<factory>, schema_version: str = '1')[source]#

Bases: object

Everything required to reproduce one workflow in another process.

dependencies: List[str]#

classmethod from_dict(value: Mapping[str, Any]) → WorkflowSpec[source]#

inputs: Mapping[str, Any]#

metadata: Mapping[str, Any]#

parameters: Mapping[str, Any]#

schema_version: str = '1'#

seed: int | None = None#

to_dict() → Dict[str, Any][source]#

validate(*, stochastic: bool = False) → None[source]#

workflow_id: str#

PyHydroGeophysX.apply_agc(data: ndarray, dt: float, window: float = 0.05, rms: float = 1.0) → ndarray[source]#

Apply gate-based automatic gain control to traces.

The implementation follows the MATLAB example’s gate/interpolation logic while handling zero-energy windows safely.

PyHydroGeophysX.calibrate_to_reference(path: str, method: str, geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], ref_resistivity: float, *, max_probe: int = 6, log: ~typing.Callable[[str], None] = <function noop>) → float[source]#

Find the data_scale that makes the recovered near-surface resistivity match a known/expected value.

The amplitude scale and the absolute resistivity level are degenerate — the EM data alone cannot fix the level (any data_scale fits the data, with the resistivity shifting to compensate). This breaks the degeneracy with EXTERNAL information: the user supplies ref_resistivity (a known background, e.g. from a borehole or regional geology), and this inverts a few probe soundings at two trial scales, fits the near-surface resistivity’s log-linear response to the scale, and solves for the scale that yields ref_resistivity. Returns the current data_scale unchanged if calibration is not possible.

PyHydroGeophysX.create_difference_gif(mesh: Any, models: Sequence[ndarray], reference: ndarray, filename: str, *, mode: str = 'difference', titles: Sequence[str] | None = None, cmap: str = 'RdBu_r', symmetric: bool = True, label: str = '', figsize: Tuple[float, float] = (8, 2.5), dpi: int = 150, duration: int = 100, coverage: ndarray | Sequence[ndarray] | None = None) → str[source]#

Create a GIF animation showing model changes relative to a reference.

Parameters:

mesh (pygimli.Mesh)
models (sequence of array-like) – Model snapshots for each frame.
reference (array-like) – Baseline model to subtract / divide.
mode ('difference' | 'ratio' | 'percent_change')
symmetric (bool) – Center colorbar at zero / one.

:param (other parameters same as create_timelapse_gif()):

Returns:: Path to the saved GIF file.
Return type:: str

PyHydroGeophysX.create_mesh_from_layers(surface: ndarray, line1: ndarray, line2: ndarray, bottom_depth: float = 30.0, quality: float = 28, area: float = 40) → Tuple[pygimli.Mesh, ndarray, ndarray][source]#

Create mesh from layer boundaries and get cell centers and markers.

Parameters:

surface – Surface coordinates [[x,z],…]
line1 – First layer boundary coordinates
line2 – Second layer boundary coordinates
bottom_depth – Depth below surface minimum for mesh bottom
quality – Mesh quality parameter
area – Maximum cell area

Returns:

PyGIMLI mesh mesh_centers: Array of cell center coordinates markers: Array of cell markers

Return type:

mesh

PyHydroGeophysX.create_timelapse_gif(mesh: Any, models: Sequence[ndarray], filename: str, *, titles: Sequence[str] | None = None, cmap: Any = None, cmin: float | None = None, cmax: float | None = None, log_scale: bool = False, label: str = '', xlabel: str = 'Distance (m)', ylabel: str = 'Elevation (m)', coverage: ndarray | Sequence[ndarray] | None = None, figsize: Tuple[float, float] = (8, 2.5), dpi: int = 150, duration: int = 100, first_frame_duration: int = 500, loop: int = 0) → str[source]#

Create a GIF animation of time-lapse model snapshots.

Each frame renders the model on the given mesh using pg.show and captures it as a PIL image. Requires Pillow.

Parameters:

mesh (pygimli.Mesh) – Mesh shared by all models.
models (sequence of array-like) – Model values for each frame.
filename (str) – Output GIF path.
titles (sequence of str, optional) – Title per frame.
cmap (str or Colormap, optional) – Colormap. Defaults to BlueDarkRed18_18_r if available.
cmin (float, optional) – Fixed color limits.
cmax (float, optional) – Fixed color limits.
log_scale (bool) – Logarithmic color scale.
label (str) – Colorbar label.
coverage (array or list of arrays, optional) – Coverage mask(s).
figsize (tuple) – Figure size per frame.
dpi (int) – Resolution.
duration (int) – Milliseconds per frame.
first_frame_duration (int) – Duration of the first frame in ms (longer for visual pause).
loop (int) – Number of loops (0 = infinite).

Returns:

Path to the saved GIF file.

Return type:

str

PyHydroGeophysX.create_timelapse_mp4(mesh: Any, models: Sequence[ndarray], filename: str, *, titles: Sequence[str] | None = None, cmap: Any = None, cmin: float | None = None, cmax: float | None = None, log_scale: bool = False, label: str = '', xlabel: str = 'Distance (m)', ylabel: str = 'Elevation (m)', coverage: ndarray | Sequence[ndarray] | None = None, figsize: Tuple[float, float] = (8, 2.5), dpi: int = 150, fps: int = 10) → str[source]#

Create an MP4 video of time-lapse model snapshots using matplotlib.

Requires ffmpeg to be available on the system path.

Parameters:

mesh (pygimli.Mesh)
models (sequence of array-like)
filename (str) – Output .mp4 path.
fps (int) – Frames per second.

:param (other parameters same as create_timelapse_gif()):

Returns:: Path to the saved MP4 file.
Return type:: str

PyHydroGeophysX.ertforandjac(fob: Any, rhomodel: Any, xr: Any) → Any[source]#

Forward model and Jacobian for ERT.

Parameters:

fob (pygimli.ERTModelling) – ERT forward operator.
rhomodel (pg.RVector) – Resistivity model.
xr (np.ndarray) – Log-transformed model parameter.

Returns:

Log-transformed forward response. J (np.ndarray): Jacobian matrix.

Return type:

dr (np.ndarray)

PyHydroGeophysX.ertforandjac2(fob: Any, xr: Any, mesh: Any) → Any[source]#

Alternative ERT forward model and Jacobian using log-resistivity values.

Parameters:

fob (pygimli.ERTModelling) – ERT forward operator.
xr (np.ndarray) – Log-transformed model parameter.
mesh (pg.Mesh) – Mesh for the forward model.

Returns:

Log-transformed forward response. J (np.ndarray): Jacobian matrix.

Return type:

dr (np.ndarray)

PyHydroGeophysX.ertforward(fob: Any, mesh: Any, rhomodel: Any, xr: Any) → Any[source]#

Forward model for ERT.

Parameters:

fob (pygimli.ERTModelling) – ERT forward operator.
mesh (pg.Mesh) – Mesh for the forward model.
rhomodel (pg.RVector) – Resistivity model vector.
xr (np.ndarray) – Log-transformed model parameter (resistivity).

Returns:

Log-transformed forward response. rhomodel (pg.RVector): Updated resistivity model.

Return type:

dr (np.ndarray)

PyHydroGeophysX.ertforward2(fob: Any, xr: Any, mesh: Any) → Any[source]#

Simplified ERT forward model.

Parameters:

fob (pygimli.ERTModelling) – ERT forward operator.
xr (np.ndarray) – Log-transformed model parameter.
mesh (pg.Mesh) – Mesh for the forward model.

Returns:

Log-transformed forward response.

Return type:

dr (np.ndarray)

PyHydroGeophysX.estimate_data_scale(path: str, method: str, geom: ~typing.Dict[str, ~typing.Any], *, max_soundings: int = 8, log: ~typing.Callable[[str], None] = <function noop>) → float[source]#

Estimate the amplitude calibration (data_scale) for normalized data.

Normalized airborne responses (e.g. moment-normalized dB/dt) differ from the studio’s 1D forward by a near-constant amplitude factor. This fits each sounding’s decay SHAPE to a grid of half-space forward responses at the current geometry and takes the geometric-mean amplitude ratio forward/observed at the best-fitting resistivity. Returns 1.0 if it cannot be estimated (so the caller can fall back to no scaling).

PyHydroGeophysX.export_first_breaks(picks: Iterable[FirstBreakPick | Dict[str, Any]], filename: str) → str[source]#: Export first-break picks to CSV.

PyHydroGeophysX.export_mesh_to_vtk(mesh: Any, filename: str, cell_data: Dict[str, ndarray] | None = None) → str[source]#

Export a PyGIMLi mesh to VTK unstructured grid format.

This uses PyGIMLi’s built-in mesh.exportVTK when available and then optionally injects extra cell data fields.

Parameters:

mesh (pygimli.Mesh) – The mesh to export.
filename (str) – Output .vtk path.
cell_data (dict of str -> array, optional) – Additional named cell-data arrays to write.

Returns:

Path to the written file.

Return type:

str

PyHydroGeophysX.export_points_to_vtk(points: ndarray, scalars: Dict[str, ndarray], filename: str) → str[source]#

Export point-cloud data (e.g. cell centers) to VTK PolyData.

Parameters:

points (array of shape (N, 3)) – XYZ coordinates.
scalars (dict of str -> array) – Named scalar values at each point.
filename (str) – Output .vtk path.

Returns:

Path to the written file.

Return type:

str

PyHydroGeophysX.export_results_to_csv(result: Any, filename: str) → str[source]#: Export inversion result arrays to a tabular CSV file.

PyHydroGeophysX.export_structured_vtk(values: ndarray, filename: str, *, dx: float = 1.0, dy: float = 1.0, dz: float = 1.0, origin: tuple = (0.0, 0.0, 0.0), scalar_name: str = 'model') → str[source]#

Export a 3-D numpy array to VTK structured-points format.

Useful for MODFLOW/ParFlow grids that map directly to a regular grid.

Parameters:

values (3-D array) – Shape (nz, ny, nx) — layer, row, column ordering.
filename (str) – Output .vtk path.
dx (float) – Cell spacing in each direction.
dy (float) – Cell spacing in each direction.
dz (float) – Cell spacing in each direction.
origin (tuple of float) – (x0, y0, z0) origin of the grid.
scalar_name (str) – Name for the scalar dataset.

Returns:

Path to the written file.

Return type:

str

PyHydroGeophysX.export_structured_vtk_multi(scalars: Dict[str, ndarray], filename: str, *, dx: float = 1.0, dy: float = 1.0, dz: float = 1.0, origin: tuple = (0.0, 0.0, 0.0)) → str[source]#

Export multiple 3-D arrays as named scalars in a single VTK file.

All arrays must have the same shape (nz, ny, nx).

Parameters:

scalars (dict of str -> 3-D array) – Named scalar datasets (e.g. {"resistivity": res, "water_content": wc}).
filename (str) – Output .vtk path.

Returns:

Path to the written file.

Return type:

str

PyHydroGeophysX.export_timelapse_structured_vtk(models: Sequence[ndarray], output_dir: str, prefix: str = 'timelapse', *, dx: float = 1.0, dy: float = 1.0, dz: float = 1.0, origin: tuple = (0.0, 0.0, 0.0), scalar_name: str = 'model') → List[str][source]#

Export time-lapse 3-D arrays as a numbered VTK series.

Parameters:

models (sequence of 3-D arrays) – Each entry has shape (nz, ny, nx).
output_dir (str) – Output directory.
prefix (str) – File prefix.
dx (float) – Grid spacing.
dy (float) – Grid spacing.
dz (float) – Grid spacing.
origin (tuple) – Grid origin.
scalar_name (str) – Scalar field name.

Returns:

Paths to generated VTK files.

Return type:

list of str

PyHydroGeophysX.export_timelapse_vtk(mesh: Any, models: Sequence[ndarray], output_dir: str, prefix: str = 'timelapse', *, scalar_name: str = 'model', extra_data: Dict[str, Sequence[ndarray]] | None = None) → List[str][source]#

Export time-lapse models as a numbered VTK series for ParaView.

Creates files prefix_0000.vtk, prefix_0001.vtk, … and a .pvd collection file that ParaView can load as a time series.

Parameters:

mesh (pygimli.Mesh) – Mesh shared by all timesteps.
models (sequence of array-like) – Model values per timestep.
output_dir (str) – Directory for output files.
prefix (str) – Filename prefix.
scalar_name (str) – Name for the primary scalar field.
extra_data (dict of str -> sequence of arrays, optional) – Additional scalar fields per timestep.

Returns:

Paths to all generated VTK files.

Return type:

list of str

PyHydroGeophysX.export_to_vtk(result: Any, mesh: Any, filename: str) → str[source]#

Export result/model values to a lightweight VTK PolyData file.

The export writes cell-center points with a single scalar field.

PyHydroGeophysX.extract_profile(grid: Dict[str, ndarray], p1: Sequence[float], p2: Sequence[float], n: int = 200) → Dict[str, ndarray][source]#: Sample a gridded field along the line p1 -> p2 (bilinear).

PyHydroGeophysX.fdem_invert(data: ~typing.Dict[str, ~typing.Any], geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#: Invert an FDEM sounding for a layered resistivity model (Occam 1D).

PyHydroGeophysX.first_breaks_to_traveltime(picks: Iterable[FirstBreakPick | Dict[str, Any]], filename: str, receiver_spacing: float = 1.0, shot_spacing: float | None = None) → str[source]#: Export first breaks to a PyGIMLi/BERT travel-time .dat file.

PyHydroGeophysX.forward_bodies(xobs: ~numpy.ndarray, yobs: ~numpy.ndarray, kind: str, bodies: ~typing.List[~typing.Dict[str, ~typing.Any]], field: ~typing.Dict[str, ~typing.Any] | None = None, log: ~typing.Callable[[str], None] = <function noop>) → ndarray[source]#: Sum the anomaly of a list of bodies. kind = ‘gravity’ or ‘magnetics’.

PyHydroGeophysX.generate_python(spec: WorkflowSpec, path: str | Path) → Path[source]#: Generate a standalone, editable script using only the public workflow API.

PyHydroGeophysX.get_joint_capabilities(include_planned: bool = True) → List[JointPairCapability][source]#: List implemented capabilities and, optionally, planned method pairs.

PyHydroGeophysX.get_workflow(workflow_id: str) → WorkflowDescriptor[source]#

PyHydroGeophysX.gravity_prism(xobs: ndarray, yobs: ndarray, body: Dict[str, Any]) → ndarray[source]#: Vertical gravity (mGal) of a right rectangular prism (Nagy 1966). z down.

PyHydroGeophysX.gravity_sphere(xobs: ndarray, yobs: ndarray, body: Dict[str, Any]) → ndarray[source]#: Vertical gravity (mGal) of a buried sphere. z positive down, obs at z=0.

PyHydroGeophysX.grid_data(x: ndarray, y: ndarray, value: ndarray, nx: int = 120, ny: int = 120, method: str = 'linear') → Dict[str, ndarray][source]#: Grid scattered station values onto a regular map. Returns xx, yy, zz.

PyHydroGeophysX.interpolate_structure_to_profile(structure_data: List[ndarray], X_grid: ndarray, Y_grid: ndarray, X_pro: ndarray, Y_pro: ndarray) → ndarray[source]#

Interpolate multiple structure layers onto profile

Parameters:

structure_data – List of 2D arrays for each layer
X_grid – X coordinates of original grid
Y_grid – Y coordinates of original grid
X_pro – X coordinates of profile points
Y_pro – Y coordinates of profile points

Returns:

Array of interpolated values with shape (n_layers, n_points)

PyHydroGeophysX.interpolate_to_mesh(property_values: ndarray, profile_distance: ndarray, depth_values: ndarray, mesh_x: ndarray, mesh_y: ndarray, mesh_markers: ndarray, ID: Any, layer_markers: list = [3, 0, 2]) → ndarray[source]#

Interpolate property values from profile to mesh with layer-specific handling.

Parameters:

property_values – Property values array (n_points)
profile_distance – Distance along profile (n_points)
depth_values – Depth values array (n_layers, n_points)
mesh_x – X coordinates of mesh cells
mesh_y – Y coordinates of mesh cells
mesh_markers – Markers indicating different layers in mesh
layer_markers – List of marker values for each layer

Returns:

Interpolated values for mesh cells

PyHydroGeophysX.interpolate_to_profile(data: ndarray, X_grid: ndarray, Y_grid: ndarray, X_pro: ndarray, Y_pro: ndarray, method: str = 'linear') → ndarray[source]#

Interpolate 2D data onto a profile line

Parameters:

data – 2D array of values to interpolate
X_grid – X coordinates of original grid (meshgrid)
Y_grid – Y coordinates of original grid (meshgrid)
X_pro – X coordinates of profile points
Y_pro – Y coordinates of profile points
method – Interpolation method (‘linear’ or ‘nearest’)

Returns:

Interpolated values along profile

PyHydroGeophysX.invert_gravmag(x, y, value, kind: str, *, z: ~numpy.ndarray | None = None, field: ~typing.Dict[str, ~typing.Any] | None = None, detrend: int = 0, n_xy: int = 22, n_z: int = 12, max_iterations: int = 20, beta0_ratio: float = 1.0, max_stations: int = 600, relative_error: float = 0.03, noise_floor: float | None = None, solver: str = 'simpeg', auto_beta: bool = True, target_chi2: float = 1.0, chi2_tolerance: float = 0.2, max_beta_trials: int = 6, sensitivity_power: float = 1.0, out_dir: str | None = None, random_seed: int | None = 42, log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#

Run a SimPEG 3D potential-field inversion under the survey.

gravity recovers a density-contrast model (g/cc); magnetics recovers a susceptibility model (SI) and needs field = {inclination, declination, strength_nT}. z is optional per-station elevation (m, positive upward); a missing value falls back to 1 m. detrend (0..3) removes a polynomial regional trend before inversion. The returned grid uses elevation increasing upward. random_seed makes SimPEG’s eigenvalue-based beta estimate reproducible. Raises InversionBackendUnavailable if SimPEG is missing.

PyHydroGeophysX.invert_line(path: str, method: str, geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], *, spacing: float = 50.0, positions: ~numpy.ndarray | None = None, heights: ~numpy.ndarray | None = None, max_soundings: int = 12, lines: ~typing.Sequence[int] | None = None, doi_blank: bool = True, doi_factor: float = 0.5, ref_resistivity: float = 0.0, out_dir: ~pathlib.Path | None = None, initial_models: ~numpy.ndarray | None = None, log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#

Invert a line on a shared fixed-layer grid.

inv["lci_mode"] selects how the soundings are coupled:

simultaneous (the default whenever lateral_smoothness is positive): Solves the whole line as one system, with the lateral constraint part of what is being minimized. See PyHydroGeophysX.inversion.em1d_lci.
sequential: The older block-coordinate passes: each station is re-inverted on its own against the distance-weighted model its neighbours had at the end of the previous pass. Kept because it needs no analytic Jacobian, so it still runs against a forward operator that cannot supply one.
off: Independent 1D inversion per sounding, no lateral coupling.

The models are laid side by side to form a resistivity(position, depth) section ready for Model3DView.show_model(): edges = (ex, ey, ez) (ez is elevation, increasing upward) and model3d of shape (n_pos, 1, n_depth).

positions gives the along-line distance of each sounding (the section x-axis); heights overrides the sensor height per sounding. When doi_blank is set, cells below a per-sounding depth of investigation (a diffusion-depth estimate scaled by doi_factor) are blanked (NaN) so the unconstrained deep part of an early-time sounding is not shown as railed.

inv["robust_errors"] retains all imported gates and iteratively inflates effective errors for large residuals. It overrides hard rejection. The main chi2 uses ORIGINAL errors; result["robust"] records effective errors and a separate effective chi2. Import-time flags and QC still apply.

inv["auto_lambda"] re-solves the line at other smoothness weights to reach target_chi2. inv["reject_outliers"] drops the gates the converged model cannot explain (beyond outlier_threshold sigma, over outlier_passes cycles, never below min_data_fraction of the gates) and solves again; what it removed is reported under result["outliers"]. They address different causes, so they can be used together: relaxing the smoothness helps when the model is too stiff for the data, rejection helps when a minority of gates are simply wrong.

lines restricts the run to the named survey lines, so a line whose data is thinner than the rest can be given its own settings instead of one set having to suit every line. Passing None runs from the first station, as before. The lateral constraint already groups by line, so a line inverted on its own is tied exactly as it would be inside a whole-survey run; what changes is which settings reach it, and that the other lines are not re-solved. max_soundings then counts within the selection.

Stations arrive ordered by line, so a selection is a contiguous block. A set of lines that is not contiguous is refused rather than quietly widened to the span that encloses it, which would invert the lines in between under settings chosen for their neighbours.

PyHydroGeophysX.list_workflows() → Tuple[WorkflowDescriptor, ...][source]#

PyHydroGeophysX.load_line_geometry(path: str) → Dict[str, Any][source]#: Load per-sounding line geometry: along-line positions (m), optional sensor heights (m), and optional map coordinates x/y (e.g. easting/northing) for plan-view depth slices. Recognizes header names (distance/position for the position; alt/height for the height; easting/northing for the map coordinates, which also derive the distance when no distance column is present). A header-less file is read by column order (1 column = position; 2+ = position, height). positions is shifted to start at 0.

PyHydroGeophysX.load_recipe(path: str | Path) → WorkflowSpec[source]#

PyHydroGeophysX.load_sounding(path: str, method: str, sounding: int = 0, *, moment: str = 'HM', use_flags: bool = True, max_relative_std: float | None = None, gate_rejection: str = 'truncate', reject_negative: bool = False, min_gates_per_moment: Mapping[str, int] | None = None, ttem_loop_area: float | None = None, ttem_gex_path: str | None = None, ttem_tfi_path: str | None = None) → Dict[str, Any][source]#

Load one sounding from a sounding file.

The first column is the abscissa (FDEM: frequency Hz; TDEM: time s). The remaining columns hold the response(s) — a single sounding, or several stacked side by side so one file can carry a whole survey line (common for airborne EM exports). sounding picks which one (0-based):

TDEM: each extra column is one sounding’s response → column 1 + sounding.
FDEM: response columns come in (real, imag) pairs → one sounding is the pair starting at 1 + 2*sounding; a lone trailing real column gives imag = 0.

The returned dict also reports n_soundings so the caller can offer a picker. use_flags applies only to TEMcompany project databases; see load_temcompany_sounding().

PyHydroGeophysX.magnetic_dipole(xobs: ndarray, yobs: ndarray, body: Dict[str, Any], field: Dict[str, Any]) → ndarray[source]#: Total-field magnetic anomaly (nT) of an induced/magnetized sphere (a dipole).

PyHydroGeophysX.normalize_traces(data: ndarray, trace_axis: int = 1, eps: float = 1e-12) → ndarray[source]#: Normalize each seismic trace by its maximum absolute amplitude.

PyHydroGeophysX.pick_first_breaks(data: SeismicDataset | SeismicShotGather | ndarray, dt: float | None = None, headers: Sequence[SeismicTraceHeader] | None = None, threshold: float = 0.2, noise_multiplier: float = 5.0, min_time: float = 0.0, max_time: float | None = None, polarity: float = 1.0) → List[FirstBreakPick][source]#

Pick first arrivals using a simple amplitude/noise threshold.

This assisted picker is intended as a starting point for GUI review rather than a final scientific picking algorithm.

PyHydroGeophysX.plot_convergence(chi2_history: Sequence[float], *, ax: Any = None, target_chi2: float = 1.0, ylabel: str = '$\\chi^2$', title: str = 'Inversion Convergence') → Tuple[source]#

Plot chi-squared convergence curve.

Parameters:

chi2_history (sequence of float) – Chi-squared value per iteration.
target_chi2 (float) – Target misfit (plotted as a dashed line).

Return type:

fig, ax

PyHydroGeophysX.plot_coverage(mesh: Any, coverage: ndarray, *, ax: Any = None, cmap: str = 'YlGn', threshold: float | None = None, title: str = 'Data Coverage') → Tuple[source]#

Plot a coverage / sensitivity map.

Parameters:

mesh (pygimli.Mesh)
coverage (array-like) – Coverage values per cell.
threshold (float, optional) – If given, overlay a contour at this level.

Return type:

fig, ax, cbar

PyHydroGeophysX.plot_difference_map(mesh: Any, model_a: ndarray, model_b: ndarray, *, mode: str = 'difference', ax: Any = None, cmap: str = 'RdBu_r', symmetric: bool = True, label: str = '', title: str = '', coverage: ndarray | None = None) → Tuple[source]#

Plot the difference or ratio between two models.

Parameters:

mesh (pygimli.Mesh)
model_a (array-like) – Two model arrays. result = model_b - model_a (difference) or model_b / model_a (ratio).
model_b (array-like) – Two model arrays. result = model_b - model_a (difference) or model_b / model_a (ratio).
mode ('difference' | 'ratio' | 'percent_change')
symmetric (bool) – If True, center the colorbar at zero (difference) or one (ratio).
coverage (array-like, optional) – Coverage mask.

Return type:

fig, ax, cbar

PyHydroGeophysX.plot_electrode_layout(positions: Dict[str, ndarray], *, ax: Any = None, color_by: str = 'z', cmap: str = 'terrain', title: str = 'Electrode Layout') → Tuple[source]#

Scatter-plot electrode positions colored by elevation.

Parameters:: positions (dict) – Must contain 'x' and 'y' keys; optionally 'z'.
Return type:: fig, ax

PyHydroGeophysX.plot_model_section(mesh: Any, values: ndarray, *, ax: Any = None, cmap: Any = None, cmin: float | None = None, cmax: float | None = None, log_scale: bool = False, label: str = '', xlabel: str = 'Distance (m)', ylabel: str = 'Elevation (m)', title: str = '', coverage: ndarray | None = None, orientation: str = 'vertical') → Tuple[source]#

Plot a 2D model cross-section on a PyGIMLi mesh.

Parameters:

mesh (pygimli.Mesh) – The mesh to plot on.
values (array-like) – Cell values (resistivity, velocity, water content, etc.).
ax (matplotlib.axes.Axes, optional) – Axes to draw on. Created if None.
cmap (str or Colormap, optional) – Colormap. Defaults to BlueDarkRed18_18_r if available.
cmin (float, optional) – Color limits.
cmax (float, optional) – Color limits.
log_scale (bool) – Use logarithmic color scaling.
label (str) – Colorbar label.
coverage (array-like, optional) – Coverage array for masking low-sensitivity cells.
orientation (str) – Colorbar orientation ('vertical' or 'horizontal').

Return type:

fig, ax, cbar

PyHydroGeophysX.plot_monitoring_timeseries(times: ndarray, series: Dict[str, ndarray], *, true_series: Dict[str, ndarray] | None = None, uncertainties: Dict[str, Tuple[ndarray, ndarray]] | None = None, ax: Any = None, ylabel: str = 'Value', title: str = 'Monitoring Point Time Series') → Tuple[source]#

Plot estimated (and optionally true) time-series at monitoring points.

Parameters:

times (array-like) – Time axis.
series (dict of str -> array) – Estimated values keyed by point name.
true_series (dict of str -> array, optional) – True / reference values for comparison (dashed lines).
uncertainties (dict of str -> (lower, upper), optional) – Uncertainty bounds per point for shading.

Return type:

fig, ax

PyHydroGeophysX.plot_pseudosection_matrix(data_matrix: ndarray, *, ax: Any = None, cmap: Any = None, vmin: float | None = None, vmax: float | None = None, xlabel: str = 'Time', ylabel: str = 'Measurement #', label: str = 'Apparent resistivity ($\\Omega\\cdot$m)', title: str = '') → Tuple[source]#

Plot a time-lapse apparent resistivity matrix as a heatmap.

Parameters:: data_matrix (2-D array) – Shape (n_times, n_measurements) or similar.
Return type:: fig, ax, im

PyHydroGeophysX.plot_timelapse_snapshots(mesh: Any, models: Sequence[ndarray], *, titles: Sequence[str] | None = None, ncols: int = 4, cmap: Any = None, cmin: float | None = None, cmax: float | None = None, log_scale: bool = False, label: str = '', coverage: ndarray | Sequence[ndarray] | None = None, figsize_per_panel: Tuple[float, float] = (4.0, 2.5)) → Tuple[source]#

Plot a grid of time-lapse model snapshots.

Parameters:

mesh (pygimli.Mesh) – Mesh shared by all snapshots.
models (sequence of array-like) – Model arrays for each timestep.
titles (sequence of str, optional) – Panel titles. Defaults to 'Timestep 1', 'Timestep 2', …
ncols (int) – Number of columns.
cmap – Passed to pg.show.
cmin – Passed to pg.show.
cmax – Passed to pg.show.
log_scale – Passed to pg.show.
label – Passed to pg.show.
coverage (array or sequence of arrays, optional) – Coverage mask(s). If a single 1-D array it is reused for all panels. If 2-D, coverage[i] is used for panel i.
figsize_per_panel (tuple) – (width, height) per subplot panel.

Return type:

fig, axes

PyHydroGeophysX.plot_topography(topo_grid: ndarray, *, profile_endpoints: List[Tuple[float, float]] | None = None, ax: Any = None, cmap: str = 'terrain', title: str = 'Surface Topography') → Tuple[source]#

Plot a 2-D topography grid with optional profile line overlay.

Parameters:

topo_grid (2-D array) – Elevation raster.
profile_endpoints (list of (row, col) tuples, optional) – If two points are given, draw the profile line.

Return type:

fig, ax

PyHydroGeophysX.prepare_2D_profile_data(data: ndarray, XX: ndarray, YY: ndarray, X_pro: ndarray, Y_pro: ndarray) → ndarray[source]#

Interpolate multiple 2D gridded data layers onto a profile line.

Parameters:

data – 3D array of gridded data (n_layers, ny, nx)
XX – Coordinate grids from meshgrid
YY – Coordinate grids from meshgrid
X_pro – Profile line coordinates
Y_pro – Profile line coordinates

Returns:

Interpolated values along profile (n_layers, n_profile_points)

PyHydroGeophysX.qc_products(x: ndarray, y: ndarray, value: ndarray, *, detrend: int = 1, nx: int = 120, ny: int = 120) → Dict[str, Any][source]#: Calculate observed, regional and residual products for map/profile QC.

PyHydroGeophysX.read_seg2_seismic(file: str) → Dict[str, ndarray][source]#

Read SEG-2 seismic file into arrays.

If ObsPy is available, uses ObsPy SEG-2 reader. Otherwise falls back to text loading for already exported arrays.

PyHydroGeophysX.read_segy(file: str, max_traces: int | None = None, load_traces: bool = True, prefer_obspy: bool = True) → SeismicDataset[source]#

Read a SEG-Y file into traces, headers, and metadata.

Reader order is segyio -> ObsPy -> built-in: the most robust library available is used, and the built-in conservative reader is the final fallback so reads keep working with no third-party SEG-Y dependency installed.

Parameters:

file (str) – SEG-Y path.
max_traces (int, optional) – Maximum number of traces to read. Useful for responsive GUI previews.
load_traces (bool, optional) – If False, parse headers but skip sample arrays (built-in reader only).
prefer_obspy (bool, optional) – When True (default), try ObsPy if segyio is unavailable or fails.

Returns:

Parsed seismic dataset.

Return type:

SeismicDataset

PyHydroGeophysX.read_tem_fast(file: str) → Dict[str, ndarray][source]#: Read TEM-FAST style text data into SimPEG-ready arrays.

PyHydroGeophysX.regional_residual(x: ndarray, y: ndarray, value: ndarray, degree: int = 1) → Tuple[ndarray, ndarray][source]#: Fit a polynomial regional trend of degree (1..3); return (regional, residual).

PyHydroGeophysX.resistivity_to_saturation(resistivity: Any, porosity: Any, m: Any, rho_fluid: Any, n: Any, sigma_sur: Any = 0, a: Any = 1.0) → Any[source]#

Convert resistivity to saturation using Waxman-Smits model.

The function calculates saturated resistivity using Archie’s law: rhos = a * rho_fluid * porosity^(-m)

Then solves the Waxman-Smits equation: 1/rho = sigma_sat * S^n + sigma_sur * S^(n-1) where sigma_sat = 1/rhos

Parameters:

resistivity (array) – Resistivity values (ohm-m)
porosity (array) – Porosity values (fraction, 0-1)
m (float) – Cementation exponent (typically 1.3-2.5)
rho_fluid (float) – Fluid resistivity (ohm-m)
n (float) – Saturation exponent (typically 1.8-2.2)
sigma_sur (float) – Surface conductivity (S/m). Default is 0 (no surface effects)
a (float) – Tortuosity factor. Default is 1.0

Returns:

Saturation values (fraction, 0-1)

Return type:

array

PyHydroGeophysX.resistivity_to_water_content(resistivity: Any, rhos: Any, n: Any, porosity: Any, sigma_sur: Any = 0) → Any[source]#

Convert resistivity to water content using Waxman-Smits model.

Parameters:

resistivity (array) – Resistivity values
rhos (float) – Saturated resistivity without surface effects
n (float) – Saturation exponent
porosity (array) – Porosity values
sigma_sur (float) – Surface conductivity. Default is 0 (no surface effects).

Returns:

Volumetric water content values

Return type:

array

PyHydroGeophysX.run_joint_inversion(*args, **kwargs)[source]#: Lazily load and run the backend selected by a joint request.

PyHydroGeophysX.run_workflow(spec: WorkflowSpec | Mapping[str, Any], context: RunContext | None = None) → WorkflowRunResult[source]#: Validate and execute one registered workflow.

PyHydroGeophysX.satK(Keff: Any, Km: Any, phi: Any, Sat: Any) → Any[source]#

Calculate the saturated bulk modulus (K_sat) based on Brie’s equation.

Parameters: Keff (float): Effective bulk modulus of the dry rock (GPa). Km (float): Bulk modulus of the matrix (GPa). phi (float): Porosity of the rock. Sat (float): Saturation level of the fluid in the pores.

Returns: float: Saturated bulk modulus (GPa).

PyHydroGeophysX.save_recipe(spec: WorkflowSpec, path: str | Path) → Path[source]#: Validate and save a workflow recipe as UTF-8 JSON.

PyHydroGeophysX.setup_profile_coordinates(point1: List[int], point2: List[int], surface_data: ndarray, origin_x: float = 0.0, origin_y: float = 0.0, pixel_width: float = 1.0, pixel_height: float = -1.0, num_points: int = 200) → Tuple[ndarray, ndarray, ndarray, ndarray, ndarray][source]#

Set up profile coordinates based on surface elevation data between two points

Parameters:

point1 – Starting point indices [col, row]
point2 – Ending point indices [col, row]
surface_data – 2D array of surface elevation data
origin_x – X coordinate of origin
origin_y – Y coordinate of origin
pixel_width – Width of each pixel
pixel_height – Height of each pixel (negative for top-down)
num_points – Number of points along profile

Returns:

X coordinates along profile Y_pro: Y coordinates along profile L_profile: Distances along profile XX: X coordinate grid YY: Y coordinate grid

Return type:

X_pro

PyHydroGeophysX.tdem_invert(data: ~typing.Dict[str, ~typing.Any], geom: ~typing.Dict[str, ~typing.Any], inv: ~typing.Dict[str, ~typing.Any], log: ~typing.Callable[[str], None] = <function noop>) → Dict[str, Any][source]#: Invert a TDEM sounding for a layered resistivity model (Occam 1D).

PyHydroGeophysX.velDEM(phi: Any, Km: Any, Gm: Any, rho_b: Any, Sat: Any, alpha: Any) → Any[source]#

Calculate effective bulk modulus (Keff), shear modulus (Geff), and P-wave velocity (Vp) for a rock with varying porosity (phi) based on the DEM model, taking into account the saturation (Sat) and the crack aspect ratio (alpha).

Parameters: phi (np.array): Array of porosities. Km (float): Initial bulk modulus of the material (GPa). Gm (float): Initial shear modulus of the material (GPa). rho_b (float): Density of the solid phase (kg/m^3). Sat (float): Saturation level of the fluid in the cracks (0 to 1, where 1 is fully saturated). alpha (float): Crack aspect ratio.

Returns: Keff1 (np.array): Effective bulk modulus for each porosity value (GPa). Geff1 (np.array): Effective shear modulus for each porosity value (GPa). Vp (np.array): P-wave velocity for each porosity value (m/s).

PyHydroGeophysX.vel_porous(phi: Any, Km: Any, Gm: Any, rho_b: Any, Sat: Any, depth: Any = 1) → Any[source]#

Calculate P-wave velocity (Vp) for a rock with varying porosity (phi) based on the Hertz-Mindlin model and Hashin-Shtrikman bounds, taking into account the saturation (Sat).

Parameters: phi (np.array): Array of porosities. Km (float): Bulk modulus of the solid phase (GPa). Gm (float): Shear modulus of the solid phase (GPa). rho_b (float): Density of the solid phase (kg/m^3). Sat (float): Saturation level of the fluid in the pores (0 to 1, where 1 is fully saturated). depth (float): depth for pressure estimation (m)

Returns: Vp_h (np.array): P-wave velocity for each porosity value (upper Hashin-Shtrikman bound) (m/s). Vp_l (np.array): P-wave velocity for each porosity value (lower Hashin-Shtrikman bound) (m/s).

PyHydroGeophysX.water_content_to_resistivity(water_content: Any, rhos: Any, n: Any, porosity: Any, sigma_sur: Any = 0) → Any[source]#

Convert water content to resistivity using Waxman-Smits model.

Parameters:

water_content (array) – Volumetric water content (θ)
rhos (float) – Saturated resistivity without surface effects
n (float) – Saturation exponent
porosity (array) – Porosity values (φ)
sigma_sur (float) – Surface conductivity. Default is 0 (no surface effects).

Returns:

Resistivity values

Return type:

array