geovector-3d-mcp

Downhole 3-D geochemical vectoring with honest uncertainty, plus an MCP server so AI agents can call the maths instead of improvising it.

You have six holes through a pegmatite, all sub-economic. The fractionation indices change from hole to hole: K/Rb falls, Nb/Ta falls, Cs rises. The question every exploration team asks next is which way, and how sure are we? This package answers it from collar, survey, assay and lithology tables:

1. Desurvey every hole by minimum curvature, so sample positions are true 3-D points, not collar-plus-depth guesses.

2. Clean the geochemistry. It parses below-detection values, uses log-ratios instead of raw ratios, and computes indices only in the rock they are defined for.

3. Fit a local 3-D gradient of the index, aligned to the host structure, with a 95 % confidence cone on its direction and a significance test that it exists.

4. Propose a step-out hole (collar, azimuth, dip, target depth), or decline to when the direction is not resolved.

The benchmark data in this repository are synthetic, generated by scripts/make_benchmark.py with a planted answer. Recovering that answer shows the code is correct. It does not show that real pegmatites zone this simply.

Quick start

git clone https://github.com/YASHUJOSHI/geovector-3d-mcp.git
cd geovector-3d-mcp
pip install -e ".[dev]"          # core + MCP + Streamlit + pyproj + pytest
pytest                           # 43 tests
streamlit run app/streamlit_viewer.py

From Python:

from geovector3d import run_vectoring

r = run_vectoring("data/benchmark_stuck_block", index="ln_K_Rb",
                  host_codes=("PEG",), strike=300, dip=50)
r["gradient"]["direction_trend_deg"], r["gradient"]["cone_half_angle_deg"]
r["proposal"]          # collar, azimuth, dip, target depth, confidence, caveats

Input tables

Related MCP server: GeoAgent MCP Server

Results on the benchmark

Six inclined holes cut a 20 m pegmatite dyke (strike 300°, dip 50° NE). K/Rb falls from 140 to 38 and Nb/Ta from 12 to 2.5 down a 34° plunge toward 335°. Planted truth: trend 335.0° / plunge 34.4°.

The last two rows are the point of the package. Without structure, and without keeping wall rock out of the fit, the same data give a wrong or empty answer. A tool that reported a direction regardless would have sent a hole the wrong way.

Method

METHODOLOGY.md has the derivations. In short:

Minimum-curvature desurveying

Each survey station becomes a unit tangent t = (cos δ sin α, cos δ cos α, sin δ), with azimuth α and dip δ negative down. Between stations i and i+1, separated by ΔMD:

• dogleg angle β = 2 arcsin(|tᵢ₊₁ − tᵢ| / 2)

• ratio factor RF = (2/β) tan(β/2)

• displacement Δx = (ΔMD/2) · RF · (tᵢ + tᵢ₊₁)

β is computed from the chord rather than arccos(tᵢ·tᵢ₊₁), which loses most of its precision for the nearly straight courses that dominate real surveys. RF is replaced by its limit 1 + β²/12 when β < 10⁻⁶. A sample midpoint between stations is placed on the arc, not the chord. The sub-arc to a fraction f ends at the slerp of the two tangents and is itself a minimum-curvature course of length f·ΔMD. The tests check this against a constant-build circular arc to 10⁻⁸ m.

Because everything is carried as vectors, a vertical hole is not a special case. At δ = −90° the azimuth multiplies zero and drops out, so there is no gimbal lock.

Compositional data

• Non-detects are kept as censored values and replaced multiplicatively (Martín-Fernández et al., 2003) at 0.65 × LLD. This leaves ratios among detected parts unchanged. Analytes censored in more than 30 % of samples are reported as unusable, and a ratio of two censored values is set to missing.

• Vectoring runs on log-ratios (ln K/Rb, ln Nb/Ta, …). A pairwise log-ratio is scale-invariant and subcompositionally coherent. Raw ratios are skewed, and their gradients are dominated by the high end. CLR and ILR are provided for multi-element work but are computed on a fixed, closed composition, because CLR scores change whenever the analyte suite changes.

