mod cartesian; mod contains; mod clip; mod clip_antimeridian; mod path_sink; mod resample; pub use self::cartesian::Cartesian; pub use self::contains::polygon_contains_south; pub use self::clip::{Clip, ClipControl, Clipper, PathClipper}; pub use self::clip_antimeridian::ClipAntimeridian; pub use self::path_sink::PathSink; pub use self::resample::ResamplePath; use super::*; const F32_PRECISION: f32 = 1e-6; use crate::api::world_geo::{Feature, FeatureCollection, Geometry, PolygonData, Position}; #[derive(Clone, Copy, Debug)] pub struct RadianPoint { pub lambda: f32, // "longitude" pub phi: f32, // "latitude" } impl RadianPoint { pub fn wrap(self) -> Self { use std::f32::consts::{PI, FRAC_PI_2}; let lambda = if self.lambda.abs() <= PI { self.lambda } else { self.lambda.signum() * ((self.lambda.abs() + PI) % (2.0 * PI) - PI) }; let phi = if self.phi.abs() <= FRAC_PI_2 { self.phi } else { self.phi.signum() * ((self.phi.abs() + FRAC_PI_2) % PI - FRAC_PI_2) }; Self { lambda, phi } } } impl PartialEq for RadianPoint { fn eq(&self, other: &RadianPoint) -> bool { (self.lambda - other.lambda).abs() < F32_PRECISION && (self.phi - other.phi).abs() < F32_PRECISION } } impl Eq for RadianPoint {} impl From for RadianPoint { fn from(pos: Position) -> Self { Self { lambda: pos.longitude.to_radians(), // + std::f32::consts::FRAC_PI_4, // some test rotation phi: pos.latitude.to_radians(), }.wrap() } } #[derive(Clone, Copy, Debug)] pub struct ProjectedPoint { pub x: f32, pub y: f32, } pub trait Projection { fn project(&self, p: RadianPoint) -> ProjectedPoint; } impl Projection for &P { fn project(&self, p: RadianPoint) -> ProjectedPoint { (**self).project(p) } } pub trait DrawPath { type Point; // only valid during ring or path are "open" fn line(&mut self, to: Self::Point, stroke: bool); // closed polygon fn ring_start(&mut self) {} fn ring_end(&mut self) {} // open path, not filled fn path_start(&mut self) {} fn path_end(&mut self) {} } impl PathTransformer for &mut DP { type Point = DP::Point; type Sink = DP; fn sink(&mut self) -> &mut Self::Sink { &mut **self } fn transform_line(&mut self, to: Self::Point, stroke: bool) { self.line(to, stroke); } } pub trait PathTransformer { type Point; type Sink: DrawPath; fn sink(&mut self) -> &mut Self::Sink; fn transform_line(&mut self, to: Self::Point, stroke: bool); fn transform_ring_start(&mut self) { self.sink().ring_start(); } fn transform_ring_end(&mut self) { self.sink().ring_end(); } fn transform_path_start(&mut self) { self.sink().path_start(); } fn transform_path_end(&mut self) { self.sink().path_end(); } } impl DrawPath for PT { type Point = ::Point; fn line(&mut self, to: Self::Point, stroke: bool) { self.transform_line(to, stroke); } fn ring_start(&mut self) { self.transform_ring_start(); } fn ring_end(&mut self) { self.transform_ring_end(); } fn path_start(&mut self) { self.transform_path_start(); } fn path_end(&mut self) { self.transform_path_end(); } } pub struct ApplyProjection { pub projection: P, pub sink: S, } impl> ApplyProjection { #[allow(unused)] pub fn new(projection: P, sink: S) -> Self { Self { projection, sink, } } } impl> PathTransformer for ApplyProjection { type Sink = S; type Point = RadianPoint; fn sink(&mut self) -> &mut Self::Sink { &mut self.sink } fn transform_line(&mut self, to: Self::Point, stroke: bool) { self.sink.line(self.projection.project(to), stroke); } } #[allow(unused)] #[derive(Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Debug)] pub enum MapProjection { EquiRectangular, EqualEarth, Fahey, } impl Projection for MapProjection { fn project(&self, p: RadianPoint) -> ProjectedPoint { use std::f32::consts::{PI, FRAC_PI_2}; match self { MapProjection::EquiRectangular => { ProjectedPoint { x: p.lambda.min(PI).max(-PI) * (180.0 / PI), y: -p.phi.min(FRAC_PI_2).max(-FRAC_PI_2) * (180.0 / PI), } }, MapProjection::EqualEarth => { const A1: f32 = 1.340264; const A2: f32 = -0.081106; const A3: f32 = 0.000893; const A4: f32 = 0.003796; let m = (3.0f32).sqrt() / 2.0; let l = (m * p.phi.sin()).asin(); let l2 = l * l; let l6 = l2 * l2 * l2; ProjectedPoint { x: 70.0 * p.lambda * l.cos() / (m * (A1 + 3.0 * A2 * l2 + l6 * (7.0 * A3 + 9.0 * A4 * l2))), y: -70.0 * l * (A1 + A2 * l2 + l6 * (A3 + A4 * l2)), } }, MapProjection::Fahey => { let fahey_k = 35.0f32.to_radians().cos(); let t = (p.phi / 2.0).tan(); ProjectedPoint { x: 50.0 * p.lambda * fahey_k * (1.0 - t*t).sqrt(), y: -50.0 * (1.0 + fahey_k) * t, } } } } } pub trait DrawOnMap { fn draw_on_map(&self, projection: MapProjection) -> Html; } impl DrawOnMap for PolygonData { /* // wrapping output stream: transformRadians(transformRotate(rotate)(preclip(projectResample(postclip(stream))))); -> first map input points to radian (*pi/180) -> apply rotation (zero by default) -> clip (pre) (clipAntimeridian by default) -> resample and project -> clip (post) (identity by default) */ fn draw_on_map(&self, projection: MapProjection) -> Html { let mut sink = PathSink::new(); let mut resample_sink = ClipAntimeridian::new(ResamplePath::new(projection, sink.draw_path(), 0.5)); // let mut resample_sink = ClipAntimeridian::new(ApplyProjection::new(projection, sink.draw_path())); for line in &self.0 { let line = line.iter().cloned().map(RadianPoint::from); resample_sink.ring_start(); for point in line { resample_sink.line(point, true); } resample_sink.ring_end(); } drop(resample_sink); if let Some(stroke_path) = sink.stroke_path { html!{ } } else { html!{ } } } } impl DrawOnMap for Geometry { fn draw_on_map(&self, projection: MapProjection) -> Html { match self { Geometry::Polygon { coordinates } => coordinates.draw_on_map(projection), Geometry::MultiPolygon { coordinates } => html!{ { for coordinates.iter().map(|p| p.draw_on_map(projection)) } }, } } } impl DrawOnMap for Feature { fn draw_on_map(&self, projection: MapProjection) -> Html { html!{ { self.geometry.draw_on_map(projection) } } } } impl DrawOnMap for FeatureCollection { fn draw_on_map(&self, projection: MapProjection) -> Html { html!{ { for self.features.iter().map(|c| c.draw_on_map(projection)) } } } } fn build_grid(projection: MapProjection) -> Html
{ let mut sink = PathSink::new(); let mut resample_sink = ResamplePath::new(projection, sink.draw_path(), 0.5); for longitude in (-180..=180).step_by(10) { let height = if longitude % 90 == 0 { 90 } else { 80 }; let lambda = (longitude as f32).to_radians(); let phi = (height as f32).to_radians(); resample_sink.path_start(); resample_sink.line(RadianPoint { lambda, phi: -phi }, true); resample_sink.line(RadianPoint { lambda, phi: 0.0 }, true); resample_sink.line(RadianPoint { lambda, phi: phi }, true); resample_sink.path_end(); } for latitude in (-90..=90).step_by(10) { let phi = (latitude as f32).to_radians(); resample_sink.path_start(); for longitude in (-180..=180).step_by(5) { resample_sink.line(RadianPoint { lambda: (longitude as f32).to_radians(), phi }, true); } resample_sink.path_end(); } drop(resample_sink); html!{ } } impl Main { pub fn view_map(&self, world_geo: &FeatureCollection) -> Html { html!{ { world_geo.draw_on_map(MapProjection::Fahey) } { build_grid(MapProjection::Fahey) } } } }