// Headless test for the Live World ocean/sky stage (moons + tides). No display needed.
//
//   g++ -std=c++17 -O2 -Isrc/sim test_ocean.cpp src/sim/IcoSphere.cpp
//       src/sim/Planet.cpp src/sim/PlanetTectonics.cpp src/sim/PlanetDrift.cpp
//       src/sim/PlanetErosion.cpp src/sim/PlanetHydrology.cpp
//       src/sim/PlanetBiomes.cpp src/sim/PlanetClimate.cpp src/sim/PlanetLive.cpp
//       src/sim/PlanetOcean.cpp src/sim/PlanetBiota.cpp src/sim/PlanetFloraGen.cpp
//       src/sim/PlanetFaunaGen.cpp src/sim/PlanetFungiGen.cpp src/sim/PlanetIO.cpp
//       -o /tmp/to && /tmp/to
//
// Verifies: 1-3 moons generated deterministically; moon/sun directions are unit vectors that
// sweep with the clock; the equilibrium tide is a zero-mean two-bulge field (high under a body
// and its antipode, low at 90 deg); and save v9 round-trips the moons.

#include "Planet.hpp"
#include "Projection.hpp"
#include <cstdio>
#include <cmath>
#include <algorithm>
#include <sstream>

static int failures = 0;
static void check(bool cond, const char* what) {
    std::printf("  [%s] %s\n", cond ? "PASS" : "FAIL", what);
    if (!cond) ++failures;
}
static int nearestCell(const Planet& p, const Vec3& dir) {
    int best = 0; double bd = -2.0;
    for (int i = 0; i < (int)p.cells.size(); ++i) {
        double d = p.cells[i].unit.dot(dir);
        if (d > bd) { bd = d; best = i; }
    }
    return best;
}