• Indices are host-restricted. LCT ratios are computed only in pegmatite and granite codes, and the Ishikawa AI and CCPI only in felsic volcanics.

Local gradient

A weighted least-squares fit of f(x) = a + g·(x − x₀) around x₀, with three safeguards for the way drillholes sample space:

The cone is arctan(√(χ²₂,₀.₉₅ · λ⊥) / |g|), where λ⊥ is the largest variance perpendicular to g. It is a better guide than the design-matrix condition number, which is reported too, because it scales with the signal.

Step-out design

The target is placed along the in-plane fractionation direction, starting from the most-evolved decile of host samples. The default step is the data-support distance (median spacing of host intercepts). The hole is aimed down the host pole to cut the body at a high angle, with dip clipped to [−85°, −45°]. Longer steps are allowed but flagged beyond_support. A local gradient cannot tell you whether a trend continues for another 120 m: that is the hypothesis the hole tests.

Is a barren surface sample barren?

diagnose_false_surface_anomaly checks that immobile ratios (Zr/Hf, Nb/Ta, Ti/Zr) match the fresh parent, which rules out transported cover. It then computes tau mass-transfer coefficients (Brimhall & Dietrich, 1987) for mobile elements. It returns evidence and caveats, not a verdict, because Zr and Ti are only relatively immobile.

MCP server for AI agents

Agents should not do drillhole trigonometry or log-ratio algebra in their heads. The server exposes the engine as four deterministic tools. Every response states its sign conventions, and bad input returns an explicit error.

Built on the MCP Python SDK 2.x (MCPServer, formerly FastMCP). Add it to Claude Desktop's claude_desktop_config.json:

{
  "mcpServers": {
    "geovector3d": { "command": "geovector3d-mcp" }
  }
}

If the console script is not on your PATH, use "command": "python", "args": ["-m", "geovector3d.mcp_server"].

Repository layout

src/geovector3d/
  desurvey.py        minimum curvature, validation, interval mid-points
  compositional.py   censoring, units, CLR/ILR, LCT and alteration indices
  vector_field.py    structural frame, local gradient, cone, step-out design
  anomaly.py         immobile ratios and tau mass transfer
  workflow.py        load -> prepare -> gradient -> proposal
  export.py          CSV drill plan, GeoJSON (needs an EPSG code)
  mcp_server.py      MCP tools
app/streamlit_viewer.py
data/benchmark_stuck_block/   synthetic six-hole dataset + truth.json
scripts/                      benchmark and figure generators
tests/                        43 tests

Limitations

• Linear local model: zoning that curves within the search radius is averaged.

• One structural orientation per run. Folded hosts need domaining; unfolding is out of scope.

• Fractionation indices assume a single magmatic system. Overprinting, alteration or multiple intrusive phases break the simple gradient picture.

• GeoJSON export needs a real CRS; local mine grids export as CSV.

• Validated only on the synthetic benchmark. Real-data case studies are the next step.

How this was built

I wrote this code with an AI coding assistant. The geological framing is mine: the choice of indices, host restriction, structure-aligned search, refusing to recommend without a significant gradient, and flagging step-outs beyond data support.

References

• Brimhall, G.H. & Dietrich, W.E. (1987). Geochim. Cosmochim. Acta 51, 567–587.

• Černý, P. (1991). Rare-element granitic pegmatites. Geoscience Canada 18, 49–67.

• Egozcue, J.J. et al. (2003). Isometric logratio transformations. Math. Geol. 35, 279–300.

• Ishikawa, Y. et al. (1976). Mining Geology 26, 105–117. Alteration index.

• Large, R.R. et al. (2001). Econ. Geol. 96, 957–971. AI and CCPI box plot.

• London, D. (2008). Pegmatites. Can. Mineral. Spec. Publ. 10.

• Martín-Fernández, J.A. et al. (2003). Math. Geol. 35, 253–278. Zero replacement.

• Sawaryn, S.J. & Thorogood, J.L. (2005). SPE Drill. Complet. 20, 77–86. Minimum curvature.

MIT licence.