diff options
Diffstat (limited to 'server')
| -rw-r--r-- | server/src/api.cpp | 86 | ||||
| -rw-r--r-- | server/src/api.cppm | 19 | ||||
| -rw-r--r-- | server/src/datex2.cpp | 62 | ||||
| -rw-r--r-- | server/src/datex2.cppm | 12 | ||||
| -rw-r--r-- | server/src/geo.cpp | 42 | ||||
| -rw-r--r-- | server/src/geo.cppm | 23 | ||||
| -rw-r--r-- | server/src/geo/geo.cppm | 7 | ||||
| -rw-r--r-- | server/src/geo/multizonal.cppm | 240 | ||||
| -rw-r--r-- | server/src/geo/utm.cppm | 199 | ||||
| -rw-r--r-- | server/src/geo/utm_zone_local.cppm | 74 | ||||
| -rw-r--r-- | server/src/geo/wgs84.cppm | 115 | ||||
| -rw-r--r-- | server/src/gpx.cpp | 16 | ||||
| -rw-r--r-- | server/src/gpx.cppm | 4 | ||||
| -rw-r--r-- | server/src/log.cpp | 2 |
14 files changed, 738 insertions, 163 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; | |||
| 17 | 17 | ||
| 18 | namespace routemon::api { | 18 | namespace routemon::api { |
| 19 | 19 | ||
| 20 | auto json_value_from_point(geo::point const& p) -> json::value | 20 | auto json_value_from_point(geo::utm::zonable_wgs84_point const& p) |
| 21 | -> json::value | ||
| 21 | { | 22 | { |
| 22 | return json::array{bgeo::get<1>(p), bgeo::get<0>(p)}; | 23 | return json::array{p.lat(), p.lon()}; |
| 23 | } | 24 | } |
| 24 | 25 | ||
| 25 | auto json_value_from_linestring(geo::linestring const& ls) -> json::value | 26 | auto json_value_from_linestring(geo::utm::zonable_wgs84_linestring const& ls) |
| 27 | -> json::value | ||
| 26 | { | 28 | { |
| 27 | json::array a; | 29 | json::array a; |
| 28 | for (auto const& p : ls) | 30 | for (auto const& p : ls) |
| @@ -30,8 +32,8 @@ auto json_value_from_linestring(geo::linestring const& ls) -> json::value | |||
| 30 | return a; | 32 | return a; |
| 31 | } | 33 | } |
| 32 | 34 | ||
| 33 | auto json_value_from_linestrings(std::vector<geo::linestring> const& lss) | 35 | auto json_value_from_linestrings( |
| 34 | -> json::value | 36 | std::vector<geo::utm::zonable_wgs84_linestring> const& lss) -> json::value |
| 35 | { | 37 | { |
| 36 | json::array a; | 38 | json::array a; |
| 37 | for (auto const& ls : lss) | 39 | for (auto const& ls : lss) |
| @@ -108,24 +110,7 @@ handler::handler(log::logger const& l, datex2::situation_publication pub) | |||
| 108 | static_cast<float>(i) / static_cast<float>(nsit) * 100.0f); | 110 | static_cast<float>(i) / static_cast<float>(nsit) * 100.0f); |
| 109 | }; | 111 | }; |
| 110 | 112 | ||
| 111 | l_.info("Building indices for {} situations", nsit); | 113 | l_.info("Building indexes for {} situations", nsit); |
| 112 | |||
| 113 | auto const buffer_distance = min_distance_; | ||
| 114 | auto const points_per_circle = 8; | ||
| 115 | // Note: thomas strategy does not work for geographic_join_round; | ||
| 116 | // need to use andoyer for that. | ||
| 117 | using formula = bgeo::strategy::andoyer; | ||
| 118 | bgeo::strategy::buffer::distance_symmetric<double> distance_strategy{ | ||
| 119 | buffer_distance | ||
| 120 | }; | ||
| 121 | bgeo::strategy::buffer::geographic_join_round<formula> join_strategy{ | ||
| 122 | points_per_circle | ||
| 123 | }; | ||
| 124 | bgeo::strategy::buffer::geographic_end_round<formula> end_strategy{4}; | ||
| 125 | bgeo::strategy::buffer::geographic_point_circle<formula> circle_strategy{ | ||
| 126 | points_per_circle | ||
| 127 | }; | ||
| 128 | bgeo::strategy::buffer::geographic_side_straight<formula> side_strategy; | ||
| 129 | 114 | ||
| 130 | l_.info("Progress: 0%"); | 115 | l_.info("Progress: 0%"); |
| 131 | auto const before_build = chrono::steady_clock::now(); | 116 | auto const before_build = chrono::steady_clock::now(); |
| @@ -135,17 +120,11 @@ handler::handler(log::logger const& l, datex2::situation_publication pub) | |||
| 135 | { | 120 | { |
| 136 | for (auto const& ls : rc->relevant_line_strings) | 121 | for (auto const& ls : rc->relevant_line_strings) |
| 137 | { | 122 | { |
| 138 | auto buffered_ls = geo::multi_polygon{}; | 123 | blse_index_.insert(*ls, rc); |
| 139 | bgeo::buffer( | ||
| 140 | *ls, buffered_ls, distance_strategy, side_strategy, join_strategy, | ||
| 141 | end_strategy, circle_strategy); | ||
| 142 | auto blse = geo::box{}; | ||
| 143 | bgeo::envelope(buffered_ls, blse); | ||
| 144 | blse_index_.insert(std::make_tuple(blse, buffered_ls, rc)); | ||
| 145 | } | 124 | } |
| 146 | for (auto p : rc->relevant_points) | 125 | for (auto p : rc->relevant_points) |
| 147 | { | 126 | { |
| 148 | p_index_.insert(std::make_pair(p, rc)); | 127 | p_index_.insert(p_index_value{p, rc}); |
| 149 | } | 128 | } |
| 150 | } | 129 | } |
| 151 | 130 | ||
| @@ -160,7 +139,7 @@ handler::handler(log::logger const& l, datex2::situation_publication pub) | |||
| 160 | 139 | ||
| 161 | auto const after_build = chrono::steady_clock::now(); | 140 | auto const after_build = chrono::steady_clock::now(); |
| 162 | auto const dur_build_s = chrono::duration<double>{after_build - before_build}; | 141 | auto const dur_build_s = chrono::duration<double>{after_build - before_build}; |
| 163 | l_.info("Indices built in {}", dur_build_s); | 142 | l_.info("Indexes built in {}", dur_build_s); |
| 164 | l_.info("BLSE index size: {}", blse_index_.size()); | 143 | l_.info("BLSE index size: {}", blse_index_.size()); |
| 165 | l_.info("Point index size: {}", p_index_.size()); | 144 | l_.info("Point index size: {}", p_index_.size()); |
| 166 | } | 145 | } |
| @@ -175,16 +154,13 @@ auto handler::process_gpx(gpx::file&& gpx_file) | |||
| 175 | auto const check_periods = time::period_seq{relevant.begin(), relevant.end()}; | 154 | auto const check_periods = time::period_seq{relevant.begin(), relevant.end()}; |
| 176 | 155 | ||
| 177 | // TODO: eliminate use of overlap segments | 156 | // TODO: eliminate use of overlap segments |
| 178 | auto splits_with_overlap_segments = std::vector<geo::linestring>{}; | 157 | auto segments = std::vector<geo::utm::zonable_wgs84_linestring>{}; |
| 179 | for (auto const& track : gpx_file.tracks) | 158 | for (auto const& track : gpx_file.tracks) |
| 180 | for (auto const& seg : track.segments) | 159 | for (auto const& seg : track.segments) |
| 181 | geo::split_linestring_with_overlap_segments( | 160 | segments.push_back(seg.waypoints); |
| 182 | seg.waypoints, | ||
| 183 | 5000 /* meters max total dist until a new split is forced */, | ||
| 184 | splits_with_overlap_segments); | ||
| 185 | auto const before_query = chrono::steady_clock::now(); | 161 | auto const before_query = chrono::steady_clock::now(); |
| 186 | 162 | ||
| 187 | auto vincenty_strategy = geo::vincenty_strategy{}; | 163 | auto vincenty_strategy = geo::wgs84::vincenty_strategy{}; |
| 188 | 164 | ||
| 189 | l_.debug("Querying for relevant situations"); | 165 | l_.debug("Querying for relevant situations"); |
| 190 | auto relevant_road_closures = | 166 | auto relevant_road_closures = |
| @@ -192,28 +168,17 @@ auto handler::process_gpx(gpx::file&& gpx_file) | |||
| 192 | auto ls_checked = 0uz; | 168 | auto ls_checked = 0uz; |
| 193 | auto p_checked = 0uz; | 169 | auto p_checked = 0uz; |
| 194 | auto i = 0; | 170 | auto i = 0; |
| 195 | for (geo::linestring const& part : splits_with_overlap_segments) | 171 | for (auto const& part : segments) |
| 196 | { | 172 | { |
| 197 | l_.debug("Checking part [{}/{}]", ++i, splits_with_overlap_segments.size()); | 173 | l_.debug("Checking part [{}/{}]", ++i, segments.size()); |
| 198 | 174 | ||
| 175 | auto part_zone_lss = geo::multizonal::zoned_linestring_seg_seq{part}; | ||
| 199 | // TODO: consider buffering with min_distance_ | 176 | // TODO: consider buffering with min_distance_ |
| 200 | auto part_box = geo::box{}; | 177 | auto part_box = geo::wgs84::box{}; |
| 201 | bgeo::envelope(part, part_box); | 178 | bgeo::envelope(part, part_box); |
| 202 | 179 | ||
| 203 | for (auto it = blse_index_.qbegin(bgeo::index::intersects(part_box)); | 180 | for (auto const& rc : blse_index_.intersection(part_zone_lss)) |
| 204 | it != blse_index_.qend(); it++) | ||
| 205 | { | 181 | { |
| 206 | // Cannot use structured bindings here, as boost::geometry::get | ||
| 207 | // interferes with ADL. It is a candidate as the namespace | ||
| 208 | // boost::geometry is part of the associated namespace set, which | ||
| 209 | // happens because geo::linestring ≡ | ||
| 210 | // boost::geometry::model::linestring<geo::point> is part of the | ||
| 211 | // whole tuple type (blse_index_value) that is the value_type of the | ||
| 212 | // iterator. | ||
| 213 | |||
| 214 | geo::box const& blse = std::get<0>(*it); | ||
| 215 | geo::multi_polygon const& bls = std::get<1>(*it); | ||
| 216 | std::shared_ptr<datex2::road_closure> const& rc = std::get<2>(*it); | ||
| 217 | if (rc->validity | 182 | if (rc->validity |
| 218 | && rc->validity->intersect(check_periods).periods().empty()) | 183 | && rc->validity->intersect(check_periods).periods().empty()) |
| 219 | continue; | 184 | continue; |
| @@ -221,8 +186,7 @@ auto handler::process_gpx(gpx::file&& gpx_file) | |||
| 221 | // auto small_parts = bgeo::model::multi_linestring<geo::linestring>{}; | 186 | // auto small_parts = bgeo::model::multi_linestring<geo::linestring>{}; |
| 222 | // bgeo::intersection(part, lse, small_parts); | 187 | // bgeo::intersection(part, lse, small_parts); |
| 223 | 188 | ||
| 224 | if (bgeo::intersects(blse, part) && bgeo::intersects(bls, part)) | 189 | relevant_road_closures.emplace(rc); |
| 225 | relevant_road_closures.emplace(rc); | ||
| 226 | 190 | ||
| 227 | // if (!bgeo::is_empty(small_parts)) | 191 | // if (!bgeo::is_empty(small_parts)) |
| 228 | // if (bgeo::distance(*ls, small_parts, vincenty_strategy) < | 192 | // if (bgeo::distance(*ls, small_parts, vincenty_strategy) < |
| @@ -234,7 +198,7 @@ auto handler::process_gpx(gpx::file&& gpx_file) | |||
| 234 | it != p_index_.qend(); it++) | 198 | it != p_index_.qend(); it++) |
| 235 | { | 199 | { |
| 236 | // Cannot use structured bindings here for the same reason as above. | 200 | // Cannot use structured bindings here for the same reason as above. |
| 237 | geo::point const& p = std::get<0>(*it); | 201 | geo::utm::zonable_wgs84_point const& p = std::get<0>(*it); |
| 238 | std::shared_ptr<datex2::road_closure> const& rc = std::get<1>(*it); | 202 | std::shared_ptr<datex2::road_closure> const& rc = std::get<1>(*it); |
| 239 | if (rc->validity | 203 | if (rc->validity |
| 240 | && rc->validity->intersect(check_periods).periods().empty()) | 204 | && rc->validity->intersect(check_periods).periods().empty()) |
| @@ -301,10 +265,12 @@ auto handler::process_gpx(gpx::file&& gpx_file) | |||
| 301 | .relevant_lss = | 265 | .relevant_lss = |
| 302 | rc->relevant_line_strings | 266 | rc->relevant_line_strings |
| 303 | | views::transform( | 267 | | views::transform( |
| 304 | [](auto const& lsp) -> geo::linestring | 268 | [](auto const& lsp) |
| 269 | -> geo::utm::zonable_wgs84_linestring | ||
| 305 | { return *lsp; }) | 270 | { return *lsp; }) |
| 306 | | std::ranges:: | 271 | | std::ranges::to<std::vector< |
| 307 | to<std::vector<geo::linestring>>(), | 272 | geo::utm::zonable_wgs84_linestring |
| 273 | >>(), | ||
| 308 | }; | 274 | }; |
| 309 | }) | 275 | }) |
| 310 | | std::ranges::to<std::vector<relevant_road_closure>>(), | 276 | | std::ranges::to<std::vector<relevant_road_closure>>(), |
diff --git a/server/src/api.cppm b/server/src/api.cppm index b0f3944..77e0f85 100644 --- a/server/src/api.cppm +++ b/server/src/api.cppm | |||
| @@ -7,7 +7,8 @@ export module routemon:api; | |||
| 7 | 7 | ||
| 8 | import std; | 8 | import std; |
| 9 | import :datex2; | 9 | import :datex2; |
| 10 | import :geo; | 10 | import :geo.wgs84; |
| 11 | import :geo.multizonal; | ||
| 11 | import :gpx; | 12 | import :gpx; |
| 12 | import :log; | 13 | import :log; |
| 13 | import :time; | 14 | import :time; |
| @@ -25,20 +26,20 @@ export namespace routemon::api { | |||
| 25 | 26 | ||
| 26 | struct relevant_road_closure | 27 | struct relevant_road_closure |
| 27 | { | 28 | { |
| 28 | std::vector<geo::linestring> relevant_lss; | 29 | std::vector<geo::utm::zonable_wgs84_linestring> relevant_lss; |
| 29 | }; | 30 | }; |
| 30 | 31 | ||
| 31 | struct relevant_situation | 32 | struct relevant_situation |
| 32 | { | 33 | { |
| 33 | std::string id; | 34 | std::string id; |
| 34 | std::optional<geo::point> location; | 35 | std::optional<geo::utm::zonable_wgs84_point> location; |
| 35 | std::vector<std::string> comments; | 36 | std::vector<std::string> comments; |
| 36 | std::vector<relevant_road_closure> relevant_road_closures; | 37 | std::vector<relevant_road_closure> relevant_road_closures; |
| 37 | }; | 38 | }; |
| 38 | 39 | ||
| 39 | struct track_segment | 40 | struct track_segment |
| 40 | { | 41 | { |
| 41 | geo::linestring points; | 42 | geo::utm::zonable_wgs84_linestring points; |
| 42 | }; | 43 | }; |
| 43 | 44 | ||
| 44 | struct track | 45 | struct track |
| @@ -79,13 +80,11 @@ class handler | |||
| 79 | { | 80 | { |
| 80 | constexpr static auto const min_distance_ = 5.0; | 81 | constexpr static auto const min_distance_ = 5.0; |
| 81 | 82 | ||
| 82 | using blse_index_value = std::tuple< | 83 | using blse_index_value = std::shared_ptr<datex2::road_closure>; |
| 83 | geo::box, geo::multi_polygon, std::shared_ptr<datex2::road_closure> | 84 | using blse_index = geo::multizonal::linestring_rtree<blse_index_value>; |
| 85 | using p_index_value = std::pair< | ||
| 86 | geo::utm::zonable_wgs84_point, std::shared_ptr<datex2::road_closure> | ||
| 84 | >; | 87 | >; |
| 85 | using p_index_value = | ||
| 86 | std::pair<geo::point, std::shared_ptr<datex2::road_closure>>; | ||
| 87 | using blse_index = | ||
| 88 | bgeo::index::rtree<blse_index_value, bgeo::index::quadratic<16>>; | ||
| 89 | using p_index = bgeo::index::rtree<p_index_value, bgeo::index::quadratic<16>>; | 88 | using p_index = bgeo::index::rtree<p_index_value, bgeo::index::quadratic<16>>; |
| 90 | 89 | ||
| 91 | log::logger l_; | 90 | log::logger l_; |
diff --git a/server/src/datex2.cpp b/server/src/datex2.cpp index 4c8000b..70c5f94 100644 --- a/server/src/datex2.cpp +++ b/server/src/datex2.cpp | |||
| @@ -51,32 +51,41 @@ auto loader::add_location_from_xml( | |||
| 51 | } | 51 | } |
| 52 | auto const pos_list_str = | 52 | auto const pos_list_str = |
| 53 | std::string_view{loc_gml_xml.child_value("loc:posList")}; | 53 | std::string_view{loc_gml_xml.child_value("loc:posList")}; |
| 54 | // lat1 long1 lat2 long2 ... lat(n-1) long(n-1) latn longn | 54 | // lat1 lon1 lat2 lon2 ... lat(n-1) lon(n-1) latn lonn |
| 55 | 55 | ||
| 56 | auto ls = std::make_shared<geo::linestring>(); | 56 | auto ls = std::make_shared<geo::utm::zonable_wgs84_linestring>(); |
| 57 | 57 | ||
| 58 | auto lat_set = false; | 58 | auto lat_set = false; |
| 59 | auto lat = 0.0; | 59 | auto lat = 0.0; |
| 60 | for (auto const lat_or_long_str : std::views::split(pos_list_str, " "sv)) | 60 | for (auto const lat_or_lon_str : std::views::split(pos_list_str, " "sv)) |
| 61 | { | 61 | { |
| 62 | auto mlat_or_long = util::parse_double(std::string_view{lat_or_long_str}); | 62 | auto mlat_or_lon = util::parse_double(std::string_view{lat_or_lon_str}); |
| 63 | if (!mlat_or_long) | 63 | if (!mlat_or_lon) |
| 64 | { | 64 | { |
| 65 | warnings_.insert( | 65 | warnings_.insert( |
| 66 | std::format( | 66 | std::format( |
| 67 | "failed to parse coordinate {:?}", | 67 | "failed to parse coordinate {:?}", |
| 68 | std::string_view{lat_or_long_str})); | 68 | std::string_view{lat_or_lon_str})); |
| 69 | return; | 69 | return; |
| 70 | } | 70 | } |
| 71 | 71 | ||
| 72 | if (!lat_set) | 72 | if (!lat_set) |
| 73 | { | 73 | { |
| 74 | lat = *mlat_or_long; | 74 | lat = *mlat_or_lon; |
| 75 | lat_set = true; | 75 | lat_set = true; |
| 76 | } | 76 | } |
| 77 | else | 77 | else |
| 78 | { | 78 | { |
| 79 | bgeo::append(*ls, geo::point{*mlat_or_long, lat}); | 79 | auto mp = geo::utm::zonable_wgs84_point::from( |
| 80 | geo::wgs84::normalized_point{geo::wgs84::from_lat_lon( | ||
| 81 | lat, *mlat_or_lon)}); | ||
| 82 | if (!mp) | ||
| 83 | warnings_.insert( | ||
| 84 | std::format( | ||
| 85 | "point ({}, {}) outside of UTM range; ignoring it", lat, | ||
| 86 | *mlat_or_lon)); | ||
| 87 | else | ||
| 88 | bgeo::append(*ls, *mp); | ||
| 80 | lat = 0; | 89 | lat = 0; |
| 81 | lat_set = false; | 90 | lat_set = false; |
| 82 | } | 91 | } |
| @@ -110,11 +119,24 @@ auto loader::add_location_from_xml( | |||
| 110 | // maar heeft het UML-model van DATEX II v3 het over ETRS 89: | 119 | // maar heeft het UML-model van DATEX II v3 het over ETRS 89: |
| 111 | // https://docs.datex2.eu/_static/data/v3.7/umlmodel/html/EARoot/EA3/EA3/EA5/EA676.htm | 120 | // https://docs.datex2.eu/_static/data/v3.7/umlmodel/html/EARoot/EA3/EA3/EA5/EA676.htm |
| 112 | 121 | ||
| 113 | auto const coords_etrs89 = geo::point{*mlon, *mlat}; | 122 | auto const coords_etrs89 = geo::wgs84::point{*mlon, *mlat}; |
| 114 | auto coords_wgs84 = geo::point{}; | 123 | auto coords_wgs84 = geo::wgs84::point{}; |
| 115 | etrs89_to_wgs84_.forward(coords_etrs89, coords_wgs84); | 124 | etrs89_to_wgs84_.forward(coords_etrs89, coords_wgs84); |
| 116 | 125 | ||
| 117 | rc.relevant_points.push_back(coords_wgs84); | 126 | auto mp = geo::utm::zonable_wgs84_point::from( |
| 127 | geo::wgs84::normalized_point{coords_wgs84}); | ||
| 128 | if (!mp) | ||
| 129 | { | ||
| 130 | warnings_.insert( | ||
| 131 | std::format( | ||
| 132 | "(ETRS89) point ({}, {}) outside of UTM range after " | ||
| 133 | "transformation to WGS 84; ignoring it", | ||
| 134 | *mlat, *mlon)); | ||
| 135 | } | ||
| 136 | else | ||
| 137 | { | ||
| 138 | rc.relevant_points.push_back(*mp); | ||
| 139 | } | ||
| 118 | } | 140 | } |
| 119 | else | 141 | else |
| 120 | { | 142 | { |
| @@ -347,13 +369,23 @@ auto loader::load_situation_publication(std::string const& filename) | |||
| 347 | // het over ETRS 89: | 369 | // het over ETRS 89: |
| 348 | // https://docs.datex2.eu/_static/data/v3.7/umlmodel/html/EARoot/EA3/EA3/EA5/EA676.htm | 370 | // https://docs.datex2.eu/_static/data/v3.7/umlmodel/html/EARoot/EA3/EA3/EA5/EA676.htm |
| 349 | 371 | ||
| 350 | auto const coords_etrs89 = geo::point{*mlon, *mlat}; | 372 | auto const coords_etrs89 = geo::wgs84::point{*mlon, *mlat}; |
| 351 | auto coords_wgs84 = geo::point{}; | 373 | auto coords_wgs84 = geo::wgs84::point{}; |
| 352 | etrs89_to_wgs84_.forward(coords_etrs89, coords_wgs84); | 374 | etrs89_to_wgs84_.forward(coords_etrs89, coords_wgs84); |
| 353 | 375 | ||
| 354 | if (!sit->location) | 376 | auto mp = geo::utm::zonable_wgs84_point::from( |
| 377 | geo::wgs84::normalized_point{coords_wgs84}); | ||
| 378 | if (!mp) | ||
| 379 | { | ||
| 380 | warnings_.insert( | ||
| 381 | std::format( | ||
| 382 | "(ETRS89) point ({}, {}) outside of UTM range after " | ||
| 383 | "transformation to WGS 84; ignoring it", | ||
| 384 | *mlat, *mlon)); | ||
| 385 | } | ||
| 386 | else if (!sit->location) | ||
| 355 | { | 387 | { |
| 356 | sit->location = coords_wgs84; | 388 | sit->location = *mp; |
| 357 | } | 389 | } |
| 358 | } | 390 | } |
| 359 | } | 391 | } |
diff --git a/server/src/datex2.cppm b/server/src/datex2.cppm index 9c75013..ae66631 100644 --- a/server/src/datex2.cppm +++ b/server/src/datex2.cppm | |||
| @@ -1,6 +1,5 @@ | |||
| 1 | module; | 1 | module; |
| 2 | 2 | ||
| 3 | #include <boost/geometry/algorithms/is_empty.hpp> | ||
| 4 | #include <boost/geometry/srs/epsg.hpp> | 3 | #include <boost/geometry/srs/epsg.hpp> |
| 5 | #include <boost/geometry/srs/transformation.hpp> | 4 | #include <boost/geometry/srs/transformation.hpp> |
| 6 | 5 | ||
| @@ -9,7 +8,7 @@ module; | |||
| 9 | export module routemon:datex2; | 8 | export module routemon:datex2; |
| 10 | 9 | ||
| 11 | import std; | 10 | import std; |
| 12 | import :geo; | 11 | import :geo.utm; |
| 13 | import :time; | 12 | import :time; |
| 14 | import :util; | 13 | import :util; |
| 15 | 14 | ||
| @@ -23,14 +22,15 @@ export struct road_closure | |||
| 23 | { | 22 | { |
| 24 | std::weak_ptr<situation> parent; | 23 | std::weak_ptr<situation> parent; |
| 25 | std::optional<time::period_seq> validity; | 24 | std::optional<time::period_seq> validity; |
| 26 | std::vector<geo::point> relevant_points = {}; | 25 | std::vector<geo::utm::zonable_wgs84_point> relevant_points = {}; |
| 27 | std::vector<std::shared_ptr<geo::linestring>> relevant_line_strings = {}; | 26 | std::vector<std::shared_ptr<geo::utm::zonable_wgs84_linestring>> |
| 27 | relevant_line_strings = {}; | ||
| 28 | }; | 28 | }; |
| 29 | 29 | ||
| 30 | export struct situation | 30 | export struct situation |
| 31 | { | 31 | { |
| 32 | std::string id; | 32 | std::string id; |
| 33 | std::optional<geo::point> location = | 33 | std::optional<geo::utm::zonable_wgs84_point> location = |
| 34 | std::nullopt; // as shown on the map, not used for querying | 34 | std::nullopt; // as shown on the map, not used for querying |
| 35 | std::vector<std::string> comments = {}; | 35 | std::vector<std::string> comments = {}; |
| 36 | std::vector<std::shared_ptr<road_closure>> road_closures = {}; | 36 | std::vector<std::shared_ptr<road_closure>> road_closures = {}; |
| @@ -47,7 +47,7 @@ export class loader | |||
| 47 | // ETRS 89 (EPSG:4258) -> WGS 84 (EPSG:4326) | 47 | // ETRS 89 (EPSG:4258) -> WGS 84 (EPSG:4326) |
| 48 | bgeo::srs:: | 48 | bgeo::srs:: |
| 49 | transformation<bgeo::srs::static_epsg<4258>, bgeo::srs::static_epsg<4326>> | 49 | transformation<bgeo::srs::static_epsg<4258>, bgeo::srs::static_epsg<4326>> |
| 50 | etrs89_to_wgs84_{}; | 50 | etrs89_to_wgs84_ = {}; |
| 51 | 51 | ||
| 52 | std::multiset<std::string> warnings_; | 52 | std::multiset<std::string> warnings_; |
| 53 | 53 | ||
diff --git a/server/src/geo.cpp b/server/src/geo.cpp deleted file mode 100644 index 1776f3c..0000000 --- a/server/src/geo.cpp +++ /dev/null | |||
| @@ -1,42 +0,0 @@ | |||
| 1 | module; | ||
| 2 | |||
| 3 | #include <boost/geometry.hpp> | ||
| 4 | |||
| 5 | module routemon:geo$impl; | ||
| 6 | |||
| 7 | import :geo; | ||
| 8 | |||
| 9 | namespace routemon::geo { | ||
| 10 | |||
| 11 | auto split_linestring_with_overlap_segments( | ||
| 12 | linestring const& ls, double max_split_distance_m, | ||
| 13 | std::vector<linestring>& append_to) -> void | ||
| 14 | { | ||
| 15 | if (bgeo::is_empty(ls)) | ||
| 16 | return; | ||
| 17 | |||
| 18 | auto current_ls = linestring{}; | ||
| 19 | auto current_ls_length = 0.0; | ||
| 20 | auto previous = std::optional<point>{}; | ||
| 21 | bgeo::for_each_point( | ||
| 22 | ls, | ||
| 23 | [&](point p) -> void | ||
| 24 | { | ||
| 25 | bgeo::append(current_ls, p); | ||
| 26 | if (previous) | ||
| 27 | { | ||
| 28 | auto d = bgeo::distance(*previous, p, vincenty_strategy()); | ||
| 29 | current_ls_length += d; | ||
| 30 | if (current_ls_length > max_split_distance_m) | ||
| 31 | { | ||
| 32 | append_to.push_back(std::move(current_ls)); | ||
| 33 | current_ls = linestring{*previous, p}; | ||
| 34 | current_ls_length = d; | ||
| 35 | } | ||
| 36 | } | ||
| 37 | previous = p; | ||
| 38 | }); | ||
| 39 | append_to.emplace_back(std::move(current_ls)); | ||
| 40 | } | ||
| 41 | |||
| 42 | } // namespace routemon::geo | ||
diff --git a/server/src/geo.cppm b/server/src/geo.cppm deleted file mode 100644 index b3513c8..0000000 --- a/server/src/geo.cppm +++ /dev/null | |||
| @@ -1,23 +0,0 @@ | |||
| 1 | module; | ||
| 2 | |||
| 3 | #include <boost/geometry.hpp> | ||
| 4 | |||
| 5 | export module routemon:geo; | ||
| 6 | |||
| 7 | export namespace bgeo = boost::geometry; | ||
| 8 | |||
| 9 | export namespace routemon::geo { | ||
| 10 | |||
| 11 | using point = bgeo::model::point<double, 2, bgeo::cs::geographic<bgeo::degree>>; | ||
| 12 | using linestring = bgeo::model::linestring<point>; | ||
| 13 | using box = bgeo::model::box<point>; | ||
| 14 | using polygon = bgeo::model::polygon<geo::point>; | ||
| 15 | using multi_polygon = bgeo::model::multi_polygon<polygon>; | ||
| 16 | using stype = bgeo::srs::spheroid<double>; | ||
| 17 | using vincenty_strategy = bgeo::strategy::distance::vincenty<stype>; | ||
| 18 | |||
| 19 | auto split_linestring_with_overlap_segments( | ||
| 20 | linestring const& ls, double max_split_distance_m, | ||
| 21 | std::vector<linestring>& append_to) -> void; | ||
| 22 | |||
| 23 | } // namespace routemon::geo | ||
diff --git a/server/src/geo/geo.cppm b/server/src/geo/geo.cppm new file mode 100644 index 0000000..446f1eb --- /dev/null +++ b/server/src/geo/geo.cppm | |||
| @@ -0,0 +1,7 @@ | |||
| 1 | module; | ||
| 2 | |||
| 3 | #include <boost/geometry.hpp> | ||
| 4 | |||
| 5 | export module routemon:geo; | ||
| 6 | |||
| 7 | namespace bgeo = boost::geometry; | ||
diff --git a/server/src/geo/multizonal.cppm b/server/src/geo/multizonal.cppm new file mode 100644 index 0000000..f732994 --- /dev/null +++ b/server/src/geo/multizonal.cppm | |||
| @@ -0,0 +1,240 @@ | |||
| 1 | module; | ||
| 2 | |||
| 3 | #include <boost/geometry.hpp> | ||
| 4 | #include <boost/geometry/srs/epsg.hpp> | ||
| 5 | #include <boost/geometry/srs/transformation.hpp> | ||
| 6 | |||
| 7 | export module routemon:geo.multizonal; | ||
| 8 | |||
| 9 | import std; | ||
| 10 | import :geo; | ||
| 11 | import :geo.utm; | ||
| 12 | import :geo.utm.zone_local; | ||
| 13 | import :geo.wgs84; | ||
| 14 | |||
| 15 | namespace routemon::geo::multizonal { | ||
| 16 | |||
| 17 | template<std::default_initializable T> | ||
| 18 | class multi_zone | ||
| 19 | { | ||
| 20 | std::array<T, utm::zone::max().as_index() + 1> zones_; | ||
| 21 | |||
| 22 | public: | ||
| 23 | auto operator[](utm::zone z) -> T& | ||
| 24 | { | ||
| 25 | return zones_[z.as_index()]; | ||
| 26 | } | ||
| 27 | |||
| 28 | auto operator[](utm::zone z) const -> T const& | ||
| 29 | { | ||
| 30 | return zones_[z.as_index()]; | ||
| 31 | } | ||
| 32 | }; | ||
| 33 | |||
| 34 | namespace wgs84 | ||
| 35 | { | ||
| 36 | |||
| 37 | // Transformations from WGS 84 (EPSG:4326) | ||
| 38 | using transform_from_t = bgeo::srs::transformation<bgeo::srs::static_epsg<4326>>; | ||
| 39 | |||
| 40 | class utm_transform : public transform_from_t | ||
| 41 | { | ||
| 42 | utm::zone to_zone_; | ||
| 43 | |||
| 44 | public: | ||
| 45 | explicit utm_transform(utm::zone to_zone) | ||
| 46 | : transform_from_t{{}, bgeo::srs::epsg{to_zone.wgs84_proj_epsg()}}, | ||
| 47 | to_zone_{to_zone} | ||
| 48 | { | ||
| 49 | } | ||
| 50 | |||
| 51 | utm_transform() | ||
| 52 | : utm_transform{utm::zone::min()} | ||
| 53 | {} | ||
| 54 | |||
| 55 | auto apply(utm::zonable_wgs84_point p) -> utm::zone_local::point | ||
| 56 | { | ||
| 57 | auto local_p = utm::zone_local::point{to_zone_, 0.0, 0.0}; | ||
| 58 | forward(p, local_p); | ||
| 59 | return local_p; | ||
| 60 | } | ||
| 61 | }; | ||
| 62 | |||
| 63 | class utm_transforms : public multi_zone<utm_transform> | ||
| 64 | { | ||
| 65 | utm_transforms() | ||
| 66 | { | ||
| 67 | for (auto z = utm::zone::min(); z != utm::zone::max(); z = z.next()) | ||
| 68 | { | ||
| 69 | (*this)[z] = utm_transform{z}; | ||
| 70 | } | ||
| 71 | } | ||
| 72 | |||
| 73 | public: | ||
| 74 | static utm_transforms const& instance() | ||
| 75 | { | ||
| 76 | static utm_transforms inst; | ||
| 77 | return inst; | ||
| 78 | } | ||
| 79 | }; | ||
| 80 | |||
| 81 | auto neighbor_utm_zone(utm::zonable_wgs84_point p) -> std::pair<utm::zone, double> | ||
| 82 | { | ||
| 83 | auto separating_meridian_lon = std::round(p.lon() / 6.0) * 6.0; | ||
| 84 | auto closest_zone_middle = separating_meridian_lon < p.lon() | ||
| 85 | ? separating_meridian_lon - 3.0 | ||
| 86 | : separating_meridian_lon + 3.0; | ||
| 87 | auto closest_zone = utm::zone::for_wgs84_point(utm::zonable_wgs84_point::from(geo::wgs84::normalized_point{geo::wgs84::from_lat_lon(p.lat(), closest_zone_middle)}).value()); | ||
| 88 | |||
| 89 | auto separating_meridian_ls = geo::wgs84::linestring{ | ||
| 90 | geo::wgs84::from_lat_lon(90.0, 0.0), | ||
| 91 | geo::wgs84::from_lat_lon(0.0, separating_meridian_lon), | ||
| 92 | geo::wgs84::from_lat_lon(-90.0, 0.0), | ||
| 93 | }; | ||
| 94 | auto closest_zone_dist = bgeo::distance(p, separating_meridian_ls, geo::wgs84::vincenty_strategy{}); | ||
| 95 | |||
| 96 | return std::make_pair(closest_zone, closest_zone_dist); | ||
| 97 | } | ||
| 98 | |||
| 99 | } // namespace wgs84 | ||
| 100 | |||
| 101 | struct multizone_linestring | ||
| 102 | { | ||
| 103 | multi_zone<std::vector<utm::zone_local::prim::linestring>> segments; | ||
| 104 | |||
| 105 | // TODO: consider taking an input range instead | ||
| 106 | explicit multizone_linestring(utm::zonable_wgs84_linestring const& ls) | ||
| 107 | { | ||
| 108 | auto to_utm = wgs84::utm_transforms::instance(); | ||
| 109 | auto working = multi_zone<utm::zone_local::prim::linestring>{}; | ||
| 110 | auto mprev_p_zone = std::optional<utm::zone>{}; | ||
| 111 | auto mprev_p_alt_zone = std::optional<utm::zone>{}; | ||
| 112 | |||
| 113 | auto push = [&](utm::zonable_wgs84_point p, utm::zone z) | ||
| 114 | { working[z].push_back(to_utm[z].apply(p)); }; | ||
| 115 | auto flush = [&](utm::zone z) | ||
| 116 | { | ||
| 117 | segments[z].push_back(std::move(working[z])); | ||
| 118 | working[z] = {}; | ||
| 119 | }; | ||
| 120 | |||
| 121 | for (auto const& p : ls) | ||
| 122 | { | ||
| 123 | auto p_zone = utm::zone::for_wgs84_point(p); | ||
| 124 | auto mp_alt_zone = std::optional<utm::zone>{}; | ||
| 125 | if (auto [neighbor_zone, neighbor_zone_dist] = wgs84::neighbor_utm_zone(p); | ||
| 126 | neighbor_zone_dist < 30.0 /* m */) | ||
| 127 | mp_alt_zone = neighbor_zone; | ||
| 128 | assert(!mp_alt_zone || *mp_alt_zone != p_zone); | ||
| 129 | |||
| 130 | push(p, p_zone); | ||
| 131 | if (mp_alt_zone) | ||
| 132 | push(p, *mp_alt_zone); | ||
| 133 | if (mprev_p_zone && *mprev_p_zone != p_zone && mprev_p_zone != mp_alt_zone) | ||
| 134 | { | ||
| 135 | push(p, *mprev_p_zone); | ||
| 136 | flush(*mprev_p_zone); | ||
| 137 | } | ||
| 138 | if (mprev_p_alt_zone && *mprev_p_alt_zone != p_zone && mprev_p_alt_zone != mp_alt_zone) | ||
| 139 | { | ||
| 140 | push(p, *mprev_p_alt_zone); | ||
| 141 | flush(*mprev_p_alt_zone); | ||
| 142 | } | ||
| 143 | |||
| 144 | mprev_p_zone = p_zone; | ||
| 145 | mprev_p_alt_zone = mp_alt_zone; | ||
| 146 | } | ||
| 147 | |||
| 148 | if (mprev_p_zone && !working[*mprev_p_zone].empty()) | ||
| 149 | { | ||
| 150 | flush(*mprev_p_zone); | ||
| 151 | } | ||
| 152 | if (mprev_p_alt_zone && !working[*mprev_p_alt_zone].empty()) | ||
| 153 | { | ||
| 154 | flush(*mprev_p_alt_zone); | ||
| 155 | } | ||
| 156 | } | ||
| 157 | }; | ||
| 158 | |||
| 159 | struct zoned_linestring_seg_seq | ||
| 160 | { | ||
| 161 | std::vector<utm::zone_local::linestring> segments; | ||
| 162 | |||
| 163 | explicit zoned_linestring_seg_seq(utm::zonable_wgs84_linestring ls) | ||
| 164 | { | ||
| 165 | auto to_utm = wgs84::utm_transforms::instance(); | ||
| 166 | auto mworking_seg = std::optional<utm::zone_local::linestring>{}; | ||
| 167 | |||
| 168 | for (auto const& p : ls) | ||
| 169 | { | ||
| 170 | auto p_zone = utm::zone::for_wgs84_point(p); | ||
| 171 | if (mworking_seg && mworking_seg->zone != p_zone) | ||
| 172 | { | ||
| 173 | mworking_seg->push_back(to_utm[mworking_seg->zone].apply(p)); | ||
| 174 | segments.push_back(std::move(*mworking_seg)); | ||
| 175 | mworking_seg = std::nullopt; | ||
| 176 | } | ||
| 177 | if (!mworking_seg) | ||
| 178 | mworking_seg = utm::zone_local::linestring{p_zone}; | ||
| 179 | mworking_seg->push_back(to_utm[p_zone].apply(p)); | ||
| 180 | } | ||
| 181 | |||
| 182 | if (mworking_seg && !mworking_seg->empty()) | ||
| 183 | { | ||
| 184 | segments.push_back(std::move(*mworking_seg)); | ||
| 185 | } | ||
| 186 | } | ||
| 187 | }; | ||
| 188 | |||
| 189 | template<class T> | ||
| 190 | class linestring_rtree | ||
| 191 | { | ||
| 192 | public: | ||
| 193 | using index_value = std::tuple<utm::zone_local::prim::box, utm::zone_local::prim::linestring, T>; | ||
| 194 | |||
| 195 | private: | ||
| 196 | multi_zone<bgeo::index::rtree<index_value, bgeo::index::quadratic<16>>> local_rtrees_; | ||
| 197 | |||
| 198 | public: | ||
| 199 | template<std::convertible_to<T> U> | ||
| 200 | auto insert(utm::zonable_wgs84_linestring const& ls, U&& arg) -> void | ||
| 201 | { | ||
| 202 | auto zone_segments = multizone_linestring{ls}.segments; | ||
| 203 | for (auto z = utm::zone::min(); z != utm::zone::max(); z = z.next()) | ||
| 204 | { | ||
| 205 | for (auto ls : zone_segments[z]) | ||
| 206 | { | ||
| 207 | utm::zone_local::prim::box blse = | ||
| 208 | bgeo::return_buffer<utm::zone_local::prim::box>(bgeo::return_envelope<utm::zone_local::prim::box>(ls), 30.0 /* m */); | ||
| 209 | local_rtrees_[z].insert(index_value{blse, std::move(ls), std::forward<U>(arg)}); | ||
| 210 | } | ||
| 211 | } | ||
| 212 | } | ||
| 213 | |||
| 214 | auto size() const -> std::size_t | ||
| 215 | { | ||
| 216 | auto total_size = 0uz; | ||
| 217 | for (auto z = utm::zone::min(); z != utm::zone::max(); z = z.next()) | ||
| 218 | total_size += local_rtrees_[z].size(); | ||
| 219 | return total_size; | ||
| 220 | } | ||
| 221 | |||
| 222 | // Note: the same T may be generated more than once! | ||
| 223 | auto intersection(zoned_linestring_seg_seq const& lss) const -> std::generator<T const&> | ||
| 224 | { | ||
| 225 | for (auto const& ls : lss.segments) | ||
| 226 | { | ||
| 227 | // auto ls_box = bgeo::return_envelope<utm::zone_local::prim::box>(ls); | ||
| 228 | for (auto it = local_rtrees_[ls.zone].qbegin(bgeo::index::intersects(ls)); | ||
| 229 | it != local_rtrees_[ls.zone].qend(); it++) | ||
| 230 | { | ||
| 231 | if (bgeo::intersects(ls, std::get<1>(*it))) | ||
| 232 | { | ||
| 233 | co_yield std::get<2>(*it); | ||
| 234 | } | ||
| 235 | } | ||
| 236 | } | ||
| 237 | } | ||
| 238 | }; | ||
| 239 | |||
| 240 | } // namespace routemon::geo::multizonal | ||
diff --git a/server/src/geo/utm.cppm b/server/src/geo/utm.cppm new file mode 100644 index 0000000..83da222 --- /dev/null +++ b/server/src/geo/utm.cppm | |||
| @@ -0,0 +1,199 @@ | |||
| 1 | module; | ||
| 2 | |||
| 3 | #include <cassert> | ||
| 4 | #include <boost/geometry/geometries/geometries.hpp> | ||
| 5 | |||
| 6 | export module routemon:geo.utm; | ||
| 7 | |||
| 8 | import std; | ||
| 9 | import :geo; | ||
| 10 | import :geo.wgs84; | ||
| 11 | |||
| 12 | namespace routemon::geo::utm | ||
| 13 | { | ||
| 14 | |||
| 15 | class zonable_wgs84_point | ||
| 16 | { | ||
| 17 | wgs84::normalized_point p_; | ||
| 18 | |||
| 19 | explicit zonable_wgs84_point(wgs84::normalized_point p) | ||
| 20 | : p_{p} | ||
| 21 | {} | ||
| 22 | |||
| 23 | public: | ||
| 24 | zonable_wgs84_point() | ||
| 25 | : p_{} | ||
| 26 | { | ||
| 27 | } | ||
| 28 | |||
| 29 | static auto from(wgs84::normalized_point p) -> std::optional<zonable_wgs84_point> | ||
| 30 | { | ||
| 31 | if (p.lat() < -80 || p.lat() > 84) | ||
| 32 | return std::nullopt; | ||
| 33 | return zonable_wgs84_point{p}; | ||
| 34 | } | ||
| 35 | |||
| 36 | auto lat() const -> double | ||
| 37 | { | ||
| 38 | return p_.lat(); | ||
| 39 | } | ||
| 40 | |||
| 41 | auto lon() const -> double | ||
| 42 | { | ||
| 43 | return p_.lon(); | ||
| 44 | } | ||
| 45 | |||
| 46 | auto lon(double lon) -> void | ||
| 47 | { | ||
| 48 | p_.lon(lon); | ||
| 49 | } | ||
| 50 | |||
| 51 | auto lat(double lat) -> void | ||
| 52 | { | ||
| 53 | if (lat < -80 || lat > 84) | ||
| 54 | throw std::range_error{"latitude out of range for UTM"}; | ||
| 55 | p_.lat(lat); | ||
| 56 | } | ||
| 57 | }; | ||
| 58 | |||
| 59 | } | ||
| 60 | |||
| 61 | namespace boost::geometry::traits { | ||
| 62 | |||
| 63 | namespace { | ||
| 64 | |||
| 65 | using zonable_wgs84_point = routemon::geo::utm::zonable_wgs84_point; | ||
| 66 | |||
| 67 | } // namespace <anonymous> | ||
| 68 | |||
| 69 | template<> struct tag<zonable_wgs84_point> { using type = point_tag; }; | ||
| 70 | template<> struct dimension<zonable_wgs84_point> : boost::mpl::int_<2> {}; | ||
| 71 | template<> struct coordinate_type<zonable_wgs84_point> { using type = double; }; | ||
| 72 | template<> struct coordinate_system<zonable_wgs84_point> { using type = routemon::geo::wgs84::cs; }; | ||
| 73 | |||
| 74 | template<std::size_t index> | ||
| 75 | struct access<zonable_wgs84_point, index> { | ||
| 76 | static inline auto get(zonable_wgs84_point const& p) -> double | ||
| 77 | { | ||
| 78 | if constexpr (index == 0) | ||
| 79 | return p.lon(); | ||
| 80 | else if constexpr (index == 1) | ||
| 81 | return p.lat(); | ||
| 82 | else static_assert(false, "Out of range"); | ||
| 83 | } | ||
| 84 | |||
| 85 | static inline auto set(zonable_wgs84_point& p, double v) -> void | ||
| 86 | { | ||
| 87 | if constexpr (index == 0) | ||
| 88 | p.lon(v); | ||
| 89 | else if constexpr (index == 1) | ||
| 90 | p.lat(v); | ||
| 91 | else static_assert(false, "Out of range"); | ||
| 92 | } | ||
| 93 | }; | ||
| 94 | |||
| 95 | } // namespace boost::geometry::traits | ||
| 96 | |||
| 97 | namespace routemon::geo::utm | ||
| 98 | { | ||
| 99 | |||
| 100 | using zonable_wgs84_linestring = bgeo::model::linestring<zonable_wgs84_point, std::vector>; | ||
| 101 | |||
| 102 | // In the sense of the common WGS84 subdivisions by simple northing | ||
| 103 | // and easting (so no Norway/Svalbard exceptions). | ||
| 104 | class zone | ||
| 105 | { | ||
| 106 | // Note: calculations heavily depend on the values of the variants. | ||
| 107 | enum class hemisphere : std::uint8_t | ||
| 108 | { | ||
| 109 | northern = 0, | ||
| 110 | southern = 1, | ||
| 111 | }; | ||
| 112 | |||
| 113 | // Negative if in the southern hemisphere | ||
| 114 | // Equal to (zone_no - 1) * 2 + hemisphere | ||
| 115 | // In range [0, 119] | ||
| 116 | std::uint8_t zone_; | ||
| 117 | |||
| 118 | static constexpr auto zone_no_valid(std::uint8_t zone_no) -> bool | ||
| 119 | { | ||
| 120 | return 0 < zone_no && zone_no <= 60; | ||
| 121 | } | ||
| 122 | |||
| 123 | static constexpr auto from(std::uint8_t zone_no, hemisphere h) -> std::uint8_t | ||
| 124 | { | ||
| 125 | if (!zone_no_valid(zone_no)) | ||
| 126 | throw std::invalid_argument{"invalid UTM zone number"}; | ||
| 127 | return (zone_no - 1) * 2 + static_cast<std::uint8_t>(h); | ||
| 128 | } | ||
| 129 | |||
| 130 | public: | ||
| 131 | explicit constexpr zone(std::uint8_t zone_no, hemisphere h) | ||
| 132 | : zone_{from(zone_no, h)} | ||
| 133 | {} | ||
| 134 | |||
| 135 | auto hemisphere() const -> enum hemisphere | ||
| 136 | { | ||
| 137 | return static_cast<enum hemisphere>(zone_ % 2); | ||
| 138 | } | ||
| 139 | |||
| 140 | // Return value in range [1, 60] | ||
| 141 | auto zone_no() const -> std::uint8_t | ||
| 142 | { | ||
| 143 | return 1 + zone_ / 2; | ||
| 144 | } | ||
| 145 | |||
| 146 | // Return value in range [min(), max()] | ||
| 147 | constexpr auto as_index() const -> std::uint8_t | ||
| 148 | { | ||
| 149 | return zone_; | ||
| 150 | } | ||
| 151 | |||
| 152 | static auto for_wgs84_point(zonable_wgs84_point p) noexcept -> zone | ||
| 153 | { | ||
| 154 | auto zone_no = static_cast<std::uint8_t>(1 + (p.lon() + 180.0) / 6.0); | ||
| 155 | auto northern = p.lat() >= 0.0; | ||
| 156 | return zone{zone_no, northern ? hemisphere::northern : hemisphere::southern}; | ||
| 157 | } | ||
| 158 | |||
| 159 | auto wgs84_proj_epsg() const -> int | ||
| 160 | { | ||
| 161 | switch (hemisphere()) | ||
| 162 | { | ||
| 163 | case hemisphere::northern: | ||
| 164 | return 32600 + zone_no(); | ||
| 165 | case hemisphere::southern: | ||
| 166 | return 32700 + zone_no(); | ||
| 167 | } | ||
| 168 | } | ||
| 169 | |||
| 170 | // Guarantee: min().as_index() == 0. | ||
| 171 | static constexpr auto min() noexcept -> zone | ||
| 172 | { | ||
| 173 | return zone{1, hemisphere::northern}; | ||
| 174 | } | ||
| 175 | |||
| 176 | static constexpr auto max() noexcept -> zone | ||
| 177 | { | ||
| 178 | return zone{60, hemisphere::southern}; | ||
| 179 | } | ||
| 180 | |||
| 181 | auto next() const -> zone | ||
| 182 | { | ||
| 183 | assert(zone_ <= max().as_index()); | ||
| 184 | if (zone_ == max().as_index()) | ||
| 185 | throw std::range_error{"cannot take next of greatest UTM zone"}; | ||
| 186 | auto copy = zone{*this}; | ||
| 187 | copy.zone_++; | ||
| 188 | return copy; | ||
| 189 | } | ||
| 190 | |||
| 191 | auto operator==(zone rhs) const -> bool | ||
| 192 | { | ||
| 193 | return zone_ == rhs.zone_; | ||
| 194 | } | ||
| 195 | }; | ||
| 196 | static_assert(zone::min().as_index() == 0); | ||
| 197 | static_assert(zone::max().as_index() == 119); | ||
| 198 | |||
| 199 | } // namespace routemon::geo::utm | ||
diff --git a/server/src/geo/utm_zone_local.cppm b/server/src/geo/utm_zone_local.cppm new file mode 100644 index 0000000..98a1344 --- /dev/null +++ b/server/src/geo/utm_zone_local.cppm | |||
| @@ -0,0 +1,74 @@ | |||
| 1 | module; | ||
| 2 | |||
| 3 | #include <boost/geometry.hpp> | ||
| 4 | #include <boost/geometry/geometries/register/linestring.hpp> | ||
| 5 | |||
| 6 | export module routemon:geo.utm.zone_local; | ||
| 7 | |||
| 8 | import :geo; | ||
| 9 | import :geo.utm; | ||
| 10 | |||
| 11 | // Coordinates that are local to some known UTM zone. | ||
| 12 | namespace routemon::geo::utm::zone_local { | ||
| 13 | |||
| 14 | using cs = bgeo::cs::cartesian; | ||
| 15 | |||
| 16 | namespace prim { | ||
| 17 | |||
| 18 | using point = bgeo::model::d2::point_xy<double>; | ||
| 19 | using linestring = bgeo::model::linestring<point>; | ||
| 20 | using box = bgeo::model::box<point>; | ||
| 21 | |||
| 22 | } // namespace prim | ||
| 23 | |||
| 24 | struct point : public prim::point | ||
| 25 | { | ||
| 26 | zone zone; | ||
| 27 | |||
| 28 | explicit point(class zone zone, double x, double y) | ||
| 29 | : prim::point{x, y}, zone{zone} | ||
| 30 | { | ||
| 31 | } | ||
| 32 | }; | ||
| 33 | |||
| 34 | struct linestring : public prim::linestring | ||
| 35 | { | ||
| 36 | zone zone; | ||
| 37 | |||
| 38 | explicit linestring(class zone zone) | ||
| 39 | : zone{zone} | ||
| 40 | { | ||
| 41 | } | ||
| 42 | }; | ||
| 43 | |||
| 44 | } // namespace routemon::geo::utm::zone_local | ||
| 45 | |||
| 46 | BOOST_GEOMETRY_REGISTER_LINESTRING(routemon::geo::utm::zone_local::linestring) | ||
| 47 | |||
| 48 | namespace boost::geometry::traits { | ||
| 49 | |||
| 50 | namespace { | ||
| 51 | |||
| 52 | using zone_local_point = routemon::geo::utm::zone_local::point; | ||
| 53 | |||
| 54 | } // namespace <anonymous> | ||
| 55 | |||
| 56 | template<> struct tag<zone_local_point> { using type = point_tag; }; | ||
| 57 | template<> struct dimension<zone_local_point> : boost::mpl::int_<2> {}; | ||
| 58 | template<> struct coordinate_type<zone_local_point> { using type = double; }; | ||
| 59 | template<> struct coordinate_system<zone_local_point> { using type = routemon::geo::utm::zone_local::cs; }; | ||
| 60 | |||
| 61 | template<std::size_t index> | ||
| 62 | struct access<zone_local_point, index> { | ||
| 63 | static inline auto get(zone_local_point const& p) -> double | ||
| 64 | { | ||
| 65 | return bgeo::get<index>(static_cast<routemon::geo::utm::zone_local::prim::point const&>(p)); | ||
| 66 | } | ||
| 67 | |||
| 68 | static inline auto set(zone_local_point& p, double v) -> void | ||
| 69 | { | ||
| 70 | return bgeo::set<index>(static_cast<routemon::geo::utm::zone_local::prim::point&>(p), v); | ||
| 71 | } | ||
| 72 | }; | ||
| 73 | |||
| 74 | } // namespace boost::geometry::traits | ||
diff --git a/server/src/geo/wgs84.cppm b/server/src/geo/wgs84.cppm new file mode 100644 index 0000000..d574d39 --- /dev/null +++ b/server/src/geo/wgs84.cppm | |||
| @@ -0,0 +1,115 @@ | |||
| 1 | module; | ||
| 2 | |||
| 3 | #include <boost/geometry.hpp> | ||
| 4 | |||
| 5 | export module routemon:geo.wgs84; | ||
| 6 | |||
| 7 | import :geo; | ||
| 8 | |||
| 9 | namespace routemon::geo::wgs84 { | ||
| 10 | |||
| 11 | using cs = bgeo::cs::geographic<bgeo::degree>; | ||
| 12 | using point = bgeo::model::point<double, 2, cs>; | ||
| 13 | // TODO: guarantee that linestring provides non-static member | ||
| 14 | // auto reserve(std::size_t) -> void | ||
| 15 | // Perhaps better yet: guarantee that the backing container is a std::vector. | ||
| 16 | using linestring = bgeo::model::linestring<point, std::vector /* the default */>; | ||
| 17 | using box = bgeo::model::box<point>; | ||
| 18 | using stype = bgeo::srs::spheroid<double>; | ||
| 19 | using vincenty_strategy = bgeo::strategy::distance::vincenty<stype>; | ||
| 20 | |||
| 21 | auto from_lat_lon(double lat, double lon) -> point | ||
| 22 | { | ||
| 23 | return point{lon, lat}; | ||
| 24 | } | ||
| 25 | |||
| 26 | auto lat(point p) -> double | ||
| 27 | { | ||
| 28 | return bgeo::get<1>(p); | ||
| 29 | } | ||
| 30 | |||
| 31 | auto lon(point p) -> double | ||
| 32 | { | ||
| 33 | return bgeo::get<0>(p); | ||
| 34 | } | ||
| 35 | |||
| 36 | // Point p with latitude in [-90, 90] and longitude in [-180, 180) | ||
| 37 | auto normalize(point p) -> point | ||
| 38 | { | ||
| 39 | // Equivalent degrees in the range [0, 360) | ||
| 40 | auto nonneg_mod360 = [](double t) -> double | ||
| 41 | { return std::remainder(std::remainder(t, 360.0) + 360.0, 360.0); }; | ||
| 42 | |||
| 43 | auto p_lon = lon(p); | ||
| 44 | auto p_lat = nonneg_mod360(lat(p)) - 90.0; | ||
| 45 | if (90.0 < p_lat) | ||
| 46 | { | ||
| 47 | assert(p_lat < 270.0); // by nonneg_mod360 | ||
| 48 | p_lon += 180.0; | ||
| 49 | p_lat -= 180.0; | ||
| 50 | } | ||
| 51 | p_lon = nonneg_mod360(p_lon + 180.0) - 180.0; | ||
| 52 | |||
| 53 | return from_lat_lon(p_lat, p_lon); | ||
| 54 | } | ||
| 55 | |||
| 56 | auto is_normalized(point p) -> bool | ||
| 57 | { | ||
| 58 | return -90 <= lat(p) && lat(p) <= 90 && -180 <= lon(p) && lon(p) < 180; | ||
| 59 | } | ||
| 60 | |||
| 61 | class normalized_point | ||
| 62 | { | ||
| 63 | point p_; | ||
| 64 | |||
| 65 | public: | ||
| 66 | normalized_point(point p) | ||
| 67 | : p_{is_normalized(p) ? p : normalize(p)} | ||
| 68 | {} | ||
| 69 | |||
| 70 | normalized_point() | ||
| 71 | : p_{} | ||
| 72 | {} | ||
| 73 | |||
| 74 | auto lat() const -> double | ||
| 75 | { | ||
| 76 | return bgeo::get<1>(p_); | ||
| 77 | } | ||
| 78 | |||
| 79 | auto lon() const -> double | ||
| 80 | { | ||
| 81 | return bgeo::get<0>(p_); | ||
| 82 | } | ||
| 83 | |||
| 84 | auto lat(double lat) -> void | ||
| 85 | { | ||
| 86 | if (lat < -90 || lat > 90) | ||
| 87 | throw std::range_error{"latitude out of range"}; | ||
| 88 | bgeo::set<1>(p_, lat); | ||
| 89 | } | ||
| 90 | |||
| 91 | auto lon(double lon) -> void | ||
| 92 | { | ||
| 93 | if (lon < -180 || lon > 180) | ||
| 94 | throw std::range_error{"longitude out of range"}; | ||
| 95 | bgeo::set<0>(p_, lon); | ||
| 96 | } | ||
| 97 | }; | ||
| 98 | |||
| 99 | class normalized_linestring | ||
| 100 | { | ||
| 101 | std::vector<normalized_point> ls_; | ||
| 102 | |||
| 103 | public: | ||
| 104 | auto empty() const -> bool | ||
| 105 | { | ||
| 106 | return ls_.empty(); | ||
| 107 | } | ||
| 108 | |||
| 109 | auto push_back(normalized_point p) -> void | ||
| 110 | { | ||
| 111 | ls_.push_back(p); | ||
| 112 | } | ||
| 113 | }; | ||
| 114 | |||
| 115 | } // namespace routemon::geo::wgs84 | ||
diff --git a/server/src/gpx.cpp b/server/src/gpx.cpp index 078248f..2193e6a 100644 --- a/server/src/gpx.cpp +++ b/server/src/gpx.cpp | |||
| @@ -22,7 +22,7 @@ auto qname(std::string_view local) -> xml::qname_view | |||
| 22 | } | 22 | } |
| 23 | 23 | ||
| 24 | auto parse_wpt(xml::executor_ref e, xml::attribute_view attrs) | 24 | auto parse_wpt(xml::executor_ref e, xml::attribute_view attrs) |
| 25 | -> xml::parser<geo::point> | 25 | -> xml::parser<geo::utm::zonable_wgs84_point> |
| 26 | { | 26 | { |
| 27 | auto parse_xml_double = [](std::string_view sv) -> std::optional<double> | 27 | auto parse_xml_double = [](std::string_view sv) -> std::optional<double> |
| 28 | { return util::parse_double(sv, std::chars_format::fixed); }; | 28 | { return util::parse_double(sv, std::chars_format::fixed); }; |
| @@ -43,8 +43,12 @@ auto parse_wpt(xml::executor_ref e, xml::attribute_view attrs) | |||
| 43 | throw std::runtime_error{ | 43 | throw std::runtime_error{ |
| 44 | "expected valid latitude and longitude for waypoint" | 44 | "expected valid latitude and longitude for waypoint" |
| 45 | }; | 45 | }; |
| 46 | auto mp = geo::utm::zonable_wgs84_point::from( | ||
| 47 | geo::wgs84::normalized_point{geo::wgs84::from_lat_lon(*mlat, *mlon)}); | ||
| 48 | if (!mp) | ||
| 49 | throw std::runtime_error{"expected waypoint to be within UTM range"}; | ||
| 46 | co_await xml::ignore_contents(e); | 50 | co_await xml::ignore_contents(e); |
| 47 | co_return geo::point{*mlon, *mlat}; | 51 | co_return *mp; |
| 48 | } | 52 | } |
| 49 | 53 | ||
| 50 | auto parse_trkseg(xml::executor_ref e, xml::attribute_view) | 54 | auto parse_trkseg(xml::executor_ref e, xml::attribute_view) |
| @@ -121,7 +125,7 @@ auto parse_metadata(xml::executor_ref e, xml::attribute_view) | |||
| 121 | } | 125 | } |
| 122 | 126 | ||
| 123 | auto parse_wpt(xml::executor_ref e, xml::attribute_view attrs) | 127 | auto parse_wpt(xml::executor_ref e, xml::attribute_view attrs) |
| 124 | -> xml::parser<geo::point> | 128 | -> xml::parser<geo::utm::zonable_wgs84_point> |
| 125 | { | 129 | { |
| 126 | auto parse_xml_double = [](std::string_view sv) -> std::optional<double> | 130 | auto parse_xml_double = [](std::string_view sv) -> std::optional<double> |
| 127 | { return util::parse_double(sv, std::chars_format::fixed); }; | 131 | { return util::parse_double(sv, std::chars_format::fixed); }; |
| @@ -142,8 +146,12 @@ auto parse_wpt(xml::executor_ref e, xml::attribute_view attrs) | |||
| 142 | throw std::runtime_error{ | 146 | throw std::runtime_error{ |
| 143 | "expected valid latitude and longitude for waypoint" | 147 | "expected valid latitude and longitude for waypoint" |
| 144 | }; | 148 | }; |
| 149 | auto mp = geo::utm::zonable_wgs84_point::from( | ||
| 150 | geo::wgs84::normalized_point{geo::wgs84::from_lat_lon(*mlat, *mlon)}); | ||
| 151 | if (!mp) | ||
| 152 | throw std::runtime_error{"expected waypoint to be within UTM range"}; | ||
| 145 | co_await xml::ignore_contents(e); | 153 | co_await xml::ignore_contents(e); |
| 146 | co_return geo::point{*mlon, *mlat}; | 154 | co_return *mp; |
| 147 | } | 155 | } |
| 148 | 156 | ||
| 149 | auto parse_trkseg(xml::executor_ref e, xml::attribute_view) | 157 | auto parse_trkseg(xml::executor_ref e, xml::attribute_view) |
diff --git a/server/src/gpx.cppm b/server/src/gpx.cppm index 5872805..6267443 100644 --- a/server/src/gpx.cppm +++ b/server/src/gpx.cppm | |||
| @@ -1,7 +1,7 @@ | |||
| 1 | export module routemon:gpx; | 1 | export module routemon:gpx; |
| 2 | 2 | ||
| 3 | import std; | 3 | import std; |
| 4 | import :geo; | 4 | import :geo.utm; |
| 5 | import :util; | 5 | import :util; |
| 6 | import :xml; | 6 | import :xml; |
| 7 | 7 | ||
| @@ -15,7 +15,7 @@ struct metadata | |||
| 15 | 15 | ||
| 16 | struct track_segment | 16 | struct track_segment |
| 17 | { | 17 | { |
| 18 | geo::linestring waypoints; | 18 | geo::utm::zonable_wgs84_linestring waypoints; |
| 19 | }; | 19 | }; |
| 20 | 20 | ||
| 21 | struct track | 21 | struct track |
diff --git a/server/src/log.cpp b/server/src/log.cpp index 7fe725f..bb9fde2 100644 --- a/server/src/log.cpp +++ b/server/src/log.cpp | |||
| @@ -39,7 +39,7 @@ auto sink::write(tmp_message msg) -> void | |||
| 39 | sos << "] " << msg.txt; | 39 | sos << "] " << msg.txt; |
| 40 | for (auto const& [k, v] : msg.attrs) | 40 | for (auto const& [k, v] : msg.attrs) |
| 41 | sos << " " << k << "=" << std::quoted(v); | 41 | sos << " " << k << "=" << std::quoted(v); |
| 42 | sos << '\n'; | 42 | sos << std::endl; // we want the flush! |
| 43 | } | 43 | } |
| 44 | 44 | ||
| 45 | auto make_sink(level lvl) -> std::shared_ptr<sink> | 45 | auto make_sink(level lvl) -> std::shared_ptr<sink> |