// =========================================================================
// SHORE PROFILE GATE — is there a beach to put waves on?
//
// An FFT ocean's entire shore treatment -- displacement attenuation, shoaling,
// breaker foam -- is driven by DEPTH. All of it lives in the band between the
// waterline and roughly 20 m down, where waves feel the bottom. If the hypso
// curve's continental slope carries the seabed through that band in a few
// horizontal metres, none of it is visible at any wave scale and the coastal
// buffer has to come first.
//
// So measure the band before building anything on it: walk transects across
// the real mountain field, find the waterline crossing, and report how far you
// have to walk seaward to reach 5 m and 20 m of depth.
//
// ⚠ MOUNTAIN + HYPSO ONLY, NO EROSION FILTER. This reuses the same chain
// test_hypso_curve measures land fraction with. Erosion perturbs the shore on
// top of this; the large-scale profile that decides whether a surf zone exists
// at all is set here. Treat the numbers as the envelope, not the coastline.
//
// Build:
// g++ -O2 -std=c++17 -I ../proc-nodes/src/shared \
// -o /tmp/shore_gate tests/test_shore_profile.cpp
// =========================================================================
#include <cstdio>
#include <cmath>
#include <vector>
#include <algorithm>
#include "erosion_cpp_shim.h"
namespace dcerosion {
#include "hypso_curve.inc"
}
using dcerosion::vec2;
using dcerosion::vec3;
using dcerosion::fractal_multifractal_2d_lod;
using dcerosion::hypso_curve;
// map_dc_erosion.gd's constants.
static const float SEA = 0.00f;
static const float SHELF_W = 0.12f;
static const float SLOPE_W = 0.05f;
static const float ABYSS = -0.62f;
static const float PLAT_SLOPE = 0.13f;
static const float MTN_ONSET = 0.35f;
static const float MTN_WIDTH = 0.10f;
static const float MTN_SLOPE = 1.30f;
static const float PEAK_ONSET = 0.60f;
static const float PEAK_GAIN = 0.30f;
static const float TH = 60000.0f, MSCALE = 0.1f, LAC = 1.87f, MGAIN = 0.68f;
static const int MOCT = 8;
static const float VARIANCE = 0.0f;
static const float M_PER_UNIT = 7790.0f; // curve units -> world metres
// Surface elevation in metres relative to the waterline. Negative = underwater.
static float depth_at(float x, float z, float waterline) {
vec3 m = fractal_multifractal_2d_lod(
vec2(x / TH * MSCALE, z / TH * MSCALE), 1.0f,
MOCT, LAC, MGAIN, VARIANCE, 999.0f);
float o = hypso_curve(m.x, SEA, SHELF_W, SLOPE_W, ABYSS, PLAT_SLOPE,
MTN_ONSET, MTN_WIDTH, MTN_SLOPE,
PEAK_ONSET, PEAK_GAIN).x;
return (waterline - o) * M_PER_UNIT; // +ve = water depth
}
int main() {
float waterline = hypso_curve(SEA, SEA, SHELF_W, SLOPE_W, ABYSS, PLAT_SLOPE,
MTN_ONSET, MTN_WIDTH, MTN_SLOPE,
PEAK_ONSET, PEAK_GAIN).x;
printf("SHORE PROFILE GATE — width of the band where waves feel the bottom\n");
printf(" waterline at curve %+.4f, %.0f m per curve unit\n", waterline, M_PER_UNIT);
const float HALF = 458800.0f;
const int TRANSECTS = 64;
const float COARSE = 200.0f; // m, to bracket the crossing
const float FINE = 2.0f; // m, to walk the band
std::vector<float> w5, w20;
int found = 0;
for (int t = 0; t < TRANSECTS; t++) {
// Transects along +x at assorted z, which cross the diagonal shore.
float z = -HALF + 2.0f * HALF * (t + 0.5f) / TRANSECTS;
float prev = depth_at(-HALF, z, waterline);
for (float x = -HALF + COARSE; x < HALF; x += COARSE) {
float d = depth_at(x, z, waterline);
if ((prev < 0.0f) != (d < 0.0f)) {
// Bracketed a waterline crossing. Bisect to ~1 m.
float lo = x - COARSE, hi = x;
for (int k = 0; k < 24; k++) {
float mid = 0.5f * (lo + hi);
if ((depth_at(lo, z, waterline) < 0.0f) != (depth_at(mid, z, waterline) < 0.0f)) hi = mid;
else lo = mid;
}
float shore = 0.5f * (lo + hi);
// Walk seaward (whichever direction deepens) to 5 m and 20 m.
float dir = (depth_at(shore + FINE, z, waterline) > depth_at(shore - FINE, z, waterline)) ? 1.0f : -1.0f;
float d5 = -1.0f, d20 = -1.0f;
for (float s = 0.0f; s < 20000.0f; s += FINE) {
float dd = depth_at(shore + dir * s, z, waterline);
if (d5 < 0.0f && dd >= 5.0f) d5 = s;
if (d20 < 0.0f && dd >= 20.0f) { d20 = s; break; }
}
if (d5 >= 0.0f && d20 >= 0.0f) { w5.push_back(d5); w20.push_back(d20); found++; }
break;
}
prev = d;
}
}
if (found == 0) { printf(" FAIL — no waterline crossing found on any transect\n"); return 1; }
std::sort(w5.begin(), w5.end());
std::sort(w20.begin(), w20.end());
auto med = [](std::vector<float>& v) { return v[v.size() / 2]; };
printf(" %d transects crossed the shore\n", found);
printf(" horizontal run from waterline to 5 m depth: median %8.1f m (min %7.1f, max %8.1f)\n",
med(w5), w5.front(), w5.back());
printf(" horizontal run from waterline to 20 m depth: median %8.1f m (min %7.1f, max %8.1f)\n",
med(w20), w20.front(), w20.back());
float m20 = med(w20);
printf("\n");
if (m20 < 50.0f) {
printf(" VERDICT: NO SURF ZONE. The seabed reaches 20 m within %.0f m of\n", m20);
printf(" the shore, so shoaling, breakers and displacement attenuation all\n");
printf(" collapse into a band thinner than a single wave. A coastal buffer\n");
printf(" must flatten the nearshore BEFORE any of the ocean shore work.\n");
return 1;
}
printf(" VERDICT: surf zone is %.0f m wide at the median — wide enough to\n", m20);
printf(" carry shoaling and breakers at ocean wavelengths (tens of metres).\n");
return 0;
}