spatialrust_interchange/tiles3d/
copc.rs

1//! Bounded COPC → 3D Tiles 1.1 tileset export.
2//!
3//! Opens a COPC file once, walks its octree hierarchy in deterministic order,
4//! and writes one `pnts` tile per COPC node plus a `tileset.json` that mirrors
5//! the octree parent/child structure. The cloud is never materialized as a
6//! whole: each node is decoded, re-centered to its own `RTC_CENTER`, encoded,
7//! and dropped before the next node is processed.
8
9use std::path::Path;
10
11use spatialrust_core::HasPositions3;
12use spatialrust_io::CopcNode;
13
14use crate::tiles3d::builder::{BuiltTile, BuiltTileset, TilesetWriteReceipt};
15use crate::tiles3d::pnts::{encode_pnts, PntsFeatureTable};
16use crate::tiles3d::tileset::{
17    serialize_tileset_json, BoundingVolume, Refinement, Tile, TileContent, Tileset,
18};
19use crate::{InterchangeError, InterchangeResult};
20
21/// Options controlling COPC → 3D Tiles export.
22#[derive(Clone, Debug)]
23pub struct CopcTilesetOptions {
24    /// Maximum octree depth to export; `None` exports the whole hierarchy.
25    pub max_level: Option<i32>,
26    /// Refinement policy recorded on the root tile.
27    pub refine: Refinement,
28    /// Root geometric error; defaults to the root bounds diagonal when `None`.
29    pub root_geometric_error: Option<f64>,
30    /// Multiplier applied to the parent geometric error for each child level.
31    pub geometric_error_scale: f64,
32}
33
34impl Default for CopcTilesetOptions {
35    fn default() -> Self {
36        Self {
37            max_level: None,
38            refine: Refinement::Replace,
39            root_geometric_error: None,
40            geometric_error_scale: 0.5,
41        }
42    }
43}
44
45impl CopcTilesetOptions {
46    fn validate(&self) -> InterchangeResult<()> {
47        if let Some(level) = self.max_level {
48            if level < 0 {
49                return Err(InterchangeError::InvalidConfiguration(
50                    "max_level must be non-negative".into(),
51                ));
52            }
53        }
54        if let Some(error) = self.root_geometric_error {
55            if !error.is_finite() || error < 0.0 {
56                return Err(InterchangeError::InvalidConfiguration(
57                    "root_geometric_error must be finite and non-negative".into(),
58                ));
59            }
60        }
61        if !self.geometric_error_scale.is_finite()
62            || !(0.0..=1.0).contains(&self.geometric_error_scale)
63        {
64            return Err(InterchangeError::InvalidConfiguration(
65                "geometric_error_scale must be within [0, 1]".into(),
66            ));
67        }
68        Ok(())
69    }
70}
71
72/// Exports a COPC file into a 3D Tiles 1.1 tileset without materializing the
73/// whole cloud. Returns the byte/point/tile receipt.
74pub fn export_copc_tileset(
75    copc_path: impl AsRef<Path>,
76    out_dir: impl AsRef<Path>,
77    options: &CopcTilesetOptions,
78) -> InterchangeResult<TilesetWriteReceipt> {
79    options.validate()?;
80    let mut reader = spatialrust_io::CopcNodeReader::open(copc_path.as_ref())
81        .map_err(|error| InterchangeError::InvalidConfiguration(format!("COPC open: {error}")))?;
82
83    let nodes: Vec<CopcNode> = reader
84        .nodes()
85        .iter()
86        .filter(|node| options.max_level.map_or(true, |level| node.level <= level))
87        .copied()
88        .collect();
89    if nodes.is_empty() {
90        return Err(InterchangeError::InvalidConfiguration(
91            "COPC file contains no nodes at or above max_level".into(),
92        ));
93    }
94
95    let root_error = match options.root_geometric_error {
96        Some(error) => error,
97        None => bounds_diagonal(nodes[0].bounds.min, nodes[0].bounds.max),
98    };
99
100    let index_of = |level: i32, x: i32, y: i32, z: i32| -> Option<usize> {
101        nodes
102            .iter()
103            .position(|node| node.level == level && node.x == x && node.y == y && node.z == z)
104    };
105
106    let mut tiles: Vec<BuiltTile> = Vec::with_capacity(nodes.len());
107    let mut children_by_parent: Vec<Vec<usize>> = vec![Vec::new(); nodes.len()];
108    for (index, node) in nodes.iter().enumerate() {
109        if node.level == 0 {
110            continue;
111        }
112        if let Some(parent) = index_of(node.level - 1, node.x >> 1, node.y >> 1, node.z >> 1) {
113            children_by_parent[parent].push(index);
114        }
115    }
116    for children in &mut children_by_parent {
117        children.sort_unstable();
118    }
119
120    // Recursive materialization in deterministic parent-before-child order.
121    let mut root_tile_index = None;
122    for (index, node) in nodes.iter().enumerate() {
123        if node.level == 0 {
124            root_tile_index = Some(index);
125            break;
126        }
127    }
128    let root_index = root_tile_index.ok_or_else(|| {
129        InterchangeError::InvalidConfiguration("COPC file has no level-0 root node".into())
130    })?;
131
132    fn materialize(
133        reader: &mut spatialrust_io::CopcNodeReader,
134        nodes: &[CopcNode],
135        children: &[Vec<usize>],
136        tiles: &mut Vec<BuiltTile>,
137        index: usize,
138        parent_error: f64,
139        scale: f64,
140    ) -> InterchangeResult<Tile> {
141        let node = &nodes[index];
142        let bounds = node.bounds;
143        let center = [
144            (bounds.min[0] + bounds.max[0]) * 0.5,
145            (bounds.min[1] + bounds.max[1]) * 0.5,
146            (bounds.min[2] + bounds.max[2]) * 0.5,
147        ];
148        let geometric_error = if node.level == 0 { parent_error } else { parent_error * scale };
149
150        let cloud = reader.read_node(index).map_err(|error| {
151            InterchangeError::InvalidConfiguration(format!("COPC node: {error}"))
152        })?;
153        let (x, y, z) = cloud.positions3().map_err(|error| {
154            InterchangeError::InvalidConfiguration(format!("COPC schema: {error}"))
155        })?;
156        let rgb = extract_rgb(&cloud);
157        let mut local_positions = Vec::with_capacity(cloud.len() * 3);
158        for point_index in 0..cloud.len() {
159            local_positions.push(x[point_index] - center[0] as f32);
160            local_positions.push(y[point_index] - center[1] as f32);
161            local_positions.push(z[point_index] - center[2] as f32);
162        }
163        let table = PntsFeatureTable { positions: local_positions, rgb, rtc_center: Some(center) };
164        let pnts = encode_pnts(&table)?;
165        let uri = format!("{}.pnts", tiles.len());
166        let point_count = table.point_count();
167        tiles.push(BuiltTile { uri: uri.clone(), pnts, point_count });
168
169        let mut child_tiles = Vec::new();
170        for &child in &children[index] {
171            child_tiles.push(materialize(
172                reader,
173                nodes,
174                children,
175                tiles,
176                child,
177                geometric_error,
178                scale,
179            )?);
180        }
181
182        Ok(Tile {
183            bounding_volume: BoundingVolume::box_from_bounds(bounds.min, bounds.max)?,
184            geometric_error,
185            refine: None,
186            content: Some(TileContent { uri }),
187            children: child_tiles,
188        })
189    }
190
191    let root_tile = materialize(
192        &mut reader,
193        &nodes,
194        &children_by_parent,
195        &mut tiles,
196        root_index,
197        root_error,
198        options.geometric_error_scale,
199    )?;
200    let mut root_tile = root_tile;
201    root_tile.refine = Some(options.refine);
202
203    let built =
204        BuiltTileset { tileset: Tileset { geometric_error: root_error, root: root_tile }, tiles };
205    write_built(out_dir.as_ref(), &built)
206}
207
208/// Writes a built tileset to `dir` and returns the receipt.
209fn write_built(dir: &Path, built: &BuiltTileset) -> InterchangeResult<TilesetWriteReceipt> {
210    std::fs::create_dir_all(dir).map_err(io_error)?;
211    let document = serialize_tileset_json(&built.tileset)?;
212    std::fs::write(dir.join("tileset.json"), document.as_bytes()).map_err(io_error)?;
213
214    let mut tile_count = 0u64;
215    let mut point_count = 0u64;
216    let mut pnts_bytes = 0u64;
217    for tile in &built.tiles {
218        std::fs::write(dir.join(&tile.uri), &tile.pnts).map_err(io_error)?;
219        tile_count += 1;
220        point_count = point_count
221            .checked_add(tile.point_count as u64)
222            .ok_or_else(|| InterchangeError::InvalidConfiguration("point count overflow".into()))?;
223        pnts_bytes = pnts_bytes.checked_add(tile.pnts.len() as u64).ok_or_else(|| {
224            InterchangeError::InvalidConfiguration("pnts byte count overflow".into())
225        })?;
226    }
227    Ok(TilesetWriteReceipt {
228        tileset_json_bytes: document.len() as u64,
229        tile_count,
230        point_count,
231        pnts_bytes,
232    })
233}
234
235fn bounds_diagonal(min: [f64; 3], max: [f64; 3]) -> f64 {
236    let dx = max[0] - min[0];
237    let dy = max[1] - min[1];
238    let dz = max[2] - min[2];
239    (dx * dx + dy * dy + dz * dz).sqrt()
240}
241
242/// Extracts interleaved 8-bit RGB from a cloud's color fields when present.
243///
244/// LAS/COPC color is stored as `u16` (0–65535); 3D Tiles `pnts` RGB uses
245/// `u8` (0–255), so each channel is shifted right by eight bits. Returns
246/// `None` when the cloud has no color fields.
247fn extract_rgb(cloud: &spatialrust_core::PointCloud) -> Option<Vec<u8>> {
248    use spatialrust_core::{FieldSemantic, PointBuffer};
249
250    let field = |semantic: FieldSemantic| {
251        let name = match semantic {
252            FieldSemantic::ColorR => "red",
253            FieldSemantic::ColorG => "green",
254            FieldSemantic::ColorB => "blue",
255            _ => return None,
256        };
257        let buffer = cloud.field(name).ok()?;
258        match buffer {
259            PointBuffer::U16(values) => Some(values.as_slice()),
260            _ => None,
261        }
262    };
263
264    let (r, g, b) = match (
265        field(FieldSemantic::ColorR),
266        field(FieldSemantic::ColorG),
267        field(FieldSemantic::ColorB),
268    ) {
269        (Some(r), Some(g), Some(b)) => (r, g, b),
270        _ => return None,
271    };
272    if r.len() != cloud.len() || g.len() != cloud.len() || b.len() != cloud.len() {
273        return None;
274    }
275    let mut out = Vec::with_capacity(cloud.len() * 3);
276    for index in 0..cloud.len() {
277        out.push((r[index] >> 8) as u8);
278        out.push((g[index] >> 8) as u8);
279        out.push((b[index] >> 8) as u8);
280    }
281    Some(out)
282}
283
284fn io_error(error: std::io::Error) -> InterchangeError {
285    InterchangeError::InvalidConfiguration(format!("tileset IO failure: {error}"))
286}
287
288#[cfg(test)]
289mod tests {
290    use super::{export_copc_tileset, CopcTilesetOptions};
291    use crate::tiles3d::pnts::decode_pnts;
292    use crate::tiles3d::tileset::parse_tileset_json;
293    use spatialrust_core::PointCloudBuilder;
294    use spatialrust_io::{write_copc_file, write_copc_file_with_params, CopcWriterParams};
295
296    fn dense_grid_cloud(count: usize) -> spatialrust_core::PointCloud {
297        let mut builder = PointCloudBuilder::xyz();
298        for index in 0..count {
299            let x = (index % 31) as f32 - 15.0;
300            let y = ((index / 31) % 29) as f32 - 14.0;
301            let z = ((index / (31 * 29)) % 23) as f32 - 11.0;
302            builder.push_point([x, y, z]).unwrap();
303        }
304        builder.build().unwrap()
305    }
306
307    #[test]
308    fn exports_copc_without_full_materialization() {
309        let cloud = dense_grid_cloud(7_000);
310        let copc_path = std::env::temp_dir()
311            .join(format!("spatialrust_tiles3d_copc_{}.copc.laz", std::process::id()));
312        write_copc_file_with_params(
313            &copc_path,
314            &cloud,
315            &CopcWriterParams { max_points_per_node: 96, max_depth: 8 },
316        )
317        .unwrap();
318
319        let out_dir = std::env::temp_dir()
320            .join(format!("spatialrust_tiles3d_copc_out_{}", std::process::id()));
321        let receipt =
322            export_copc_tileset(&copc_path, &out_dir, &CopcTilesetOptions::default()).unwrap();
323        assert_eq!(receipt.point_count, cloud.len() as u64);
324        assert!(receipt.tile_count > 1);
325
326        let document = std::fs::read_to_string(out_dir.join("tileset.json")).unwrap();
327        let tileset = parse_tileset_json(&document).unwrap();
328        assert_eq!(tileset.root.content.as_ref().unwrap().uri, "0.pnts");
329
330        for tile in 0..receipt.tile_count {
331            let pnts = std::fs::read(out_dir.join(format!("{tile}.pnts"))).unwrap();
332            assert!(decode_pnts(&pnts).is_ok());
333        }
334
335        let _ = std::fs::remove_dir_all(&out_dir);
336        let _ = std::fs::remove_file(&copc_path);
337    }
338
339    #[test]
340    fn max_level_bounds_export() {
341        let cloud = dense_grid_cloud(7_000);
342        let copc_path = std::env::temp_dir()
343            .join(format!("spatialrust_tiles3d_copc_lvl_{}.copc.laz", std::process::id()));
344        write_copc_file_with_params(
345            &copc_path,
346            &cloud,
347            &CopcWriterParams { max_points_per_node: 96, max_depth: 8 },
348        )
349        .unwrap();
350
351        let out_dir = std::env::temp_dir()
352            .join(format!("spatialrust_tiles3d_copc_lvl_out_{}", std::process::id()));
353        let receipt = export_copc_tileset(
354            &copc_path,
355            &out_dir,
356            &CopcTilesetOptions { max_level: Some(0), ..Default::default() },
357        )
358        .unwrap();
359        assert_eq!(receipt.tile_count, 1);
360        assert!(
361            receipt.point_count < cloud.len() as u64,
362            "level-0 export must expose only the root node chunk"
363        );
364
365        let _ = std::fs::remove_dir_all(&out_dir);
366        let _ = std::fs::remove_file(&copc_path);
367    }
368
369    #[test]
370    fn preserves_rgb_from_las_color() {
371        use spatialrust_core::{PointCloudBuilder, StandardSchemas};
372
373        let mut builder = PointCloudBuilder::new(StandardSchemas::point_xyzrgb());
374        for index in 0..2_000usize {
375            let x = (index % 31) as f32 - 15.0;
376            let y = ((index / 31) % 29) as f32 - 14.0;
377            let z = ((index / (31 * 29)) % 23) as f32 - 11.0;
378            let r = ((index % 256) << 8) as f32;
379            let g = (((index * 7) % 256) << 8) as f32;
380            let b = (((index * 13) % 256) << 8) as f32;
381            builder.push_point([x, y, z, r, g, b]).unwrap();
382        }
383        let cloud = builder.build().unwrap();
384
385        let copc_path = std::env::temp_dir()
386            .join(format!("spatialrust_tiles3d_copc_rgb_{}.copc.laz", std::process::id()));
387        write_copc_file(&copc_path, &cloud).unwrap();
388
389        let out_dir = std::env::temp_dir()
390            .join(format!("spatialrust_tiles3d_copc_rgb_out_{}", std::process::id()));
391        let receipt =
392            export_copc_tileset(&copc_path, &out_dir, &CopcTilesetOptions::default()).unwrap();
393        assert_eq!(receipt.point_count, cloud.len() as u64);
394
395        let mut rgb_points = 0usize;
396        for tile in 0..receipt.tile_count {
397            let pnts = std::fs::read(out_dir.join(format!("{tile}.pnts"))).unwrap();
398            let decoded = decode_pnts(&pnts).unwrap();
399            let rgb = decoded.rgb.as_ref().expect("color-bearing COPC must write RGB");
400            assert_eq!(rgb.len(), decoded.point_count() * 3);
401            rgb_points += decoded.point_count();
402        }
403        assert_eq!(rgb_points, cloud.len());
404
405        let _ = std::fs::remove_dir_all(&out_dir);
406        let _ = std::fs::remove_file(&copc_path);
407    }
408}
spatialrust_interchange/tiles3d/copc.rs

spatialrust_interchange/tiles3d/
copc.rs
