2.22. pycsamt.models #

Raises:

FileNotFoundError – Raised when no executable is found and compilation is disabled or does not produce a binary.
RuntimeError – Propagated from compile() when the compiler is missing or make fails.

See also

OccamRunner.compile: Compiles the bundled Fortran source.
OccamRunner.run: Calls this method before launching the solver.

Examples

>>> from pycsamt.models.occam2d import OccamRunner
>>> runner = OccamRunner("occam_run")
>>> binary = runner.discover_binary(False)

compile(fc='gfortran', flags='-O2')#

Compile the bundled Occam2D Fortran source.

Compilation is performed in the package _source directory by invoking make with FC90 and FCFLAGS variables. The resulting executable is expected to be named Occam2D there. This method does not copy the binary into workdir; discover_binary() uses that path directly.

Parameters:

fc (str, default "gfortran") – Fortran compiler command passed to make as FC90. Use this to select another compiler that understands the bundled source.
flags (str, default "-O2") – Compiler flags passed to make as FCFLAGS. Optimization flags are usually sufficient; debug builds can pass flags such as "-g".

Returns:

Path to the compiled binary inside _source.

Return type:

Raises:

FileNotFoundError – Raised when the source directory is absent.
RuntimeError – Raised when the requested compiler is unavailable, make fails, or no executable is produced.

Examples

>>> from pycsamt.models.occam2d import OccamRunner
>>> runner = OccamRunner("occam_run")
>>> binary = runner.compile("gfortran", "-O2")

run(max_iter=None, target_misfit=None, auto_compile=True, timeout=None)#

Run Occam2D synchronously.

This method blocks until the executable exits. It resolves the binary, optionally patches the startup file, launches Occam2D <startup_file> inside workdir, and writes process streams to occam_stdout.log and occam_stderr.log.

Parameters:

max_iter (int, optional) – Temporary override for the Iterations to run field in the startup file. The override is written in place before launch, so the file must be writable.
target_misfit (float, optional) – Temporary override for the Target Misfit field in the startup file. This changes the run-control file before launch.
auto_compile (bool, default True) – Passed to discover_binary(). If True, missing binaries may trigger compilation.
timeout (float, optional) – Maximum wall-clock seconds to wait before killing the process. None (the default) waits indefinitely, matching the solver’s own Iterations to run and internal step-search limits. Occam2D’s outer iteration count does not bound the inner Lagrange-multiplier line search, so a numerically pathological model can in principle run far longer than any observed typical case; pass an explicit bound for unattended batch runs. On expiry, the process (and its full process group, on POSIX) is killed and exit_code is set to -9.

Returns:

Process exit code. A value of 0 indicates that the executable returned successfully. A value of -9 indicates the process was killed after exceeding timeout.

Return type:

int

Raises:

FileNotFoundError – Raised when the binary or patched startup file cannot be found.
RuntimeError – Propagated from compilation if automatic compilation fails.

See also

OccamRunner.run_async: Starts the same executable without blocking.
OccamRunner._patch_startup: Applies max_iter and target_misfit.
InversionResult: Loads output files after a successful run.

Examples

>>> from pycsamt.models.occam2d import OccamRunner
>>> runner = OccamRunner(workdir="occam_run")
>>> code = runner.run(max_iter=100, target_misfit=1.0)

Bound unattended batch runs against a pathological case:

>>> code = runner.run(timeout=3600)

run_forward(output_root='Forward', *, auto_compile=True)#

Run the native full-2-D forward solver and return its output.

The bundled Occam2D executable supports a forward-only -F mode. It evaluates the parameter vector in startup_file and writes an OCCAM2MTDATA_1.0 file containing the modeled data and the errors from the input data file. No inversion iteration is performed.

Parameters:

output_root (str, default "Forward") – Local root used for the generated .fwd filename. It must be non-empty and contain no path separators or dots.
auto_compile (bool, default True) – If True, allow binary discovery to compile the bundled Fortran source when no executable can be found.

Returns:

Path to the generated forward-data file.

Return type:

Raises:

FileNotFoundError – Raised when the startup file or solver binary is missing.
RuntimeError – Raised when the solver exits unsuccessfully or reports success without creating the expected output.
ValueError – Raised when output_root is not a local filename root.

See also

OccamRunner.run: Executes the iterative Occam inversion.
OccamData.read: Reads the generated forward-data file.

Examples

Evaluate a known model stored in TruthStartup:

>>> runner = OccamRunner(
...     "synthetic_case",
...     startup_file="TruthStartup",
... )
>>> forward_file = runner.run_forward("TruthForward")

run_async(auto_compile=True)#

Start Occam2D in a background process.

The method resolves the binary and launches the solver with subprocess.Popen. It returns immediately with a process handle. Standard output and standard error are redirected to stdout_log and stderr_log.

Parameters:

auto_compile (bool, default True) – Passed to discover_binary(). If True, missing binaries may trigger compilation.

Returns:

Live process handle for the background run.

Return type:

subprocess.Popen

Raises:

FileNotFoundError – Raised when no executable can be found.
RuntimeError – Propagated from automatic compilation failures.

See also

OccamRunner.wait: Blocks until the background process completes.
OccamRunner.is_running: Reports whether the process is still active.

Examples

>>> from pycsamt.models.occam2d import OccamRunner
>>> runner = OccamRunner(workdir="occam_run")
>>> process = runner.run_async()
>>> runner.is_running

wait()#

Block until the asynchronous run finishes.

Returns:: Exit code returned by the background process.
Return type:: int
Raises:: RuntimeError – Raised when no process has been started with run_async().

property is_running: bool#: Whether an asynchronous Occam2D process is active.

class pycsamt.models.occam2d.InversionResult(workdir='.', iteration=None, **kwargs)#

Bases: OccamBase

Load and summarize a completed Occam2D inversion run.

InversionResult is the post-processing access layer for an Occam2D working directory. It scans the directory for input and output files, loads the selected iteration, matches the response file, reads the log, and builds a two-dimensional log10-resistivity grid on the mesh.

The reconstruction maps the iteration parameter vector to model layers and mesh cells. If \(m_j\) is the log10-resistivity value assigned to model parameter \(j\), then the physical resistivity represented by a cell in that parameter group is

\[\rho_j = 10^{m_j}.\]

The stored rho_2d array keeps \(m_j\), not \(\rho_j\), because Occam iteration files store log10-resistivity values.

Parameters:

workdir (path-like, default ".") – Directory produced by an Occam2D run. It should contain an Occam2DMesh file, an Occam2DModel file, a data file, a log file, .iter files, and matching .resp files.
iteration (int or None, default None) – Iteration number to load. If None, the highest numbered .iter file is selected. If the requested iteration is unavailable, the loader falls back to the last available iteration file.
verbose (int or bool, default 0) – Verbosity level inherited from OccamBase. Positive values enable progress messages through the instance logger.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created.

Variables:

workdir (pathlib.Path) – Directory scanned by the loader.
log (OccamLog or None) – Parsed convergence log when a log file is found.
mesh (OccamMesh or None) – Parsed finite-element mesh.
model (OccamModel or None) – Parsed model-parameter definition.
best_iter (OccamIter or None) – Selected iteration file. The name is kept for backward compatibility; it may be the requested iteration or last available iteration.
response (OccamResponse or None) – Response file matching the selected iteration when available. If no exact match exists, the loader falls back to the last response file.
data (OccamData or None) – Parsed observed-data file when available.
iter_files (list of pathlib.Path) – All .iter files found in workdir, sorted by embedded iteration number.
resp_files (list of pathlib.Path) – All .resp files found in workdir, sorted by embedded iteration number.
rho_2d (numpy.ndarray of float or None) – Log10-resistivity grid with shape (mesh.n_zcells, mesh.n_xcells). Cells outside the model domain are stored as nan.

Notes

The loader is deliberately tolerant. Missing optional files leave corresponding attributes as None instead of failing immediately. A missing working directory still raises NotADirectoryError because there is no useful scan to perform.

See also

OccamRunner: Runs the executable that produces result files.
OccamLog: Parses convergence information loaded here.
OccamResponse: Parses modeled responses and weighted residuals.
PlotModel: Visualizes the reconstructed model grid.

Examples

Load the latest available iteration:

>>> from pycsamt.models.occam2d import InversionResult
>>> result = InversionResult(workdir="occam_run")
>>> result.final_rms

Load a specific iteration and export the model grid:

>>> result = InversionResult("occam_run", iteration=17)
>>> result.iter2dat("occam_run/final_model.dat")

Access loaded components directly:

>>> result.mesh.n_xcells
>>> result.response.rms

References