int main() {
    PlanetConfig cfg; cfg.subdivisions = 5; cfg.seed = 1337;
    Planet planet; planet.generate(cfg);
    const int n = (int)planet.cells.size();

    std::printf("Ocean: moons\n");
    int nm = (int)planet.getMoons().size();
    std::printf("    moon count: %d\n", nm);
    check(nm >= 1 && nm <= 3, "1-3 moons generated");

    // Determinism: same seed -> identical moons (and tectonic state unchanged by moon RNG).
    Planet p2; p2.generate(cfg);
    bool sameMoons = ((int)p2.getMoons().size() == nm);
    for (int i = 0; i < nm && sameMoons; ++i) {
        const Moon& a = planet.getMoons()[i]; const Moon& b = p2.getMoons()[i];
        if (a.orbitRadius != b.orbitRadius || a.periodDays != b.periodDays || a.phase != b.phase
            || a.inclination != b.inclination || a.tideWeight != b.tideWeight) sameMoons = false;
    }
    check(sameMoons, "moon generation is deterministic for a seed");
    bool sameCells = true;
    for (int i = 0; i < n; ++i) if (planet.cells[i].elevation != p2.cells[i].elevation) sameCells = false;
    check(sameCells, "moon RNG does not perturb tectonic determinism");

    std::printf("Ocean: sky directions\n");
    Vec3 s0 = planet.sunDirection(0.25, 0.10), s1 = planet.sunDirection(0.25, 0.60);
    check(std::fabs(s0.length() - 1.0) < 1e-9, "sun direction is a unit vector");
    check((s0 - s1).length() > 0.1, "sun sweeps as the day advances (time of day)");
    if (nm > 0) {
        Vec3 m0 = planet.moonDirection(0, 0.10, 3.0), m1 = planet.moonDirection(0, 0.60, 3.0);
        check(std::fabs(m0.length() - 1.0) < 1e-9, "moon direction is a unit vector");
        check((m0 - m1).length() > 0.1, "moon sweeps across the sky with the clock");
    }

    std::printf("Ocean: tides (single equatorial moon, sun off)\n");
    // Isolate one equatorial moon so the two-bulge shape is unambiguous.
    Planet t; t.generate(cfg);
    t.cfg.tideSunFactor = 0.0;
    t.moons.clear();
    Moon m; m.orbitRadius = 20; m.periodDays = 27; m.phase = 0; m.inclination = 0; m.tideWeight = 1; m.dispRadius = 0.1;
    t.moons.push_back(m);
    t.computeTides(0.0, 0.25, 5.0);
    const std::vector<double> tide = t.tide();   // copy (t.tide() is overwritten on recompute)
    double sum = 0.0, mx = -1e30, mn = 1e30;
    for (int i = 0; i < n; ++i) { sum += tide[i]; mx = std::max(mx, tide[i]); mn = std::min(mn, tide[i]); }
    double mean = sum / n;
    check(std::fabs(mean) < 0.02, "tide field is ~zero-mean");
    check(mx > 0.0 && mn < 0.0, "tide has highs and lows (bulges + troughs)");
    Vec3 md = t.moonDirection(0, 0.25, 5.0);
    int sub = nearestCell(t, md);
    int anti = nearestCell(t, md * -1.0);
    int pole = nearestCell(t, Vec3{0, 1, 0});          // 90 deg from an equatorial moon
    check(tide[sub] > 0.0 && tide[anti] > 0.0, "high tide under the moon AND its antipode");
    check(tide[pole] < 0.0, "low tide at 90 degrees from the moon");

    // The bulge moves with the clock. (Use a quarter-ish step, not 0.5: half a day apart the
    // moon is at the antipode and the tide -- being cos^2 symmetric -- is correctly identical.)
    t.computeTides(0.0, 0.40, 5.0);
    bool moved = false;
    for (int i = 0; i < n; ++i) if (std::fabs(t.tide()[i] - tide[i]) > 1e-6) { moved = true; break; }
    check(moved, "the tidal bulge moves as time of day advances");

    std::printf("Ocean: currents\n");
    const std::vector<Vec3>& cur = planet.current();   // populated by generate()'s computeClimate
    check((int)cur.size() == n, "current field sized to the grid");
    bool tangent = true, landZero = true; int flowing = 0;
    for (int i = 0; i < n; ++i) {
        if (planet.cells[i].elevation <= cfg.seaLevel) {
            if (std::fabs(cur[i].dot(planet.cells[i].unit)) > 1e-6) tangent = false;
            if (cur[i].length() > 1e-9) ++flowing;
        } else if (cur[i].length() > 1e-12) landZero = false;
    }
    check(tangent, "ocean currents are tangent to the surface");
    check(landZero, "currents are zero on land");
    check(flowing > 50, "many ocean cells carry a current");

    // Climate feedback: identical terrain, current factor on vs off.
    PlanetConfig cfgOff = cfg; cfgOff.climateCurrentFactor = 0.0;
    Planet conOn;  conOn.generate(cfg);
    Planet conOff; conOff.generate(cfgOff);
    double maxAbs = 0.0, maxWarm = -1e9, maxCold = 1e9;
    for (int i = 0; i < n; ++i) {
        double d = conOn.temperature()[i] - conOff.temperature()[i];
        maxAbs = std::max(maxAbs, std::fabs(d));
        maxWarm = std::max(maxWarm, d); maxCold = std::min(maxCold, d);
    }
    check(maxAbs <= cfg.climateCurrentFactor + 1e-6, "current temp feedback bounded by climateCurrentFactor");
    check(maxWarm > 0.1 && maxCold < -0.1, "currents both warm and cool coasts");
    Planet conOn2; conOn2.generate(cfg);
    bool det = true;
    for (int i = 0; i < n; ++i) {
        if (conOn2.temperature()[i] != conOn.temperature()[i]) det = false;
        if ((conOn2.current()[i] - conOn.current()[i]).length() != 0.0) det = false;
    }
    check(det, "currents + feedback deterministic for a seed");

    std::printf("Ocean: enclosed-sea tide cap\n");
    {
        // A single isolated ocean cell (closed-off sea) must have a tiny tidal range.
        Planet e; e.generate(cfg);
        for (auto& cc : e.cells) cc.elevation = 100.0;       // all land
        e.cells[0].elevation = -100.0;                       // one-cell sea
        e.computeTides(0.0, 0.3, 4.0);
        check(std::fabs(e.tide()[0]) <= 0.04 + 1e-9, "a one-cell sea is capped to ~0.04 m");

        // A fully open ocean (one giant body) keeps the full equilibrium tidal range.
        Planet g; g.generate(cfg);
        for (auto& cc : g.cells) cc.elevation = -100.0;       // all ocean
        g.computeTides(0.0, 0.3, 4.0);
        double gmax = 0.0; for (int i = 0; i < n; ++i) gmax = std::max(gmax, std::fabs(g.tide()[i]));
        check(gmax > 0.1, "a large open ocean keeps the full tidal range");
    }

    std::printf("Ocean: save v9\n");
    std::stringstream ss(std::ios::in | std::ios::out | std::ios::binary);
    planet.writeState(ss);
    Planet q;
    bool ok = q.readState(ss, true, true, true);
    bool moonsRT = ok && (int)q.getMoons().size() == nm;
    for (int i = 0; i < nm && moonsRT; ++i)
        if (q.getMoons()[i].periodDays != planet.getMoons()[i].periodDays
            || q.getMoons()[i].tideWeight != planet.getMoons()[i].tideWeight) moonsRT = false;
    check(moonsRT, "save v9 round-trips the moons");

    std::printf(failures ? "\nSOME OCEAN CHECKS FAILED (%d)\n" : "\nALL OCEAN CHECKS PASSED\n", failures);
    return failures ? 1 : 0;
}