mpadge.tngl.sh
Distances on Directed Graphs in R
1#include "dodgr-to-sf.h"
2
3//' Make unordered_set of all new edge names
4//' @noRd
5size_t dodgr_sf::make_edge_name_set (std::unordered_set <std::string> &new_edge_name_set,
6 const Rcpp::CharacterVector &new_edges)
7{
8 new_edge_name_set.clear ();
9 // Rcpp::CharacterVector requires long int index type
10 for (long int i = 0; i < new_edges.size (); i++)
11 {
12 new_edge_name_set.emplace (static_cast <std::string> (new_edges [i]));
13 }
14 return new_edge_name_set.size ();
15}
16
17//' Make vector of all new edge names
18//' @noRd
19void dodgr_sf::make_edge_name_vec (const size_t n,
20 const Rcpp::CharacterVector &new_edges,
21 std::vector <std::string> &new_edge_name_vec)
22{
23 new_edge_name_vec.clear ();
24 new_edge_name_vec.resize (n);
25 new_edge_name_vec [0] = static_cast <std::string> (new_edges [0]);
26 size_t count = 0;
27 for (long int i = 1; i < new_edges.size (); i++)
28 {
29 std::string new_edge_i = static_cast <std::string> (new_edges [i]);
30 if (new_edge_i != new_edge_name_vec [count])
31 new_edge_name_vec [++count] = new_edge_i;
32 }
33}
34
35
36//' Get numbers of old edges corresponding to each contracted edge
37//' @noRd
38size_t dodgr_sf::get_edgevec_sizes (const size_t nedges,
39 const Rcpp::CharacterVector &new_edges,
40 std::vector <size_t> &edgevec_sizes)
41{
42 edgevec_sizes.clear ();
43 edgevec_sizes.resize (nedges);
44 size_t count = 1, edgenum = 0;
45 for (long int i = 1; i < new_edges.size (); i++)
46 {
47 if (new_edges [i] == new_edges [i - 1])
48 count++;
49 else
50 {
51 edgevec_sizes [edgenum++] = count;
52 count = 1;
53 }
54 }
55 // And last count too:
56 edgevec_sizes [edgenum++] = count;
57 return edgenum;
58}
59
60//' Collect sets of all from and to vertices for each set of edges corresponding
61//' to each contracted edge. These sets aren't in any particular order, but the
62//' two sets may be used to match the from and to vertices. These sets are then
63//' arranged into sequences in the subsequent function,
64//' \code{order_vert_sequences}.
65//' @noRd
66void dodgr_sf::get_edge_to_vert_maps (const std::vector <size_t> &edgevec_sizes,
67 const Rcpp::DataFrame &graph_full,
68 const Rcpp::CharacterVector &old_edges,
69 const Rcpp::CharacterVector &new_edges,
70 const std::vector <std::string> &new_edge_names,
71 std::unordered_map <std::string,
72 std::vector <std::string> > &full_from_edge_map,
73 std::unordered_map <std::string,
74 std::vector <std::string> > &full_to_edge_map)
75{
76 Rcpp::CharacterVector idf_r = graph_full ["from_id"],
77 idt_r = graph_full ["to_id"];
78
79 size_t count = 1, edgenum = 0;
80 std::vector <std::string> from_node, to_node;
81 from_node.resize (edgevec_sizes [edgenum]);
82 to_node.resize (edgevec_sizes [edgenum]);
83 // old_edges is 1-indexed!
84 long int the_edge = static_cast <long int> (atoi (old_edges [0])) - 1;
85 from_node [0] = static_cast <std::string> (idf_r [the_edge]);
86 to_node [0] = static_cast <std::string> (idt_r [the_edge]);
87 for (long int i = 1; i < new_edges.size (); i++)
88 {
89 if (new_edges [i] != new_edges [i - 1])
90 {
91 full_from_edge_map.emplace (new_edge_names [edgenum], from_node);
92 full_to_edge_map.emplace (new_edge_names [edgenum], to_node);
93 from_node.clear ();
94 to_node.clear ();
95 edgenum++;
96 from_node.resize (edgevec_sizes [edgenum]);
97 to_node.resize (edgevec_sizes [edgenum]);
98 count = 0;
99 }
100 the_edge = atoi (old_edges [i]) - 1; // it's 1-indexed!
101 from_node [count] = static_cast <std::string> (idf_r [the_edge]);
102 to_node [count++] = static_cast <std::string> (idt_r [the_edge]);
103 }
104 full_from_edge_map.emplace (new_edge_names [edgenum], from_node);
105 full_to_edge_map.emplace (new_edge_names [edgenum], to_node);
106
107 from_node.clear ();
108 to_node.clear ();
109}
110
111void dodgr_sf::order_vert_sequences (Rcpp::List &edge_sequences,
112 std::vector <std::string> &new_edge_names,
113 std::unordered_map <std::string,
114 std::vector <std::string> > &full_from_edge_map,
115 std::unordered_map <std::string,
116 std::vector <std::string> > &full_to_edge_map)
117{
118 const size_t nedges = static_cast <size_t> (edge_sequences.size ());
119 for (size_t i = 0; i < nedges; i++)
120 {
121 Rcpp::checkUserInterrupt ();
122 std::map <std::string, std::string> idmap, idmap_rev;
123 std::string this_edge = new_edge_names [i];
124 std::vector <std::string> from_edges = full_from_edge_map [this_edge],
125 to_edges = full_to_edge_map [this_edge];
126 for (size_t j = 0; j < from_edges.size (); j++)
127 {
128 idmap.emplace (from_edges [j], to_edges [j]);
129 idmap_rev.emplace (to_edges [j], from_edges [j]);
130 }
131 size_t nnodes = idmap.size ();
132
133 // Find the front of the sequence of which the first idmap pair is
134 // part by stepping backwards through idmap_rev
135 std::string front_node = idmap.begin ()->first;
136 std::string back_node = idmap.begin ()->second;
137 while (idmap_rev.find (front_node) != idmap_rev.end ())
138 front_node = idmap_rev.find (front_node)->second;
139 back_node = idmap.find (front_node)->second;
140
141 std::vector <std::string> id (nnodes + 1);
142 id [0] = front_node;
143 id [1] = back_node;
144 size_t count = 2;
145 idmap.erase (front_node);
146
147 while (idmap.size () > 0)
148 {
149 std::map <std::string, std::string>::iterator it =
150 idmap.find (back_node);
151 back_node = it->second;
152 id [count++] = back_node;
153 idmap.erase (it);
154 }
155 idmap_rev.clear ();
156
157 edge_sequences [static_cast <long int> (i)] = id;
158 }
159}
160
161// Count number of non-contracted edges in the contracted graph (non-contracted
162// because they can not be simplified).
163size_t dodgr_sf::count_non_contracted_edges (const Rcpp::CharacterVector &contr_edges,
164 std::unordered_set <std::string> &new_edge_name_set)
165{
166 size_t edge_count = 0;
167 for (long int i = 0; i < contr_edges.size (); i++)
168 {
169 if (new_edge_name_set.find (static_cast <std::string>
170 (contr_edges [i])) == new_edge_name_set.end ())
171 {
172 edge_count++;
173 }
174 }
175 return edge_count;
176}
177
178// Append non-contracted edges to contracted edges in all edge-to-vertex maps
179void dodgr_sf::append_nc_edges (const size_t nc_edge_count,
180 const Rcpp::DataFrame &graph_contr,
181 std::unordered_set <std::string> &new_edge_name_set,
182 std::vector <std::string> &new_edge_name_vec,
183 const Rcpp::List &edge_sequences_contr,
184 std::vector <std::string> &all_edge_names,
185 Rcpp::List &edge_sequences_all)
186{
187 Rcpp::List edge_sequences_new (nc_edge_count);
188 std::vector <std::string> old_edge_names (nc_edge_count);
189 size_t count = 0;
190 Rcpp::CharacterVector idf_r_c = graph_contr ["from_id"],
191 idt_r_c = graph_contr ["to_id"],
192 contr_edges = graph_contr ["edge_id"];
193 for (long int i = 0; i < graph_contr.nrow (); i++)
194 {
195 if (new_edge_name_set.find (static_cast <std::string>
196 (contr_edges [i])) == new_edge_name_set.end ())
197 {
198 old_edge_names [count] = contr_edges [i];
199 Rcpp::CharacterVector idvec (2);
200 idvec [0] = idf_r_c [i];
201 idvec [1] = idt_r_c [i];
202 edge_sequences_new [static_cast <long int> (count++)] = idvec;
203 }
204 }
205 // Then just join the two edge_sequence Lists together, along with vectors
206 // of edge names
207 const size_t total_edges = static_cast <size_t> (edge_sequences_contr.size ()) +
208 nc_edge_count;
209 all_edge_names.resize (total_edges);
210 // These two have different indexing types:
211 for (size_t i = 0; i < static_cast <size_t> (edge_sequences_contr.size ()); i++)
212 all_edge_names [i] = new_edge_name_vec [i];
213 for (long int i = 0; i < edge_sequences_contr.size (); i++)
214 edge_sequences_all [i] = edge_sequences_contr [i];
215 for (size_t i = 0; i < static_cast <size_t> (edge_sequences_new.size ()); i++)
216 all_edge_names [static_cast <size_t> (edge_sequences_contr.size ()) + i] =
217 old_edge_names [i];
218 for (long int i = 0; i < edge_sequences_new.size (); i++)
219 edge_sequences_all [edge_sequences_contr.size () + i] = edge_sequences_new [i];
220}
221
222// from osmdata/src/get-bbox.cpp
223Rcpp::NumericVector rcpp_get_bbox_sf (double xmin, double xmax, double ymin, double ymax)
224{
225 std::vector <std::string> names;
226 names.push_back ("xmin");
227 names.push_back ("ymin");
228 names.push_back ("xmax");
229 names.push_back ("ymax");
230
231 Rcpp::NumericVector bbox (4, NA_REAL);
232 bbox (0) = xmin;
233 bbox (1) = xmax;
234 bbox (2) = ymin;
235 bbox (3) = ymax;
236
237 bbox.attr ("names") = names;
238 bbox.attr ("class") = "bbox";
239
240 return bbox;
241}
242
243// Convert the Rcpp::List of vertex names to corresponding coordinates,
244// construct matrices of sequential coordinates for each contracted edge, and
245// convert these to sf geometry objects
246void dodgr_sf::xy_to_sf (const Rcpp::DataFrame &graph_full,
247 const Rcpp::List &edge_sequences,
248 const std::vector <std::string> &all_edge_names,
249 Rcpp::List &res)
250{
251 const size_t total_edges = static_cast <size_t> (res.size ());
252
253 Rcpp::CharacterVector idf_r = graph_full ["from_id"],
254 idt_r = graph_full ["to_id"];
255
256 Rcpp::NumericVector xf = graph_full ["from_lon"],
257 yf = graph_full ["from_lat"],
258 xt = graph_full ["to_lon"],
259 yt = graph_full ["to_lat"];
260
261 std::unordered_map <std::string, size_t> edge_num_map;
262 for (long int i = 0; i < idf_r.size (); i++)
263 {
264 std::string idft = static_cast <std::string> (idf_r [i]) + "-" +
265 static_cast <std::string> (idt_r [i]);
266 edge_num_map.emplace (idft, i);
267 }
268
269 double xmin = INFINITE_DOUBLE, xmax = -INFINITE_DOUBLE,
270 ymin = INFINITE_DOUBLE, ymax = -INFINITE_DOUBLE;
271 for (size_t i = 0; i < total_edges; i++)
272 {
273 Rcpp::CharacterVector idvec = edge_sequences [static_cast <long int> (i)];
274 Rcpp::NumericVector x (idvec.size ()), y (idvec.size ());
275 // Fill first edge
276 std::string id0 = static_cast <std::string> (idvec [0]),
277 id1 = static_cast <std::string> (idvec [1]);
278 std::string id01 = id0 + "-" + id1;
279 long int j = static_cast <long int> (edge_num_map.find (id01)->second);
280 x [0] = xf [j];
281 y [0] = yf [j];
282 x [1] = xt [j];
283 y [1] = yt [j];
284 // Then just remaining "to" values
285 idvec.erase (0);
286 idvec.erase (0);
287 long int count = 2;
288 while (idvec.size () > 0)
289 {
290 id0 = id1;
291 id1 = static_cast <std::string> (idvec [0]);
292 id01 = id0 + "-" + id1;
293 j = static_cast <long int> (edge_num_map.find (id01)->second);
294 x [count] = xt [j];
295 y [count] = yt [j];
296 count++;
297 idvec.erase (0);
298 }
299 Rcpp::NumericMatrix mat (static_cast <const int> (x.size ()), 2);
300 mat (Rcpp::_, 0) = x;
301 mat (Rcpp::_, 1) = y;
302 xmin = std::min (xmin, static_cast <double> (Rcpp::min (x)));
303 xmax = std::max (xmax, static_cast <double> (Rcpp::max (x)));
304 ymin = std::min (ymin, static_cast <double> (Rcpp::min (y)));
305 ymax = std::max (ymax, static_cast <double> (Rcpp::max (y)));
306
307 // Then some sf necessities:
308 mat.attr ("class") =
309 Rcpp::CharacterVector::create ("XY", "LINESTRING", "sfg");
310 res (i) = mat;
311 }
312
313 res.attr ("names") = all_edge_names;
314 res.attr ("n_empty") = 0;
315 res.attr ("class") = Rcpp::CharacterVector::create ("sfc_LINESTRING", "sfc");
316 res.attr ("precision") = 0.0;
317 res.attr ("bbox") = rcpp_get_bbox_sf (xmin, ymin, xmax, ymax);
318 Rcpp::List crs = Rcpp::List::create (static_cast <int> (4326), osm_p4s);
319 crs.attr ("names") = Rcpp::CharacterVector::create ("epsg", "proj4string");
320 crs.attr ("class") = "crs";
321 res.attr ("crs") = crs;
322}
323
324//' rcpp_aggregate_to_sf
325//'
326//' Aggregate a dodgr network data.frame to an sf LINESTRING data.frame
327//'
328//' @param graph_full Rcpp::DataFrame containing the **full** graph
329//' @param graph_contr Rcpp::DataFrame containing the **contracted** graph
330//' @param edge_map Rcpp::DataFrame containing the edge map returned from
331//' \code{dodgr_contract_graph}
332//'
333//' @return Rcpp::List object of `sf::LINESTRING` geoms
334//'
335//' @noRd
336// [[Rcpp::export]]
337Rcpp::List rcpp_aggregate_to_sf (const Rcpp::DataFrame &graph_full,
338 const Rcpp::DataFrame &graph_contr, const Rcpp::DataFrame &edge_map)
339{
340 Rcpp::CharacterVector old_edges, new_edges,
341 contr_edges = graph_contr ["edge_id"];
342
343 if (edge_map.nrow () > 0) {
344 old_edges = edge_map ["edge_old"];
345 new_edges = edge_map ["edge_new"];
346 } else {
347 // If a previously contracted graph is submitted, the edge_map will be
348 // empty and graph_full == graph_contr, so just use the latter only:
349 old_edges = contr_edges;
350 new_edges = contr_edges;
351 }
352
353 // store vector of names used to name the resultant sf-objects. A
354 // corresponding set is also used below to insert non-contracted edges in
355 // final result
356 std::vector <std::string> new_edge_names;
357 std::unordered_set <std::string> new_edge_set;
358 const size_t nedges = dodgr_sf::make_edge_name_set (new_edge_set, new_edges);
359 dodgr_sf::make_edge_name_vec (nedges, new_edges, new_edge_names);
360
361 // Then get numbers of original edges for each contracted edge
362 std::vector <size_t> edgevec_sizes;
363 size_t check = dodgr_sf::get_edgevec_sizes (nedges, new_edges, edgevec_sizes);
364 if (check != nedges)
365 Rcpp::stop ("number of new edges in contracted graph not the right size"); // # nocov
366
367 /* Map contracted edge ids onto corresponding pairs of from and to vertices.
368 * The first vertex of a sequence won't be necessarily at the start of a
369 * sequence, so two maps are made, one from each node to the next, and the
370 * other holding the same values in reverse. The latter enables sequences to
371 * be traced in reverse by matching end points.
372 */
373 std::unordered_map <std::string,
374 std::vector <std::string> > full_from_edge_map, full_to_edge_map;
375 dodgr_sf::get_edge_to_vert_maps (edgevec_sizes, graph_full, old_edges, new_edges,
376 new_edge_names, full_from_edge_map, full_to_edge_map);
377 edgevec_sizes.clear ();
378
379 Rcpp::List edge_sequences (nedges); // holds final sequences
380 // The main work done in this whole file:
381 dodgr_sf::order_vert_sequences (edge_sequences, new_edge_names, full_from_edge_map,
382 full_to_edge_map);
383 full_from_edge_map.clear ();
384 full_to_edge_map.clear ();
385
386 // Then append the non-contracted edges that are in the contracted graph
387 size_t nc_edge_count = dodgr_sf::count_non_contracted_edges (contr_edges,
388 new_edge_set);
389 std::vector <std::string> all_edge_names;
390 const size_t total_edges = nedges + nc_edge_count;
391 Rcpp::List edge_sequences_all (total_edges);
392 dodgr_sf::append_nc_edges (nc_edge_count, graph_contr,
393 new_edge_set, new_edge_names, edge_sequences,
394 all_edge_names, edge_sequences_all);
395 new_edge_names.clear ();
396 new_edge_set.clear ();
397
398 Rcpp::List res (total_edges);
399 dodgr_sf::xy_to_sf (graph_full, edge_sequences_all, all_edge_names, res);
400
401 return res;
402}