import_bag: seed CUBE predicted surface from Chart prior to reject false-deep detections (#89)

cube_bathymetry/include/cube_bathymetry/store_import.h

@@ -146,9 +146,23 @@ namespace cube
 /// uncertainty squared, floored at a small positive epsilon (Node::setPredictedDepth
 /// requires a finite positive variance; a single-sample CUBE cell can persist a
 /// zero or non-finite uncertainty, so every finite-depth cell is still seeded).
-/// Does NOT insert a CUBE hypothesis and does NOT mark the sheet dirty.
+/// Does NOT mark the sheet dirty.
+///
+/// @param seed_settled When true (default), ALSO reseed each cell's SETTLED depth
+///   as a CUBE hypothesis (`setSettledDepthAt`) so the cell round-trips through
+///   `values()` and survives the next whole-tile save -- the lossless warm-start
+///   reload of a persisted draft tile (ADR-0001, #21/#70). When false, seed ONLY
+///   the predicted surface: the cell carries a **blunder-rejection prior** (it
+///   turns on `Node::insert`'s predicted-surface gate so a false-deep sounding
+///   below `target - blunder_scalar*sqrt(var)` is rejected) but does NOT
+///   fill/settle the cell -- no hypothesis, no `values()` output, the sheet stays
+///   clean. This is the `Chart` (contour) prior path (cube#89): settling coarse
+///   contour depths would contaminate the survey layer (and its co-estimated
+///   backscatter) with non-measured fill, so the prior only gates, it does not
+///   fill; survey-falls-through-to-chart gap-filling stays a query-time concern.
   void primeFromTile(
-    const marine_bathymetry_store::BathymetryTile & tile, GeoMapSheet & map_sheet);
+    const marine_bathymetry_store::BathymetryTile & tile, GeoMapSheet & map_sheet,
+    bool seed_settled = true);
 
 /// @brief Load every tile of @p layer from @p store into @p map_sheet,
 ///        priming predicted depths cell-by-cell.
@@ -158,10 +172,15 @@ namespace cube
 /// are seeded. Warm-starts slope correction from persisted draft tiles on
 /// restart; does not reconstruct CUBE hypothesis state (#21). A no-op if @p layer
 /// holds no tiles.
+///
+/// @param seed_settled Forwarded to @ref primeFromTile (default true =
+///   settled+predicted warm-start reload; false = predicted-only blunder-rejection
+///   prior, e.g. priming the predicted surface from a `Chart` layer, cube#89).
   void loadIntoSheet(
     const marine_bathymetry_store::BathymetryStore & store,
     marine_bathymetry_store::SourceLayer layer,
-    GeoMapSheet & map_sheet);
+    GeoMapSheet & map_sheet,
+    bool seed_settled = true);
 
 }  // namespace cube
 

cube_bathymetry/src/import_bag_main.cpp

@@ -73,6 +73,11 @@
     "node writes draft. #85)\n";
   std::cout << "  -o <store_dir>: Output bathymetry-store directory (created if "
     "needed)\n";
+  std::cout << "  --prior <store_dir>: seed the CUBE predicted surface from this "
+    "store's Chart (contour) layer so blunder rejection drops false-deep "
+    "detections (#89). The chart must be at the survey GGGS level. NOTE: a "
+    "coarse/shallow-biased chart can also reject LEGITIMATE deeper-than-charted "
+    "returns; the rejection margin is tunable via the blunder_* params\n";
   std::cout << "  --bs-store <dir>: Also write an MBES backscatter store layer "
     "(Processed) from the same CUBE pass (optional; surfaces the co-estimated "
     "intensity -- uncorrected by default, angle-corrected when "
@@ -278,6 +283,7 @@ int main(int argc, char * argv[])
 
   std::vector<std::string> bagfile_names;
   std::string store_dir;
+  std::string prior_dir;  // optional: Chart-prior store to seed the predicted surface (#89)
   std::string bs_store_dir;  // optional: MBES backscatter store output (#80)
   std::string bathy_layer_str = "processed";  // bathy target layer (#85): draft|processed
   std::string detections_topic;  // required
@@ -318,6 +324,8 @@ int main(int argc, char * argv[])
       usage();
     } else if (*arg == "-o") {
       store_dir = next_value("-o");
+    } else if (*arg == "--prior") {
+      prior_dir = next_value("--prior");
     } else if (*arg == "--bathy-layer") {
       bathy_layer_str = next_value("--bathy-layer");
     } else if (*arg == "--bs-store") {
@@ -489,6 +497,81 @@ int main(int argc, char * argv[])
     backscatter_mode, std::move(backscatter_curve.points),
     backscatter_curve.tl_removed, backscatter_curve.absorption_db_per_m);
 
