UTM projection

1 files changed, 26 insertions, 60 deletions

diff --git a/server/src/api.cpp b/server/src/api.cpp
index 5dd3fc5..137666b 100644
--- a/server/src/api.cpp
+++ b/server/src/api.cpp
@@ -17,12 +17,14 @@ namespace views = std::views;
 
 namespace routemon::api {
 
-auto json_value_from_point(geo::point const& p) -> json::value
+auto json_value_from_point(geo::utm::zonable_wgs84_point const& p)
+    -> json::value
 {
-  return json::array{bgeo::get<1>(p), bgeo::get<0>(p)};
+  return json::array{p.lat(), p.lon()};
 }
 
-auto json_value_from_linestring(geo::linestring const& ls) -> json::value
+auto json_value_from_linestring(geo::utm::zonable_wgs84_linestring const& ls)
+    -> json::value
 {
   json::array a;
   for (auto const& p : ls)
@@ -30,8 +32,8 @@ auto json_value_from_linestring(geo::linestring const& ls) -> json::value
   return a;
 }
 
-auto json_value_from_linestrings(std::vector<geo::linestring> const& lss)
-    -> json::value
+auto json_value_from_linestrings(
+    std::vector<geo::utm::zonable_wgs84_linestring> const& lss) -> json::value
 {
   json::array a;
   for (auto const& ls : lss)
@@ -108,24 +110,7 @@ handler::handler(log::logger const& l, datex2::situation_publication pub)
         static_cast<float>(i) / static_cast<float>(nsit) * 100.0f);
   };
 
-  l_.info("Building indices for {} situations", nsit);
-
-  auto const buffer_distance = min_distance_;
-  auto const points_per_circle = 8;
-  // Note: thomas strategy does not work for geographic_join_round;
-  // need to use andoyer for that.
-  using formula = bgeo::strategy::andoyer;
-  bgeo::strategy::buffer::distance_symmetric<double> distance_strategy{
-    buffer_distance
-  };
-  bgeo::strategy::buffer::geographic_join_round<formula> join_strategy{
-    points_per_circle
-  };
-  bgeo::strategy::buffer::geographic_end_round<formula> end_strategy{4};
-  bgeo::strategy::buffer::geographic_point_circle<formula> circle_strategy{
-    points_per_circle
-  };
-  bgeo::strategy::buffer::geographic_side_straight<formula> side_strategy;
+  l_.info("Building indexes for {} situations", nsit);
 
   l_.info("Progress: 0%");
   auto const before_build = chrono::steady_clock::now();
@@ -135,17 +120,11 @@ handler::handler(log::logger const& l, datex2::situation_publication pub)
     {
       for (auto const& ls : rc->relevant_line_strings)
       {
-        auto buffered_ls = geo::multi_polygon{};
-        bgeo::buffer(
-            *ls, buffered_ls, distance_strategy, side_strategy, join_strategy,
-            end_strategy, circle_strategy);
-        auto blse = geo::box{};
-        bgeo::envelope(buffered_ls, blse);
-        blse_index_.insert(std::make_tuple(blse, buffered_ls, rc));
+        blse_index_.insert(*ls, rc);
       }
       for (auto p : rc->relevant_points)
       {
-        p_index_.insert(std::make_pair(p, rc));
+        p_index_.insert(p_index_value{p, rc});
       }
     }
 
@@ -160,7 +139,7 @@ handler::handler(log::logger const& l, datex2::situation_publication pub)
 
   auto const after_build = chrono::steady_clock::now();
   auto const dur_build_s = chrono::duration<double>{after_build - before_build};
-  l_.info("Indices built in {}", dur_build_s);
+  l_.info("Indexes built in {}", dur_build_s);
   l_.info("BLSE index size: {}", blse_index_.size());
   l_.info("Point index size: {}", p_index_.size());
 }
@@ -175,16 +154,13 @@ auto handler::process_gpx(gpx::file&& gpx_file)
   auto const check_periods = time::period_seq{relevant.begin(), relevant.end()};
 
   // TODO: eliminate use of overlap segments
