oxigdal_3d/pointcloud/
las.rs

1//! LAS/LAZ point cloud format support
2//!
3//! Provides comprehensive LAS 1.4 format reading and writing with LAZ compression support.
4
5use crate::error::{Error, Result};
6use rstar::{AABB, RTree, RTreeObject};
7use serde::{Deserialize, Serialize};
8use std::fs::File;
9use std::io::{BufReader, BufWriter};
10use std::path::Path;
11
12/// LAS point classification codes (ASPRS Standard)
13#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Serialize, Deserialize)]
14#[repr(u8)]
15pub enum Classification {
16    /// Never classified
17    NeverClassified = 0,
18    /// Unclassified
19    Unclassified = 1,
20    /// Ground
21    Ground = 2,
22    /// Low vegetation
23    LowVegetation = 3,
24    /// Medium vegetation
25    MediumVegetation = 4,
26    /// High vegetation
27    HighVegetation = 5,
28    /// Building
29    Building = 6,
30    /// Low point (noise)
31    LowPoint = 7,
32    /// Reserved
33    Reserved = 8,
34    /// Water
35    Water = 9,
36    /// Rail
37    Rail = 10,
38    /// Road surface
39    RoadSurface = 11,
40    /// Reserved (overlap)
41    ReservedOverlap = 12,
42    /// Wire - guard (shield)
43    WireGuard = 13,
44    /// Wire - conductor (phase)
45    WireConductor = 14,
46    /// Transmission tower
47    TransmissionTower = 15,
48    /// Wire-structure connector (insulator)
49    WireConnector = 16,
50    /// Bridge deck
51    BridgeDeck = 17,
52    /// High noise
53    HighNoise = 18,
54    /// Other/custom
55    Other(u8),
56}
57
58impl From<u8> for Classification {
59    fn from(value: u8) -> Self {
60        match value {
61            0 => Classification::NeverClassified,
62            1 => Classification::Unclassified,
63            2 => Classification::Ground,
64            3 => Classification::LowVegetation,
65            4 => Classification::MediumVegetation,
66            5 => Classification::HighVegetation,
67            6 => Classification::Building,
68            7 => Classification::LowPoint,
69            8 => Classification::Reserved,
70            9 => Classification::Water,
71            10 => Classification::Rail,
72            11 => Classification::RoadSurface,
73            12 => Classification::ReservedOverlap,
74            13 => Classification::WireGuard,
75            14 => Classification::WireConductor,
76            15 => Classification::TransmissionTower,
77            16 => Classification::WireConnector,
78            17 => Classification::BridgeDeck,
79            18 => Classification::HighNoise,
80            other => Classification::Other(other),
81        }
82    }
83}
84
85impl From<Classification> for u8 {
86    fn from(class: Classification) -> Self {
87        match class {
88            Classification::NeverClassified => 0,
89            Classification::Unclassified => 1,
90            Classification::Ground => 2,
91            Classification::LowVegetation => 3,
92            Classification::MediumVegetation => 4,
93            Classification::HighVegetation => 5,
94            Classification::Building => 6,
95            Classification::LowPoint => 7,
96            Classification::Reserved => 8,
97            Classification::Water => 9,
98            Classification::Rail => 10,
99            Classification::RoadSurface => 11,
100            Classification::ReservedOverlap => 12,
101            Classification::WireGuard => 13,
102            Classification::WireConductor => 14,
103            Classification::TransmissionTower => 15,
104            Classification::WireConnector => 16,
105            Classification::BridgeDeck => 17,
106            Classification::HighNoise => 18,
107            Classification::Other(val) => val,
108        }
109    }
110}
111
112/// RGB color values
113#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
114pub struct ColorRgb {
115    /// Red channel (0-65535)
116    pub red: u16,
117    /// Green channel (0-65535)
118    pub green: u16,
119    /// Blue channel (0-65535)
120    pub blue: u16,
121}
122
123/// RGB + NIR color values
124#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
125pub struct ColorRgbNir {
126    /// Red channel (0-65535)
127    pub red: u16,
128    /// Green channel (0-65535)
129    pub green: u16,
130    /// Blue channel (0-65535)
131    pub blue: u16,
132    /// Near-infrared channel (0-65535)
133    pub nir: u16,
134}
135
136/// 3D bounding box
137#[derive(Debug, Clone, Copy, PartialEq, Serialize, Deserialize)]
138pub struct Bounds3d {
139    /// Minimum X coordinate
140    pub min_x: f64,
141    /// Maximum X coordinate
142    pub max_x: f64,
143    /// Minimum Y coordinate
144    pub min_y: f64,
145    /// Maximum Y coordinate
146    pub max_y: f64,
147    /// Minimum Z coordinate
148    pub min_z: f64,
149    /// Maximum Z coordinate
150    pub max_z: f64,
151}
152
153impl Bounds3d {
154    /// Create new bounds
155    pub fn new(min_x: f64, max_x: f64, min_y: f64, max_y: f64, min_z: f64, max_z: f64) -> Self {
156        Self {
157            min_x,
158            max_x,
159            min_y,
160            max_y,
161            min_z,
162            max_z,
163        }
164    }
165
166    /// Check if bounds contain a point
167    pub fn contains(&self, x: f64, y: f64, z: f64) -> bool {
168        x >= self.min_x
169            && x <= self.max_x
170            && y >= self.min_y
171            && y <= self.max_y
172            && z >= self.min_z
173            && z <= self.max_z
174    }
175
176    /// Check if bounds intersect
177    pub fn intersects(&self, other: &Bounds3d) -> bool {
178        self.min_x <= other.max_x
179            && self.max_x >= other.min_x
180            && self.min_y <= other.max_y
181            && self.max_y >= other.min_y
182            && self.min_z <= other.max_z
183            && self.max_z >= other.min_z
184    }
185
186    /// Get bounds center
187    pub fn center(&self) -> (f64, f64, f64) {
188        (
189            (self.min_x + self.max_x) / 2.0,
190            (self.min_y + self.max_y) / 2.0,
191            (self.min_z + self.max_z) / 2.0,
192        )
193    }
194}
195
196/// Point format enumeration
197#[derive(Debug, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
198pub enum PointFormat {
199    /// Format 0: X, Y, Z, Intensity, Return Number, etc.
200    Format0 = 0,
201    /// Format 1: Format 0 + GPS Time
202    Format1 = 1,
203    /// Format 2: Format 0 + RGB
204    Format2 = 2,
205    /// Format 3: Format 1 + RGB
206    Format3 = 3,
207    /// Format 4: Format 1 + Wave Packets
208    Format4 = 4,
209    /// Format 5: Format 3 + Wave Packets
210    Format5 = 5,
211    /// Format 6: Format 1 + Extended (LAS 1.4)
212    Format6 = 6,
213    /// Format 7: Format 6 + RGB
214    Format7 = 7,
215    /// Format 8: Format 7 + NIR
216    Format8 = 8,
217    /// Format 9: Format 6 + Wave Packets
218    Format9 = 9,
219    /// Format 10: Format 8 + Wave Packets
220    Format10 = 10,
221}
222
223impl TryFrom<u8> for PointFormat {
224    type Error = Error;
225
226    fn try_from(value: u8) -> Result<Self> {
227        match value {
228            0 => Ok(PointFormat::Format0),
229            1 => Ok(PointFormat::Format1),
230            2 => Ok(PointFormat::Format2),
231            3 => Ok(PointFormat::Format3),
232            4 => Ok(PointFormat::Format4),
233            5 => Ok(PointFormat::Format5),
234            6 => Ok(PointFormat::Format6),
235            7 => Ok(PointFormat::Format7),
236            8 => Ok(PointFormat::Format8),
237            9 => Ok(PointFormat::Format9),
238            10 => Ok(PointFormat::Format10),
239            _ => Err(Error::UnsupportedPointFormat(value)),
240        }
241    }
242}
243
244/// 3D point with attributes
245#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
246pub struct Point {
247    /// X coordinate
248    pub x: f64,
249    /// Y coordinate
250    pub y: f64,
251    /// Z coordinate
252    pub z: f64,
253    /// Intensity
254    pub intensity: u16,
255    /// Return number
256    pub return_number: u8,
257    /// Number of returns
258    pub number_of_returns: u8,
259    /// Classification
260    pub classification: Classification,
261    /// Scan angle
262    pub scan_angle: i16,
263    /// User data
264    pub user_data: u8,
265    /// Point source ID
266    pub point_source_id: u16,
267    /// GPS time (optional)
268    pub gps_time: Option<f64>,
269    /// RGB color (optional)
270    pub color: Option<ColorRgb>,
271    /// NIR value (optional)
272    pub nir: Option<u16>,
273}
274
275impl Point {
276    /// Create a new point with minimum attributes
277    pub fn new(x: f64, y: f64, z: f64) -> Self {
278        Self {
279            x,
280            y,
281            z,
282            intensity: 0,
283            return_number: 1,
284            number_of_returns: 1,
285            classification: Classification::Unclassified,
286            scan_angle: 0,
287            user_data: 0,
288            point_source_id: 0,
289            gps_time: None,
290            color: None,
291            nir: None,
292        }
293    }
294
295    /// Check if point is ground
296    pub fn is_ground(&self) -> bool {
297        self.classification == Classification::Ground
298    }
299
300    /// Check if point is vegetation
301    pub fn is_vegetation(&self) -> bool {
302        matches!(
303            self.classification,
304            Classification::LowVegetation
305                | Classification::MediumVegetation
306                | Classification::HighVegetation
307        )
308    }
309
310    /// Check if point is building
311    pub fn is_building(&self) -> bool {
312        self.classification == Classification::Building
313    }
314
315    /// Distance to another point
316    pub fn distance_to(&self, other: &Point) -> f64 {
317        let dx = self.x - other.x;
318        let dy = self.y - other.y;
319        let dz = self.z - other.z;
320        (dx * dx + dy * dy + dz * dz).sqrt()
321    }
322
323    /// 2D distance (ignoring Z)
324    pub fn distance_2d(&self, other: &Point) -> f64 {
325        let dx = self.x - other.x;
326        let dy = self.y - other.y;
327        (dx * dx + dy * dy).sqrt()
328    }
329}
330
331// Implement RTreeObject for spatial indexing
332impl RTreeObject for Point {
333    type Envelope = AABB<[f64; 3]>;
334
335    fn envelope(&self) -> Self::Envelope {
336        AABB::from_point([self.x, self.y, self.z])
337    }
338}
339
340/// LAS header information
341#[derive(Debug, Clone, Serialize, Deserialize)]
342pub struct LasHeader {
343    /// Version (e.g., "1.4")
344    pub version: String,
345    /// Point format
346    pub point_format: PointFormat,
347    /// Number of points
348    pub point_count: u64,
349    /// Bounds
350    pub bounds: Bounds3d,
351    /// Scale factors
352    pub scale: (f64, f64, f64),
353    /// Offset values
354    pub offset: (f64, f64, f64),
355    /// System identifier
356    pub system_identifier: String,
357    /// Generating software
358    pub generating_software: String,
359}
360
361/// Point cloud data structure
362#[derive(Debug, Clone, Serialize, Deserialize)]
363pub struct PointCloud {
364    /// Header information
365    pub header: LasHeader,
366    /// Points
367    pub points: Vec<Point>,
368}
369
370impl PointCloud {
371    /// Create a new point cloud
372    pub fn new(header: LasHeader, points: Vec<Point>) -> Self {
373        Self { header, points }
374    }
375
376    /// Get number of points
377    pub fn len(&self) -> usize {
378        self.points.len()
379    }
380
381    /// Check if empty
382    pub fn is_empty(&self) -> bool {
383        self.points.is_empty()
384    }
385
386    /// Get bounds
387    pub fn bounds(&self) -> Bounds3d {
388        self.header.bounds
389    }
390
391    /// Filter points by classification
392    pub fn filter_by_classification(&self, class: Classification) -> Vec<&Point> {
393        self.points
394            .iter()
395            .filter(|p| p.classification == class)
396            .collect()
397    }
398
399    /// Filter ground points
400    pub fn ground_points(&self) -> Vec<&Point> {
401        self.filter_by_classification(Classification::Ground)
402    }
403
404    /// Calculate actual bounds from points
405    pub fn calculate_bounds(&self) -> Option<Bounds3d> {
406        if self.points.is_empty() {
407            return None;
408        }
409
410        let mut min_x = f64::INFINITY;
411        let mut max_x = f64::NEG_INFINITY;
412        let mut min_y = f64::INFINITY;
413        let mut max_y = f64::NEG_INFINITY;
414        let mut min_z = f64::INFINITY;
415        let mut max_z = f64::NEG_INFINITY;
416
417        for point in &self.points {
418            min_x = min_x.min(point.x);
419            max_x = max_x.max(point.x);
420            min_y = min_y.min(point.y);
421            max_y = max_y.max(point.y);
422            min_z = min_z.min(point.z);
423            max_z = max_z.max(point.z);
424        }
425
426        Some(Bounds3d::new(min_x, max_x, min_y, max_y, min_z, max_z))
427    }
428}
429
430/// Spatial index for fast point queries
431pub struct SpatialIndex {
432    tree: RTree<Point>,
433    points: Vec<Point>,
434}
435
436impl SpatialIndex {
437    /// Create a new spatial index from points
438    pub fn new(points: Vec<Point>) -> Self {
439        let tree = RTree::bulk_load(points.clone());
440        Self { tree, points }
441    }
442
443    /// Find nearest point
444    pub fn nearest(&self, x: f64, y: f64, z: f64) -> Option<&Point> {
445        let mut nearest: Option<&Point> = None;
446        let mut min_dist_sq = f64::INFINITY;
447
448        for point in &self.points {
449            let dx = point.x - x;
450            let dy = point.y - y;
451            let dz = point.z - z;
452            let dist_sq = dx * dx + dy * dy + dz * dz;
453
454            if dist_sq < min_dist_sq {
455                min_dist_sq = dist_sq;
456                nearest = Some(point);
457            }
458        }
459
460        nearest
461    }
462
463    /// Find k nearest points
464    pub fn nearest_k(&self, x: f64, y: f64, z: f64, k: usize) -> Vec<&Point> {
465        let mut dists: Vec<(usize, f64)> = self
466            .points
467            .iter()
468            .enumerate()
469            .map(|(i, p)| {
470                let dx = p.x - x;
471                let dy = p.y - y;
472                let dz = p.z - z;
473                let dist_sq = dx * dx + dy * dy + dz * dz;
474                (i, dist_sq)
475            })
476            .collect();
477
478        dists.sort_by(|a, b| a.1.partial_cmp(&b.1).unwrap_or(std::cmp::Ordering::Equal));
479
480        dists
481            .iter()
482            .take(k)
483            .map(|(i, _)| &self.points[*i])
484            .collect()
485    }
486
487    /// Find points within radius
488    pub fn within_radius(&self, x: f64, y: f64, z: f64, radius: f64) -> Vec<&Point> {
489        let radius_sq = radius * radius;
490        self.points
491            .iter()
492            .filter(|p| {
493                let dx = p.x - x;
494                let dy = p.y - y;
495                let dz = p.z - z;
496                dx * dx + dy * dy + dz * dz <= radius_sq
497            })
498            .collect()
499    }
500
501    /// Find points within bounds
502    pub fn within_bounds(&self, bounds: &Bounds3d) -> Vec<&Point> {
503        let aabb = AABB::from_corners(
504            [bounds.min_x, bounds.min_y, bounds.min_z],
505            [bounds.max_x, bounds.max_y, bounds.max_z],
506        );
507        self.tree.locate_in_envelope(&aabb).collect()
508    }
509}
510
511/// Point record for internal processing
512#[derive(Debug, Clone)]
513pub struct PointRecord {
514    /// Raw LAS point
515    pub raw: las::Point,
516}
517
518/// LAS file reader
519pub struct LasReader {
520    reader: las::Reader,
521    header: LasHeader,
522}
523
524impl LasReader {
525    /// Open a LAS/LAZ file
526    pub fn open<P: AsRef<Path>>(path: P) -> Result<Self> {
527        let file = File::open(path)?;
528        let buf_reader = BufReader::new(file);
529        let reader = las::Reader::new(buf_reader)?;
530
531        let las_header = reader.header();
532        let version = format!(
533            "{}.{}",
534            las_header.version().major,
535            las_header.version().minor
536        );
537        let fmt_u8 = las_header
538            .point_format()
539            .to_u8()
540            .map_err(|_| Error::UnsupportedPointFormat(0))?;
541        let point_format = PointFormat::try_from(fmt_u8)?;
542
543        let bounds = Bounds3d::new(
544            las_header.bounds().min.x,
545            las_header.bounds().max.x,
546            las_header.bounds().min.y,
547            las_header.bounds().max.y,
548            las_header.bounds().min.z,
549            las_header.bounds().max.z,
550        );
551
552        let header = LasHeader {
553            version,
554            point_format,
555            point_count: las_header.number_of_points(),
556            bounds,
557            scale: (
558                las_header.transforms().x.scale,
559                las_header.transforms().y.scale,
560                las_header.transforms().z.scale,
561            ),
562            offset: (
563                las_header.transforms().x.offset,
564                las_header.transforms().y.offset,
565                las_header.transforms().z.offset,
566            ),
567            system_identifier: las_header.system_identifier().to_string(),
568            generating_software: las_header.generating_software().to_string(),
569        };
570
571        Ok(Self { reader, header })
572    }
573
574    /// Get header
575    pub fn header(&self) -> &LasHeader {
576        &self.header
577    }
578
579    /// Read all points
580    pub fn read_all(&mut self) -> Result<PointCloud> {
581        let mut points = Vec::new();
582
583        for result in self.reader.points() {
584            let las_point = result?;
585            let point = Self::convert_point(&las_point)?;
586            points.push(point);
587        }
588
589        Ok(PointCloud::new(self.header.clone(), points))
590    }
591
592    /// Read points with limit
593    pub fn read_n(&mut self, n: usize) -> Result<Vec<Point>> {
594        let mut points = Vec::with_capacity(n);
595
596        for (i, result) in self.reader.points().enumerate() {
597            if i >= n {
598                break;
599            }
600            let las_point = result?;
601            let point = Self::convert_point(&las_point)?;
602            points.push(point);
603        }
604
605        Ok(points)
606    }
607
608    /// Convert LAS point to our Point structure
609    fn convert_point(las_point: &las::Point) -> Result<Point> {
610        let class_u8: u8 = las_point.classification.into();
611        let classification = Classification::from(class_u8);
612
613        let color = las_point.color.map(|c| ColorRgb {
614            red: c.red,
615            green: c.green,
616            blue: c.blue,
617        });
618
619        let nir = las_point.nir;
620
621        Ok(Point {
622            x: las_point.x,
623            y: las_point.y,
624            z: las_point.z,
625            intensity: las_point.intensity,
626            return_number: las_point.return_number,
627            number_of_returns: las_point.number_of_returns,
628            classification,
629            scan_angle: las_point.scan_angle as i16,
630            user_data: las_point.user_data,
631            point_source_id: las_point.point_source_id,
632            gps_time: las_point.gps_time,
633            color,
634            nir,
635        })
636    }
637}
638
639/// LAS file writer
640pub struct LasWriter {
641    writer: las::Writer<BufWriter<File>>,
642}
643
644impl LasWriter {
645    /// Create a new LAS file
646    pub fn create<P: AsRef<Path>>(path: P, header: &LasHeader) -> Result<Self> {
647        let mut builder = las::Builder::default();
648
649        // Set version
650        let version_parts: Vec<&str> = header.version.split('.').collect();
651        if version_parts.len() == 2 {
652            if let (Ok(major), Ok(minor)) = (
653                version_parts[0].parse::<u8>(),
654                version_parts[1].parse::<u8>(),
655            ) {
656                builder.version = las::Version::new(major, minor);
657            }
658        }
659
660        // Set point format
661        builder.point_format = las::point::Format::new(header.point_format as u8)
662            .map_err(|e| Error::Las(e.to_string()))?;
663
664        // Set transforms
665        builder.transforms = las::Vector {
666            x: las::Transform {
667                scale: header.scale.0,
668                offset: header.offset.0,
669            },
670            y: las::Transform {
671                scale: header.scale.1,
672                offset: header.offset.1,
673            },
674            z: las::Transform {
675                scale: header.scale.2,
676                offset: header.offset.2,
677            },
678        };
679
680        // Set system identifier and generating software
681        builder.system_identifier = header.system_identifier.clone();
682        builder.generating_software = header.generating_software.clone();
683
684        let file = File::create(path)?;
685        let buf_writer = BufWriter::new(file);
686        let writer = las::Writer::new(
687            buf_writer,
688            builder
689                .into_header()
690                .map_err(|e| Error::Las(e.to_string()))?,
691        )?;
692
693        Ok(Self { writer })
694    }
695
696    /// Write a point
697    pub fn write_point(&mut self, point: &Point) -> Result<()> {
698        let las_point = las::Point {
699            x: point.x,
700            y: point.y,
701            z: point.z,
702            intensity: point.intensity,
703            return_number: point.return_number,
704            number_of_returns: point.number_of_returns,
705            classification: las::point::Classification::new(point.classification.into())
706                .map_err(|e| Error::Las(e.to_string()))?,
707            scan_angle: point.scan_angle as f32,
708            user_data: point.user_data,
709            point_source_id: point.point_source_id,
710            gps_time: point.gps_time,
711            color: point.color.map(|c| las::Color {
712                red: c.red,
713                green: c.green,
714                blue: c.blue,
715            }),
716            nir: point.nir,
717            ..Default::default()
718        };
719
720        self.writer.write_point(las_point)?;
721        Ok(())
722    }
723
724    /// Write multiple points
725    pub fn write_points(&mut self, points: &[Point]) -> Result<()> {
726        for point in points {
727            self.write_point(point)?;
728        }
729        Ok(())
730    }
731
732    /// Finish writing and close file
733    pub fn close(mut self) -> Result<()> {
734        self.writer.close()?;
735        Ok(())
736    }
737}
738
739#[cfg(test)]
740mod tests {
741    use super::*;
742    use approx::assert_relative_eq;
743
744    #[test]
745    fn test_classification_conversion() {
746        assert_eq!(Classification::from(2u8), Classification::Ground);
747        assert_eq!(u8::from(Classification::Ground), 2);
748    }
749
750    #[test]
751    fn test_point_creation() {
752        let point = Point::new(1.0, 2.0, 3.0);
753        assert_relative_eq!(point.x, 1.0);
754        assert_relative_eq!(point.y, 2.0);
755        assert_relative_eq!(point.z, 3.0);
756        assert_eq!(point.classification, Classification::Unclassified);
757    }
758
759    #[test]
760    fn test_point_classification_checks() {
761        let mut point = Point::new(0.0, 0.0, 0.0);
762
763        point.classification = Classification::Ground;
764        assert!(point.is_ground());
765        assert!(!point.is_building());
766
767        point.classification = Classification::Building;
768        assert!(point.is_building());
769        assert!(!point.is_ground());
770
771        point.classification = Classification::HighVegetation;
772        assert!(point.is_vegetation());
773    }
774
775    #[test]
776    fn test_point_distance() {
777        let p1 = Point::new(0.0, 0.0, 0.0);
778        let p2 = Point::new(3.0, 4.0, 0.0);
779
780        assert_relative_eq!(p1.distance_to(&p2), 5.0);
781        assert_relative_eq!(p1.distance_2d(&p2), 5.0);
782    }
783
784    #[test]
785    fn test_bounds_contains() {
786        let bounds = Bounds3d::new(0.0, 10.0, 0.0, 10.0, 0.0, 10.0);
787
788        assert!(bounds.contains(5.0, 5.0, 5.0));
789        assert!(!bounds.contains(15.0, 5.0, 5.0));
790    }
791
792    #[test]
793    fn test_bounds_intersects() {
794        let b1 = Bounds3d::new(0.0, 10.0, 0.0, 10.0, 0.0, 10.0);
795        let b2 = Bounds3d::new(5.0, 15.0, 5.0, 15.0, 5.0, 15.0);
796        let b3 = Bounds3d::new(20.0, 30.0, 20.0, 30.0, 20.0, 30.0);
797
798        assert!(b1.intersects(&b2));
799        assert!(!b1.intersects(&b3));
800    }
801
802    #[test]
803    fn test_bounds_center() {
804        let bounds = Bounds3d::new(0.0, 10.0, 0.0, 10.0, 0.0, 10.0);
805        let (cx, cy, cz) = bounds.center();
806
807        assert_relative_eq!(cx, 5.0);
808        assert_relative_eq!(cy, 5.0);
809        assert_relative_eq!(cz, 5.0);
810    }
811
812    #[test]
813    fn test_point_cloud_creation() {
814        let header = LasHeader {
815            version: "1.4".to_string(),
816            point_format: PointFormat::Format0,
817            point_count: 0,
818            bounds: Bounds3d::new(0.0, 0.0, 0.0, 0.0, 0.0, 0.0),
819            scale: (0.01, 0.01, 0.01),
820            offset: (0.0, 0.0, 0.0),
821            system_identifier: "OxiGDAL".to_string(),
822            generating_software: "oxigdal-3d".to_string(),
823        };
824
825        let points = vec![Point::new(1.0, 2.0, 3.0), Point::new(4.0, 5.0, 6.0)];
826
827        let cloud = PointCloud::new(header, points);
828        assert_eq!(cloud.len(), 2);
829        assert!(!cloud.is_empty());
830    }
831
832    #[test]
833    fn test_point_cloud_bounds_calculation() {
834        let header = LasHeader {
835            version: "1.4".to_string(),
836            point_format: PointFormat::Format0,
837            point_count: 0,
838            bounds: Bounds3d::new(0.0, 0.0, 0.0, 0.0, 0.0, 0.0),
839            scale: (0.01, 0.01, 0.01),
840            offset: (0.0, 0.0, 0.0),
841            system_identifier: "OxiGDAL".to_string(),
842            generating_software: "oxigdal-3d".to_string(),
843        };
844
845        let points = vec![Point::new(0.0, 0.0, 0.0), Point::new(10.0, 10.0, 10.0)];
846
847        let cloud = PointCloud::new(header, points);
848        let bounds = cloud
849            .calculate_bounds()
850            .expect("Failed to calculate bounds");
851
852        assert_relative_eq!(bounds.min_x, 0.0);
853        assert_relative_eq!(bounds.max_x, 10.0);
854    }
855
856    #[test]
857    fn test_spatial_index() {
858        let points = vec![
859            Point::new(0.0, 0.0, 0.0),
860            Point::new(1.0, 1.0, 1.0),
861            Point::new(2.0, 2.0, 2.0),
862        ];
863
864        let index = SpatialIndex::new(points);
865
866        let nearest = index.nearest(0.5, 0.5, 0.5);
867        assert!(nearest.is_some());
868
869        let within = index.within_radius(0.0, 0.0, 0.0, 2.0);
870        assert!(!within.is_empty());
871    }
872}
oxigdal_3d/pointcloud/las.rs

oxigdal_3d/pointcloud/
las.rs
