diff --git a/Cargo.toml b/Cargo.toml index 23bfdc4..c06d494 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -14,17 +14,19 @@ license = "MIT OR Apache-2.0" [dependencies] geojson = { version = ">=0.16, <=0.24", optional = true } -geo-types= { version = "0.7" } +geo-types = { version = "0.7" } lazy_static = "1.0" serde_json = { version = "^1.0", optional = true } +serde = { version = "1.0", optional = true } rustc-hash = "1.0" slab = "0.4" +num-traits = "0.2" [dev-dependencies] serde_json = "^1.0" [features] -geojson = ["dep:geojson", "dep:serde_json"] +geojson = ["dep:geojson", "dep:serde_json", "dep:serde"] f32 = [] [package.metadata.docs.rs] diff --git a/benches/bench.rs b/benches/bench.rs index 0e623df..ac86ad8 100644 --- a/benches/bench.rs +++ b/benches/bench.rs @@ -2,8 +2,7 @@ extern crate contour; extern crate test; -use contour::contour_rings; -use contour::ContourBuilder; +use contour::{contour_rings, ContourBuilder}; use test::{black_box, Bencher}; #[rustfmt::skip] diff --git a/examples/ex.rs b/examples/ex.rs index f8274c0..03adc0f 100644 --- a/examples/ex.rs +++ b/examples/ex.rs @@ -1,7 +1,9 @@ use contour::{ContourBuilder, Float}; use geojson::{FeatureCollection, GeoJson}; -use std::fs::File; -use std::io::{BufWriter, Write}; +use std::{ + fs::File, + io::{BufWriter, Write}, +}; fn main() { let pot_pop_fr = include_str!("../tests/fixtures/pot_pop_fr.json"); @@ -46,7 +48,8 @@ fn main() { let features = contours .iter() .map(|contour| contour.to_geojson()) - .collect::>(); + .collect::, _>>() + .unwrap(); let geojson_str = GeoJson::from(FeatureCollection { bbox: None, @@ -91,7 +94,8 @@ fn main() { let features = contours .iter() .map(|contour| contour.to_geojson()) - .collect::>(); + .collect::, _>>() + .unwrap(); let geojson_str = GeoJson::from(FeatureCollection { bbox: None, diff --git a/src/band.rs b/src/band.rs index b41249b..6e760f2 100644 --- a/src/band.rs +++ b/src/band.rs @@ -1,36 +1,38 @@ -use crate::Float; +use crate::{Float, GridValue}; use geo_types::MultiPolygon; /// An isoband has the geometry and min / max values of a contour ring, built by [`ContourBuilder`]. #[derive(Debug, Clone)] -pub struct Band { +pub struct Band { pub(crate) geometry: MultiPolygon, - pub(crate) min_v: Float, - pub(crate) max_v: Float, + pub(crate) min_v: V, + pub(crate) max_v: V, } -impl Band { +impl Band { /// Borrow the [`MultiPolygon`](geo_types::MultiPolygon) geometry of this contour. pub fn geometry(&self) -> &MultiPolygon { &self.geometry } /// Get the owned polygons and thresholds (min and max) of this band. - pub fn into_inner(self) -> (MultiPolygon, Float, Float) { + pub fn into_inner(self) -> (MultiPolygon, V, V) { (self.geometry, self.min_v, self.max_v) } /// Get the minimum value used to construct this band. - pub fn min_v(&self) -> Float { + pub fn min_v(&self) -> V { self.min_v } /// Get the maximum value used to construct this band. - pub fn max_v(&self) -> Float { + pub fn max_v(&self) -> V { self.max_v } +} - #[cfg(feature = "geojson")] +#[cfg(feature = "geojson")] +impl Band { /// Convert the band to a struct from the `geojson` crate. /// /// To get a string representation, call to_geojson().to_string(). @@ -53,21 +55,21 @@ impl Band { /// # 0., 0., 0., 0., 0., 0., 0., 0., 0., 0. /// ], &[0.5, 1.5, 2.5]).unwrap(); /// - /// let geojson_string = contours[0].to_geojson().to_string(); + /// let geojson_string = contours[0].to_geojson().unwrap().to_string(); /// /// assert_eq!(&geojson_string[0..27], r#"{"geometry":{"coordinates":"#); /// ``` - pub fn to_geojson(&self) -> geojson::Feature { + pub fn to_geojson(&self) -> serde_json::Result { let mut properties = geojson::JsonObject::with_capacity(2); - properties.insert("min_v".to_string(), self.min_v.into()); - properties.insert("max_v".to_string(), self.max_v.into()); + properties.insert("min_v".to_string(), serde_json::to_value(self.min_v)?); + properties.insert("max_v".to_string(), serde_json::to_value(self.max_v)?); - geojson::Feature { + Ok(geojson::Feature { bbox: None, geometry: Some(geojson::Geometry::from(self.geometry())), id: None, properties: Some(properties), foreign_members: None, - } + }) } } diff --git a/src/contour.rs b/src/contour.rs index 37ed75e..f013a3d 100644 --- a/src/contour.rs +++ b/src/contour.rs @@ -1,30 +1,32 @@ -use crate::Float; +use crate::{Float, GridValue}; use geo_types::MultiPolygon; /// A contour has the geometry and threshold of a contour ring, built by [`ContourBuilder`]. #[derive(Debug, Clone)] -pub struct Contour { +pub struct Contour { pub(crate) geometry: MultiPolygon, - pub(crate) threshold: Float, + pub(crate) threshold: V, } -impl Contour { +impl Contour { /// Borrow the [`MultiPolygon`](geo_types::MultiPolygon) geometry of this contour. pub fn geometry(&self) -> &MultiPolygon { &self.geometry } /// Get the owned polygons and threshold of this contour. - pub fn into_inner(self) -> (MultiPolygon, Float) { + pub fn into_inner(self) -> (MultiPolygon, V) { (self.geometry, self.threshold) } /// Get the threshold used to construct this contour. - pub fn threshold(&self) -> Float { + pub fn threshold(&self) -> V { self.threshold } +} - #[cfg(feature = "geojson")] +#[cfg(feature = "geojson")] +impl Contour { /// Convert the contour to a struct from the `geojson` crate. /// /// To get a string representation, call to_geojson().to_string(). @@ -47,20 +49,23 @@ impl Contour { /// # 0., 0., 0., 0., 0., 0., 0., 0., 0., 0. /// ], &[0.5]).unwrap(); /// - /// let geojson_string = contours[0].to_geojson().to_string(); + /// let geojson_string = contours[0].to_geojson().unwrap().to_string(); /// /// assert_eq!(&geojson_string[0..27], r#"{"geometry":{"coordinates":"#); /// ``` - pub fn to_geojson(&self) -> geojson::Feature { + pub fn to_geojson(&self) -> serde_json::Result { let mut properties = geojson::JsonObject::with_capacity(1); - properties.insert("threshold".to_string(), self.threshold.into()); + properties.insert( + "threshold".to_string(), + serde_json::to_value(self.threshold)?, + ); - geojson::Feature { + Ok(geojson::Feature { bbox: None, geometry: Some(geojson::Geometry::from(self.geometry())), id: None, properties: Some(properties), foreign_members: None, - } + }) } } diff --git a/src/contourbuilder.rs b/src/contourbuilder.rs index 830a548..388ee97 100644 --- a/src/contourbuilder.rs +++ b/src/contourbuilder.rs @@ -1,7 +1,9 @@ -use crate::area::{area, contains}; -use crate::error::{new_error, ErrorKind, Result}; -use crate::isoringbuilder::IsoRingBuilder; -use crate::{Band, Contour, Float, Line, Ring}; +use crate::{ + area::{area, contains}, + error::{new_error, ErrorKind, Result}, + isoringbuilder::IsoRingBuilder, + Band, Contour, Error, Float, GridValue, Line, Ring, +}; use geo_types::{LineString, MultiLineString, MultiPolygon, Polygon}; use rustc_hash::FxHashMap; @@ -43,10 +45,10 @@ impl ContourBuilder { dx, dy, smooth, - x_origin: 0., - y_origin: 0., - x_step: 1., - y_step: 1., + x_origin: 0.0, + y_origin: 0.0, + x_step: 1.0, + y_step: 1.0, } } @@ -74,32 +76,41 @@ impl ContourBuilder { self } - fn smooth_linear(&self, ring: &mut Ring, values: &[Float], value: Float) { + fn smooth_linear(&self, ring: &mut Ring, values: &[V], value: V) -> Result<()> { let dx = self.dx; let dy = self.dy; let len_values = values.len(); - ring.iter_mut() - .map(|point| { - let x = point.x; - let y = point.y; - let xt = x.trunc() as usize; - let yt = y.trunc() as usize; - let mut v0; - let ix = yt * dx + xt; - if ix < len_values { - let v1 = values[ix]; - if x > 0.0 && x < (dx as Float) && (xt as Float - x).abs() < Float::EPSILON { - v0 = values[yt * dx + xt - 1]; - point.x = x + (value - v0) / (v1 - v0) - 0.5; - } - if y > 0.0 && y < (dy as Float) && (yt as Float - y).abs() < Float::EPSILON { - v0 = values[(yt - 1) * dx + xt]; - point.y = y + (value - v0) / (v1 - v0) - 0.5; - } + ring.iter_mut().try_for_each(|point| { + let x = point.x; + let y = point.y; + let xt = x.trunc() as usize; + let yt = y.trunc() as usize; + let mut v0; + let ix = yt * dx + xt; + if ix < len_values { + let v1 = values[ix]; + if x > 0.0 && x < (dx as Float) && (xt as Float - x).abs() < Float::EPSILON { + v0 = values[yt * dx + xt - 1]; + point.x += num_traits::cast::<_, Float>(value - v0) + .ok_or_else(|| new_error(ErrorKind::BadCast))? + / num_traits::cast::<_, Float>(v1 - v0) + .ok_or_else(|| new_error(ErrorKind::BadCast))? + - 0.5; } - }) - .for_each(drop); + if y > 0.0 && y < (dy as Float) && (yt as Float - y).abs() < Float::EPSILON { + v0 = values[(yt - 1) * dx + xt]; + point.y += num_traits::cast::<_, Float>(value - v0) + .ok_or_else(|| new_error(ErrorKind::BadCast))? + / num_traits::cast::<_, Float>(v1 - v0) + .ok_or_else(|| new_error(ErrorKind::BadCast))? + - 0.5; + } + } + Ok::<_, Error>(()) + })?; + + Ok(()) } /// Computes isolines according the given input `values` and the given `thresholds`. @@ -111,7 +122,7 @@ impl ContourBuilder { /// /// * `values` - The slice of values to be used. /// * `thresholds` - The slice of thresholds values to be used. - pub fn lines(&self, values: &[Float], thresholds: &[Float]) -> Result> { + pub fn lines(&self, values: &[V], thresholds: &[V]) -> Result>> { if values.len() != self.dx * self.dy { return Err(new_error(ErrorKind::BadDimension)); } @@ -122,19 +133,19 @@ impl ContourBuilder { .collect() } - fn line( + fn line( &self, - values: &[Float], - threshold: Float, + values: &[V], + threshold: V, isoring: &mut IsoRingBuilder, - ) -> Result { + ) -> Result> { let mut result = isoring.compute(values, threshold)?; let mut linestrings = Vec::new(); - result.drain(..).for_each(|mut ring| { + for mut ring in result.drain(..) { // Smooth the ring if needed if self.smooth { - self.smooth_linear(&mut ring, values, threshold); + self.smooth_linear(&mut ring, values, threshold)?; } // Compute the polygon coordinates according to the grid properties if needed if (self.x_origin, self.y_origin) != (0.0, 0.0) @@ -146,9 +157,10 @@ impl ContourBuilder { }); } linestrings.push(LineString(ring)); - }); + } + Ok(Line { - geometry: MultiLineString::(linestrings), + geometry: MultiLineString(linestrings), threshold, }) } @@ -162,7 +174,11 @@ impl ContourBuilder { /// /// * `values` - The slice of values to be used. /// * `thresholds` - The slice of thresholds values to be used. - pub fn contours(&self, values: &[Float], thresholds: &[Float]) -> Result> { + pub fn contours( + &self, + values: &[V], + thresholds: &[V], + ) -> Result>> { if values.len() != self.dx * self.dy { return Err(new_error(ErrorKind::BadDimension)); } @@ -173,19 +189,19 @@ impl ContourBuilder { .collect() } - fn contour( + fn contour( &self, - values: &[Float], - threshold: Float, + values: &[V], + threshold: V, isoring: &mut IsoRingBuilder, - ) -> Result { + ) -> Result> { let (mut polygons, mut holes) = (Vec::new(), Vec::new()); let mut result = isoring.compute(values, threshold)?; - result.drain(..).for_each(|mut ring| { + for mut ring in result.drain(..) { // Smooth the ring if needed if self.smooth { - self.smooth_linear(&mut ring, values, threshold); + self.smooth_linear(&mut ring, values, threshold)?; } // Compute the polygon coordinates according to the grid properties if needed if (self.x_origin, self.y_origin) != (0.0, 0.0) @@ -197,11 +213,11 @@ impl ContourBuilder { }); } if area(&ring) > 0.0 { - polygons.push(Polygon::::new(LineString::new(ring), vec![])) + polygons.push(Polygon::new(LineString::new(ring), vec![])) } else { holes.push(LineString::new(ring)); } - }); + } holes.drain(..).for_each(|hole| { for polygon in &mut polygons { @@ -213,7 +229,7 @@ impl ContourBuilder { }); Ok(Contour { - geometry: MultiPolygon::(polygons), + geometry: MultiPolygon(polygons), threshold, }) } @@ -228,7 +244,7 @@ impl ContourBuilder { /// * `values` - The slice of values to be used. /// * `thresholds` - The slice of thresholds values to be used /// (have to be equal to or greater than 2). - pub fn isobands(&self, values: &[Float], thresholds: &[Float]) -> Result> { + pub fn isobands(&self, values: &[V], thresholds: &[V]) -> Result>> { // We will compute rings as previously, but we will // iterate over the contours in pairs and use the paths from the lower threshold // and the path from the upper threshold to create the isoband. @@ -250,7 +266,7 @@ impl ContourBuilder { .map(|mut ring| { // Smooth the ring if needed if self.smooth { - self.smooth_linear(&mut ring, values, *threshold); + self.smooth_linear(&mut ring, values, *threshold)?; } ring.dedup(); // Compute the polygon coordinates according to the grid properties if needed @@ -262,13 +278,13 @@ impl ContourBuilder { point.y = point.y * self.y_step + self.y_origin; }); } - ring + Ok(ring) }) - .filter(|ring| ring.len() > 3) - .collect::>(); + .filter(|ring| ring.as_ref().map(|v| v.len() > 3).unwrap_or(true)) + .collect::>>()?; Ok((rings, *threshold)) }) - .collect::, Float)>>>()?; + .collect::, V)>>>()?; // We now have the rings for each isolines for all the given thresholds, // we can iterate over them in pairs to compute the isobands. @@ -281,7 +297,7 @@ impl ContourBuilder { }) .collect::>(); - let mut bands: Vec = Vec::new(); + let mut bands: Vec> = Vec::new(); // Reconstruction of the polygons b.into_iter().for_each(|(rings, min_v, max_v)| { let mut rings_and_area = rings @@ -314,7 +330,7 @@ impl ContourBuilder { for (i, (ring, _)) in rings_and_area.into_iter().enumerate() { if *enclosed_by_n.get(&i).unwrap() % 2 == 0 { - polygons.push(Polygon::::new(ring.into(), vec![])); + polygons.push(Polygon::new(ring.into(), vec![])); } else { interior_rings.push(ring.into()); } @@ -331,7 +347,7 @@ impl ContourBuilder { polygons.reverse(); bands.push(Band { - geometry: MultiPolygon::(polygons), + geometry: MultiPolygon(polygons), min_v: *min_v, max_v: *max_v, }); diff --git a/src/error.rs b/src/error.rs index 37b6f7c..173d5ef 100644 --- a/src/error.rs +++ b/src/error.rs @@ -1,6 +1,4 @@ -use std::error::Error as StdError; -use std::fmt; -use std::result; +use std::{error::Error as StdError, fmt, result}; /// A crate private constructor for `Error`. pub(crate) fn new_error(kind: ErrorKind) -> Error { @@ -32,6 +30,7 @@ impl Error { pub enum ErrorKind { BadDimension, Unexpected, + BadCast, #[cfg(feature = "geojson")] JsonError(serde_json::error::Error), } @@ -48,6 +47,7 @@ impl StdError for Error { match *self.0 { ErrorKind::BadDimension => None, ErrorKind::Unexpected => None, + ErrorKind::BadCast => None, #[cfg(feature = "geojson")] ErrorKind::JsonError(ref err) => Some(err), } @@ -62,6 +62,7 @@ impl fmt::Display for Error { "The length of provided values doesn't match the (dx, dy) dimensions of the grid" ), ErrorKind::Unexpected => write!(f, "Unexpected error while computing contours"), + ErrorKind::BadCast => write!(f, "Failed to cast grid value to Float"), #[cfg(feature = "geojson")] ErrorKind::JsonError(ref err) => err.fmt(f), } diff --git a/src/isoringbuilder.rs b/src/isoringbuilder.rs index 4024028..3429c78 100644 --- a/src/isoringbuilder.rs +++ b/src/isoringbuilder.rs @@ -1,5 +1,8 @@ -use crate::error::{new_error, ErrorKind, Result}; -use crate::{Float, Pt, Ring}; +use crate::{ + error::{new_error, ErrorKind, Result}, + Float, GridValue, Pt, Ring, +}; +use geo_types::Coord; use lazy_static::lazy_static; use rustc_hash::FxHashMap; use slab::Slab; @@ -49,9 +52,10 @@ struct Fragment { /// * `threshold` - The threshold value. /// * `dx` - The number of columns in the grid. /// * `dy` - The number of rows in the grid. -pub fn contour_rings( - values: &[Float], - threshold: Float, + +pub fn contour_rings( + values: &[V], + threshold: V, dx: usize, dy: usize, ) -> Result> { @@ -94,7 +98,7 @@ impl IsoRingBuilder { /// /// * `values` - The slice of values to be used. /// * `threshold` - The threshold value to use. - pub fn compute(&mut self, values: &[Float], threshold: Float) -> Result> { + pub fn compute(&mut self, values: &[V], threshold: V) -> Result> { macro_rules! case_stitch { ($ix:expr, $x:ident, $y:ident, $result:expr) => { CASES[$ix] @@ -169,7 +173,7 @@ impl IsoRingBuilder { (point.x as usize) * 2 + (point.y as usize) * (self.dx + 1usize) * 4 } - // Stitchs segments to rings. + // Stitches segments to rings. fn stitch( &mut self, line: &[Vec], @@ -177,11 +181,11 @@ impl IsoRingBuilder { y: i64, result: &mut Vec, ) -> Result<()> { - let start = Pt { + let start = Coord { x: line[0][0] + x as Float, y: line[0][1] + y as Float, }; - let end = Pt { + let end = Coord { x: line[1][0] + x as Float, y: line[1][1] + y as Float, }; diff --git a/src/lib.rs b/src/lib.rs index 551ef19..ba14918 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -49,7 +49,7 @@ //! "properties": {"threshold": 0.5}, //! }); //! -//! assert_eq!(res[0].to_geojson(), std::convert::TryFrom::try_from(output).unwrap()); +//! assert_eq!(res[0].to_geojson().unwrap(), std::convert::TryFrom::try_from(output).unwrap()); //! ``` //! //! [`contour_rings`]: fn.contour_rings.html @@ -63,23 +63,26 @@ mod error; mod isoringbuilder; mod line; +pub trait GridValue: PartialOrd + Copy + Num + NumCast {} +impl GridValue for T where T: PartialOrd + Copy + Num + NumCast {} + #[cfg(feature = "f32")] pub type Float = f32; #[cfg(not(feature = "f32"))] pub type Float = f64; -#[cfg(feature = "f32")] -pub type Pt = geo_types::Coord; -#[cfg(not(feature = "f32"))] -pub type Pt = geo_types::Coord; +pub type Pt = geo_types::Coord; pub type Ring = Vec; -pub use crate::band::Band; -pub use crate::contour::Contour; -pub use crate::contourbuilder::ContourBuilder; -pub use crate::error::{Error, ErrorKind, Result}; -pub use crate::isoringbuilder::contour_rings; -pub use crate::line::Line; +pub use crate::{ + band::Band, + contour::Contour, + contourbuilder::ContourBuilder, + error::{Error, ErrorKind, Result}, + isoringbuilder::contour_rings, + line::Line, +}; +use num_traits::{Num, NumCast}; #[cfg(test)] mod tests { @@ -559,7 +562,7 @@ mod tests { 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0. ], &[0.5]).unwrap(); - match res[0].to_geojson().geometry.unwrap().value { + match res[0].to_geojson().unwrap().geometry.unwrap().value { geojson::Value::MultiPolygon(p) => { assert_eq!( p, diff --git a/src/line.rs b/src/line.rs index 967534c..c3a0cd9 100644 --- a/src/line.rs +++ b/src/line.rs @@ -1,30 +1,32 @@ -use crate::Float; +use crate::{Float, GridValue}; use geo_types::MultiLineString; /// A line has the geometry and threshold of a contour ring, built by [`ContourBuilder`]. #[derive(Debug, Clone)] -pub struct Line { +pub struct Line { pub(crate) geometry: MultiLineString, - pub(crate) threshold: Float, + pub(crate) threshold: V, } -impl Line { +impl Line { /// Borrow the [`MultiLineString`](geo_types::MultiLineString) geometry of this contour. pub fn geometry(&self) -> &MultiLineString { &self.geometry } /// Get the owned lines and threshold of this contour. - pub fn into_inner(self) -> (MultiLineString, Float) { + pub fn into_inner(self) -> (MultiLineString, V) { (self.geometry, self.threshold) } /// Get the threshold used to construct this isoline. - pub fn threshold(&self) -> Float { + pub fn threshold(&self) -> V { self.threshold } +} - #[cfg(feature = "geojson")] +#[cfg(feature = "geojson")] +impl Line { /// Convert the line to a struct from the `geojson` crate. /// /// To get a string representation, call to_geojson().to_string(). @@ -47,20 +49,23 @@ impl Line { /// # 0., 0., 0., 0., 0., 0., 0., 0., 0., 0. /// ], &[0.5]).unwrap(); /// - /// let geojson_string = contours[0].to_geojson().to_string(); + /// let geojson_string = contours[0].to_geojson().unwrap().to_string(); /// /// assert_eq!(&geojson_string[0..27], r#"{"geometry":{"coordinates":"#); /// ``` - pub fn to_geojson(&self) -> geojson::Feature { + pub fn to_geojson(&self) -> serde_json::Result { let mut properties = geojson::JsonObject::with_capacity(1); - properties.insert("threshold".to_string(), self.threshold.into()); + properties.insert( + "threshold".to_string(), + serde_json::to_value(self.threshold)?, + ); - geojson::Feature { + Ok(geojson::Feature { bbox: None, geometry: Some(geojson::Geometry::from(self.geometry())), id: None, properties: Some(properties), foreign_members: None, - } + }) } }