diff --git a/Cargo.toml b/Cargo.toml index bfb638a..a600bd8 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -14,11 +14,19 @@ path = "tools/dla/main.rs" name = "point-cloud" path = "tools/point-cloud/main.rs" +[[bin]] +name = "transform" +path = "tools/transform/main.rs" + [dependencies] +clap = {version="4.0", features=["derive"]} +hex = "0.4" +geo = "0.23" kdtree = "0.6" log = "0.4" petgraph = "0.6" rand = "0.8" rand_distr = "0.4" stderrlog = "0.5" -clap = {version="4.0", features=["derive"]} +wkb = "0.7" +wkt = "0.10" diff --git a/generative/lib.rs b/generative/lib.rs index 6042870..fe4b163 100644 --- a/generative/lib.rs +++ b/generative/lib.rs @@ -1,2 +1,3 @@ pub mod dla; pub mod stdio; +pub mod wkio; diff --git a/generative/stdio.rs b/generative/stdio.rs index 84db36a..80c6490 100644 --- a/generative/stdio.rs +++ b/generative/stdio.rs @@ -2,9 +2,9 @@ use std::fs::File; use std::io::{BufReader, BufWriter, Read, Write}; use std::path::PathBuf; -pub fn get_output_writer(output: Option) -> Result>, String> { +pub fn get_output_writer(output: &Option) -> Result>, String> { match output { - Some(path) => match File::create(&path) { + Some(path) => match File::create(path) { Err(why) => Err(format!( "Couldn't create: '{}' because: '{}'", path.display(), @@ -16,9 +16,9 @@ pub fn get_output_writer(output: Option) -> Result) -> Result>, String> { +pub fn get_input_reader(input: &Option) -> Result>, String> { match input { - Some(path) => match File::open(&path) { + Some(path) => match File::open(path) { Err(why) => Err(format!( "Couldn't open: '{}' because: '{}'", path.display(), diff --git a/generative/wkio.rs b/generative/wkio.rs new file mode 100644 index 0000000..26f3733 --- /dev/null +++ b/generative/wkio.rs @@ -0,0 +1,377 @@ +use clap::ValueEnum; +use geo::Geometry; +use hex::{decode, encode_upper}; +use log::warn; +use std::io::{BufRead, BufReader, Lines, Read, Write}; +use std::str::FromStr; +use wkb::{geom_to_wkb, wkb_to_geom, write_geom_to_wkb}; +use wkt::{ToWkt, Wkt}; + +#[derive(Debug, Clone, ValueEnum)] +pub enum GeometryFormat { + /// One WKT geometry per line. Ignores trailing garbage; does not skip over leading garbage. + Wkt, + /// Stringified hex encoded WKB, one geometry per line + WkbHex, + /// Raw WKB bytes with no separator between geometries + WkbRaw, + // TODO: Flat? + // TODO: Splines? +} + +impl std::fmt::Display for GeometryFormat { + fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result { + match self { + // important: Should match clap::ValueEnum format + GeometryFormat::Wkt => write!(f, "wkt"), + GeometryFormat::WkbHex => write!(f, "wkb-hex"), + GeometryFormat::WkbRaw => write!(f, "wkb-raw"), + } + } +} + +pub fn read_geometries( + reader: R, + format: &GeometryFormat, +) -> Box>> +where + R: Read + 'static, +{ + match format { + GeometryFormat::Wkt => Box::new(read_wkt_geometries(reader)), + GeometryFormat::WkbHex => Box::new(read_wkbhex_geometries(reader)), + GeometryFormat::WkbRaw => Box::new(read_wkbraw_geometries(reader)), + } +} + +pub fn write_geometries(writer: W, geometries: G, format: &GeometryFormat) +where + W: Write, + G: IntoIterator>, +{ + match format { + GeometryFormat::Wkt => write_wkt_geometries(writer, geometries), + GeometryFormat::WkbHex => write_wkbhex_geometries(writer, geometries), + GeometryFormat::WkbRaw => write_wkbraw_geometries(writer, geometries), + } +} + +pub struct WktGeometries +where + R: Read, +{ + lines: Lines>, +} + +pub struct WkbHexGeometries +where + R: Read, +{ + lines: Lines>, +} + +pub struct WkbRawGeometries +where + R: Read, +{ + reader: BufReader, +} + +impl Iterator for WktGeometries +where + R: Read, +{ + type Item = Geometry; + + fn next(&mut self) -> Option { + match self.lines.next() { + Some(Ok(line)) => match Wkt::::from_str(line.as_str()) { + Ok(geometry) => match geometry.try_into() { + Ok(geometry) => Some(geometry), + Err(e) => { + warn!("Failed to convert '{}' to geo geometry: {:?}", line, e); + None + } + }, + Err(e) => { + warn!("Failed to parse '{}' as WKT: {:?}", line, e); + None + } + }, + Some(Err(e)) => { + warn!("Failed to read line: {:?}", e); + None + } + None => None, + } + } +} + +impl Iterator for WkbRawGeometries +where + R: Read, +{ + type Item = Geometry; + + fn next(&mut self) -> Option { + // This is the only way to tell if a BufRead is exhausted without using the nightly only + // unstable has_data_left() API. + match self.reader.fill_buf() { + Ok(buf) => { + if buf.is_empty() { + return None; + } + } + Err(e) => { + warn!("Failed to read WKB: {:?}", e); + return None; + } + } + + match wkb_to_geom(&mut self.reader) { + Ok(geom) => Some(geom), + Err(e) => { + warn!("Failed to parse WKB: {:?}", e); + None + } + } + } +} + +impl Iterator for WkbHexGeometries +where + R: Read, +{ + type Item = Geometry; + + fn next(&mut self) -> Option { + match self.lines.next() { + Some(line) => match line { + Ok(line) => match decode(line) { + Ok(buf) => match wkb_to_geom(&mut &buf[..]) { + Ok(geom) => Some(geom), + Err(e) => { + warn!("Failed to parse WKB(hex): {:?}", e); + None + } + }, + Err(e) => { + warn!("Failed to decode WKB(hex): {:?}", e); + None + } + }, + Err(e) => { + warn!("Failed to read WKB(hex) from line: {:?}", e); + None + } + }, + None => None, + } + } +} + +/// Return an iterator to the WKT geometries passed in through the given BufReader +/// +/// Expects one geometry per line (LF or CRLF). Parsing any given line ends after either the first +/// failure or the first geometry yielded, whichever comes first. That is, a line can have trailing +/// garbage, but not leading garbage. +fn read_wkt_geometries(reader: R) -> WktGeometries +where + R: Read, +{ + WktGeometries { + // TODO: Is there a nice way to implement whitespace-separated geometries? + lines: BufReader::new(reader).lines(), + } +} + +fn read_wkbhex_geometries(reader: R) -> WkbHexGeometries +where + R: Read, +{ + WkbHexGeometries { + lines: BufReader::new(reader).lines(), + } +} + +fn read_wkbraw_geometries(reader: R) -> WkbRawGeometries +where + R: Read, +{ + WkbRawGeometries { + reader: BufReader::new(reader), + } +} + +/// Write the given geometries with the given Writer in WKT format +/// +/// Each geometry will be written on its own line. +fn write_wkt_geometries(mut writer: W, geometries: G) +where + W: Write, + G: IntoIterator>, +{ + for geometry in geometries { + let wkt_geom = geometry.to_wkt(); + writeln!(writer, "{}", wkt_geom.item).expect("Writing failed"); + } +} + +fn write_wkbhex_geometries(mut writer: W, geometries: G) +where + W: Write, + G: IntoIterator>, +{ + for geom in geometries { + match geom_to_wkb(&geom) { + Ok(buffer) => { + writeln!(writer, "{}", encode_upper(buffer)).unwrap(); + } + Err(e) => { + warn!("Failed to serialize geometry to WKB: {:?}", e); + } + } + } +} + +fn write_wkbraw_geometries(mut writer: W, geometries: G) +where + W: Write, + G: IntoIterator>, +{ + for geom in geometries { + // TODO: What's this about the endianity byte? + if let Err(e) = write_geom_to_wkb(&geom, &mut writer) { + warn!("Failed to write geometry: {:?}", e); + } + } +} + +#[cfg(test)] +mod tests { + use super::*; + use geo::{Geometry, Point}; + + #[test] + fn test_read_simple_point() { + let input = b"POINT(1 1)"; + let mut geometries = read_wkt_geometries(&input[..]); + let geometry = geometries.next(); + assert_ne!(geometry, None); + + let geometry = geometry.unwrap(); + let point: Result, _> = geometry.try_into(); + assert!(point.is_ok()); + let point = point.unwrap(); + + let expected = Point::new(1.0, 1.0); + assert_eq!(point, expected); + } + + #[test] + fn test_empty() { + let input = b""; + let mut geometries = read_wkt_geometries(&input[..]); + assert_eq!(geometries.next(), None); + } + + #[test] + fn test_nothing_but_garbage() { + let input = b"garbage"; + let mut geometries = read_wkt_geometries(&input[..]); + assert_eq!(geometries.next(), None); + } + + #[test] + fn test_each_geometry_must_be_on_its_own_line() { + let input = b"POINT(1 1)\nPOINT(2 2)\rPOINT(3 3)\r\nPOINT(4 4)\nPOINT(5 5) POINT(6 6)\nPOINT(7 7)\tPOINT(8 8)"; + let geometries = read_wkt_geometries(&input[..]); + let actual: Vec> = geometries.collect(); + let expected = vec![ + Geometry::Point(Point::new(1.0, 1.0)), + Geometry::Point(Point::new(2.0, 2.0)), // fails to grab point 3 because it's separated by a single \r + Geometry::Point(Point::new(4.0, 4.0)), + Geometry::Point(Point::new(5.0, 5.0)), // fails to grab point 6 because it's separated by a space + Geometry::Point(Point::new(7.0, 7.0)), // fails to grab point 8 because it's separated by a tab + ]; + + assert_eq!(actual, expected); + } + + #[test] + fn test_wkb_single_input() { + let input_wkt = b"POINT(2 3.5)"; + let input_wkbhex = b"010100000000000000000000400000000000000C40"; + let input_wkbraw = decode(input_wkbhex).unwrap(); + + let expected = Geometry::Point(Point::new(2.0, 3.5)); + + let mut wkt_geometries = read_wkt_geometries(&input_wkt[..]); + assert_eq!(wkt_geometries.next().unwrap(), expected); + + let mut wkbraw_geometries = read_wkbraw_geometries(input_wkbraw.as_slice()); + assert_eq!(wkbraw_geometries.next().unwrap(), expected); + + let mut wkbhex_geometries = read_wkbhex_geometries(&input_wkbhex[..]); + assert_eq!(wkbhex_geometries.next().unwrap(), expected); + } + + #[test] + fn test_wkb_multi_input() { + let input_wkt = b"POINT(1 1)\nPOINT(2 3.5)"; + let input_wkbhex = b"0101000000000000000000F03F000000000000F03F\n010100000000000000000000400000000000000C40"; + let buffer: Vec = input_wkbhex[..] + .lines() + .flat_map(|l| decode(l.unwrap()).unwrap()) + .collect(); + + let expected1 = Geometry::Point(Point::new(1.0, 1.0)); + let expected2 = Geometry::Point(Point::new(2.0, 3.5)); + + let mut wkt_geometries = read_wkt_geometries(&input_wkt[..]); + assert_eq!(wkt_geometries.next().unwrap(), expected1); + assert_eq!(wkt_geometries.next().unwrap(), expected2); + + let mut wkbhex_geometries = read_wkbhex_geometries(&input_wkbhex[..]); + assert_eq!(wkbhex_geometries.next().unwrap(), expected1); + assert_eq!(wkbhex_geometries.next().unwrap(), expected2); + + let mut wkbraw_geometries = read_wkbraw_geometries(buffer.as_slice()); + assert_eq!(wkbraw_geometries.next().unwrap(), expected1); + assert_eq!(wkbraw_geometries.next().unwrap(), expected2); + } + + #[test] + fn test_wkbhex_output() { + let input_wkbhex = b"0101000000000000000000F03F000000000000F03F\n010100000000000000000000400000000000000C40\n"; + let geometries = read_wkbhex_geometries(&input_wkbhex[..]); + + let mut output_buffer = Vec::::new(); + write_wkbhex_geometries(&mut output_buffer, geometries); + + assert_eq!(output_buffer, input_wkbhex); + } + + #[test] + fn test_wkbraw_output() { + let input_wkbhex = b"0101000000000000000000F03F000000000000F03F\n010100000000000000000000400000000000000C40\n"; + let buffer: Vec = input_wkbhex[..] + .lines() + .flat_map(|l| decode(l.unwrap()).unwrap()) + .collect(); + let geometries = read_wkbraw_geometries(buffer.as_slice()); + + let mut output_buffer = Vec::::new(); + write_wkbraw_geometries(&mut output_buffer, geometries); + + assert_eq!(output_buffer, buffer); + } + + #[test] + fn test_cant_parse_3d_sad_face() { + let wkt = b"POINT Z(1 2 3)"; + let geometries = read_wkt_geometries(&wkt[..]); + + assert_eq!(geometries.count(), 0); + } +} diff --git a/tools/point-cloud/main.rs b/tools/point-cloud/main.rs index 186bb14..fb7d88e 100644 --- a/tools/point-cloud/main.rs +++ b/tools/point-cloud/main.rs @@ -85,7 +85,7 @@ fn main() { eprintln!("Generating {} points with seed {}", num_points, seed); let points = generate(num_points as usize, args.domain, &mut rng); - let mut writer = get_output_writer(args.output).unwrap(); + let mut writer = get_output_writer(&args.output).unwrap(); for point in points { writeln!( writer, diff --git a/tools/transform/cmdline.rs b/tools/transform/cmdline.rs new file mode 100644 index 0000000..d256ad7 --- /dev/null +++ b/tools/transform/cmdline.rs @@ -0,0 +1,72 @@ +use clap::{Parser, ValueEnum}; +use std::path::PathBuf; + +use generative::wkio::GeometryFormat; + +#[derive(Debug, Clone, ValueEnum)] +pub enum TransformCenter { + /// Center the affine transform on (0, 0) + Origin, + /// Center the transform on the center of each geometry's bounding box + EachGeometry, + /// Center the transform on the center of the entire collection's bounding box + WholeCollection, +} + +/// Perform transformations on 2D geometries +/// +/// Transformations are applied in the order: +/// +/// 1. rotation +/// 2. scale +/// 3. offset +/// 4. skew +/// +/// If you want to apply transformations in any other order, you can chain invocations of this +/// command, specifying only one transformation per invocation. +/// +/// If you want to apply transformations to 3D geometries, they must first be projected to 2D using +/// the project.py tool. +#[derive(Debug, Parser)] +#[clap(name = "transform", verbatim_doc_comment)] +pub struct CmdlineOptions { + /// Increase logging verbosity. Defaults to ERROR level. + #[clap(short, long, action = clap::ArgAction::Count)] + pub verbosity: u8, + + /// Output file to write result to. Defaults to stdout. + #[clap(short, long)] + pub output: Option, + + /// Output geometry format. + #[clap(short = 'O', long, default_value_t = GeometryFormat::Wkt)] + pub output_format: GeometryFormat, + + /// Input file to read input from. Defaults to stdin. + #[clap(short, long)] + pub input: Option, + + /// Input geometry format. + #[clap(short = 'I', long, default_value_t = GeometryFormat::Wkt)] + pub input_format: GeometryFormat, + + /// How to center the affine transformation + #[clap(long, default_value = "origin")] + pub center: TransformCenter, + + /// Degrees CCW rotation, applied before any other transformation + #[clap(short, long, default_value_t = 0.0)] + pub rotation: f64, + + /// The (x, y) multiplicative scale, applied after rotation + #[clap(short = 's', long, number_of_values = 2)] + pub scale: Option>, + + /// The (x, y) additive offset, applied after scale + #[clap(short = 't', long, number_of_values = 2)] + pub offset: Option>, + + /// Degrees (x, y) skew, applied after offset + #[clap(short = 'S', long, number_of_values = 2)] + pub skew: Option>, +} diff --git a/tools/transform/main.rs b/tools/transform/main.rs new file mode 100644 index 0000000..d5a3cf5 --- /dev/null +++ b/tools/transform/main.rs @@ -0,0 +1,103 @@ +mod cmdline; +use cmdline::{CmdlineOptions, TransformCenter}; + +use clap::Parser; +use generative::stdio::{get_input_reader, get_output_writer}; +use generative::wkio::{read_geometries, write_geometries}; +use geo::{coord, AffineOps, AffineTransform, BoundingRect, Coord, Geometry, Rect}; +use stderrlog::ColorChoice; +use wkt::ToWkt; + +fn build_transform(args: &CmdlineOptions, center: Coord) -> AffineTransform { + let mut transform = AffineTransform::rotate(args.rotation, center); + if let Some(scale) = &args.scale { + // The clap parser guarantees that the Vec used for scale, offset, and skew have + // exactly 2 values. + transform = transform.scaled(scale[0], scale[1], center); + } + if let Some(offset) = &args.offset { + transform = transform.translated(offset[0], offset[1]); + } + if let Some(skew) = &args.skew { + transform = transform.skewed(skew[0], skew[1], center); + } + + transform +} + +fn main() { + let args = cmdline::CmdlineOptions::parse(); + + stderrlog::new() + .verbosity(args.verbosity as usize + 1) // Default to WARN level. + .color(ColorChoice::Auto) + .init() + .expect("Failed to initialize stderrlog"); + + let reader = get_input_reader(&args.input).unwrap(); + let writer = get_output_writer(&args.output).unwrap(); + let geometries = read_geometries(reader, &args.input_format); // lazily loaded + + match args.center { + TransformCenter::Origin => { + let center = coord! {x:0.0, y: 0.0}; + let transform = build_transform(&args, center); + let transformed = geometries.map(|geom| geom.affine_transform(&transform)); + write_geometries(writer, transformed, &args.output_format); + } + TransformCenter::EachGeometry => { + let transformed = geometries.map(|geom| { + let center = geom + .bounding_rect() + .unwrap_or_else(|| { + panic!( + "Geometry '{}' didn't have a bounding rectangle", + geom.to_wkt() + ) + }) + .center(); + let transform = build_transform(&args, center); + geom.affine_transform(&transform) + }); + write_geometries(writer, transformed, &args.output_format); + } + // more expensive for large numbers of geometries (has to load all of them into RAM before + // performing the transformations) + TransformCenter::WholeCollection => { + // Read geometries into memory so we can loop over them twice + let geometries: Vec> = geometries.collect(); + // Calculate the center of the bounding box; needed to build the AffineTransform + let mut min_x = f64::MAX; + let mut min_y = f64::MAX; + let mut max_x = f64::MIN; + let mut max_y = f64::MIN; + for geom in geometries.iter() { + let temp = geom.bounding_rect().unwrap_or_else(|| { + panic!( + "Geometry '{}' didn't have a bounding rectangle", + geom.to_wkt() + ) + }); + + let min = temp.min(); + let max = temp.max(); + + min_x = min_x.min(min.x); + min_y = min_y.min(min.y); + max_x = max_x.max(max.x); + max_y = max_y.max(max.y); + } + let rect = Rect::new(coord! {x:min_x, y:min_y}, coord! {x:max_x, y:max_y}); + let center = rect.center(); + let transform = build_transform(&args, center); + + // Instead of applying the transformation in-place all at once _and then_ writing the + // results, we lazily perform the transformation so that we can pipeline the + // transformation and the serialization. + let transformed = geometries + .into_iter() + .map(|geom| geom.affine_transform(&transform)); + write_geometries(writer, transformed, &args.output_format); + } + } +}