+  // Chart-prior prime (#89): seed the CUBE predicted surface from the prior store's
+  // Chart (contour) layer BEFORE any soundings are added. primeFromTile lazy-creates
+  // a node per chart cell carrying a predicted depth, which turns ON CUBE's
+  // predicted-surface blunder gate (Node::insert is a pass-through when
+  // predicted_depth_ is NaN); the survey soundings then hit nodes that already carry
+  // a predicted depth, so a false-deep detection below
+  // `target - blunder_scalar*sqrt(predicted_var)` is rejected. seed_settled=false:
+  // seed ONLY the predicted surface, never settle the coarse contour into the survey
+  // layer (that would contaminate both the bathy and co-estimated backscatter).
+  //
+  // ALIGNMENT: load() restores each Chart tile at the GGGS level encoded in its
+  // filename (the store is multi-level; the fromCellSize() arg below only sets the
+  // store's default cellIndex() level, it does NOT pin the level tiles load at). A
+  // chart cell only coincides with -- and thus gates -- a survey node when its level
+  // matches the survey level, i.e. the Chart layer was imported at the survey's
+  // resolution. A chart tile at a different level primes a node the survey soundings
+  // never land on, so it cannot gate. Rather than silently build an un-gated store
+  // (the operator passed --prior to reject outliers), we validate below and refuse
+  // when nothing would gate. A future refinement could resample the chart.
+  if (!prior_dir.empty()) {
+    marine_bathymetry_store::BathymetryStore prior =
+      marine_bathymetry_store::BathymetryStore::fromCellSize(
+      static_cast<float>(geo_map_sheet.nominalCellSizeMeters()));
+    marine_bathymetry_store::SourceRegistry prior_reg;
+    try {
+      marine_bathymetry_store::load(prior, prior_dir, &prior_reg);
+    } catch (const std::exception & e) {
+      std::cerr << "ERROR: --prior store '" << prior_dir << "' failed to load: "
+                << e.what() << std::endl;
+      return 1;
+    }
+    const auto & chart_tiles =
+      prior.tiles(marine_bathymetry_store::SourceLayer::Chart);
+    const std::size_t prior_tiles = chart_tiles.size();
+    if (prior_tiles == 0) {
+      std::cerr << "ERROR: --prior store '" << prior_dir << "' has no Chart-layer "
+                << "tiles; nothing would seed the predicted surface and blunder "
+                << "rejection would be INACTIVE. Check the path and that the store "
+                << "has a chart/ layer (refusing to build an un-gated store)."
+                << std::endl;
+      return 1;
+    }
+    // Only tiles at the survey level can coincide with a survey node and gate
+    // (see ALIGNMENT above). The prior store's default level is fromCellSize() of
+    // the survey cell size, so it equals the survey sheet's level.
+    const uint8_t survey_level = prior.level().level();
+    std::size_t mismatched = 0;
+    for (const auto & grid_tile : chart_tiles) {
+      if (grid_tile.second.index().level() != survey_level) {
+        ++mismatched;
+      }
+    }
+    if (mismatched == prior_tiles) {
+      std::cerr << "ERROR: all " << prior_tiles << " Chart tile(s) in '" << prior_dir
+                << "' are at a GGGS level other than the survey level "
+                << static_cast<int>(survey_level) << "; none would coincide with a "
+                << "survey node, so blunder rejection would be INACTIVE. Re-import "
+                << "the chart at the survey resolution (refusing to build an "
+                << "un-gated store)." << std::endl;
+      return 1;
+    }
+    if (mismatched > 0) {
+      std::cerr << "WARNING: " << mismatched << " of " << prior_tiles << " Chart "
+                << "tile(s) are not at the survey level "
+                << static_cast<int>(survey_level) << " and will not gate; only the "
+                << "level-matched tiles seed the predicted surface." << std::endl;
+    }
+    cube::loadIntoSheet(
+      prior, marine_bathymetry_store::SourceLayer::Chart, geo_map_sheet,
+      /*seed_settled=*/false);
+    std::cout << "Primed CUBE predicted surface from " << (prior_tiles - mismatched)
+              << " of " << prior_tiles << " Chart tile(s) in " << prior_dir
+              << " (blunder rejection active, #89)." << std::endl;
+  }
+
   std::cout << "reading messages..." << std::endl;
 
   int ping_count = 0;

cube_bathymetry/src/store_import.cpp

