diff --git a/Cargo.lock b/Cargo.lock index 6004c23..1cba836 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -99,6 +99,7 @@ name = "app" version = "0.0.1" dependencies = [ "bitcode", + "bytemuck", "cesiumtiles", "cesiumtiles-gltf", "cesiumtiles-gltf-json", @@ -106,6 +107,8 @@ dependencies = [ "clap", "env_logger", "glob", + "itertools", + "kv-extsort", "log", "pcd-core", "pcd-exporter", @@ -117,6 +120,7 @@ dependencies = [ "serde_json", "tempfile", "thiserror", + "tinymvt", ] [[package]] @@ -172,6 +176,20 @@ name = "bytemuck" version = "1.20.0" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "8b37c88a63ffd85d15b406896cc343916d7cf57838a847b3a6f2ca5d39a5695a" +dependencies = [ + "bytemuck_derive", +] + +[[package]] +name = "bytemuck_derive" +version = "1.8.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "bcfcc3cd946cb52f0bbfdbbcfa2f4e24f75ebb6c0e1002f7c25904fada18b9ec" +dependencies = [ + "proc-macro2", + "quote", + "syn", +] [[package]] name = "byteorder" @@ -301,6 +319,15 @@ version = "0.8.7" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "773648b94d0e5d620f64f280777445740e61fe701025087ec8b57f45c791888b" +[[package]] +name = "crossbeam-channel" +version = "0.5.13" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "33480d6946193aa8033910124896ca395333cae7e2d1113d1fef6c3272217df2" +dependencies = [ + "crossbeam-utils", +] + [[package]] name = "crossbeam-deque" version = "0.8.5" @@ -518,6 +545,19 @@ dependencies = [ "wasm-bindgen", ] +[[package]] +name = "kv-extsort" +version = "0.1.0" +source = "git+https://github.com/MIERUNE/kv-extsort-rs.git#6ac60e183cca6525876b24137cb1b822236e6e0d" +dependencies = [ + "bytemuck", + "crossbeam-channel", + "log", + "rayon", + "tempfile", + "thiserror", +] + [[package]] name = "las" version = "0.9.1" diff --git a/app/Cargo.toml b/app/Cargo.toml index 9f11d19..1ed95c9 100644 --- a/app/Cargo.toml +++ b/app/Cargo.toml @@ -24,3 +24,7 @@ serde_json = "1.0.133" glob = "0.3.1" tempfile = "3.14.0" bitcode = "0.6.3" +bytemuck = "1.20.0" +tinymvt = "0.0.1" +kv-extsort = { git = "https://github.com/MIERUNE/kv-extsort-rs.git" } +itertools = "0.13.0" diff --git a/app/src/main.rs b/app/src/main.rs index 89466d4..b603cf5 100644 --- a/app/src/main.rs +++ b/app/src/main.rs @@ -1,17 +1,28 @@ use std::collections::HashMap; +use std::convert::Infallible; use std::ffi::OsStr; use std::fs::File; -use std::io::{BufWriter, Write}; +use std::io::{BufWriter, Read as _, Write}; +use std::sync::{mpsc, Arc}; +use std::thread; use std::{ fs, path::{Path, PathBuf}, }; +use bitcode::{Decode, Encode}; +use bytemuck::{Pod, Zeroable}; use chrono::Local; use clap::Parser; use env_logger::Builder; use glob::glob; +use itertools::Itertools as _; use log::LevelFilter; +use pcd_parser::reader::csv::CsvPointReader; +use pcd_parser::reader::las::LasPointReader; +use pcd_parser::reader::PointReader; +use pcd_transformer::projection::transform_point; +use projection_transform::vshift::Jgd2011ToWgs84; use rayon::iter::{IntoParallelRefIterator as _, ParallelIterator as _}; use tempfile::tempdir; @@ -25,13 +36,9 @@ use pcd_exporter::{ gltf::generate_quantized_glb, tiling::{geometric_error, TileContent, TileTree}, }; -use pcd_parser::parsers::csv::CsvParserProvider; -use pcd_parser::parsers::{get_extension, Extension}; -use pcd_parser::parsers::{las::LasParserProvider, ParserProvider as _}; -use pcd_transformer::{ - builder::PointCloudTransformBuilder, runner::PointCloudTransformer, Transformer, -}; +use pcd_parser::parser::{get_extension, Extension}; use projection_transform::cartesian::geodetic_to_geocentric; +use tinymvt::tileid::hilbert; #[derive(Parser, Debug)] #[command( @@ -244,6 +251,78 @@ fn export_tiles_to_glb( Ok(tile_contents) } +#[derive(Debug, Copy, Clone, Ord, PartialOrd, Eq, PartialEq, Pod, Zeroable, Encode, Decode)] +#[repr(C)] +pub struct SortKey { + pub tile_id: u64, +} + +#[derive(Clone, Copy, Debug)] +pub enum TileIdMethod { + Hilbert, +} + +impl TileIdMethod { + pub fn zxy_to_id(&self, z: u8, x: u32, y: u32) -> u64 { + match self { + TileIdMethod::Hilbert => hilbert::zxy_to_id(z, x, y), + } + } + + pub fn id_to_zxy(&self, tile_id: u64) -> (u8, u32, u32) { + match self { + TileIdMethod::Hilbert => hilbert::id_to_zxy(tile_id), + } + } +} + +struct RunFileIterator { + files: std::vec::IntoIter, + current: Option>, +} + +impl RunFileIterator { + fn new(files: Vec) -> Self { + RunFileIterator { + files: files.into_iter(), + current: None, + } + } + + fn read_run_file(path: PathBuf) -> Result, Infallible> { + let file = File::open(path).unwrap(); + let mut buf_reader = std::io::BufReader::new(file); + let mut buffer = Vec::new(); + buf_reader.read_to_end(&mut buffer).unwrap(); + let data: Vec<(SortKey, Point)> = bitcode::decode(&buffer[..]).unwrap(); + Ok(data) + } +} + +impl Iterator for RunFileIterator { + type Item = (SortKey, Point); + + fn next(&mut self) -> Option { + loop { + if let Some(ref mut iter) = self.current { + if let Some(item) = iter.next() { + return Some(item); + } + } + + match self.files.next() { + Some(file) => { + let data = RunFileIterator::read_run_file(file).unwrap(); + self.current = Some(data.into_iter()); + } + None => { + return None; + } + } + } + } +} + fn main() { Builder::new() .format(|buf, record| { @@ -275,93 +354,129 @@ fn main() { let output_path = PathBuf::from(args.output); std::fs::create_dir_all(&output_path).unwrap(); - log::info!("start parsing..."); + let tmp_dir_path = tempdir().unwrap(); + + let min_zoom = args.min; + let max_zoom = args.max; + + log::info!( + "start parse and transform and tiling at max_zoom ({})...", + max_zoom + ); let start_local = std::time::Instant::now(); - let extension = check_and_get_extension(&input_files).unwrap(); - let parser = match extension { - Extension::Las => { - let las_parser_provider = LasParserProvider { - filenames: input_files, - epsg: args.epsg, - }; - let provider = las_parser_provider; - provider.get_parser() - } - Extension::Laz => { - let las_parser_provider = LasParserProvider { - filenames: input_files, - epsg: args.epsg, - }; - let provider = las_parser_provider; - provider.get_parser() - } - Extension::Csv => { - let csv_parser_provider = CsvParserProvider { - filenames: input_files, - epsg: args.epsg, - }; - let provider = csv_parser_provider; - provider.get_parser() - } - Extension::Txt => { - let csv_parser_provider = CsvParserProvider { - filenames: input_files, - epsg: args.epsg, - }; - let provider = csv_parser_provider; - provider.get_parser() + let jgd2wgs = Arc::new(Jgd2011ToWgs84::default()); + + // Number of points processed at once + let chunk_size = 10_000_000; + + let (tx, rx) = mpsc::channel::>(); + let handle = thread::spawn(move || { + let mut buffer = Vec::with_capacity(chunk_size); + + let extension = check_and_get_extension(&input_files).unwrap(); + let mut reader: Box = match extension { + Extension::Las => Box::new(LasPointReader::new(input_files).unwrap()), + Extension::Laz => Box::new(LasPointReader::new(input_files).unwrap()), + Extension::Csv => Box::new(CsvPointReader::new(input_files).unwrap()), + Extension::Txt => Box::new(CsvPointReader::new(input_files).unwrap()), + }; + + while let Ok(Some(p)) = reader.next_point() { + buffer.push(p); + if buffer.len() >= chunk_size { + if tx.send(buffer.clone()).is_err() { + break; + } + buffer = Vec::with_capacity(chunk_size); + } } - }; - // TODO: Allow each chunk to be retrieved - let point_cloud = match parser.parse() { - Ok(point_cloud) => point_cloud, - Err(e) => { - log::error!("Failed to parse point cloud: {:?}", e); - return; + if !buffer.is_empty() { + let _ = tx.send(buffer.clone()); + buffer.clear(); } - }; - log::info!("finish parsing in {:?}", start_local.elapsed()); + }); - log::info!("start transforming..."); - let output_epsg = 4979; - let start_local = std::time::Instant::now(); - let transform_builder = PointCloudTransformBuilder::new(output_epsg); - let transformer = PointCloudTransformer::new(Box::new(transform_builder)); + for (current_run_index, chunk) in rx.into_iter().enumerate() { + let mut keyed_points: Vec<(SortKey, Point)> = chunk + .into_iter() + .map(|p| { + let transformed = transform_point(p, args.epsg, &jgd2wgs); - let transformed = transformer.execute(point_cloud.clone()); - log::info!("Finish transforming in {:?}", start_local.elapsed()); + let tile_coords = + tiling::scheme::zxy_from_lng_lat(max_zoom, transformed.x, transformed.y); + let tile_id = + TileIdMethod::Hilbert.zxy_to_id(tile_coords.0, tile_coords.1, tile_coords.2); - let min_zoom = args.min; - let max_zoom = args.max; + (SortKey { tile_id }, transformed) + }) + .collect(); - let tmp_dir_path = tempdir().unwrap(); + keyed_points.sort_by_key(|(k, _)| k.tile_id); - log::info!("start tiling at max_zoom ({})...", max_zoom); + let run_file_path = tmp_dir_path + .path() + .join(format!("run_{}.bin", current_run_index)); + let file = fs::File::create(run_file_path).unwrap(); + let mut writer = std::io::BufWriter::new(file); + + let encoded = bitcode::encode(&keyed_points); + writer.write_all(&encoded).unwrap(); + } + + handle.join().expect("Reading thread panicked"); + + log::info!( + "Finish transforming and tiling in {:?}", + start_local.elapsed() + ); + + log::info!("start sorting..."); let start_local = std::time::Instant::now(); + let pattern = tmp_dir_path.path().join("run_*.bin"); + let run_files = glob::glob(pattern.to_str().unwrap()) + .unwrap() + .map(|r| r.unwrap()) + .collect::>(); + + let tile_id_iter = RunFileIterator::new(run_files); + + let config = kv_extsort::SortConfig::default().max_chunk_bytes(1024 * 1024 * 1024); + let sorted_iter = kv_extsort::sort( + tile_id_iter.map(|(key, point)| { + let encoded_point = bitcode::encode(&point); + std::result::Result::<_, Infallible>::Ok((key, encoded_point)) + }), + config, + ); - // calc tile coords for max_zoom - let tile_pairs: Vec<((u8, u32, u32), Point)> = transformed - .points - .par_iter() - .map(|p| { - let tile_coords = tiling::scheme::zxy_from_lng_lat(max_zoom, p.x, p.y); - (tile_coords, p.clone()) - }) - .collect(); + let grouped_iter = sorted_iter.chunk_by(|res| match res { + Ok((key, _)) => (false, *key), + Err(_) => (true, SortKey { tile_id: 0 }), + }); - // grouping by tile - let mut tile_map: HashMap<(u8, u32, u32), Vec> = HashMap::new(); - for (tile, pt) in tile_pairs { - tile_map.entry(tile).or_default().push(pt); - } + for ((_, key), group) in &grouped_iter { + let points = group + .into_iter() + .map(|r| r.map(|(_, p)| bitcode::decode::(&p).unwrap())) + .collect::, _>>() + .unwrap(); - // export to tile files - for (tile, pts) in tile_map { - write_points_to_tile(tmp_dir_path.path(), tile, &pts).unwrap(); - } + let tile_id = key.tile_id; + let tile = TileIdMethod::Hilbert.id_to_zxy(tile_id); - log::info!("Finish tiling at max_zoom in {:?}", start_local.elapsed()); + let (z, x, y) = tile; + let tile_path = tmp_dir_path.path().join(format!("{}/{}/{}.bin", z, x, y)); + + fs::create_dir_all(tile_path.parent().unwrap()).unwrap(); + + let file = fs::File::create(tile_path).unwrap(); + let mut writer = BufWriter::new(file); + + let encoded = bitcode::encode(&points); + writer.write_all(&encoded).unwrap(); + } + log::info!("Finish sorting in {:?}", start_local.elapsed()); log::info!("start zoom aggregation..."); let start_local = std::time::Instant::now(); diff --git a/pcd-exporter/examples/las_to_tiles.rs b/pcd-exporter/examples/las_to_tiles.rs index 10b8d6b..4ac0caf 100644 --- a/pcd-exporter/examples/las_to_tiles.rs +++ b/pcd-exporter/examples/las_to_tiles.rs @@ -12,7 +12,7 @@ use pcd_exporter::{ gltf::generate_quantized_glb, tiling::{geometric_error, TileContent, TileTree}, }; -use pcd_parser::parsers::{las::LasParserProvider, ParserProvider as _}; +use pcd_parser::parser::{las::LasParserProvider, ParserProvider as _}; use pcd_transformer::{ builder::PointCloudTransformBuilder, runner::PointCloudTransformer, Transformer, }; diff --git a/pcd-parser/examples/read_csv.rs b/pcd-parser/examples/read_csv.rs index 8175b16..e7dcdf7 100644 --- a/pcd-parser/examples/read_csv.rs +++ b/pcd-parser/examples/read_csv.rs @@ -1,6 +1,6 @@ use std::path::PathBuf; -use pcd_parser::parsers::{csv::CsvParserProvider, ParserProvider as _}; +use pcd_parser::parser::{csv::CsvParserProvider, ParserProvider as _}; fn main() { let provider = CsvParserProvider { diff --git a/pcd-parser/examples/read_las.rs b/pcd-parser/examples/read_las.rs index 363089d..613feb4 100644 --- a/pcd-parser/examples/read_las.rs +++ b/pcd-parser/examples/read_las.rs @@ -1,6 +1,6 @@ use std::path::PathBuf; -use pcd_parser::parsers::{las::LasParserProvider, ParserProvider as _}; +use pcd_parser::parser::{las::LasParserProvider, ParserProvider as _}; fn main() { let provider = LasParserProvider { diff --git a/pcd-parser/src/lib.rs b/pcd-parser/src/lib.rs index be756a0..ca52771 100644 --- a/pcd-parser/src/lib.rs +++ b/pcd-parser/src/lib.rs @@ -1 +1,2 @@ -pub mod parsers; +pub mod parser; +pub mod reader; diff --git a/pcd-parser/src/parsers/csv/mod.rs b/pcd-parser/src/parser/csv/mod.rs similarity index 100% rename from pcd-parser/src/parsers/csv/mod.rs rename to pcd-parser/src/parser/csv/mod.rs diff --git a/pcd-parser/src/parser/las/mod.rs b/pcd-parser/src/parser/las/mod.rs new file mode 100644 index 0000000..8149787 --- /dev/null +++ b/pcd-parser/src/parser/las/mod.rs @@ -0,0 +1,74 @@ +use std::{error::Error, path::PathBuf}; + +use las::Reader; + +use pcd_core::pointcloud::point::{Color, Point, PointAttributes, PointCloud}; +use projection_transform::crs::EpsgCode; + +use super::{Parser, ParserProvider}; + +pub struct LasParserProvider { + pub filenames: Vec, + pub epsg: EpsgCode, +} + +impl ParserProvider for LasParserProvider { + fn get_parser(&self) -> Box { + Box::new(LasParser { + filenames: self.filenames.clone(), + epsg: self.epsg, + }) + } +} + +pub struct LasParser { + pub filenames: Vec, + pub epsg: EpsgCode, +} + +impl Parser for LasParser { + fn parse(&self) -> Result> { + let mut points = Vec::new(); + + for f in self.filenames.iter() { + let mut reader = Reader::from_path(f.to_str().unwrap()).unwrap(); + for las_point in reader.points() { + let las_point = las_point.unwrap(); + + let color = las_point.color.map(|c| Color { + r: c.red, + g: c.green, + b: c.blue, + }); + + let attributes = PointAttributes { + intensity: Some(las_point.intensity), + return_number: Some(las_point.return_number), + classification: Some(format!("{:?}", las_point.classification)), + scanner_channel: Some(las_point.user_data), + scan_angle: Some(las_point.scan_angle), + user_data: Some(las_point.user_data), + point_source_id: Some(las_point.point_source_id), + gps_time: Some(las_point.gps_time.unwrap_or(0.0)), + }; + + let point = Point { + x: las_point.x, + y: las_point.y, + z: las_point.z, + color: color.unwrap_or(Color { + r: 65535, + g: 65535, + b: 65535, + }), + attributes, + }; + + points.push(point); + } + } + + let point_cloud = PointCloud::new(points, self.epsg); + Ok(point_cloud) + } +} diff --git a/pcd-parser/src/parsers/mod.rs b/pcd-parser/src/parser/mod.rs similarity index 100% rename from pcd-parser/src/parsers/mod.rs rename to pcd-parser/src/parser/mod.rs diff --git a/pcd-parser/src/parsers/las/mod.rs b/pcd-parser/src/parsers/las/mod.rs deleted file mode 100644 index 9026e01..0000000 --- a/pcd-parser/src/parsers/las/mod.rs +++ /dev/null @@ -1,95 +0,0 @@ -use std::{error::Error, path::PathBuf, sync::mpsc::channel, thread}; - -use las::Reader; - -use pcd_core::pointcloud::point::{Color, Point, PointAttributes, PointCloud}; -use projection_transform::crs::EpsgCode; - -use super::{Parser, ParserProvider}; - -pub struct LasParserProvider { - pub filenames: Vec, - pub epsg: EpsgCode, -} - -impl ParserProvider for LasParserProvider { - fn get_parser(&self) -> Box { - Box::new(LasParser { - filenames: self.filenames.clone(), - epsg: self.epsg, - }) - } -} - -pub struct LasParser { - pub filenames: Vec, - pub epsg: EpsgCode, -} - -impl Parser for LasParser { - fn parse(&self) -> Result> { - let mut points = Vec::new(); - - let (tx, rx) = channel(); - - let handles: Vec<_> = self - .filenames - .iter() - .cloned() - .map(|filename| { - let tx = tx.clone(); - thread::spawn(move || { - let mut reader = Reader::from_path(filename).unwrap(); - for las_point in reader.points() { - let las_point = las_point.unwrap(); - - let color = las_point.color.map(|c| Color { - r: c.red, - g: c.green, - b: c.blue, - }); - - let attributes = PointAttributes { - intensity: Some(las_point.intensity), - return_number: Some(las_point.return_number), - classification: Some(format!("{:?}", las_point.classification)), - scanner_channel: Some(las_point.user_data), - scan_angle: Some(las_point.scan_angle), - user_data: Some(las_point.user_data), - point_source_id: Some(las_point.point_source_id), - gps_time: Some(las_point.gps_time.unwrap_or(0.0)), - }; - - let point = Point { - x: las_point.x, - y: las_point.y, - z: las_point.z, - color: color.unwrap_or(Color { - r: 65535, - g: 65535, - b: 65535, - }), - attributes, - }; - - tx.send(point).unwrap(); - } - drop(tx); - }) - }) - .collect(); - - drop(tx); - - for point in rx { - points.push(point); - } - - for handle in handles { - handle.join().unwrap(); - } - - let point_cloud = PointCloud::new(points, self.epsg); - Ok(point_cloud) - } -} diff --git a/pcd-parser/src/reader/csv/mod.rs b/pcd-parser/src/reader/csv/mod.rs new file mode 100644 index 0000000..1971bc3 --- /dev/null +++ b/pcd-parser/src/reader/csv/mod.rs @@ -0,0 +1,244 @@ +use std::{collections::HashMap, error::Error, fs::File, io, path::PathBuf}; + +use csv::ReaderBuilder; +use pcd_core::pointcloud::point::{Color, Point, PointAttributes}; + +use super::PointReader; + +fn create_field_mapping( + headers: &csv::StringRecord, + has_headers: bool, +) -> Result, Box> { + let mut mapping = HashMap::new(); + + let attribute_names = vec![ + "x", + "y", + "z", + "intensity", + "return_number", + "classification", + "scanner_channel", + "scan_angle", + "user_data", + "point_source_id", + "gps_time", + "r", + "g", + "b", + "red", + "green", + "blue", + ]; + + if has_headers { + for (index, header) in headers.iter().enumerate() { + let normalized_header = header.to_lowercase().replace(['_', '-'], ""); + + for attr_name in &attribute_names { + let normalized_attr = attr_name.to_lowercase().replace(['_', '-'], ""); + if normalized_header == normalized_attr { + mapping.insert(attr_name.to_string(), index); + break; + } + } + } + } else { + for (index, attr_name) in attribute_names.iter().enumerate() { + mapping.insert(attr_name.to_string(), index); + } + } + + for attr_name in &["x", "y", "z"] { + if !mapping.contains_key(*attr_name) { + return Err(format!( + "Required attribute '{}' is missing in CSV headers or mapping.", + attr_name + ) + .into()); + } + } + + Ok(mapping) +} + +fn get_field_value<'a>( + record: &'a csv::StringRecord, + field_mapping: &HashMap, + field_name: &str, +) -> Option<&'a str> { + if let Some(&index) = field_mapping.get(field_name) { + let value = record.get(index); + value + } else { + None + } +} + +fn parse_optional_field( + record: &csv::StringRecord, + field_mapping: &HashMap, + field_name: &str, +) -> Option { + if let Some(value_str) = get_field_value(record, field_mapping, field_name) { + if value_str.trim().is_empty() { + None + } else { + Some(value_str.to_string()) + } + } else { + None + } +} + +pub struct CsvPointReader { + pub files: Vec, + pub current_file_index: usize, + pub current_reader: Option>, + pub field_mapping: HashMap, +} + +impl CsvPointReader { + pub fn new(files: Vec) -> io::Result { + let mut reader = CsvPointReader { + files, + current_file_index: 0, + current_reader: None, + field_mapping: HashMap::new(), + }; + + reader.open_next_file()?; + Ok(reader) + } + + fn open_next_file(&mut self) -> io::Result<()> { + if self.current_file_index < self.files.len() { + let path = &self.files[self.current_file_index]; + self.current_file_index += 1; + + let mut rdr = ReaderBuilder::new().has_headers(true).from_path(path)?; + + let headers = rdr.headers()?.clone(); + let has_headers = !headers.iter().all(|h| h.trim().is_empty()); + + let mapping = create_field_mapping(&headers, has_headers) + .map_err(|e| io::Error::new(io::ErrorKind::Other, e))?; + + self.field_mapping = mapping; + self.current_reader = Some(rdr); + Ok(()) + } else { + self.current_reader = None; + Ok(()) + } + } + + fn parse_point(&self, record: &csv::StringRecord) -> Result> { + let x_str = get_field_value(record, &self.field_mapping, "x").ok_or("Missing 'x' field")?; + let y_str = get_field_value(record, &self.field_mapping, "y").ok_or("Missing 'y' field")?; + let z_str = get_field_value(record, &self.field_mapping, "z").ok_or("Missing 'z' field")?; + + let x: f64 = x_str.parse()?; + let y: f64 = y_str.parse()?; + let z: f64 = z_str.parse()?; + + let r = parse_optional_field(record, &self.field_mapping, "r") + .or_else(|| parse_optional_field(record, &self.field_mapping, "red")) + .unwrap_or("65535".to_string()) + .parse::()? + .floor() as u16; + + let g = parse_optional_field(record, &self.field_mapping, "g") + .or_else(|| parse_optional_field(record, &self.field_mapping, "green")) + .unwrap_or("65535".to_string()) + .parse::()? + .floor() as u16; + + let b = parse_optional_field(record, &self.field_mapping, "b") + .or_else(|| parse_optional_field(record, &self.field_mapping, "blue")) + .unwrap_or("65535".to_string()) + .parse::()? + .floor() as u16; + + let color = Color { r, g, b }; + + let attributes = PointAttributes { + intensity: None, + return_number: None, + classification: None, + scanner_channel: None, + scan_angle: None, + user_data: None, + point_source_id: None, + gps_time: None, + }; + // TODO: To be implemented in the future + // let attributes = PointAttributes { + // intensity: parse_optional_field(&record, &field_mapping, "intensity") + // .unwrap_or(None), + // return_number: parse_optional_field(&record, &field_mapping, "return_number") + // .unwrap_or(None), + // classification: get_field_value(&record, &field_mapping, "classification") + // .map(|v| v.to_string()), + // scanner_channel: parse_optional_field( + // &record, + // &field_mapping, + // "scanner_channel", + // ) + // .unwrap_or(None), + // scan_angle: parse_optional_field(&record, &field_mapping, "scan_angle") + // .unwrap_or(None), + // user_data: parse_optional_field(&record, &field_mapping, "user_data") + // .unwrap_or(None), + // point_source_id: parse_optional_field( + // &record, + // &field_mapping, + // "point_source_id", + // ) + // .unwrap_or(None), + // gps_time: parse_optional_field(&record, &field_mapping, "gps_time") + // .unwrap_or(None), + // }; + + Ok(Point { + x, + y, + z, + color, + attributes, + }) + } +} + +impl PointReader for CsvPointReader { + fn next_point(&mut self) -> io::Result> { + loop { + if self.current_reader.is_none() { + self.open_next_file()?; + if self.current_reader.is_none() { + return Ok(None); + } + } + + let reader = self.current_reader.as_mut().unwrap(); + let mut record = csv::StringRecord::new(); + + match reader.read_record(&mut record) { + Ok(true) => match self.parse_point(&record) { + Ok(p) => return Ok(Some(p)), + Err(e) => { + eprintln!("Error parsing CSV point: {}", e); + return Err(io::Error::new(io::ErrorKind::Other, format!("{}", e))); + } + }, + Ok(false) => { + self.current_reader = None; + } + Err(e) => { + eprintln!("Error reading CSV record: {}", e); + return Err(io::Error::new(io::ErrorKind::Other, e)); + } + } + } + } +} diff --git a/pcd-parser/src/reader/las/mod.rs b/pcd-parser/src/reader/las/mod.rs new file mode 100644 index 0000000..7b27b89 --- /dev/null +++ b/pcd-parser/src/reader/las/mod.rs @@ -0,0 +1,138 @@ +use std::{ + fs::File, + io::{self, BufReader}, + path::PathBuf, +}; + +use las::Reader; +use pcd_core::pointcloud::point::{Color, Point, PointAttributes}; + +use super::PointReader; + +pub struct LasPointReader { + pub files: Vec, + pub current_file_index: usize, + pub current_reader: Option, +} + +impl LasPointReader { + pub fn new(files: Vec) -> io::Result { + Ok(Self { + files, + current_file_index: 0, + current_reader: None, + }) + } + + pub fn open_next_file(&mut self) -> io::Result<()> { + if self.current_file_index < self.files.len() { + let path = &self.files[self.current_file_index]; + let file = File::open(path).unwrap(); + let reader = Reader::new(BufReader::new(file)).unwrap(); + self.current_reader = Some(reader); + self.current_file_index += 1; + Ok(()) + } else { + self.current_reader = None; + Ok(()) + } + } + + pub fn convert_las_point(las_point: las::Point) -> Point { + let color = las_point + .color + .map(|c| Color { + r: c.red, + g: c.green, + b: c.blue, + }) + .unwrap_or(Color { + r: 65535, + g: 65535, + b: 65535, + }); + + let attributes = PointAttributes { + intensity: Some(las_point.intensity), + return_number: Some(las_point.return_number), + classification: Some(format!("{:?}", las_point.classification)), + scanner_channel: Some(las_point.user_data), + scan_angle: Some(las_point.scan_angle), + user_data: Some(las_point.user_data), + point_source_id: Some(las_point.point_source_id), + gps_time: Some(las_point.gps_time.unwrap_or(0.0)), + }; + + Point { + x: las_point.x, + y: las_point.y, + z: las_point.z, + color, + attributes, + } + } +} + +impl PointReader for LasPointReader { + fn next_point(&mut self) -> io::Result> { + loop { + if self.current_reader.is_none() { + self.open_next_file()?; + if self.current_reader.is_none() { + return Ok(None); + } + } + + let reader = self.current_reader.as_mut().unwrap(); + match reader.points().next() { + Some(Ok(las_point)) => { + let p = Self::convert_las_point(las_point); + return Ok(Some(p)); + } + Some(Err(e)) => { + eprintln!("Error reading LAS point: {}", e); + return Err(io::Error::new(io::ErrorKind::Other, e)); + } + None => { + self.current_reader = None; + } + } + } + } +} + +pub struct PointIterator { + pub reader: R, + pub chunk_size: usize, +} + +impl PointIterator { + pub fn new(reader: R, chunk_size: usize) -> Self { + Self { reader, chunk_size } + } +} + +impl Iterator for PointIterator { + type Item = Vec; + + fn next(&mut self) -> Option { + let mut buffer = Vec::with_capacity(self.chunk_size); + + for _ in 0..self.chunk_size { + match self.reader.next_point() { + Ok(Some(p)) => buffer.push(p), + Ok(None) => break, + Err(e) => { + eprintln!("Error reading point: {}", e); + break; + } + } + } + + if buffer.is_empty() { + None + } else { + Some(buffer) + } + } +} diff --git a/pcd-parser/src/reader/mod.rs b/pcd-parser/src/reader/mod.rs new file mode 100644 index 0000000..4bc7c83 --- /dev/null +++ b/pcd-parser/src/reader/mod.rs @@ -0,0 +1,9 @@ +pub mod csv; +pub mod las; + +use pcd_core::pointcloud::point::Point; +use std::io; + +pub trait PointReader { + fn next_point(&mut self) -> io::Result>; +} diff --git a/pcd-transformer/examples/pointcloud_transform.rs b/pcd-transformer/examples/pointcloud_transform.rs index dcfac32..8700e9e 100644 --- a/pcd-transformer/examples/pointcloud_transform.rs +++ b/pcd-transformer/examples/pointcloud_transform.rs @@ -2,7 +2,7 @@ use std::path::PathBuf; use pcd_core::pointcloud::point::{Point, PointCloud}; use pcd_exporter::gltf::generate_quantized_glb; -use pcd_parser::parsers::{las::LasParserProvider, ParserProvider as _}; +use pcd_parser::parser::{las::LasParserProvider, ParserProvider as _}; use pcd_transformer::{ builder::PointCloudTransformBuilder, runner::{PointCloudTransformer, Transformer}, diff --git a/pcd-transformer/src/lib.rs b/pcd-transformer/src/lib.rs index 8803851..2eda587 100644 --- a/pcd-transformer/src/lib.rs +++ b/pcd-transformer/src/lib.rs @@ -1,4 +1,5 @@ pub mod builder; +pub mod projection; pub mod runner; pub mod transform; diff --git a/pcd-transformer/src/projection.rs b/pcd-transformer/src/projection.rs new file mode 100644 index 0000000..23bddeb --- /dev/null +++ b/pcd-transformer/src/projection.rs @@ -0,0 +1,86 @@ +use pcd_core::pointcloud::point::Point; +use projection_transform::{crs::*, jprect::JPRZone, vshift::Jgd2011ToWgs84}; + +pub fn transform_point(point: Point, input_epsg: EpsgCode, jgd2wgs: &Jgd2011ToWgs84) -> Point { + match input_epsg { + EPSG_JGD2011_JPRECT_I + | EPSG_JGD2011_JPRECT_II + | EPSG_JGD2011_JPRECT_III + | EPSG_JGD2011_JPRECT_IV + | EPSG_JGD2011_JPRECT_V + | EPSG_JGD2011_JPRECT_VI + | EPSG_JGD2011_JPRECT_VII + | EPSG_JGD2011_JPRECT_VIII + | EPSG_JGD2011_JPRECT_IX + | EPSG_JGD2011_JPRECT_X + | EPSG_JGD2011_JPRECT_XI + | EPSG_JGD2011_JPRECT_XII + | EPSG_JGD2011_JPRECT_XIII + | EPSG_JGD2011_JPRECT_XIV + | EPSG_JGD2011_JPRECT_XV + | EPSG_JGD2011_JPRECT_XVI + | EPSG_JGD2011_JPRECT_XVII + | EPSG_JGD2011_JPRECT_XVIII + | EPSG_JGD2011_JPRECT_XIX + | EPSG_JGD2011_JPRECT_I_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_II_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_III_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_IV_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_V_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_VI_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_VII_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_VIII_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_IX_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_X_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_XI_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_XII_JGD2011_HEIGHT + | EPSG_JGD2011_JPRECT_XIII_JGD2011_HEIGHT => { + transform_from_jgd2011(point, Some(input_epsg), jgd2wgs) + } + _ => { + panic!("Unsupported input CRS: {}", input_epsg); + } + } +} + +fn rectangular_to_lnglat(x: f64, y: f64, height: f64, input_epsg: EpsgCode) -> (f64, f64, f64) { + let zone = JPRZone::from_epsg(input_epsg).unwrap(); + let proj = zone.projection(); + let (lng, lat, height) = proj.project_inverse(x, y, height).unwrap(); + (lng, lat, height) +} + +fn transform_from_jgd2011( + point: Point, + rectangular: Option, + jgd2wgs: &Jgd2011ToWgs84, +) -> Point { + let output_epsg = EPSG_WGS84_GEOGRAPHIC_3D; + + match output_epsg { + EPSG_WGS84_GEOGRAPHIC_3D => { + let x = point.x; + let y = point.y; + let z = point.z; + + let (lng, lat, height) = if let Some(input_epsg) = rectangular { + rectangular_to_lnglat(x, y, z, input_epsg) + } else { + (x, y, z) + }; + + let (lng, lat, height) = jgd2wgs.convert(lng, lat, height); + + Point { + x: lng, + y: lat, + z: height, + color: point.color.clone(), + attributes: point.attributes.clone(), + } + } + _ => { + panic!("Unsupported output CRS: {}", output_epsg); + } + } +}