-  auto splits_with_overlap_segments = std::vector<geo::linestring>{};
+  auto segments = std::vector<geo::utm::zonable_wgs84_linestring>{};
   for (auto const& track : gpx_file.tracks)
     for (auto const& seg : track.segments)
-      geo::split_linestring_with_overlap_segments(
-          seg.waypoints,
-          5000 /* meters max total dist until a new split is forced */,
-          splits_with_overlap_segments);
+      segments.push_back(seg.waypoints);
   auto const before_query = chrono::steady_clock::now();
 
-  auto vincenty_strategy = geo::vincenty_strategy{};
+  auto vincenty_strategy = geo::wgs84::vincenty_strategy{};
 
   l_.debug("Querying for relevant situations");
   auto relevant_road_closures =
@@ -192,28 +168,17 @@ auto handler::process_gpx(gpx::file&& gpx_file)
   auto ls_checked = 0uz;
   auto p_checked = 0uz;
   auto i = 0;
-  for (geo::linestring const& part : splits_with_overlap_segments)
+  for (auto const& part : segments)
   {
-    l_.debug("Checking part [{}/{}]", ++i, splits_with_overlap_segments.size());
+    l_.debug("Checking part [{}/{}]", ++i, segments.size());
 
+    auto part_zone_lss = geo::multizonal::zoned_linestring_seg_seq{part};
     // TODO: consider buffering with min_distance_
-    auto part_box = geo::box{};
+    auto part_box = geo::wgs84::box{};
     bgeo::envelope(part, part_box);
 
-    for (auto it = blse_index_.qbegin(bgeo::index::intersects(part_box));
-         it != blse_index_.qend(); it++)
+    for (auto const& rc : blse_index_.intersection(part_zone_lss))
     {
-      // Cannot use structured bindings here, as boost::geometry::get
-      // interferes with ADL. It is a candidate as the namespace
-      // boost::geometry is part of the associated namespace set, which
-      // happens because geo::linestring ≡
-      // boost::geometry::model::linestring<geo::point> is part of the
-      // whole tuple type (blse_index_value) that is the value_type of the
-      // iterator.
-
-      geo::box const& blse = std::get<0>(*it);
-      geo::multi_polygon const& bls = std::get<1>(*it);
-      std::shared_ptr<datex2::road_closure> const& rc = std::get<2>(*it);
       if (rc->validity
           && rc->validity->intersect(check_periods).periods().empty())
         continue;
@@ -221,8 +186,7 @@ auto handler::process_gpx(gpx::file&& gpx_file)
       // auto small_parts = bgeo::model::multi_linestring<geo::linestring>{};
       // bgeo::intersection(part, lse, small_parts);
 
-      if (bgeo::intersects(blse, part) && bgeo::intersects(bls, part))
-        relevant_road_closures.emplace(rc);
+      relevant_road_closures.emplace(rc);
 
       // if (!bgeo::is_empty(small_parts))
       //   if (bgeo::distance(*ls, small_parts, vincenty_strategy) <
@@ -234,7 +198,7 @@ auto handler::process_gpx(gpx::file&& gpx_file)
          it != p_index_.qend(); it++)
     {
       // Cannot use structured bindings here for the same reason as above.
-      geo::point const& p = std::get<0>(*it);
+      geo::utm::zonable_wgs84_point const& p = std::get<0>(*it);
       std::shared_ptr<datex2::road_closure> const& rc = std::get<1>(*it);
       if (rc->validity
           && rc->validity->intersect(check_periods).periods().empty())
@@ -301,10 +265,12 @@ auto handler::process_gpx(gpx::file&& gpx_file)
                             .relevant_lss =
                                 rc->relevant_line_strings
                                 | views::transform(
-                                    [](auto const& lsp) -> geo::linestring
+                                    [](auto const& lsp)
+                                        -> geo::utm::zonable_wgs84_linestring
                                     { return *lsp; })
-                                | std::ranges::
-                                    to<std::vector<geo::linestring>>(),
+                                | std::ranges::to<std::vector<
+                                    geo::utm::zonable_wgs84_linestring
+                                >>(),
                           };
                         })
                     | std::ranges::to<std::vector<relevant_road_closure>>(),