@@ -149,7 +149,8 @@ mapSheetToBackscatterCells(
 }
 
 void primeFromTile(
-  const marine_bathymetry_store::BathymetryTile & tile, GeoMapSheet & map_sheet)
+  const marine_bathymetry_store::BathymetryTile & tile, GeoMapSheet & map_sheet,
+  bool seed_settled)
 {
   const gggs::GridIndex & grid = tile.index();
   const std::vector<double> & depth = tile.depthBand();
@@ -178,6 +179,16 @@ void primeFromTile(
     map_sheet.setPredictedDepthAt(
       *it, static_cast<float>(d), static_cast<float>(variance));
 
+    if (!seed_settled) {
+      // Predicted-only prime (cube#89): seed the slope/blunder-rejection prior
+      // but do NOT settle the cell. A `Chart` (contour) prior must turn the
+      // blunder gate on WITHOUT filling the survey layer with coarse contour
+      // depths (which would also contaminate the co-estimated backscatter with
+      // non-measured cells). The survey accumulates real data on top; gap-filling
+      // survey->chart is a query-time overlay, not a settled fill.
+      continue;
+    }
+
     // Lossless reload (ADR-0001): also reseed the SETTLED depth as a CUBE
     // hypothesis so the primed cell round-trips through values() and survives the
     // next whole-tile save. Without this, a tile primed at startup, partially
@@ -191,12 +202,13 @@ void primeFromTile(
 void loadIntoSheet(
   const marine_bathymetry_store::BathymetryStore & store,
   marine_bathymetry_store::SourceLayer layer,
-  GeoMapSheet & map_sheet)
+  GeoMapSheet & map_sheet,
+  bool seed_settled)
 {
   // Single fused grid per layer (unh_marine_autonomy#221): iterate the layer's
   // tiles directly. Empty map -> no-op.
   for (const auto & grid_tile : store.tiles(layer)) {
-    primeFromTile(grid_tile.second, map_sheet);
+    primeFromTile(grid_tile.second, map_sheet, seed_settled);
   }
 }
 

cube_bathymetry/test/test_store_import.cpp

@@ -76,6 +76,25 @@ std::vector<GeoSounding> makeSoundingsWithIntensity(float intensity)
   return soundings;
 }
 
+// The SAME lat/lon footprint as makeSoundings() (so the touchdown cells coincide
+// cell-for-cell) but every sounding sits at a single, much deeper depth. Used to
+// stage a false-deep "blunder" against a shallow predicted surface (#89).
+std::vector<GeoSounding> makeDeepSoundings(float depth)
+{
+  std::vector<GeoSounding> soundings;
+  const double base_lat = 43.07;
+  const double base_lon = -70.76;
+  for (int i = 0; i < 5; ++i) {
+    gz4d::GeoPointLatLongDegrees point(
+      base_lat + i * 1e-5, base_lon + i * 1e-5, depth);
+    GeoSounding s(point);
+    s.sounding.vertical_error = 0.5f;
+    s.sounding.horizontal_error = 0.1f;
+    soundings.push_back(s);
+  }
+  return soundings;
+}
+
 constexpr int64_t kStamp = 1234567890123456789LL;
 constexpr uint16_t kSource = 7;
 constexpr float kIntensity = 42.5f;
@@ -370,4 +389,129 @@ TEST(StoreImport, BackscatterNaNPropagation)
     << "soundings without intensity must surface no backscatter cells";
 }
 
