From ec0134cbf4ac1379f82fd66c71871ff89a1c9647 Mon Sep 17 00:00:00 2001 From: Jaiveer <32986920+CallMeJai@users.noreply.github.com> Date: Fri, 22 Dec 2023 14:20:05 -0800 Subject: [PATCH 1/6] solution to percent geodesic without formal tests --- Cargo.toml | 4 ++++ src/percent_distance.rs | 18 ++++++++++++++++++ 2 files changed, 22 insertions(+) create mode 100644 src/percent_distance.rs diff --git a/Cargo.toml b/Cargo.toml index 445a66e..42042ba 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -51,6 +51,10 @@ path = "src/geodesic.rs" name = "geodesic2" path = "src/geodesic2.rs" +[[bin]] +name = "percent_distance" +path = "src/percent_distance.rs" + [[bin]] name = "test" path = "src/test.rs" diff --git a/src/percent_distance.rs b/src/percent_distance.rs new file mode 100644 index 0000000..b69cf3c --- /dev/null +++ b/src/percent_distance.rs @@ -0,0 +1,18 @@ + +use geographiclib_rs::{Geodesic, DirectGeodesic, InverseGeodesic}; +// given two points (A, B) and a proportion P, +// find the point X along the geodesic P of the length of the geodesic from A +fn percent_geodesic(a: (f64, f64), b: (f64, f64), p: f64) -> (f64, f64) { + let geod = Geodesic::wgs84(); + let (s_ab, azi_a, _, _) = geod.inverse(a.0, a.1, b.0, b.1); + geod.direct(a.0, a.1, azi_a, p * s_ab) +} + +fn main() { + println!("24 km case:"); + println!("{:?}", percent_geodesic((52.0, 5.0), (51.4, 6.0), 0.25)); + println!("1000 km case:"); + println!("{:?}", percent_geodesic((42.0, 29.0), (39.0, -77.0), 0.5)); + println!("12200 km case:"); + println!("{:?}", percent_geodesic((42.0, 29.0), (-35.0, -70.0), 0.75)); +} \ No newline at end of file From d55f276225cfa051fd21357f48aeb6aac4d8a111 Mon Sep 17 00:00:00 2001 From: Jaiveer <32986920+CallMeJai@users.noreply.github.com> Date: Sat, 23 Dec 2023 14:23:27 -0800 Subject: [PATCH 2/6] implement percent_linestring --- Cargo.toml | 2 ++ src/percent_distance.rs | 22 ++++++++++++++++++++++ 2 files changed, 24 insertions(+) diff --git a/Cargo.toml b/Cargo.toml index b94af7b..4519ba1 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -31,6 +31,8 @@ tokio-postgres = "*" tokio-pg-mapper = "*" tokio-pg-mapper-derive = "*" actix_block_ai_crawling = "0.2.8" +geo = "0.27.0" +itertools = "0.12.0" [[bin]] name = "actix" diff --git a/src/percent_distance.rs b/src/percent_distance.rs index b69cf3c..211ef7c 100644 --- a/src/percent_distance.rs +++ b/src/percent_distance.rs @@ -1,5 +1,7 @@ use geographiclib_rs::{Geodesic, DirectGeodesic, InverseGeodesic}; +use geo::LineString; +use itertools::Itertools; // given two points (A, B) and a proportion P, // find the point X along the geodesic P of the length of the geodesic from A fn percent_geodesic(a: (f64, f64), b: (f64, f64), p: f64) -> (f64, f64) { @@ -8,6 +10,26 @@ fn percent_geodesic(a: (f64, f64), b: (f64, f64), p: f64) -> (f64, f64) { geod.direct(a.0, a.1, azi_a, p * s_ab) } +fn percent_linestring(l: LineString, p: f64) -> (f64, f64) { + let geod = Geodesic::wgs84(); + let mut distances: Vec = Vec::new(); + for (a, b) in l.0.iter().tuple_windows() { + distances.push(geod.inverse(a.x, a.y, b.x, b.y)); + } + let distance: f64 = p * distances.iter().sum::(); + let mut sum = 0.0; + let mut i = 0; + for _ in 0..distances.len() { + sum += distances[i]; + if sum > distance { + break; + } + i += 1; + } + let (azi_a, _, _) = geod.inverse(l[i].x, l[i].y, l[i + 1].x, l[i + 1].y); + geod.direct(l[i].x, l[i].y, azi_a, distance - sum + distances[i]) +} + fn main() { println!("24 km case:"); println!("{:?}", percent_geodesic((52.0, 5.0), (51.4, 6.0), 0.25)); From 445716589ac7bc021e5195e26d2f387cd04e527d Mon Sep 17 00:00:00 2001 From: Jaiveer <32986920+CallMeJai@users.noreply.github.com> Date: Mon, 8 Jan 2024 10:38:52 -0800 Subject: [PATCH 3/6] Return geo::Coord instead of (f64, f64) Makes testing easier since geo::Coord has default implementations such as relative_eq. If the original tuple is needed, Coord implements x_y() which gives the original tuple. --- src/percent_distance.rs | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/src/percent_distance.rs b/src/percent_distance.rs index 211ef7c..c9a6db1 100644 --- a/src/percent_distance.rs +++ b/src/percent_distance.rs @@ -1,16 +1,17 @@ use geographiclib_rs::{Geodesic, DirectGeodesic, InverseGeodesic}; -use geo::LineString; +use geo::{coord, Coord, LineString}; use itertools::Itertools; // given two points (A, B) and a proportion P, // find the point X along the geodesic P of the length of the geodesic from A -fn percent_geodesic(a: (f64, f64), b: (f64, f64), p: f64) -> (f64, f64) { +fn percent_geodesic(a: (f64, f64), b: (f64, f64), p: f64) -> Coord { let geod = Geodesic::wgs84(); let (s_ab, azi_a, _, _) = geod.inverse(a.0, a.1, b.0, b.1); - geod.direct(a.0, a.1, azi_a, p * s_ab) + let (x, y) = geod.direct(a.0, a.1, azi_a, p * s_ab); + coord! {x: x, y: y,} } -fn percent_linestring(l: LineString, p: f64) -> (f64, f64) { +fn percent_linestring(l: LineString, p: f64) -> Coord { let geod = Geodesic::wgs84(); let mut distances: Vec = Vec::new(); for (a, b) in l.0.iter().tuple_windows() { @@ -27,7 +28,8 @@ fn percent_linestring(l: LineString, p: f64) -> (f64, f64) { i += 1; } let (azi_a, _, _) = geod.inverse(l[i].x, l[i].y, l[i + 1].x, l[i + 1].y); - geod.direct(l[i].x, l[i].y, azi_a, distance - sum + distances[i]) + let (x, y) = geod.direct(l[i].x, l[i].y, azi_a, distance - sum + distances[i]); + coord! {x: x, y: y,} } fn main() { From a7fbc80eff5583e6a0f5417a5cefd581a39af9f8 Mon Sep 17 00:00:00 2001 From: Jaiveer <32986920+CallMeJai@users.noreply.github.com> Date: Mon, 8 Jan 2024 12:20:28 -0800 Subject: [PATCH 4/6] Write basic unit tests for percent_geodesic --- src/percent_distance.rs | 32 ++++++++++++++++++++++++++++++++ 1 file changed, 32 insertions(+) diff --git a/src/percent_distance.rs b/src/percent_distance.rs index c9a6db1..3fc5343 100644 --- a/src/percent_distance.rs +++ b/src/percent_distance.rs @@ -39,4 +39,36 @@ fn main() { println!("{:?}", percent_geodesic((42.0, 29.0), (39.0, -77.0), 0.5)); println!("12200 km case:"); println!("{:?}", percent_geodesic((42.0, 29.0), (-35.0, -70.0), 0.75)); +} + +#[cfg(test)] +mod tests { + use super::percent_geodesic; + use approx::assert_relative_eq; + #[test] + fn test_percent_geodesic_from_both_ends_short() { + let a = (52.0, 5.0); + let b = (51.4, 6.0); + assert_relative_eq!(percent_geodesic(a, b, 0.5), percent_geodesic(b, a, 0.5), epsilon=1e-8); + assert_relative_eq!(percent_geodesic(a, b, 0.75), percent_geodesic(b, a, 0.25), epsilon=1e-8); + assert_relative_eq!(percent_geodesic(a, b, 0.125), percent_geodesic(b, a, 0.875), epsilon=1e-8); + } + + #[test] + fn test_percent_geodesic_from_both_ends_long() { + let a = (42.0, 29.0); + let b = (39.0, -77.0); + assert_relative_eq!(percent_geodesic(a, b, 0.5), percent_geodesic(b, a, 0.5), epsilon=1e-8); + assert_relative_eq!(percent_geodesic(a, b, 0.75), percent_geodesic(b, a, 0.25), epsilon=1e-8); + assert_relative_eq!(percent_geodesic(a, b, 0.125), percent_geodesic(b, a, 0.875), epsilon=1e-8); + } + + #[test] + fn test_percent_geodesic_from_both_ends_very_long() { + let a = (42.0, 29.0); + let b = (-35.0, -70.0); + assert_relative_eq!(percent_geodesic(a, b, 0.5), percent_geodesic(b, a, 0.5), epsilon=1e-8); + assert_relative_eq!(percent_geodesic(a, b, 0.75), percent_geodesic(b, a, 0.25), epsilon=1e-8); + assert_relative_eq!(percent_geodesic(a, b, 0.125), percent_geodesic(b, a, 0.875), epsilon=1e-8); + } } \ No newline at end of file From 05ebf1804abb09fa52e945e7af4a6183c10132f6 Mon Sep 17 00:00:00 2001 From: Jaiveer <32986920+CallMeJai@users.noreply.github.com> Date: Mon, 8 Jan 2024 12:22:17 -0800 Subject: [PATCH 5/6] Make library functions public --- src/percent_distance.rs | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/percent_distance.rs b/src/percent_distance.rs index 3fc5343..2d9120b 100644 --- a/src/percent_distance.rs +++ b/src/percent_distance.rs @@ -4,14 +4,14 @@ use geo::{coord, Coord, LineString}; use itertools::Itertools; // given two points (A, B) and a proportion P, // find the point X along the geodesic P of the length of the geodesic from A -fn percent_geodesic(a: (f64, f64), b: (f64, f64), p: f64) -> Coord { +pub fn percent_geodesic(a: (f64, f64), b: (f64, f64), p: f64) -> Coord { let geod = Geodesic::wgs84(); let (s_ab, azi_a, _, _) = geod.inverse(a.0, a.1, b.0, b.1); let (x, y) = geod.direct(a.0, a.1, azi_a, p * s_ab); coord! {x: x, y: y,} } -fn percent_linestring(l: LineString, p: f64) -> Coord { +pub fn percent_linestring(l: LineString, p: f64) -> Coord { let geod = Geodesic::wgs84(); let mut distances: Vec = Vec::new(); for (a, b) in l.0.iter().tuple_windows() { From 297dee114708153ae4db0e9ecd5d5f4f4abcc68e Mon Sep 17 00:00:00 2001 From: Jaiveer <32986920+CallMeJai@users.noreply.github.com> Date: Mon, 8 Jan 2024 23:59:48 -0800 Subject: [PATCH 6/6] Write basic unit test for percent_linestring --- src/percent_distance.rs | 32 +++++++++++++++++++++++++++++++- 1 file changed, 31 insertions(+), 1 deletion(-) diff --git a/src/percent_distance.rs b/src/percent_distance.rs index 2d9120b..db17ff5 100644 --- a/src/percent_distance.rs +++ b/src/percent_distance.rs @@ -43,8 +43,10 @@ fn main() { #[cfg(test)] mod tests { - use super::percent_geodesic; + use super::{percent_geodesic, percent_linestring}; use approx::assert_relative_eq; + use gtfs_structures::{GtfsReader, Shape}; + use geo::{coord, LineString}; #[test] fn test_percent_geodesic_from_both_ends_short() { let a = (52.0, 5.0); @@ -71,4 +73,32 @@ mod tests { assert_relative_eq!(percent_geodesic(a, b, 0.75), percent_geodesic(b, a, 0.25), epsilon=1e-8); assert_relative_eq!(percent_geodesic(a, b, 0.125), percent_geodesic(b, a, 0.875), epsilon=1e-8); } + + fn shape_to_line_string(s: &Vec) -> LineString { + s.iter().map(|p| coord! {x: p.latitude, y: p.longitude}).collect() + } + + fn shape_to_reversed_line_string(s: &Vec) -> LineString { + s.iter().rev().map(|p| coord! {x: p.latitude, y: p.longitude}).collect() + } + + #[test] + fn test_percent_linestring_from_both_ends() { + let gtfs = GtfsReader::default() + .read_stop_times(false) + .read_shapes(true) + .unkown_enum_as_default(false) + .read("gtfs_rail").unwrap(); + for (_, shape) in gtfs.shapes { + assert_relative_eq!(percent_linestring(shape_to_line_string(&shape), 0.5), + percent_linestring(shape_to_reversed_line_string(&shape), 0.5), + epsilon=1e-8); + assert_relative_eq!(percent_linestring(shape_to_line_string(&shape), 0.75), + percent_linestring(shape_to_reversed_line_string(&shape), 0.25), + epsilon=1e-8); + assert_relative_eq!(percent_linestring(shape_to_line_string(&shape), 0.125), + percent_linestring(shape_to_reversed_line_string(&shape), 0.875), + epsilon=1e-8); + } + } } \ No newline at end of file