diff --git a/.gitignore b/.gitignore index 008b35d..e11e9ad 100644 --- a/.gitignore +++ b/.gitignore @@ -2,5 +2,6 @@ /cover /deps /doc +/priv erl_crash.dump *.ez diff --git a/lib/exshape.ex b/lib/exshape.ex index 6423a3f..aef4776 100644 --- a/lib/exshape.ex +++ b/lib/exshape.ex @@ -8,10 +8,10 @@ defmodule Exshape do defp open_file(c, size), do: File.stream!(c, [], size) - defp zip(nil, nil), do: [] - defp zip(nil, d), do: Dbf.read(d) - defp zip(s, nil), do: Shp.read(s) - defp zip(s, d), do: Stream.zip(Shp.read(s), Dbf.read(d)) + defp zip(nil, nil, _), do: [] + defp zip(nil, d, _), do: Dbf.read(d) + defp zip(s, nil, opts), do: Shp.read(s, opts) + defp zip(s, d, opts), do: Stream.zip(Shp.read(s, opts), Dbf.read(d)) defp unzip!(path, cwd, false), do: :zip.extract(to_charlist(path), cwd: cwd) defp unzip!(path, cwd, true) do @@ -78,12 +78,13 @@ defmodule Exshape do end open_file(Path.join(cwd, file), size) end - }) + }, + opts) end end @spec from_filesystem(Filesystem.t) :: [layer] - def from_filesystem(fs) do + def from_filesystem(fs, opts \\ []) do fs.list.() |> Enum.filter(&keep_file?/1) |> Enum.map(fn {:zip_file, filename, _, _, _, _} -> filename end) @@ -110,7 +111,8 @@ defmodule Exshape do # zip up the unzipped shp and dbf components stream = zip( shp && fs.stream.(shp), - dbf && fs.stream.(dbf) + dbf && fs.stream.(dbf), + opts ) {Path.basename(root), prj_contents, stream} diff --git a/lib/exshape/shp.ex b/lib/exshape/shp.ex index 329c636..19e3bfa 100644 --- a/lib/exshape/shp.ex +++ b/lib/exshape/shp.ex @@ -1,5 +1,9 @@ defmodule Exshape.Shp do + require Rustler + use Rustler, otp_app: :exshape, crate: :exshape_shape, mode: :release + defmodule State do + @enforce_keys [:nest_polygon] defstruct mode: :header, shape_type: nil, emit: [], @@ -7,7 +11,8 @@ defmodule Exshape.Shp do item: nil, part_index: 0, measures: [], - z_values: [] + z_values: [], + nest_polygon: nil end @magic_nodata_num :math.pow(10, 38) * -1 @@ -111,8 +116,8 @@ defmodule Exshape.Shp do %{s | measures: [], z_values: []} end - defp emit(s, %Polygon{} = p) do - %{s | mode: :record_header, emit: [%{p | points: nest_polygon(p)} | s.emit]} + defp emit(%State{nest_polygon: nest_polygon} = s, %Polygon{} = p) do + %{s | mode: :record_header, emit: [%{p | points: nest_polygon.(p)} | s.emit]} end defp emit(s, %Polyline{} = p) do @@ -134,9 +139,9 @@ defmodule Exshape.Shp do %{s | mode: :record_header, emit: [polylinem | s.emit]} |> reset_unzipped end - defp emit(s, %PolygonM{} = pm) do + defp emit(%State{nest_polygon: nest_polygon} = s, %PolygonM{} = pm) do p = zip_measures(pm, s) - polylinem = %{p | points: nest_polygon(p)} + polylinem = %{p | points: nest_polygon.(p)} %{s | mode: :record_header, emit: [polylinem | s.emit]} |> reset_unzipped end @@ -158,12 +163,12 @@ defmodule Exshape.Shp do %{s | mode: :record_header, emit: [polylinez | s.emit]} |> reset_unzipped end - defp emit(s, %PolygonZ{} = pz) do + defp emit(%State{nest_polygon: nest_polygon} = s, %PolygonZ{} = pz) do p = pz |> zip_measures(s) |> zip_zvals(s) - polygonz = %{p | points: nest_polygon(p)} + polygonz = %{p | points: nest_polygon.(p)} %{s | mode: :record_header, emit: [polygonz | s.emit]} |> reset_unzipped end @@ -218,8 +223,13 @@ defmodule Exshape.Shp do parts end + def native_nest_polygon(p) do + {:ok, r} = native_nest_polygon_impl(unflatten_parts(p)) + r + end + defp native_nest_polygon_impl(_p), do: throw :nif_not_loaded - def nest_polygon(p) do + def beam_nest_polygon(p) do {polys, holes} = unflatten_parts(p) |> Enum.split_with(&is_clockwise?/1) Enum.reduce(holes, Enum.map(polys, fn p -> [p] end), fn hole, polys -> @@ -754,8 +764,13 @@ defmodule Exshape.Shp do |> Stream.run ``` """ - def read(byte_stream) do - Stream.transform(byte_stream, {<<>>, %State{}}, fn bin, {buf, state} -> + def read(byte_stream, opts \\ []) do + native = Keyword.get(opts, :native, true) + + state = %State{ + nest_polygon: if(native, do: &native_nest_polygon/1, else: &beam_nest_polygon/1) + } + Stream.transform(byte_stream, {<<>>, state}, fn bin, {buf, state} -> case do_read(state, buf <> bin) do {_, %State{mode: :done}} = s -> {:halt, s} {buf, %State{emit: emit} = s} -> {Enum.reverse(emit), {buf, %{s | emit: []}}} diff --git a/mix.exs b/mix.exs index ab5df76..89c4c38 100644 --- a/mix.exs +++ b/mix.exs @@ -33,7 +33,7 @@ defmodule Exshape.Mixfile do # Type "mix help compile.app" for more information def application do # Specify extra applications you'll use from Erlang/Elixir - [extra_applications: [:logger]] + [extra_applications: [:logger, :rustler]] end # Dependencies can be Hex packages: @@ -49,7 +49,8 @@ defmodule Exshape.Mixfile do [ {:elixir_uuid, "~> 1.2"}, {:ex_doc, ">= 0.0.0", only: :dev}, - {:poison, "~> 3.1", only: :test} + {:poison, "~> 3.1", only: :test}, + {:rustler, "~> 0.22.0"}, ] end end diff --git a/mix.lock b/mix.lock index 646fa60..bbf56eb 100644 --- a/mix.lock +++ b/mix.lock @@ -1,10 +1,12 @@ %{ - "earmark": {:hex, :earmark, "1.2.5", "4d21980d5d2862a2e13ec3c49ad9ad783ffc7ca5769cf6ff891a4553fbaae761", [:mix], [], "hexpm"}, - "elixir_uuid": {:hex, :elixir_uuid, "1.2.0", "ff26e938f95830b1db152cb6e594d711c10c02c6391236900ddd070a6b01271d", [:mix], [], "hexpm"}, - "ex_doc": {:hex, :ex_doc, "0.19.1", "519bb9c19526ca51d326c060cb1778d4a9056b190086a8c6c115828eaccea6cf", [:mix], [{:earmark, "~> 1.1", [hex: :earmark, repo: "hexpm", optional: false]}, {:makeup_elixir, "~> 0.7", [hex: :makeup_elixir, repo: "hexpm", optional: false]}], "hexpm"}, - "makeup": {:hex, :makeup, "0.5.1", "966c5c2296da272d42f1de178c1d135e432662eca795d6dc12e5e8787514edf7", [:mix], [{:nimble_parsec, "~> 0.2.2", [hex: :nimble_parsec, repo: "hexpm", optional: false]}], "hexpm"}, - "makeup_elixir": {:hex, :makeup_elixir, "0.8.0", "1204a2f5b4f181775a0e456154830524cf2207cf4f9112215c05e0b76e4eca8b", [:mix], [{:makeup, "~> 0.5.0", [hex: :makeup, repo: "hexpm", optional: false]}, {:nimble_parsec, "~> 0.2.2", [hex: :nimble_parsec, repo: "hexpm", optional: false]}], "hexpm"}, - "nimble_parsec": {:hex, :nimble_parsec, "0.2.2", "d526b23bdceb04c7ad15b33c57c4526bf5f50aaa70c7c141b4b4624555c68259", [:mix], [], "hexpm"}, - "poison": {:hex, :poison, "3.1.0", "d9eb636610e096f86f25d9a46f35a9facac35609a7591b3be3326e99a0484665", [:mix], [], "hexpm"}, + "earmark": {:hex, :earmark, "1.2.5", "4d21980d5d2862a2e13ec3c49ad9ad783ffc7ca5769cf6ff891a4553fbaae761", [:mix], [], "hexpm", "c57508ddad47dfb8038ca6de1e616e66e9b87313220ac5d9817bc4a4dc2257b9"}, + "elixir_uuid": {:hex, :elixir_uuid, "1.2.0", "ff26e938f95830b1db152cb6e594d711c10c02c6391236900ddd070a6b01271d", [:mix], [], "hexpm", "e4d6e26434471761ed45a3545239da87af7b70904dd4442a55f87d06b137c56b"}, + "ex_doc": {:hex, :ex_doc, "0.19.1", "519bb9c19526ca51d326c060cb1778d4a9056b190086a8c6c115828eaccea6cf", [:mix], [{:earmark, "~> 1.1", [hex: :earmark, repo: "hexpm", optional: false]}, {:makeup_elixir, "~> 0.7", [hex: :makeup_elixir, repo: "hexpm", optional: false]}], "hexpm", "dc87f778d8260da0189a622f62790f6202af72f2f3dee6e78d91a18dd2fcd137"}, + "makeup": {:hex, :makeup, "0.5.1", "966c5c2296da272d42f1de178c1d135e432662eca795d6dc12e5e8787514edf7", [:mix], [{:nimble_parsec, "~> 0.2.2", [hex: :nimble_parsec, repo: "hexpm", optional: false]}], "hexpm", "259748a45dfcf5f49765a7c29c9594791c82de23e22d7a3e6e59533fe8e8935b"}, + "makeup_elixir": {:hex, :makeup_elixir, "0.8.0", "1204a2f5b4f181775a0e456154830524cf2207cf4f9112215c05e0b76e4eca8b", [:mix], [{:makeup, "~> 0.5.0", [hex: :makeup, repo: "hexpm", optional: false]}, {:nimble_parsec, "~> 0.2.2", [hex: :nimble_parsec, repo: "hexpm", optional: false]}], "hexpm", "393d17c5a648e3b30522b2a4743bd1dc3533e1227c8c2823ebe8c3a8e5be5913"}, + "nimble_parsec": {:hex, :nimble_parsec, "0.2.2", "d526b23bdceb04c7ad15b33c57c4526bf5f50aaa70c7c141b4b4624555c68259", [:mix], [], "hexpm", "4ababf5c44164f161872704e1cfbecab3935fdebec66c72905abaad0e6e5cef6"}, + "poison": {:hex, :poison, "3.1.0", "d9eb636610e096f86f25d9a46f35a9facac35609a7591b3be3326e99a0484665", [:mix], [], "hexpm", "fec8660eb7733ee4117b85f55799fd3833eb769a6df71ccf8903e8dc5447cfce"}, + "rustler": {:hex, :rustler, "0.22.0", "e2930f9d6933e910f87526bb0a7f904e32b62a7e838a3ca4a884ee7fdfb957ed", [:mix], [{:toml, "~> 0.5.2", [hex: :toml, repo: "hexpm", optional: false]}], "hexpm", "01f5989dd511ebec09be481e07d3c59773d5373c5061e09d3ebc3ef61811b49d"}, + "toml": {:hex, :toml, "0.5.2", "e471388a8726d1ce51a6b32f864b8228a1eb8edc907a0edf2bb50eab9321b526", [:mix], [], "hexpm", "f1e3dabef71fb510d015fad18c0e05e7c57281001141504c6b69d94e99750a07"}, "uuid": {:hex, :uuid, "1.1.8", "e22fc04499de0de3ed1116b770c7737779f226ceefa0badb3592e64d5cfb4eb9", [:mix], [], "hexpm"}, } diff --git a/native/exshape_shape/.gitignore b/native/exshape_shape/.gitignore new file mode 100644 index 0000000..b83d222 --- /dev/null +++ b/native/exshape_shape/.gitignore @@ -0,0 +1 @@ +/target/ diff --git a/native/exshape_shape/Cargo.lock b/native/exshape_shape/Cargo.lock new file mode 100644 index 0000000..977025d --- /dev/null +++ b/native/exshape_shape/Cargo.lock @@ -0,0 +1,158 @@ +# This file is automatically @generated by Cargo. +# It is not intended for manual editing. +version = 3 + +[[package]] +name = "derivative" +version = "2.1.3" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "eaed5874effa6cde088c644ddcdcb4ffd1511391c5be4fdd7a5ccd02c7e4a183" +dependencies = [ + "proc-macro2", + "quote", + "syn", +] + +[[package]] +name = "either" +version = "1.6.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e78d4f1cc4ae33bbfc157ed5d5a5ef3bc29227303d595861deb238fcec4e9457" + +[[package]] +name = "exshape_shape" +version = "0.1.0" +dependencies = [ + "derivative", + "float_extras", + "itertools", + "rustler", + "rustler_codegen", +] + +[[package]] +name = "float_extras" +version = "0.1.6" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b22b70f8649ea2315955f1a36d964b0e4da482dfaa5f0d04df0d1fb7c338ab7a" +dependencies = [ + "libc", +] + +[[package]] +name = "heck" +version = "0.3.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "87cbf45460356b7deeb5e3415b5563308c0a9b057c85e12b06ad551f98d0a6ac" +dependencies = [ + "unicode-segmentation", +] + +[[package]] +name = "itertools" +version = "0.10.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "37d572918e350e82412fe766d24b15e6682fb2ed2bbe018280caa810397cb319" +dependencies = [ + "either", +] + +[[package]] +name = "lazy_static" +version = "1.4.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e2abad23fbc42b3700f2f279844dc832adb2b2eb069b2df918f455c4e18cc646" + +[[package]] +name = "libc" +version = "0.2.81" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "1482821306169ec4d07f6aca392a4681f66c75c9918aa49641a2595db64053cb" + +[[package]] +name = "proc-macro2" +version = "1.0.24" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "1e0704ee1a7e00d7bb417d0770ea303c1bccbabf0ef1667dae92b5967f5f8a71" +dependencies = [ + "unicode-xid", +] + +[[package]] +name = "quote" +version = "1.0.8" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "991431c3519a3f36861882da93630ce66b52918dcf1b8e2fd66b397fc96f28df" +dependencies = [ + "proc-macro2", +] + +[[package]] +name = "rustler" +version = "0.22.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b787d3b2a80007f41cd4c0c310cdeb3936192768159585f65ecc7e96faf97fc3" +dependencies = [ + "lazy_static", + "rustler_codegen", + "rustler_sys", +] + +[[package]] +name = "rustler_codegen" +version = "0.22.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b5a1f867002b6f0130f47abf215cac4405646db6f5d7b009b21c890980490aa4" +dependencies = [ + "heck", + "proc-macro2", + "quote", + "syn", +] + +[[package]] +name = "rustler_sys" +version = "2.1.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9fb96034ff33723615fd19223d58c987c1f6476342e83557a6e467ef95f83bda" +dependencies = [ + "unreachable", +] + +[[package]] +name = "syn" +version = "1.0.57" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "4211ce9909eb971f111059df92c45640aad50a619cf55cd76476be803c4c68e6" +dependencies = [ + "proc-macro2", + "quote", + "unicode-xid", +] + +[[package]] +name = "unicode-segmentation" +version = "1.7.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "bb0d2e7be6ae3a5fa87eed5fb451aff96f2573d2694942e40543ae0bbe19c796" + +[[package]] +name = "unicode-xid" +version = "0.2.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f7fe0bb3479651439c9112f72b6c505038574c9fbb575ed1bf3b797fa39dd564" + +[[package]] +name = "unreachable" +version = "1.0.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "382810877fe448991dfc7f0dd6e3ae5d58088fd0ea5e35189655f84e6814fa56" +dependencies = [ + "void", +] + +[[package]] +name = "void" +version = "1.0.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "6a02e4885ed3bc0f2de90ea6dd45ebcbb66dacffe03547fadbb0eeae2770887d" diff --git a/native/exshape_shape/Cargo.toml b/native/exshape_shape/Cargo.toml new file mode 100644 index 0000000..c8ab39e --- /dev/null +++ b/native/exshape_shape/Cargo.toml @@ -0,0 +1,17 @@ +[package] +name = "exshape_shape" +version = "0.1.0" +authors = [] +edition = "2018" + +[lib] +name = "exshape_shape" +path = "src/lib.rs" +crate-type = ["dylib"] + +[dependencies] +rustler = "0.22.0" +rustler_codegen = "0.22.0" +itertools = "0.10" +derivative = "2.1" +float_extras = "0.1" diff --git a/native/exshape_shape/src/lib.rs b/native/exshape_shape/src/lib.rs new file mode 100644 index 0000000..9c01ebe --- /dev/null +++ b/native/exshape_shape/src/lib.rs @@ -0,0 +1,80 @@ +use rustler::{self, Encoder, Env, Term}; +use itertools::{Itertools, Either}; + +mod atoms { + use ::rustler; + pub use rustler::types::atom::*; + rustler::atoms! { x, y } +} + +rustler::init!("Elixir.Exshape.Shp", [native_nest_polygon_impl]); + +mod point; +mod lineseg; +mod ring; +mod poly; + +use ring::Ring; +use poly::Poly; +use point::Point; + +struct Yes(T); +impl Encoder for Yes { + fn encode<'a>(&self, env: Env<'a>) -> Term<'a> { + (atoms::ok(), &self.0).encode(env) + } +} + +#[rustler::nif(schedule = "DirtyCpu")] +fn native_nest_polygon_impl<'a>(rings: Vec>) -> Yes>> { + let (mut polys, holes) = rings.into_iter().partition_map(|ring| { + if ring.is_clockwise() { + Either::Left(ring.into()) + } else { + Either::Right(ring) + } + }); + nest_holes(&mut polys, holes); + + Yes(polys) +} + +fn nest_holes<'a>(polys: &mut Vec>, holes: Vec>) { + if holes.len() == 1 { + // if there's only a single hole, we won't bother slicing the + // polygons, since we'd just throw away all that work anyway. + let hole = holes.into_iter().next().unwrap(); + process(polys, hole, Ring::contains_unsliced) + } else { + for hole in holes { + process(polys, hole, Ring::contains); + } + } +} + +fn process<'a>(polys: &mut Vec>, hole: Ring<'a>, contain: fn(&Ring<'a>, &Point) -> bool) { + match polys.len() { + 0 => { + polys.push(hole.into()); + } + 1 => { + polys[0].push(hole); + } + _ => { + // in the original, this is recursive, but we'll do it + // iteratively. What we want to do is find the first poly + // which contains the first point of the ring and push the + // hole onto it. If it doesn't fit in any poly, just smash + // it onlo the last. + let pt = hole.first_point(); + match polys.iter_mut().find(|poly| contain(poly.first_ring(), pt)) { + Some(poly) => { + poly.push(hole); + } + None => { + polys.last_mut().unwrap().push(hole); + } + } + } + } +} diff --git a/native/exshape_shape/src/lineseg.rs b/native/exshape_shape/src/lineseg.rs new file mode 100644 index 0000000..2de56cd --- /dev/null +++ b/native/exshape_shape/src/lineseg.rs @@ -0,0 +1,7 @@ +use crate::point::Point; + +#[derive(Clone, Copy, Debug)] +pub struct LineSeg { + pub a: Point, + pub b: Point +} diff --git a/native/exshape_shape/src/point.rs b/native/exshape_shape/src/point.rs new file mode 100644 index 0000000..986f62a --- /dev/null +++ b/native/exshape_shape/src/point.rs @@ -0,0 +1,5 @@ +#[derive(Clone, Copy, Debug)] +pub struct Point { + pub x: f64, + pub y: f64 +} diff --git a/native/exshape_shape/src/poly.rs b/native/exshape_shape/src/poly.rs new file mode 100644 index 0000000..a6309ef --- /dev/null +++ b/native/exshape_shape/src/poly.rs @@ -0,0 +1,39 @@ +use rustler::{Decoder, Encoder, Error, NifResult, Term}; + +use crate::ring::Ring; + +pub struct Poly<'a> { + rings: Vec> +} + +impl <'a> Poly<'a> { + pub fn first_ring(&self) -> &Ring<'a> { + &self.rings[0] + } + + pub fn push(&mut self, ring: Ring<'a>) { + self.rings.push(ring) + } +} + +impl <'a> From> for Poly<'a> { + fn from(ring: Ring<'a>) -> Self { + Self { rings: vec![ring] } + } +} + +impl <'a> Decoder<'a> for Poly<'a> { + fn decode(term: Term<'a>) -> NifResult { + let rings = term.decode::>()?; + if rings.is_empty() { + return Err(Error::BadArg); + } + Ok(Poly { rings }) + } +} + +impl <'a> Encoder for Poly<'a> { + fn encode<'b>(&self, env: rustler::Env<'b>) -> Term<'b> { + self.rings.encode(env) + } +} diff --git a/native/exshape_shape/src/ring.rs b/native/exshape_shape/src/ring.rs new file mode 100644 index 0000000..1313d1a --- /dev/null +++ b/native/exshape_shape/src/ring.rs @@ -0,0 +1,253 @@ +use itertools::Itertools; +use std::cell::{Ref, RefCell}; +use derivative::Derivative; + +use rustler::{Decoder, Encoder, Env, ListIterator, NifResult, Term, Error}; + +use crate::point::Point; +use crate::lineseg::LineSeg; +use crate::atoms; + +#[derive(Debug)] +struct Slices { + segments: Vec>, + y_min: f64, + y_max: f64 +} + +#[derive(Derivative)] +#[derivative(Debug)] +pub struct Ring<'a> { + #[derivative(Debug = "ignore")] + term: Term<'a>, + points: Vec, + slices: RefCell> +} + +// The terms that we get from shapefiles will only ever be floats, but +// that's annoying to write by hand in test code, so we'll accept +// integers too. +fn decode_floatish<'a>(term: Term<'a>) -> NifResult { + term.decode().or_else(|e| { + match e { + Error::BadArg => term.decode::().map(|i| i as f64), + other => Err(other) + } + }) +} + +impl <'a> Decoder<'a> for Ring<'a> { + fn decode(term: Term<'a>) -> NifResult { + // could define a Decoder for Point and just use Vec's Decoder + // impl, but this way we can look up the atoms just once per + // ring instead of once per point... + + let env = term.get_env(); + let x = atoms::x().to_term(env); + let y = atoms::y().to_term(env); + let points = + term.decode::>()?.map(|pt| { + Ok(Point { x : decode_floatish(pt.map_get(x)?)?, + y : decode_floatish(pt.map_get(y)?)? }) + }).collect::>>()?; + + if points.is_empty() { + return Err(Error::BadArg); + } + + Ok( + Ring { + term, + points, + slices: RefCell::new(None) + } + ) + } +} + +impl <'a> Encoder for Ring<'a> { + fn encode<'b>(&self, env: Env<'b>) -> Term<'b> { + self.term.in_env(env) + } +} + +impl <'a> Ring<'a> { + pub fn first_point(&self) -> &Point { + &self.points[0] // guaranteed to exist because the decoder requires non-emptiness + } + + pub fn is_clockwise(&self) -> bool { + let (_, area) = + self.points[1..].iter().fold((&self.points[0], 0.0), |(prev_pt, s), pt| { + (pt, s + (pt.x - prev_pt.x) * (pt.y + prev_pt.y)) + }); + area >= 0.0 + } + + fn slices(&self) -> Ref { + let mut slices = self.slices.borrow(); + if slices.is_none() { + drop(slices); + *self.slices.borrow_mut() = Some(slice(&self.points)); + slices = self.slices.borrow(); + } + Ref::map(slices, |t| t.as_ref().unwrap()) + } + + fn slice_for(&self, pt: &Point) -> Option>> { + let slices = self.slices(); + if pt.y < slices.y_min || pt.y > slices.y_max { + None + } else { + Some(Ref::map(slices, |slices| { + &slices.segments[band_for(slices.y_min, slices.y_max, pt.y, slices.segments.len())] + })) + } + } + + pub fn contains(&self, pt: &Point) -> bool { + match self.slice_for(pt) { + None => { + false + } + Some(vec) => { + let Point { x, y } = *pt; + vec.iter().fold(false, move |c, lineseg| { + if ((lineseg.a.y > y) != (lineseg.b.y > y)) && (x < ((((lineseg.b.x - lineseg.a.x) * (y - lineseg.a.y)) / (lineseg.b.y - lineseg.a.y)) + lineseg.a.x)) { + !c + } else { + c + } + }) + } + } + } + + pub fn contains_unsliced(&self, pt: &Point) -> bool { + let Point { x, y } = *pt; + self.points.iter().fold((false, self.points.last().unwrap()), move |(c, j), i| { + let c = + if ((i.y > y) != (j.y > y)) && (x < ((((j.x - i.x) * (y - i.y)) / (j.y - i.y)) + i.x)) { + !c + } else { + c + }; + (c, i) + }).0 + } +} + +fn band_for(y_min: f64, y_max: f64, y: f64, bands: usize) -> usize { + let range = y_max - y_min; + let frac = (y - y_min) / range; + (bands as f64 * frac).floor() as usize +} + +fn slice(points: &Vec) -> Slices { + let (y_min, y_max) = + points.iter().fold((f64::INFINITY, f64::NEG_INFINITY), |(min, max), pt| { + (min.min(pt.y), max.max(pt.y)) + }); + + // poke out the range a tiny bit to avoid edge cases + let y_min = float_extras::f64::nextafter(y_min, f64::NEG_INFINITY); + let y_max = float_extras::f64::nextafter(y_max, f64::INFINITY); + + let mut segments = vec![Vec::new(); 10]; + + for (&a, &b) in points.last().into_iter().chain(points.into_iter()).tuple_windows() { + let a_seg = band_for(y_min, y_max, a.y.min(y_max).max(y_min), segments.len()); + let b_seg = band_for(y_min, y_max, b.y.min(y_max).max(y_min), segments.len()); + + let min_seg = a_seg.min(b_seg); + let max_seg = a_seg.max(b_seg); + + for seg in min_seg..=max_seg { + segments[seg].push(LineSeg { a, b }); + } + } + + Slices { + y_min, + y_max, + segments + } +} + +#[cfg(test)] +mod test { + use std::cell::RefCell; + use rustler::{Term, Env}; + use super::Ring; + use crate::point::Point; + + fn fake_term() -> Term<'static> { + // SAFETY: this in fact isn't safe :) + // But trying to do anything with the term will crash anyway, + // because the nif dynamic library won't be loaded while + // running the tests. + unsafe { + Term::new(Env::new(&(), std::ptr::null_mut()), 0) + } + } + + fn unit_square() -> Ring<'static> { + Ring { + term: fake_term(), + points: vec![Point { x: -0.5, y: -0.5 }, + Point { x: -0.5, y: 0.5 }, + Point { x: 0.5, y: 0.5 }, + Point { x: 0.5, y: -0.5 }], + slices: RefCell::new(None) + } + } + + fn u_shape() -> Ring<'static> { + Ring { + term: fake_term(), + points: vec![Point { x: -0.5, y: -0.5 }, + Point { x: -0.5, y: 0.5 }, + Point { x: -0.4, y: 0.5 }, + Point { x: -0.4, y: -0.4 }, + Point { x: 0.4, y: -0.4}, + Point { x: 0.4, y: 0.5 }, + Point { x: 0.5, y: 0.5 }, + Point { x: 0.5, y: -0.5 }], + slices: RefCell::new(None) + } + } + + fn contains(ring: &Ring, pt: Point) { + assert!(ring.contains(&pt)); + assert!(ring.contains_unsliced(&pt)); + } + + fn doesnt_contain(ring: &Ring, pt: Point) { + assert!(!ring.contains(&pt)); + assert!(!ring.contains_unsliced(&pt)); + } + + #[test] + fn basic_sanity_check() { + let sq = unit_square(); + contains(&sq, Point { x: 0.0, y: 0.0 }); + doesnt_contain(&sq, Point { x: 10.0, y: 0.0 }); + doesnt_contain(&sq, Point { x: -10.0, y: 0.0 }); + doesnt_contain(&sq, Point { x: 0.0, y: 10.0 }); + doesnt_contain(&sq, Point { x: 0.0, y: -10.0 }); + } + + #[test] + fn more_complex_sanity_check() { + let shape = u_shape(); + contains(&shape, Point { x: -0.45, y: 0.0 }); + contains(&shape, Point { x: 0.45, y: 0.0 }); + contains(&shape, Point { x: 0.0, y: -0.45 }); + doesnt_contain(&shape, Point { x: 0.0, y: 0.45 }); + doesnt_contain(&shape, Point { x: 0.0, y: 0.0 }); + doesnt_contain(&shape, Point { x: 10.0, y: 0.0 }); + doesnt_contain(&shape, Point { x: -10.0, y: 0.0 }); + doesnt_contain(&shape, Point { x: 0.0, y: 10.0 }); + doesnt_contain(&shape, Point { x: 0.0, y: -10.0 }); + } +} diff --git a/test/exshape_test.exs b/test/exshape_test.exs index 4d1af9f..3491214 100644 --- a/test/exshape_test.exs +++ b/test/exshape_test.exs @@ -115,4 +115,28 @@ defmodule ExshapeTest do [_, {shape, _}] = Enum.into(stream, []) assert (shape.points |> List.flatten |> length) == 174045 end + + test "native and beam produce the same results" do + ["archive", "chicago_zoning", "co-parcels", "hoods", "howard-beach", "row_181", "seattle_basketball_points", "speed_enforcement", "zillow"] + |> Enum.each(fn path -> + [{_, _, beam_stream}] = zip(path) |> Exshape.from_zip(native: false) + beam = Enum.into(beam_stream, []) + [{_, _, native_stream}] = zip(path) |> Exshape.from_zip(native: true) + native = Enum.into(native_stream, []) + assert beam == native + end) + + ["multipatch", "multipointm", "multipoint", "multipointz", "pointm", "point", "pointz", "polygonm", "polygon", "polygons", "polygonz", "polylinem", "polyline", "polylinez", "Neighborhoods/neighborhoods_orleans"] + |> Enum.each(fn path -> + beam = shp(path) + |> File.stream!([], 2048) + |> Exshape.Shp.read(native: false) + |> Enum.into([]) + native = shp(path) + |> File.stream!([], 2048) + |> Exshape.Shp.read(native: true) + |> Enum.into([]) + assert beam == native + end) + end end diff --git a/test/shp_test.exs b/test/shp_test.exs index 7a5de4a..1735711 100644 --- a/test/shp_test.exs +++ b/test/shp_test.exs @@ -10,6 +10,13 @@ defmodule ShpTest do } doctest Exshape + def nest_polygon(p) do + beam = Shp.beam_nest_polygon(p) + native = Shp.native_nest_polygon(p) + assert beam == native + beam + end + describe "regular geoms" do test "can read points" do [_header | points] = fixture("point.shp") @@ -304,7 +311,7 @@ defmodule ShpTest do describe "nesting" do test "can nest holes" do - assert Shp.nest_polygon(%Polygon{ + assert nest_polygon(%Polygon{ parts: [0, 5], points: Enum.reverse([ %Point{x: 0, y: 4}, @@ -340,7 +347,7 @@ defmodule ShpTest do end test "appends a part to the polygon when the part is clockwise" do - assert Shp.nest_polygon(%Polygon{ + assert nest_polygon(%Polygon{ parts: [0, 5], points: Enum.reverse([ %Point{x: 0, y: 4}, @@ -422,7 +429,7 @@ defmodule ShpTest do end test "can nest many holes" do - assert Shp.nest_polygon(%Polygon{ + assert nest_polygon(%Polygon{ parts: [0, 5, 10], points: Enum.reverse([ %Point{x: 0, y: 5}, @@ -471,7 +478,7 @@ defmodule ShpTest do end test "can nest holes and rings" do - assert Shp.nest_polygon(%Polygon{ + assert nest_polygon(%Polygon{ parts: [0, 5, 10], points: Enum.reverse([ %Point{x: 0, y: 5}, diff --git a/test/test_helper.exs b/test/test_helper.exs index 9b1256e..cf59491 100644 --- a/test/test_helper.exs +++ b/test/test_helper.exs @@ -4,5 +4,6 @@ defmodule TestHelper do end def zip(name), do: "#{__DIR__}/fixtures/#{name}.zip" + def shp(name), do: "#{__DIR__}/fixtures/#{name}.shp" end ExUnit.start()