commit c554d2c95bfd8d2b3902f2baac65955e674b4fb4
parent bc1454ab9f9ecc1060fe1377e0d921d92a40cadf
Author: Anders Damsgaard <anders@adamsgaard.dk>
Date: Sat, 27 Jun 2026 08:21:50 +0200
feat: add 3D IDW voxel interpolation of TEM layers
Adds a "3D interpolation…" menu item that interpolates loaded TEM
LineStringZ layers into a styled voxel PointZ GeoPackage layer via the
new numpy-backed tem_loader/interpolate.py core. The dialog derives cell
counts from extent and cell size (read-only, live-updated) with a Z cell
size default of 10 m. Bumps version to 0.3.0 and adds numpy to CI/tests.
Diffstat:
11 files changed, 1369 insertions(+), 82 deletions(-)
diff --git a/.gitignore b/.gitignore
@@ -1,5 +1,10 @@
.artifacts/
+# Python and macOS junk
+__pycache__/
+*.pyc
+.DS_Store
+
# Generated plugin/package outputs
*.gpkg
*.zip
diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
@@ -9,7 +9,7 @@ test:
stage: test
before_script:
- python -m pip install --upgrade pip
- - python -m pip install flake8
+ - python -m pip install flake8 numpy
script:
- make verify
diff --git a/AGENTS.md b/AGENTS.md
@@ -1,6 +1,6 @@
# TEM Loader Plugin
-QGIS 3.40+ plugin for loading Time-Domain Electromagnetic (TEM) inversion XYZ files as styled 3D vector layers. Python, no external dependencies beyond QGIS itself.
+QGIS 3.40+ plugin for loading Time-Domain Electromagnetic (TEM) inversion XYZ files as styled 3D vector layers. Python; the only third-party dependency is numpy, which ships with QGIS. Tests and CI must install numpy.
## Commands
@@ -33,8 +33,9 @@ CI runs `make test` on Python 3.12; the GitLab pipeline also produces `tem_loade
## Architecture
- `tem_loader/core.py` — XYZ parsing, CSV writing, resistivity color mapping. No QGIS dependency; fully unit-testable.
-- `tem_loader/tem_loader.py` — QGIS plugin glue: dialog, layer creation, CRS resolution, styling, group organization.
-- `tem_loader/styles/*.qml` — QGIS layer style files for points, DOI, and layers (supports 3D renderer with dynamic color).
+- `tem_loader/interpolate.py` — Vectorized (numpy) IDW interpolation core: `GridSpec`, `axis_count`/`axis_cell_size`, `samples_from_layers`, `idw`, `interpolate_voxels`. Depends on numpy but has no QGIS dependency; fully unit-testable (numpy must be installed for tests).
+- `tem_loader/tem_loader.py` — QGIS plugin glue: dialogs, layer creation, CRS resolution, styling, group organization.
+- `tem_loader/styles/*.qml` — QGIS layer style files for points, DOI, layers, and voxels (supports 3D renderer with dynamic color).
- `test/data/` — fixture XYZ files for the three supported export families: TEMImage, Aarhus Workbench, and SCI Workbench.
## Conventions
diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -1,5 +1,10 @@
# Changelog
+## [0.3.0] - 2026-06-27
+
+### Added
+- "3D interpolation…" menu item: IDW interpolation of loaded TEM layer features into a voxel PointZ layer written to a GeoPackage. The interpolation dialog lets the user select input LineStringZ layers and optional DOI PointZ layers, define a 3D grid extent with cell sizes and a live derived cell-count preview, and control IDW parameters (power, log10 mode, vertical anisotropy, per-sounding DOI clipping). Output voxels are styled with the resistivity color ramp and rendered as cubes in the 3D Map View.
+
## [0.2.2] - 2026-05-28
### Added
diff --git a/README.md b/README.md
@@ -25,6 +25,8 @@ Each file gets its own layer group. Layers are styled with pre-built QML styles
## Usage
+### Load TEM XYZ files
+
1. Open a QGIS project.
2. Go to **Plugins > TEM Loader > Load TEM XYZ files…**.
3. Select one or more `.xyz` inversion files.
@@ -33,6 +35,20 @@ Each file gets its own layer group. Layers are styled with pre-built QML styles
6. A GeoPackage (`.gpkg`) is written beside each source file with `points`, optional `doi`, and `layers` layers, and the corresponding layers are added to the project with `points` above `doi` above `layers`.
7. Invalid or malformed `.xyz` inputs are rejected with explicit validation errors, and when multiple files are selected the loader continues with the remaining files while reporting the failing input filename.
+### 3D interpolation
+
+After loading one or more TEM datasets, the plugin can interpolate their resistivity model into a regular 3D voxel grid:
+
+1. Load at least one dataset with **Load TEM XYZ files…** so that `layers` (LineStringZ) and optionally `doi` (PointZ) layers are present in the project.
+2. Go to **Plugins > TEM Loader > 3D interpolation…**.
+3. Select the input `layers` from the list (multi-select). Select any `doi` layers to enable per-sounding DOI clipping.
+4. Set the 3D grid extent (X, Y, Z) and cell size. The cell counts are derived from the extent and cell size and shown as a live `nx × ny × nz = N cells` preview that updates as you edit. Use **Use extent of selected input layers** to fill X/Y from the layer bounding box.
+5. Adjust IDW parameters (power, log10 mode, vertical anisotropy, search radius, and min/max neighbors) and optionally disable DOI clipping. Leave the search radius at **Auto** to derive it from the median sounding spacing, or **Max neighbors** at **Unlimited** to weight every sample in range.
+6. Click **OK** and choose an output `.gpkg` path.
+7. The resulting `voxels` PointZ layer is added under an **Interpolation** group, styled with the resistivity color ramp, and visible as colored cubes in the 3D Map View.
+
+The interpolation is pure Python (no external dependencies). For large grids a progress bar is shown and can be cancelled. The grid is capped at 5 000 000 cells to prevent accidental out-of-memory runs; a warning is shown when the cap is exceeded.
+
## XYZ File Format
Input files are whitespace-delimited with optional `/`-prefixed metadata lines. The loader currently supports these export variants:
diff --git a/tem_loader/interpolate.py b/tem_loader/interpolate.py
@@ -0,0 +1,340 @@
+import math
+from collections import namedtuple
+
+import numpy as np
+
+from .core import resistivity_color
+
+VOXEL_SAFETY_CAP = 5_000_000
+
+Sample = namedtuple('Sample', ['x', 'y', 'z', 'resistivity', 'doi_elev'])
+
+
+def axis_count(extent, cell_size):
+ if cell_size <= 0:
+ raise ValueError('cell size must be positive')
+ if extent <= 0:
+ raise ValueError('extent must be positive')
+ return max(1, math.ceil(extent / cell_size))
+
+
+def axis_cell_size(extent, count):
+ if count < 1:
+ raise ValueError('cell count must be >= 1')
+ if extent <= 0:
+ raise ValueError('extent must be positive')
+ return extent / count
+
+
+class GridSpec:
+ def __init__(self, xmin, ymin, zmin, dx, dy, dz, nx, ny, nz):
+ self.xmin = xmin
+ self.ymin = ymin
+ self.zmin = zmin
+ self.dx = dx
+ self.dy = dy
+ self.dz = dz
+ self.nx = nx
+ self.ny = ny
+ self.nz = nz
+
+ @classmethod
+ def from_cell_size(cls, bounds, dx, dy, dz):
+ xmin, xmax, ymin, ymax, zmin, zmax = bounds
+ nx = axis_count(xmax - xmin, dx)
+ ny = axis_count(ymax - ymin, dy)
+ nz = axis_count(zmax - zmin, dz)
+ return cls(xmin, ymin, zmin, dx, dy, dz, nx, ny, nz)
+
+ @classmethod
+ def from_cell_count(cls, bounds, nx, ny, nz):
+ xmin, xmax, ymin, ymax, zmin, zmax = bounds
+ dx = axis_cell_size(xmax - xmin, nx)
+ dy = axis_cell_size(ymax - ymin, ny)
+ dz = axis_cell_size(zmax - zmin, nz)
+ return cls(xmin, ymin, zmin, dx, dy, dz, nx, ny, nz)
+
+ def cell_count(self):
+ return self.nx * self.ny * self.nz
+
+ def centers(self):
+ for iz in range(self.nz):
+ z = self.zmin + (iz + 0.5) * self.dz
+ for iy in range(self.ny):
+ y = self.ymin + (iy + 0.5) * self.dy
+ for ix in range(self.nx):
+ x = self.xmin + (ix + 0.5) * self.dx
+ yield ix, iy, iz, x, y, z
+
+
+def doi_by_xy_from_points(doi_points):
+ result = {}
+ for pt in doi_points:
+ x = pt.get('X') if isinstance(pt, dict) else pt[0]
+ y = pt.get('Y') if isinstance(pt, dict) else pt[1]
+ zdoi = pt.get('ZDOI') if isinstance(pt, dict) else pt[2]
+ if x is not None and y is not None and zdoi is not None:
+ result[(float(x), float(y))] = float(zdoi)
+ return result
+
+
+def samples_from_layers(layer_rows, doi_by_xy=None):
+ samples = []
+ for row in layer_rows:
+ try:
+ x = float(row['X'])
+ y = float(row['Y'])
+ z = float(row['ZMid'])
+ except (KeyError, ValueError, TypeError):
+ continue
+ if math.isnan(x) or math.isnan(y) or math.isnan(z):
+ continue
+ res = row.get('Resistivity')
+ if res is None:
+ continue
+ try:
+ res = float(res)
+ except (ValueError, TypeError):
+ continue
+ if math.isnan(res):
+ continue
+ doi_elev = None
+ if doi_by_xy is not None:
+ # Joined by exact (x, y); both layers come from the same loader
+ # parse, so their coordinates are bit-identical.
+ doi_elev = doi_by_xy.get((x, y))
+ samples.append(Sample(x, y, z, res, doi_elev))
+ return samples
+
+
+class _SampleArrays:
+ """Columnar view of the sample set for vectorized interpolation.
+
+ Missing DOI elevations are stored as NaN so the gate test reduces to a
+ single masked comparison."""
+
+ __slots__ = ('x', 'y', 'z', 'res', 'doi')
+
+ def __init__(self, samples):
+ n = len(samples)
+ self.x = np.fromiter((s.x for s in samples), dtype=float, count=n)
+ self.y = np.fromiter((s.y for s in samples), dtype=float, count=n)
+ self.z = np.fromiter((s.z for s in samples), dtype=float, count=n)
+ self.res = np.fromiter((s.resistivity for s in samples), dtype=float,
+ count=n)
+ self.doi = np.fromiter(
+ (np.nan if s.doi_elev is None else s.doi_elev for s in samples),
+ dtype=float, count=n,
+ )
+
+ def __len__(self):
+ return self.x.shape[0]
+
+
+def _default_radius(samples):
+ """Three times the median nearest-neighbor spacing of the distinct sounding
+ locations. Operating on unique (x, y) positions makes the estimate
+ independent of feature ordering and of how many depth layers each sounding
+ has."""
+ locations = sorted({(s.x, s.y) for s in samples})
+ if len(locations) < 2:
+ return None
+ # The nearest-neighbor scan is O(n**2); cap the working set with an
+ # evenly-spaced subsample so dense surveys stay fast (the subsample slightly
+ # over-estimates spacing, which the ample 3x multiplier absorbs).
+ cap = 500
+ if len(locations) > cap:
+ step = len(locations) / cap
+ locations = [locations[int(i * step)] for i in range(cap)]
+ pts = np.asarray(locations, dtype=float)
+ diff = pts[:, None, :] - pts[None, :, :]
+ dist = np.hypot(diff[..., 0], diff[..., 1])
+ np.fill_diagonal(dist, np.inf)
+ nearest = dist.min(axis=1)
+ nearest = nearest[nearest > 0]
+ if nearest.size == 0:
+ return None
+ return float(np.median(nearest) * 3.0)
+
+
+def _bucket_index(arrays, radius):
+ """Map integer (bx, by) cells of side ``radius`` to the sample indices that
+ fall inside them, so a query only weighs samples in its 3x3 neighborhood."""
+ if radius is None or radius <= 0:
+ return None, None
+ bx = np.floor(arrays.x / radius).astype(np.int64)
+ by = np.floor(arrays.y / radius).astype(np.int64)
+ buckets = {}
+ for i in range(len(arrays)):
+ key = (int(bx[i]), int(by[i]))
+ buckets.setdefault(key, []).append(i)
+ return {k: np.asarray(v) for k, v in buckets.items()}, radius
+
+
+def _candidate_indices(buckets, bucket_size, x, y):
+ if buckets is None:
+ return None
+ bx = int(math.floor(x / bucket_size))
+ by = int(math.floor(y / bucket_size))
+ parts = []
+ for dbx in (-1, 0, 1):
+ for dby in (-1, 0, 1):
+ arr = buckets.get((bx + dbx, by + dby))
+ if arr is not None:
+ parts.append(arr)
+ if not parts:
+ return np.empty(0, dtype=np.int64)
+ return np.concatenate(parts)
+
+
+def _interpolate_at(arrays, idx, qx, qy, qz, *, power, radius, max_neighbors,
+ min_neighbors, anisotropy_z, log10, gate_doi):
+ """Inverse-distance interpolation at a single query point.
+
+ Returns ``(value, num_valid)``. ``value`` is ``None`` when fewer than
+ ``min_neighbors`` usable samples remain; ``num_valid`` is always the count of
+ samples that passed the DOI gate and radius cull. Shared by :func:`idw` and
+ :func:`interpolate_voxels` so both honor identical weighting rules."""
+ if idx is None:
+ sx, sy, sz, sres, sdoi = (arrays.x, arrays.y, arrays.z, arrays.res,
+ arrays.doi)
+ else:
+ sx, sy, sz, sres, sdoi = (arrays.x[idx], arrays.y[idx], arrays.z[idx],
+ arrays.res[idx], arrays.doi[idx])
+
+ dx = sx - qx
+ dy = sy - qy
+ horizontal2 = dx * dx + dy * dy
+
+ keep = np.ones(horizontal2.shape, dtype=bool)
+ if gate_doi:
+ # NaN DOI (no DOI for that sounding) never gates; only drop samples
+ # whose DOI elevation sits above the query depth.
+ keep &= ~((~np.isnan(sdoi)) & (qz < sdoi))
+ if radius is not None:
+ keep &= horizontal2 <= radius * radius
+
+ horizontal2 = horizontal2[keep]
+ sz = sz[keep]
+ res = sres[keep]
+
+ if anisotropy_z != 0:
+ vdz = (sz - qz) / anisotropy_z
+ d2 = horizontal2 + vdz * vdz
+ else:
+ d2 = horizontal2
+
+ if max_neighbors is not None and d2.shape[0] > max_neighbors:
+ nearest = np.argpartition(d2, max_neighbors)[:max_neighbors]
+ d2 = d2[nearest]
+ res = res[nearest]
+
+ num_valid = int(d2.shape[0])
+ if num_valid < min_neighbors:
+ return None, num_valid
+
+ # Coincident samples: average all exact hits rather than picking the first.
+ exact = res[d2 == 0.0]
+ if exact.size:
+ usable = exact[exact > 0] if log10 else exact
+ if usable.size:
+ if log10:
+ return float(10.0 ** np.log10(usable).mean()), num_valid
+ return float(usable.mean()), num_valid
+
+ if log10:
+ # In log space a non-positive resistivity has no defined log; drop it
+ # entirely rather than letting it contribute weight at log10(res)==0.
+ positive = res > 0
+ d2 = d2[positive]
+ res = res[positive]
+
+ if d2.size == 0:
+ return None, num_valid
+
+ weights = 1.0 / (d2 ** (power / 2.0))
+ total_w = float(weights.sum())
+ if total_w == 0.0:
+ return None, num_valid
+
+ values = np.log10(res) if log10 else res
+ avg = float((weights * values).sum() / total_w)
+ return (10.0 ** avg if log10 else avg), num_valid
+
+
+def idw(samples, queries, power=2.0, radius=None, max_neighbors=None,
+ min_neighbors=1, anisotropy_z=1.0, log10=True, gate_doi=True):
+ if radius is None:
+ radius = _default_radius(samples)
+
+ arrays = _SampleArrays(samples)
+ buckets, bucket_size = _bucket_index(arrays, radius)
+
+ results = []
+ for qx, qy, qz in queries:
+ idx = _candidate_indices(buckets, bucket_size, qx, qy)
+ value, _ = _interpolate_at(
+ arrays, idx, qx, qy, qz,
+ power=power, radius=radius, max_neighbors=max_neighbors,
+ min_neighbors=min_neighbors, anisotropy_z=anisotropy_z,
+ log10=log10, gate_doi=gate_doi,
+ )
+ results.append(value)
+ return results
+
+
+def interpolate_voxels(layer_rows, grid, *, doi_points=None,
+ clip_below_doi=True, power=2.0, radius=None,
+ max_neighbors=None, min_neighbors=1, anisotropy_z=1.0,
+ log10=True, progress=None):
+ if grid.cell_count() > VOXEL_SAFETY_CAP:
+ raise ValueError(
+ f'grid has {grid.cell_count()} cells, exceeding the safety cap of '
+ f'{VOXEL_SAFETY_CAP}; reduce the cell count or increase cell size'
+ )
+
+ doi_by_xy = None
+ if doi_points is not None:
+ doi_by_xy = doi_by_xy_from_points(doi_points)
+
+ samples = samples_from_layers(layer_rows, doi_by_xy)
+
+ if radius is None:
+ radius = _default_radius(samples)
+
+ arrays = _SampleArrays(samples)
+ buckets, bucket_size = _bucket_index(arrays, radius)
+
+ voxels = []
+ total = grid.cell_count()
+ done = 0
+
+ for ix, iy, iz, cx, cy, cz in grid.centers():
+ idx = _candidate_indices(buckets, bucket_size, cx, cy)
+ res_val, num_valid = _interpolate_at(
+ arrays, idx, cx, cy, cz,
+ power=power, radius=radius, max_neighbors=max_neighbors,
+ min_neighbors=min_neighbors, anisotropy_z=anisotropy_z,
+ log10=log10, gate_doi=clip_below_doi,
+ )
+ done += 1
+
+ if res_val is not None:
+ voxels.append({
+ 'X': cx,
+ 'Y': cy,
+ 'Z': cz,
+ 'Resistivity': res_val,
+ 'Color': resistivity_color(res_val),
+ 'CellX': ix,
+ 'CellY': iy,
+ 'CellZ': iz,
+ 'NumSamples': num_valid,
+ 'Geometry': f'POINT Z ({cx} {cy} {cz})',
+ })
+
+ # A truthy return from progress requests cancellation.
+ if progress is not None and progress(done, total):
+ break
+
+ return voxels
diff --git a/tem_loader/metadata.txt b/tem_loader/metadata.txt
@@ -3,11 +3,13 @@ name=TEM Loader
qgisMinimumVersion=3.40.0
qgisMaximumVersion=4.99.0
description=Load TEM inversion XYZ files as styled 3D vector layers
-version=0.2.2
+version=0.3.0
author=Anders Damsgaard
email=andam@geus.dk
about=QGIS plugin for loading geophysical inversion models. Adds a menu item to the "Plugin" menu. The models are loaded as 3D points for the terrain surface with line numbers as labels, 3D points for depth of investigation (DOI), and model layers as 3D lines. The plugin supports profile visualization using the "Elevation Profile" tool in QGIS and 3D Map Views.
changelog=
+ 0.3.0
+ * Add "3D interpolation…" menu item: IDW voxel interpolation of loaded TEM layers into a styled PointZ GeoPackage layer
0.2.2
* Hide layers below DOI by default with new import option that clips layer geometry at the DOI
* Exclude DOI layer from Elevation Profile by default
diff --git a/tem_loader/styles/voxels.qml b/tem_loader/styles/voxels.qml
@@ -1,9 +1,9 @@
<!DOCTYPE qgis PUBLIC 'http://mrcc.com/qgis.dtd' 'SYSTEM'>
-<qgis autoRefreshMode="Disabled" autoRefreshTime="0" hasScaleBasedVisibilityFlag="0" labelsEnabled="0" layerType="Vector" maxScale="0" minScale="100000000" readOnly="0" simplifyAlgorithm="0" simplifyDrawingHints="1" simplifyDrawingTol="1" simplifyLocal="1" simplifyMaxScale="1" styleCategories="AllStyleCategories" symbologyReferenceScale="-1" version="4.0.2-Norrköping">
- <renderer-3d layer="" type="vector">
+<qgis autoRefreshMode="Disabled" autoRefreshTime="0" hasScaleBasedVisibilityFlag="0" labelsEnabled="0" layerType="Vector" maxScale="0" minScale="100000000" readOnly="0" simplifyAlgorithm="0" simplifyDrawingHints="0" simplifyDrawingTol="1" simplifyLocal="1" simplifyMaxScale="1" styleCategories="AllStyleCategories" symbologyReferenceScale="-1" version="4.0.3-Norrköping">
+ <renderer-3d layer="voxels_be42d7fc_9810_40a8_972c_a8a80fbaa732" type="vector">
<vector-layer-3d-tiling max-chunk-features="50000" show-bounding-boxes="0" zoom-levels-count="3"/>
- <symbol material_type="phong" shape="box" type="point">
- <data alt-binding="centroid" alt-clamping="absolute"/>
+ <symbol material_type="phong" shape="cylinder" type="point">
+ <data alt-clamping="absolute"/>
<material ambient="26,26,26,255,rgb:0.1000076,0.1000076,0.1000076,1" diffuse="128,128,128,255,rgb:0.5019608,0.5019608,0.5019608,1" ka="1" kd="1" ks="1" opacity="1" shininess="0" specular="0,0,0,255,rgb:0,0,0,1">
<data-defined-properties>
<Option type="Map">
@@ -21,10 +21,51 @@
</material>
<shape-properties>
<Option type="Map">
- <Option name="size" type="double" value="50"/>
+ <Option name="length" type="double" value="10"/>
+ <Option name="model" type="QString" value=""/>
+ <Option name="radius" type="double" value="10"/>
</Option>
</shape-properties>
<transform matrix="1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1"/>
+ <symbol alpha="1" clip_to_extent="1" force_rhr="0" frame_rate="10" is_animated="0" name="symbol" type="marker">
+ <data_defined_properties>
+ <Option type="Map">
+ <Option name="name" type="QString" value=""/>
+ <Option name="properties"/>
+ <Option name="type" type="QString" value="collection"/>
+ </Option>
+ </data_defined_properties>
+ <layer class="SimpleMarker" enabled="1" id="{7f8264d1-d36e-4bdd-b028-1506179f893b}" locked="0" pass="0">
+ <Option type="Map">
+ <Option name="angle" type="QString" value="0"/>
+ <Option name="cap_style" type="QString" value="square"/>
+ <Option name="color" type="QString" value="231,113,72,255,rgb:0.9058824,0.4431373,0.282353,1"/>
+ <Option name="horizontal_anchor_point" type="QString" value="1"/>
+ <Option name="joinstyle" type="QString" value="bevel"/>
+ <Option name="name" type="QString" value="circle"/>
+ <Option name="offset" type="QString" value="0,0"/>
+ <Option name="offset_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="offset_unit" type="QString" value="MM"/>
+ <Option name="outline_color" type="QString" value="35,35,35,255,rgb:0.1372549,0.1372549,0.1372549,1"/>
+ <Option name="outline_style" type="QString" value="solid"/>
+ <Option name="outline_width" type="QString" value="0"/>
+ <Option name="outline_width_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="outline_width_unit" type="QString" value="MM"/>
+ <Option name="scale_method" type="QString" value="diameter"/>
+ <Option name="size" type="QString" value="2"/>
+ <Option name="size_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="size_unit" type="QString" value="MM"/>
+ <Option name="vertical_anchor_point" type="QString" value="1"/>
+ </Option>
+ <data_defined_properties>
+ <Option type="Map">
+ <Option name="name" type="QString" value=""/>
+ <Option name="properties"/>
+ <Option name="type" type="QString" value="collection"/>
+ </Option>
+ </data_defined_properties>
+ </layer>
+ </symbol>
</symbol>
</renderer-3d>
<flags>
@@ -33,7 +74,13 @@
<Searchable>1</Searchable>
<Private>0</Private>
</flags>
- <elevation binding="Centroid" clamping="Absolute" customToleranceEnabled="0" extrusion="0" extrusionEnabled="0" respectLayerSymbol="1" showMarkerSymbolInSurfacePlots="0" symbology="Point" type="IndividualFeatures" zoffset="0" zscale="1">
+ <temporal accumulate="0" durationField="fid" durationUnit="min" enabled="0" endExpression="" endField="" fixedDuration="0" limitMode="0" mode="0" startExpression="" startField="">
+ <fixedRange>
+ <start></start>
+ <end></end>
+ </fixedRange>
+ </temporal>
+ <elevation binding="Centroid" clamping="Absolute" customToleranceEnabled="0" extrusion="0" extrusionEnabled="0" respectLayerSymbol="1" showMarkerSymbolInSurfacePlots="0" symbology="Line" type="IndividualFeatures" zoffset="0" zscale="1">
<data-defined-properties>
<Option type="Map">
<Option name="name" type="QString" value=""/>
@@ -50,12 +97,35 @@
<Option name="type" type="QString" value="collection"/>
</Option>
</data_defined_properties>
- <layer class="SimpleLine" enabled="1" id="" locked="0" pass="0">
+ <layer class="SimpleLine" enabled="1" id="{94ddf49e-8495-49e5-917a-b395b2af0384}" locked="0" pass="0">
<Option type="Map">
+ <Option name="align_dash_pattern" type="QString" value="0"/>
+ <Option name="capstyle" type="QString" value="square"/>
+ <Option name="customdash" type="QString" value="5;2"/>
+ <Option name="customdash_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="customdash_unit" type="QString" value="MM"/>
+ <Option name="dash_pattern_offset" type="QString" value="0"/>
+ <Option name="dash_pattern_offset_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="dash_pattern_offset_unit" type="QString" value="MM"/>
+ <Option name="draw_inside_polygon" type="QString" value="0"/>
+ <Option name="joinstyle" type="QString" value="bevel"/>
<Option name="line_color" type="QString" value="114,155,111,255,rgb:0.4470588,0.6078432,0.4352941,1"/>
<Option name="line_style" type="QString" value="solid"/>
<Option name="line_width" type="QString" value="0.6"/>
<Option name="line_width_unit" type="QString" value="MM"/>
+ <Option name="offset" type="QString" value="0"/>
+ <Option name="offset_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="offset_unit" type="QString" value="MM"/>
+ <Option name="ring_filter" type="QString" value="0"/>
+ <Option name="trim_distance_end" type="QString" value="0"/>
+ <Option name="trim_distance_end_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="trim_distance_end_unit" type="QString" value="MM"/>
+ <Option name="trim_distance_start" type="QString" value="0"/>
+ <Option name="trim_distance_start_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="trim_distance_start_unit" type="QString" value="MM"/>
+ <Option name="tweak_dash_pattern_on_corners" type="QString" value="0"/>
+ <Option name="use_custom_dash" type="QString" value="0"/>
+ <Option name="width_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
</Option>
<data_defined_properties>
<Option type="Map">
@@ -76,9 +146,18 @@
<Option name="type" type="QString" value="collection"/>
</Option>
</data_defined_properties>
- <layer class="SimpleFill" enabled="1" id="" locked="0" pass="0">
+ <layer class="SimpleFill" enabled="1" id="{061c705d-9984-4474-be8a-6e0e9bc68bce}" locked="0" pass="0">
<Option type="Map">
+ <Option name="border_width_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
<Option name="color" type="QString" value="114,155,111,255,rgb:0.4470588,0.6078432,0.4352941,1"/>
+ <Option name="joinstyle" type="QString" value="bevel"/>
+ <Option name="offset" type="QString" value="0,0"/>
+ <Option name="offset_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="offset_unit" type="QString" value="MM"/>
+ <Option name="outline_color" type="QString" value="35,35,35,255,rgb:0.1372549,0.1372549,0.1372549,1"/>
+ <Option name="outline_style" type="QString" value="solid"/>
+ <Option name="outline_width" type="QString" value="0.26"/>
+ <Option name="outline_width_unit" type="QString" value="MM"/>
<Option name="style" type="QString" value="solid"/>
</Option>
<data_defined_properties>
@@ -100,12 +179,27 @@
<Option name="type" type="QString" value="collection"/>
</Option>
</data_defined_properties>
- <layer class="SimpleMarker" enabled="1" id="" locked="0" pass="0">
+ <layer class="SimpleMarker" enabled="1" id="{f39c2212-b7bb-4795-b5c2-7b37ec8a5d35}" locked="0" pass="0">
<Option type="Map">
+ <Option name="angle" type="QString" value="0"/>
+ <Option name="cap_style" type="QString" value="square"/>
<Option name="color" type="QString" value="114,155,111,255,rgb:0.4470588,0.6078432,0.4352941,1"/>
+ <Option name="horizontal_anchor_point" type="QString" value="1"/>
+ <Option name="joinstyle" type="QString" value="bevel"/>
<Option name="name" type="QString" value="square"/>
+ <Option name="offset" type="QString" value="0,0"/>
+ <Option name="offset_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="offset_unit" type="QString" value="MM"/>
+ <Option name="outline_color" type="QString" value="35,35,35,255,rgb:0.1372549,0.1372549,0.1372549,1"/>
+ <Option name="outline_style" type="QString" value="solid"/>
+ <Option name="outline_width" type="QString" value="0"/>
+ <Option name="outline_width_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="outline_width_unit" type="QString" value="MM"/>
+ <Option name="scale_method" type="QString" value="diameter"/>
<Option name="size" type="QString" value="3"/>
+ <Option name="size_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
<Option name="size_unit" type="QString" value="MM"/>
+ <Option name="vertical_anchor_point" type="QString" value="1"/>
</Option>
<data_defined_properties>
<Option type="Map">
@@ -118,18 +212,15 @@
</symbol>
</profileMarkerSymbol>
</elevation>
- <renderer-v2 attr="Color" enableorderby="0" forceraster="0" referencescale="-1" symbollevels="0" type="categorizedSymbol">
- <categories>
- <category label="Resistivity color" render="true" symbol="0" uuid="{voxel-color}" value=""/>
- </categories>
+ <renderer-v2 enableorderby="0" forceraster="0" referencescale="-1" symbollevels="0" type="singleSymbol">
<symbols>
<symbol alpha="1" clip_to_extent="1" force_rhr="0" frame_rate="10" is_animated="0" name="0" type="marker">
<data_defined_properties>
<Option type="Map">
<Option name="name" type="QString" value=""/>
<Option name="properties" type="Map">
- <Option name="color" type="Map">
- <Option name="active" type="bool" value="true"/>
+ <Option name="alpha" type="Map">
+ <Option name="active" type="bool" value="false"/>
<Option name="field" type="QString" value="Color"/>
<Option name="type" type="int" value="2"/>
</Option>
@@ -137,17 +228,38 @@
<Option name="type" type="QString" value="collection"/>
</Option>
</data_defined_properties>
- <layer class="SimpleMarker" enabled="1" id="" locked="0" pass="0">
+ <layer class="SimpleMarker" enabled="1" id="{7d4774e8-68fa-4179-ba02-6f537e0b1400}" locked="0" pass="0">
<Option type="Map">
+ <Option name="angle" type="QString" value="0"/>
+ <Option name="cap_style" type="QString" value="square"/>
<Option name="color" type="QString" value="0,0,0,255,rgb:0,0,0,1"/>
+ <Option name="horizontal_anchor_point" type="QString" value="1"/>
+ <Option name="joinstyle" type="QString" value="bevel"/>
<Option name="name" type="QString" value="square"/>
- <Option name="size" type="QString" value="2"/>
+ <Option name="offset" type="QString" value="0,0"/>
+ <Option name="offset_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="offset_unit" type="QString" value="MM"/>
+ <Option name="outline_color" type="QString" value="35,35,35,255,rgb:0.1372549,0.1372549,0.1372549,1"/>
+ <Option name="outline_style" type="QString" value="no"/>
+ <Option name="outline_width" type="QString" value="0"/>
+ <Option name="outline_width_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="outline_width_unit" type="QString" value="MM"/>
+ <Option name="scale_method" type="QString" value="diameter"/>
+ <Option name="size" type="QString" value="3"/>
+ <Option name="size_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
<Option name="size_unit" type="QString" value="MM"/>
+ <Option name="vertical_anchor_point" type="QString" value="1"/>
</Option>
<data_defined_properties>
<Option type="Map">
<Option name="name" type="QString" value=""/>
- <Option name="properties"/>
+ <Option name="properties" type="Map">
+ <Option name="fillColor" type="Map">
+ <Option name="active" type="bool" value="true"/>
+ <Option name="field" type="QString" value="Color"/>
+ <Option name="type" type="int" value="2"/>
+ </Option>
+ </Option>
<Option name="type" type="QString" value="collection"/>
</Option>
</data_defined_properties>
@@ -164,4 +276,273 @@
</Option>
</data-defined-properties>
</renderer-v2>
+ <selection mode="Default">
+ <selectionColor invalid="1"/>
+ <selectionSymbol>
+ <symbol alpha="1" clip_to_extent="1" force_rhr="0" frame_rate="10" is_animated="0" name="" type="marker">
+ <data_defined_properties>
+ <Option type="Map">
+ <Option name="name" type="QString" value=""/>
+ <Option name="properties"/>
+ <Option name="type" type="QString" value="collection"/>
+ </Option>
+ </data_defined_properties>
+ <layer class="SimpleMarker" enabled="1" id="{69d96a7c-eb66-4934-9ce3-de72887750eb}" locked="0" pass="0">
+ <Option type="Map">
+ <Option name="angle" type="QString" value="0"/>
+ <Option name="cap_style" type="QString" value="square"/>
+ <Option name="color" type="QString" value="255,0,0,255,rgb:1,0,0,1"/>
+ <Option name="horizontal_anchor_point" type="QString" value="1"/>
+ <Option name="joinstyle" type="QString" value="bevel"/>
+ <Option name="name" type="QString" value="circle"/>
+ <Option name="offset" type="QString" value="0,0"/>
+ <Option name="offset_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="offset_unit" type="QString" value="MM"/>
+ <Option name="outline_color" type="QString" value="35,35,35,255,rgb:0.1372549,0.1372549,0.1372549,1"/>
+ <Option name="outline_style" type="QString" value="solid"/>
+ <Option name="outline_width" type="QString" value="0"/>
+ <Option name="outline_width_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="outline_width_unit" type="QString" value="MM"/>
+ <Option name="scale_method" type="QString" value="diameter"/>
+ <Option name="size" type="QString" value="2"/>
+ <Option name="size_map_unit_scale" type="QString" value="3x:0,0,0,0,0,0"/>
+ <Option name="size_unit" type="QString" value="MM"/>
+ <Option name="vertical_anchor_point" type="QString" value="1"/>
+ </Option>
+ <data_defined_properties>
+ <Option type="Map">
+ <Option name="name" type="QString" value=""/>
+ <Option name="properties"/>
+ <Option name="type" type="QString" value="collection"/>
+ </Option>
+ </data_defined_properties>
+ </layer>
+ </symbol>
+ </selectionSymbol>
+ </selection>
+ <customproperties>
+ <Option type="Map">
+ <Option name="embeddedWidgets/count" type="int" value="0"/>
+ <Option name="variableNames"/>
+ <Option name="variableValues"/>
+ </Option>
+ </customproperties>
+ <blendMode>0</blendMode>
+ <featureBlendMode>0</featureBlendMode>
+ <layerOpacity>1</layerOpacity>
+ <geometryOptions geometryPrecision="0" removeDuplicateNodes="0">
+ <activeChecks/>
+ <checkConfiguration/>
+ </geometryOptions>
+ <legend showLabelLegend="0" type="default-vector"/>
+ <referencedLayers/>
+ <referencingLayers/>
+ <fieldConfiguration>
+ <field configurationFlags="NoFlag" name="fid">
+ <editWidget type="TextEdit">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="X">
+ <editWidget type="TextEdit">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="Y">
+ <editWidget type="TextEdit">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="Z">
+ <editWidget type="TextEdit">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="Resistivity">
+ <editWidget type="TextEdit">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="Color">
+ <editWidget type="TextEdit">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="CellX">
+ <editWidget type="Range">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="CellY">
+ <editWidget type="Range">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="CellZ">
+ <editWidget type="Range">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ <field configurationFlags="NoFlag" name="NumSamples">
+ <editWidget type="Range">
+ <config>
+ <Option/>
+ </config>
+ </editWidget>
+ </field>
+ </fieldConfiguration>
+ <aliases>
+ <alias field="fid" index="0" name=""/>
+ <alias field="X" index="1" name=""/>
+ <alias field="Y" index="2" name=""/>
+ <alias field="Z" index="3" name=""/>
+ <alias field="Resistivity" index="4" name=""/>
+ <alias field="Color" index="5" name=""/>
+ <alias field="CellX" index="6" name=""/>
+ <alias field="CellY" index="7" name=""/>
+ <alias field="CellZ" index="8" name=""/>
+ <alias field="NumSamples" index="9" name=""/>
+ </aliases>
+ <defaults>
+ <default applyOnUpdate="0" expression="" field="fid"/>
+ <default applyOnUpdate="0" expression="" field="X"/>
+ <default applyOnUpdate="0" expression="" field="Y"/>
+ <default applyOnUpdate="0" expression="" field="Z"/>
+ <default applyOnUpdate="0" expression="" field="Resistivity"/>
+ <default applyOnUpdate="0" expression="" field="Color"/>
+ <default applyOnUpdate="0" expression="" field="CellX"/>
+ <default applyOnUpdate="0" expression="" field="CellY"/>
+ <default applyOnUpdate="0" expression="" field="CellZ"/>
+ <default applyOnUpdate="0" expression="" field="NumSamples"/>
+ </defaults>
+ <constraints>
+ <constraint constraints="3" exp_strength="0" field="fid" notnull_strength="1" unique_strength="1"/>
+ <constraint constraints="0" exp_strength="0" field="X" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="Y" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="Z" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="Resistivity" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="Color" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="CellX" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="CellY" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="CellZ" notnull_strength="0" unique_strength="0"/>
+ <constraint constraints="0" exp_strength="0" field="NumSamples" notnull_strength="0" unique_strength="0"/>
+ </constraints>
+ <constraintExpressions>
+ <constraint desc="" exp="" field="fid"/>
+ <constraint desc="" exp="" field="X"/>
+ <constraint desc="" exp="" field="Y"/>
+ <constraint desc="" exp="" field="Z"/>
+ <constraint desc="" exp="" field="Resistivity"/>
+ <constraint desc="" exp="" field="Color"/>
+ <constraint desc="" exp="" field="CellX"/>
+ <constraint desc="" exp="" field="CellY"/>
+ <constraint desc="" exp="" field="CellZ"/>
+ <constraint desc="" exp="" field="NumSamples"/>
+ </constraintExpressions>
+ <expressionfields/>
+ <attributeactions>
+ <defaultAction key="Canvas" value="{00000000-0000-0000-0000-000000000000}"/>
+ </attributeactions>
+ <attributetableconfig actionWidgetStyle="dropDown" sortExpression="" sortOrder="0">
+ <columns>
+ <column hidden="0" name="fid" type="field" width="-1"/>
+ <column hidden="0" name="X" type="field" width="-1"/>
+ <column hidden="0" name="Y" type="field" width="-1"/>
+ <column hidden="0" name="Z" type="field" width="-1"/>
+ <column hidden="0" name="Resistivity" type="field" width="-1"/>
+ <column hidden="0" name="Color" type="field" width="-1"/>
+ <column hidden="0" name="CellX" type="field" width="-1"/>
+ <column hidden="0" name="CellY" type="field" width="-1"/>
+ <column hidden="0" name="CellZ" type="field" width="-1"/>
+ <column hidden="0" name="NumSamples" type="field" width="-1"/>
+ <column hidden="1" type="actions" width="-1"/>
+ </columns>
+ </attributetableconfig>
+ <conditionalstyles>
+ <rowstyles/>
+ <fieldstyles/>
+ </conditionalstyles>
+ <storedexpressions/>
+ <editform tolerant="1"></editform>
+ <editforminit/>
+ <editforminitcodesource>0</editforminitcodesource>
+ <editforminitfilepath></editforminitfilepath>
+ <editforminitcode><![CDATA[# -*- coding: utf-8 -*-
+"""
+QGIS forms can have a Python function that is called when the form is
+opened.
+
+Use this function to add extra logic to your forms.
+
+Enter the name of the function in the "Python Init function"
+field.
+An example follows:
+"""
+from qgis.PyQt.QtWidgets import QWidget
+
+def my_form_open(dialog, layer, feature):
+ geom = feature.geometry()
+ control = dialog.findChild(QWidget, "MyLineEdit")
+]]></editforminitcode>
+ <featformsuppress>0</featformsuppress>
+ <editorlayout>generatedlayout</editorlayout>
+ <editable>
+ <field editable="1" name="CellX"/>
+ <field editable="1" name="CellY"/>
+ <field editable="1" name="CellZ"/>
+ <field editable="1" name="Color"/>
+ <field editable="1" name="NumSamples"/>
+ <field editable="1" name="Resistivity"/>
+ <field editable="1" name="X"/>
+ <field editable="1" name="Y"/>
+ <field editable="1" name="Z"/>
+ <field editable="1" name="fid"/>
+ </editable>
+ <labelOnTop>
+ <field labelOnTop="0" name="CellX"/>
+ <field labelOnTop="0" name="CellY"/>
+ <field labelOnTop="0" name="CellZ"/>
+ <field labelOnTop="0" name="Color"/>
+ <field labelOnTop="0" name="NumSamples"/>
+ <field labelOnTop="0" name="Resistivity"/>
+ <field labelOnTop="0" name="X"/>
+ <field labelOnTop="0" name="Y"/>
+ <field labelOnTop="0" name="Z"/>
+ <field labelOnTop="0" name="fid"/>
+ </labelOnTop>
+ <reuseLastValuePolicy>
+ <field name="CellX" reuseLastValuePolicy="NotAllowed"/>
+ <field name="CellY" reuseLastValuePolicy="NotAllowed"/>
+ <field name="CellZ" reuseLastValuePolicy="NotAllowed"/>
+ <field name="Color" reuseLastValuePolicy="NotAllowed"/>
+ <field name="NumSamples" reuseLastValuePolicy="NotAllowed"/>
+ <field name="Resistivity" reuseLastValuePolicy="NotAllowed"/>
+ <field name="X" reuseLastValuePolicy="NotAllowed"/>
+ <field name="Y" reuseLastValuePolicy="NotAllowed"/>
+ <field name="Z" reuseLastValuePolicy="NotAllowed"/>
+ <field name="fid" reuseLastValuePolicy="NotAllowed"/>
+ </reuseLastValuePolicy>
+ <dataDefinedFieldProperties/>
+ <widgets/>
+ <previewExpression>"Color"</previewExpression>
+ <mapTip enabled="1"></mapTip>
+ <layerGeometryType>0</layerGeometryType>
</qgis>
diff --git a/tem_loader/tem_loader.py b/tem_loader/tem_loader.py
@@ -361,16 +361,16 @@ class _InterpolationDialog(QDialog):
self._dx = QDoubleSpinBox()
self._dy = QDoubleSpinBox()
self._dz = QDoubleSpinBox()
- self._nx = QSpinBox()
- self._ny = QSpinBox()
- self._nz = QSpinBox()
for sb in (self._dx, self._dy, self._dz):
sb.setRange(0.001, 1e6)
sb.setDecimals(3)
sb.setValue(100.0)
- for sb in (self._nx, self._ny, self._nz):
- sb.setRange(1, 10000)
- sb.setValue(10)
+ self._dz.setValue(10.0)
+ # Cell counts are derived from extent and cell size, not editable.
+ self._nx = QLabel('—')
+ self._ny = QLabel('—')
+ self._nz = QLabel('—')
+ self._count_x = self._count_y = self._count_z = None
self._summary_label = QLabel()
@@ -408,21 +408,13 @@ class _InterpolationDialog(QDialog):
use_extent_btn = QPushButton('Use extent of selected input layers')
use_extent_btn.clicked.connect(self._fill_extent_from_inputs)
- self._dx.valueChanged.connect(self._on_cell_size_changed)
- self._dy.valueChanged.connect(self._on_cell_size_changed)
- self._dz.valueChanged.connect(self._on_cell_size_changed)
- self._nx.valueChanged.connect(self._on_cell_count_changed)
- self._ny.valueChanged.connect(self._on_cell_count_changed)
- self._nz.valueChanged.connect(self._on_cell_count_changed)
- self._xmin.valueChanged.connect(self._update_summary)
- self._xmax.valueChanged.connect(self._update_summary)
- self._ymin.valueChanged.connect(self._update_summary)
- self._ymax.valueChanged.connect(self._update_summary)
- self._zmin.valueChanged.connect(self._update_summary)
- self._zmax.valueChanged.connect(self._update_summary)
-
self._updating = False
+ for sb in (self._dx, self._dy, self._dz,
+ self._xmin, self._xmax, self._ymin, self._ymax,
+ self._zmin, self._zmax):
+ sb.valueChanged.connect(self._recompute_counts)
+
buttons = QDialogButtonBox(_dialog_buttons())
buttons.accepted.connect(self.accept)
buttons.rejected.connect(self.reject)
@@ -465,7 +457,7 @@ class _InterpolationDialog(QDialog):
layout.addWidget(buttons)
self.setLayout(layout)
- self._update_summary()
+ self._recompute_counts()
def _fill_extent_from_inputs(self):
layers = self._selected_input_layers()
@@ -493,7 +485,7 @@ class _InterpolationDialog(QDialog):
self._zmin.setValue(min(zs))
self._zmax.setValue(max(zs))
self._updating = False
- self._update_summary()
+ self._recompute_counts()
def _selected_input_layers(self):
layers = []
@@ -515,48 +507,34 @@ class _InterpolationDialog(QDialog):
layers.append(lyr)
return layers
- def _on_cell_size_changed(self):
- if self._updating:
- return
- self._updating = True
+ @staticmethod
+ def _axis_count(extent, cell_size):
try:
- xext = self._xmax.value() - self._xmin.value()
- yext = self._ymax.value() - self._ymin.value()
- zext = self._zmax.value() - self._zmin.value()
- if xext > 0:
- self._nx.setValue(interp.axis_count(xext, self._dx.value()))
- if yext > 0:
- self._ny.setValue(interp.axis_count(yext, self._dy.value()))
- if zext > 0:
- self._nz.setValue(interp.axis_count(zext, self._dz.value()))
+ return interp.axis_count(extent, cell_size)
except (ValueError, ZeroDivisionError):
- pass
- finally:
- self._updating = False
- self._update_summary()
+ return None
- def _on_cell_count_changed(self):
+ def _recompute_counts(self):
if self._updating:
return
- self._updating = True
- try:
- xext = self._xmax.value() - self._xmin.value()
- yext = self._ymax.value() - self._ymin.value()
- zext = self._zmax.value() - self._zmin.value()
- if xext > 0:
- self._dx.setValue(interp.axis_cell_size(xext, self._nx.value()))
- if yext > 0:
- self._dy.setValue(interp.axis_cell_size(yext, self._ny.value()))
- if zext > 0:
- self._dz.setValue(interp.axis_cell_size(zext, self._nz.value()))
- except (ValueError, ZeroDivisionError):
- pass
- finally:
- self._updating = False
+ self._count_x = self._axis_count(
+ self._xmax.value() - self._xmin.value(), self._dx.value())
+ self._count_y = self._axis_count(
+ self._ymax.value() - self._ymin.value(), self._dy.value())
+ self._count_z = self._axis_count(
+ self._zmax.value() - self._zmin.value(), self._dz.value())
+ self._nx.setText('—' if self._count_x is None else str(self._count_x))
+ self._ny.setText('—' if self._count_y is None else str(self._count_y))
+ self._nz.setText('—' if self._count_z is None else str(self._count_z))
self._update_summary()
def _update_summary(self):
- nx, ny, nz = self._nx.value(), self._ny.value(), self._nz.value()
+ if None in (self._count_x, self._count_y, self._count_z):
+ self._summary_label.setText(
+ 'Set a non-empty extent to compute the grid.'
+ )
+ return
+ nx, ny, nz = self._count_x, self._count_y, self._count_z
total = nx * ny * nz
self._summary_label.setText(
f'{nx} × {ny} × {nz} = {total:,} cells'
@@ -568,11 +546,11 @@ class _InterpolationDialog(QDialog):
self._ymin.value(), self._ymax.value(),
self._zmin.value(), self._zmax.value(),
)
- grid = interp.GridSpec.from_cell_count(
+ grid = interp.GridSpec.from_cell_size(
bounds,
- self._nx.value(),
- self._ny.value(),
- self._nz.value(),
+ self._dx.value(),
+ self._dy.value(),
+ self._dz.value(),
)
radius = self._radius.value()
max_neighbors = self._max_neighbors.value()
diff --git a/test/test_core.py b/test/test_core.py
@@ -678,6 +678,123 @@ class PluginTests(unittest.TestCase):
def addLayout(self, layout):
self.items.append(layout)
+ class FakeHBoxLayout:
+ def __init__(self):
+ self.items = []
+
+ def addWidget(self, widget):
+ self.items.append(widget)
+
+ def addLayout(self, layout):
+ self.items.append(layout)
+
+ class FakeDoubleSpinBox:
+ def __init__(self):
+ self.minimum = None
+ self.maximum = None
+ self.decimals = None
+ self._value = 0.0
+ self.enabled = True
+ self.valueChanged = FakeSignal()
+
+ def setRange(self, minimum, maximum):
+ self.minimum = minimum
+ self.maximum = maximum
+
+ def setDecimals(self, d):
+ self.decimals = d
+
+ def setValue(self, value):
+ old = self._value
+ self._value = float(value)
+ if old != self._value:
+ self.valueChanged.emit(self._value)
+
+ def setEnabled(self, enabled):
+ self.enabled = bool(enabled)
+
+ def value(self):
+ return self._value
+
+ class FakeLabel:
+ def __init__(self, text=''):
+ self._text = text
+
+ def setText(self, text):
+ self._text = text
+
+ def text(self):
+ return self._text
+
+ class FakeListWidgetItem:
+ class _RoleData:
+ pass
+
+ def __init__(self, text=''):
+ self._text = text
+ self._selected = False
+ self._data = {}
+
+ def setSelected(self, selected):
+ self._selected = bool(selected)
+
+ def isSelected(self):
+ return self._selected
+
+ def setData(self, role, value):
+ self._data[role] = value
+
+ def data(self, role):
+ return self._data.get(role)
+
+ def text(self):
+ return self._text
+
+ class FakeListWidget:
+ MultiSelection = 'MultiSelection'
+
+ class SelectionMode:
+ MultiSelection = 'MultiSelection'
+
+ def __init__(self):
+ self._items = []
+
+ def setSelectionMode(self, mode):
+ pass
+
+ def addItem(self, item):
+ self._items.append(item)
+
+ def count(self):
+ return len(self._items)
+
+ def item(self, index):
+ return self._items[index]
+
+ class FakeProgressDialog:
+ def __init__(self, *args, **kwargs):
+ self._value = 0
+
+ def setWindowTitle(self, title):
+ pass
+
+ def setMinimumDuration(self, ms):
+ pass
+
+ def setValue(self, value):
+ self._value = value
+
+ def value(self):
+ return self._value
+
+ def wasCanceled(self):
+ return False
+
+ class FakePushButton:
+ def __init__(self, text=''):
+ self.text = text
+ self.clicked = FakeSignal()
+
class FakeQMetaType:
class Type:
QString = "QString"
@@ -692,9 +809,16 @@ class PluginTests(unittest.TestCase):
qtwidgets.QCheckBox = FakeCheckBox
qtwidgets.QDialog = FakeDialog
qtwidgets.QDialogButtonBox = FakeDialogButtonBox
+ qtwidgets.QDoubleSpinBox = FakeDoubleSpinBox
qtwidgets.QFileDialog = FakeFileDialog
qtwidgets.QFormLayout = FakeFormLayout
+ qtwidgets.QHBoxLayout = FakeHBoxLayout
+ qtwidgets.QLabel = FakeLabel
+ qtwidgets.QListWidget = FakeListWidget
+ qtwidgets.QListWidgetItem = FakeListWidgetItem
qtwidgets.QMessageBox = FakeMessageBox
+ qtwidgets.QProgressDialog = FakeProgressDialog
+ qtwidgets.QPushButton = FakePushButton
qtwidgets.QSpinBox = FakeSpinBox
qtwidgets.QVBoxLayout = FakeVBoxLayout
@@ -2028,3 +2152,46 @@ class PluginTests(unittest.TestCase):
"model.xyz: failed to load layers: doi",
message_box.warnings[0][2],
)
+
+ def test_voxels_3d_renderer_uses_color_field(self):
+ import xml.etree.ElementTree as ET
+ tree = ET.parse(STYLES_DIR / "voxels.qml")
+ renderer = tree.getroot().find("./renderer-3d")
+ self.assertIsNotNone(renderer)
+ ambient = renderer.find(
+ ".//material/data-defined-properties/Option/Option[@name='properties']"
+ "/Option[@name='ambient']"
+ )
+ self.assertIsNotNone(ambient)
+ active = ambient.find("./Option[@name='active']")
+ expr = ambient.find("./Option[@name='expression']")
+ self.assertEqual(active.attrib["value"], "true")
+ self.assertEqual(expr.attrib["value"], '"Color"')
+
+ def test_integer_fields_includes_voxel_attrs(self):
+ module, _, _ = self._import_plugin_module()
+ for name in ('NumSamples', 'CellX', 'CellY', 'CellZ'):
+ self.assertIn(name, module.INTEGER_FIELDS)
+
+ def test_build_fields_types_for_voxel_attrs(self):
+ module, _, _ = self._import_plugin_module()
+ rows = [{
+ 'X': 1.0,
+ 'Y': 2.0,
+ 'Z': -5.0,
+ 'Resistivity': 50.0,
+ 'Color': '#00ff00',
+ 'CellX': 0,
+ 'CellY': 0,
+ 'CellZ': 1,
+ 'NumSamples': 3,
+ 'Geometry': 'POINT Z (1 2 -5)',
+ }]
+ fields = module._build_fields(rows)
+ type_map = {field.name(): field.field_type for field in fields}
+ self.assertEqual(type_map['CellX'], module.QMetaType.Type.Int)
+ self.assertEqual(type_map['CellY'], module.QMetaType.Type.Int)
+ self.assertEqual(type_map['CellZ'], module.QMetaType.Type.Int)
+ self.assertEqual(type_map['NumSamples'], module.QMetaType.Type.Int)
+ self.assertEqual(type_map['Resistivity'], module.QMetaType.Type.Double)
+ self.assertEqual(type_map['Color'], module.QMetaType.Type.QString)
diff --git a/test/test_interpolate.py b/test/test_interpolate.py
@@ -0,0 +1,392 @@
+import math
+import unittest
+
+from tem_loader.interpolate import (
+ GridSpec,
+ Sample,
+ _default_radius,
+ axis_cell_size,
+ axis_count,
+ doi_by_xy_from_points,
+ idw,
+ interpolate_voxels,
+ samples_from_layers,
+)
+
+
+class AxisMathTests(unittest.TestCase):
+ def test_axis_count_ceil(self):
+ self.assertEqual(axis_count(10.0, 3.0), 4)
+ self.assertEqual(axis_count(9.0, 3.0), 3)
+ self.assertEqual(axis_count(0.1, 3.0), 1)
+
+ def test_axis_count_exact(self):
+ self.assertEqual(axis_count(100.0, 10.0), 10)
+
+ def test_axis_count_rejects_non_positive_cell_size(self):
+ with self.assertRaises(ValueError):
+ axis_count(10.0, 0.0)
+ with self.assertRaises(ValueError):
+ axis_count(10.0, -1.0)
+
+ def test_axis_count_rejects_non_positive_extent(self):
+ with self.assertRaises(ValueError):
+ axis_count(0.0, 1.0)
+ with self.assertRaises(ValueError):
+ axis_count(-5.0, 1.0)
+
+ def test_axis_cell_size_round_trip(self):
+ extent = 100.0
+ cell_size = 7.0
+ n = axis_count(extent, cell_size)
+ cs = axis_cell_size(extent, n)
+ self.assertLessEqual(cs, cell_size)
+ self.assertGreaterEqual(n * cs, extent)
+
+ def test_axis_cell_size_rejects_zero_count(self):
+ with self.assertRaises(ValueError):
+ axis_cell_size(10.0, 0)
+
+ def test_axis_cell_size_rejects_non_positive_extent(self):
+ with self.assertRaises(ValueError):
+ axis_cell_size(0.0, 5)
+
+
+class GridSpecTests(unittest.TestCase):
+ def test_from_cell_size_cell_count(self):
+ bounds = (0.0, 100.0, 0.0, 50.0, -10.0, 0.0)
+ g = GridSpec.from_cell_size(bounds, 10.0, 10.0, 2.0)
+ self.assertEqual(g.nx, 10)
+ self.assertEqual(g.ny, 5)
+ self.assertEqual(g.nz, 5)
+ self.assertEqual(g.cell_count(), 250)
+
+ def test_from_cell_count_consistency(self):
+ bounds = (0.0, 100.0, 0.0, 50.0, -10.0, 0.0)
+ g = GridSpec.from_cell_count(bounds, 5, 5, 2)
+ self.assertEqual(g.nx, 5)
+ self.assertEqual(g.ny, 5)
+ self.assertEqual(g.nz, 2)
+ self.assertAlmostEqual(g.dx, 20.0)
+ self.assertAlmostEqual(g.dy, 10.0)
+ self.assertAlmostEqual(g.dz, 5.0)
+
+ def test_centers_half_cell_offset(self):
+ g = GridSpec(0.0, 0.0, 0.0, 10.0, 10.0, 10.0, 2, 2, 1)
+ centers = list(g.centers())
+ self.assertEqual(len(centers), 4)
+ xs = sorted({c[3] for c in centers})
+ ys = sorted({c[4] for c in centers})
+ zs = sorted({c[5] for c in centers})
+ self.assertAlmostEqual(xs[0], 5.0)
+ self.assertAlmostEqual(xs[1], 15.0)
+ self.assertAlmostEqual(ys[0], 5.0)
+ self.assertAlmostEqual(ys[1], 15.0)
+ self.assertAlmostEqual(zs[0], 5.0)
+
+ def test_covered_extent_at_least_requested(self):
+ bounds = (0.0, 7.0, 0.0, 5.0, 0.0, 3.0)
+ g = GridSpec.from_cell_size(bounds, 3.0, 3.0, 3.0)
+ self.assertGreaterEqual(g.nx * g.dx, 7.0)
+ self.assertGreaterEqual(g.ny * g.dy, 5.0)
+ self.assertGreaterEqual(g.nz * g.dz, 3.0)
+
+
+class IDWTests(unittest.TestCase):
+ def _s(self, x, y, z, res, doi=None):
+ return Sample(x, y, z, res, doi)
+
+ def test_single_sample_returns_its_value(self):
+ samples = [self._s(0.0, 0.0, 0.0, 100.0)]
+ result = idw(samples, [(1.0, 0.0, 0.0)], radius=None)
+ self.assertAlmostEqual(result[0], 100.0, places=3)
+
+ def test_equidistant_power2_is_mean_in_log_space(self):
+ samples = [
+ self._s(-1.0, 0.0, 0.0, 10.0),
+ self._s(1.0, 0.0, 0.0, 1000.0),
+ ]
+ result = idw(samples, [(0.0, 0.0, 0.0)], power=2.0, log10=True, radius=None)
+ expected = 10.0 ** ((math.log10(10.0) + math.log10(1000.0)) / 2.0)
+ self.assertAlmostEqual(result[0], expected, places=5)
+
+ def test_exact_hit_returns_exact_value(self):
+ samples = [
+ self._s(0.0, 0.0, 0.0, 50.0),
+ self._s(1.0, 0.0, 0.0, 200.0),
+ ]
+ result = idw(samples, [(0.0, 0.0, 0.0)], radius=None)
+ self.assertAlmostEqual(result[0], 50.0)
+
+ def test_radius_excludes_far_points(self):
+ samples = [
+ self._s(0.0, 0.0, 0.0, 50.0),
+ self._s(100.0, 0.0, 0.0, 9999.0),
+ ]
+ result = idw(samples, [(1.0, 0.0, 0.0)], radius=10.0, min_neighbors=1)
+ self.assertIsNotNone(result[0])
+ self.assertAlmostEqual(result[0], 50.0, places=1)
+
+ def test_min_neighbors_unmet_returns_none(self):
+ samples = [self._s(0.0, 0.0, 0.0, 50.0)]
+ result = idw(samples, [(100.0, 0.0, 0.0)], radius=5.0, min_neighbors=2)
+ self.assertIsNone(result[0])
+
+ def test_anisotropy_z_changes_vertical_weighting(self):
+ # Sample A is 1 unit away horizontally; sample B is 10 units away
+ # vertically. With anisotropy_z=1 the vertical distance dominates so A
+ # wins; with anisotropy_z=10 the vertical gap shrinks so B gains weight.
+ samples = [
+ self._s(1.0, 0.0, 0.0, 10.0),
+ self._s(0.0, 0.0, 10.0, 100.0),
+ ]
+ result_tight = idw(samples, [(0.0, 0.0, 0.0)], anisotropy_z=1.0,
+ log10=False, radius=None)
+ result_loose = idw(samples, [(0.0, 0.0, 0.0)], anisotropy_z=10.0,
+ log10=False, radius=None)
+ self.assertIsNotNone(result_tight[0])
+ self.assertIsNotNone(result_loose[0])
+ self.assertNotAlmostEqual(result_tight[0], result_loose[0], places=1)
+
+ def test_log10_two_samples_geometric_mean(self):
+ samples = [
+ self._s(-1.0, 0.0, 0.0, 10.0),
+ self._s(1.0, 0.0, 0.0, 1000.0),
+ ]
+ result = idw(samples, [(0.0, 0.0, 0.0)], log10=True, radius=None)
+ expected = math.sqrt(10.0 * 1000.0)
+ self.assertAlmostEqual(result[0], expected, places=4)
+
+ def test_doi_gate_drops_sounding_below_doi(self):
+ samples = [
+ self._s(0.0, 0.0, -50.0, 100.0, doi=-30.0),
+ self._s(1.0, 0.0, -10.0, 200.0, doi=None),
+ ]
+ result = idw(samples, [(0.5, 0.0, -40.0)], gate_doi=True, radius=None,
+ log10=False)
+ result_no_gate = idw(samples, [(0.5, 0.0, -40.0)], gate_doi=False,
+ radius=None, log10=False)
+ self.assertIsNotNone(result[0])
+ self.assertIsNotNone(result_no_gate[0])
+ self.assertNotAlmostEqual(result[0], result_no_gate[0], places=1)
+
+ def test_doi_gate_none_doi_never_gated(self):
+ samples = [self._s(0.0, 0.0, -100.0, 50.0, doi=None)]
+ result = idw(samples, [(0.0, 0.0, -200.0)], gate_doi=True, radius=None)
+ self.assertIsNotNone(result[0])
+
+ def test_log10_skips_non_positive_resistivity(self):
+ # A zero/negative resistivity has no log; it must not pull the average
+ # toward log10==0 (i.e. 1 ohm-m). Result should equal the positive-only
+ # interpolation.
+ with_bad = [
+ self._s(-1.0, 0.0, 0.0, 100.0),
+ self._s(1.0, 0.0, 0.0, 0.0),
+ ]
+ positive_only = [self._s(-1.0, 0.0, 0.0, 100.0)]
+ result = idw(with_bad, [(0.0, 0.0, 0.0)], log10=True, radius=None)
+ expected = idw(positive_only, [(0.0, 0.0, 0.0)], log10=True, radius=None)
+ self.assertAlmostEqual(result[0], expected[0], places=6)
+
+ def test_multiple_exact_hits_are_averaged(self):
+ # Two samples coincide with the query; both must contribute (geometric
+ # mean in log space), not just whichever appears first.
+ samples = [
+ self._s(0.0, 0.0, 0.0, 10.0),
+ self._s(0.0, 0.0, 0.0, 1000.0),
+ ]
+ result = idw(samples, [(0.0, 0.0, 0.0)], log10=True, radius=None)
+ self.assertAlmostEqual(result[0], math.sqrt(10.0 * 1000.0), places=4)
+
+
+class DefaultRadiusTests(unittest.TestCase):
+ def _s(self, x, y):
+ return Sample(x, y, 0.0, 100.0, None)
+
+ def test_independent_of_sample_order(self):
+ locations = [(0.0, 0.0), (10.0, 0.0), (0.0, 10.0), (10.0, 10.0)]
+ forward = [self._s(x, y) for x, y in locations]
+ reversed_samples = list(reversed(forward))
+ self.assertAlmostEqual(
+ _default_radius(forward), _default_radius(reversed_samples)
+ )
+
+ def test_ignores_duplicate_depth_layers(self):
+ # Many depth layers at the same (x, y) must not collapse the radius to 0.
+ stacked = [Sample(0.0, 0.0, z, 100.0, None) for z in range(20)]
+ stacked += [Sample(30.0, 0.0, z, 100.0, None) for z in range(20)]
+ # Median nearest-neighbor spacing is 30; radius is 3x that.
+ self.assertAlmostEqual(_default_radius(stacked), 90.0)
+
+
+class SamplesFromLayersTests(unittest.TestCase):
+ def _row(self, x, y, zmid, res):
+ return {
+ 'X': x, 'Y': y, 'ZMid': zmid,
+ 'Resistivity': res,
+ 'Geometry': f'LINESTRING Z ({x} {y} {zmid + 1}, {x} {y} {zmid - 1})',
+ }
+
+ def test_basic_extraction(self):
+ rows = [self._row(1.0, 2.0, -5.0, 100.0)]
+ result = samples_from_layers(rows)
+ self.assertEqual(len(result), 1)
+ self.assertEqual(result[0].x, 1.0)
+ self.assertEqual(result[0].y, 2.0)
+ self.assertEqual(result[0].z, -5.0)
+ self.assertEqual(result[0].resistivity, 100.0)
+ self.assertIsNone(result[0].doi_elev)
+
+ def test_doi_join_by_xy(self):
+ rows = [self._row(1.0, 2.0, -5.0, 100.0)]
+ doi_by_xy = {(1.0, 2.0): -30.0}
+ result = samples_from_layers(rows, doi_by_xy)
+ self.assertEqual(result[0].doi_elev, -30.0)
+
+ def test_skips_nan_resistivity(self):
+ rows = [
+ self._row(1.0, 2.0, -5.0, float('nan')),
+ self._row(3.0, 4.0, -5.0, 50.0),
+ ]
+ result = samples_from_layers(rows)
+ self.assertEqual(len(result), 1)
+ self.assertEqual(result[0].x, 3.0)
+
+ def test_skips_none_resistivity(self):
+ rows = [{'X': 1.0, 'Y': 2.0, 'ZMid': -5.0, 'Resistivity': None,
+ 'Geometry': ''}]
+ result = samples_from_layers(rows)
+ self.assertEqual(len(result), 0)
+
+
+class DoiByXYTests(unittest.TestCase):
+ def test_dict_rows(self):
+ doi_points = [{'X': 1.0, 'Y': 2.0, 'ZDOI': -30.0}]
+ result = doi_by_xy_from_points(doi_points)
+ self.assertEqual(result, {(1.0, 2.0): -30.0})
+
+ def test_multiple_points(self):
+ doi_points = [
+ {'X': 0.0, 'Y': 0.0, 'ZDOI': -10.0},
+ {'X': 1.0, 'Y': 1.0, 'ZDOI': -20.0},
+ ]
+ result = doi_by_xy_from_points(doi_points)
+ self.assertEqual(result[(0.0, 0.0)], -10.0)
+ self.assertEqual(result[(1.0, 1.0)], -20.0)
+
+
+class InterpolateVoxelsTests(unittest.TestCase):
+ def _make_layer_rows(self, x, y, zmid, res, doi_elev=None):
+ return {
+ 'X': x, 'Y': y, 'ZMid': zmid,
+ 'Resistivity': res,
+ 'Geometry': f'LINESTRING Z ({x} {y} {zmid + 1}, {x} {y} {zmid - 1})',
+ }
+
+ def test_basic_output_keys(self):
+ rows = [self._make_layer_rows(0.0, 0.0, 0.0, 100.0)]
+ g = GridSpec(0.0, 0.0, 0.0, 10.0, 10.0, 10.0, 1, 1, 1)
+ voxels = interpolate_voxels(rows, g)
+ self.assertEqual(len(voxels), 1)
+ v = voxels[0]
+ for key in ('X', 'Y', 'Z', 'Resistivity', 'Color', 'CellX', 'CellY',
+ 'CellZ', 'NumSamples', 'Geometry'):
+ self.assertIn(key, v)
+
+ def test_color_set_from_ramp(self):
+ from tem_loader.core import resistivity_color as rc
+ rows = [self._make_layer_rows(0.0, 0.0, 0.0, 31.0)]
+ g = GridSpec(0.0, 0.0, 0.0, 10.0, 10.0, 10.0, 1, 1, 1)
+ voxels = interpolate_voxels(rows, g)
+ self.assertEqual(voxels[0]['Color'], rc(voxels[0]['Resistivity']))
+
+ def test_geometry_is_point_z_wkt(self):
+ rows = [self._make_layer_rows(0.0, 0.0, 0.0, 50.0)]
+ g = GridSpec(0.0, 0.0, 0.0, 10.0, 10.0, 10.0, 1, 1, 1)
+ voxels = interpolate_voxels(rows, g)
+ self.assertTrue(voxels[0]['Geometry'].startswith('POINT Z'))
+
+ def test_num_samples_reported(self):
+ rows = [
+ self._make_layer_rows(0.0, 0.0, 0.0, 50.0),
+ self._make_layer_rows(1.0, 1.0, 0.0, 100.0),
+ ]
+ # Cell center is at (1, 1, 1), close to both samples.
+ g = GridSpec(0.0, 0.0, 0.0, 2.0, 2.0, 2.0, 1, 1, 1)
+ voxels = interpolate_voxels(rows, g)
+ self.assertTrue(len(voxels) > 0, 'no voxels produced')
+ self.assertGreaterEqual(voxels[0]['NumSamples'], 1)
+
+ def test_per_sample_doi_gating(self):
+ # Sounding A (0,0) has DOI at z=-60; sounding B (5,0) has DOI at z=-80.
+ # Cells: zmin=-90, dz=30, nz=3 → centers at z=-75, -45, -15.
+ # CellZ=0 (z=-75): A gated (−75 < −60), B usable (−75 > −80). B-only.
+ # CellZ=1 (z=-45): both usable.
+ # CellZ=2 (z=-15): both usable.
+ rows = [
+ self._make_layer_rows(0.0, 0.0, -10.0, 50.0),
+ self._make_layer_rows(5.0, 0.0, -10.0, 200.0),
+ ]
+ doi_dict = [
+ {'X': 0.0, 'Y': 0.0, 'ZDOI': -60.0},
+ {'X': 5.0, 'Y': 0.0, 'ZDOI': -80.0},
+ ]
+ g = GridSpec(-5.0, -5.0, -90.0, 20.0, 20.0, 30.0, 1, 1, 3)
+
+ voxels_by_z = {v['CellZ']: v for v in interpolate_voxels(
+ rows, g, doi_points=doi_dict, clip_below_doi=True, min_neighbors=1
+ )}
+ # CellZ=0 (z=-75): only B should contribute
+ b_only = voxels_by_z.get(0)
+ self.assertIsNotNone(b_only, 'expected a voxel in B-only zone at z=-75')
+ self.assertEqual(b_only['NumSamples'], 1)
+ # CellZ=1 (z=-45): both should contribute
+ both = voxels_by_z.get(1)
+ self.assertIsNotNone(both, 'expected a voxel in both-usable zone at z=-45')
+ self.assertGreaterEqual(both['NumSamples'], 2)
+
+ def test_safety_cap_raises(self):
+ bounds = (0.0, 200.0, 0.0, 200.0, 0.0, 200.0)
+ g = GridSpec.from_cell_size(bounds, 1.0, 1.0, 1.0)
+ rows = [self._make_layer_rows(100.0, 100.0, 100.0, 50.0)]
+ with self.assertRaises(ValueError):
+ interpolate_voxels(rows, g)
+
+ def test_progress_callback(self):
+ rows = [self._make_layer_rows(0.0, 0.0, 0.0, 50.0)]
+ g = GridSpec(0.0, 0.0, 0.0, 10.0, 10.0, 10.0, 2, 2, 2)
+ calls = []
+ interpolate_voxels(rows, g, progress=lambda done, total: calls.append((done, total)))
+ self.assertTrue(all(total == g.cell_count() for _, total in calls))
+ self.assertEqual(calls[-1][0], g.cell_count())
+
+ def test_progress_truthy_return_cancels(self):
+ rows = [self._make_layer_rows(0.0, 0.0, 0.0, 50.0)]
+ g = GridSpec(0.0, 0.0, 0.0, 10.0, 10.0, 10.0, 4, 4, 4)
+ calls = []
+
+ def cancel_after_first(done, total):
+ calls.append(done)
+ return done >= 1
+
+ interpolate_voxels(rows, g, progress=cancel_after_first)
+ self.assertEqual(calls, [1])
+
+ def test_matches_idw_values(self):
+ # interpolate_voxels and idw share one core, so a voxel's resistivity
+ # must equal idw evaluated at the same cell center.
+ rows = [
+ self._make_layer_rows(0.0, 0.0, 0.0, 50.0),
+ self._make_layer_rows(20.0, 0.0, -5.0, 500.0),
+ ]
+ g = GridSpec(-10.0, -10.0, -10.0, 10.0, 10.0, 10.0, 3, 1, 2)
+ voxels = interpolate_voxels(rows, g, power=2.0, log10=True)
+ samples = samples_from_layers(rows)
+ for v in voxels:
+ expected = idw(samples, [(v['X'], v['Y'], v['Z'])],
+ power=2.0, log10=True)[0]
+ self.assertAlmostEqual(v['Resistivity'], expected, places=6)
+
+
+if __name__ == '__main__':
+ unittest.main()