+// Predicted-only prime (seed_settled=false, #89): seeds the predicted surface
+// (turns the blunder gate on) but must NOT settle the cell -- no CUBE hypothesis,
+// no values() output, sheet stays clean. Contrasted in the same test against the
+// default seed_settled=true, which DOES settle the cell.
+TEST(StoreImport, PredictedOnlyPrimeSeedsNoSettledHypothesis)
+{
+  // Build a tile with finite depths from the synthetic soundings.
+  GeoMapSheet source(1.0f);
+  source.addSoundings(makeSoundings());
+  auto grids = source.grids();
+  ASSERT_FALSE(grids.empty());
+
+  const std::vector<DepthAndUncertainty> values = grids.front()->values();
+  const marine_bathymetry_store::BathymetryTile tile =
+    geoGridToTile(*grids.front(), kStamp, kSource);
+
+  // Predicted-only prime: seeds the predicted surface but creates no hypothesis.
+  GeoMapSheet predicted_only(1.0f);
+  primeFromTile(tile, predicted_only, /*seed_settled=*/false);
+
+  // Contrast: the default prime (seed_settled=true) DOES settle each cell.
+  GeoMapSheet settled(1.0f);
+  primeFromTile(tile, settled, /*seed_settled=*/true);
+
+  auto po_grid = predicted_only.gridAt(tile.index());
+  auto st_grid = settled.gridAt(tile.index());
+  ASSERT_NE(po_grid, nullptr);
+  ASSERT_NE(st_grid, nullptr);
+
+  // values() is a const method (it flushes the median pre-filter), so it is safe
+  // to call through the const grid handle.
+  const std::vector<DepthAndUncertainty> po_values = po_grid->values();
+  const std::vector<DepthAndUncertainty> st_values = st_grid->values();
+
+  gggs::CellAreaIterator it(tile.index());
+  std::size_t k = 0;
+  std::size_t checked = 0;
+  for (; it.valid() && k < values.size(); it.next(), ++k) {
+    if (std::isnan(values[k].depth)) {
+      continue;
+    }
+    ++checked;
+    // The predicted surface IS seeded for the predicted-only prime.
+    EXPECT_FLOAT_EQ(po_grid->predictedDepthAt(*it), values[k].depth);
+    // Predicted-only: NO settled hypothesis -> values() carries no estimate.
+    EXPECT_TRUE(std::isnan(po_values[k].depth))
+      << "predicted-only prime must not settle a hypothesis";
+    // Default prime: the SAME cell carries a settled estimate (the distinguishing
+    // contrast -- predicted-only vs predicted+settled).
+    EXPECT_FALSE(std::isnan(st_values[k].depth))
+      << "default prime (seed_settled=true) must settle the cell";
+  }
+  EXPECT_GT(checked, 0u) << "tile should carry some finite cells to prime";
+
+  // Priming reproduces persisted data -- it must never mark the sheet dirty.
+  EXPECT_TRUE(predicted_only.dirtyGrids().empty());
+}
+
+// The load-bearing behavior (#89): a seeded (shallow) predicted surface turns on
+// CUBE's blunder gate so a false-deep sounding is REJECTED, whereas with no prior
+// the same deep sounding is accepted. Demonstrates the bathy store is protected
+// from false-deep detections by the Chart-prior prime.
+TEST(StoreImport, SeededPredictedSurfaceRejectsDeepBlunder)
+{
+  // A shallow predicted surface (~ -10 m) from synthetic soundings -> the Chart-
+  // like prior tile we seed from.
+  GeoMapSheet shallow_src(1.0f);
+  shallow_src.addSoundings(makeSoundings());
+  auto shallow_tiles = mapSheetToTiles(shallow_src, kStamp, kSource);
+  ASSERT_FALSE(shallow_tiles.empty());
+
+  // A clearly-too-deep blunder at the SAME locations (identical touchdown cells),
+  // far below any shallow-surface blunder limit. With target ~ -10 m and the
+  // defaults (blunder_minimum 10, blunder_percent 0.25, blunder_scalar 3), the
+  // limit is ~ -20 m, so -150 m is unambiguously a blunder.
+  constexpr float kDeep = -150.0f;
+  const std::vector<GeoSounding> deep = makeDeepSoundings(kDeep);
+
+  // Baseline: WITHOUT a prior, predicted_depth_ stays INVALID_DATA, so Node::insert
+  // skips the blunder gate and the deep sounding is accepted -> finite deep cells.
+  GeoMapSheet no_prior(1.0f);
+  no_prior.addSoundings(deep);
+  std::size_t accepted_deep = 0;
+  for (const auto & grid : no_prior.grids()) {
+    for (const auto & v : grid->values()) {
+      if (!std::isnan(v.depth) && v.depth < -100.0f) {
+        ++accepted_deep;
+      }
+    }
+  }
+  EXPECT_GT(accepted_deep, 0u)
+    << "without a predicted surface the deep blunder must be accepted";
+
+  // WITH the shallow predicted surface primed (predicted-only): the blunder gate is
+  // active at every primed cell, so the deep sounding that touches it is rejected.
+  GeoMapSheet primed(1.0f);
+  for (const auto & gt : shallow_tiles) {
+    primeFromTile(gt.second, primed, /*seed_settled=*/false);
+  }
+  primed.addSoundings(deep);
+
+  // Every cell that carries a (shallow) predicted depth must end with NO settled
+  // estimate: the deep sounding it received was blunder-rejected, and the
+  // predicted-only prime left no hypothesis of its own.
+  std::size_t gated_cells = 0;
+  for (const auto & gt : shallow_tiles) {
+    const auto & tile = gt.second;
+    auto grid = primed.gridAt(tile.index());
+    ASSERT_NE(grid, nullptr);
+    const std::vector<DepthAndUncertainty> vals = grid->values();
+    const std::vector<double> & depth = tile.depthBand();
+    gggs::CellAreaIterator it(tile.index());
+    std::size_t k = 0;
+    for (; it.valid() && k < depth.size() && k < vals.size(); it.next(), ++k) {
+      if (std::isnan(depth[k])) {
+        continue;  // not a primed cell
+      }
+      ++gated_cells;
+      EXPECT_TRUE(std::isnan(vals[k].depth))
+        << "a deep blunder must be rejected where a shallow predicted surface gates";
+    }
+  }
+  EXPECT_GT(gated_cells, 0u) << "the shallow tile should prime some cells to gate";
+}
+
 }  // namespace cube