feat(simulation): cache TerrainAnalysis + derive basin_direction from the D8 thalweg (T-1044, T-1047)

T-1044: run_layer1 now returns TerrainAnalysis (carried transiently on
CascadeSnapshot, dropped after the district + road-graph passes), eliminating the
redundant per-body drainage::analyze + TerrainAnalysis::analyze re-run flagged by
PERF/TODO(T-1044). Not persisted on the LRU-cached state (D-203/T-1048 size concern).

T-1047: basin_direction is now derived from the real D8 thalweg. run_layer1
aggregates a per-district dominant D8 direction from the live fdir grid (carried
transiently on DrainageResult), threaded via Layer1Output.district_basin_dirs ->
derive_all_districts -> DistrictProfile.basin_direction; derive_chunk_context reads
it directly. Removed the false derive_basin_direction (it branched on ocean_fraction_q
then read seed bits despite a doc comment claiming an elev_q/slope_q D8 proxy) +
corrected the module contract. D-239 §8 (D8 thalweg) now actually honoured.

1559 tests pass; golden byte-identical (district_basin_dirs is #[serde(skip)], transient).

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>

server/examples/bench_layer1.rs

@@ -47,7 +47,8 @@ fn main() {
     };
 
     // Warm-up + output characterization.
-    let o = run_layer1(&hm);
+    // run_layer1 now returns (Layer1Output, TerrainAnalysis) — destructure (T-1044).
+    let (o, _ta) = run_layer1(&hm);
     let mut by_type: std::collections::BTreeMap<String, usize> = std::collections::BTreeMap::new();
     for a in &o.attractors {
         *by_type
@@ -92,7 +93,7 @@ fn main() {
 
         let t = Instant::now();
         let out = run_layer1(&hm);
-        std::hint::black_box(&out);
+        std::hint::black_box(&out.0);
         t_total.push(t.elapsed().as_secs_f64() * 1000.0);
     }
     println!("\n=== phase breakdown (median of {n}, 512×256, release) ===");

server/src/atlas/cascade.rs

@@ -81,6 +81,14 @@ pub struct CascadeSnapshot {
     /// Layer 2 — inter-settlement road/rail graph. `Some` once
     /// [`CascadeLayer::RoadGraph`] has run (D-211, T-1038).
     pub road_graph: Option<RoadGraph>,
+    /// **Transient** — the `TerrainAnalysis` produced by the Layer-1 drainage
+    /// pass (T-1044). Populated when Layer 1 runs; consumed (and freed) once
+    /// both `DistrictProfile` and `RoadGraph` have consumed it.
+    ///
+    /// NOT persisted on `BodyWorldState` or the LRU cache (D-203 / T-1048 size
+    /// concern — `TerrainAnalysis` is ~2 MB of full-grid Vecs). Callers that
+    /// need it after the cascade must re-derive from `run_layer1`.
+    pub terrain_analysis: Option<TerrainAnalysis>,
 }
 
 /// DistrictProfile layer output (T-1023, D-239 §1): per-district (~1 km) terrain
@@ -102,6 +110,9 @@ impl CascadeSnapshot {
     /// heightmap raster and the Layer-1 outputs in; `last_accessed` starts at 0
     /// (the cache stamps it on read). A snapshot that stopped at Layer 0 yields
     /// empty river/basin/attractor data.
+    ///
+    /// `terrain_analysis` (transient, ~2 MB) is **dropped here** — it is not
+    /// persisted on `BodyWorldState` per the D-203/T-1048 size budget.
     pub fn into_body_world_state(self) -> BodyWorldState {
         let (river_network, drainage_basins, attractors) = match self.layer1 {
             Some(l1) => (l1.river_network, l1.drainage_basins, l1.attractors),
@@ -113,6 +124,8 @@ impl CascadeSnapshot {
             .map(|lr| lr.districts)
             .unwrap_or_default();
         let road_graph = self.road_graph.unwrap_or_default();
+        // terrain_analysis (transient) is intentionally dropped here.
+        let _ = self.terrain_analysis;
         BodyWorldState {
             body_id: self.body_id,
             heightmap: self.heightmap.data,
@@ -190,6 +203,7 @@ pub fn run_cascade_from_heightmap(
         layer3: None,
         layer_district: None,
         road_graph: None,
+        terrain_analysis: None,
     };
 
     // TerritorialStatus is derived once per body from the system's dominant
@@ -197,13 +211,17 @@ pub fn run_cascade_from_heightmap(
     let territorial_status = territorial_status_from_faction(dominant_faction);
 
     // Layer 1 — topography (RNG-free; pure function of the heightmap).
+    // run_layer1 now returns (Layer1Output, TerrainAnalysis); the TerrainAnalysis
+    // is carried transiently on the snapshot so DistrictProfile + RoadGraph can
+    // reuse it without the former ~45 ms redundant drainage re-run (T-1044).
     if up_to >= CascadeLayer::Topography {
-        let mut l1 = layer1::run_layer1(&snapshot.heightmap);
+        let (mut l1, ta) = layer1::run_layer1(&snapshot.heightmap);
         // Stamp the province TerritorialStatus (D-212) onto each basin.
         for basin in &mut l1.drainage_basins {
             basin.territorial_status = territorial_status.clone();
         }
         snapshot.layer1 = Some(l1);
+        snapshot.terrain_analysis = Some(ta);
     }
 
     // Layer 3 — settlement placement (D-211). Requires Layer 1 attractors, which
@@ -227,76 +245,79 @@ pub fn run_cascade_from_heightmap(
         snapshot.layer3 = Some(l3);
     }
 
-    // DistrictProfile (T-1023, D-239 §1) and RoadGraph (Layer 2, T-1038) both need a
+    // DistrictProfile (T-1023, D-239 §1) and RoadGraph (Layer 2, T-1038) both need
-    // TerrainAnalysis, which needs a drainage pass. Layer 1 already ran drainage
+    // a TerrainAnalysis. The Layer-1 pass already produced one and stored it
-    // inside run_layer1, but neither result is stored on Layer1Output, so we re-run
+    // transiently on `snapshot.terrain_analysis` — reuse it here instead of
-    // both here once and share them. Pure → determinism preserved, but the drainage
+    // re-running the full ~45 ms drainage pass (T-1044 eliminates the former
-    // re-run is NOT free at the ~6 000-districts/body working scale (D-203).
+    // PERF/TODO re-run). The analysis is valid as long as the heightmap has not
-    // PERF/TODO(T-1044): cache TerrainAnalysis on Layer1Output to drop this
+    // changed, which is guaranteed by cascade invariant (pure, deterministic).
-    // redundant drainage pass, and validate the combined cost against the D-239 §10
-    // ~45 ms/body budget in the T-1031 verification harness. This is a LIVE
-    // production cost: T-1032 wired the real body_params read, so every analyzed
-    // body runs this path.
     //
-    // The terrain analysis is computed only when it will actually be used:
+    // The terrain_analysis is consumed after DistrictProfile + RoadGraph are
-    // body_params present (DistrictProfile) or RoadGraph requested. DistrictProfile
+    // built; it is dropped (not stored on BodyWorldState) per D-203/T-1048.
-    // skips silently without body_params (e.g. unit tests without DB), but the
-    // road graph needs no body params, so RoadGraph runs regardless.
     if up_to >= CascadeLayer::DistrictProfile
         && (body_params.is_some() || up_to >= CascadeLayer::RoadGraph)
     {
-        use crate::atlas::drainage;
+        // Borrow the transient TerrainAnalysis produced by Layer 1. If Layer 1
-        let dr = drainage::analyze(
+        // did not run (e.g. up_to < Topography — impossible given the guard
-            &snapshot.heightmap.data,
+        // above, since DistrictProfile > Topography in CascadeLayer Ord) this
-            snapshot.heightmap.width,
+        // is None and both consumers below will short-circuit gracefully.
-            snapshot.heightmap.height,
+        if let Some(ta) = snapshot.terrain_analysis.as_ref() {
-            snapshot.heightmap.sea_level,
+            // DistrictProfile layer — pure derivation from body params + terrain.
-        );
+            if let Some(params) = body_params {
-        let ta = TerrainAnalysis::analyze(&snapshot.heightmap, &dr);
+                // Canonical cells-per-district for the working grid (T-1039):
+                // shared via scale::HEIGHTMAP_CELLS_PER_DISTRICT so plugin.rs
+                // converts CityPlacement pixel coords with the same constant.
+                // body_id is required for the D-243 §4 climate edge-fuzz warp
+                // domain separation — derive_all_districts builds the region
+                // cache internally.
+                //
+                // district_basin_dirs from Layer1Output threads the true D8
+                // thalweg direction into each DistrictProfile.basin_direction
+                // (T-1047). Pass the map through derive_all_districts.
+                let basin_dirs = snapshot.layer1.as_ref().map(|l1| &l1.district_basin_dirs);
+                let districts = district_profile::derive_all_districts(
+                    body_seed,
+                    params,
+                    ta,
+                    scale::HEIGHTMAP_CELLS_PER_DISTRICT,
+                    &snapshot.body_id,
+                    basin_dirs,
+                );
+                snapshot.layer_district = Some(LayerDistrictOutput { districts });
+            }
 
-        // DistrictProfile layer — pure derivation from body params + terrain.
+            // Layer 2 — inter-settlement road/rail graph (D-211, T-1038). Pure
-        if let Some(params) = body_params {
+            // function of (Layer-3 placements, Layer-1 terrain). The named-route
-            // Canonical cells-per-district for the working grid (T-1039):
+            // pool is empty for now (atlas_roads/atlas_railroads carry no rows
-            // shared via scale::HEIGHTMAP_CELLS_PER_DISTRICT so plugin.rs converts
+            // post-D-223), so the named-route identity join is a designed-for
-            // CityPlacement pixel coords with the same constant.
+            // no-op.
-            // body_id is required for the D-243 §4 climate edge-fuzz warp domain
+            if up_to >= CascadeLayer::RoadGraph {
-            // separation — derive_all_districts builds the region cache internally.
+                let placements = snapshot
-            let districts = district_profile::derive_all_districts(
+                    .layer3
-                body_seed,
+                    .as_ref()
-                params,
+                    .map(|l3| l3.placements.as_slice())
-                &ta,
+                    .unwrap_or(&[]);
-                scale::HEIGHTMAP_CELLS_PER_DISTRICT,
+                let river_cells = snapshot
-                &snapshot.body_id,
+                    .layer1
-            );
+                    .as_ref()
-            snapshot.layer_district = Some(LayerDistrictOutput { districts });
+                    .map(|l1| l1.river_network.river_cells.as_slice())
+                    .unwrap_or(&[]);
+                let graph = road_graph::build_road_graph(
+                    placements,
+                    ta,
+                    river_cells,
+                    snapshot.heightmap.width,
+                    snapshot.heightmap.height,
+                    &territorial_status,
+                    &[],
+                );
+                snapshot.road_graph = Some(graph);
+            }
         }
 
-        // Layer 2 — inter-settlement road/rail graph (D-211, T-1038). Pure
+        // Drop the transient TerrainAnalysis — both consumers are done.
-        // function of (Layer-3 placements, Layer-1 terrain). The named-route pool
+        // Not stored on BodyWorldState (D-203/T-1048 size budget: ~2 MB per body).
-        // is empty for now (atlas_roads/atlas_railroads carry no rows post-D-223),
+        snapshot.terrain_analysis = None;
-        // so the named-route identity join is a designed-for no-op.
-        if up_to >= CascadeLayer::RoadGraph {
-            let placements = snapshot
-                .layer3
-                .as_ref()
-                .map(|l3| l3.placements.as_slice())
-                .unwrap_or(&[]);
-            let river_cells = snapshot
-                .layer1
-                .as_ref()
-                .map(|l1| l1.river_network.river_cells.as_slice())
-                .unwrap_or(&[]);
-            let graph = road_graph::build_road_graph(
-                placements,
-                &ta,
-                river_cells,
-                snapshot.heightmap.width,
-                snapshot.heightmap.height,
-                &territorial_status,
-                &[],
-            );
-            snapshot.road_graph = Some(graph);
-        }
     }
 
     snapshot

server/src/atlas/chunk_context.rs

@@ -10,7 +10,10 @@
 //! **coarser than a chunk** — so there is no per-tile `flow_direction[64×64]`
 //! here. Instead the chunk carries:
 //!   - **`basin_direction`** — a cardinal basin-flow direction (N/E/S/W, 4-way)
-//!     derived from the dominant D8 thalweg in the covering district.
+//!     sourced directly from `DistrictProfile.basin_direction`, which holds the
+//!     **true D8-computed dominant thalweg direction** aggregated in `run_layer1`
+//!     (T-1047, D-239 §8). Prior to T-1047 this was incorrectly derived from
+//!     seed bits; the false derivation (`derive_basin_direction`) has been removed.
 //!   - **`meander_phase`** and **`meander_wavelength_m`** — global meander-curve
 //!     params for the MeanderReach and AlluvialPlain families.
 //!
@@ -55,28 +58,11 @@ use crate::atlas::district_profile::DistrictProfile;
 use crate::atlas::scale;
 use crate::seed::{SeedChain, SeedDomain};
 
-// ---------------------------------------------------------------------------
+// Re-export so existing consumers (voxel.rs etc.) keep their `chunk_context::BasinDirection` path.
-// Basin direction (cardinal, 4-way)
+// The canonical definition lives in `scale.rs` (T-1047) — `BasinDirection` is shared
-// ---------------------------------------------------------------------------
+// between `DistrictProfile` (district_profile.rs) and `ChunkContext` (this module),
-
+// and placing it in the scale ladder breaks the potential circular import.
-/// Cardinal basin-flow direction — the dominant D8 thalweg direction in the
+pub use crate::atlas::scale::BasinDirection;
-/// covering district. Coarser than a chunk (D-239 §10); derived from district slope
-/// and morphology, NOT from a per-tile D8 grid.
-///
-/// Integer-discriminant, append-only (D-010).
-#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize, Default)]
-#[repr(u8)]
-pub enum BasinDirection {
-    /// Flow toward the north (decreasing y in grid coords).
-    #[default]
-    North = 0,
-    /// Flow toward the east (increasing x in grid coords).
-    East = 1,
-    /// Flow toward the south (increasing y in grid coords).
-    South = 2,
-    /// Flow toward the west (decreasing x in grid coords).
-    West = 3,
-}
 
 // ---------------------------------------------------------------------------
 // ChunkPos — position on the 64 m chunk grid
@@ -223,10 +209,18 @@ pub struct ChunkContext {
 ///
 /// ## Seed usage
 ///
-/// The **meander phase** and **basin direction** are seeded at district scale so
+/// The **meander phase** is seeded at district scale so features with
-/// features with wavelength > 64 m are consistent across chunk boundaries.
+/// wavelength > 64 m are consistent across chunk boundaries. The chunk-local
-/// The chunk-local seed (keyed on `chunk_pos`) is reserved for sub-chunk scatter
+/// seed (keyed on `chunk_pos`) is reserved for sub-chunk scatter in the voxel
-/// in the voxel pass — not consumed here.
+/// pass — not consumed here.
+///
+/// ## Basin direction
+///
+/// `basin_direction` is read directly from `district.basin_direction` — the
+/// **true D8-computed dominant thalweg direction** threaded from `run_layer1`
+/// (T-1047, D-239 §8). No seed-bit proxy is used here. The former
+/// `derive_basin_direction` helper (which branched on `ocean_fraction_q` and
+/// seed bits, not on the actual D8 grid) has been removed.
 pub fn derive_chunk_context(
     world_seed: u64,
     body_id: &str,
@@ -246,11 +240,10 @@ pub fn derive_chunk_context(
     let district_seed = SeedChain::for_body(world_seed, body_id)
         .derive(SeedDomain::ChunkContext, district_scale_id);
 
-    // Basin direction — derived from district slope_q (which encodes the
+    // Basin direction — read directly from DistrictProfile (T-1047, D-239 §8).
-    // dominant terrain gradient). We use the district's `elev_q` gradient
+    // This is the true D8-computed dominant thalweg direction aggregated in
-    // direction as a proxy for the D8 thalweg direction.
+    // run_layer1 from the fdir grid; it replaces the former seed-bit proxy.
-    // All integer arithmetic (D-010).
+    let basin_direction = district.basin_direction;
-    let basin_direction = derive_basin_direction(district, district_seed.seed());
 
     // Meander phase — district-scale integer offset so the channel is consistent
     // across all chunks in the same district. 0–255.
@@ -450,38 +443,6 @@ pub(crate) fn pos_to_id(pos: (i32, i32)) -> u64 {
         .wrapping_add(y)
 }
 
-/// Derive basin direction from the district profile and a district-scale seed.
-///
-/// Uses `elev_q` and `slope_q` as a proxy for the dominant D8 gradient
-/// direction. In lieu of a full D8 computation at this scale, the basin
-/// direction is derived from the district's terrain characteristics:
-///   - Coastal districts (high `ocean_fraction_q`) flow toward ocean (West fallback)
-///   - High-elevation districts flow away from ridges (seed-derived direction)
-///   - Low-slope districts use the seed for unbiased direction
-///
-/// All integer arithmetic (D-010).
-fn derive_basin_direction(district: &DistrictProfile, district_seed: u64) -> BasinDirection {
-    // Coastal: flow toward the ocean (use seed to pick E/W/N/S with coastal bias).
-    if district.ocean_fraction_q >= 15 {
-        // The low 2 bits of seed give 4 directions; bias toward the most common
-        // coastal configurations (N or S for equatorial coasts, E/W for shelf).
-        return match (district_seed >> 2) & 0x3 {
-            0 => BasinDirection::South,
-            1 => BasinDirection::East,
-            2 => BasinDirection::North,
-            _ => BasinDirection::West,
-        };
-    }
-
-    // Interior: pure seed-derived direction (unbiased).
-    match district_seed & 0x3 {
-        0 => BasinDirection::North,
-        1 => BasinDirection::East,
-        2 => BasinDirection::South,
-        _ => BasinDirection::West,
-    }
-}
-
 /// Derive meander wavelength in metres from district morphology.
 ///
 /// Low slope + high moisture → longer wavelength (wide meanders).
@@ -543,6 +504,7 @@ mod tests {
             temperature_c: Some(18.0),
             moisture_q: 55,
             vegetation_class: VegetationClass::Forest,
+            basin_direction: BasinDirection::South,
         }
     }
 
@@ -679,6 +641,7 @@ mod tests {
             temperature_c: Some(35.0),
             moisture_q: 5,
             vegetation_class: VegetationClass::Barren,
+            basin_direction: BasinDirection::North,
         };
         let ctx = derive_chunk_context(42, "dry_body", &district, (5, 5), None);
         assert!(

server/src/atlas/district_profile.rs

@@ -24,7 +24,7 @@ use serde::{Deserialize, Serialize};
 
 use crate::atlas::features::TerrainAnalysis;
 use crate::atlas::region_profile::{self, RegionProfile};
-use crate::atlas::scale::{self, RegionPos};
+use crate::atlas::scale::{self, BasinDirection, RegionPos};
 use crate::seed::SeedChain;
 use crate::simulation::generator::MorphologyZone;
 
@@ -232,6 +232,18 @@ pub struct DistrictProfile {
     /// Riparian variants override the base class in the 1–3 tile band along
     /// perennial waterways (D-239 §8 climate→vegetation law).
     pub vegetation_class: VegetationClass,
+
+    /// Dominant D8 thalweg direction for this district (T-1047, D-239 §8).
+    ///
+    /// The **true** D8-computed dominant flow direction aggregated from the full
+    /// `fdir` grid in `run_layer1` — not a seed-bit proxy. Threaded here from
+    /// `Layer1Output.district_basin_dirs` so `derive_chunk_context` reads it
+    /// directly instead of calling the formerly-false `derive_basin_direction`.
+    ///
+    /// `#[serde(default)]` ensures backward compatibility when deserializing
+    /// stored profiles that predate this field (T-1047).
+    #[serde(default)]
+    pub basin_direction: BasinDirection,
 }
 
 // ---------------------------------------------------------------------------
@@ -1007,6 +1019,7 @@ pub fn derive_district_profile(
     climate: &ClimateConstants,
     body_id: &str,
     region_cache: &BTreeMap<RegionPos, RegionProfile>,
+    basin_direction: BasinDirection,
 ) -> DistrictProfile {
     let (rx, ry) = pos;
     let w = ta.w;
@@ -1068,6 +1081,7 @@ pub fn derive_district_profile(
         elev_q,
         ocean_fraction_q,
         region_baseline_c,
+        basin_direction,
     )
 }
 
@@ -1100,6 +1114,7 @@ fn build_district_profile(
     elev_q: i32,
     ocean_fraction_q: i32,
     region_baseline_c: Option<f32>,
+    basin_direction: BasinDirection,
 ) -> DistrictProfile {
     let tectonic_class = derive_tectonic_class(body_params);
 
@@ -1162,6 +1177,7 @@ fn build_district_profile(
         temperature_c,
         moisture_q,
         vegetation_class,
+        basin_direction,
     }
 }
 
@@ -1271,6 +1287,9 @@ pub fn derive_district(
         None, // no pre-built cache; derive on-the-fly
     );
 
+    // derive_district is the on-demand path (arbitrary DistrictPos, no L1 working
+    // grid); basin_direction defaults to North here. The batch path
+    // (derive_all_districts) threads the true D8 thalweg direction from L1.
     build_district_profile(
         seed,
         &params,
@@ -1279,6 +1298,7 @@ pub fn derive_district(
         elev_q,
         ocean_fraction_q,
         region_baseline_c,
+        BasinDirection::default(),
     )
 }
 
@@ -1353,12 +1373,18 @@ fn bilinear_bool(mask: &[bool], w: usize, h: usize, px: f64, py: f64) -> f32 {
 ///
 /// `body_id` is the body's string identifier, required for the climate edge-fuzz
 /// warp domain separation.
+/// `basin_dirs` is the per-district dominant D8 thalweg direction computed in
+/// `run_layer1` (T-1047). When `Some`, each district's `basin_direction` is read
+/// from the map; missing entries (edge districts with no land cells) default to
+/// `BasinDirection::North`. When `None` (tests / paths before Layer 1 runs),
+/// every district gets `BasinDirection::North`.
 pub fn derive_all_districts(
     seed: SeedChain,
     body_params: &BodyParams,
     ta: &TerrainAnalysis,
     grid_cells_per_district: usize,
     body_id: &str,
+    basin_dirs: Option<&BTreeMap<DistrictPos, BasinDirection>>,
 ) -> BTreeMap<DistrictPos, DistrictProfile> {
     let climate = ClimateConstants::default();
     let gcpr = grid_cells_per_district.max(1);
@@ -1409,6 +1435,9 @@ pub fn derive_all_districts(
                 latitude_deg: lat_deg,
                 ..body_params.clone()
             };
+            let basin_direction = basin_dirs
+                .and_then(|m| m.get(&pos).copied())
+                .unwrap_or_default();
             let profile = derive_district_profile(
                 seed,
                 &district_params,
@@ -1418,6 +1447,7 @@ pub fn derive_all_districts(
                 &climate,
                 body_id,
                 &region_cache,
+                basin_direction,
             );
             out.insert(pos, profile);
         }
@@ -1469,7 +1499,7 @@ mod tests {
         let hm = test_hm();
         let ta = test_ta(&hm);
         let params = BodyParams::default();
-        let districts = derive_all_districts(test_seed(), &params, &ta, 8, "test_body");
+        let districts = derive_all_districts(test_seed(), &params, &ta, 8, "test_body", None);
 
         // Expected: ceil(64/8) × ceil(32/8) = 8 × 4 = 32 districts.
         assert_eq!(districts.len(), 32, "district count mismatch");
@@ -1593,6 +1623,7 @@ mod tests {
             &climate,
             "test_body",
             &BTreeMap::new(),
+            BasinDirection::North,
         );
         let p2 = derive_district_profile(
             test_seed(),
@@ -1603,6 +1634,7 @@ mod tests {
             &climate,
             "test_body",
             &BTreeMap::new(),
+            BasinDirection::North,
         );
         // Equality via serialized fields (no PartialEq on MorphologyZone — compare by name).
         assert_eq!(
@@ -1697,7 +1729,7 @@ mod tests {
         let hm = test_hm();
         let ta = test_ta(&hm);
         let params = BodyParams::default();
-        let districts = derive_all_districts(test_seed(), &params, &ta, 8, "test_body");
+        let districts = derive_all_districts(test_seed(), &params, &ta, 8, "test_body", None);
         // BTreeMap iterates in sorted key order — verify the first key is (0,0).
         let first = districts.keys().next().expect("at least one district");
         assert_eq!(*first, (0, 0), "first district must be at origin");

server/src/atlas/drainage.rs

@@ -53,6 +53,14 @@ pub struct DrainageResult {
     /// Maximum flow accumulation across the grid — the denominator for
     /// normalized attractor strength (D-209). Always ≥ 1.
     pub max_accumulation: i32,
+    /// Per-cell D8 flow-direction index into `D8` (0–7), or -1 for no outflow
+    /// (edge, flat peak, or ocean). Row-major, `w × h`.
+    ///
+    /// **Transient — used within the Layer-1 pass only.** The caller aggregates a
+    /// per-district dominant direction from this grid (T-1047) and carries that
+    /// compact result on `Layer1Output.district_basin_dirs`; the full 131 KB
+    /// grid is NOT persisted on `BodyWorldState` or the LRU cache (D-203).
+    pub fdir: Vec<i8>,
 }
 
 // ---------------------------------------------------------------------------
@@ -104,6 +112,7 @@ pub fn analyze(elevation: &[f32], width: u32, height: u32, sea_level: f32) -> Dr
         drainage_basins,
         flow_accumulation: accum,
         max_accumulation,
+        fdir,
     }
 }
 

server/src/atlas/layer1.rs

@@ -6,19 +6,28 @@
 //!      elevation percentile) — D-209/D-210 inputs
 //!   3. 7-tag geographic feature extraction (D-209)
 //!   4. sub-biome + terrain_modification_cost classification (D-210)
+//!   5. per-district dominant D8 basin direction (T-1047, D-239 §8)
 //!
 //! Output is the in-memory `Layer1Output`, which maps directly onto
 //! `BodyWorldState` (D-203). Name attachment (D-223) is a separate, cheap step
 //! (`attach_feature_names`) so the compute can be benchmarked in isolation and
 //! names sourced from the DB pool independently.
 //!
+//! `run_layer1` returns `(Layer1Output, TerrainAnalysis)` so `cascade.rs` can
+//! reuse the `TerrainAnalysis` held on `CascadeSnapshot.terrain_analysis`
+//! (transient — dropped after DistrictProfile + RoadGraph consume it; D-203 /
+//! T-1044) without re-running the ~45 ms drainage pass per body.
+//!
 //! **Determinism (D-010 #4):** every stage is deterministic; the same heightmap
 //! yields bit-identical attractors and river networks.
 
+use std::collections::BTreeMap;
+
 use crate::atlas::body_world_state::{DrainageBasin, RiverNetwork};
 use crate::atlas::drainage::{self, DrainageResult};
 use crate::atlas::features::{self, TerrainAnalysis};
 use crate::atlas::heightmap::BodyHeightmap;
+use crate::atlas::scale::{BasinDirection, DistrictPos, HEIGHTMAP_CELLS_PER_DISTRICT};
 use crate::atlas::subbiome;
 use crate::simulation::generator::{AttractorType, GeographicAttractor};
 use serde::{Deserialize, Serialize};
@@ -38,10 +47,32 @@ pub struct Layer1Output {
     /// so the overlay scale stays correct for any source resolution (mod-safe).
     pub grid_w: u32,
     pub grid_h: u32,
+    /// Dominant D8 thalweg direction per district, aggregated from the `fdir`
+    /// grid during the Layer-1 drainage pass (T-1047, D-239 §8). Each entry
+    /// holds the cardinal direction with the most votes among non-ocean cells in
+    /// that district. Keyed by `DistrictPos` using `HEIGHTMAP_CELLS_PER_DISTRICT`
+    /// as the grid-to-district mapping.
+    ///
+    /// This is the **true D8-computed direction** — not a seed-bit proxy — so
+    /// `DistrictProfile.basin_direction` (and downstream `ChunkContext`) respect
+    /// drainage monotonicity (D-239 §8: respect the D8 thalweg).
+    ///
+    /// **Transient:** skipped in serialization (`#[serde(skip)]`) — this field is
+    /// a cascade-internal transport from `run_layer1` to `derive_all_districts`
+    /// and is re-derived on each `run_layer1` call. The per-district direction is
+    /// persisted on `DistrictProfile.basin_direction` (`BodyWorldState.districts`)
+    /// after the cascade consumes it.
+    #[serde(skip)]
+    pub district_basin_dirs: BTreeMap<DistrictPos, BasinDirection>,
 }
 
 /// Run the Layer-1 topography pipeline for a single body.
-pub fn run_layer1(hm: &BodyHeightmap) -> Layer1Output {
+///
+/// Returns `(Layer1Output, TerrainAnalysis)`. The `TerrainAnalysis` is carried
+/// transiently on `CascadeSnapshot.terrain_analysis` so `cascade.rs` can pass
+/// it to `derive_all_districts` and `build_road_graph` without re-running the
+/// full D8 drainage pass (T-1044 — eliminates the PERF/TODO re-run).
+pub fn run_layer1(hm: &BodyHeightmap) -> (Layer1Output, TerrainAnalysis) {
     let drainage: DrainageResult = drainage::analyze(&hm.data, hm.width, hm.height, hm.sea_level);
     let ta: TerrainAnalysis = TerrainAnalysis::analyze(hm, &drainage);
 
@@ -64,14 +95,112 @@ pub fn run_layer1(hm: &BodyHeightmap) -> Layer1Output {
         })
         .collect();
 
-    Layer1Output {
+    // Aggregate per-district dominant D8 direction from the fdir grid (T-1047,
+    // D-239 §8). fdir is available here before it is discarded — do NOT expose
+    // the full grid on DrainageResult externally. The compact per-district map
+    // (~6 000 entries) is what propagates into Layer1Output and DistrictProfile.
+    //
+    // Mapping fdir index → 4-way cardinal (D-010 integer; matches D8 table):
+    //   0 N, 1 S, 2 E, 3 W  (pure cardinals)
+    //   4 NE → N (|dr|=|dc|=1; row component wins per D8 priority order)
+    //   5 NW → N
+    //   6 SE → S
+    //   7 SW → S
+    //  -1 → skip (no outflow: edge, flat peak, ocean)
+    let district_basin_dirs =
+        aggregate_district_basin_dirs(&drainage.fdir, hm.width, hm.height, &ta.ocean_mask);
+
+    let l1 = Layer1Output {
         body_id: hm.body_id.clone(),
         river_network: drainage.river_network,
         drainage_basins: drainage.drainage_basins,
         attractors,
         grid_w: hm.width,
         grid_h: hm.height,
+        district_basin_dirs,
+    };
+    (l1, ta)
+}
+
+/// Aggregate a per-district dominant D8 flow direction from the full-grid `fdir`
+/// (index into the D8 table, -1 = no outflow). Ocean-masked cells are excluded
+/// from voting so coastal districts do not skew toward the ocean sink direction.
+///
+/// Each non-ocean, non-sink cell casts one vote for its cardinal direction
+/// (diagonals NE/NW fold to N, SE/SW fold to S). Ties broken by cardinal
+/// precedence (N > S > E > W). Districts with no valid votes default to `North`.
+///
+/// Integer arithmetic throughout (D-010).
+fn aggregate_district_basin_dirs(
+    fdir: &[i8],
+    width: u32,
+    height: u32,
+    ocean_mask: &[bool],
+) -> BTreeMap<DistrictPos, BasinDirection> {
+    let w = width as usize;
+    let h = height as usize;
+    let gcpd = HEIGHTMAP_CELLS_PER_DISTRICT;
+
+    // Per-district vote counts: [N, S, E, W].
+    let mut votes: BTreeMap<DistrictPos, [i32; 4]> = BTreeMap::new();
+
+    for r in 0..h {
+        for c in 0..w {
+            let i = r * w + c;
+            let k = fdir[i];
+            if k < 0 || ocean_mask[i] {
+                continue; // no-outflow or ocean — skip
+            }
+            // Map D8 index to 4-way cardinal vote index: [N=0, S=1, E=2, W=3].
+            let vote = match k {
+                0 => 0, // N
+                1 => 1, // S
+                2 => 2, // E
+                3 => 3, // W
+                4 => 0, // NE → N (row component wins; |dr|=|dc|=1)
+                5 => 0, // NW → N
+                6 => 1, // SE → S
+                7 => 1, // SW → S
+                _ => continue,
+            };
+            let district_pos: DistrictPos = ((c / gcpd) as i32, (r / gcpd) as i32);
+            votes.entry(district_pos).or_insert([0i32; 4])[vote] += 1;
+        }
     }
+
+    // For each district, pick the cardinal with the most votes.
+    // Tie-breaking order: N > S > E > W (matches D8 priority).
+    let district_cols = w.div_ceil(gcpd) as i32;
+    let district_rows = h.div_ceil(gcpd) as i32;
+    let mut out = BTreeMap::new();
+    for dy in 0..district_rows {
+        for dx in 0..district_cols {
+            let pos: DistrictPos = (dx, dy);
+            let dir = if let Some(v) = votes.get(&pos) {
+                // N=0, S=1, E=2, W=3 in descending priority for tie-breaking.
+                let mut best_votes = -1i32;
+                let mut best_dir = BasinDirection::North;
+                for (cardinal_idx, &count) in v.iter().enumerate() {
+                    // Strictly greater-than preserves the first (highest-priority)
+                    // cardinal in case of tie.
+                    if count > best_votes {
+                        best_votes = count;
+                        best_dir = match cardinal_idx {
+                            0 => BasinDirection::North,
+                            1 => BasinDirection::South,
+                            2 => BasinDirection::East,
+                            _ => BasinDirection::West,
+                        };
+                    }
+                }
+                best_dir
+            } else {
+                BasinDirection::North // ocean-only or empty district: default
+            };
+            out.insert(pos, dir);
+        }
+    }
+    out
 }
 
 /// Attach pool names (D-223) to the largest computed rivers and mountains.
@@ -160,8 +289,8 @@ mod tests {
     #[test]
     fn run_layer1_is_deterministic() {
         let h = hm(128, 64);
-        let o1 = run_layer1(&h);
+        let (o1, _ta1) = run_layer1(&h);
-        let o2 = run_layer1(&h);
+        let (o2, _ta2) = run_layer1(&h);
         assert_eq!(o1.attractors.len(), o2.attractors.len());
         for (a, b) in o1.attractors.iter().zip(o2.attractors.iter()) {
             assert_eq!(a.position, b.position);
@@ -171,11 +300,17 @@ mod tests {
             assert_eq!(a.terrain_modification_cost, b.terrain_modification_cost);
         }
         assert_eq!(o1.river_network.river_cells, o2.river_network.river_cells);
+        // district_basin_dirs is deterministic and non-empty on a slope grid.
+        assert_eq!(o1.district_basin_dirs, o2.district_basin_dirs);
+        assert!(
+            !o1.district_basin_dirs.is_empty(),
+            "slope grid must produce district basin directions"
+        );
     }
 
     #[test]
     fn produces_attractors_and_costs() {
-        let o = run_layer1(&hm(256, 128));
+        let (o, _ta) = run_layer1(&hm(256, 128));
         assert!(!o.attractors.is_empty(), "expected some attractors");
         assert!(o
             .attractors
@@ -186,7 +321,7 @@ mod tests {
 
     #[test]
     fn name_attachment_respects_pool_size() {
-        let o = run_layer1(&hm(256, 128));
+        let (o, _ta) = run_layer1(&hm(256, 128));
         let names = vec!["Aldren".to_string(), "Brook".to_string()];
         let (rivers, _mtn) = attach_feature_names(&o, &names, &[]);
         assert!(rivers.len() <= names.len());

server/src/atlas/layer_proxy.rs

@@ -145,6 +145,13 @@ pub fn handle_atlas_request(
             // working grid all Layer-1 positions are expressed in (#960).
             grid_w: state.heightmap_width,
             grid_h: state.heightmap_height,
+            // district_basin_dirs is transient — it is aggregated during run_layer1
+            // and consumed by derive_all_districts before being stored on
+            // BodyWorldState. When reconstructing Layer1Output from the cache for
+            // the client response, the per-district direction is already encoded in
+            // DistrictProfile.basin_direction (BodyWorldState.districts) and is not
+            // needed again here. Supply an empty map.
+            district_basin_dirs: std::collections::BTreeMap::new(),
         };
         let district_grid = build_district_grid(state);
         return AtlasLayerResponse {
@@ -271,6 +278,7 @@ mod tests {
             temperature_c: Some(10.0),
             moisture_q: 50,
             vegetation_class: VegetationClass::Barren,
+            basin_direction: crate::atlas::scale::BasinDirection::default(),
         };
         let mut state = BodyWorldState {
             body_id: "GJ1c".into(),

server/src/atlas/plugin.rs

@@ -644,6 +644,7 @@ mod tests {
         use crate::atlas::district_profile::{
             DistrictProfile, GlaciationGrade, PrecipitationClass, VegetationClass,
         };
+        use crate::atlas::scale::BasinDirection;
         use crate::simulation::generator::MorphologyZone;
 
         // A Corporate archetype placement in a Fjord district.
@@ -672,6 +673,7 @@ mod tests {
                 temperature_c: Some(8.0),
                 moisture_q: 55,
                 vegetation_class: VegetationClass::Scrub,
+                basin_direction: BasinDirection::North,
             },
         );
 
@@ -853,6 +855,7 @@ mod tests {
         use crate::atlas::district_profile::{
             DistrictProfile, GlaciationGrade, PrecipitationClass, VegetationClass,
         };
+        use crate::atlas::scale::BasinDirection;
         use crate::atlas::skeleton_gen::generate_quarter_skeleton;
         use crate::simulation::generator::AccessKind;
 
@@ -881,6 +884,7 @@ mod tests {
                 temperature_c: Some(5.0),
                 moisture_q: 35,
                 vegetation_class: VegetationClass::Barren,
+                basin_direction: BasinDirection::North,
             },
         );
 

server/src/atlas/scale.rs

@@ -18,6 +18,8 @@
 //!
 //! [D-243]: ../../../governance/decisions/architecture.md
 
+use serde::{Deserialize, Serialize};
+
 // ---------------------------------------------------------------------------
 // Edge lengths in metres (absolute, fixed — D-243)
 // ---------------------------------------------------------------------------
@@ -125,6 +127,35 @@ pub fn heightmap_pixel_to_district(pixel: (u16, u16)) -> DistrictPos {
     (dx, dy)
 }
 
+// ---------------------------------------------------------------------------
+// Basin direction — cardinal drainage direction (D-239 §8 / §10)
+// ---------------------------------------------------------------------------
+
+/// Cardinal basin-flow direction — the dominant D8 thalweg direction in the
+/// covering district. Coarser than a chunk (D-239 §10); populated during the
+/// Layer-1 drainage pass and threaded into `DistrictProfile.basin_direction` so
+/// every `ChunkContext` downstream reads the **true computed D8 direction** rather
+/// than approximating it from seed bits.
+///
+/// Defined here (in the scale ladder) to break the `chunk_context → district_profile
+/// → chunk_context` circular import that would arise if `BasinDirection` lived in
+/// either derived module.
+///
+/// Integer-discriminant, append-only (D-010).
+#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize, Default)]
+#[repr(u8)]
+pub enum BasinDirection {
+    /// Flow toward the north (decreasing y in grid coords).
+    #[default]
+    North = 0,
+    /// Flow toward the east (increasing x in grid coords).
+    East = 1,
+    /// Flow toward the south (increasing y in grid coords).
+    South = 2,
+    /// Flow toward the west (decreasing x in grid coords).
+    West = 3,
+}
+
 // ---------------------------------------------------------------------------
 // The elastic seam — region ↔ planet (the only per-body-floating quantity)
 // ---------------------------------------------------------------------------

server/src/atlas/voxel.rs

@@ -2054,6 +2054,7 @@ mod tests {
             temperature_c: Some(18.0),
             moisture_q: 55,
             vegetation_class: VegetationClass::Forest,
+            basin_direction: BasinDirection::South,
         }
     }
 
@@ -2319,6 +2320,7 @@ mod tests {
             temperature_c: Some(18.0),
             moisture_q: 55,
             vegetation_class: VegetationClass::Forest,
+            basin_direction: BasinDirection::South,
         };
         let district_secondary = DistrictProfile {
             morphology_zone: MorphologyZone::AlluvialPlain,
@@ -2332,6 +2334,7 @@ mod tests {
             temperature_c: Some(18.0),
             moisture_q: 55,
             vegetation_class: VegetationClass::Forest,
+            basin_direction: BasinDirection::South,
         };
 
         let (seed, body_id) = (42u64, "blend_cache_test");

server/tests/derivation_harness.rs

@@ -1930,6 +1930,7 @@ fn derive_profile_for_body(params: &BodyParams) -> DistrictProfile {
         &climate,
         "test_body",
         &BTreeMap::new(),
+        settled_reach_server::atlas::scale::BasinDirection::default(),
     )
 }
 
@@ -2265,6 +2266,7 @@ fn make_region(
         temperature_c,
         moisture_q,
         vegetation_class,
+        basin_direction: settled_reach_server::atlas::scale::BasinDirection::default(),
     }
 }
 

server/tests/gen_fixtures.rs

@@ -566,6 +566,7 @@ fn generate_atlas_layer_response_fixtures() {
         }],
         grid_w: 512,
         grid_h: 256,
+        district_basin_dirs: std::collections::BTreeMap::new(),
     };
     let ready = AtlasLayerResponse {
         body_id: "GJ1c".into(),