Skip to main content

ras_commander_fast_map/
sloping_session.rs

1//! Explicit sloping batch path with independently selected native rendering controls. Geometry, terrain and ownership are shared
2//! across profiles; hydraulic connectivity is reconstructed for every profile.
3use crate::retained::Retained;
4use crate::{
5    MapError,
6    model::{self, CrsIdentity, NodataPolicy, QueryConvention},
7    sink,
8    sloping::{self, Cell, Corner, Face, Mesh},
9    terrain::NativeTerrain,
10};
11use rayon::prelude::*;
12use serde::Deserialize;
13use serde_json::{Value, json};
14use std::{
15    collections::HashSet,
16    fs,
17    io::Write,
18    path::Path,
19    sync::{Arc, Mutex},
20    time::Instant,
21};
22
23type CachedGeometry = (Vec<Mesh>, Vec<Vec<Option<sloping::CellStencil>>>);
24type WindowKey = (
25    String,
26    String,
27    [usize; 4],
28    Vec<[usize; 2]>,
29    Vec<u8>,
30    Option<Vec<u8>>,
31);
32struct WindowData {
33    owner: Vec<u32>,
34    ground: Vec<f32>,
35    weights: Option<sloping::PixelWeights>,
36}
37
38/// One explicitly owned, serialized native session. One complete geometry and
39/// terrain share bounded LRU window and hydraulic caches. Dropping releases all.
40#[derive(Default)]
41pub struct PreparedSession {
42    terrain: Option<(String, Arc<NativeTerrain>)>,
43    geometry: Option<(String, Arc<CachedGeometry>)>,
44    windows: Retained<WindowKey, WindowData>,
45    hydraulics: Mutex<Retained<String, Vec<sloping::PreparedProfile>>>,
46    pool: Option<(usize, Arc<rayon::ThreadPool>)>,
47}
48
49fn identity(path: &Path) -> Result<String, MapError> {
50    let path = path.canonicalize()?;
51    let metadata = path.metadata()?;
52    Ok(format!(
53        "{:?}:{}:{:?}",
54        path,
55        metadata.len(),
56        metadata.modified()?
57    ))
58}
59fn geometry_key(req: &Request, base: &Path) -> Result<String, MapError> {
60    let mut key = String::new();
61    for filename in &req.mesh_files {
62        let path = model::resolve(base, filename);
63        key.push_str(&identity(&path)?);
64        let document: Value = serde_json::from_reader(fs::File::open(&path)?)
65            .map_err(|e| MapError::Invalid(e.to_string()))?;
66        if let Some(arrays) = document["arrays"].as_object() {
67            for descriptor in arrays.values() {
68                if let Some(name) = descriptor["path"].as_str() {
69                    key.push_str(&identity(&model::resolve(path.parent().unwrap(), name))?);
70                }
71            }
72        }
73    }
74    Ok(key)
75}
76
77pub const SCHEMA: &str = "ras-commander-sloping-request/1";
78
79#[derive(Deserialize)]
80#[serde(deny_unknown_fields)]
81pub struct Request {
82    pub schema: String,
83    #[serde(default)]
84    pub options: sloping::Options,
85    pub terrain: String,
86    pub epsg: u32,
87    pub crs_identity: CrsIdentity,
88    pub terrain_transform: [f64; 6],
89    pub window: [usize; 4],
90    pub geometry_file: String,
91    #[serde(default)]
92    pub aoi_file: Option<String>,
93    pub mesh_files: Vec<String>,
94    pub owner_cells: Vec<[usize; 2]>,
95    pub profiles: Vec<Profile>,
96    pub threads: usize,
97    #[serde(default = "one_profile")]
98    pub max_inflight_profiles: usize,
99    pub max_working_bytes: u64,
100    #[serde(default)]
101    pub stream_outputs: bool,
102    /// Retained multi-window ownership/terrain/weight ceiling. Zero disables it.
103    #[serde(default = "window_budget")]
104    pub window_cache_bytes: u64,
105    /// Retained reconstructed timestep ceiling. Zero disables it.
106    #[serde(default)]
107    pub hydraulic_cache_bytes: u64,
108    /// Write native LZW GeoTIFF products and stream only their completion records.
109    #[serde(default)]
110    pub geotiff_outputs: bool,
111    /// Per mesh: HDF-ordered float32 (x, y, polyline length) face normals.
112    /// With it, every profile supplies `face_velocity_files` and a velocity
113    /// magnitude product is rendered (method `rasmapper-sloping-velocity/0`).
114    #[serde(default)]
115    pub velocity_normals_files: Option<Vec<String>>,
116
117    /// Declares `terrain` as the top-priority source of a registered terrain.
118    #[serde(default)]
119    pub terrain_scope: Option<crate::terrain_scope::TerrainScope>,
120}
121
122fn one_profile() -> usize {
123    1
124}
125fn window_budget() -> u64 {
126    64 << 20
127}
128
129fn admit_profiles(
130    requested: usize,
131    threads: usize,
132    profiles: usize,
133    budget: u64,
134    serial: u64,
135    additional: u64,
136) -> Result<(usize, u64), MapError> {
137    if serial > budget || additional == 0 || requested == 0 || threads == 0 || profiles == 0 {
138        return Err(MapError::Budget(
139            "sloping decoder plan exceeds budget".into(),
140        ));
141    }
142    let limit = requested.min(threads).min(profiles);
143    let extra = ((budget - serial) / additional).min((limit - 1) as u64);
144    Ok((extra as usize + 1, serial + extra * additional))
145}
146
147fn map_batch<T: Sync, R: Send>(
148    batch: &[T],
149    pool: &rayon::ThreadPool,
150    compute: impl Fn(&T) -> Result<R, MapError> + Sync + Send,
151) -> Result<Vec<R>, MapError> {
152    pool.install(|| batch.par_iter().map(compute).collect())
153}
154
155struct RenderedProfile {
156    depth: Vec<f32>,
157    wse: Vec<f32>,
158    velocity: Option<Vec<f32>>,
159    hydraulic_seconds: f64,
160    render_seconds: f64,
161    hydraulic_reused: bool,
162    finite_pixels: usize,
163    flooded_pixels: usize,
164}
165
166#[allow(clippy::too_many_arguments)] // Explicit immutable geometry/window/cache borrows.
167fn render(
168    profile: &Profile,
169    base: &Path,
170    meshes: &[Mesh],
171    stencils: &[Vec<Option<sloping::CellStencil>>],
172    window: &sloping::RasterWindow<'_>,
173    pool: &rayon::ThreadPool,
174    options: sloping::Options,
175    weights: Option<&sloping::PixelWeights>,
176    normals: Option<&[Vec<crate::velocity::FaceNormal>]>,
177    foreign: &[usize],
178    hydraulic_cache: Option<&Mutex<Retained<String, Vec<sloping::PreparedProfile>>>>,
179) -> Result<RenderedProfile, MapError> {
180    let tick = Instant::now();
181    let cache_key = if hydraulic_cache.is_some() && normals.is_none() {
182        let mut key =
183            serde_json::to_string(&options).map_err(|e| MapError::Invalid(e.to_string()))?;
184        for file in &profile.values_files {
185            key.push_str(&identity(&model::resolve(base, file))?);
186        }
187        Some(key)
188    } else {
189        None
190    };
191    let cached = match (&cache_key, hydraulic_cache) {
192        (Some(key), Some(cache)) => cache.lock().unwrap().get(key),
193        _ => None,
194    };
195    let hydraulic_reused = cached.is_some();
196    let mut new_prepared = Vec::new();
197    let mut cell_values = Vec::new();
198    if !hydraulic_reused {
199        for ((mesh, stencil), filename) in meshes.iter().zip(stencils).zip(&profile.values_files) {
200            let wse = values(&model::resolve(base, filename), mesh.cells.len())?;
201            new_prepared.push(sloping::prepare_profile_with_options(
202                mesh, stencil, &wse, 0.001, pool, options,
203            )?);
204            if normals.is_some() {
205                cell_values.push(wse);
206            }
207        }
208    }
209    let prepared = cached.unwrap_or_else(|| Arc::new(new_prepared));
210    if !hydraulic_reused && let (Some(key), Some(cache)) = (cache_key, hydraulic_cache) {
211        let bytes = prepared.iter().map(|p| p.retained_bytes()).sum::<u64>()
212            + (prepared.capacity() * std::mem::size_of::<sloping::PreparedProfile>()) as u64
213            + key.capacity() as u64
214            + (2 * std::mem::size_of::<(String, Arc<Vec<sloping::PreparedProfile>>, u64)>()
215                + std::mem::size_of::<Vec<sloping::PreparedProfile>>()
216                + 2 * std::mem::size_of::<usize>()) as u64;
217        cache.lock().unwrap().insert(key, prepared.clone(), bytes);
218    }
219    crate::terrain_scope::deny_wet_foreign(foreign, window.owner, window.owner_cells, |m, c| {
220        prepared[m].is_wet(c)
221    })?;
222    let hydraulic_seconds = tick.elapsed().as_secs_f64();
223    let tick = Instant::now();
224    let (mut depth, mut wse) =
225        sloping::render_profile_cached(stencils, &prepared, window, 0.001, pool, weights)?;
226    let render_seconds = tick.elapsed().as_secs_f64();
227    let (finite_pixels, flooded_pixels) = crate::finalize::apply(&mut depth, &mut wse, pool);
228    let velocity = match (normals, &profile.face_velocity_files) {
229        (Some(normals), Some(files)) => {
230            let mut face_velocity = Vec::new();
231            for (mesh, file) in meshes.iter().zip(files) {
232                face_velocity.push(crate::velocity::read_face_velocity(
233                    &model::resolve(base, file),
234                    mesh.faces.len(),
235                )?);
236            }
237            let inputs: Vec<_> = meshes
238                .iter()
239                .zip(stencils)
240                .zip(&cell_values)
241                .zip(normals)
242                .zip(&face_velocity)
243                .map(|((((mesh, stencil), wse), normal), velocity)| {
244                    (
245                        mesh,
246                        stencil.as_slice(),
247                        wse.as_slice(),
248                        normal.as_slice(),
249                        velocity.as_slice(),
250                    )
251                })
252                .collect();
253            Some(crate::velocity::render(&inputs, window, &depth, pool)?)
254        }
255        _ => None,
256    };
257    Ok(RenderedProfile {
258        depth,
259        wse,
260        velocity,
261        hydraulic_seconds,
262        render_seconds,
263        hydraulic_reused,
264        finite_pixels,
265        flooded_pixels,
266    })
267}
268
269#[derive(Deserialize)]
270#[serde(deny_unknown_fields)]
271pub struct Profile {
272    pub id: String,
273    pub values_files: Vec<String>,
274    /// Per mesh: HDF-ordered float32 face-normal velocities of this profile.
275    #[serde(default)]
276    pub face_velocity_files: Option<Vec<String>>,
277}
278
279#[derive(Deserialize)]
280pub struct Geometry {
281    #[serde(default)]
282    pub vertex_elevations: Option<Vec<f32>>,
283    pub cells: Vec<InputCell>,
284    pub faces: Vec<InputFace>,
285    pub vertices: Vec<InputVertex>,
286}
287#[derive(Deserialize)]
288pub struct InputCell {
289    pub point: [f64; 2],
290    pub minimum_elevation: f32,
291    #[serde(rename = "virtual")]
292    pub virtual_cell: bool,
293    pub faces: Vec<usize>,
294}
295#[derive(Deserialize)]
296pub struct InputFace {
297    pub cells: [usize; 2],
298    pub vertices: [usize; 2],
299    pub minimum_elevation: f32,
300    pub perimeter: bool,
301    pub polyline: Vec<[f64; 2]>,
302}
303#[derive(Deserialize)]
304pub struct InputVertex {
305    pub point: [f64; 2],
306    pub faces: Vec<usize>,
307    pub ccw: Vec<usize>,
308    pub cw: Vec<usize>,
309}
310
311impl Geometry {
312    pub fn prepare(&self) -> Result<(Mesh, Vec<Option<sloping::CellStencil>>), MapError> {
313        if self
314            .cells
315            .iter()
316            .any(|c| c.virtual_cell != (c.faces.len() == 1))
317        {
318            return Err(MapError::Invalid(
319                "sloping virtual-cell flag contradicts topology".into(),
320            ));
321        }
322        let faces =
323            self.faces
324                .iter()
325                .map(|f| {
326                    let a = self.cells.get(f.cells[0]).ok_or_else(|| {
327                        MapError::Invalid("sloping cell address out of range".into())
328                    })?;
329                    let b = self.cells.get(f.cells[1]).ok_or_else(|| {
330                        MapError::Invalid("sloping cell address out of range".into())
331                    })?;
332                    Ok(Face {
333                        cells: f.cells,
334                        vertices: f.vertices,
335                        minimum: f.minimum_elevation,
336                        perimeter: f.perimeter,
337                        application: sloping::face_application(&f.polyline, [a.point, b.point])?,
338                    })
339                })
340                .collect::<Result<Vec<_>, MapError>>()?;
341        let mesh = Mesh {
342            cells: self
343                .cells
344                .iter()
345                .map(|c| Cell {
346                    minimum: c.minimum_elevation,
347                    virtual_cell: c.virtual_cell,
348                })
349                .collect(),
350            faces,
351            corners: self
352                .vertices
353                .iter()
354                .map(|v| Corner {
355                    point: v.point,
356                    faces: v.faces.clone(),
357                    ccw: v.ccw.clone(),
358                    cw: v.cw.clone(),
359                })
360                .collect(),
361        };
362        let cell_faces = self
363            .cells
364            .iter()
365            .map(|c| c.faces.clone())
366            .collect::<Vec<_>>();
367        let mut stencils = sloping::prepare_stencils(&mesh, &cell_faces)?;
368        for (cell, stencil) in stencils.iter_mut().enumerate() {
369            if let Some(stencil) = stencil {
370                for &face in &self.cells[cell].faces {
371                    let f = &self.faces[face];
372                    if f.cells[0] == cell {
373                        stencil.physical_ring.extend(f.polyline.iter().copied());
374                    } else {
375                        stencil
376                            .physical_ring
377                            .extend(f.polyline.iter().rev().copied());
378                    }
379                }
380                stencil.physical_ring.dedup();
381            }
382        }
383        if let Some(elevations) = &self.vertex_elevations {
384            sloping::set_vertex_elevations(&mut stencils, elevations, mesh.corners.len())?;
385        }
386        Ok((mesh, stencils))
387    }
388}
389
390pub fn method_record() -> Value {
391    json!({"method_id":"rasmapper-plain-sloping/1", "backend":"rust", "experimental":true,
392        "native_version_qualified":"7.0.1", "stream_outputs":true, "native_ownership":true,"aoi_mask":crate::aoi::METHOD,
393        "enhanced_options":true,"enhanced_method_id":"rasmapper-sloping-with-faces/1",
394        "depth_weighted_faces":true,
395        "native_donation_history":true,
396        "bounded_profile_concurrency":true,
397        "persistent_session":true,"multi_window_cache":true,"hydraulic_cache":true,"native_geotiff_outputs":true,"bounded_pixel_weights":true,
398        "velocity":crate::velocity::METHOD_ID,
399
400        "terrain_scope":crate::terrain_scope::METHOD_ID,
401        "options":{"CellRenderMode":"Sloping","CellStencilMode":"JustFacepoints",
402            "FacepointAdjustmentMode":"None","ShallowBehaviorMode":"None",
403            "UseDepthWeightedFaces":false,"FaceWSMode":"Adjusted"},
404        "native_plot_tolerance":0.001_f32,"extra_depth_cutoff":null,
405        "scope":"2D cells; plain cell-corner sloping; no 1D, slopingPretty or additional rendering options",
406        "reuse":"geometry, terrio ownership and terrain prepared once; bounded Rayon pixel and reconstruction workers"})
407}
408
409fn selected_method(options: sloping::Options) -> Value {
410    let mut record = method_record();
411    if options.render_mode == sloping::RenderMode::WithFaces {
412        record["method_id"] = "rasmapper-sloping-with-faces/1".into();
413        record["options"]["CellStencilMode"] = "WithFaces".into();
414        record["options"]["ShallowBehaviorMode"] = if options.reduce_shallow_to_horizontal {
415            "ReduceToHorizontal"
416        } else {
417            "None"
418        }
419        .into();
420        record["options"]["UseDepthWeightedFaces"] = options.use_depth_weighted_faces.into();
421        record["scope"] = "physical 2D cells; cell-corner and face sloping; no 1D or standalone storage-area surfaces".into();
422    }
423    record
424}
425
426fn values(path: &Path, count: usize) -> Result<Vec<f32>, MapError> {
427    let expected = count
428        .checked_mul(4)
429        .ok_or_else(|| MapError::Budget("sloping values size overflow".into()))?;
430    if fs::metadata(path)?.len() != expected as u64 {
431        return Err(MapError::Invalid("sloping profile length mismatch".into()));
432    }
433    let raw = fs::read(path)?;
434    if raw.len() != expected {
435        return Err(MapError::Invalid("sloping profile changed size".into()));
436    }
437    let values: Vec<_> = raw
438        .as_chunks::<4>()
439        .0
440        .iter()
441        .map(|v| f32::from_le_bytes(*v))
442        .collect();
443    if values.iter().any(|v| !v.is_finite()) {
444        return Err(MapError::Invalid(
445            "sloping profiles require finite cell WSE including virtual cells".into(),
446        ));
447    }
448    Ok(values)
449}
450
451pub fn run(req: &Request, base: &Path, out: &Path) -> Result<Value, MapError> {
452    run_with_session(req, base, out, None)
453}
454
455pub fn run_with_session(
456    req: &Request,
457    base: &Path,
458    out: &Path,
459    mut session: Option<&mut PreparedSession>,
460) -> Result<Value, MapError> {
461    req.options.validate()?;
462    if req.geotiff_outputs && (!req.stream_outputs || req.velocity_normals_files.is_some()) {
463        return Err(MapError::Unsupported(
464            "native GeoTIFF requires stream_outputs and no velocity".into(),
465        ));
466    }
467    if let Some(scope) = &req.terrain_scope {
468        scope.validate()?;
469        if req.options.render_mode == sloping::RenderMode::WithFaces {
470            return Err(MapError::Unsupported(
471                "terrain_scope is qualified for plain sloping only".into(),
472            ));
473        }
474    }
475    let started = Instant::now();
476    let [col, row, width, height] = req.window;
477    let count = width
478        .checked_mul(height)
479        .filter(|&n| n > 0 && n <= model::MAX_PROFILE_PIXELS as usize)
480        .ok_or_else(|| MapError::Budget("sloping window empty or too large".into()))?;
481    if req.schema != SCHEMA
482        || !(1..=64).contains(&req.threads)
483        || !(1..=64).contains(&req.max_inflight_profiles)
484        || req.epsg == 0
485        || req.mesh_files.is_empty()
486        || req.owner_cells.is_empty()
487        || req.owner_cells.len() >= u32::MAX as usize
488        || req.profiles.is_empty()
489    {
490        return Err(MapError::Invalid("invalid sloping request".into()));
491    }
492    if let Some(files) = &req.velocity_normals_files
493        && files.len() != req.mesh_files.len()
494    {
495        return Err(MapError::Invalid(
496            "velocity_normals_files needs one file per mesh".into(),
497        ));
498    }
499    if req.profiles.iter().any(|p| {
500        p.face_velocity_files.is_some() != req.velocity_normals_files.is_some()
501            || p.face_velocity_files
502                .as_ref()
503                .is_some_and(|f| f.len() != req.mesh_files.len())
504    }) {
505        return Err(MapError::Invalid(
506            "velocity_normals_files and every profile's face_velocity_files (one per mesh) go together".into(),
507        ));
508    }
509    let mut ids = HashSet::new();
510    for profile in &req.profiles {
511        model::validate_id(&profile.id)?;
512        if !ids.insert(model::id_key(&profile.id))
513            || profile.values_files.len() != req.mesh_files.len()
514        {
515            return Err(MapError::Invalid(
516                "duplicate sloping profile id or mesh value count mismatch".into(),
517            ));
518        }
519    }
520    let geometry_bytes = req
521        .mesh_files
522        .iter()
523        .chain(std::iter::once(&req.geometry_file))
524        .chain(req.aoi_file.iter())
525        .try_fold(0_u64, |sum, p| {
526            sum.checked_add(fs::metadata(model::resolve(base, p))?.len())
527                .ok_or_else(|| MapError::Budget("sloping geometry size overflow".into()))
528        })?;
529    let binary_bytes = req.mesh_files.iter().try_fold(0_u64, |sum, name| {
530        sum.checked_add(crate::sloping_binary::input_bytes(&model::resolve(
531            base, name,
532        ))?)
533        .ok_or_else(|| MapError::Budget("sloping binary geometry size overflow".into()))
534    })?;
535    let geometry_bytes = geometry_bytes
536        .checked_add(binary_bytes)
537        .ok_or_else(|| MapError::Budget("sloping geometry size overflow".into()))?;
538    let cache = 64_u64 << 20;
539    // JSON parsing, topology, stencils, one hydraulic profile and IO copies are
540    // conservatively admitted, rather than claiming an OS-enforced RSS limit.
541    let persistent_bytes = if session.is_some() {
542        req.window_cache_bytes
543    } else {
544        0
545    };
546    let hydraulic_bytes = if session.is_some() {
547        req.hydraulic_cache_bytes
548    } else {
549        0
550    };
551    if let Some(s) = session.as_mut() {
552        s.windows.resize(persistent_bytes);
553        s.hydraulics.lock().unwrap().resize(hydraulic_bytes);
554    }
555    let aoi_scratch = crate::aoi::scratch_bytes(
556        req.aoi_file
557            .as_ref()
558            .map(|p| model::resolve(base, p))
559            .as_deref(),
560        req.threads,
561    )?;
562    let planned = (count as u64)
563        .checked_mul(
564            if req.options.render_mode == sloping::RenderMode::WithFaces {
565                56
566            } else {
567                40
568            },
569        )
570        .and_then(|v| {
571            geometry_bytes
572                .checked_mul(
573                    if req.options.render_mode == sloping::RenderMode::WithFaces {
574                        20
575                    } else {
576                        16
577                    },
578                )
579                .and_then(|g| v.checked_add(g))
580        })
581        .and_then(|v| v.checked_add(cache))
582        .and_then(|v| v.checked_add(persistent_bytes))
583        .and_then(|v| v.checked_add(hydraulic_bytes))
584        .and_then(|v| v.checked_add(aoi_scratch))
585        // Terrain scope: one validity byte per expanded-window pixel.
586        .and_then(|v| {
587            v.checked_add(if req.terrain_scope.is_some() {
588                (width as u64 + 2) * (height as u64 + 2)
589            } else {
590                0
591            })
592        })
593        .ok_or_else(|| MapError::Budget("sloping array plan overflow".into()))?;
594    if planned > req.max_working_bytes {
595        return Err(MapError::Budget("sloping array plan exceeds budget".into()));
596    }
597    let key = if session.is_some() {
598        geometry_key(req, base)?
599    } else {
600        String::new()
601    };
602    // Release an incompatible previous model before opening/loading this one.
603    if let Some(s) = session.as_mut()
604        && s.geometry.as_ref().is_some_and(|(old, _)| old != &key)
605    {
606        s.geometry = None;
607        s.hydraulics.lock().unwrap().clear();
608        s.windows.clear();
609    }
610    let terrain_path = model::resolve(base, &req.terrain);
611    let terrain_key = if session.is_some() {
612        format!(
613            "{}:{}:{}",
614            identity(&terrain_path)?,
615            req.epsg,
616            serde_json::to_string(&req.crs_identity)
617                .map_err(|e| MapError::Invalid(e.to_string()))?
618        )
619    } else {
620        String::new()
621    };
622    let terrain_reused = session.as_ref().is_some_and(|s| {
623        s.terrain
624            .as_ref()
625            .is_some_and(|(key, _)| key == &terrain_key)
626    });
627    if !terrain_reused && let Some(s) = session.as_mut() {
628        s.terrain = None;
629        s.windows.clear();
630    }
631    let terrain = if terrain_reused {
632        session
633            .as_ref()
634            .unwrap()
635            .terrain
636            .as_ref()
637            .unwrap()
638            .1
639            .clone()
640    } else {
641        Arc::new(NativeTerrain::open(
642            model::resolve(base, &req.terrain)
643                .to_str()
644                .ok_or_else(|| MapError::Invalid("terrain path is not UTF8".into()))?,
645            req.epsg,
646            &req.crs_identity,
647            QueryConvention::PixelCenter,
648            NodataPolicy::Dry,
649            cache,
650        )?)
651    };
652    if let Some(s) = session.as_mut() {
653        s.terrain = Some((terrain_key.clone(), terrain.clone()));
654    }
655    if terrain.nodata.is_none()
656        || terrain.gt != req.terrain_transform
657        || (terrain.gt[1] + terrain.gt[5]).abs() > terrain.gt[1].abs() * 1e-12
658        || col.checked_add(width).is_none_or(|v| v > terrain.width)
659        || row.checked_add(height).is_none_or(|v| v > terrain.height)
660        || terrain.max_block_bytes() > cache
661    {
662        return Err(MapError::Unsupported(
663            "sloping terrain identity, square grid, block or window mismatch".into(),
664        ));
665    }
666    let admitted = planned
667        .checked_add(
668            if req.options.render_mode == sloping::RenderMode::WithFaces {
669                (terrain.block.0 as u64)
670                    .checked_mul(terrain.block.1 as u64)
671                    .and_then(|v| v.checked_mul(req.threads as u64))
672                    .and_then(|v| v.checked_mul(64))
673                    .ok_or_else(|| {
674                        MapError::Budget("native donation replay tile plan overflow".into())
675                    })?
676            } else {
677                0
678            },
679        )
680        .ok_or_else(|| MapError::Budget("native donation replay plan overflow".into()))?
681        .checked_add(
682            terrain
683                .decoder_transient_bytes_per_load()
684                .checked_mul(req.threads as u64)
685                .ok_or_else(|| MapError::Budget("sloping parallel decoder plan overflow".into()))?,
686        )
687        .ok_or_else(|| MapError::Budget("sloping decoder plan overflow".into()))?;
688    // Retain the serial plan. Each additional slot conservatively reserves
689    // output/work arrays and the existing full geometry/hydraulic allowance.
690    // This is planned admission, not an OS-enforced RSS limit.
691    let additional = (count as u64)
692        .checked_mul(
693            if req.options.render_mode == sloping::RenderMode::WithFaces {
694                32
695            } else {
696                16
697            },
698        )
699        .and_then(|v| {
700            geometry_bytes
701                .checked_mul(
702                    if req.options.render_mode == sloping::RenderMode::WithFaces {
703                        20
704                    } else {
705                        16
706                    },
707                )
708                .and_then(|g| v.checked_add(g))
709        })
710        .ok_or_else(|| MapError::Budget("sloping concurrent profile plan overflow".into()))?;
711    let (effective, admitted) = admit_profiles(
712        req.max_inflight_profiles,
713        req.threads,
714        req.profiles.len(),
715        req.max_working_bytes,
716        admitted,
717        additional,
718    )?;
719    fs::create_dir(out)?;
720    let mut report = json!({"schema":"ras-commander-sloping-report/1","status":"in_progress", "method":selected_method(req.options),
721        "threads":req.threads,"window":req.window,"planned_bytes":admitted,
722        "requested_inflight_profiles":req.max_inflight_profiles,"effective_inflight_profiles":effective,
723        "additional_profile_planned_bytes":additional,"profiles":[]});
724    sink::write_json(&out.join("report.json"), &report)?;
725    let execution = std::panic::catch_unwind(std::panic::AssertUnwindSafe(
726        || -> Result<(), MapError> {
727            if let Some(s) = session.as_mut()
728                && s.pool
729                    .as_ref()
730                    .is_some_and(|(threads, _)| *threads != req.threads)
731            {
732                s.pool = None;
733            }
734            let pool = if let Some((threads, pool)) = session
735                .as_ref()
736                .and_then(|s| s.pool.as_ref())
737                .filter(|(threads, _)| *threads == req.threads)
738            {
739                let _ = threads;
740                pool.clone()
741            } else {
742                Arc::new(
743                    rayon::ThreadPoolBuilder::new()
744                        .num_threads(req.threads)
745                        .build()
746                        .map_err(|e| MapError::Internal(e.to_string()))?,
747                )
748            };
749            if let Some(s) = session.as_mut() {
750                s.pool = Some((req.threads, pool.clone()));
751            }
752            let tick = Instant::now();
753            let geometry_reused = session
754                .as_ref()
755                .is_some_and(|s| s.geometry.as_ref().is_some_and(|(old, _)| old == &key));
756            if !geometry_reused && let Some(s) = session.as_mut() {
757                s.geometry = None;
758                s.hydraulics.lock().unwrap().clear();
759                s.windows.clear();
760            }
761            let geometry = if geometry_reused {
762                session
763                    .as_ref()
764                    .unwrap()
765                    .geometry
766                    .as_ref()
767                    .unwrap()
768                    .1
769                    .clone()
770            } else {
771                let mut meshes = Vec::new();
772                let mut stencils = Vec::new();
773                for filename in &req.mesh_files {
774                    let source = crate::sloping_binary::load(&model::resolve(base, filename))?;
775                    let (mesh, stencil) = source.prepare()?;
776                    meshes.push(mesh);
777                    stencils.push(stencil);
778                }
779                Arc::new((meshes, stencils))
780            };
781            if let Some(s) = session.as_mut() {
782                s.geometry = Some((key.clone(), geometry.clone()));
783            }
784            let (meshes, stencils) = &*geometry;
785            let velocity_normals = req
786                .velocity_normals_files
787                .as_ref()
788                .map(|files| -> Result<Vec<_>, MapError> {
789                    files
790                        .iter()
791                        .zip(meshes)
792                        .map(|(file, mesh)| {
793                            crate::velocity::read_normals(
794                                &model::resolve(base, file),
795                                mesh.faces.len(),
796                            )
797                        })
798                        .collect()
799                })
800                .transpose()?;
801            if velocity_normals.is_some() {
802                report["velocity"] = json!({"method": crate::velocity::method_record()});
803            }
804            report["session"] = json!({"terrain_reused":terrain_reused,"geometry_reused":geometry_reused,"retained_window_budget_bytes":persistent_bytes});
805            for &[mesh, cell] in &req.owner_cells {
806                if mesh >= meshes.len()
807                    || cell >= meshes[mesh].cells.len()
808                    || meshes[mesh].cells[cell].virtual_cell
809                {
810                    return Err(MapError::Invalid(
811                        "sloping ownership requires physical mesh cells".into(),
812                    ));
813                }
814            }
815            report["geometry_prepare_seconds"] = tick.elapsed().as_secs_f64().into();
816            report["ownership_prepare_seconds"] = 0.0.into();
817            report["terrain_prepare_seconds"] = 0.0.into();
818            let transform = [
819                terrain.gt[0] + col as f64 * terrain.gt[1],
820                terrain.gt[1],
821                0.,
822                terrain.gt[3] + row as f64 * terrain.gt[5],
823                0.,
824                terrain.gt[5],
825            ];
826            let grid = terrio::ops::rasterize::Grid {
827                geo_transform: transform,
828                cols: width,
829                rows: height,
830                crs: Default::default(),
831            };
832            let tick = Instant::now();
833            let window_key = if session.is_some() {
834                (
835                    terrain_key.clone(),
836                    key.clone(),
837                    req.window,
838                    req.owner_cells.clone(),
839                    fs::read(model::resolve(base, &req.geometry_file))?,
840                    req.aoi_file
841                        .as_ref()
842                        .map(|p| fs::read(model::resolve(base, p)))
843                        .transpose()?,
844                )
845            } else {
846                (
847                    String::new(),
848                    String::new(),
849                    req.window,
850                    Vec::new(),
851                    Vec::new(),
852                    None,
853                )
854            };
855            let cached_window = session.as_mut().and_then(|s| s.windows.get(&window_key));
856            let window_reused = cached_window.is_some();
857            let mut data = if let Some(data) = cached_window {
858                data
859            } else {
860                let cells = crate::ownership::load(
861                    &model::resolve(base, &req.geometry_file),
862                    req.owner_cells.len(),
863                    false,
864                )?;
865                let mut owner = pool.install(|| crate::ownership::burn(cells, grid.clone()))?;
866                if let Some(filename) = &req.aoi_file {
867                    let mask = pool.install(|| {
868                        crate::aoi::mask(&model::resolve(base, filename), terrain.gt, req.window)
869                    })?;
870                    pool.install(|| {
871                        owner
872                            .par_iter_mut()
873                            .zip(mask.par_iter())
874                            .for_each(|(label, &m)| {
875                                if m == 0 {
876                                    *label = 0;
877                                }
878                            })
879                    });
880                }
881                report["ownership_prepare_seconds"] = tick.elapsed().as_secs_f64().into();
882                let tick = Instant::now();
883                let mut ground = vec![sloping::NODATA; count];
884                let (bw, bh) = terrain.block;
885                // Row stripes own disjoint output slices. Decodes run inside the same
886                // bounded pool, and the shared single-flight cache handles boundary
887                // blocks seen by adjacent stripes. Skip blocks outside ownership/AOI.
888                let stripe_pixels = width
889                    .checked_mul(bh)
890                    .ok_or_else(|| MapError::Budget("sloping terrain stripe overflow".into()))?;
891                pool.install(|| {
892                    ground
893                        .par_chunks_mut(stripe_pixels)
894                        .enumerate()
895                        .try_for_each(|(stripe, destination)| -> Result<(), MapError> {
896                            let first = row + stripe * bh;
897                            let last = first + destination.len() / width;
898                            for by in first / bh..=(last - 1) / bh {
899                                for bx in col / bw..=(col + width - 1) / bw {
900                                    let win = terrain.block_window(bx, by);
901                                    let r0 = first.max(win.y);
902                                    let r1 = last.min(win.y + win.height);
903                                    let c0 = col.max(win.x);
904                                    let c1 = (col + width).min(win.x + win.width);
905                                    if !(r0..r1).any(|r| {
906                                        owner[(r - row) * width + c0 - col
907                                            ..(r - row) * width + c1 - col]
908                                            .iter()
909                                            .any(|&label| label != 0)
910                                    }) {
911                                        continue;
912                                    }
913                                    let block = terrain.get_block(bx, by)?;
914                                    for r in r0..r1 {
915                                        for c in c0..c1 {
916                                            if owner[(r - row) * width + c - col] == 0 {
917                                                continue;
918                                            }
919                                            if let Some(value) = terrain.classify(
920                                                block[(r - win.y) * win.width + c - win.x],
921                                            ) {
922                                                if !value.is_finite() {
923                                                    return Err(MapError::Terrain(
924                                                        "nonfinite sloping terrain".into(),
925                                                    ));
926                                                }
927                                                destination[(r - first) * width + c - col] =
928                                                    value as f32;
929                                            }
930                                        }
931                                    }
932                                }
933                            }
934                            Ok(())
935                        })
936                })?;
937                report["terrain_prepare_seconds"] = tick.elapsed().as_secs_f64().into();
938                Arc::new(WindowData {
939                    owner,
940                    ground,
941                    weights: None,
942                })
943            };
944            let window_base_bytes = count as u64 * 8
945                + (2 * std::mem::size_of::<(WindowKey, Arc<WindowData>, u64)>()
946                    + std::mem::size_of::<WindowData>()
947                    + 2 * std::mem::size_of::<usize>()) as u64
948                + window_key.0.len() as u64
949                + window_key.1.len() as u64
950                + (window_key.3.len() * std::mem::size_of::<[usize; 2]>()) as u64
951                + window_key.4.len() as u64
952                + window_key.5.as_ref().map_or(0, |v| v.len()) as u64;
953            let weight_budget = persistent_bytes
954                .saturating_sub(window_base_bytes)
955                .min(32 << 20);
956            let (bw, bh) = terrain.block;
957            if session.is_some()
958                && !window_reused
959                && weight_budget > 0
960                && req.options.render_mode != sloping::RenderMode::WithFaces
961            {
962                let window = sloping::RasterWindow {
963                    source_transform: terrain.gt,
964                    source_tile: [bw, bh],
965                    width,
966                    height,
967                    transform,
968                    owner: &data.owner,
969                    terrain: &data.ground,
970                    owner_cells: &req.owner_cells,
971                };
972                let weights = sloping::PixelWeights::prepare(
973                    stencils,
974                    &window,
975                    weight_budget as usize,
976                    &pool,
977                )?;
978                Arc::get_mut(&mut data).unwrap().weights = Some(weights);
979            }
980            let weight_bytes = data.weights.as_ref().map_or(0, |w| w.bytes());
981            report["session"]["window_reused"] = window_reused.into();
982            report["session"]["pixel_weight_cache_bytes"] = weight_bytes.into();
983            let window_bytes = window_base_bytes + weight_bytes as u64;
984            if let Some(s) = session.as_mut() {
985                s.windows.insert(window_key, data.clone(), window_bytes);
986                report["session"]["retained_window_bytes"] = s.windows.bytes().into();
987                report["session"]["retained_window_count"] = s.windows.len().into();
988            }
989            let window = sloping::RasterWindow {
990                source_transform: terrain.gt,
991                source_tile: [bw, bh],
992                width,
993                height,
994                transform,
995                owner: &data.owner,
996                terrain: &data.ground,
997                owner_cells: &req.owner_cells,
998            };
999            let foreign = match &req.terrain_scope {
1000                Some(scope) => {
1001                    let scoped = crate::terrain_scope::check(
1002                        scope,
1003                        base,
1004                        &terrain,
1005                        req.window,
1006                        &data.owner,
1007                    )?;
1008                    report["terrain_scope"] = scoped.record;
1009                    scoped.foreign
1010                }
1011                None => Vec::new(),
1012            };
1013            for batch in req.profiles.chunks(effective) {
1014                let results = map_batch(batch, &pool, |profile| {
1015                    render(
1016                        profile,
1017                        base,
1018                        meshes,
1019                        stencils,
1020                        &window,
1021                        &pool,
1022                        req.options,
1023                        data.weights.as_ref(),
1024                        velocity_normals.as_deref(),
1025                        &foreign,
1026                        session
1027                            .as_ref()
1028                            .filter(|_| hydraulic_bytes > 0)
1029                            .map(|s| &s.hydraulics),
1030                    )
1031                })?;
1032                // Indexed collection and ordered emission preserve the pipe/report
1033                // contract even when later profiles finish first.
1034                for (profile, rendered) in batch.iter().zip(results) {
1035                    let RenderedProfile {
1036                        depth,
1037                        wse,
1038                        velocity,
1039                        hydraulic_seconds,
1040                        render_seconds,
1041                        hydraulic_reused,
1042                        finite_pixels,
1043                        flooded_pixels,
1044                    } = rendered;
1045                    let dir = out.join(&profile.id);
1046                    fs::create_dir(&dir)?;
1047                    if req.stream_outputs {
1048                        let mut stream =
1049                            std::io::BufWriter::with_capacity(65536, std::io::stdout().lock());
1050                        let mut header =
1051                            json!({"id":profile.id,"status":"complete","pixels":count});
1052                        if velocity.is_some() {
1053                            header["velocity"] = true.into();
1054                        }
1055                        if req.geotiff_outputs {
1056                            let tick = Instant::now();
1057                            header["products"] = "native-geotiff/1".into();
1058                            header["depth_finite_pixels"] = finite_pixels.into();
1059                            header["flooded_pixels"] = flooded_pixels.into();
1060                            crate::geotiff::write(
1061                                &dir.join("depth.tif"),
1062                                depth,
1063                                width,
1064                                height,
1065                                transform,
1066                                req.epsg,
1067                            )?;
1068                            crate::geotiff::write(
1069                                &dir.join("wse.tif"),
1070                                wse,
1071                                width,
1072                                height,
1073                                transform,
1074                                req.epsg,
1075                            )?;
1076                            header["write_seconds"] = tick.elapsed().as_secs_f64().into();
1077                            writeln!(stream, "{header}")?;
1078                            stream.flush()?;
1079                        } else {
1080                            writeln!(stream, "{header}")?;
1081                            crate::horizontal::write_f32(&mut stream, &depth)?;
1082                            crate::horizontal::write_f32(&mut stream, &wse)?;
1083                            if let Some(velocity) = &velocity {
1084                                crate::horizontal::write_f32(&mut stream, velocity)?;
1085                            }
1086                            stream.flush()?;
1087                        }
1088                    } else {
1089                        crate::horizontal::write_values(&dir.join("depth.bin"), &depth)?;
1090                        crate::horizontal::write_values(&dir.join("wse.bin"), &wse)?;
1091                        if let Some(velocity) = &velocity {
1092                            crate::horizontal::write_values(&dir.join("velocity.bin"), velocity)?;
1093                        }
1094                    }
1095                    let unit = json!({"id":profile.id,"status":"complete","hydraulic_seconds":hydraulic_seconds,"render_seconds":render_seconds,"hydraulic_reused":hydraulic_reused});
1096                    sink::write_json(&dir.join("unit.json"), &unit)?;
1097                    report["profiles"].as_array_mut().unwrap().push(unit);
1098                    sink::write_json(&out.join("report.json"), &report)?;
1099                }
1100            }
1101            if let Some(s) = session.as_ref() {
1102                let cache = s.hydraulics.lock().unwrap();
1103                report["session"]["retained_hydraulic_bytes"] = cache.bytes().into();
1104                report["session"]["retained_hydraulic_count"] = cache.len().into();
1105            }
1106            Ok(())
1107        },
1108    ));
1109    match execution {
1110        Ok(Ok(())) => report["status"] = "complete".into(),
1111        other => {
1112            report["status"] = "aborted".into();
1113            report["error"] = match other {
1114                Ok(Err(e)) => e.to_string(),
1115                _ => "sloping internal panic".into(),
1116            }
1117            .into();
1118            for profile in &req.profiles {
1119                let dir = out.join(&profile.id);
1120                for filename in [
1121                    "depth.bin",
1122                    "wse.bin",
1123                    "depth.bin.partial",
1124                    "wse.bin.partial",
1125                    "depth.tif",
1126                    "wse.tif",
1127                ] {
1128                    let _ = fs::remove_file(dir.join(filename));
1129                }
1130                if dir.is_dir() {
1131                    sink::write_json(
1132                        &dir.join("unit.json"),
1133                        &json!({"id":profile.id,"status":"failed","error":"session aborted"}),
1134                    )?;
1135                }
1136            }
1137            report["profiles"] = req
1138                .profiles
1139                .iter()
1140                .map(|p| json!({"id":p.id,"status":"failed","error":"session aborted"}))
1141                .collect::<Vec<_>>()
1142                .into();
1143        }
1144    }
1145    report["total_seconds"] = started.elapsed().as_secs_f64().into();
1146    report["terrain_cache"] = json!(terrain.stats());
1147    sink::write_json(&out.join("report.json"), &report)?;
1148    Ok(report)
1149}
1150
1151#[cfg(test)]
1152mod tests {
1153    use super::*;
1154    use std::sync::{
1155        Barrier,
1156        atomic::{AtomicUsize, Ordering},
1157    };
1158
1159    #[test]
1160    fn admission_caps_threads_profiles_memory_and_overflow() {
1161        assert_eq!(admit_profiles(8, 4, 6, 1000, 100, 50).unwrap(), (4, 250));
1162        assert_eq!(admit_profiles(8, 4, 2, 1000, 100, 50).unwrap(), (2, 150));
1163        assert_eq!(admit_profiles(8, 4, 6, 149, 100, 50).unwrap(), (1, 100));
1164        assert!(admit_profiles(8, 4, 6, 99, 100, 50).is_err());
1165        assert_eq!(
1166            admit_profiles(64, 64, 64, u64::MAX, u64::MAX - 1, u64::MAX).unwrap(),
1167            (1, u64::MAX - 1)
1168        );
1169        assert!(admit_profiles(0, 1, 1, 100, 50, 10).is_err());
1170    }
1171
1172    #[test]
1173    fn batches_really_overlap_without_exceeding_bound_and_keep_order() {
1174        let pool = rayon::ThreadPoolBuilder::new()
1175            .num_threads(2)
1176            .build()
1177            .unwrap();
1178        let active = AtomicUsize::new(0);
1179        let peak = AtomicUsize::new(0);
1180        let barrier = Barrier::new(2);
1181        let mut output = Vec::new();
1182        for batch in [0, 1, 2, 3].chunks(2) {
1183            output.extend(
1184                map_batch(batch, &pool, |&i| {
1185                    let n = active.fetch_add(1, Ordering::SeqCst) + 1;
1186                    peak.fetch_max(n, Ordering::SeqCst);
1187                    barrier.wait();
1188                    active.fetch_sub(1, Ordering::SeqCst);
1189                    Ok(i)
1190                })
1191                .unwrap(),
1192            );
1193        }
1194        assert_eq!(peak.load(Ordering::SeqCst), 2);
1195        assert_eq!(output, vec![0, 1, 2, 3]);
1196        assert!(
1197            map_batch(&[0, 1], &pool, |_| -> Result<(), MapError> {
1198                Err(MapError::Invalid("failed profile".into()))
1199            })
1200            .is_err()
1201        );
1202    }
1203}