[InversionResult-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

[InversionResult-2]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

iter2dat(output_file)#

Write the selected model as a three-column ASCII file.

The exported file contains one row for each finite cell in rho_2d. Columns are x_center, z_center, and log10_rho. Horizontal coordinates are centered around the profile midpoint and depths are positive downward.

This format is useful for external plotting tools and for workflows that expect Bo Yang-style iter2dat output. The values remain in log10-resistivity units:

\[m = \log_{10}(\rho).\]

Parameters:: output_file (path-like) – Destination path for the exported ASCII model. Parent directories are created when needed.
Returns:: Path to the file that was written.
Return type:: pathlib.Path
Raises:: RuntimeError – Raised when the result is not fully loaded and the mesh or reconstructed grid is unavailable.

See also

InversionResult.rho_2d: Grid used to generate the exported values.
PlotModel: Visualizes the same reconstructed model grid.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> result = InversionResult("occam_run")
>>> out = result.iter2dat("occam_run/final_model.dat")

plot_model(**kwargs)#: Plot the reconstructed 2-D resistivity model.

plot_response(**kwargs)#: Plot observed and modeled response curves.

plot_misfit(**kwargs)#: Plot RMS misfit as a function of iteration.

plot_pseudo(**kwargs)#: Plot an observed-data pseudosection.

property final_rms: float#: RMS misfit of the selected iteration.

property n_iterations: int#: Number of iteration files discovered in workdir.

summary()#

Return a short text summary of the loaded inversion.

Return type:: str

class pycsamt.models.occam2d.OccamData(title='pycsamt Occam2D data file', config=None, **kwargs)#

Bases: OccamBase

Represent an Occam2D magnetotelluric data file.

OccamData stores the station list, profile offsets, global frequency table, data-type codes, datum values, and uncertainty values written to OccamDataFile.dat. The object is both a container for parsed files and the product of EDI conversion by from_edi().

The Occam2D data file uses one row for each datum. Apparent resistivity is stored in logarithmic form, while phase is stored in degrees:

\[d_\rho = \log_{10}(\rho_a), \qquad \sigma_d = \frac{\sigma_\rho}{\ln(10)} .\]

For PyCSAMT EDI arrays, TE mode is taken from \(Z_{xy}\) and TM mode is taken from \(Z_{yx}\). TM phase is shifted by \(180^\circ\) so passive-MT \(Z_{yx}\) phases are written in the first quadrant.

Parameters:

title (str, default "pycsamt Occam2D data file") – Free-text title written into the Occam data-file header. Use it to record survey name, processing version, inversion purpose, or other provenance attached to the generated OccamDataFile.dat.
config (OccamConfig, optional) – Configuration object controlling data selection, mesh geometry, startup controls, file names, and executable discovery. It centralizes choices such as modes, error floors, frequency limits, layer counts, cell sizes, target misfit, starting resistivity, and Occam file names. If omitted, a default OccamConfig is created.
verbose (int or bool, default 0) – Verbosity level for progress reporting. 0 or False keeps the object quiet. Positive values enable progress messages through the instance logger; larger values may be used by callers to request more diagnostic detail.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created automatically. Pass an explicit logger when integrating Occam2D objects into an application-level logging setup.

Variables:

format_str (str) – Occam format tag. The current writer uses "OCCAM2MTDATA_1.0".
title (str) – Free-text title written to the data-file header.
config (OccamConfig) – Configuration used for default modes, frequency limits, and error floors during EDI conversion.
sites (list of str) – Station names ordered along the profile. The same order is used by mesh construction and all one-based site indices.
offsets (numpy.ndarray of float, shape (n_sites,)) – Station chainages in metres. Values are sorted from low to high during from_edi().
frequencies (numpy.ndarray of float, shape (n_frequencies,)) – Global frequency table in hertz, sorted from high to low as expected by Occam2D.
data_blocks (numpy.ndarray of float, shape (n_data, 5)) – Data rows with columns site_index, freq_index, type_code, datum, and error. Indices are one-based because they are written directly to Occam files.
elevations (numpy.ndarray of float, shape (n_sites,)) – Per-site elevation in metres above sea level, same order as sites. Populated by from_edi() via pycsamt.topo.extract_elevation() when the EDI/Sites source carries topography; all-zero otherwise. Occam2D’s own file formats have no elevation field of their own – see station_elevations() and has_topography.

Notes

Occam2D type codes distinguish both data kind and component. The common MT rows are 1 for RhoTE, 2 for PhsTE, 5 for RhoTM, and 6 for PhsTM. Additional impedance and tipper codes are exposed through DATA_TYPE_CODES for readers and future writers.

See also

OccamConfig: Supplies default modes, frequency bounds, and error floors.
OccamMesh.from_data: Builds mesh geometry from station offsets in OccamData.
OccamResponse: Reads modeled responses and residuals for the same rows.
InputBuilder: Coordinates writing data, mesh, model, and startup files.

Examples

Build a data file from EDI sites and write it to disk:

>>> from pycsamt.models.occam2d import OccamData
>>> from pycsamt.site import Sites
>>> sites = Sites.from_any("edi")
>>> data = OccamData.from_edi(sites, modes=["TE", "TM"])
>>> data.write("occam_run/OccamDataFile.dat")

Create a synthetic container for tests or scripted workflows:

>>> import numpy as np
>>> from pycsamt.models.occam2d import OccamData
>>> data = OccamData(title="synthetic profile")
>>> data.sites = ["S00", "S01"]
>>> data.offsets = np.array([0.0, 1000.0])
>>> data.frequencies = np.array([100.0, 10.0])
>>> data.data_blocks = np.array([[1, 1, 1, 2.0, 0.05]])

Read an existing Occam data file:

>>> from pycsamt.models.occam2d import OccamData
>>> data = OccamData.read("occam_run/OccamDataFile.dat")
>>> data.n_sites, data.n_frequencies, data.n_data

References

[OccamData-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

[OccamData-2]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

elevations: ndarray#: Per-site elevation (m a.s.l.), same order as sites. Populated by from_edi() via pycsamt.topo when the EDI/Sites source carries topography; all-zero (“no elevation known”) otherwise – Occam2D’s own file formats carry no elevation field at all.

classmethod from_edi(source, modes=None, config=None, title='pycsamt Occam2D data file', **kwargs)#

Build an Occam data object from EDI-derived stations.

This constructor normalizes the accepted input source to a list of site-like objects, estimates station chainages, merges all available frequencies into a common descending frequency table, applies frequency limits, and writes TE/TM apparent-resistivity and phase rows using Occam type codes.

The EDI arrays are interpreted with the convention

\[\mathrm{TE} = Z_{xy}, \qquad \mathrm{TM} = Z_{yx} .\]

For each accepted apparent-resistivity value \(\rho_a\), the stored datum is \(\log_{10}(\rho_a)\). The relative resistivity floor is converted to log10 uncertainty by \(\sigma_d=\sigma_\rho/\ln(10)\). Phase rows use degree errors and the TM phase is shifted by \(180^\circ\).

Parameters:

source (Sites, EDICollection, or iterable) – EDI-derived survey source used to build the Occam data file. Accepted inputs include pycsamt.site.Sites, an EDI collection, or any iterable of site-like objects. Each item must expose frequency, apparent resistivity, and phase arrays. Coordinates are strongly preferred because they allow station offsets to be ordered along profile. When coordinates are absent, fallback spacing is used.
modes (list of str, optional) – Electromagnetic modes written to the data file. Supported values are "TE" for the \(Z_{xy}\) component and "TM" for the \(Z_{yx}\) component. Both apparent resistivity and phase rows are written for each selected mode. If omitted, config.modes is used.
config (OccamConfig, optional) – Configuration object controlling data selection, mesh geometry, startup controls, file names, and executable discovery. It centralizes choices such as modes, error floors, frequency limits, layer counts, cell sizes, target misfit, starting resistivity, and Occam file names. If omitted, a default OccamConfig is created.
title (str, default "pycsamt Occam2D data file") – Free-text title written into the Occam data-file header. Use it to record survey name, processing version, inversion purpose, or other provenance attached to the generated OccamDataFile.dat.
**kwargs – Additional keyword arguments passed to the OccamData constructor. This is commonly used for verbose or logger when progress messages are desired.

Returns:

Populated data object ready to be written as an OCCAM2MTDATA_1.0 file.

Return type:

OccamData

Raises:

ValueError – Raised when the source has no sites, no frequency arrays, no frequencies remain after filtering, or no requested mode is supported.

See also

OccamConfig: Provides default modes, frequency bounds, and error floors.
OccamData.write: Serializes the returned object to OccamDataFile.dat.
OccamMesh.from_data: Uses the returned station offsets to build mesh geometry.

Examples

Build TE and TM rows from a site collection:

>>> from pycsamt.models.occam2d import OccamData
>>> from pycsamt.site import Sites
>>> sites = Sites.from_any("edi")
>>> data = OccamData.from_edi(sites, modes=["TE", "TM"])

Restrict the frequency range through OccamConfig:

>>> from pycsamt.models.occam2d import OccamConfig
>>> from pycsamt.models.occam2d import OccamData
>>> cfg = OccamConfig(freq_min=0.1, freq_max=1000.0)
>>> data = OccamData.from_edi(sites, config=cfg)

Use only TM data with stronger phase floor:

>>> cfg = OccamConfig(modes=["TM"], error_floor_phase=1.0)
>>> data = OccamData.from_edi(sites, config=cfg)

References

classmethod read(path, **kwargs)#

Read an existing OCCAM2MTDATA_1.0 file.

The parser reads the format tag, title, station names, offsets, frequency table, and numeric data block. Site and frequency indices are kept in the one-based form used by Occam2D so a read-write round trip preserves the original file structure.

Parameters:

path (path-like) – Path to an Occam2D input or output file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path. Relative paths are interpreted from the current working directory of the Python process. Readers require the file to exist, while writers create parent directories when the owning method supports output.
**kwargs – Additional keyword arguments forwarded to the OccamData constructor before parsed values are attached. Use this for config, verbose, or logger.

Returns:

Parsed data-file container with arrays populated from path.

Return type:

OccamData

Raises:

FileNotFoundError – Raised when path does not exist.
ValueError – Raised when the format tag is missing or not "OCCAM2MTDATA_1.0".

Examples

>>> from pycsamt.models.occam2d import OccamData
>>> data = OccamData.read("occam_run/OccamDataFile.dat")
>>> data.type_codes

write(path)#

Write this object as an OCCAM2MTDATA_1.0 file.

The writer serializes the current title, station list, offsets, frequency table, and data rows using the Occam2D text layout. Parent directories are created before writing. The method does not modify the object, so it can be used repeatedly for round-trip checks or alternative run directories.

Parameters:: path (path-like) – Path to an Occam2D input or output file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path. Relative paths are interpreted from the current working directory of the Python process. Readers require the file to exist, while writers create parent directories when the owning method supports output.
Returns:: Path to the file that was written.
Return type:: pathlib.Path

See also

OccamData.read: Parses a file written by this method.
InputBuilder.build: Calls this method as part of complete input generation.

Examples

>>> from pycsamt.models.occam2d import OccamData
>>> data = OccamData.read("source/OccamDataFile.dat")
>>> written = data.write("copy/OccamDataFile.dat")

property n_sites: int#

property n_frequencies: int#

property n_data: int#

property type_codes: ndarray#: Unique data-type codes present in this dataset.

property has_topography: bool#: True when elevations carries real, non-zero relief.

station_elevations()#

Return {station_name: elevation_m} for stations with topography.

Built from sites/elevations (populated by from_edi() via pycsamt.topo). Empty when the source carried no real elevation. The returned mapping matches the station_elevations parameter accepted by pycsamt.format.adapters.occam2d.occam2d_to_pcsf().

Return type:: dict of str to float

class pycsamt.models.occam2d.OccamMesh(config=None, **kwargs)#

Bases: OccamBase

Represent the Occam2D PW2D finite-element mesh.

OccamMesh stores the two-dimensional grid consumed by the Occam2D forward solver. Horizontal cell widths define the profile direction. Vertical widths define air and earth layers, and character rows mark whether cells are fixed, air, or boundary cells in the PW2D mesh format.

Node coordinates are cumulative sums of cell widths:

\[x_j = \sum_{i=0}^{j-1} \Delta x_i, \qquad z_k = \sum_{i=0}^{k-1} \Delta z_i.\]

Depth \(z\) is positive downward. Mesh construction uses station offsets from OccamData, horizontal padding on both profile ends, optional air layers, and a geometrically expanding earth-layer thickness sequence.

Parameters:

config (OccamConfig, optional) – Configuration object controlling the number of active layers, number of air layers, near-surface cell sizes, depth scaling, and horizontal padding. If omitted, a default OccamConfig is created.
verbose (int or bool, default 0) – Verbosity level inherited from OccamBase. Positive values enable progress messages through the instance logger.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created.

Variables:

comment (str) – First line of the mesh file, usually a provenance comment beginning with "MESH FILE".
x_widths (numpy.ndarray of float, shape (n_xcells,)) – Horizontal cell widths in metres.
z_widths (numpy.ndarray of float, shape (n_zcells,)) – Vertical layer thicknesses in metres.
x_nodes (numpy.ndarray of float, shape (n_xcells + 1,)) – Cumulative horizontal node positions in metres.
z_nodes (numpy.ndarray of float, shape (n_zcells + 1,)) – Cumulative depth node positions in metres, positive downward.
cell_rows (list[str]) – Raw PW2D cell-type rows. Each character encodes the cell type at one horizontal position. The "?" character marks cells that may contribute to free inversion parameters.
n_airlayers (int) – Number of rows treated as air layers.

Notes

The mesh file stores widths rather than absolute node coordinates. x_nodes and z_nodes are reconstructed by cumulative summation when reading or building a mesh. The generated mesh uses seven padding cells on each side to match the boundary-column code used by OccamModel.from_mesh().

See also

OccamData: Provides station offsets used to build the mesh.
OccamModel.from_mesh: Converts mesh cells into inversion-parameter columns.
InputBuilder: Builds data, mesh, model, and startup files together.

Examples

Build a mesh from an Occam data file:

>>> from pycsamt.models.occam2d import OccamData
>>> from pycsamt.models.occam2d import OccamMesh
>>> data = OccamData.read("occam_run/OccamDataFile.dat")
>>> mesh = OccamMesh.from_data(data)
>>> mesh.write("occam_run/Occam2DMesh")

Read an existing PW2D mesh:

>>> from pycsamt.models.occam2d import OccamMesh
>>> mesh = OccamMesh.read("occam_run/Occam2DMesh")
>>> mesh.n_xcells, mesh.n_zcells

References

[OccamMesh-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

[OccamMesh-2]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

N_PAD: int = 7#: Fixed horizontal padding-cell count on each profile end. Matches the boundary code 7 used by OccamModel.from_mesh() and the hardcoded n_pad in from_data() – an architectural invariant, not something inferred per-mesh, so it also applies to meshes rebuilt via read().

cell_centers_survey_x()#

Horizontal cell-center coordinates in survey (offset) space.

x_widths/x_nodes are zero-based at the outer edge of the left padding, not at the survey’s own zero offset (see OccamData’s offsets, which starts at 0 after normalization). Comparing or resampling a solved model against anything expressed in survey/offset coordinates – true models, AI display grids, station chainage – must shift by the total left-padding width, or the comparison silently lands entirely inside the padding zone.

Returns:: Cell-center x-coordinates, in metres, in the same zero-based convention as OccamData’s offsets.
Return type:: numpy.ndarray of float, shape (n_xcells,)

classmethod from_data(data, config=None, **kwargs)#

Build a PW2D mesh from Occam data offsets.

The method creates a finite-element mesh spanning the profile described by data.offsets. It uses seven padding cells on each side, station-zone cells near the configured horizontal cell size, optional air layers, and geometrically expanding earth layers.

The horizontal padding is chosen to match the boundary columns used by OccamModel.from_mesh(). Interior station-zone cells are adjusted so the number of cells remains compatible with the model parameter grouping.

Parameters:

data (OccamData) – Populated data object. Its offsets array must contain station chainages in metres. Offsets are sorted before mesh construction.
config (OccamConfig, optional) – Configuration object controlling mesh geometry. The builder uses cell_size_horizontal, n_airlayers, n_layers, max_depth, cell_size_vertical_top, and depth_scale. Earth layers expand geometrically from cell_size_vertical_top by depth_scale and stop – truncating the last layer if needed – once cumulative depth reaches max_depth (default 1500 m) or n_layers layers have been added, whichever comes first. If omitted, a default OccamConfig is created.
**kwargs – Additional keyword arguments forwarded to the OccamMesh constructor. Use this for verbose or logger.

Returns:

Mesh object ready to be written as Occam2DMesh or passed to OccamModel.from_mesh().

Return type:

OccamMesh

Raises:

ValueError – Raised when the data object contains no station offsets.

See also

OccamData.from_edi: Creates the offsets used by this method.
OccamModel.from_mesh: Builds the inversion-parameter mapping.

Examples

>>> from pycsamt.models.occam2d import OccamData
>>> from pycsamt.models.occam2d import OccamMesh
>>> data = OccamData.read("OccamDataFile.dat")
>>> mesh = OccamMesh.from_data(data)
>>> mesh.n_airlayers

classmethod read(path, **kwargs)#

Read an existing Occam2DMesh PW2D file.

The reader parses the comment line, control line, horizontal widths, vertical widths, air-layer count, and cell-type character rows. Node arrays are rebuilt from cumulative sums of widths.

Parameters:

path (path-like) – Path to the mesh file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
**kwargs – Additional keyword arguments forwarded to the OccamMesh constructor before parsed values are attached. Use this for config, verbose, or logger.

Returns:

Parsed mesh container with widths, nodes, and cell rows populated.

Return type:

OccamMesh

Raises:

FileNotFoundError – Raised when path does not exist.
ValueError – Raised when the file is too short or the control line cannot be parsed.

Examples

>>> from pycsamt.models.occam2d import OccamMesh
>>> mesh = OccamMesh.read("occam_run/Occam2DMesh")
>>> mesh.x_nodes.shape

write(path)#

Write this mesh in PW2D format.

The writer serializes the current comment, control values, horizontal widths, vertical widths, and cell-type rows to the native Occam2D mesh format. Parent directories are created before writing.

Parameters:: path (path-like) – Destination path for the mesh file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
Returns:: Path to the file that was written.
Return type:: pathlib.Path

See also

OccamMesh.read: Parses mesh files written by this method.
InputBuilder.build: Calls this method during input-file generation.

Examples

>>> from pycsamt.models.occam2d import OccamMesh
>>> mesh = OccamMesh.read("source/Occam2DMesh")
>>> written = mesh.write("copy/Occam2DMesh")

property n_xcells: int#: Number of horizontal cells.

property n_zcells: int#: Number of vertical layers.

property n_params: int#: Number of free model parameters (cells coded as '?').

pycsamt.models.occam2d.resample_rho_to_grid(rho_2d, mesh, x, z)#

Resample a solved Occam2D resistivity model onto a regular grid.

rho_2d is defined on the mesh’s own irregular, padding-inclusive cell grid. x/z are typically a regular display or comparison grid expressed in survey coordinates (station-offset-relative, e.g. an AI benchmark’s true-model grid or a fixed display grid) – not the mesh’s own zero-at-outer-padding coordinate system. Resampling must therefore go through OccamMesh.cell_centers_survey_x(), not the mesh’s raw node coordinates, or every query lands inside the (typically tens of kilometres of) horizontal padding and is clamped to a single, flat, constant-per-row value regardless of the model’s real structure.

Parameters:

rho_2d (numpy.ndarray, shape (n_zcells, n_xcells)) – Resistivity (or log-resistivity) on the mesh’s own grid, including any air-layer rows at the top (stripped internally using mesh.n_airlayers).
mesh (OccamMesh) – Mesh the model was solved on.
x (numpy.ndarray) – Target cell-center coordinates, in metres, in survey (offset) and depth (positive-down) coordinates respectively.
z (numpy.ndarray) – Target cell-center coordinates, in metres, in survey (offset) and depth (positive-down) coordinates respectively.

Returns:

rho_2d resampled onto the x/z grid.

Return type:

numpy.ndarray, shape (z.size, x.size)

class pycsamt.models.occam2d.OccamModel(name='MODEL MADE BY PYCSAMT', description='SMOOTH INVERSION', config=None, **kwargs)#

Bases: OccamBase

Represent the Occam2D model-parameter definition.

OccamModel links a finite-element OccamMesh to the inversion parameter vector used by Occam2D. The model file does not store resistivity values. Instead, it defines how mesh cells are grouped into free or fixed parameters. The startup and iteration files then store one value for each parameter counted by n_params.

Each model layer contains integer column codes. Boundary code 7 marks fixed edge columns tied to the binding value, while active even codes represent free inversion columns [OccamModel-1]. If \(p_j\) is the code for one model column, then the number of mesh cells represented by that column is encoded by the code value itself. The PyCSAMT builder uses 2 for interior columns and 7 for the two boundary columns:

\[\mathbf{p} = [7,\;2,\;2,\;\ldots,\;2,\;7].\]

Parameters:

name (str, default "MODEL MADE BY PYCSAMT") – Model name written to the Model Name header field. Use this for a short label that identifies the model family, processing run, or inversion setup.
description (str, default "SMOOTH INVERSION") – Description written to the Description header field. It is intended for human-readable provenance and is preserved when the model is written to disk.
config (OccamConfig, optional) – Configuration object used for default file names and related Occam2D settings. If omitted, a default OccamConfig is created.
verbose (int or bool, default 0) – Verbosity level inherited from OccamBase. Positive values enable progress messages through the instance logger.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created.

Variables:

format_str (str) – Occam model format tag. The writer uses "OCCAM2MTMOD_1.0".
name (str) – Model name written to the file header.
description (str) – Human-readable model description.
config (OccamConfig) – Configuration used by this object.
mesh_file (str) – Filename of the associated mesh used by the model file, usually "Occam2DMesh".
mesh_type (str) – Mesh type string written to the header. Occam2D PW2D meshes use "PW2D".
statics_file (str) – Optional static-shift file. The default "none" means no static correction file is referenced.
prejudice_file (str) – Optional prejudice model file. The default "none" means no prejudice file is referenced.
binding_offset (float, default 0.0) – Horizontal offset of the binding column. This is the reference value used by boundary columns.
n_layers (int) – Number of active model layers after the header. It is normally the number of non-air mesh rows.
layers (list of dict) –
Per-layer parameter specification. Each entry has the following keys:

n_mergeint
Number of mesh z-rows merged into this layer.

n_colsint
Number of model columns in this layer.

paramsnumpy.ndarray of int, shape (n_cols,)
Parameter codes for each model column. Code 7 marks a boundary column; active even values mark free inversion parameters.
n_exceptions (int) – Number of exception records parsed from an existing model file, in the file’s native (possibly negative) sign convention. Ignored by write() whenever exceptions is non-empty; the written count is then always -len(exceptions).
exceptions (list of (int, int, float)) – Roughness-penalty exceptions as one-based (brick_i, brick_j, expen) triples, in the brick-pair convention (see “Roughness-penalty exceptions” above). Empty by default. Indices should reference bricks produced by the same parameter traversal as pycsamt.ai.inversion.mapping2d.map_ai_grid_to_occam(); pycsamt.ai.inversion.mapping2d.build_occam_parameter_adjacency() enumerates the valid adjacent pairs for this model.

Notes

n_params is the sum of n_cols over all model layers. This value must match the Param Count in the startup and iteration files. n_free_params excludes boundary columns with code 7.

See also

OccamMesh: Defines finite-element cells grouped by this model.
OccamStartup.from_model: Creates an initial vector with matching size.
InputBuilder: Builds data, mesh, model, and startup files together.

Examples

Build a model definition from an existing mesh:

>>> from pycsamt.models.occam2d import OccamMesh
>>> from pycsamt.models.occam2d import OccamModel
>>> mesh = OccamMesh.read("occam_run/Occam2DMesh")
>>> model = OccamModel.from_mesh(mesh)
>>> model.n_params

Read and write an existing model file:

>>> from pycsamt.models.occam2d import OccamModel
>>> model = OccamModel.read("occam_run/Occam2DModel")
>>> model.write("copy/Occam2DModel")

Create a custom empty container for tests:

>>> from pycsamt.models.occam2d import OccamModel
>>> model = OccamModel(name="SYNTHETIC MODEL")
>>> model.n_layers, model.n_params

References

[OccamModel-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

[OccamModel-2]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

classmethod from_mesh(mesh, config=None, **kwargs)#

Build a model definition from a populated mesh.

The method converts a finite-element mesh into the column mapping required by OCCAM2MTMOD_1.0. Air rows are ignored. Each remaining earth row becomes one model layer with n_merge = 1. Horizontally, the model uses fixed boundary columns and active interior columns:

\[\mathbf{p} = [7,\;2,\;2,\;\ldots,\;2,\;7].\]

The leading and trailing 7 codes represent seven mesh cells each at the profile boundaries. Interior 2 codes represent two mesh cells per free model column. Therefore the total horizontal cell count is

\[n_x = 7 + 2n_i + 7,\]

where \(n_i\) is the number of interior model columns. Meshes built by OccamMesh.from_data() are constructed to satisfy this layout.

Parameters:

mesh (OccamMesh) – Populated mesh object defining horizontal and vertical finite-element cells. It must provide n_xcells, n_zcells, and n_airlayers. Air layers are excluded from the model; all other z-cells become inversion layers.
config (OccamConfig, optional) – Configuration object used for file names and related Occam2D defaults. If omitted, a default OccamConfig is created.
**kwargs – Additional keyword arguments forwarded to the OccamModel constructor. This is commonly used for name, description, verbose, or logger.

Returns:

Model-definition object ready to be written as an Occam2DModel file. The returned object has n_layers equal to the number of active earth rows and layers populated with parameter-code arrays.

Return type:

OccamModel

Raises:

ValueError – Raised when the mesh has no active earth layers or fewer than fourteen horizontal cells. Fourteen cells are required for the two seven-cell boundary columns.

See also

OccamMesh.from_data: Builds meshes that match this parameterization.
OccamModel.write: Serializes the returned model definition.
OccamStartup.from_model: Creates a startup vector with matching parameter count.

Examples

Build a model from a mesh read from disk:

>>> from pycsamt.models.occam2d import OccamMesh
>>> from pycsamt.models.occam2d import OccamModel
>>> mesh = OccamMesh.read("occam_run/Occam2DMesh")
>>> model = OccamModel.from_mesh(mesh)
>>> model.n_layers

Pass metadata through to the model header:

>>> model = OccamModel.from_mesh(
...     mesh,
...     name="PROFILE A MODEL",
...     description="smooth TE-TM inversion",
... )

classmethod read(path, **kwargs)#

Read an existing OCCAM2MTMOD_1.0 model file.

The reader parses the model header and the per-layer parameter-code blocks. Numeric header values are cast to int or float where appropriate. Layer params arrays are stored as numpy.int32 for comparison with generated model mappings.

Parameters:

path (path-like) – Path to an Occam2D model file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
**kwargs – Additional keyword arguments forwarded to the OccamModel constructor before parsed values are attached. Use this for config, verbose, or logger.

Returns:

Parsed model-definition container with header fields and layer mappings populated from path.

Return type:

OccamModel

Raises:

FileNotFoundError – Raised when path does not exist.
ValueError – Raised when the format tag is missing or is not "OCCAM2MTMOD_1.0".

See also

OccamModel.write: Writes model definitions in the same format.
OccamStartup.read: Reads startup or iteration files that depend on the same parameter count.

Examples

>>> from pycsamt.models.occam2d import OccamModel
>>> model = OccamModel.read("occam_run/Occam2DModel")
>>> model.n_layers, model.n_params

write(path)#

Write this model in OCCAM2MTMOD_1.0 format.

The writer serializes the current header fields and layer mappings to the Occam2D model layout. Parent directories are created before writing. The object is not modified, so the same instance can be written to multiple run directories.

Parameters:: path (path-like) – Destination path for the model file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
Returns:: Path to the file that was written.
Return type:: pathlib.Path

See also

OccamModel.read: Parses model files written by this method.
InputBuilder.build: Calls this method during input generation.

Examples

>>> from pycsamt.models.occam2d import OccamModel
>>> model = OccamModel.read("source/Occam2DModel")
>>> written = model.write("copy/Occam2DModel")

property n_params: int#: Total model cells, equal to iter Param Count.

property n_free_params: int#: Number of free (non-boundary) model cells.

class pycsamt.models.occam2d.OccamStartup(config=None, description='startup created by pycsamt', **kwargs)#

Bases: OccamBase

Represent an Occam2D startup control file.

OccamStartup stores the iteration-zero OCCAMITER_FLEX file passed to the Occam2D executable. It defines run controls, file references, inversion options, and the initial model vector. Unlike .iter files produced by the solver, a valid startup file has Iteration: 0.

The startup parameter vector is initialized as a uniform half-space:

\[m_i = \log_{10}(\rho_0), \qquad i = 1, \ldots, N_p.\]

Here \(\rho_0\) is config.initial_rho and \(N_p\) is the number of model parameters defined by OccamModel.

Parameters:

config (OccamConfig, optional) – Configuration object providing file names, inversion controls, starting resistivity, target misfit, roughness settings, and debug level. If omitted, a default OccamConfig is created.
description (str, default "startup created by pycsamt") – Description written to the Description header. Use it to record the purpose or provenance of the run.
verbose (int or bool, default 0) – Verbosity level inherited from OccamBase. Positive values enable progress messages through the instance logger.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created.

Variables:

format_str (str) – File format tag, always "OCCAMITER_FLEX".
description (str) – Human-readable startup description.
model_file (str) – Model file name referenced by the startup file.
data_file (str) – Data file name referenced by the startup file.
datetime_str (str) – Creation or file timestamp string.
max_iterations (int) – Maximum number of iterations requested from Occam2D.
target_misfit (float) – Target normalized RMS misfit.
roughness_type (int) – Roughness penalty type written to the startup file.
diagonal_penalties (int) – Flag controlling diagonal roughness penalties.
stepsize_cut_count (int) – Maximum number of step-size cuts in a line search.
debug_level (int) – Debug verbosity passed to the Fortran executable.
iteration (int) – Always 0 for a valid Startup file.
lagrange_value (float) – Initial Lagrange multiplier written as Lagrange Value.
roughness_value (float) – Initial roughness value written before the first run.
misfit_value (float) – Initial misfit value written before the first run.
misfit_reached (bool) – Whether the target misfit has already been reached. Startup files normally use False.
n_params (int) – Number of model parameters. This must match the model file and the length of param_values.
param_values (numpy.ndarray of float, shape (n_params,)) – Initial log10-resistivity values. Values are uniform after from_model() and equal to log10(config.initial_rho).

Notes

OccamStartup writes the same flexible iteration format that Occam later uses for .iter files. The distinction is semantic: startup files carry Iteration: 0 and are input to the solver, while OccamIter files carry non-zero iteration numbers and are output from the solver.

See also

OccamModel: Provides the parameter count for the startup vector.
OccamIter: Reads iteration files produced after running Occam2D.
OccamRunner: Launches the executable with this startup file.

Examples

Build a startup file from a model definition:

>>> from pycsamt.models.occam2d import OccamModel
>>> from pycsamt.models.occam2d import OccamStartup
>>> model = OccamModel.read("occam_run/Occam2DModel")
>>> startup = OccamStartup.from_model(model)
>>> startup.write("occam_run/Startup")

Read an existing startup file:

>>> from pycsamt.models.occam2d import OccamStartup
>>> startup = OccamStartup.read("occam_run/Startup")
>>> startup.n_params

References

[OccamStartup-1]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

[OccamStartup-2]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

classmethod from_model(model, config=None, **kwargs)#

Build a startup object from a model definition.

The method uses model.n_params to size the startup vector and fills every entry with the log10 value of starting half-space resistivity:

\[m_i = \log_{10}(\rho_0), \qquad i = 1,\ldots,N_p.\]

This produces the standard smooth-inversion initial model: a homogeneous half-space whose value is later updated by the Occam solver.

Parameters:

model (OccamModel) – Populated model-definition object. It must contain at least one parameter so the vector can be sized consistently with the Occam2DModel file.
config (OccamConfig, optional) – Configuration object providing initial_rho, model and data file names, iteration controls, and inversion settings. If omitted, a default OccamConfig is created.
**kwargs – Additional keyword arguments forwarded to the OccamStartup constructor. Use this for description, verbose, or logger.

Returns:

Startup object with n_params and uniform param_values populated.

Return type:

OccamStartup

Raises:

ValueError – Raised when model.n_params is not positive.

See also

OccamStartup.write: Serializes the generated startup object.
OccamModel: Supplies the parameter count used here.

Examples

>>> from pycsamt.models.occam2d import OccamModel
>>> from pycsamt.models.occam2d import OccamStartup
>>> model = OccamModel.read("occam_run/Occam2DModel")
>>> startup = OccamStartup.from_model(model)
>>> startup.param_values.shape

classmethod read(path, **kwargs)#

Read an Occam2D startup file.

The reader parses an OCCAMITER_FLEX file and then validates that its Iteration header is zero. Use OccamIter.read() for non-zero iteration files produced by the solver.

Parameters:

path (path-like) – Path to the startup file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
**kwargs – Additional keyword arguments forwarded to the OccamStartup constructor before parsed values are attached.

Returns:

Parsed startup object with header fields and parameter vector populated.

Return type:

OccamStartup

Raises:

FileNotFoundError – Raised when path does not exist.
ValueError – Raised when the file is not OCCAMITER_FLEX or has a non-zero iteration number.

Examples

>>> from pycsamt.models.occam2d import OccamStartup
>>> startup = OccamStartup.read("occam_run/Startup")
>>> startup.iteration

write(path)#

Write this startup file in OCCAMITER_FLEX format.

Parameters:: path (path-like) – Destination path for the startup file. Parent directories are created when needed.
Returns:: Path to the file that was written.
Return type:: pathlib.Path

See also

OccamStartup.read: Parses files written by this method.
OccamRunner: Passes the written startup file to the executable.

class pycsamt.models.occam2d.OccamPrejudice(parameter_indices=None, target_values=None, weights=None, config=None, **kwargs)#

Bases: OccamBase

Represent a sparse Occam2D model-prejudice file.

OccamPrejudice stores selected Occam parameter indices, their preferred log10-resistivity values, and non-negative penalty weights. The object follows the same container and I/O conventions as OccamData, OccamModel, and OccamStartup.

For a target \(m_j^{\mathrm{target}}\) and weight \(w_j\), the intended local penalty is proportional to

\[w_j^2 \left(m_j-m_j^{\mathrm{target}}\right)^2.\]

The bundled solver assembles its prejudice terms using prewt**2 in the Hessian and prewt*premod in the right-hand side. Consequently, write() encodes the native prejudice field as target_values * weights. read() reverses that encoding so users always work with physical target values.

Parameters:

parameter_indices (iterable of int, optional) – One-based Occam model-parameter indices. Values must be positive and unique. If omitted, an empty prejudice object is created.
target_values (iterable of float, optional) – Preferred model values in log10 resistivity, ordered like parameter_indices. Values must be finite. If omitted, an empty array is used.
weights (iterable of float, optional) – Non-negative prejudice weights, ordered like parameter_indices. Values must be finite. A zero weight leaves the associated target inactive.
config (OccamConfig, optional) – Configuration object associated with the Occam2D project. The prejudice format has no configuration fields of its own, but retaining the project configuration matches the other Occam2D file containers. If omitted, a default configuration is created.
verbose (int or bool, default 0) – Verbosity level inherited from OccamBase. Positive values enable progress messages through the instance logger.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created.

Variables:

format_str (str) – File-format identifier. The supported value is "OCCAM2MTPREJ_2.0".
parameter_indices (numpy.ndarray of int, shape (n_prejudiced,)) – One-based indices of prejudiced Occam model parameters.
target_values (numpy.ndarray of float, shape (n_prejudiced,)) – Decoded physical targets in log10 resistivity.
weights (numpy.ndarray of float, shape (n_prejudiced,)) – Non-negative native penalty weights.
config (OccamConfig) – Occam2D project configuration associated with the object.
path (pathlib.Path or None) – Most recent path read or written. Inherited from OccamBase.

Notes

The public targets are not the raw second column stored in the native file when a weight differs from one. Directly editing that column without applying the encoding can shift the effective penalty centre.

Zero-weight records are accepted for controlled experiments and round-trip fidelity. from_dense() omits them by default to keep production prejudice files sparse.

See also

OccamModel: References a prejudice file through prejudice_file.
OccamStartup: Stores the model vector to which prejudice penalties apply.
InputBuilder: Creates the remaining Occam2D input files.

Examples

Create and write a sparse prejudice definition:

>>> from pycsamt.models.occam2d import OccamPrejudice
>>> prejudice = OccamPrejudice(
...     parameter_indices=[2, 5],
...     target_values=[1.5, 2.3],
...     weights=[2.0, 0.5],
... )
>>> prejudice.write("occam_run/DUHIPrejudice")

Build sparse records from dense model vectors:

>>> prejudice = OccamPrejudice.from_dense(
...     target_values=[1.5, 2.0, 2.5],
...     weights=[4.0, 0.0, 1.0],
... )
>>> prejudice.parameter_indices.tolist()
[1, 3]

Read a file and inspect its decoded target values:

>>> restored = OccamPrejudice.read("occam_run/DUHIPrejudice")
>>> restored.n_prejudiced
2

References

[OccamPrejudice-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

[OccamPrejudice-2]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

classmethod from_dense(target_values, weights, *, include_zero_weight=False, config=None, **kwargs)#

Build sparse prejudice records from dense model vectors.

Dense input vectors are interpreted in Occam parameter order. Their zero-based array positions are converted to one-based solver indices. By default, entries with zero weight are omitted from the sparse result.

Parameters:

target_values (iterable of float) – Preferred log10-resistivity value for every Occam model parameter.
weights (iterable of float) – Non-negative prejudice weight for every Occam model parameter. The vector must have the same length as target_values.
include_zero_weight (bool, default False) – If True, preserve inactive zero-weight records. Otherwise, omit them from the sparse object.
config (OccamConfig, optional) – Occam2D project configuration attached to the returned object. If omitted, a default configuration is created.
**kwargs – Additional keyword arguments forwarded to the OccamPrejudice constructor. Use this for verbose or logger.

Returns:

Sparse prejudice object in one-based Occam parameter order.

Return type:

OccamPrejudice

Raises:

ValueError – Raised when the dense vectors have different lengths or contain invalid target or weight values.

See also

OccamPrejudice.write: Encodes and writes the sparse result.
OccamPrejudice.validate_parameter_count: Checks the indices against a model parameter count.

Examples

>>> from pycsamt.models.occam2d import OccamPrejudice
>>> prejudice = OccamPrejudice.from_dense(
...     [1.5, 2.0, 2.5],
...     [4.0, 0.0, 1.0],
... )
>>> prejudice.parameter_indices.tolist()
[1, 3]

validate()#

Validate the current sparse prejudice records.

Validation checks vector lengths, one-based unique indices, finite targets, and finite non-negative weights. The method does not require at least one record; empty prejudice objects are valid containers.

Returns:: The validated object. Returning self supports fluent preparation workflows.
Return type:: OccamPrejudice
Raises:: ValueError – Raised when the arrays have different lengths, an index is non-positive or duplicated, a target is non-finite, or a weight is non-finite or negative.

See also

OccamPrejudice.validate_parameter_count: Performs model-size validation after record validation.
OccamPrejudice.write: Calls this method before serialization.

Examples

>>> prejudice = OccamPrejudice([1], [2.0], [1.0])
>>> prejudice.validate() is prejudice
True

validate_parameter_count(n_params)#

Validate prejudice indices against an Occam model size.

Parameters:

n_params (int) – Total number of parameters declared by the Occam model and startup files. It must be positive.

Returns:

The validated object.

Return type:

OccamPrejudice

Raises:

TypeError – Raised when n_params is not an integer.
ValueError – Raised when n_params is not positive or a prejudice index exceeds it.

See also

OccamModel.n_params: Supplies the expected model parameter count.
OccamPrejudice.validate: Checks record-level constraints first.

Examples

>>> prejudice = OccamPrejudice([1, 3], [2.0, 2.5], [1, 1])
>>> prejudice.validate_parameter_count(3) is prejudice
True

classmethod read(path, config=None, **kwargs)#

Read and decode an OCCAM2MTPREJ_2.0 file.

The native second column is divided by its positive weight to recover the public physical target. For a zero-weight record, the raw value is retained because the record has no numerical influence and cannot be uniquely decoded.

Parameters:

path (path-like) – Path to the Occam prejudice file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
config (OccamConfig, optional) – Occam2D project configuration attached to the returned object. If omitted, a default configuration is created.
**kwargs – Additional keyword arguments forwarded to the OccamPrejudice constructor. Use this for verbose or logger.

Returns:

Parsed container with decoded target values and path set to the source file.

Return type:

OccamPrejudice

Raises:

FileNotFoundError – Raised when path does not exist.
ValueError – Raised when the header is invalid, the declared record count does not match the file, a record does not have three columns, or parsed values fail validation.

See also

OccamPrejudice.write: Applies the inverse encoding during serialization.
OccamPrejudice.from_dense: Creates sparse records from dense vectors.

Examples

>>> from pycsamt.models.occam2d import OccamPrejudice
>>> prejudice = OccamPrejudice.read("DUHIPrejudice")
>>> prejudice.path.name
'DUHIPrejudice'

write(path)#

Write this object in OCCAM2MTPREJ_2.0 format.

The writer validates the current records, converts each public target to the solver-native value target * weight, creates parent directories, and stores the destination on path. A zero-weight record retains its public target in the native column for round-trip fidelity; the solver ignores that record.

Parameters:: path (path-like) – Destination path for the prejudice file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
Returns:: Path to the file that was written.
Return type:: pathlib.Path
Raises:: ValueError – Raised when the current records fail validation.

See also

OccamPrejudice.read: Reads and decodes files written by this method.
OccamModel.prejudice_file: References the written file from the Occam model.

Examples

>>> prejudice = OccamPrejudice([1], [1.5], [2.0])
>>> written = prejudice.write("run/DUHIPrejudice")
>>> written.name
'DUHIPrejudice'

property n_prejudiced: int#

Return the number of sparse prejudice records.

Returns:: Number of entries in parameter_indices, equivalent to the native Param Count header value.
Return type:: int

Examples

>>> OccamPrejudice([1, 3], [2.0, 2.5], [1, 1]).n_prejudiced
2

property native_values: ndarray#

Return solver-encoded prejudice values without writing.

Positive-weight targets are multiplied by their weights. A zero-weight entry retains its target value for round-trip fidelity. The returned array is newly allocated and modifying it does not change this object.

Returns:: Values written in the second native file column.
Return type:: numpy.ndarray of float, shape (n_prejudiced,)

See also

OccamPrejudice.target_values: Decoded physical target values.
OccamPrejudice.write: Serializes these encoded values.

Examples

>>> prejudice = OccamPrejudice([1], [1.5], [2.0])
>>> prejudice.native_values.tolist()
[3.0]

class pycsamt.models.occam2d.OccamIter(**kwargs)#

Bases: OccamBase

Represent an Occam2D iteration file.

OccamIter reads OCCAMITER_FLEX files written by the Occam2D executable after one or more inversion iterations. These files have the same structural format as Startup but carry Iteration values greater than zero and store the accepted model vector.

The parameter vector is stored in log10-resistivity units. Physical resistivity is recovered as

\[\rho_i = 10^{m_i},\]

where \(m_i\) is the stored value for parameter \(i\).

Parameters:

verbose (int or bool, default 0) – Verbosity level inherited from OccamBase. Positive values enable progress messages through the instance logger.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created.

Variables:

format_str (str) – File format tag, usually "OCCAMITER_FLEX".
description (str) – Iteration description written by Occam.
model_file (str) – Model file referenced by this iteration.
data_file (str) – Data file referenced by this iteration.
datetime_str (str) – Date and time string written by Occam.
max_iterations (int) – Iteration limit stored in the file.
target_misfit (float) – Target normalized RMS misfit.
roughness_type (int) – Roughness penalty type.
diagonal_penalties (int) – Diagonal penalty flag.
stepsize_cut_count (int) – Maximum number of line-search step-size cuts.
debug_level (int) – Debug verbosity setting.
iteration (int) – Iteration number. Valid .iter files use values greater than zero.
lagrange_value (float) – Lagrange multiplier accepted at this iteration.
roughness_value (float) – Model roughness value reported by Occam.
misfit_value (float) – Normalized RMS misfit at this iteration.
misfit_reached (bool) – True if the target misfit was achieved.
n_params (int) – Number of model parameters.
param_values (numpy.ndarray of float, shape (n_params,)) – Accepted log10-resistivity values for this iteration.

See also

OccamStartup: Represents the corresponding iteration-zero file.
InversionResult: Selects iteration files from a run directory.
OccamResponse: Reads the response file for the same iteration.

Examples

Read an iteration file and convert to resistivity:

>>> from pycsamt.models.occam2d import OccamIter
>>> iteration = OccamIter.read("occam_run/ITER17.iter")
>>> rho = iteration.to_resistivity()

Inspect log10-resistivity statistics:

>>> iteration.log10_rho_stats

References

[OccamIter-1]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

[OccamIter-2]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

classmethod read(path, **kwargs)#

Read an Occam2D .iter file.

The reader parses an OCCAMITER_FLEX file and validates that the Iteration value is non-zero. Startup files should be loaded with OccamStartup.read().

Parameters:

path (path-like) – Path to the iteration file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
**kwargs – Additional keyword arguments forwarded to the OccamIter constructor before parsed values are attached.

Returns:

Parsed iteration object with header fields and parameter vector populated.

Return type:

OccamIter

Raises:

FileNotFoundError – Raised when path does not exist.
ValueError – Raised when the file is not OCCAMITER_FLEX or has Iteration: 0.

Examples

>>> from pycsamt.models.occam2d import OccamIter
>>> iteration = OccamIter.read("ITER17.iter")
>>> iteration.misfit_value

write(path)#

Write this iteration in OCCAMITER_FLEX format.

Parameters:: path (path-like) – Destination path for the iteration file. Parent directories are created when needed.
Returns:: Path to the file that was written.
Return type:: pathlib.Path

to_resistivity()#

Return resistivity values from log10 parameters.

Returns:: Resistivity values in ohm metres computed as \(10^m\), where \(m\) is each element of param_values.
Return type:: numpy.ndarray of float

property log10_rho_stats: dict#: Return summary statistics for log10 resistivity.

class pycsamt.models.occam2d.OccamResponse(**kwargs)#

Bases: OccamBase

Represent an Occam2D response file.

OccamResponse stores the forward response written by the Occam2D executable for one inversion iteration. The response table has one row per datum and seven columns: site index, frequency index, type code, error-floor value, observed datum, modeled datum, and weighted residual.

The response residual is already weighted by the data uncertainty used by Occam. The global RMS misfit is therefore computed directly from the residual column:

\[\mathrm{RMS} = \sqrt{ \frac{1}{N} \sum_{i=1}^{N} r_i^2 },\]

where \(r_i\) is the weighted residual for datum \(i\) and \(N\) is the number of response rows. Values near one are often consistent with data errors that are neither under-estimated nor over-estimated [OccamResponse-1].

Parameters:

verbose (int or bool, default 0) – Verbosity level inherited from OccamBase. Positive values enable progress messages through the instance logger.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created.

Variables:

data (numpy.ndarray of float, shape (n_data, 7)) – Full raw table from the .resp file. Columns are site_index, freq_index, type_code, error_floor, observed, modeled, and residual.
observed (numpy.ndarray of float, shape (n_data,)) – Observed data values from column 4. Values follow the datum convention of OccamData: log10 apparent resistivity for rho rows and degrees for phase rows.
modeled (numpy.ndarray of float, shape (n_data,)) – Forward-model predictions from column 5, ordered in the same row order as observed.
residuals (numpy.ndarray of float, shape (n_data,)) – Weighted residuals from column 6. These values are the residuals used to compute rms.
rms (float) – Root-mean-square weighted residual for all response rows. Empty objects use 0.0.

Notes

Response files do not include a header. The parser accepts any line with seven numeric columns and skips non-numeric lines. Site and frequency indices are one-based to match the Occam data file.

See also

OccamData: Defines the observed data rows and type codes.
InversionResult: Loads the response for a selected iteration.
Plot2D.response: Visualizes observed and modeled response curves.

Examples

Read a response file and inspect its global RMS:

>>> from pycsamt.models.occam2d import OccamResponse
>>> response = OccamResponse.read("occam_run/RESP17.resp")
>>> response.rms

Summarize misfit by station and frequency index:

>>> response.misfit_per_site()
>>> response.misfit_per_frequency()

Select phase rows by Occam type code:

>>> phase_tm = response.data[response.data[:, 2] == 6]
>>> phase_tm.shape[0]

References

[OccamResponse-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

[OccamResponse-2]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

classmethod read(path, data_fn=None, **kwargs)#

Read an Occam2D response file.

The reader parses seven-column numeric rows from a response file produced by the Occam2D executable. The first three columns are stored in the raw table as floats because the file itself is numeric text, but convenience properties expose site, frequency, and type codes as integers.

Parameters:

path (path-like) – Path to the response file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path.
data_fn (path-like, optional) – Optional data-file path reserved for consistency checks between observed data rows and response rows. It is currently accepted for API stability but is not used by the parser.
**kwargs – Additional keyword arguments forwarded to the OccamResponse constructor. Use this for verbose or logger.

Returns:

Parsed response container with raw data, observed values, modeled values, residuals, and global RMS populated.

Return type:

OccamResponse

Raises:

FileNotFoundError – Raised when path does not exist.
ValueError – Raised when no valid seven-column numeric response rows can be parsed.

See also

OccamResponse.misfit_per_site: Computes station-index RMS values from residuals.
OccamResponse.misfit_per_frequency: Computes frequency-index RMS values.

Examples

>>> from pycsamt.models.occam2d import OccamResponse
>>> response = OccamResponse.read("RESP17.resp")
>>> response.n_data
>>> response.type_codes

property n_data: int#

property site_indices: ndarray#: 1-based site indices (int).

property freq_indices: ndarray#: 1-based frequency indices (int).

property type_codes: ndarray#: Unique data-type codes present in this response.

misfit_per_site()#

Return RMS misfit for each site index.

The returned values are computed from the weighted residual column:

\[\mathrm{RMS}_s = \sqrt{ \frac{1}{N_s} \sum_{i \in s} r_i^2 }.\]

Returns:: Mapping from one-based site index to RMS weighted residual. Empty response objects return an empty dictionary.
Return type:: dict of int to float

Examples

>>> from pycsamt.models.occam2d import OccamResponse
>>> response = OccamResponse.read("RESP17.resp")
>>> per_site = response.misfit_per_site()
>>> per_site[1]

misfit_per_frequency()#

Return RMS misfit for each frequency index.

The returned values group residuals by Occam’s one-based frequency index:

\[\mathrm{RMS}_f = \sqrt{ \frac{1}{N_f} \sum_{i \in f} r_i^2 }.\]

Returns:: Mapping from one-based frequency index to RMS weighted residual. Empty response objects return an empty dictionary.
Return type:: dict of int to float

Examples

>>> from pycsamt.models.occam2d import OccamResponse
>>> response = OccamResponse.read("RESP17.resp")
>>> per_freq = response.misfit_per_frequency()
>>> max(per_freq.values())

class pycsamt.models.occam2d.OccamLog(**kwargs)#

Bases: OccamBase

Represent an Occam2D convergence log.

OccamLog parses the text log written by the Occam2D Fortran executable. The file records one block per inversion iteration, including accepted misfit, model roughness, Lagrange multiplier, and line-search step size. Parsed arrays align by index, so iterations[i], rms[i], roughness[i], lagrange[i], and stepsize[i] describe the same iteration.

The main convergence statistic is the normalized RMS data misfit. If \(r_i\) are weighted residuals for \(N\) data, the reported quantity is commonly interpreted as

\[\phi_d = \sqrt{\frac{1}{N}\sum_{i=1}^N r_i^2} .\]

An inversion is usually considered well weighted when \(\phi_d \approx 1\). The practical target still depends on error estimates and modeling assumptions [OccamLog-1].

Parameters:

verbose (int or bool, default 0) – Verbosity level for progress reporting. 0 or False keeps the object quiet. Positive values enable progress messages through the instance logger; larger values may be used by callers to request more diagnostic detail.
logger (logging.Logger, optional) – Logger used for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created automatically. Pass an explicit logger when integrating Occam2D objects into an application-level logging setup.

Variables:

iterations (ndarray of int, shape (n_iter,)) – One-based iteration numbers parsed from ** ITERATION blocks. The values are preserved as written by Occam.
rms (ndarray of float, shape (n_iter,)) – Accepted normalized RMS misfit for each iteration. When an iteration contains repeated search steps, the parser keeps the last AND IS = value in that block.
roughness (ndarray of float, shape (n_iter,)) – Model roughness reported by Occam. The final entry may be nan when a run stops before writing ROUGHNESS IS.
lagrange (ndarray of float, shape (n_iter,)) – Accepted Lagrange multiplier, \(\mu\), for each iteration. Values are read from MINIMUM TOL FROM or INTERCEPT IS AT MU lines.
stepsize (ndarray of float, shape (n_iter,)) – Accepted step size for each iteration. The final entry may be nan if convergence problems stop the run early.

Notes

The parser is intentionally tolerant of Occam2D log variants. It ignores intermediate TOFMU search lines and keeps the last accepted values in each iteration block. This behavior matches logs where divergence problems trigger repeated Lagrange searches before a step is accepted.

See also

OccamRunner: Produces the log file by launching the executable.
InversionResult: Loads logs with model, iteration, and response files.
Plot2D.misfit: Visualizes RMS convergence from an OccamLog object.

Examples

Read a log and inspect the best iteration:

>>> from pycsamt.models.occam2d import OccamLog
>>> log = OccamLog.read("occam_run/LogFile.logfile")
>>> log.best_iteration

Check whether the inversion reached the common RMS target:

>>> log.converged

Print a compact report for scripts:

>>> log.summary()

References

[OccamLog-1]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

[OccamLog-2]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

classmethod read(path, **kwargs)#

Read an Occam2D log file.

The reader scans the file line by line with a small state machine. A ** ITERATION line starts a new block and saves the previous block. Within each block, the last accepted RMS, roughness, Lagrange multiplier, and step size are retained.

This is useful for logs where Occam cuts the step size or repeats the Lagrange search. Intermediate trial values are not stored because they do not describe the accepted model.

Parameters:

path (path-like) – Path to an Occam2D input or output file. The value may be a string, pathlib.Path, or any object accepted by pathlib.Path. Relative paths are interpreted from the current working directory of the Python process. Readers require the file to exist, while writers create parent directories when the owning method supports output.
**kwargs – Additional keyword arguments forwarded to the OccamLog constructor. Use this for verbose or logger when integrating the parser into a larger workflow.

Returns:

Parsed convergence-log container with one array entry per completed iteration block.

Return type:

OccamLog

Raises:

FileNotFoundError – Raised when path does not exist.

See also

OccamLog.summary: Returns a short text summary of convergence.
OccamLog.best_iteration: Reports the iteration with the lowest finite RMS misfit.

Examples

>>> from pycsamt.models.occam2d import OccamLog
>>> log = OccamLog.read("occam_run/LogFile.logfile")
>>> log.n_iter
>>> log.rms[-1]

property n_iter: int#: Number of parsed iterations.

property converged: bool#: True if any iteration achieved RMS ≤ 1.0 (Occam target).

property best_iteration: int#: 1-based iteration number with the lowest finite RMS misfit.

summary()#

Return a one-paragraph convergence summary.

Return type:: str

class pycsamt.models.occam2d.PlotModel(result=None, rho_min=1.0, rho_max=1000.0, depth_max=None, show_stations=True, profile_distance_unit='km', section='inversion', **kwargs)#

Bases: _OccamPlotBase

Plot a two-dimensional Occam resistivity model.

PlotModel displays the selected iteration model from InversionResult as a depth section. It replaces plotOccam2DMT.m, and the companion profile extractor replaces ExtractOccam2DMTProfile.m.

Occam iteration files store model parameters as \(\log_{10}\) resistivity. The plot converts them back to ohm metres before drawing:

\[\rho(x, z) = 10^{m(x, z)}.\]

The mesh is centered around the profile midpoint. Depth is positive downward.

Parameters:

result (InversionResult) – Loaded result containing rho_2d and mesh. Station markers are drawn when result.data.offsets is available.
rho_min (float, default 1.0, 1000.0) – Color-scale limits in ohm metres. They are passed to matplotlib.colors.LogNorm; both values must be positive.
rho_max (float, default 1.0, 1000.0) – Color-scale limits in ohm metres. They are passed to matplotlib.colors.LogNorm; both values must be positive.
depth_max (float, optional) – Maximum display depth in metres. If omitted, the full mesh depth is shown.
show_stations (bool, default True) – If True, overlay station triangles at the surface using the offsets stored in the data file.
profile_distance_unit ({"m", "km"}, default "km") – Unit used on the horizontal axis and by extract_profile().
figsize (tuple of float, optional) – Matplotlib figure size. The default suits a profile section.
cmap (str, default "jet_r") – Colormap used for the resistivity image.
dpi (int, default 100) – Figure resolution in dots per inch.
section (str | SectionStyle)

Returns:

Figure containing the resistivity model section.

Return type:

Raises:

RuntimeError – If the result does not contain rho_2d.

plot()#: Return the model section as a Matplotlib figure.

extract_profile(x0, x1)#

Extract centered coordinates, depth centers, and the log10-resistivity grid between x0 and x1.

Parameters:

x0 (float)
x1 (float)

Return type:

tuple

See also

InversionResult: Reconstructs rho_2d from Occam output files.
PlotSounding1D: Extracts model columns at station positions.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> from pycsamt.models.occam2d import PlotModel
>>> result = InversionResult("occam_run")
>>> fig = PlotModel(result, depth_max=2000).plot()
>>> plotter = PlotModel(result)
>>> x, z, log_rho = plotter.extract_profile(-1.0, 1.0)

References

[PlotModel-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

plot()#

Return a pcolormesh Figure of the 2-D resistivity model.

Return type:: matplotlib.figure.Figure

extract_profile(x0, x1)#

Extract (x_centers, z_centers, rho_subset) between x0 and x1.

Coordinates use profile_distance_unit. x0 and x1 are in the centered profile frame.

Returns:

(x_centers, z_centers, rho_2d_subset) where rho_2d_subset has shape (n_zcells, n_cols_in_range).

Return type:

tuple

Parameters:

x0 (float)
x1 (float)

class pycsamt.models.occam2d.PlotResponse(result=None, stations=None, modes=None, period_axis=True, max_stations=9, **kwargs)#

Bases: _OccamPlotBase

Plot observed and modeled Occam response curves.

PlotResponse compares observed data from the Occam file with modeled values from an Occam .resp file. It replaces plotOccam2DMTResponse.m.

Apparent-resistivity rows are stored as log10 values. They are converted before plotting:

\[\rho_a = 10^{d_\rho}.\]

Phase rows are plotted in degrees. The frequency index in the table is mapped back to physical frequency when corresponding OccamData object is available.

Parameters:

result (InversionResult) – Loaded result containing response and ideally data. The response must expose the seven-column table.
stations (list of str, list of int, or None, default None) – Stations to plot. Strings are matched against result.data.sites. Integers are one-based site indices when names are unavailable. If None, stations are sampled from the response table.
modes (list of str, optional) – Electromagnetic modes to draw. Supported values are "TE" and "TM".
period_axis (bool, default True) – Reserved for interface clarity. The implementation uses period when frequencies are available and indices otherwise.
max_stations (int, default 9) – Maximum number of station columns when stations is None.
figsize (tuple of float, optional) – Figure size. If omitted, width scales with station count.
cmap (str, default "jet_r") – Stored for a consistent interface. It is not used by this curve plot.
dpi (int, default 100) – Figure resolution in dots per inch.

Returns:

Figure with apparent-resistivity and phase panels.

Return type:

Raises:

RuntimeError – If response data are missing, modes are absent, or no stations can be selected.

See also

PlotResponseGrid: Compact version designed for many stations.
OccamResponse: Reader for the .resp file used by this plot.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> from pycsamt.models.occam2d import PlotResponse
>>> result = InversionResult("occam_run")
>>> fig = PlotResponse(result, modes=["TM"]).plot()

References

[PlotResponse-1]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

plot()#

Return a Figure with rho_a / phase subplots per station.

Return type:: matplotlib.figure.Figure

class pycsamt.models.occam2d.PlotPseudo(result=None, mode='TM', data_type='rho', **kwargs)#

Bases: _OccamPlotBase

Plot an Occam observed-data pseudosection.

PlotPseudo displays one data component from the Occam data file as a station-period view. It is the Python replacement for plotOccam2DMTPseudo.m.

The horizontal axis is station offset in kilometres when offsets are available. The vertical axis is \(\log_{10}(T)\), where \(T = 1/f\) is period in seconds. Apparent resistivity is converted from log10 storage to ohm metres; phase data remain in degrees.

Parameters:

result (InversionResult) – Loaded result containing an OccamData object with data blocks, offsets, and frequencies.
mode ({"TE", "TM"}, default "TM") – Electromagnetic mode to display. "TE" maps to codes 1 and 2; "TM" maps to codes 5 and 6.
data_type ({"rho", "phase"}, default "rho") – Quantity to display. "rho" selects apparent resistivity; "phase" selects phase.
figsize (tuple of float, optional) – Matplotlib figure size.
cmap (str, default "jet_r") – Colormap used for the pseudosection.
dpi (int, default 100) – Figure resolution in dots per inch.

Returns:

Pseudosection figure.

Return type:

Raises:

RuntimeError – If no data blocks are available or the selected type is not present.
ValueError – If mode and data_type are unsupported.

See also

OccamData: Provides the data block for the pseudosection.
PlotSiteMisfit: Builds a residual pseudosection from response values.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> from pycsamt.models.occam2d import PlotPseudo
>>> result = InversionResult("occam_run")
>>> fig = PlotPseudo(result, mode="TE").plot()

plot()#

Return a pseudosection pcolormesh Figure.

Return type:: matplotlib.figure.Figure

class pycsamt.models.occam2d.PlotMisfit(result=None, show_roughness=True, show_lagrange=False, target_line=True, **kwargs)#

Bases: _OccamPlotBase

Plot Occam2D convergence metrics by iteration.

PlotMisfit visualizes the convergence history stored in an OccamLog attached to an InversionResult. It replaces the MATLAB plotOccamIterMisfit.m view.

The main curve is the normalized root-mean-square data misfit:

\[\mathrm{RMS} = \sqrt{\frac{1}{N}\sum_{i=1}^{N} r_i^2},\]

Here \(r_i\) is the weighted residual for datum i. A run is commonly acceptable when the RMS approaches the target value of 1.0 [PlotMisfit-1].

Parameters:

result (InversionResult) – Loaded inversion result. It must expose log with iterations, rms, roughness, lagrange, and n_iter attributes.
show_roughness (bool, default True) – If True, add a secondary y-axis for roughness. Roughness is plotted on a log scale because it can vary by several orders of magnitude.
show_lagrange (bool, default False) – If True, add a lower panel for the accepted Lagrange multiplier at each iteration.
target_line (bool, default True) – If True, draw a dashed line at RMS equal to 1.0.
figsize (tuple of float, optional) – Matplotlib figure size passed through the shared base. If omitted, the number of panels controls the size.
cmap (str, default "jet_r") – Stored for consistency with other plot classes. It is not used by this line plot.
dpi (int, default 100) – Figure resolution in dots per inch.

Returns:

Figure containing the convergence plot.

Return type:

Raises:

RuntimeError – If result.log is missing or has no iterations.

See also

OccamLog: Parses convergence values from the Occam log file.
InversionResult.plot_misfit: Convenience wrapper that instantiates this class.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> from pycsamt.models.occam2d import PlotMisfit
>>> result = InversionResult("occam_run")
>>> fig = PlotMisfit(result, show_lagrange=True).plot()

References

[PlotMisfit-1]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

plot()#

Return a convergence Figure (RMS ± roughness vs iteration).

Return type:: matplotlib.figure.Figure

class pycsamt.models.occam2d.PlotSounding1D(result=None, stations=None, max_stations=16, depth_max=None, rho_min=1.0, rho_max=1000.0, overlay=False, **kwargs)#

Bases: _OccamPlotBase

Plot station-centered 1-D profiles from a 2-D Occam model.

PlotSounding1D samples the reconstructed 2-D resistivity grid at the mesh column nearest each station. The result is a set of resistivity-depth curves to compare vertical structure below stations.

The plotted resistivity is converted from the log10 grid:

\[\rho(z) = 10^{m(z)}.\]

Air layers are omitted using result.mesh.n_airlayers.

Parameters:

result (InversionResult) – Loaded result containing rho_2d, mesh, and data.offsets.
stations (list of str or None, default None) – Station names to plot. If None, stations are sampled from all available offsets.
max_stations (int, default 16) – Maximum station profiles when stations is None.
depth_max (float, optional) – Maximum plotted depth. If omitted, mesh depth controls the lower limit.
rho_min (float, default 1.0, 1000.0) – Horizontal resistivity-axis limits in ohm metres.
rho_max (float, default 1.0, 1000.0) – Horizontal resistivity-axis limits in ohm metres.
overlay (bool, default False) – If True, draw selected profiles on one axis. If False, draw one panel per station.
figsize (tuple of float, optional) – Figure size. Defaults depend on station count.
cmap (str, default "jet_r") – Stored for interface consistency. Overlay plots use a tabular Matplotlib colormap internally.
dpi (int, default 100) – Figure resolution in dots per inch.

Returns:

Resistivity-depth profile figure.

Return type:

Raises:

RuntimeError – If rho_2d is missing, station offsets are missing, or station selection is empty.

See also

PlotModel.extract_profile: Extracts a horizontal interval from the same grid.
InversionResult: Provides the reconstructed model grid.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> from pycsamt.models.occam2d import PlotSounding1D
>>> result = InversionResult("occam_run")
>>> fig = PlotSounding1D(result, overlay=True).plot()

plot()#

Return a Figure of 1-D ρ–depth soundings.

Return type:: matplotlib.figure.Figure

class pycsamt.models.occam2d.PlotSiteMisfit(result=None, modes=None, show_residual_map=True, rms_target=1.0, **kwargs)#

Bases: _OccamPlotBase

Plot per-site Occam response misfit diagnostics.

PlotSiteMisfit summarizes the fit between observed and modeled values at each station. The top panel is a bar chart of RMS residual by station and data type. The optional lower panel is a residual pseudosection.

Residuals are normalized by the error column in the Occam response table:

\[r_i = \frac{d_i^{obs} - d_i^{pred}}{\sigma_i}.\]

Per-site RMS values use these normalized residuals. If a response error is non-positive, that residual is ignored.

Parameters:

result (InversionResult) – Loaded result containing response. Station labels and frequencies come from result.data.
modes (list of str, optional) – Modes included in the summary. Supported values are "TE" and "TM". If omitted, both are used.
show_residual_map (bool, default True) – If True, draw the normalized residual map under the bar chart.
rms_target (float, default 1.0) – Target RMS value drawn in the bar panel. Use None to omit the target line.
figsize (tuple of float, optional) – Figure size. Defaults scale with station count.
cmap (str, default "jet_r") – Stored for the shared interface. The residual map uses a diverging colormap.
dpi (int, default 100) – Figure resolution in dots per inch.

Returns:

Figure containing per-site RMS diagnostics.

Return type:

Raises:

RuntimeError – If response data are missing or requested type codes are absent.

See also

OccamResponse.misfit_per_site: Returns a simpler per-site RMS dictionary.
PlotResponseGrid: Shows observed and modeled curves for many stations.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> from pycsamt.models.occam2d import PlotSiteMisfit
>>> result = InversionResult("occam_run")
>>> fig = PlotSiteMisfit(result).plot()

plot()#

Return a per-site misfit Figure.

Return type:: matplotlib.figure.Figure

class pycsamt.models.occam2d.PlotResponseGrid(result=None, stations=None, n_cols=5, modes=None, max_stations=25, **kwargs)#

Bases: _OccamPlotBase

Plot a compact grid of observed and modeled responses.

PlotResponseGrid is designed to scan many stations. Each station uses two axes: apparent resistivity above and phase below. Observed values are drawn as points and modeled values as lines. Titles include per-site RMS when response errors are available.

Apparent-resistivity values are converted from log10 Occam storage before plotting:

\[\rho_a = 10^{d_\rho}.\]

Parameters:

result (InversionResult) – Loaded result with response and ideally data.
stations (list of str or None, default None) – Station names to include. If omitted, station indices sampled from the response table.
n_cols (int, default 5) – Maximum number of station columns in each grid row.
modes (list of str, optional) – Modes to draw. Use "TE", "TM", or both.
max_stations (int, default 25) – Maximum station count included when stations is None.
figsize (tuple of float, optional) – Figure size. Defaults scale with columns and rows.
cmap (str, default "jet_r") – Stored for interface consistency. It is unused by this curve plot.
dpi (int, default 100) – Figure resolution in dots per inch.

Returns:

Compact response-grid figure.

Return type:

Raises:

RuntimeError – If response data or station choices are missing.

See also

PlotResponse: Larger response panels for a smaller station subset.
PlotSiteMisfit: Per-site residual summary from the response table.

Examples

>>> from pycsamt.models.occam2d import InversionResult
>>> from pycsamt.models.occam2d import PlotResponseGrid
>>> result = InversionResult("occam_run")
>>> fig = PlotResponseGrid(result, n_cols=4).plot()

plot()#

Return a compact response-grid Figure.

Return type:: matplotlib.figure.Figure

class pycsamt.models.occam2d.OccamConfig(modes=<factory>, error_floor_rho=0.05, error_floor_phase=0.5, freq_min=None, freq_max=None, n_layers=30, max_depth=1500.0, n_airlayers=5, cell_size_horizontal=100.0, cell_size_vertical_top=10.0, depth_scale=1.2, n_padding_x=7, max_iterations=100, target_misfit=1.0, roughness_type=1, diagonal_penalties=0, stepsize_cut_count=8, debug_level=1, initial_rho=100.0, lagrange_start=5.0, data_file='OccamDataFile.dat', mesh_file='Occam2DMesh', model_file='Occam2DModel', startup_file='Startup', binary_name='Occam2D')#

Bases: object

Collect settings that define an Occam2D run.

OccamConfig groups the options shared by the Occam2D builder, file containers, runner, and result loaders. It is a plain dataclass, so users may set fields at construction time or mutate them before calling InputBuilder.

The configuration controls four parts of the workflow:

data selection and error floors;
mesh geometry and depth discretization;
startup and inversion-control values;
file names and executable discovery.

2.22. Data Options#

modeslist of str, default [“TE”, “TM”]: Electromagnetic modes written to the Occam data file. "TE" selects the \(Z_{xy}\) component and "TM" selects \(Z_{yx}\). Each selected mode writes apparent-resistivity and phase rows.
error_floor_rhofloat: Relative apparent-resistivity error floor. A value of 0.05 means five percent. Because Occam stores apparent resistivity as \(\log_{10}(\rho_a)\), the builder converts this floor before writing data.
error_floor_phasefloat: Absolute phase error floor in degrees. This keeps phase rows with unrealistically small source errors from dominating the normalized data misfit.
freq_minfloat or None: Lower frequency limit in hertz. Frequencies below this value are excluded when built from EDI sources. None leaves the lower bound open.
freq_maxfloat or None: Upper frequency limit in hertz. Frequencies above this value are excluded when built from EDI sources. None leaves the upper bound open.

2.22. Mesh Options#

n_layersint: Maximum number of active earth layers below the air layers. Larger values represent more vertical structure but increase the parameter count. The mesh builder may stop earlier, once max_depth is reached.
max_depthfloat: Target maximum depth in metres for the earth-layer column. OccamMesh.from_data() stops adding geometrically expanding layers – truncating the last one if needed – once cumulative depth reaches this value, or once n_layers layers have been added, whichever comes first.
n_airlayersint: Number of air layers above the earth model. These layers stabilize finite-element boundaries near the surface.
cell_size_horizontalfloat: Target horizontal cell width in metres near stations. Smaller values give finer lateral detail and bigger meshes.
cell_size_vertical_topfloat: Thickness, in metres, of the top earth layer and air layers used by the current mesh builder.
depth_scalefloat: Geometric multiplier applied to layer thickness with depth. Values greater than 1 make layers progressively thicker.
n_padding_xint: Number of horizontal padding cells added on each side. Padding moves side boundaries away from the survey profile.

2.22. Startup Options#

max_iterationsint: Maximum number of Occam iterations requested in the startup file.
target_misfitfloat: Target normalized RMS misfit. Values near 1.0 are typical when data errors are realistic.
roughness_typeint: Roughness penalty type passed to Occam. A value of 1 selects the standard gradient penalty; 2 selects curvature when supported by the executable.
diagonal_penaltiesint: Flag controlling diagonal roughness penalties in the startup file. 0 disables them.
stepsize_cut_countint: Maximum number of Lagrange step-size reductions allowed during a line-search stage.
debug_levelint: Debug verbosity passed to the Occam executable.
initial_rhofloat: Starting half-space resistivity in ohm metres. The startup vector is initialized as \(\log_{10}\) of this value.
lagrange_startfloat: Initial Lagrange multiplier written to startup.

2.22. File Options#

data_filestr: Data filename written inside the run directory.
mesh_filestr: Mesh filename written inside the run directory.
model_filestr: Model filename written inside the run directory.
startup_filestr: Startup filename passed to the Occam executable.
binary_namestr: Executable name searched by OccamRunner.

Notes

InputBuilder.build accepts one-shot overrides for common data and mesh fields. Overrides update the same OccamConfig instance stored on the builder.

2.22. Source-Of-Truth Files#

Users can generate an editable configuration file before building an Occam2D run. Python is the default template format because it supports rich inline comments and can be read safely by from_file() using literal parsing. YAML templates also keep comments. JSON templates store explanations in a "_schema" metadata block and editable values under "config" because standard JSON has no comment syntax.

The recommended workflow is:

Generate a template with write_template().
Edit values in the generated file.
Load the edited file with from_file() or read().
Pass the resulting configuration to builders and runners.

See also

InputBuilder: Consumes this configuration while writing input files.
OccamData.from_edi: Uses data options to select modes, frequencies, and errors.
OccamMesh.from_data: Uses mesh options to build finite-element geometry.
OccamStartup.from_model: Uses startup options to initialize inversion controls.
OccamRunner: Uses file and binary settings during execution.

Examples

Create a standard configuration for the builder:

>>> from pycsamt.models.occam2d import OccamConfig
>>> from pycsamt.models.occam2d import InputBuilder
>>> cfg = OccamConfig(n_layers=32, target_misfit=1.0)
>>> cfg.error_floor_rho = 0.07
>>> builder = InputBuilder([], workdir="run", config=cfg)

Restrict the frequency band and use only TM mode:

>>> cfg = OccamConfig(modes=["TM"])
>>> cfg.freq_min = 0.1
>>> cfg.freq_max = 1000.0

Configure a finer near-station mesh:

>>> cfg = OccamConfig()
>>> cfg.cell_size_horizontal = 50.0
>>> cfg.cell_size_vertical_top = 5.0
>>> cfg.depth_scale = 1.15
>>> cfg.max_depth = 800.0  # shallow near-surface target, in metres

Generate a documented source-of-truth template:

>>> path = OccamConfig.write_template("occam2d_config.py")
>>> cfg = OccamConfig.from_file(path)
>>> cfg.binary_name
'Occam2D'

Use YAML when the configuration will be edited outside Python:

>>> OccamConfig.write_template("occam2d_config.yml")
PosixPath('occam2d_config.yml')

References

[OccamConfig-1]

Constable, S. C., Parker, R. L., and Constable, C. G., “Occam’s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data”, Geophysics, 52(3), 289-300, 1987.

[OccamConfig-2]

deGroot-Hedlin, C., and Constable, S., “Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data”, Geophysics, 55(12), 1613-1624, 1990.

modes: list[str]#

error_floor_rho: float = 0.05#

error_floor_phase: float = 0.5#

freq_min: float | None = None#

freq_max: float | None = None#

n_layers: int = 30#

max_depth: float = 1500.0#

n_airlayers: int = 5#

cell_size_horizontal: float = 100.0#

cell_size_vertical_top: float = 10.0#

depth_scale: float = 1.2#

n_padding_x: int = 7#

max_iterations: int = 100#

target_misfit: float = 1.0#

roughness_type: int = 1#

diagonal_penalties: int = 0#

stepsize_cut_count: int = 8#

debug_level: int = 1#

initial_rho: float = 100.0#

lagrange_start: float = 5.0#

data_file: str = 'OccamDataFile.dat'#

mesh_file: str = 'Occam2DMesh'#

model_file: str = 'Occam2DModel'#

startup_file: str = 'Startup'#

binary_name: str = 'Occam2D'#

to_template(path='occam2d_config.py', *, fmt=None)#

Write this configuration as an editable template.

Parameters:

path (path-like, default "occam2d_config.py") – Destination file. If the path has no suffix, the suffix is inferred from fmt and defaults to .py.
fmt ({"py", "json", "yml", "yaml"}, optional) – Template format. Python and YAML templates include comments. JSON templates include a "_schema" metadata block because standard JSON does not support comments.

Returns:

Path of the generated source-of-truth file.

Return type:

classmethod write_template(path='occam2d_config.py', *, fmt=None)#

Write a default editable Occam2D configuration file.

Parameters:

path (path-like, default "occam2d_config.py") – Destination file. Suffixes .py, .json, .yml, and .yaml select the output format.
fmt ({"py", "json", "yml", "yaml"}, optional) – Explicit output format. When omitted, the suffix of path is used; paths without a suffix produce a Python template.

Returns:

Path of the generated template.

Return type:

classmethod from_file(path, *, strict=True)#

Create a configuration from a source-of-truth file.

Parameters:

path (path-like) – Python, JSON, YML, or YAML configuration file generated by write_template() or following the same structure.
strict (bool, default True) – If True, unknown editable keys raise ValueError. If False, unknown keys are ignored. Metadata keys beginning with "_" are always ignored.

Returns:

Configuration populated from edited file values.

Return type:

OccamConfig

classmethod read(path, *, strict=True)#

Create a configuration from a source-of-truth file.

Parameters:

path (path-like) – Python, JSON, YML, or YAML configuration file generated by write_template() or following the same structure.
strict (bool, default True) – If True, unknown editable keys raise ValueError. If False, unknown keys are ignored. Metadata keys beginning with "_" are always ignored.

Returns:

Configuration populated from edited file values.

Return type:

OccamConfig

Parameters:

modes (list[str])
error_floor_rho (float)
error_floor_phase (float)
freq_min (float | None)
freq_max (float | None)
n_layers (int)
max_depth (float)
n_airlayers (int)
cell_size_horizontal (float)
cell_size_vertical_top (float)
depth_scale (float)
n_padding_x (int)
max_iterations (int)
target_misfit (float)
roughness_type (int)
diagonal_penalties (int)
stepsize_cut_count (int)
debug_level (int)
initial_rho (float)
lagrange_start (float)
data_file (str)
mesh_file (str)
model_file (str)
startup_file (str)
binary_name (str)

class pycsamt.models.occam2d.SyntheticSite(name, x_m, freq, rho, phase, rho_err, phase_err)#

Bases: object

Minimal object accepted by OccamData.from_edi (via InputBuilder).

Not a real pycsamt.site.Site – there is no EDI file behind forward-modelled array data to justify building one. coords fakes a longitude so OccamData()’s lat/lon-to-metres offset recovery lands on the real along-profile chainage that was passed in.

Parameters:

name (str) – Station label.
x_m (float) – Along-profile position, metres.
freq (array_like, shape (n_freq,)) – Frequencies, Hz.
rho (ndarray, shape (n_freq, 2, 2)) – Apparent resistivity, Ω·m. Only [:, 0, 1] (TE / Zxy) and [:, 1, 0] (TM / Zyx) are read.
phase (ndarray, shape (n_freq, 2, 2)) – Phase, degrees, raw (not pre-shifted – OccamData.from_edi applies the conventional TM +180 degree normalisation into the first quadrant itself).
rho_err (ndarray, shape (n_freq, 2, 2)) – Per-datum errors, same units as rho/phase.
phase_err (ndarray, shape (n_freq, 2, 2)) – Per-datum errors, same units as rho/phase.

pycsamt.models.occam2d.sites_from_response(resp, station_x, station_names, *, rho_err_frac=0.05, phase_err_deg=1.5)#

Build one SyntheticSite per station from a real ForwardResponse2D.

Parameters:

resp (ForwardResponse2D) – Real forward-modelled response, e.g. from MT2DForward(freqs, grid).run(). resp.rho_a_te/rho_a_tm and resp.phase_te/phase_tm must have shape (n_freq, n_stations), matching station_x’s length.
station_x (sequence of float) – Along-profile position of each station, metres, same order and length as resp’s station axis.
station_names (sequence of str) – Station labels, same length as station_x.
rho_err_frac (float, default 0.05) – Assigned relative resistivity error (5%).
phase_err_deg (float, default 1.5) – Assigned absolute phase error, degrees.

Return type:

list of SyntheticSite

Examples

>>> import numpy as np
>>> from pycsamt.forward.em2d import MT2DForward
>>> from pycsamt.forward.grid2d import Grid2D
>>> from pycsamt.models.occam2d.synthetic import sites_from_response
>>> grid = Grid2D.halfspace(rho=100.0, nx=20, x_max=2000.0, n_stations=5)
>>> resp = MT2DForward(np.array([100.0, 10.0]), grid, verbose=False).run()
>>> sites = sites_from_response(resp, grid.x_stations, ["S0", "S1", "S2", "S3", "S4"])
>>> len(sites), sites[0].name
(5, 'S0')
>>> sites[0].rho.shape
(2, 2, 2)

2.22.2.2. pycsamt.models.mare2dem#

pycsamt.models.mare2dem — Python interface to the MARE2DEM EM inversion code.

MARE2DEM (Modeling with Adaptively Refined Elements for 2.5-D Electromagnetic inversion) implements finite-element MT and CSEM forward modelling and Occam-style regularized inversion. It is an MPI-parallel Fortran code developed by Kerry Key at LDEO, Columbia University.

Because the compiled binary is ~49 MB, the source tree is not bundled in the repository. Use SourceManager to download and build it before running inversions.

Ported MATLAB scripts (complete)#

I/O (a_util/io/ + a_util/data/):

.emdata / .EMResp reader + writer
.resistivity reader + writer (all anisotropy modes)
Triangle .poly PSLG reader + writer
.settings parallel-decomposition writer
.emdata_group data-group file reader + writer
Group-RMS CSV log reader
Most-recently-modified file finder

Data management:

Data-type code lookup table
High-level .emdata builder from MT / CSEM survey configs
ZMM impedance file reader + MT data-file builder
Synthetic noise addition + make_synthetic_data wrapper
Multi-file merge utility

Geometry (a_util/geom/, a_util/mapping/):

Topography interpolation + slope angles
Line-segment intersections (bounding-box pre-filtered)
Douglas-Peucker polyline simplification
PSLG polygon simplification (collinear node removal)
Area-weighted triangle region centroids
Survey-profile line orientation
Station-to-profile projection
UTM ↔ Lon/Lat conversion (pyproj with pure-Python WGS-84 fallback)
Survey area-of-interest estimator
Triangle FEM region flood-fill assignment

Model construction:

2-D resistivity grid → MARE2DEM .poly + .resistivity
Topography import + profile projection

Mesh generation:

Topography-aware PSLG construction + Triangle refinement (build_survey_mesh, run_triangle)
Triangle .node/.ele reader (read_triangulation) and conversion to the solver-neutral TriMesh contract (tri_mesh_from_poly)

Model comparison:

Log10 (or custom) difference of two .resistivity files

Plotting:

RMS convergence curve
Survey map (Rx/Tx positions in UTM)
Receiver geometry QC (6-panel)
Transmitter geometry QC
.poly PSLG mesh plot
Resistivity section as a color-filled triangular mesh (via read_triangulation), falling back to a histogram when no mesh is present

Quick start#

Step 1 — Download and compile (once per machine):

from pycsamt.models.mare2dem import SourceManager

sm = SourceManager(verbose=1)
sm.download()  # git-clone from Bitbucket
sm.build()  # requires Intel oneAPI + MKL

Step 2 — Create a data file from MT survey parameters:

import numpy as np
from pycsamt.models.mare2dem import MTSurveyConfig, make_data_file

mt = MTSurveyConfig(
    frequencies=np.logspace(-3, 3, 20),
    rx_y=np.linspace(-5000, 5000, 20),
    rx_type="marine",
    lTE=True,
    lTM=True,
)
em = make_data_file("survey.emdata", topo=-1000.0, mt=mt)

Step 3 — Prepare a full input set and run:

from pycsamt.models.mare2dem import (
    Mare2DEMConfig,
    InputBuilder,
    Mare2DEMRunner,
)

cfg = Mare2DEMConfig(initial_rho=1.0, n_procs=8)
InputBuilder(cfg).build("survey.emdata", workdir="./run")
result = Mare2DEMRunner("./run", cfg).run("mare2dem")
result.print_summary()

References

Key, K. (2016). MARE2DEM: A 2-D inversion code for controlled-source electromagnetic and magnetotelluric data. Geophysical Journal International, 207(1), 571–588. doi:10.1093/gji/ggw290.

class pycsamt.models.mare2dem.Mare2DEMConfig(source_dir=None, fc_compiler=None, cc_compiler=None, binary='MARE2DEM', use_mpi=True, n_procs=4, mpi_command='mpirun', max_iterations=150, target_rms=1.0, initial_rho=1.0, data_file='mare2dem.emdata', resistivity_file='mare2dem.resistivity', settings_file='mare2dem.settings')#

Bases: object

Collect settings that define a MARE2DEM run.

Mare2DEMConfig is the central configuration object for the MARE2DEM wrapper. It is a plain dataclass whose fields cover source management, MPI execution, inversion control, starting model, and default file names.

The configuration is shared by SourceManager, InputBuilder, Mare2DEMRunner, and InversionResult so that all components of a workflow use consistent parameters.

2.22. Source Management#

source_dirpath-like or None, default None: Explicit path to the MARE2DEM Fortran source tree. When None, SourceManager applies a four-level fallback: the PYCSAMT_MARE2DEM_SOURCE environment variable, the bundled _source/ directory inside the package (writable dev installs), and finally the platform user-data directory (PyPI installs). Set this field when the source tree lives at an unconventional location.
fc_compilerstr or None, default None: MPI-Fortran compiler used to build MARE2DEM from source. When None, SourceManager.build() auto-detects in the order mpiifort (Intel oneAPI), mpifort (generic OpenMPI/MPICH). Intel compilers are preferred because the MARE2DEM Makefile uses Intel-specific pre-processor flags and the Intel MKL is required.
cc_compilerstr or None, default None: MPI-C compiler used for the Triangle mesh library and ScaLAPACK. Auto-detects mpiicc then mpicc when None.

2.22. Binary And MPI#

binarystr, default “MARE2DEM”: Name or absolute path of the compiled MARE2DEM executable. The runner first checks whether the name is on PATH, then looks inside the resolved source directory, and finally checks the user-data directory. Use an absolute path to pin a specific build.
use_mpibool, default True: Whether to launch MARE2DEM through an MPI wrapper. MARE2DEM is an MPI-parallel code; serial execution is not supported by the official distribution. Set to False only when testing with a special single-process build.
n_procsint, default 4: Number of MPI processes requested when use_mpi is True. Each process handles a subset of the spatial wavenumbers. Typical values range from 4 to the number of physical cores available on the node.
mpi_commandstr, default “mpirun”: MPI launcher command. Common alternatives are "mpiexec" and "srun" (SLURM). The value is prepended to the MARE2DEM command line when use_mpi is True.

2.22. Inversion Control#

max_iterationsint, default 150: Maximum number of Occam inversion iterations MARE2DEM is allowed to perform before stopping regardless of the misfit target.
target_rmsfloat, default 1.0: Normalized RMS misfit target. MARE2DEM stops iterating when the data misfit falls below this value. Lower values demand a tighter fit and may produce more structured models.

2.22. Initial Model#

initial_rhofloat, default 1.0: Starting half-space resistivity in ohm-metres for the initial homogeneous resistivity model. MARE2DEM internally uses log10-resistivity, so the value must be positive.

2.22. File Names#

data_filestr, default “mare2dem.emdata”: Default MARE2DEM data filename. The file uses the .emdata extension and records source–receiver geometry, observed data, and data uncertainties for one or more EM methods (MT, CSEM, or combined).
resistivity_filestr, default “mare2dem.resistivity”: Default MARE2DEM resistivity model filename. The file describes the 2-D finite-element mesh and layer resistivities and is the primary input file stem passed to the binary.
settings_filestr, default “mare2dem.settings”: Default MARE2DEM inversion-settings filename. It controls inversion type (Occam vs NLCG), iteration limits, target misfit, output verbosity, and optional regularization parameters.

ivar resistivity_stem:: The stem of resistivity_file without its extension. MARE2DEM receives this stem as its only positional argument and derives the data and settings filenames from it.
vartype resistivity_stem:: str

Notes

MARE2DEM implements 2.5-D finite-element MT and CSEM forward modelling in the frequency domain with an Occam-style regularized inversion [Mare2DEMConfig-1]. The inversion minimizes an objective of the form

\[\begin{split}\Phi(m) = \| W_d (F(m) - d) \|_2^2 + \\lambda \| \\nabla m \|_2^2,\end{split}\]

where \(F(m)\) is the 2.5-D FEM forward operator, \(d\) is the observed data vector, and \(W_d\) is the data-weighting operator. The configuration fields target_rms and max_iterations bound the inversion iteration.

2.22. Source-Of-Truth Files#

The recommended workflow is:

Generate a template with write_template().
Edit the values in the generated file.
Load the edited file with from_file().
Pass the configuration to SourceManager, InputBuilder, Mare2DEMRunner.

See also

SourceManager: Download and compile the MARE2DEM Fortran source.
InputBuilder: Write MARE2DEM resistivity model, data, and settings files.
Mare2DEMRunner: Launch the MARE2DEM binary subprocess.
InversionResult: Load MARE2DEM inversion output files.

Examples

Create a default configuration:

>>> from pycsamt.models.mare2dem.config import Mare2DEMConfig
>>> cfg = Mare2DEMConfig()
>>> cfg.resistivity_stem
'mare2dem'

Configure a parallel run with 8 MPI processes:

>>> cfg = Mare2DEMConfig(use_mpi=True, n_procs=8)

Point to a custom source directory:

>>> cfg = Mare2DEMConfig(source_dir="/opt/mare2dem_source")

References

[Mare2DEMConfig-1]

Key, K. (2016). MARE2DEM: A 2-D inversion code for controlled-source electromagnetic and magnetotelluric data. Geophysical Journal International, 207(1), 571-588. doi:10.1093/gji/ggw290.

source_dir: str | Path | None = None#

fc_compiler: str | None = None#

cc_compiler: str | None = None#

binary: str = 'MARE2DEM'#

use_mpi: bool = True#

n_procs: int = 4#

mpi_command: str = 'mpirun'#

max_iterations: int = 150#

target_rms: float = 1.0#

initial_rho: float = 1.0#

data_file: str = 'mare2dem.emdata'#

resistivity_file: str = 'mare2dem.resistivity'#

settings_file: str = 'mare2dem.settings'#

property resistivity_stem: str#: Return the stem of resistivity_file without extension.

to_template(path='mare2dem_config.py', *, fmt=None)#

Write this configuration as an editable template file.

Parameters:

path (path-like, default "mare2dem_config.py") – Destination file.
fmt ({"py", "json", "yml", "yaml"}, optional) – Output format. Defaults to .py when the path has no recognized suffix.

Returns:

Path of the generated template.

Return type:

classmethod write_template(path='mare2dem_config.py', *, fmt=None)#

Write a default editable MARE2DEM configuration file.

Parameters:

path (path-like, default "mare2dem_config.py") – Destination file. Suffixes .py, .json, .yml, and .yaml select the output format.
fmt ({"py", "json", "yml", "yaml"}, optional) – Explicit output format.

Returns:

Path of the generated source-of-truth file.

Return type:

Examples

>>> from pycsamt.models.mare2dem.config import Mare2DEMConfig
>>> path = Mare2DEMConfig.write_template("mare2dem_config.py")
>>> path.name
'mare2dem_config.py'

classmethod from_file(path, *, strict=True)#

Create a configuration from a source-of-truth file.

Parameters:

path (path-like) – Python, JSON, YML, or YAML configuration file generated by write_template().
strict (bool, default True) – If True, unknown keys raise ValueError.

Returns:

Configuration populated from the edited file.

Return type:

Mare2DEMConfig

Examples

>>> from pycsamt.models.mare2dem.config import Mare2DEMConfig
>>> Mare2DEMConfig.write_template("mare2dem_config.json")
PosixPath('mare2dem_config.json')
>>> cfg = Mare2DEMConfig.from_file("mare2dem_config.json")
>>> cfg.binary
'MARE2DEM'

classmethod read(path, *, strict=True)#

Create a configuration from a source-of-truth file.

Parameters:

path (path-like) – Python, JSON, YML, or YAML configuration file generated by write_template().
strict (bool, default True) – If True, unknown keys raise ValueError.

Returns:

Configuration populated from the edited file.

Return type:

Mare2DEMConfig

Examples

>>> from pycsamt.models.mare2dem.config import Mare2DEMConfig
>>> Mare2DEMConfig.write_template("mare2dem_config.json")
PosixPath('mare2dem_config.json')
>>> cfg = Mare2DEMConfig.from_file("mare2dem_config.json")
>>> cfg.binary
'MARE2DEM'

Parameters:

source_dir (str | Path | None)
fc_compiler (str | None)
cc_compiler (str | None)
binary (str)
use_mpi (bool)
n_procs (int)
mpi_command (str)
max_iterations (int)
target_rms (float)
initial_rho (float)
data_file (str)
resistivity_file (str)
settings_file (str)

class pycsamt.models.mare2dem.SourceManager(config=None, source_dir=None, **kwargs)#

Bases: Mare2DEMBase

Manage the MARE2DEM Fortran source: download, build, and locate.

SourceManager is the entry point for obtaining a working MARE2DEM binary. It separates source management from inversion execution so that users who already have a compiled binary can skip directly to Mare2DEMRunner.

2.22. Source-directory resolution#

The directory where sources are stored follows this priority:

source_dir constructor argument.
config.source_dir field.
PYCSAMT_MARE2DEM_SOURCE environment variable.
Bundled _source/ inside the installed package (writable dev installs only).
Platform user-data directory — the safe PyPI fallback:
- Linux / WSL: ~/.local/share/pycsamt/mare2dem/
- macOS: ~/Library/Application Support/pycsamt/mare2dem/
- Windows: %LOCALAPPDATA%\\pycsamt\\mare2dem\\

param config:: Configuration object. Supplies source_dir, fc_compiler, cc_compiler, and binary.
type config:: Mare2DEMConfig, optional
param source_dir:: Explicit path that overrides config.source_dir and the environment variable.
type source_dir:: path-like, optional
param verbose:: Verbosity level. 0 or False keeps the object quiet. Positive values enable diagnostic messages through the instance logger. Larger values may be used to request more detailed run, parsing, or build information.
type verbose:: int or bool, default 0
param logger:: Logger for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created automatically. Provide a logger when integrating MARE2DEM workflows into an application-wide logging configuration.
type logger:: logging.Logger, optional

Notes

MARE2DEM requires Intel compilers (mpiifx/mpiicx on current oneAPI releases, or the classic mpiifort/mpiicc on older ones) and the Intel MKL. When Intel oneAPI is installed, source the setvars.sh script before calling build():

source /opt/intel/oneapi/setvars.sh
python -c "
from pycsamt.models.mare2dem import SourceManager
sm = SourceManager(verbose=1)
sm.download()
sm.build()
"

Windows native builds are not supported. Use WSL2.

Examples

Quick-start — download and build:

>>> from pycsamt.models.mare2dem import SourceManager
>>> sm = SourceManager(verbose=1)
>>> sm.download()          # clones from Bitbucket
>>> sm.build()             # compiles with auto-detected compilers

Check status without downloading:

>>> sm.print_status()

Point to sources already on disk:

>>> sm = SourceManager(source_dir="/data/mare2dem_source")
>>> sm.build(inc_file="/data/mare2dem_source/include/habanero.inc")

See also

Mare2DEMRunner: Launch the compiled MARE2DEM executable for inversion.
Mare2DEMConfig: Configuration controlling compiler selection and binary name.

References

[SourceManager-1]

Key, K. (2016). MARE2DEM: A 2-D inversion code for controlled-source electromagnetic and magnetotelluric data. Geophysical Journal International, 207(1), 571-588. doi:10.1093/gji/ggw290.

REPO_URL: str = 'https://bitbucket.org/mare2dem/mare2dem_source'#

ARCHIVE_URL: str = 'https://bitbucket.org/mare2dem/mare2dem_source/get/master.tar.gz'#

resolve_source_dir()#

Return the directory where MARE2DEM sources should live.

The resolution order is:

Explicit source_dir argument passed to the constructor.
config.source_dir field.
PYCSAMT_MARE2DEM_SOURCE environment variable.
Bundled _source/ inside the package (only when writable — i.e. editable / development installs).
Platform user-data directory (always writable).

Returns:: Resolved directory (created if it does not yet exist).
Return type:: pathlib.Path

resolve_binary()#

Return the path to the compiled MARE2DEM binary or None.

Resolution order:

Binary name found on PATH.
<source_dir>/MARE2DEM.
Platform user-data directory MARE2DEM.

Return type:: Path | None

resolve_triangle_binary(name=None)#

Return the path to a Triangle mesh-generator executable or None.

Triangle meshing is a lighter-weight, independent concern from the MARE2DEM solver itself (no MPI/Fortran/LAPACK build required), so this resolves separately from resolve_binary() rather than assuming a full MARE2DEM build has happened.

Parameters:

name (str, optional) – Explicit executable name to look for. Tries both "triangle" and "Triangle" when omitted.

Returns:

pathlib.Path or None – Resolved executable, or None if not found anywhere.
Resolution order (per candidate name)
1. PATH lookup via shutil.which().
1. <source_dir>/<name> – Triangle is typically a byproduct of – build(), since MARE2DEM’s own Makefile compiles it.
1. Platform user-data directory <name>.
Each directory candidate is checked with ``shutil.which(name,
path=directory)`` before a literal-file fallback, so a bare name
like "triangle" also resolves to "triangle.exe" via
PATHEXT on Windows – the same fix already required for
pycsamt.forward.maxwell.external.resolve_executable()’s
equivalent search-path loop.

Return type:

Path | None

is_downloaded()#

Return True when the source tree appears populated.

Return type:: bool

is_built()#

Return True when the MARE2DEM binary exists.

Return type:: bool

download(*, method='auto', force=False)#

Download the MARE2DEM source tree.

Parameters:

method ({"auto", "git", "archive"}, default "auto") – Download strategy. "auto" tries "git" first and falls back to "archive".
force (bool, default False) – Re-download even if the source tree is already present.

Returns:

Source directory containing the downloaded tree.

Return type:

Raises:

RuntimeError – When neither git nor requests is available, or when the download fails.

build(*, inc_file=None, fc=None, cc=None, clean_first=False)#

Compile MARE2DEM from the downloaded source.

Parameters:

inc_file (path-like or None, default None) – Explicit Make include file. When None, one is auto-generated from the detected compilers and MKL.
fc (str or None, default None) – MPI-Fortran compiler override. Falls back to config.fc_compiler then auto-detection.
cc (str or None, default None) – MPI-C compiler override. Falls back to config.cc_compiler then auto-detection.
clean_first (bool, default False) – Run make clean_all before compiling to start fresh.

Returns:

Path to the compiled MARE2DEM binary.

Return type:

Raises:

FileNotFoundError – When the source tree is not present. Call download() first.
RuntimeError – When the build fails or the binary is not found after the build.

status()#

Return a summary dictionary of source and build status.

Returns:: Keys: source_dir, downloaded, binary_path, built, fc, cc, mklroot.
Return type:: dict

print_status()#

Print a human-readable source and build status report.

Return type:: None

Parameters:

config (Mare2DEMConfig | None)
source_dir (str | Path | None)

class pycsamt.models.mare2dem.Mare2DEMFileType(*values)#

Bases: Enum

Enumeration of recognized MARE2DEM file types.

EMDATA = 1#

RESISTIVITY = 2#

POLY = 3#

SETTINGS = 4#

LOG = 5#

RESPONSE = 6#

SENSITIVITY = 7#

UNKNOWN = 8#

pycsamt.models.mare2dem.detect_file_type(path)#

Detect the MARE2DEM file type from the file extension and name.

Parameters:: path (path-like) – Path to the file to classify.
Returns:: Detected file type, or Mare2DEMFileType.UNKNOWN.
Return type:: Mare2DEMFileType

pycsamt.models.mare2dem.is_emdata_file(path)#

Return True when path is a MARE2DEM observed-data file.

Parameters:: path (str | Path)
Return type:: bool

pycsamt.models.mare2dem.is_resistivity_file(path)#

Return True when path is a MARE2DEM resistivity model file.

Parameters:: path (str | Path)
Return type:: bool

pycsamt.models.mare2dem.is_settings_file(path)#

Return True when path is a MARE2DEM settings file.

Parameters:: path (str | Path)
Return type:: bool

pycsamt.models.mare2dem.is_log_file(path)#

Return True when path is a MARE2DEM log file.

Parameters:: path (str | Path)
Return type:: bool

pycsamt.models.mare2dem.is_response_file(path)#

Return True when path is a MARE2DEM predicted-response file.

Parameters:: path (str | Path)
Return type:: bool

pycsamt.models.mare2dem.code_label(code)#

Return "Component — Representation" for an integer data code.

Parameters:: code (int) – Integer data-type code from the DATA block of a .emdata file.
Returns:: Human-readable label, e.g. "Zxy (TE) — Phase" or "Unknown (code=7)".
Return type:: str

Examples

>>> code_label(104)
'Zxy (TE) — Phase'
>>> code_label(27)
'Ex — Log10 Amplitude'

pycsamt.models.mare2dem.code_component(code)#

Return the EM component name for code, or "" if unknown.

Parameters:: code (int)
Return type:: str

pycsamt.models.mare2dem.code_representation(code)#

Return the representation name for code, or "" if unknown.

Parameters:: code (int)
Return type:: str

pycsamt.models.mare2dem.is_mt_code(code)#

Return True when code belongs to the MT data type range.

Parameters:: code (int)
Return type:: bool

pycsamt.models.mare2dem.is_csem_code(code)#

Return True when code belongs to the CSEM data type range.

Parameters:: code (int)
Return type:: bool

class pycsamt.models.mare2dem.EMDataFile(path=None, format='EMData_2.3', is_response=False, comment='', utm=<factory>, csem=None, mt=None, dc=None, data=<factory>)#

Bases: object

Complete contents of one MARE2DEM .emdata or .EMResp file.

Variables:

path (pathlib.Path or None) – Source file path, if loaded from disk.
format (str) – Format string from the file header (e.g. "EMData_2.3").
is_response (bool) – True when the file is a response file (8-column DATA block).
comment (str) – Optional free-text comment from the file header.
utm (UTMOrigin) – UTM origin metadata.
csem (CSEMConfig or None) – CSEM configuration section, or None if absent.
mt (MTConfig or None) – MT configuration section, or None if absent.
dc (DCConfig or None) – DC configuration section, or None if absent.
data (numpy.ndarray, shape (n_data, 6) or (n_data, 8)) – DATA block. Columns for observed data: [type, freq#, tx#, rx#, data, std_err]. Response files add [response, residual].

Parameters:

path (Path | None)
format (str)
is_response (bool)
comment (str)
utm (UTMOrigin)
csem (CSEMConfig | None)
mt (MTConfig | None)
dc (DCConfig | None)
data (ndarray)

path: Path | None = None#

format: str = 'EMData_2.3'#

is_response: bool = False#

comment: str = ''#

utm: UTMOrigin#

csem: CSEMConfig | None = None#

mt: MTConfig | None = None#

dc: DCConfig | None = None#

data: ndarray#

property n_data: int#: Total number of data rows.

property n_mt_frequencies: int#: Number of MT frequencies.

property n_mt_receivers: int#: Number of MT receivers.

property n_csem_transmitters: int#: Number of CSEM transmitters.

property n_csem_receivers: int#: Number of CSEM receivers.

class pycsamt.models.mare2dem.UTMOrigin(grid=0, hemi='N', north0=0.0, east0=0.0, theta=0.0)#

Bases: object

UTM mapping metadata for the 2-D profile.

Variables:

grid (int) – UTM zone number.
hemi (str) – Hemisphere letter ('N' or 'S').
north0 (float) – Northing of profile origin in metres.
east0 (float) – Easting of profile origin in metres.
theta (float) – Profile strike direction (degrees).

Parameters:

grid (int)
hemi (str)
north0 (float)
east0 (float)
theta (float)

grid: int = 0#

hemi: str = 'N'#

north0: float = 0.0#

east0: float = 0.0#

theta: float = 0.0#

class pycsamt.models.mare2dem.CSEMConfig(phase_convention='lag', reciprocity_used='', frequencies=<factory>, time_offsets=None, tdem_waveform=None, transmitters=<factory>, transmitter_type=<factory>, transmitter_name=<factory>, receivers=<factory>, receiver_name=<factory>)#

Bases: object

CSEM source–receiver configuration and data parameters.

Variables:

phase_convention (str) – 'lag' or 'lead'.
reciprocity_used (str) – 'yes', 'no', or ''.
frequencies (numpy.ndarray, shape (n_freq,)) – CSEM frequencies in Hz.
time_offsets (numpy.ndarray or None) – Time offsets for TDEM modeling.
tdem_waveform (numpy.ndarray or None) – TDEM waveform, shape (n_pts, 2).
transmitters (numpy.ndarray, shape (n_tx, 7)) – Columns: x, y, z, azimuth, dip, length, solve_corr.
transmitter_type (list of str) – Dipole type per transmitter ('edipole' or 'bdipole').
transmitter_name (list of str) – Transmitter labels.
receivers (numpy.ndarray, shape (n_rx, 8)) – Columns: x, y, z, theta, alpha, beta, length, solve_corr.
receiver_name (list of str) – Receiver labels.

Parameters:

phase_convention (str)
reciprocity_used (str)
frequencies (ndarray)
time_offsets (ndarray | None)
tdem_waveform (ndarray | None)
transmitters (ndarray)
transmitter_type (list[str])
transmitter_name (list[str])
receivers (ndarray)
receiver_name (list[str])

phase_convention: str = 'lag'#

reciprocity_used: str = ''#

frequencies: ndarray#

time_offsets: ndarray | None = None#

tdem_waveform: ndarray | None = None#

transmitters: ndarray#

transmitter_type: list[str]#

transmitter_name: list[str]#

receivers: ndarray#

receiver_name: list[str]#

class pycsamt.models.mare2dem.MTConfig(frequencies=<factory>, receivers=<factory>, receiver_name=<factory>)#

Bases: object

MT receiver configuration and frequency list.

Variables:

frequencies (numpy.ndarray, shape (n_freq,)) – MT frequencies in Hz.
receivers (numpy.ndarray, shape (n_rx, 8)) – Columns: x, y, z, theta, alpha, beta, length, solve_static.
receiver_name (list of str) – Receiver labels (station names).

Parameters:

frequencies (ndarray)
receivers (ndarray)
receiver_name (list[str])

frequencies: ndarray#

receivers: ndarray#

receiver_name: list[str]#

class pycsamt.models.mare2dem.DCConfig(tx_electrodes=<factory>, rx_electrodes=<factory>, transmitters=<factory>, receivers=<factory>, transmitter_name=<factory>, receiver_name=<factory>)#

Bases: object

DC resistivity electrode and receiver configuration.

Variables:

tx_electrodes (numpy.ndarray, shape (n_tx_el, 3)) – Transmitter electrode positions (x, y, z).
rx_electrodes (numpy.ndarray, shape (n_rx_el, 3)) – Receiver electrode positions (x, y, z).
transmitters (numpy.ndarray, shape (n_tx, 2)) – Integer electrode-index pairs (A, B) for each transmitter.
receivers (numpy.ndarray, shape (n_rx, 2)) – Integer electrode-index pairs (M, N) for each receiver.
transmitter_name (list of str)
receiver_name (list of str)

Parameters:

tx_electrodes (ndarray)
rx_electrodes (ndarray)
transmitters (ndarray)
receivers (ndarray)
transmitter_name (list[str])
receiver_name (list[str])

tx_electrodes: ndarray#

rx_electrodes: ndarray#

transmitters: ndarray#

receivers: ndarray#

transmitter_name: list[str]#

receiver_name: list[str]#

pycsamt.models.mare2dem.read_emdata(path, *, silent=False)#

Read a MARE2DEM .emdata or .EMResp file.

Port of m2d_readEMData2DFile.m.

Parameters:

path (path-like) – File to read.
silent (bool, default False) – Suppress warnings on unrecognised format strings.

Returns:

Parsed file contents.

Return type:

Raises:

FileNotFoundError – When path does not exist.
ValueError – When the format string in the header is not recognised and silent=False.

Examples

>>> from pycsamt.models.mare2dem.iotools.emdata import read_emdata
>>> em = read_emdata("survey.emdata")
>>> em.n_mt_receivers
12

pycsamt.models.mare2dem.write_emdata(em, path)#

Write an EMDataFile to path.

Port of m2d_writeEMData2DFile.m.

Parameters:

em (EMDataFile) – Data to write.
path (path-like) – Destination .emdata file.

Returns:

Path of the written file.

Return type:

Examples

>>> from pycsamt.models.mare2dem.iotools.emdata import write_emdata
>>> write_emdata(em, "survey_out.emdata")
PosixPath('survey_out.emdata')

class pycsamt.models.mare2dem.ResistivityFile(resistivity_file='mare2dem.resistivity', poly_file='mare2dem.poly', data_file='mare2dem.emdata', settings_file='mare2dem.settings', version='mare2dem_1.1', anisotropy='isotropic', target_misfit=1.0, max_iterations=100, iteration=0, log10_lagrange=5.0, roughness=None, misfit=None, date_time='', bounds_transform='bandpass', global_bounds=<factory>, roughness_penalty_method='gradient', yz_penalty_weights=<factory>, penalty_cut_weight=0.1, roughness_with_prejudice=False, beta_mgs=0.0, anisotropy_penalty_weight=None, anisotropy_ratio_roughness_weight=None, debug_level=1, inversion_method='occam', rms_threshold=0.85, converge_slowly='no', resistivity=<factory>, free_parameter=<factory>, bounds=<factory>, prejudice=<factory>, data_group_file='', joint_inv_weight_type='', penalty_file='', fixed_mu_cut=None)#

Bases: object

Contents of one MARE2DEM .resistivity file.

Variables:

resistivity_file (str) – Output filename (used when writing).
poly_file (str) – Triangle mesh (.poly) file stem referenced by this model.
data_file (str) – .emdata file referenced by this model.
settings_file (str) – .settings file referenced by this model.
version (str) – Format version string (e.g. "mare2dem_1.1").
anisotropy (str) – Anisotropy mode: "isotropic" (default), "triaxial", "tix", "tiy", "tiz", "tiz_ratio", "isotropic_ip", "isotropic_complex".
target_misfit (float) – Target normalized RMS misfit.
max_iterations (int) – Maximum inversion iterations.
iteration (int) – Current iteration number (output from inversion).
log10_lagrange (float) – Log10 Lagrange (regularization) trade-off value.
roughness (float or None) – Model roughness at this iteration.
misfit (float or None) – Normalized misfit at this iteration.
date_time (str) – Date/time stamp written by MARE2DEM.
bounds_transform (str) – Bounds transform type ("bandpass").
global_bounds (numpy.ndarray, shape (2,)) – Lower and upper log10-resistivity bounds.
roughness_penalty_method (str) – Smoothing type ("gradient" or "first_difference").
yz_penalty_weights (numpy.ndarray, shape (2,)) – Smoothing weights (y, z).
penalty_cut_weight (float) – Penalty cut weight.
roughness_with_prejudice (bool) – Use prejudice in regularization.
beta_mgs (float) – Minimum gradient support weight.
anisotropy_penalty_weight (float or None) – Anisotropic penalty weight.
anisotropy_ratio_roughness_weight (float or None) – Anisotropic ratio roughness weight.
debug_level (int) – Verbosity / print level.
inversion_method (str) – Inversion method string.
rms_threshold (float) – Misfit decrease threshold.
converge_slowly (str) – "yes" or "no".
resistivity (numpy.ndarray, shape (n_regions, nrho)) – Linear resistivity (ohm-m) per region per component. Despite the “Global Bounds”/log10_lagrange fields elsewhere in this file being log10-scaled, the region table itself is linear – confirmed against a real compiled MARE2DEM binary: its reader (mare2dem_io.f90, comment "! rho is linear") applies log10() itself only internally, for its own inversion parameterization. Writing log10(rho) here instead produces a self-consistent but physically wrong forward response (magnitude off by a constant factor at every frequency; see pycsamt.forward.maxwell.mare2dem’s module docstring for the full story).
free_parameter (numpy.ndarray, shape (n_regions, nrho)) – Free-parameter index (0 = fixed, >0 = free parameter number).
bounds (numpy.ndarray, shape (n_regions, 2*nrho)) – Lower / upper bounds per region component, linear ohm-m (same convention as resistivity above, not log10).
prejudice (numpy.ndarray, shape (n_regions, 2*nrho)) – Prejudice value and weight per region component.

Parameters:

resistivity_file (str)
poly_file (str)
data_file (str)
settings_file (str)
version (str)
anisotropy (str)
target_misfit (float)
max_iterations (int)
iteration (int)
log10_lagrange (float)
roughness (float | None)
misfit (float | None)
date_time (str)
bounds_transform (str)
global_bounds (ndarray)
roughness_penalty_method (str)
yz_penalty_weights (ndarray)
penalty_cut_weight (float)
roughness_with_prejudice (bool)
beta_mgs (float)
anisotropy_penalty_weight (float | None)
anisotropy_ratio_roughness_weight (float | None)
debug_level (int)
inversion_method (str)
rms_threshold (float)
converge_slowly (str)
resistivity (ndarray)
free_parameter (ndarray)
bounds (ndarray)
prejudice (ndarray)
data_group_file (str)
joint_inv_weight_type (str)
penalty_file (str)
fixed_mu_cut (float | None)

resistivity_file: str = 'mare2dem.resistivity'#

poly_file: str = 'mare2dem.poly'#

data_file: str = 'mare2dem.emdata'#

settings_file: str = 'mare2dem.settings'#

version: str = 'mare2dem_1.1'#

anisotropy: str = 'isotropic'#

target_misfit: float = 1.0#

max_iterations: int = 100#

iteration: int = 0#

log10_lagrange: float = 5.0#

roughness: float | None = None#

misfit: float | None = None#

date_time: str = ''#

bounds_transform: str = 'bandpass'#

global_bounds: ndarray#

roughness_penalty_method: str = 'gradient'#

yz_penalty_weights: ndarray#

penalty_cut_weight: float = 0.1#

roughness_with_prejudice: bool = False#

beta_mgs: float = 0.0#

anisotropy_penalty_weight: float | None = None#

anisotropy_ratio_roughness_weight: float | None = None#

debug_level: int = 1#

inversion_method: str = 'occam'#

rms_threshold: float = 0.85#

converge_slowly: str = 'no'#

resistivity: ndarray#

free_parameter: ndarray#

bounds: ndarray#

prejudice: ndarray#

data_group_file: str = ''#

joint_inv_weight_type: str = ''#

penalty_file: str = ''#

fixed_mu_cut: float | None = None#

property num_regions: int#: Number of resistivity regions.

pycsamt.models.mare2dem.read_resistivity(path, *, no_data=False)#

Read a MARE2DEM .resistivity file.

Port of m2d_readResistivity.m.

Parameters:

path (path-like) – File to read.
no_data (bool, default False) – When True, stop reading before the region resistivity table. Useful for quickly inspecting header metadata in large iteration output files.

Returns:

Parsed file contents.

Return type:

ResistivityFile

Examples

>>> from pycsamt.models.mare2dem.iotools.resistivity import (
...     read_resistivity,
... )
>>> rf = read_resistivity("mare2dem.resistivity")
>>> rf.num_regions
1234

pycsamt.models.mare2dem.write_resistivity(rf, path=None)#

Write a ResistivityFile to path.

Port of m2d_writeResistivity.m.

Parameters:

rf (ResistivityFile) – Model to write.
path (path-like or None) – Destination file. When None, uses rf.resistivity_file.

Returns:

Path of the written file.

Return type:

Examples

>>> from pycsamt.models.mare2dem.iotools.resistivity import (
...     write_resistivity,
... )
>>> write_resistivity(rf, "mare2dem_iter10.resistivity")
PosixPath('mare2dem_iter10.resistivity')

class pycsamt.models.mare2dem.PolyFile#

Bases: object

Contents of one Triangle .poly PSLG file.

Variables:

nodes (numpy.ndarray, shape (n_nodes, 2)) – Node (x, y) coordinates. A MARE2DEM .poly file uses a right-handed coordinate system where the second coordinate is depth (z positive down).
node_attributes (numpy.ndarray, shape (n_nodes, n_attr) or shape (0,)) – Optional per-node attribute values.
node_boundary_markers (numpy.ndarray, shape (n_nodes,) or shape (0,)) – Optional per-node boundary marker integers.
segments (numpy.ndarray, shape (n_segs, 2) with int dtype) – Segment endpoint node-index pairs (1-based).
segment_markers (numpy.ndarray, shape (n_segs,) or shape (0,)) – Optional per-segment boundary marker integers.
holes (numpy.ndarray, shape (n_holes, 2) or shape (0, 2)) – Hole point coordinates.
regions (numpy.ndarray, shape (n_reg, 4) or shape (0, 4)) – Region attribute and area-constraint points. Columns: x, y, attribute, max_area.

property n_nodes: int#

property n_segments: int#

property n_holes: int#

property n_regions: int#

pycsamt.models.mare2dem.read_poly(path)#

Read a Triangle .poly PSLG file.

Port of m2d_readPoly.m.

Parameters:: path (path-like) – File to read. If the node count is zero (triangulation already done), the function automatically looks for the companion .node and .ele files in the same directory.
Returns:: Parsed PSLG contents.
Return type:: PolyFile

Examples

>>> from pycsamt.models.mare2dem.iotools.poly import read_poly
>>> poly = read_poly("mare2dem.poly")
>>> poly.n_nodes
4812

pycsamt.models.mare2dem.read_triangulation(node_path)#

Read a Triangle .node/.ele pair into plain arrays.

Shared by tri_mesh_from_poly() and PlotModel so there is one parser for Triangle’s element-file format, not two.

Parameters:

node_path (path-like) – The .node file. Its companion .ele file (same stem, same directory) is read alongside it.

Returns:

nodes (numpy.ndarray, shape (n_nodes, 2)) – Node (x, y) coordinates.
triangles (numpy.ndarray of int, shape (n_triangles, 3)) – 0-based node-index connectivity, normalized to 0-based regardless of whether the source files used Triangle’s -z (0-based) or default (1-based) numbering – detected from the .node file’s own first index.
region_attrs (numpy.ndarray of int, shape (n_triangles,)) – Triangle’s per-element region attribute (its 4th .ele column, present for -A/region-constrained runs), or all zeros when the .ele file carries no region attribute column.

Raises:

FileNotFoundError – If node_path or its companion .ele file is missing.

Return type:

tuple[ndarray, ndarray, ndarray]

Examples

Normally called on files written by run_triangle(), not constructed by hand.

pycsamt.models.mare2dem.write_poly(pf, path)#

Write a PolyFile to path.

Port of m2d_writePoly.m.

Parameters:

pf (PolyFile) – PSLG data to write.
path (path-like) – Destination .poly file.

Returns:

Path of the written file.

Return type:

Examples

>>> from pycsamt.models.mare2dem.iotools.poly import write_poly
>>> write_poly(pf, "mare2dem.poly")
PosixPath('mare2dem.poly')

pycsamt.models.mare2dem.write_triangulation(nodes, triangles, region_attrs, node_path)#

Write a .node/.ele/.neigh triple plus a companion stub .poly.

Inverse of read_triangulation(): hands an already-triangulated mesh (e.g. a TriMesh) to MARE2DEM directly, without re-invoking Triangle. Node/element indices are written 1-based, matching the convention MARE2DEM itself reads (Triangle’s default, distinct from the -z 0-based output run_triangle() requests for its own refined mesh).

The .neigh file (triangle-triangle adjacency, normally produced by Triangle’s own -n switch) is computed here directly, since we are bypassing Triangle entirely. Confirmed against a real compiled MARE2DEM binary that this file is required, not optional – MARE2DEM’s readPoly reads it unconditionally and stops with “Did you remember to include the ‘n’ flag when calling Triangle?” when it is missing.

Parameters:

nodes (array-like, shape (n_nodes, 2)) – Node (x, y) coordinates.
triangles (array-like of int, shape (n_triangles, 3)) – 0-based node-index connectivity.
region_attrs (array-like of int, shape (n_triangles,), optional) – Per-triangle region attribute. Defaults to all 1 when omitted.
node_path (path-like) – Destination .node file. The companion .ele/.neigh/ .poly files are written alongside it with the same stem.

Returns:

The written .node path.

Return type:

Raises:

ValueError – If the connectivity is non-manifold (an edge shared by more than two triangles), since no valid .neigh file can describe that.

Examples

>>> nodes = [[0, 0], [1, 0], [0, 1]]
>>> triangles = [[0, 1, 2]]
>>> path = write_triangulation(nodes, triangles, None, "mesh.node")
>>> read_triangulation(path)[1].tolist()
[[0, 1, 2]]

class pycsamt.models.mare2dem.SettingsFile(tolerance=1.0, tx_per_group=10, csem_rx_per_group=40, csem_freq_per_group=1, mt_rx_per_group=40, mt_freq_per_group=1, use_mesh_coarsening=True, use_mt_scattered_field=False, print_adaptive=True, print_decomposition=True)#

Bases: object

Parameters for one MARE2DEM .settings file.

Variables:

tolerance (float) – Target solution accuracy in percent.
tx_per_group (int) – Maximum transmitters per parallel group (≤ 10 recommended).
csem_rx_per_group (int) – CSEM receivers per parallel group.
csem_freq_per_group (int) – CSEM frequencies per group (usually 1).
mt_rx_per_group (int) – MT receivers per parallel group.
mt_freq_per_group (int) – MT frequencies per group (usually 1).
use_mesh_coarsening (bool) – Enable moving-window mesh coarsening for long profiles.
use_mt_scattered_field (bool) – Use scattered-field MT formulation (for deep-water scenarios).
print_adaptive (bool) – Print adaptive refinement iteration stats.
print_decomposition (bool) – Print parallel decomposition settings.

Parameters:

tolerance (float)
tx_per_group (int)
csem_rx_per_group (int)
csem_freq_per_group (int)
mt_rx_per_group (int)
mt_freq_per_group (int)
use_mesh_coarsening (bool)
use_mt_scattered_field (bool)
print_adaptive (bool)
print_decomposition (bool)

tolerance: float = 1.0#

tx_per_group: int = 10#

csem_rx_per_group: int = 40#

csem_freq_per_group: int = 1#

mt_rx_per_group: int = 40#

mt_freq_per_group: int = 1#

use_mesh_coarsening: bool = True#

use_mt_scattered_field: bool = False#

print_adaptive: bool = True#

print_decomposition: bool = True#

pycsamt.models.mare2dem.write_settings(sf, path, *, overwrite=True)#

Write a MARE2DEM .settings file.

Port of m2d_writeSettingsFile.m.

Parameters:

sf (SettingsFile) – Settings to write.
path (path-like) – Destination file.
overwrite (bool, default True) – If False and the file already exists, return the existing path without writing.

Returns:

Path of the written (or existing) file.

Return type:

Examples

>>> from pycsamt.models.mare2dem.iotools.settings import (
...     SettingsFile,
...     write_settings,
... )
>>> sf = SettingsFile(tx_per_group=5, csem_rx_per_group=20)
>>> write_settings(sf, "mare2dem.settings")
PosixPath('mare2dem.settings')

class pycsamt.models.mare2dem.GroupRMSLog(path=None, headers=<factory>, rms_log=<factory>)#

Bases: object

Contents of one MARE2DEM group-RMS log file.

Variables:

path (pathlib.Path or None) – Source file.
headers (list of str) – Column names from the first line of the file.
rms_log (numpy.ndarray, shape (n_iterations, n_groups)) – RMS values, one row per iteration.

Parameters:

path (Path | None)
headers (list[str])
rms_log (ndarray)

path: Path | None = None#

headers: list[str]#

rms_log: ndarray#

property n_iterations: int#

property n_groups: int#

pycsamt.models.mare2dem.read_group_rms_log(path)#

Read a MARE2DEM group-level RMS log file.

Port of m2d_read_group_rms_log.m.

Parameters:: path (path-like) – CSV-like log file written by MARE2DEM.
Returns:: Parsed log with headers and numeric RMS table.
Return type:: GroupRMSLog

Examples

>>> from pycsamt.models.mare2dem.iotools.group_rms import (
...     read_group_rms_log,
... )
>>> log = read_group_rms_log("mare2dem_group_rms.log")
>>> log.n_iterations
42

class pycsamt.models.mare2dem.DataGroupFile(path=None, comment='', group_names=<factory>, group_indices=<factory>)#

Bases: object

Contents of one MARE2DEM .emdata_group file.

Variables:

path (pathlib.Path or None) – Source file path.
comment (str) – Optional free-text comment from the file header.
group_names (list of str) – Ordered list of group name strings.
group_indices (numpy.ndarray, shape (n_data,)) – 1-based group index for each datum. The value at position i selects the group name at group_names[group_indices[i] - 1].

Parameters:

path (Path | None)
comment (str)
group_names (list[str])
group_indices (ndarray)

path: Path | None = None#

comment: str = ''#

group_names: list[str]#

group_indices: ndarray#

property n_groups: int#

property n_data: int#

pycsamt.models.mare2dem.read_data_group(path)#

Read a MARE2DEM .emdata_group file.

Port of m2d_readDataGroupFile.m.

Parameters:

path (path-like) – File to read (Format: EMDataGroup_1.0).

Returns:

Parsed data-group file.

Return type:

DataGroupFile

Raises:

FileNotFoundError – When path does not exist.
ValueError – When the format string is not EMDataGroup_1.0.

Examples

>>> from pycsamt.models.mare2dem.iotools.data_group import read_data_group
>>> dg = read_data_group("survey.emdata_group")
>>> dg.group_names
['MT', 'Seafloor CSEM', 'Towed CSEM']

pycsamt.models.mare2dem.write_data_group(dg, path)#

Write a DataGroupFile to path.

Port of m2d_writeDataGroupFile.m.

Parameters:

dg (DataGroupFile) – Data to write.
path (path-like) – Destination file.

Returns:

Path of the written file.

Return type:

Raises:

ValueError – When DataGroupFile.group_names is empty or DataGroupFile.group_indices are out of range.

Examples

>>> from pycsamt.models.mare2dem.iotools.data_group import (
...     DataGroupFile,
...     write_data_group,
... )
>>> import numpy as np
>>> dg = DataGroupFile(
...     group_names=["MT", "CSEM"], group_indices=np.array([1, 1, 2, 2])
... )
>>> write_data_group(dg, "survey.emdata_group")
PosixPath('survey.emdata_group')

pycsamt.models.mare2dem.get_most_recent(file_or_keyword, pattern='*.resistivity', *, search_dir='.')#

Return the most recently modified file matching pattern.

Port of m2d_getMostRecent.m.

Parameters:

file_or_keyword (str or path-like) – Either a literal file path, or one of the special keywords 'lastiter', 'last', or 'newest'. When a keyword is given the function scans search_dir for files matching pattern and returns the most recently modified one.
pattern (str, default “*.resistivity”) – Glob pattern used when a keyword is given.
search_dir (path-like, default ".") – Directory to search when a keyword is given.

Returns:

Resolved file path, or None when the keyword is used but no matching files are found.

Return type:

pathlib.Path or None

Examples

Load the latest inversion model:

>>> from pycsamt.models.mare2dem.iotools.most_recent import get_most_recent
>>> path = get_most_recent("newest", "*.resistivity", search_dir="./run")
>>> path
PosixPath('run/mare2dem.0020.resistivity')

Use a literal path (pass-through):

>>> get_most_recent("mare2dem.0010.resistivity")
PosixPath('mare2dem.0010.resistivity')

pycsamt.models.mare2dem.parse_topo(topo, y)#

Interpolate topography to profile positions y.

Port of m2d_parseTopo.m.

Parameters:

topo (float or array-like, shape (n_pts, 2)) – Topography input. A single scalar gives a flat surface at constant depth. An (n_pts, 2) array provides piecewise- linear topography with columns [y, z].
y (array-like, shape (n,)) – Profile positions at which to interpolate.

Returns:

z (numpy.ndarray, shape (n,)) – Depth of the topographic surface at each y position.
slope_angle (numpy.ndarray, shape (n,)) – Surface slope angle in degrees, positive clockwise from the +y axis towards +z (down). Zero for flat topography.
on_node (numpy.ndarray of bool, shape (n,)) – True where a y position falls exactly on a topography node. The slope is not well-defined at these points.

Return type:

tuple[ndarray, ndarray, ndarray]

Examples

Flat topography at depth 1000 m (e.g. flat seafloor):

>>> import numpy as np
>>> from pycsamt.models.mare2dem.geom.topo import parse_topo
>>> y = np.array([-5000.0, 0.0, 5000.0])
>>> z, slope, on_node = parse_topo(1000.0, y)
>>> z
array([1000., 1000., 1000.])

Variable topography:

>>> topo_xy = np.array([[0.0, 500.0], [5000.0, 1000.0], [10000.0, 800.0]])
>>> z, slope, on_node = parse_topo(topo_xy, np.array([2500.0, 5000.0]))

pycsamt.models.mare2dem.topo_depth(topo, y)#

Return interpolated topographic depth at positions y.

Convenience wrapper around parse_topo().

Parameters:

topo (float | ndarray)
y (ndarray)

Return type:

ndarray

pycsamt.models.mare2dem.topo_slope(topo, y)#

Return topographic slope angle (degrees) at positions y.

Convenience wrapper around parse_topo().

Parameters:

topo (float | ndarray)
y (ndarray)

Return type:

ndarray

pycsamt.models.mare2dem.get_intersections(xya, xyb, *, tol=np.float64(2.220446049250313e-13))#

Find intersections between segments in xya and segment xyb.

Port of m2d_getIntersections.m.

Parameters:

xya (array-like, shape (n_a, 4)) – Set of line segments. Each row is [x0, x1, y0, y1].
xyb (array-like, shape (1, 4) or (4,)) – Single query segment [x0, x1, y0, y1].
tol (float) – Interior-point tolerance (intersection on exact endpoint is excluded).

Returns:

intersect (numpy.ndarray of int) – Original row indices in xya where an interior intersection with xyb exists.
xi (numpy.ndarray of float) – x coordinates of intersection points.
yi (numpy.ndarray of float) – y coordinates of intersection points.
pa (numpy.ndarray of float, shape (n_a,)) – Parametric position along each segment in xya (-1 for segments not tested).
pb (numpy.ndarray of float, shape (n_a,)) – Parametric position along xyb (-1 for untested).

Return type:

tuple[ndarray, ndarray, ndarray, ndarray, ndarray]

Examples

>>> import numpy as np
>>> from pycsamt.models.mare2dem.geom.intersections import (
...     get_intersections,
... )
>>> segs_a = np.array([[0.0, 2.0, 1.0, 1.0]])  # horizontal segment
>>> seg_b = np.array([1.0, 1.0, 0.0, 2.0])  # vertical segment
>>> inter, xi, yi, pa, pb = get_intersections(segs_a, seg_b)
>>> xi
array([1.])
>>> yi
array([1.])

pycsamt.models.mare2dem.do_rects_overlap(a_rect, many_rects, *, tol=np.float64(2.220446049250313e-14))#

Return indices of rows in many_rects that overlap a_rect.

Port of the inner doRectsOverlap function in m2d_getIntersections.m.

Parameters:

a_rect (array-like, shape (4,)) – Single bounding rectangle [x0, x1, y0, y1].
many_rects (array-like, shape (n, 4)) – Array of bounding rectangles, same format.
tol (float) – Tolerance for the overlap test.

Returns:

Row indices of many_rects whose bounding boxes overlap a_rect.

Return type:

numpy.ndarray of int

pycsamt.models.mare2dem.dp_simplify(points, tolerance)#

Simplify a polyline using the Douglas-Peucker algorithm.

Port of m2d_dpsimplify.m.

Parameters:

points (array-like, shape (n, 2)) – Input polyline vertices.
tolerance (float) – Maximum allowable perpendicular deviation. Points whose distance from the simplified segment exceeds this value are retained.

Returns:

Simplified polyline with m <= n vertices.

Return type:

numpy.ndarray, shape (m, 2)

Examples

>>> import numpy as np
>>> from pycsamt.models.mare2dem.geom.simplify import dp_simplify
>>> pts = np.column_stack(
...     [np.linspace(0, 10, 100), np.sin(np.linspace(0, np.pi, 100))]
... )
>>> simplified = dp_simplify(pts, tolerance=0.05)
>>> len(simplified) < 100
True

pycsamt.models.mare2dem.simplify_poly(nodes, adjacency, *, tol=1e-10)#

Remove redundant interior nodes from a PSLG polygon.

Port of m2d_simplify_poly.m.

Parameters:

nodes (array-like, shape (n_nodes, 2)) – Node (y, z) coordinates.
adjacency (array-like, shape (n_nodes, n_nodes)) – Dense or sparse symmetric adjacency matrix. A non-zero adjacency[i, j] means node i and node j are connected by a segment. Diagonal entries are ignored (self-loops).
tol (float, default 1e-10) – Collinearity tolerance used to identify redundant nodes.

Returns:

nodes_out (numpy.ndarray, shape (m_nodes, 2)) – Pruned node array with redundant interior nodes removed.
adjacency_out (numpy.ndarray, shape (m_nodes, m_nodes)) – Updated adjacency matrix.

Return type:

Examples

>>> import numpy as np
>>> from pycsamt.models.mare2dem.geom.simplify_poly import simplify_poly
>>> nodes = np.array([[0.0, 0.0], [1.0, 0.0], [2.0, 0.0], [3.0, 0.0]])
>>> adj = np.zeros((4, 4))
>>> adj[0, 1] = adj[1, 0] = 1
... adj[1, 2] = adj[2, 1] = 1
... adj[2, 3] = adj[3, 2] = 1
>>> n_out, a_out = simplify_poly(nodes, adj)
>>> len(n_out)  # interior collinear nodes 1 and 2 removed
2

pycsamt.models.mare2dem.get_centroids(nodes, elements, tri_index)#

Compute area-weighted centroids of mesh regions.

Port of m2d_getCentroids.m.

Parameters:

nodes (array-like, shape (n_nodes, 2)) – Node (y, z) coordinates.
elements (array-like, shape (n_elements, 3)) – Element connectivity (1-based indices accepted).
tri_index (array-like, shape (n_elements,)) – Region index for each triangle (1-based).

Returns:

Columns: y_centroid, z_centroid, total_area of each region.

Return type:

numpy.ndarray, shape (n_regions, 3)

Examples

>>> import numpy as np
>>> from pycsamt.models.mare2dem.geom.centroids import get_centroids
>>> nodes = np.array([[0.0, 0.0], [1.0, 0.0], [0.0, 1.0], [1.0, 1.0]])
>>> elems = np.array([[1, 2, 3], [2, 4, 3]])
>>> tri_idx = np.array([1, 1])
>>> get_centroids(nodes, elems, tri_idx)
array([[0.5, 0.5, 1. ]])

pycsamt.models.mare2dem.triangle_centroids(nodes, elements)#

Return the centroid of each triangle element.

Parameters:

nodes (array-like, shape (n_nodes, 2)) – Node (y, z) coordinates.
elements (array-like, shape (n_elements, 3)) – Element connectivity, 0-based or 1-based node indices. The function detects 1-based indices automatically.

Returns:

Per-element centroid coordinates.

Return type:

numpy.ndarray, shape (n_elements, 2)

pycsamt.models.mare2dem.triangle_areas(nodes, elements)#

Return the signed area of each triangle element.

Parameters:

nodes (array-like, shape (n_nodes, 2))
elements (array-like, shape (n_elements, 3))

Returns:

Area of each triangle (absolute value).

Return type:

numpy.ndarray, shape (n_elements,)

pycsamt.models.mare2dem.lonlat_to_utm(lon, lat, *, zone=None, south_hemi=None, ellipsoid='wgs84')#

Convert geographic longitude / latitude to UTM.

Port of LonLat2UTM.m. When pyproj is installed the conversion is delegated to it; otherwise the pure-Python Snyder formulas are used.

Parameters:

lon (float or array-like) – Longitude in decimal degrees (−180 to +180).
lat (float or array-like) – Latitude in decimal degrees.
zone (int, optional) – Force a specific UTM zone. Auto-computed from median longitude when omitted.
south_hemi (bool, optional) – Force southern-hemisphere false northing. Auto-detected from median latitude when omitted.
ellipsoid (str, default "wgs84") – Ellipsoid name (see ELLIPSOIDS).

Returns:

easting (numpy.ndarray) – UTM easting in metres.
northing (numpy.ndarray) – UTM northing in metres.
zone (int) – UTM zone used.
south_hemi (bool) – Hemisphere flag used.

Return type:

tuple[ndarray, ndarray, int, bool]

Examples

>>> from pycsamt.models.mare2dem.geom.utm import lonlat_to_utm
>>> e, n, zone, sh = lonlat_to_utm(-70.0, 42.0)
>>> zone
19

pycsamt.models.mare2dem.utm_to_lonlat(easting, northing, zone, south_hemi=False, *, ellipsoid='wgs84')#

Convert UTM coordinates to geographic longitude / latitude.

Port of UTM2LonLat.m.

Parameters:

easting (float or array-like) – UTM easting in metres.
northing (float or array-like) – UTM northing in metres.
zone (int) – UTM zone number.
south_hemi (bool or str, default False) – True / 'S' for southern hemisphere.
ellipsoid (str, default "wgs84") – Ellipsoid name (see ELLIPSOIDS).

Returns:

lon (numpy.ndarray) – Longitude in decimal degrees.
lat (numpy.ndarray) – Latitude in decimal degrees.

Return type:

Examples

>>> from pycsamt.models.mare2dem.geom.utm import utm_to_lonlat
>>> lon, lat = utm_to_lonlat(330000.0, 4650000.0, 19, False)

pycsamt.models.mare2dem.estimate_area_of_interest(em)#

Estimate the y and z limits for display or mesh generation.

Port of m2d_estimateAreaOfInterest.m.

Parameters:

em (EMDataFile) – Survey data file supplying MT, CSEM, and DC receiver / transmitter positions.

Returns:

ylim (numpy.ndarray, shape (2,) or None) – [y_min, y_max] recommended profile extent in metres. None when no survey geometry is available.
zlim (numpy.ndarray, shape (2,) or None) – [z_min, z_max] recommended depth extent in metres. None when no survey geometry is available.

Return type:

tuple[ndarray | None, ndarray | None]

Examples

>>> from pycsamt.models.mare2dem import read_emdata
>>> from pycsamt.models.mare2dem.geom.area_of_interest import (
...     estimate_area_of_interest,
... )
>>> em = read_emdata("survey.emdata")
>>> ylim, zlim = estimate_area_of_interest(em)
>>> ylim
array([-5500.,  5500.])

pycsamt.models.mare2dem.survey_points(em)#

Return every receiver/transmitter/electrode (y, z) position.

Shared by estimate_area_of_interest() and build_survey_mesh() so there is one place that knows how to walk an EMDataFile’s MT/CSEM/DC geometry, not two.

Parameters:: em (EMDataFile) – Survey data file supplying MT, CSEM, and DC receiver/ transmitter positions.
Returns:: Stacked (y, z) positions, empty (shape (0, 2)) when em carries no geometry.
Return type:: numpy.ndarray, shape (n, 2)

Examples

>>> from pycsamt.models.mare2dem import read_emdata
>>> from pycsamt.models.mare2dem.geom.area_of_interest import (
...     survey_points,
... )
>>> em = read_emdata("survey.emdata")
>>> survey_points(em).shape[1]
2

pycsamt.models.mare2dem.get_triangle_regions(points, triangles, segments, region_seeds=None)#

Assign each triangle to a numbered region.

Port of m2d_getTriangleRegions.m.

Parameters:

points (array-like, shape (n_nodes, 2)) – Node (y, z) coordinates.
triangles (array-like, shape (n_tri, 3)) – Triangle connectivity, 1-based node indices.
segments (array-like, shape (n_segs, 2)) – Boundary segment node-index pairs (1-based) that separate regions.
region_seeds (array-like, shape (n_seeds, 2) or None) – One seed point per pre-defined region. Each seed is a (y, z) coordinate inside a region. Pass None if no pre-defined seeds are given; all regions are then discovered automatically.

Returns:

tri_index (numpy.ndarray, shape (n_tri,)) – Region index for each triangle (1-based).
region_map (numpy.ndarray, shape (n_new_regions,)) – Maps new region indices back to input seed indices (1-based); 0 for regions with no corresponding seed.

Return type:

Notes

The algorithm builds a sparse adjacency matrix from the boundary segments, zeros out neighbour links that cross boundaries, then performs a BFS flood fill from each seed (or from untouched triangles after the seed pass).

Examples

>>> import numpy as np
>>> from pycsamt.models.mare2dem.geom.triangle_regions import (
...     get_triangle_regions,
... )
>>> pts = np.array([[0.0, 0.0], [1.0, 0.0], [0.5, 1.0]])
>>> tris = np.array([[1, 2, 3]])
>>> segs = np.array([[1, 2]])
>>> tri_index, region_map = get_triangle_regions(pts, tris, segs)
>>> tri_index
array([1])

pycsamt.models.mare2dem.get_line_orientation(northings, eastings)#

Estimate the survey profile orientation from station UTM positions.

Port of m2d_getLineOrientation.m.

The function fits a line to the input (northing, easting) pairs and returns its bearing. The result is in the range 0–180°:

0° / 180° → N–S profile
90° → E–W profile

Parameters:

northings (array-like, shape (n,)) – UTM northing coordinates in metres.
eastings (array-like, shape (n,)) – UTM easting coordinates in metres.

Returns:

Line orientation in degrees clockwise from geographic north (0 ≤ result ≤ 180).

Return type:

float

Examples

>>> import numpy as np
>>> from pycsamt.models.mare2dem.geom.line_orientation import (
...     get_line_orientation,
... )
>>> northings = np.array([0.0, 1000.0, 2000.0])
>>> eastings = np.array([0.0, 0.0, 0.0])
>>> get_line_orientation(northings, eastings)  # N-S profile
0.0

pycsamt.models.mare2dem.project_onto_line(northings, eastings, utm0_north, utm0_east, line_orientation)#

Project UTM (northing, easting) positions onto a survey profile.

Parameters:

northings (array-like) – Station UTM northings in metres.
eastings (array-like) – Station UTM eastings in metres.
utm0_north (float) – Profile origin northing in metres.
utm0_east (float) – Profile origin easting in metres.
line_orientation (float) – Profile orientation in degrees clockwise from north.

Returns:

x (numpy.ndarray) – Cross-profile (x) offset in metres.
y (numpy.ndarray) – Along-profile (y) offset in metres.

Return type:

pycsamt.forward.maxwell.contracts_tri.TriMesh

pycsamt.models.mare2dem.build_pslg(boundary, *, interior_points=None, region_seeds=None)#

Assemble a closed-polygon PSLG from boundary and interior points.

Parameters:

boundary (array-like, shape (n, 2)) – Ordered (y, z) vertices of a single closed outer boundary polygon (do not repeat the first point at the end). Consecutive vertices become PSLG boundary segments, including the closing edge from the last vertex back to the first.
interior_points (array-like, shape (m, 2), optional) – Extra points Triangle must place a mesh node at exactly (e.g. receiver locations), added with no connecting segments.
region_seeds (array-like, shape (k, 2), optional) – One interior point per resistivity region, written as Triangle region-constraint rows (attribute i + 1, no area cap). Omit for a single uniform region (region_ids will be all zero).

Returns:

PSLG ready for write_poly().

Return type:

PolyFile

Examples

>>> box = [[0, 0], [1000, 0], [1000, 500], [0, 500]]
>>> pf = build_pslg(box)
>>> pf.n_nodes, pf.n_segments
(4, 4)

pycsamt.models.mare2dem.build_survey_mesh(em, *, topo=0.0, depth_below_topo_m=None, simplify_tolerance_m=None, min_angle=30.0, max_area=None, workdir='.', stem='mesh', executable=None, config=None)#

Build a topography-aware triangular mesh for a survey.

Composes estimate_area_of_interest(), parse_topo(), dp_simplify(), build_pslg(), and run_triangle(), then converts the refined mesh via tri_mesh_from_poly().

Parameters:

em (EMDataFile) – Survey data supplying receiver/transmitter positions, used both for the area-of-interest estimate and as required PSLG interior points so Triangle places a node exactly at every receiver.
topo (float or array-like, shape (n, 2), default=0.0) – Topography, forwarded to parse_topo(). A scalar gives a flat surface at that depth.
depth_below_topo_m (float, optional) – Bottom-boundary depth. Defaults to the estimated area-of-interest zlim[1].
simplify_tolerance_m (float, optional) – Douglas-Peucker tolerance applied to the sampled topography polyline before meshing. Left unsimplified when None.
min_angle (float, optional) – Forwarded to run_triangle().
max_area (float, optional) – Forwarded to run_triangle().
workdir (path-like, default=".") – Directory the .poly/.node/.ele files are written to.
stem (str, default="mesh") – Base filename (without extension) for the written mesh files.
executable (PathLike | None) – Forwarded to run_triangle().
config (Mare2DEMConfig | None) – Forwarded to run_triangle().

Returns:

Triangular mesh spanning the survey’s area of interest, with a node at every receiver position.

Return type:

Raises:

ValueError – If em has no receiver/transmitter geometry to mesh.
pycsamt.models.mare2dem.triangle_exec.TriangleRunError – If no Triangle executable is found or refinement fails.

Examples

Requires a real Triangle executable, so this is illustrative only:

>>> mesh = build_survey_mesh(em, topo=0.0)

pycsamt.models.mare2dem.tri_mesh_from_poly(refined_poly_path)#

Convert a Triangle-refined .poly mesh to a TriMesh.

Parameters:: refined_poly_path (path-like) – The .1.poly file produced by run_triangle() (or any .poly/.node path with a companion .node/.ele pair of the same stem).
Returns:: Validated triangular mesh; region_ids are Triangle’s own per-element region attribute (all zero when the mesh was built without region constraints).
Return type:: pycsamt.forward.maxwell.contracts_tri.TriMesh
Raises:: FileNotFoundError – If the companion .node/.ele files are missing.

Examples

Normally called on run_triangle()’s own return value:

>>> from pycsamt.models.mare2dem.triangle_exec import run_triangle
>>> refined = run_triangle("mesh.poly")
>>> mesh = tri_mesh_from_poly(refined)

pycsamt.models.mare2dem.run_triangle(poly_path, *, executable=None, min_angle=30.0, max_area=None, extra_args=(), config=None, timeout=None)#

Refine a Triangle .poly PSLG into a quality FEM mesh.

Parameters:

poly_path (path-like) – Input .poly PSLG file, e.g. from write_poly().
executable (path-like or sequence of path-like, optional) – Explicit Triangle executable, or an [interpreter, script] sequence (e.g. [sys.executable, "fake_triangle.py"]) for a scripted test double standing in for a real binary. Resolved via resolve_triangle_binary() when omitted.
min_angle (float or None, default=30.0) – Minimum triangle angle in degrees, Triangle’s -q quality constraint. Pass None to disable it.
max_area (float or None, optional) – Maximum triangle area, Triangle’s -a constraint.
extra_args (sequence of str, optional) – Additional flags appended verbatim (e.g. ("-Y",) to forbid Steiner points on PSLG boundary segments).
config (Mare2DEMConfig, optional) – Forwarded to SourceManager when resolving the executable.
timeout (float, optional) – Subprocess timeout in seconds.

Returns:

Path to the refined <stem>.1.poly file – Triangle’s own output-naming convention when run with -p on <stem>.poly. Its companion <stem>.1.node/<stem>.1.ele files sit alongside it; read all three together with read_triangulation().

Return type:

Raises:

FileNotFoundError – If poly_path does not exist.
TriangleRunError – If no executable is found, the process exits non-zero, or the expected .1.poly output is missing afterward.

Examples

>>> run_triangle("missing.poly")
Traceback (most recent call last):
...
FileNotFoundError: poly_path not found: missing.poly

exception pycsamt.models.mare2dem.TriangleRunError#

Bases: RuntimeError

Raised when Triangle cannot be found, fails, or produces no output.

Examples

>>> isinstance(TriangleRunError("no triangle"), RuntimeError)
True

class pycsamt.models.mare2dem.MTSurveyConfig(frequencies=<factory>, rx_y=<factory>, rx_type='marine', rx_z=None, rx_z_offset=None, rx_beta=None, rx_name=None, lTE=True, lTM=True, lZDet=False, lTipper=False, lTipperRealImag=False, lMTFields=False)#

Bases: object

MT receiver and data-selection configuration.

Parameters:

frequencies (array-like) – MT frequencies in Hz.
rx_y (array-like) – Receiver y positions along the 2-D profile.
rx_type (str) – Receiver placement mode: 'land', 'marine', or 'amphibious'.
rx_z (array-like or None) – Override receiver depth positions. When supplied, slope- based tilt angles are not computed.
rx_z_offset (float or array-like or None) – Offset from topography (< 0 above, > 0 below). Overrides the rx_type depth rule.
rx_beta (float or array-like or None) – Override receiver beta (y-tilt) angles in degrees.
rx_name (list of str or None) – Station names.
lTE (bool) – Include log10 apparent resistivity and phase (TE mode).
lTM (bool) – Include log10 apparent resistivity and phase (TM mode).
lZDet (bool) – Include impedance-determinant apparent resistivity and phase.
lTipper (bool) – Include TE tipper amplitude and phase.
lTipperRealImag (bool) – Include TE tipper real and imaginary parts.
lMTFields (bool) – Include TE and TM field vectors (Ex, Ey, Ez, Hx, Hy, Hz).

frequencies: ndarray#

rx_y: ndarray#

rx_type: str = 'marine'#

rx_z: ndarray | None = None#

rx_z_offset: float | ndarray | None = None#

rx_beta: float | ndarray | None = None#

rx_name: list[str] | None = None#

lTE: bool = True#

lTM: bool = True#

lZDet: bool = False#

lTipper: bool = False#

lTipperRealImag: bool = False#

lMTFields: bool = False#

class pycsamt.models.mare2dem.CSEMSurveyConfig(frequencies=<factory>, tx_y=<factory>, rx_y=None, rx_r=None, tx_z=None, tx_z_offset=None, rx_z=None, rx_z_offset=None, rx_type='marine', tx_type='edipole', tx_azimuth=None, tx_dip=None, tx_length=None, rx_length=None, rx_beta=None, tx_name=None, rx_name=None, phase_convention='lag', min_range=None, max_range=None, lEx=True, lEy=False, lEz=False, lBx=False, lBy=False, lBz=False)#

Bases: object

CSEM source–receiver and data-selection configuration.

Parameters:

frequencies (array-like) – CSEM frequencies in Hz.
tx_y (array-like) – Transmitter y positions.
rx_y (array-like or None) – Receiver y positions. Mutually exclusive with rx_r.
rx_r (array-like or None) – Towed receiver offsets from transmitter positions.
tx_z (array-like or None) – Override transmitter depth positions.
tx_z_offset (float or array-like or None) – Transmitter offset relative to topography.
rx_z (array-like or None) – Override receiver depth positions.
rx_z_offset (float or array-like or None) – Receiver offset relative to topography.
rx_type (str) – Receiver placement mode.
tx_type (str) – Dipole type: 'edipole' or 'bdipole'.
tx_azimuth (float or array-like or None) – Transmitter azimuth (degrees from x towards y).
tx_dip (float or array-like or None) – Transmitter dip (degrees positive down).
tx_length (float or array-like or None) – Electric dipole length in metres.
rx_length (float or array-like or None) – Receiver dipole length in metres.
rx_beta (float or array-like or None) – Receiver beta (y-tilt) angles in degrees.
tx_name (list of str or None)
rx_name (list of str or None)
phase_convention (str) – 'lag' (default) or 'lead'.
min_range (float or None) – Minimum Tx–Rx range to include in data file.
max_range (float or None) – Maximum Tx–Rx range.
lEx (bool) – Include Ex log10 amplitude and phase.
lEy (bool)
lEz (bool)
lBx (bool)
lBy (bool)
lBz (bool)

frequencies: ndarray#

tx_y: ndarray#

rx_y: ndarray | None = None#

rx_r: ndarray | None = None#

tx_z: ndarray | None = None#

tx_z_offset: float | ndarray | None = None#

rx_z: ndarray | None = None#

rx_z_offset: float | ndarray | None = None#

rx_type: str = 'marine'#

tx_type: str = 'edipole'#

tx_azimuth: float | ndarray | None = None#

tx_dip: float | ndarray | None = None#

tx_length: float | ndarray | None = None#

rx_length: float | ndarray | None = None#

rx_beta: float | ndarray | None = None#

tx_name: list[str] | None = None#

rx_name: list[str] | None = None#

phase_convention: str = 'lag'#

min_range: float | None = None#

max_range: float | None = None#

lEx: bool = True#

lEy: bool = False#

lEz: bool = False#

lBx: bool = False#

lBy: bool = False#

lBz: bool = False#

pycsamt.models.mare2dem.make_data_file(out_file, topo, *, mt=None, csem=None, file_type='forward', comment='Created by pycsamt make_data_file', utm=None)#

Create a MARE2DEM .emdata file from survey parameters.

Port of m2d_makeDataFile.m.

Parameters:

out_file (path-like) – Output .emdata filename.
topo (float or array-like (n, 2)) – Topography: scalar depth (flat) or [y, z] table.
mt (MTSurveyConfig, optional) – MT receiver/data configuration.
csem (CSEMSurveyConfig, optional) – CSEM transmitter/receiver/data configuration.
file_type (str, default "forward") – "forward" or "inversion" (inversion data not yet implemented).
comment (str) – Optional comment written on the second line of the file.
utm (UTMOrigin, optional) – UTM mapping metadata.

Returns:

The constructed data object (also written to out_file).

Return type:

Examples

MT forward dataset at 10 frequencies, 20 land stations:

>>> import numpy as np
>>> from pycsamt.models.mare2dem.survey import (
...     MTSurveyConfig,
...     make_data_file,
... )
>>> freqs = np.logspace(-3, 3, 10)
>>> rx_y = np.linspace(-10000, 10000, 20)
>>> mt_cfg = MTSurveyConfig(
...     frequencies=freqs, rx_y=rx_y, rx_type="land", lTE=True, lTM=True
... )
>>> em = make_data_file("survey.emdata", topo=0.0, mt=mt_cfg)
>>> em.n_mt_receivers
20

class pycsamt.models.mare2dem.ZMMStation(name='', latitude=0.0, longitude=0.0, declination=0.0, periods=<factory>, apres_te=<factory>, phase_te=<factory>, apres_te_se=<factory>, phase_te_se=<factory>, apres_tm=<factory>, phase_tm=<factory>, apres_tm_se=<factory>, phase_tm_se=<factory>, tipper_zy=None, tipper_zy_se=None, x_profile=0.0, y_profile=0.0, z_profile=0.0)#

Bases: object

Impedance tensor data from one .zmm file.

Variables:

name (str) – Station name.
latitude (float) – Geographic latitude in decimal degrees.
longitude (float) – Geographic longitude in decimal degrees.
declination (float) – Geomagnetic declination to apply (degrees).
periods (numpy.ndarray) – Period samples in seconds.
apres_te (numpy.ndarray) – TE-mode apparent resistivity (Ω·m).
phase_te (numpy.ndarray) – TE-mode phase (degrees).
apres_te_se (numpy.ndarray) – TE apparent resistivity standard error.
phase_te_se (numpy.ndarray) – TE phase standard error (degrees).
apres_tm (numpy.ndarray) – TM-mode apparent resistivity.
phase_tm (numpy.ndarray) – TM-mode phase.
apres_tm_se (numpy.ndarray)
phase_tm_se (numpy.ndarray)
tipper_zy (numpy.ndarray or None) – Complex TE tipper values.
tipper_zy_se (numpy.ndarray or None) – Tipper standard errors.

Parameters:

name (str)
latitude (float)
longitude (float)
declination (float)
periods (ndarray)
apres_te (ndarray)
phase_te (ndarray)
apres_te_se (ndarray)
phase_te_se (ndarray)
apres_tm (ndarray)
phase_tm (ndarray)
apres_tm_se (ndarray)
phase_tm_se (ndarray)
tipper_zy (ndarray | None)
tipper_zy_se (ndarray | None)
x_profile (float)
y_profile (float)
z_profile (float)

name: str = ''#

latitude: float = 0.0#

longitude: float = 0.0#

declination: float = 0.0#

periods: ndarray#

apres_te: ndarray#

phase_te: ndarray#

apres_te_se: ndarray#

phase_te_se: ndarray#

apres_tm: ndarray#

phase_tm: ndarray#

apres_tm_se: ndarray#

phase_tm_se: ndarray#

tipper_zy: ndarray | None = None#

tipper_zy_se: ndarray | None = None#

x_profile: float = 0.0#

y_profile: float = 0.0#

z_profile: float = 0.0#

pycsamt.models.mare2dem.read_zmm(path)#

Read one EMTF .zmm impedance file.

Parameters:: path (path-like) – File to read.
Returns:: Parsed station data.
Return type:: ZMMStation
Raises:: FileNotFoundError – When path does not exist.

Notes

The reader parses the ZMM 2-column-per-component format:

Period ZXX.r ZXX.i ZXY.r ZXY.i ZYX.r ZYX.i ZYY.r ZYY.i
TZX.r TZX.i TZY.r TZY.i Coh_XY Coh_YX

Apparent resistivities and phases are derived from the complex impedance components ZXY (TE) and ZYX (TM) using:

ρ_a = |Z|² / (ω · μ₀) φ = atan2(Z.imag, Z.real)

where ω = 2π / T and μ₀ = 4π × 10⁻⁷ H/m.

pycsamt.models.mare2dem.make_mt_data_from_zmm(zmm_files, out_file, *, output_modes='all', error_floor_te=0.0, error_floor_tm=0.0, error_floor_tipper=0.0, omit_periods=None, line_orientation=None, declination=0.0, utm0=None, utm_zone='', topo=0.0, rx_z_offset=-0.1)#

Create a MARE2DEM MT data file from a list of .zmm files.

Port of m2d_makeMTDataFromZmm.m.

Parameters:

zmm_files (list of path-like) – .zmm impedance files in profile order.
out_file (path-like) – Output .emdata filename.
output_modes (str, default "all") – Comma-separated data types to include: "TE", "TM", "tipper", "TE+tipper", "all impedance", or "all" (TE + TM + tipper).
error_floor_te (float) – Relative error floor for TE (e.g. 0.05 = 5 %).
error_floor_tm (float) – Relative error floor for TM.
error_floor_tipper (float) – Absolute tipper error floor.
omit_periods (array-like, shape (n, 2) or None) – Period bands to omit: each row is [T_min, T_max] in seconds.
line_orientation (float or None) – Profile orientation in degrees clockwise from North. Auto-fit from station positions when None.
declination (float, default 0.0) – Geomagnetic declination correction in degrees.
utm0 ((float, float) or None) – Profile UTM origin (northing, easting) in metres.
utm_zone (str, default "") – UTM zone string, e.g. "19N". Auto-detected when empty.
topo (float or array-like) – Topography for receiver placement (see parse_topo()).
rx_z_offset (float, default -0.1) – Vertical offset of receivers relative to topography (m). Negative → above (marine), positive → below (land).

Returns:

Constructed data file (also written to out_file).

Return type:

Examples

>>> from pycsamt.models.mare2dem.zmm import make_mt_data_from_zmm
>>> em = make_mt_data_from_zmm(
...     ["S001.zmm", "S002.zmm"],
...     "line1_mt.emdata",
...     output_modes="TE+tipper",
...     error_floor_te=0.05,
... )
>>> em.n_mt_receivers
2

pycsamt.models.mare2dem.make_mt_data_from_stations(stations, out_file, *, output_modes='all', error_floor_te=0.0, error_floor_tm=0.0, error_floor_tipper=0.0, omit_periods=None, line_orientation=None, declination=0.0, utm0=None, utm_zone='', topo=0.0, rx_z_offset=-0.1)#

Create a MARE2DEM MT data file from prepared ZMMStation objects.

Backend shared by make_mt_data_from_zmm() (stations parsed from .zmm files) and pycsamt.models.mare2dem.edi.make_mt_data_from_edi() (stations converted from EDI impedances). See make_mt_data_from_zmm() for the parameter documentation; the only difference is that stations are supplied directly.

Parameters:

stations (list[ZMMStation])
out_file (str | Path)
output_modes (str)
error_floor_te (float)
error_floor_tm (float)
error_floor_tipper (float)
omit_periods (ndarray | None)
line_orientation (float | None)
declination (float)
utm0 (tuple[float, float] | None)
utm_zone (str)
topo (Any)
rx_z_offset (float)

Return type:

EMDataFile

pycsamt.models.mare2dem.stations_from_edi(source, *, default_rel_error=0.05, confidence_weighting=False, confidence_method='composite', confidence_weights=None, confidence_min=0.05, confidence_power=1.0)#

Convert an EDI source into ZMMStation objects.

Parameters:

source (path-like, Sites, or EDI collection) – Anything accepted by pycsamt.emtools._core.ensure_sites().
default_rel_error (float, default 0.05) – Relative impedance error assumed when the EDIs carry no error block (5 %). Error floors applied later can only raise it.
confidence_weighting (bool, default False) – If True, inflate the relative impedance errors by the frequency-level confidence ratio (CR) before apparent-resistivity and phase standard errors are computed.
confidence_method (str) – Passed to pycsamt.emtools.qc.frequency_confidence_table().
confidence_weights (Mapping[str, float] | None) – Passed to pycsamt.emtools.qc.frequency_confidence_table().
confidence_min (float, default 0.05) – Lower bound applied to CR before inverting it. This prevents one severely degraded datum from producing infinite uncertainty.
confidence_power (float, default 1.0) – Exponent in the uncertainty multiplier (1 / max(CR, confidence_min)) ** confidence_power.

Returns:

One station per EDI with valid impedance, coordinates, and a frequency table matching the first station’s. TE is Zxy, TM is Zyx; apparent resistivity uses the field-units convention ρ_a = 0.2 |Z|² / f and TM phase is wrapped by +180° into the first quadrant (same conventions as read_zmm()).

Return type:

list of ZMMStation

Raises:

ValueError – When no station yields usable impedance + coordinates.

pycsamt.models.mare2dem.make_mt_data_from_edi(source, out_file, *, output_modes='all', error_floor_te=0.05, error_floor_tm=0.05, error_floor_tipper=0.0, default_rel_error=0.05, confidence_weighting=False, confidence_method='composite', confidence_weights=None, confidence_min=0.05, confidence_power=1.0, **kwargs)#

Create a MARE2DEM MT .emdata file from EDI data.

Parameters:

source (path-like, Sites, or EDI collection) – Anything accepted by pycsamt.emtools._core.ensure_sites().
out_file (path-like) – Output .emdata filename.
output_modes (str) – See make_mt_data_from_zmm().
error_floor_te (float) – See make_mt_data_from_zmm().
error_floor_tm (float) – See make_mt_data_from_zmm().
error_floor_tipper (float) – See make_mt_data_from_zmm().
default_rel_error (float, default 0.05) – Assumed relative impedance error when the EDIs carry none.
confidence_weighting (bool, default False) – If enabled, CR-derived uncertainty inflation is applied before MARE2DEM data weights are written. The propagated errors are sigma_rho = 2 rho sigma_Z/|Z| and sigma_phi = degrees(sigma_Z/|Z|) after multiplying the relative impedance error by (1 / max(CR, confidence_min)) ** confidence_power.
**kwargs (Any) – Remaining keyword arguments of make_mt_data_from_stations() (omit_periods, line_orientation, utm0, utm_zone, topo, rx_z_offset, declination).
confidence_method (str)
confidence_weights (Mapping[str, float] | None)
confidence_min (float)
confidence_power (float)
**kwargs

Returns:

Constructed data file (also written to out_file).

Return type:

Examples

>>> from pycsamt.models.mare2dem.edi import make_mt_data_from_edi
>>> em = make_mt_data_from_edi(
...     "data/AMT/WILLY_DATA/L22PLT",
...     "run/mare2dem.emdata",
...     error_floor_te=0.05,
...     error_floor_tm=0.05,
... )
>>> em.n_mt_receivers
25

class pycsamt.models.mare2dem.NoiseConfig(mt_rel_noise=0.05, mt_rel_noise_tipper=0.01, mt_abs_noise_tipper=0.01, csem_rel_noise_e=0.05, csem_rel_noise_b=0.05, max_sigma_factor=2.0)#

Bases: object

Per-data-type noise specifications.

Variables:

mt_rel_noise (float) – Relative noise for MT apparent resistivity (e.g. 0.05 = 5 %). Phase noise is mt_rel_noise / 2.
mt_rel_noise_tipper (float) – Relative noise for tipper real/imaginary components.
mt_abs_noise_tipper (float) – Absolute noise floor for tipper. Tipper data with amplitude below this value are dropped (NaN).
csem_rel_noise_e (float) – Relative noise for CSEM electric-field responses.
csem_rel_noise_b (float) – Relative noise for CSEM magnetic-field responses.
max_sigma_factor (float, default 2.0) – Clip noise at this many standard deviations to avoid unrealistically large random draws (matches MATLAB abs(randn) > 2 clipping).

Parameters:

mt_rel_noise (float)
mt_rel_noise_tipper (float)
mt_abs_noise_tipper (float)
csem_rel_noise_e (float)
csem_rel_noise_b (float)
max_sigma_factor (float)

mt_rel_noise: float = 0.05#

mt_rel_noise_tipper: float = 0.01#

mt_abs_noise_tipper: float = 0.01#

csem_rel_noise_e: float = 0.05#

csem_rel_noise_b: float = 0.05#

max_sigma_factor: float = 2.0#

pycsamt.models.mare2dem.add_synthetic_noise(em, noise, *, seed=None)#

Add synthetic Gaussian noise to MARE2DEM response data.

Port of m2d_addSyntheticNoise.m.

The function reads predicted responses from column 7 of the 8-column DATA block (is_response=True), adds noise according to noise, and stores the result in columns 5 and 6 (data, std_err). The output DATA block has 6 columns (observed-data format).

Parameters:

em (EMDataFile) – Response file with 8-column DATA block (Format: EMResp_2.3).
noise (NoiseConfig) – Noise specification.
seed (int or None, default None) – Random-number generator seed for reproducibility.

Returns:

New EMDataFile with 6-column DATA block containing the noisy synthetic data and standard errors.

Return type:

Raises:

ValueError – When em.data does not have 8 columns (not a response file).

Examples

>>> from pycsamt.models.mare2dem import read_emdata
>>> from pycsamt.models.mare2dem.noise import (
...     NoiseConfig,
...     add_synthetic_noise,
... )
>>> resp = read_emdata("forward.EMResp")
>>> nc = NoiseConfig(mt_rel_noise=0.05)
>>> noisy = add_synthetic_noise(resp, nc, seed=42)
>>> noisy.n_data
200

pycsamt.models.mare2dem.make_synthetic_data(in_file, out_file, noise, *, seed=None)#

Read a MARE2DEM response file, add noise, and write a data file.

Port of m2d_makeSyntheticData.m.

Parameters:

in_file (path-like) – MARE2DEM response file (.EMResp) with 8-column DATA block.
out_file (path-like) – Output .emdata file with 6-column noisy DATA block.
noise (NoiseConfig) – Noise specification.
seed (int or None, default None) – RNG seed for reproducibility.

Returns:

The noisy data object (also written to out_file).

Return type:

Examples

>>> from pycsamt.models.mare2dem.noise import (
...     NoiseConfig,
...     make_synthetic_data,
... )
>>> nc = NoiseConfig(mt_rel_noise=0.05, csem_rel_noise_e=0.05)
>>> em = make_synthetic_data(
...     "forward.EMResp", "synthetic.emdata", nc, seed=0
... )

pycsamt.models.mare2dem.merge_data_files(files_to_merge, out_file, *, keep_duplicate_rx=False)#

Merge MARE2DEM .emdata files and write the result.

Port of m2d_mergeDataFiles.m.

Parameters:

files_to_merge (list of path-like) – At least two .emdata files to merge. Include the full path if files are not in the current directory.
out_file (path-like) – Output .emdata filename.
keep_duplicate_rx (bool, default False) – Keep identical receiver locations (needed for towed CSEM arrays).

Returns:

Merged file contents (also written to out_file).

Return type:

Examples

>>> from pycsamt.models.mare2dem.merge import merge_data_files
>>> em = merge_data_files(["mt.emdata", "csem.emdata"], "joint.emdata")
>>> em.n_data
3200

pycsamt.models.mare2dem.merge_emdata(files, *, keep_duplicate_rx=False, comment='')#

Merge a list of EMDataFile objects into one.

Parameters:

files (list of EMDataFile) – At least two files to merge.
keep_duplicate_rx (bool, default False) – When True, identical receiver locations are kept separate (needed for towed CSEM arrays with identical offsets).
comment (str) – Comment written into the merged output file header.

Returns:

Merged data file (not yet written to disk).

Return type:

Raises:

ValueError – When fewer than two files are given, or when UTM origins differ.

pycsamt.models.mare2dem.grid_to_mare2dem(Y, Z, Rho, *, padding_y=50000.0, padding_z=50000.0, out_dir='.', model_name='mare2dem', data_file='mare2dem.emdata', target_misfit=1.0, max_iterations=100)#

Create a MARE2DEM model from a regular 2-D resistivity grid.

Port of m2d_gridToM2D.m.

Parameters:

Y (array-like, shape (nz, ny)) – y-coordinates (along-profile) of grid cell centres in metres. Values must vary along columns; all rows for the same column share the same y value.
Z (array-like, shape (nz, ny)) – Depth coordinates of grid cell centres in metres (positive down). Values must vary along rows; all columns for the same row share the same z value.
Rho (array-like, shape (nz, ny)) – Resistivity at each cell centre in Ω·m.
padding_y (float, default 50000.0) – Lateral padding in metres added outside the grid.
padding_z (float, default 50000.0) – Vertical padding in metres added below and above the grid.
out_dir (path-like, default ".") – Output directory.
model_name (str, default "mare2dem") – Stem name used for all output files.
data_file (str, default "mare2dem.emdata") – Name of the associated data file written into the model header.
target_misfit (float, default 1.0) – Target normalized RMS misfit.
max_iterations (int, default 100) – Maximum inversion iterations.

Returns:

Keys "resistivity", "poly", "settings".

Return type:

dict[str, pathlib.Path]

Examples

>>> import numpy as np
>>> from pycsamt.models.mare2dem.grid_to_m2d import grid_to_mare2dem
>>> y1d = np.linspace(-5000, 5000, 21)
>>> z1d = np.linspace(0, 3000, 11)
>>> Y, Z = np.meshgrid(y1d, z1d)
>>> Rho = np.ones_like(Y) * 10.0  # 10 Ω·m half-space
>>> files = grid_to_mare2dem(Y, Z, Rho, out_dir="/tmp/m2d_grid_test")
>>> files["resistivity"].exists()
True

class pycsamt.models.mare2dem.TopoConfig(topo_file='', col_longitude=None, col_latitude=None, col_elevation_m=None, col_depth_m=None, col_distance_km=None, col_distance_m=None)#

Bases: object

Configuration for loading one topography file.

Variables:

topo_file (str or path-like) – Path to the topography file (whitespace-separated columns).
col_longitude (int or None) – 1-based column index for longitude. Required when the file has geographic coordinates.
col_latitude (int or None) – 1-based column index for latitude.
col_elevation_m (int or None) – 1-based column for elevation in metres (positive up). Mutually exclusive with col_depth_m.
col_depth_m (int or None) – 1-based column for depth in metres (positive down).
col_distance_km (int or None) – 1-based column for along-profile distance in km.
col_distance_m (int or None) – 1-based column for along-profile distance in metres.

Parameters:

topo_file (str | Path)
col_longitude (int | None)
col_latitude (int | None)
col_elevation_m (int | None)
col_depth_m (int | None)
col_distance_km (int | None)
col_distance_m (int | None)

topo_file: str | Path = ''#

col_longitude: int | None = None#

col_latitude: int | None = None#

col_elevation_m: int | None = None#

col_depth_m: int | None = None#

col_distance_km: int | None = None#

col_distance_m: int | None = None#

class pycsamt.models.mare2dem.TopoProfile(y_topo, z_topo, northings=None, eastings=None)#

Bases: object

Loaded and projected topography profile.

Variables:

y_topo (numpy.ndarray) – Along-profile position in metres (MARE2DEM y axis).
z_topo (numpy.ndarray) – Depth in metres, positive down (MARE2DEM z axis).
northings (numpy.ndarray or None) – UTM northing (m) — present when converted from Lon/Lat.
eastings (numpy.ndarray or None) – UTM easting (m) — present when converted from Lon/Lat.

Parameters:

y_topo (ndarray)
z_topo (ndarray)
northings (ndarray | None)
eastings (ndarray | None)

y_topo: ndarray#

z_topo: ndarray#

northings: ndarray | None = None#

eastings: ndarray | None = None#

pycsamt.models.mare2dem.import_topo(cfg, *, utm_north0=0.0, utm_east0=0.0, utm_theta=0.0, utm_zone=None, south_hemi=False, orientation_tol=5.0)#

Import a topography file and project it onto the survey profile.

Port of m2d_importTopo.m.

Parameters:

cfg (TopoConfig) – Column layout and file path.
utm_north0 (float) – Profile UTM origin northing (metres).
utm_east0 (float) – Profile UTM origin easting (metres).
utm_theta (float) – Profile UTM theta angle (degrees) — the stUTM.theta field from the .emdata UTM block. Note: the survey-line direction is theta + 90°.
utm_zone (int or None) – UTM zone number. Required only when Lon/Lat input is used.
south_hemi (bool, default False) – Southern-hemisphere flag for UTM conversion.
orientation_tol (float, default 5.0) – Maximum allowed angle (degrees) between the topography profile orientation and the survey line orientation.

Returns:

Projected topography (y, z) in MARE2DEM coordinates.

Return type:

TopoProfile

Raises:

FileNotFoundError – When the topography file does not exist.
ValueError – When the topography and survey orientations differ by more than orientation_tol degrees.

Examples

>>> from pycsamt.models.mare2dem.import_topo import TopoConfig, import_topo
>>> cfg = TopoConfig(
...     topo_file="topo.txt", col_distance_km=1, col_depth_m=2
... )
>>> prof = import_topo(cfg, utm_north0=0.0, utm_east0=0.0)
>>> prof.y_topo
array([0., ...])

pycsamt.models.mare2dem.diff_resistivity(file1, file2, out_file, *, diff_fn=None)#

Difference two MARE2DEM resistivity files and write the result.

Port of diffMARE2DEM_Resistivity.m.

Parameters:

file1 (path-like) – First .resistivity file (minuend, or reference model).
file2 (path-like) – Second .resistivity file (subtrahend, or inverted model).
out_file (path-like) – Destination .resistivity file for the difference model.
diff_fn (callable or None, default None) –
Function (A, B) -> C applied element-wise to the resistivity arrays. Both A and B have shape (n_regions, nrho). The default is:
```
lambda A, B: np.log10(A) - np.log10(B)
```
Custom alternatives:
- Absolute percentage difference: lambda A, B: np.abs((A - B) / A * 100)
- Linear difference: lambda A, B: A - B

Returns:

The difference resistivity model (also written to out_file).

Return type:

ResistivityFile

Raises:

ValueError – When the two files have different numbers of regions.

Examples

Default log10 difference:

>>> from pycsamt.models.mare2dem.diff import diff_resistivity
>>> dm = diff_resistivity(
...     "mare2dem_iter00.resistivity",
...     "mare2dem_iter20.resistivity",
...     "mare2dem_diff.resistivity",
... )
>>> dm.num_regions
4812

Percentage change:

>>> dm = diff_resistivity(
...     "iter00.resistivity",
...     "iter20.resistivity",
...     "pct_change.resistivity",
...     diff_fn=lambda A, B: np.abs((A - B) / A * 100),
... )

class pycsamt.models.mare2dem.EMData(path=None, **kwargs)#

Bases: object

Thin wrapper around EMDataFile for backwards compatibility.

Provides the same path, header, data, and write() interface as the original stub while delegating to the full parser.

Parameters:: path (path-like, optional) – Source file. Read immediately when provided.

property header: dict[str, Any]#

property data: ndarray | None#

property n_data: int#

write(path)#

Write the stored .emdata to path.

Parameters:: path (str | Path)
Return type:: Path

class pycsamt.models.mare2dem.ResistivityModel(path=None, **kwargs)#

Bases: object

Thin wrapper around ResistivityFile for backwards compatibility.

Parameters:: path (path-like, optional) – Source .resistivity file.

property header: dict[str, Any]#

property n_elements: int#

property n_nodes: int#

write(path)#

Write the model to path.

Parameters:: path (str | Path)
Return type:: Path

classmethod halfspace(log10_rho=0.0, *, n_nodes=0)#

Return a homogeneous half-space resistivity model stub.

Parameters:

log10_rho (float)
n_nodes (int)

Return type:

ResistivityModel

class pycsamt.models.mare2dem.PolyMesh(path=None, **kwargs)#

Bases: object

Thin wrapper around PolyFile for backwards compatibility.

Parameters:: path (path-like, optional) – Source .poly file.

property vertices: ndarray | None#

property segments: ndarray | None#

property holes: ndarray | None#

write(path)#

Parameters:: path (str | Path)
Return type:: Path

class pycsamt.models.mare2dem.Mare2DEMLog(path)#

Bases: object

Parse the MARE2DEM per-iteration OccamLog.2012.0 log file.

MARE2DEM writes one block per completed iteration containing Model Misfit, Roughness, and Optimal Mu lines. This parser extracts those values and exposes them as IterationRecord objects.

Parameters:

path (path-like) – Path to the log file (usually *.logfile or *.log).

Variables:

path (pathlib.Path)
iterations (list of IterationRecord)
converged (bool)

property final_rms: float | None#: RMS at the last logged iteration.

property n_iterations: int#: Number of completed iterations in the log.

rms_history()#

Return per-iteration RMS values in order.

Return type:: list[float]

class pycsamt.models.mare2dem.IterationRecord(iteration, rms, roughness, lambda_)#

Bases: object

One completed iteration from a MARE2DEM log file.

Variables:

iteration (int) – Iteration number.
rms (float) – Normalized RMS misfit at this iteration.
roughness (float) – Model roughness.
lambda (float) – Log10 of the Lagrange (regularization) multiplier.

Parameters:

iteration (int)
rms (float)
roughness (float)
lambda_ (float)

iteration: int#

rms: float#

roughness: float#

lambda_: float#

class pycsamt.models.mare2dem.InversionResult(workdir, config=None, **kwargs)#

Bases: Mare2DEMBase

Load and expose MARE2DEM inversion output files.

InversionResult scans a MARE2DEM run directory after the binary has finished and loads the iteration log, final resistivity model, observed-data file, and predicted-response file.

Parameters:

workdir (path-like) – Directory to scan for MARE2DEM output files.
config (Mare2DEMConfig, optional) – Configuration providing default file stems. When omitted, the scanner looks for any .log, .resistivity, .emdata, and *_MARE2DEM.emdata files.
**kwargs – Forwarded to Mare2DEMBase.

Variables:

workdir (pathlib.Path) – Absolute path of the scanned run directory.
config (Mare2DEMConfig) – Configuration used for file-name hints.
log (Mare2DEMLog or None) – Parsed iteration log.
model (ResistivityModel or None) – Final inverted resistivity model.
data (EMData or None) – Observed data file.
response (EMData or None) – Predicted-response file (*_MARE2DEM.emdata).

Examples

>>> from pycsamt.models.mare2dem import InversionResult
>>> result = InversionResult("./mare2dem_run")
>>> result.log.final_rms
0.98
>>> result.log.converged
True

property converged: bool#: True when the log reports successful convergence.

property final_rms: float | None#: Final normalized RMS from the log, or None.

property n_iterations: int#: Number of completed inversion iterations.

summary()#

Return a human-readable summary string.

Return type:: str

print_summary()#

Print summary() to stdout.

Return type:: None

class pycsamt.models.mare2dem.Mare2DEMRunner(workdir, config=None, **kwargs)#

Bases: Mare2DEMBase

Launch MARE2DEM inversion subprocesses.

Mare2DEMRunner is the execution layer of the MARE2DEM wrapper. It receives the stem of a .resistivity file prepared by InputBuilder, selects the configured MARE2DEM executable, and launches the MPI process from workdir. After the run it optionally loads output into an InversionResult.

The command has the logical form:

mpirun -np 8 MARE2DEM mare2dem

where mare2dem is the stem of mare2dem.resistivity.

Parameters:

workdir (path-like, default ".") – Directory that contains, or will receive, a MARE2DEM run. The builder writes the resistivity model, data, and settings files here. The runner executes the MARE2DEM binary from this directory so the input file stem is resolved relative to the run folder. The directory is created before output is written.
config (Mare2DEMConfig, optional) – Configuration object controlling the resistivity model, data component selection, inversion parameters, source management, executable name, and MPI settings. If omitted, a default Mare2DEMConfig is created. Pass an explicit configuration when several MARE2DEM objects must share exactly the same run parameters.
verbose (int or bool, default 0) – Verbosity level. 0 or False keeps the object quiet. Positive values enable diagnostic messages through the instance logger. Larger values may be used to request more detailed run, parsing, or build information.
logger (logging.Logger, optional) – Logger for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created automatically. Provide a logger when integrating MARE2DEM workflows into an application-wide logging configuration.

Variables:

workdir (pathlib.Path) – Directory from which the subprocess is launched.
config (Mare2DEMConfig) – Configuration used for binary name, MPI, and source management.

Notes

Binary resolution is performed by _resolve_binary(), which first checks PATH, then delegates to SourceManager.resolve_binary() for locally compiled binaries.

See also

SourceManager: Download and compile the MARE2DEM binary.
InputBuilder: Write the resistivity model, data, and settings files.
InversionResult: Load MARE2DEM output after the run.

References

[Mare2DEMRunner-1]

Key, K. (2016). MARE2DEM: A 2-D inversion code for controlled-source electromagnetic and magnetotelluric data. Geophysical Journal International, 207(1), 571-588. doi:10.1093/gji/ggw290.

run(resistivity_stem, *, use_mpi=None, n_procs=None, extra_args=None, timeout=None, load_result=True)#

Run a MARE2DEM inversion subprocess.

MARE2DEM receives one positional argument: the stem of the .resistivity file. It derives the data filename by replacing the extension with .emdata and the settings filename with .settings.

Parameters:

resistivity_stem (path-like) – Stem or full path to the .resistivity file. MARE2DEM strips the extension itself; you may pass either "run" or "run.resistivity". Relative paths are interpreted from workdir.
use_mpi (bool, optional) – MPI override. Falls back to config.use_mpi.
n_procs (int, optional) – Number of MPI processes. Falls back to config.n_procs.
extra_args (sequence of str, optional) – Additional command-line arguments appended to the MARE2DEM invocation.
timeout (float, optional) – Maximum run time in seconds. None means no timeout.
load_result (bool, default True) – Whether to scan workdir and return an InversionResult after the run completes.

Returns:

Parsed result object when load_result is True. None otherwise.

Return type:

InversionResult or None

Raises:

FileNotFoundError – When the MARE2DEM binary cannot be located. Build it first with SourceManager.
subprocess.CalledProcessError – When MARE2DEM exits with a non-zero return code.
subprocess.TimeoutExpired – When timeout is set and the process exceeds it.

Examples

Serial run (special single-process build):

>>> from pycsamt.models.mare2dem import Mare2DEMConfig, Mare2DEMRunner
>>> cfg = Mare2DEMConfig(use_mpi=False)
>>> runner = Mare2DEMRunner("./mare2dem_run", config=cfg)
>>> result = runner.run("mare2dem")

MPI run with 8 processes:

>>> cfg = Mare2DEMConfig(use_mpi=True, n_procs=8)
>>> runner = Mare2DEMRunner("./mare2dem_run", config=cfg)
>>> result = runner.run("mare2dem")

command(resistivity_stem, *, use_mpi=None, n_procs=None)#

Return the MARE2DEM command string without executing it.

Parameters:

resistivity_stem (path-like) – Resistivity file stem passed to MARE2DEM.
use_mpi (bool, optional) – MPI override.
n_procs (int, optional) – Process-count override.

Returns:

Shell-quoted command string for display or logging.

Return type:

str

Examples

>>> from pycsamt.models.mare2dem import Mare2DEMConfig, Mare2DEMRunner
>>> cfg = Mare2DEMConfig(use_mpi=True, n_procs=4)
>>> runner = Mare2DEMRunner("./run", config=cfg)
>>> "mpirun" in runner.command("mare2dem")
True

class pycsamt.models.mare2dem.InputBuilder(config=None, **kwargs)#

Bases: Mare2DEMBase

Prepare a MARE2DEM working directory from survey parameters.

InputBuilder produces the three required input files for a MARE2DEM inversion run:

the .emdata observed-data file;
the starting .resistivity model;
the .settings parallel-decomposition control file.

Parameters:

config (Mare2DEMConfig, optional) – Configuration for inversion parameters and file names.
verbose (int or bool, default 0) – Verbosity level. 0 or False keeps the object quiet. Positive values enable diagnostic messages through the instance logger. Larger values may be used to request more detailed run, parsing, or build information.
logger (logging.Logger, optional) – Logger for progress and diagnostic messages. If omitted, a class-specific PyCSAMT logger is created automatically. Provide a logger when integrating MARE2DEM workflows into an application-wide logging configuration.

Examples

Build from an existing .emdata file:

>>> from pycsamt.models.mare2dem import Mare2DEMConfig, InputBuilder
>>> cfg = Mare2DEMConfig(initial_rho=1.0, max_iterations=100)
>>> builder = InputBuilder(config=cfg)
>>> files = builder.build("survey.emdata", workdir="./run")

Build from MTSurveyConfig:

>>> import numpy as np
>>> from pycsamt.models.mare2dem.survey import MTSurveyConfig
>>> mt = MTSurveyConfig(
...     frequencies=np.logspace(-3, 3, 10),
...     rx_y=np.linspace(-5000, 5000, 20),
...     rx_type="marine", lTE=True, lTM=True,
... )
>>> files = builder.build(None, workdir="./run", mt=mt)

See also

Mare2DEMRunner: Launch MARE2DEM on the written input files.
make_data_file: Low-level data file generator.

write_settings(workdir='.', *, filename=None, **sf_kwargs)#

Write the MARE2DEM .settings file.

Parameters:

workdir (path-like, default ".") – Target directory.
filename (str, optional) – Override the settings filename.
**sf_kwargs – Extra keyword arguments forwarded to SettingsFile.

Returns:

Path of the written file.

Return type:

write_resistivity(workdir='.', *, filename=None, poly_file=None)#

Write a homogeneous half-space .resistivity file.

Parameters:

workdir (path-like, default ".") – Target directory.
filename (str, optional) – Override the resistivity filename.
poly_file (str, optional) – Poly file reference written into the model header.

Returns:

Path of the written file.

Return type:

build(source, workdir='.', *, topo=0.0, mt=None, csem=None, data_filename=None, model_filename=None, settings_filename=None)#

Write a MARE2DEM input set to workdir.

Parameters:

source (path-like, EMDataFile, or None) – Existing .emdata file (copied to workdir) or None (generate from mt/csem config objects).
workdir (path-like, default ".") – Target directory.
topo (float or array-like) – Topography for receiver/transmitter placement (used only when source is None).
mt (MTSurveyConfig, optional) – MT survey config (used when source is None).
csem (CSEMSurveyConfig, optional) – CSEM survey config (used when source is None).
data_filename (str, optional) – Override data file name.
model_filename (str, optional) – Override resistivity model file name.
settings_filename (str, optional) – Override settings file name.

Returns:

Keys: "data", "model", "settings".

Return type:

dict[str, pathlib.Path]

class pycsamt.models.mare2dem.PlotConvergence(log_or_result, **kwargs)#

Bases: Mare2DEMBase

Plot RMS misfit convergence from a MARE2DEM log.

Parameters:: log_or_result (Mare2DEMLog or InversionResult) – Source of iteration records.

Examples

>>> from pycsamt.models.mare2dem import InversionResult, PlotConvergence
>>> result = InversionResult("./run")
>>> pc = PlotConvergence(result)
>>> fig = pc.plot()

plot(ax=None, *, savefig=None, dpi=150, target_rms=None)#

Draw the RMS-vs-iteration convergence curve.

Parameters:

ax (matplotlib.axes.Axes | None)
savefig (str | Path | None)
dpi (int)
target_rms (float | None)

Return type:

class pycsamt.models.mare2dem.PlotSurveyLayout(em, **kwargs)#

Bases: Mare2DEMBase

Plot survey receiver and transmitter positions on a map.

Port of plotMARE2DEM_SurveyLayout.m (map view only; GUI removed).

Parameters:

em (EMDataFile) – Data file supplying receiver / transmitter positions and the UTM origin metadata.
**kwargs – Forwarded to Mare2DEMBase.

Examples

>>> from pycsamt.models.mare2dem import read_emdata
>>> from pycsamt.models.mare2dem.plot import PlotSurveyLayout
>>> em = read_emdata("survey.emdata")
>>> fig = PlotSurveyLayout(em).plot()

plot(ax=None, *, savefig=None, dpi=150, units='km')#

Draw the survey map: receivers and transmitters in UTM.

Parameters:

ax (matplotlib.axes.Axes, optional) – Axes to draw on.
savefig (path-like, optional) – Save figure to this path.
dpi (int, default 150) – DPI for the saved figure.
units ({"m", "km"}, default "km") – Display units for axes labels.

Return type:

class pycsamt.models.mare2dem.PlotRxParams(em, **kwargs)#

Bases: Mare2DEMBase

Plot receiver geometry parameters (x, y, z, θ, α, β).

Port of the plotRxParams sub-function in plotMARE2DEM_SurveyLayout.m.

Parameters:

em (EMDataFile) – Survey data file.
**kwargs – Forwarded to Mare2DEMBase.

plot(*, fig=None, savefig=None, dpi=150, units='km')#

Draw 6-panel receiver parameter overview.

Parameters:

savefig (str | Path | None)
dpi (int)
units (str)

Return type:

class pycsamt.models.mare2dem.PlotTxParams(em, **kwargs)#

Bases: Mare2DEMBase

Plot CSEM transmitter geometry parameters (x, y, z, azimuth, dip).

Port of the plotTxParams sub-function in plotMARE2DEM_SurveyLayout.m.

Parameters:: em (Any)

plot(*, fig=None, savefig=None, dpi=150, units='km')#

Draw transmitter parameter overview.

Parameters:

savefig (str | Path | None)
dpi (int)
units (str)

Return type:

pycsamt.models.mare2dem.plot_poly(poly_file, ax=None, *, linewidth=1.0, color='k', savefig=None, dpi=150)#

Plot a Triangle .poly PSLG mesh file.

Port of m2d_plot_poly.m.

Parameters:

poly_file (path-like) – Path to the .poly file.
ax (matplotlib.axes.Axes, optional) – Axes to draw on.
linewidth (float, default 1.0)
color (str, default "k")
savefig (path-like, optional)
dpi (int, default 150)

Returns:

The axes with the PSLG drawn.

Return type:

matplotlib.axes.Axes

Examples

>>> from pycsamt.models.mare2dem.plot import plot_poly
>>> ax = plot_poly("mare2dem.poly")

class pycsamt.models.mare2dem.PlotModel(model_or_result, **kwargs)#

Bases: Mare2DEMBase

Plot the log10-resistivity 2-D section.

Renders a colour-filled triangular mesh via matplotlib.axes.Axes.tripcolor() when Triangle mesh output files (.node + .ele) are found next to the .resistivity file. Falls back to a histogram of region resistivity values when no mesh is present.

Parameters:

model_or_result (ResistivityModel or InversionResult) – Source of resistivity values and mesh location.
**kwargs – Forwarded to Mare2DEMBase.

plot(ax=None, *, savefig=None, dpi=150, cmap='turbo_r', vmin=None, vmax=None)#

Plot the 2-D resistivity section.

When Triangle mesh output files (.node + .ele) are found next to the .resistivity file, a colour-filled triangular section is rendered. Otherwise a histogram of resistivity values is shown.

Parameters:

ax (matplotlib.axes.Axes, optional) – Target axes. When None a new figure is created.
savefig (path-like, optional) – Save the figure to this path.
dpi (int, default 150) – Output resolution in dots per inch.
cmap (str, default "turbo_r") – Colour map for the resistivity section.
vmin (float, optional) – Lower log10-rho bound for the colour axis.
vmax (float, optional) – Upper log10-rho bound for the colour axis.

Return type:

class pycsamt.models.mare2dem.PlotResponse(result, **kwargs)#

Bases: Mare2DEMBase

Compare observed and predicted MT responses.

Generates a per-receiver grid of subplots overlaying TE and TM apparent resistivity and phase versus period.

Parameters:

result (InversionResult) – Inversion output containing observed data and MARE2DEM predicted response.
**kwargs – Forwarded to Mare2DEMBase.

plot(ax=None, *, savefig=None, dpi=150, station=None, max_rx=4, figsize=None)#

Plot observed vs predicted MT responses.

One figure row per receiver, two columns: apparent resistivity on the left and phase on the right.

Parameters:

ax (matplotlib.axes.Axes, optional) – Ignored – the method creates its own figure.
savefig (path-like, optional) – Save the figure to this path.
dpi (int, default 150) – Output resolution in dots per inch.
station (str, optional) – Plot only this receiver by name. When None, the first max_rx receivers are plotted.
max_rx (int, default 4) – Maximum number of receivers when station is None.
figsize ((float, float), optional) – Figure size in inches.

Return type:

Raises:

ValueError – When the result contains no MT data.

2.22.2.3. pycsamt.models.modem#

pycsamt.models.modem — Python interface to the ModEM MT inversion code.

ModEM (Modular EM inversion system) supports 2-D and 3-D magnetotelluric inversion using a nonlinear conjugate-gradient (NLCG) algorithm.

References

Egbert, G.D. & Kelbert, A. (2012). Computational recipes for electromagnetic inverse problems. Geophysical Journal International, 189(1), 251-267.

Kelbert, A., Meqbel, N., Egbert, G.D., Tandon, K. (2014). ModEM: A modular system for inversion of electromagnetic geophysical data. Computers & Geosciences, 66, 40-53.

Quick start — 3D#

>>> from pycsamt.models.modem import ModEmConfig, InputBuilder
>>> cfg = ModEmConfig(mode="3d", initial_rho=100.0)
>>> builder = InputBuilder(config=cfg)
>>> builder.build(edi_source, workdir="./modem_run")
>>> from pycsamt.models.modem import ModEmRunner
>>> ModEmRunner(workdir="./modem_run", config=cfg).run()
>>> from pycsamt.models.modem import InversionResult
>>> result = InversionResult(workdir="./modem_run")
>>> result.plot_model().savefig("model.png")

class pycsamt.models.modem.ModEmConfig(mode='3d', component_type='Full_Impedance', sign_convention='exp(+i\\omega t)', units='[mV/km]/[nT]', error_floor_z=0.05, error_floor_z_floor=0.0, freq_min=None, freq_max=None, nx_2d=100, nz_2d=50, n_airlayers_2d=5, cell_size_h_2d=100.0, cell_size_v_top_2d=10.0, depth_scale_2d=1.2, n_padding_x_2d=7, nx=20, ny=20, nz=30, n_airlayers=5, cell_size_h=500.0, cell_size_v_top=10.0, depth_scale=1.2, n_padding_xy=7, smooth_x=0.1, smooth_y=0.1, smooth_z=0.1, n_smooth_iter=2, qmr_iters_per_divcor=40, max_divcor=20, max_iter_divcor=100, tol_em_fwd=1e-07, tol_em_adj=1e-07, tol_divcor=1e-05, max_iterations=100, target_rms=1.05, initial_lambda=10.0, lambda_divisor=100.0, initial_alpha=10.0, rms_diff_tol=0.0005, lambda_exit=0.0001, initial_rho=100.0, data_file='ModEMData.dat', model_file='ModEM_Model.rho', covariance_file='ModEM.cov', control_file='ModEM.inv', fwd_control_file='ModEM_fwd.ctrl', log_file='Modular_NLCG.log', output_stem='ModEM_out', binary_2d='Mod2DMT', binary_3d='Mod3DMT', use_mpi=False, n_procs=4, mpi_command='mpirun')#

Bases: object

Collect settings that define a ModEM run.

ModEmConfig is the central configuration object for the ModEM v2 subpackage. It is a plain dataclass, so values may be set at construction time, changed before building files, or shared across data, model, covariance, control, runner, and result objects. The class does not perform file I/O or launch the ModEM executable; it records the choices used by those objects.

The configuration covers six workflow areas:

dimensionality and data-component selection;
impedance units, sign convention, and error floors;
2-D and 3-D starting-model geometry;
3-D covariance smoothing and active-cell behavior;
nonlinear inversion-control values;
file names, executable names, and MPI launch settings.

The starting-model classes store resistivity in logarithmic form. For a positive half-space resistivity \(\rho_0\), the initial model value is commonly

\[m_0 = \ln(\rho_0).\]

ModEM then updates this model while balancing data fit and regularization [ModEmConfig-1], [ModEmConfig-2].

2.22. Dimensionality#

mode{“2d”, “3d”}, default “3d”: Dimensionality of the ModEM workflow. "2d" selects two-dimensional input formats and the Mod2DMT binary. "3d" selects three-dimensional formats, writes a covariance file, and uses the Mod3DMT binary by default. The value also controls which model class is built from survey data.

2.22. Data Options#

component_typestr, default “Full_Impedance”: ModEM data component family written to the data file. Typical 2-D choices are "TE_Impedance" and "TM_Impedance". Common 3-D choices include "Full_Impedance", "Off_Diagonal_Impedance", "Determinant_Impedance", and "Full_Vertical_Components". The selected component controls which impedance tensor entries are exported.
sign_conventionstr, default “exp(+i\omega t)”: Time-harmonic sign convention recorded in the ModEM data header. Use "exp(+i\\omega t)" or "exp(-i\\omega t)" to match the convention used by the impedance estimates. A mismatch changes the sign of imaginary components and can lead to inconsistent phase responses.
unitsstr, default “[mV/km]/[nT]”: Impedance units written to the data-file header. The default corresponds to common magnetotelluric field units. Use "[V/m]/[T]" when the impedance tensors are stored in SI units. Unit consistency is important because ModEM interprets data values directly from the file.
error_floor_zfloat, default 0.05: Relative impedance-error floor expressed as a fraction of \(|Z|\). A value of 0.05 enforces a five percent minimum uncertainty on each impedance component. The floor prevents very small formal errors from dominating the objective function and stabilizes inversion weighting.
error_floor_z_floorfloat, default 0.0: Absolute lower bound applied to impedance errors after the relative floor. Use this when some components have very small amplitudes and \(|Z|\)-scaled errors alone would still be too small for a stable inversion.
freq_minfloat, optional: Lower frequency limit in hertz. Frequencies below this value are excluded when building ModEM data from EDI-like sources. Use it to remove low-frequency samples that are noisy, sparsely sampled, or outside the desired depth range.
freq_maxfloat, optional: Upper frequency limit in hertz. Frequencies above this value are excluded when building ModEM data from EDI-like sources. Use it to remove high-frequency samples affected by near-surface noise, instrument limits, or processing artefacts.

2.22. 2-D Grid Options#

nx_2dint, default 100: Number of core horizontal cells in a 2-D ModEM model. The cells describe the inversion region along profile before lateral padding is added. Larger values can represent more lateral structure but increase model size and run time.
nz_2dint, default 50: Number of active earth layers in a 2-D ModEM model. This count excludes air layers. More layers allow finer depth variation but require enough period coverage to constrain the additional parameters.
n_airlayers_2dint, default 5: Number of air layers placed above the earth in 2-D models. Air layers allow the forward solver to represent the air-earth boundary. Their resistivity is normally fixed to a very high value and is not interpreted geologically.
cell_size_h_2dfloat, default 100.0: Nominal horizontal cell width in metres near the 2-D station zone. Smaller values increase near-station resolution and file size. The value should be chosen with station spacing and shortest useful period in mind.
cell_size_v_top_2dfloat, default 10.0: Thickness in metres of the shallowest earth layer in a 2-D model. Subsequent layers grow geometrically according to depth_scale_2d. Choose a value fine enough to represent shallow sensitivity without over-refining the mesh.
depth_scale_2dfloat, default 1.2: Geometric growth factor for 2-D earth-layer thicknesses. Values slightly greater than one create gradually thicker cells with depth. Larger values reach great depths with fewer layers but reduce vertical resolution.
n_padding_x_2dint, default 7: Number of lateral padding cells added to each side of the 2-D station zone. Padding moves artificial boundaries away from the survey line and helps reduce edge effects in forward responses.

2.22. 3-D Grid Options#

nxint, default 20: Number of core cells along the ModEM 3-D x direction, commonly local northing. Padding cells are added outside this core region. Increase the value when station coverage or expected structure requires finer north-south detail.
nyint, default 20: Number of core cells along the ModEM 3-D y direction, commonly local easting. Padding cells are added outside this core region. Increase the value when the survey has dense east-west coverage or strong lateral gradients.
nzint, default 30: Number of active earth layers in a 3-D ModEM model. The count excludes air layers. The depth range and vertical resolution are controlled jointly by cell_size_v_top, depth_scale, and this layer count.
n_airlayersint, default 5: Number of air layers placed above the earth in 3-D models. These layers help represent the air-earth boundary in the forward solver and are usually assigned very high resistivity values.
cell_size_hfloat, default 500.0: Nominal horizontal cell width in metres for 3-D core cells in both x and y directions. It should be consistent with station spacing, expected target size, and the shortest periods retained in the data file.
cell_size_v_topfloat, default 10.0: Thickness in metres of the shallowest 3-D earth layer. Deeper layers grow according to depth_scale. Smaller values improve shallow resolution but increase the number of cells needed to reach the same depth.
depth_scalefloat, default 1.2: Geometric growth factor for 3-D earth-layer thicknesses. The thickness of deeper layers approximately follows \(h_k=h_0 s^k\), where \(h_0\) is the top-layer thickness and \(s\) is this scale factor.
n_padding_xyint, default 7: Number of horizontal padding cells added on each side of the 3-D core grid in both x and y directions. Padding expands the computational domain so boundary conditions are farther from the survey area.

2.22. Covariance Options#

smooth_xfloat, default 0.1: Smoothing coefficient applied between neighbouring cells in the ModEM x direction. Larger values impose stronger model smoothness along this direction. The value should be balanced with data coverage and expected geological strike.
smooth_yfloat, default 0.1: Smoothing coefficient applied between neighbouring cells in the ModEM y direction. For 3-D inversions this controls lateral regularization perpendicular to smooth_x.
smooth_zfloat, default 0.1: Smoothing coefficient applied between neighbouring cells in depth. Increasing this value discourages rapid vertical resistivity changes, while smaller values allow sharper layering when supported by the data.
n_smooth_iterint, default 2: Number of times ModEM applies the covariance smoothing operator. 0 disables repeated smoothing. Higher values strengthen regularization and can produce smoother, more conservative model updates.

2.22. Inversion Control#

max_iterationsint, default 100: Maximum number of nonlinear conjugate-gradient iterations requested in the ModEM control file. The run may stop earlier when the target RMS, lambda-exit threshold, or convergence criteria are reached.
target_rmsfloat, default 1.05: Target normalized root-mean-square misfit. When data errors are realistic, values near \(1\) indicate a fit comparable to the assigned uncertainties. Smaller targets demand tighter data fit and may increase model roughness.
initial_lambdafloat, default 10.0: Initial damping or regularization trade-off parameter used by ModEM. Larger values favour smoother, smaller model updates at the start of inversion. The value is reduced during successful iterations according to lambda_divisor.
lambda_divisorfloat, default 100.0: Factor by which ModEM reduces lambda after successful steps. Larger divisors make lambda decrease faster, while smaller divisors keep stronger regularization for more iterations.
initial_alphafloat, default 10.0: Initial line-search step length used by the inversion control file. It controls the first trial update along the search direction and may be adjusted internally by ModEM during the nonlinear solve.
rms_diff_tolfloat, default 5.0e-4: RMS-change tolerance used as a convergence or restart criterion. When successive RMS values differ by less than this threshold, the inversion is considered to have made little progress.
lambda_exitfloat, default 1.0e-4: Lower lambda threshold used to stop the inversion. Once the regularization trade-off parameter is smaller than this value, additional reductions are unlikely to improve the solution in a stable way.

2.22. Initial Model#

initial_rhofloat, default 100.0: Starting half-space resistivity in ohm metres. Builders fill active earth cells with this value before inversion. A positive value is required because model writers store resistivity in logarithmic form for ModEM-compatible files.

2.22. File Names#

data_filestr, default “ModEMData.dat”: Default observed-data filename used by workflows that consume configuration file names. The path is interpreted relative to the run working directory unless callers provide an absolute path.
model_filestr, default “ModEM_Model.rho”: Default model filename used by configured runner or result workflows. Builders may override this with mode-specific names such as "m0.ws" or "m0.rho".
covariance_filestr, default “ModEM.cov”: Default 3-D covariance filename. It identifies the file containing smoothing coefficients and active-cell masks.
control_filestr, default “ModEM.inv”: Default inversion-control filename. It stores nonlinear solver settings and output naming controls.
log_filestr, default “Modular_NLCG.log”: Default ModEM nonlinear conjugate-gradient log filename. Result loaders use this name when scanning run output.
output_stemstr, default “ModEM_out”: Stem used by ModEM for generated model, response, and log outputs. It corresponds to the model and data output-name field in the control file. Use a distinctive stem when several inversions share one directory.

2.22. Binary And MPI#

binary_2dstr, default “Mod2DMT”: Name or path of the ModEM 2-D executable. If only a name is supplied, the runner searches the working directory and the system PATH. Use an absolute path when the binary is installed in a non-standard location.
binary_3dstr, default “Mod3DMT”: Name or path of the ModEM 3-D executable. If only a name is supplied, the runner searches the working directory and the system PATH. MPI-enabled builds can be launched with use_mpi and n_procs.
use_mpibool, default False: Whether to launch the 3-D executable through an MPI command. This requires a ModEM binary compiled with MPI support. When False, the runner starts the executable directly as a serial process.
n_procsint, default 4: Number of MPI processes requested when use_mpi is true. The value is passed to the MPI launcher and should be compatible with the machine, scheduler allocation, and ModEM build.
mpi_commandstr, default “mpirun”: MPI launcher used when use_mpi is true. Common values include "mpirun" and "mpiexec". Cluster environments may require a site-specific wrapper.

ivar is_3d:: Derived property that returns True when mode is "3d" after whitespace stripping and lower-casing.
vartype is_3d:: bool
ivar binary_name:: Derived property returning binary_3d for 3-D mode and binary_2d otherwise.
vartype binary_name:: str

Notes

The dataclass is intentionally permissive. It records user choices but does not validate every field at construction time. Downstream builders, readers, writers, and runners check the subset of fields they need. This keeps quick configuration experiments lightweight while preserving clear failure points near the operation that requires a value.

For 2-D workflows, choose component_type from TE/TM impedance families and use the *_2d grid fields. For 3-D workflows, use full, off-diagonal, determinant, or vertical component families and the 3-D grid/covariance fields.

2.22. Source-Of-Truth Files#

Users can generate an editable configuration file before building or running a model. Python is the default format because it supports rich inline comments and can still be read safely by from_file() using literal parsing. YAML files also keep comments. JSON files cannot contain comments, so the generated JSON template stores explanations in a "_schema" metadata block and editable values under "config".

The recommended workflow is:

Generate a template with write_template().
Edit the values in the generated file.
Load the edited file with from_file() or read().
Pass the resulting configuration to builders and runners.

The same reusable configuration I/O layer is designed for other model subpackages, including Occam2D.

See also

InputBuilder: Consumes this configuration while writing ModEM input files.
ModEmData.from_edi: Uses data options to select components, units, error floors, and frequency limits.
ModEmModel2D.halfspace: Uses 2-D grid and initial-resistivity settings.
ModEmModel3D.halfspace: Uses 3-D grid and initial-resistivity settings.
ModEmCovariance.from_model: Uses covariance smoothing settings for 3-D runs.
ModEmControl.from_config: Converts inversion-control fields into a ModEM control object.
ModEmRunner: Uses binary, MPI, mode, and file-name settings.

Examples

Create a default 3-D configuration:

>>> from pycsamt.models.modem.config import ModEmConfig
>>> cfg = ModEmConfig()
>>> cfg.is_3d
True
>>> cfg.binary_name
'Mod3DMT'

Configure a 2-D TE workflow:

>>> cfg = ModEmConfig(
...     mode="2d",
...     component_type="TE_Impedance",
...     nx_2d=120,
...     nz_2d=60,
... )
>>> cfg.binary_name
'Mod2DMT'

Set data weighting and frequency limits:

>>> cfg = ModEmConfig(
...     error_floor_z=0.07,
...     error_floor_z_floor=1e-6,
...     freq_min=0.01,
...     freq_max=1000.0,
... )

Configure an MPI-enabled 3-D run:

>>> cfg = ModEmConfig(
...     mode="3d",
...     use_mpi=True,
...     n_procs=16,
...     binary_3d="Mod3DMT_MPI",
... )

Generate a documented source-of-truth file and read it back:

>>> path = ModEmConfig.write_template("modem_config.py")
>>> cfg = ModEmConfig.from_file(path)
>>> cfg.mode
'3d'

Generate a JSON template for environments where JSON is easier to exchange:

>>> ModEmConfig.write_template("modem_config.json")
PosixPath('modem_config.json')

References

[ModEmConfig-1]

Egbert, G. D., and Kelbert, A., “Computational recipes for electromagnetic inverse problems”, Geophysical Journal International, 189(1), 251-267, 2012, doi:10.1111/j.1365-246X.2011.05347.x.

[ModEmConfig-2]

Kelbert, A., Meqbel, N., Egbert, G. D., and Tandon, K., “ModEM: A modular system for inversion of electromagnetic geophysical data”, Computers and Geosciences, 66, 40-53, 2014, doi:10.1016/j.cageo.2014.01.010.

mode: str = '3d'#

component_type: str = 'Full_Impedance'#

sign_convention: str = 'exp(+i\\omega t)'#

units: str = '[mV/km]/[nT]'#

error_floor_z: float = 0.05#

error_floor_z_floor: float = 0.0#

freq_min: float | None = None#

freq_max: float | None = None#

nx_2d: int = 100#

nz_2d: int = 50#

n_airlayers_2d: int = 5#

cell_size_h_2d: float = 100.0#

cell_size_v_top_2d: float = 10.0#

depth_scale_2d: float = 1.2#

n_padding_x_2d: int = 7#

nx: int = 20#

ny: int = 20#

nz: int = 30#

n_airlayers: int = 5#

cell_size_h: float = 500.0#

cell_size_v_top: float = 10.0#

depth_scale: float = 1.2#

n_padding_xy: int = 7#

smooth_x: float = 0.1#

smooth_y: float = 0.1#

smooth_z: float = 0.1#

n_smooth_iter: int = 2#

qmr_iters_per_divcor: int = 40#

max_divcor: int = 20#

max_iter_divcor: int = 100#

tol_em_fwd: float = 1e-07#

tol_em_adj: float = 1e-07#

tol_divcor: float = 1e-05#

max_iterations: int = 100#

target_rms: float = 1.05#

initial_lambda: float = 10.0#

lambda_divisor: float = 100.0#

initial_alpha: float = 10.0#

rms_diff_tol: float = 0.0005#

lambda_exit: float = 0.0001#

initial_rho: float = 100.0#

data_file: str = 'ModEMData.dat'#

model_file: str = 'ModEM_Model.rho'#

covariance_file: str = 'ModEM.cov'#

control_file: str = 'ModEM.inv'#

fwd_control_file: str = 'ModEM_fwd.ctrl'#

log_file: str = 'Modular_NLCG.log'#

output_stem: str = 'ModEM_out'#

binary_2d: str = 'Mod2DMT'#

binary_3d: str = 'Mod3DMT'#

use_mpi: bool = False#

n_procs: int = 4#

mpi_command: str = 'mpirun'#

property is_3d: bool#: Return True when mode selects 3-D ModEM.

property binary_name: str#: Return the executable name implied by mode.

to_template(path='modem_config.py', *, fmt=None)#

Write this configuration as an editable template.

Parameters:

path (path-like, default "modem_config.py") – Destination file. If the path has no suffix, the suffix is inferred from fmt and defaults to .py.
fmt ({"py", "json", "yml", "yaml"}, optional) – Template format. Python and YAML templates include comments. JSON templates include a "_schema" documentation block because standard JSON does not support comments.

Returns:

Path of the generated template.

Return type:

classmethod write_template(path='modem_config.py', *, fmt=None)#

Write a default editable ModEM configuration file.

Parameters:

path (path-like, default "modem_config.py") – Destination file. Suffixes .py, .json, .yml, and .yaml select the output format.
fmt ({"py", "json", "yml", "yaml"}, optional) – Explicit output format. When omitted, the suffix of path is used; paths without a suffix produce a Python template.

Returns:

Path of the generated source-of-truth file.

Return type: