mpadge.tngl.sh
Distances on Directed Graphs in R
1#include "turn_penalty.h"
2
3/* This code adds time penalties for turning across oncoming traffic in
4 * insersections. Works by sorting all outgoing edges in a clockwise direction,
5 * and applying penalties to edges with sorted indices >= 2. (Right-side travel
6 * simply reverses indices to effectively be anti-clockwise. Incoming edges are
7 * also sorted, but sort order is not used.) Implements the following steps:
8 *
9 * 1. The `fill_edges` function creates unordered_map between centre vertices of
10 * junctions and a `std::pair <RTEdgeSet, RTEdgeSet>` of (incoming, outgoing)
11 * vertices, where `RTEdgeSet = std::set <OneEdge, clockwise_sort>`. This
12 * function also fills an unordered_set of `junction_vertices`.
13 * 2. The `replace_junctions` function scans through that map and creates new
14 * edges connecting the incoming vertex to the outgoing vertex in each
15 * direction, and applies penalties to all edges with sorted order > 2.
16 */
17
18
19void routetimes::fill_edges (const Rcpp::DataFrame &graph,
20 std::unordered_map <std::string,
21 std::pair <RTEdgeSet, RTEdgeSet> > &the_edges,
22 std::unordered_set <std::string> &junction_vertices)
23{
24 std::vector <std::string> vx0 = graph [".vx0"],
25 vx1 = graph [".vx1"],
26 edge_ = graph ["edge_"];
27 std::vector <double> vx0_x = graph [".vx0_x"],
28 vx0_y = graph [".vx0_y"],
29 vx1_x = graph [".vx1_x"],
30 vx1_y = graph [".vx1_y"];
31
32 const size_t n = static_cast <size_t> (graph.nrow ());
33
34 for (size_t i = 0; i < n; i++)
35 {
36 OneEdge edge;
37 edge.v0 = vx0 [i];
38 edge.v1 = vx1 [i];
39 edge.edge = edge_ [i];
40 edge.x = vx1_x [i] - vx0_x [i];
41 edge.y = vx1_y [i] - vx0_y [i];
42
43 // the incoming edge: vx0->vx1 with key = vx1
44 routetimes::replace_one_map_edge (the_edges, vx1 [i], edge, true);
45 // the outgoing edge: vx0->vx1 with key = vx0
46 routetimes::replace_one_map_edge (the_edges, vx0 [i], edge, false);
47 }
48
49 erase_non_junctions (the_edges, junction_vertices);
50}
51
52void routetimes::replace_one_map_edge (
53 std::unordered_map <std::string,
54 std::pair <RTEdgeSet, RTEdgeSet> > &the_edges,
55 std::string key, OneEdge edge, bool incoming)
56{
57 std::pair <RTEdgeSet, RTEdgeSet> edge_set;
58 if (the_edges.find (key) != the_edges.end ())
59 {
60 edge_set = the_edges.at (key);
61 the_edges.erase (key);
62 }
63 if (incoming)
64 edge_set.first.emplace (edge);
65 else
66 edge_set.second.emplace (edge);
67
68 the_edges.emplace (key, edge_set);
69}
70
71/* Remove all edge items with < 3 outgoing edges. A junction has 2 or more
72 * out_edges / in_edges, but presume that flow in and out of 1->2 junctions is
73 * regulated by lights, and that only proper cross-junctions (1->3) incur
74 * additional waiting times. The result: min_edges = 4
75 */
76void routetimes::erase_non_junctions (
77 std::unordered_map <std::string,
78 std::pair <RTEdgeSet, RTEdgeSet> > &the_edges,
79 std::unordered_set <std::string> &junction_vertices)
80{
81 const int min_edges = 4;
82
83 std::unordered_set <std::string> removes;
84 for (auto e: the_edges)
85 {
86 if (e.second.second.size () < min_edges) // set of outgoing edges
87 removes.emplace (e.first);
88 else
89 {
90 for (auto r: e.second.second) // RTEdgeSet
91 junction_vertices.emplace (r.v0);
92 }
93 }
94 for (auto r: removes)
95 the_edges.erase (r);
96}
97
98/* Replace junction mid-points with direct connections from in -> 0 -> out with
99 * (in, out) pairs and a binary flag for turning penalty.
100 */
101void routetimes::replace_junctions (
102 const std::unordered_map <std::string,
103 std::pair <RTEdgeSet, RTEdgeSet> > &the_edges,
104 std::vector <OneCompoundEdge> &junctions,
105 bool left_side)
106{
107 // Estimate size of resultant vector - this is not simply n * (n - 1),
108 // because one way streeets mean in-edges aren't always the same as out, so
109 // numbers have to be explicitly counted.
110 size_t n_junctions = 0;
111 for (auto e: the_edges)
112 for (auto ei: e.second.first) // incoming edge
113 {
114 std::unordered_set <std::string> out_edges;
115 for (auto ej: e.second.second) // outgoing edge
116 if (ej.edge != ei.edge)
117 out_edges.emplace (ej.edge);
118
119 for (auto ej: e.second.second)
120 if (out_edges.find (ej.edge) != out_edges.end ())
121 n_junctions++;
122 }
123 junctions.resize (n_junctions);
124
125 size_t count = 0;
126 for (auto e: the_edges)
127 {
128 // iterate over each incoming edge, and establish penalties for all
129 // other ongoing.
130 for (auto ei: e.second.first) // incoming edge
131 {
132 std::unordered_map <std::string, size_t> out_edges;
133 size_t i = 0;
134 for (auto ej: e.second.second) // outgoing edge
135 if (ej.edge != ei.edge)
136 {
137 out_edges.emplace (ej.edge, i++);
138 }
139 // For right-side travel, simply reverse the clockwise sequence of
140 // values, to give anti-clockwise sequential counts:
141 if (!left_side)
142 for (auto oe: out_edges)
143 {
144 out_edges [oe.first] = out_edges.size () - oe.second;
145 }
146
147 // Then use those out_edges to construct new edge lists:
148 for (auto ej: e.second.second) // outgoing edge
149 if (out_edges.find (ej.edge) != out_edges.end ())
150 {
151 OneCompoundEdge edge;
152 edge.v0 = ei.v0;
153 edge.v1 = ej.v1;
154 edge.edge0 = ei.edge;
155 edge.edge1 = ej.edge;
156 edge.penalty = false;
157 if (out_edges.at (ej.edge) > 2)
158 edge.penalty = true;
159 junctions [count++] = edge;
160 }
161 }
162 }
163}
164
165// This is hard-coded for weight_streetnet.sc structure
166Rcpp::DataFrame routetimes::expand_edges (const Rcpp::DataFrame &graph,
167 std::vector <OneCompoundEdge> &junctions, int turn_penalty)
168{
169 const size_t hash_len = 10; // for new edge IDs
170
171 std::unordered_map <std::string, double> map_x, map_y, map_d, map_dw,
172 map_time, map_timew;
173 std::unordered_map <std::string, std::string> map_object, map_highway;
174
175 std::vector <std::string> vx0 = graph [".vx0"],
176 vx1 = graph [".vx1"],
177 edge_ = graph ["edge_"],
178 object = graph ["object_"],
179 highway = graph ["highway"];
180 std::vector <double> vx0_x = graph [".vx0_x"],
181 vx0_y = graph [".vx0_y"],
182 vx1_x = graph [".vx1_x"],
183 vx1_y = graph [".vx1_y"],
184 d = graph ["d"],
185 dw = graph ["d_weighted"],
186 time = graph ["time"],
187 timew = graph ["time_weighted"];
188
189 std::unordered_set <std::string> edge_set;
190 for (auto j: junctions)
191 {
192 edge_set.emplace (j.edge0);
193 edge_set.emplace (j.edge1);
194 }
195
196 for (size_t i = 0; i < static_cast <size_t> (graph.nrow ()); i++)
197 {
198 if (edge_set.find (edge_ [i]) != edge_set.end ())
199 {
200 map_x.emplace (vx0 [i], vx0_x [i]);
201 map_y.emplace (vx0 [i], vx0_y [i]);
202 map_x.emplace (vx1 [i], vx1_x [i]);
203 map_y.emplace (vx1 [i], vx1_y [i]);
204
205 map_d.emplace (edge_ [i], d [i]);
206 map_dw.emplace (edge_ [i], dw [i]);
207 map_time.emplace (edge_ [i], time [i]);
208 map_timew.emplace (edge_ [i], timew [i]);
209
210 map_object.emplace (edge_ [i], object [i]);
211 map_highway.emplace (edge_ [i], highway [i]);
212 }
213 }
214
215 const size_t n = junctions.size ();
216 std::vector <std::string> vx0_out (n),
217 vx1_out (n),
218 edge_out (n),
219 object_out (n),
220 highway_out (n),
221 old_edge_in (n),
222 old_edge_out (n);
223 std::vector <double> vx0_x_out (n),
224 vx0_y_out (n),
225 vx1_x_out (n),
226 vx1_y_out (n),
227 d_out (n),
228 dw_out (n),
229 time_out (n),
230 timew_out (n);
231
232 for (size_t i = 0; i < n; i++)
233 {
234 vx0_out [i] = junctions [i].v0;
235 vx0_x_out [i] = map_x.at (junctions [i].v0);
236 vx0_y_out [i] = map_y.at (junctions [i].v0);
237 vx1_out [i] = junctions [i].v1;
238 vx1_x_out [i] = map_x.at (junctions [i].v1);
239 vx1_y_out [i] = map_y.at (junctions [i].v1);
240 edge_out [i] = "j_" + sc::random_id (hash_len);
241 old_edge_in [i] = junctions [i].edge0;
242 old_edge_out [i] = junctions [i].edge1;
243
244 // Map all others to properties of out edge:
245 object_out [i] = map_object.at (junctions [i].edge1);
246 highway_out [i] = map_highway.at (junctions [i].edge1);
247
248 d_out [i] = map_d.at (junctions [i].edge0) +
249 map_d.at (junctions [i].edge1);
250 dw_out [i] = map_dw.at (junctions [i].edge0) +
251 map_dw.at (junctions [i].edge1);
252 time_out [i] = map_time.at (junctions [i].edge0) +
253 map_time.at (junctions [i].edge1);
254 timew_out [i] = map_timew.at (junctions [i].edge0) +
255 map_timew.at (junctions [i].edge1);
256 if (junctions [i].penalty)
257 {
258 time_out [i] += turn_penalty;
259 timew_out [i] += turn_penalty;
260 }
261 }
262
263 Rcpp::DataFrame res = Rcpp::DataFrame::create (
264 Rcpp::Named (".vx0") = vx0_out,
265 Rcpp::Named (".vx1") = vx1_out,
266 Rcpp::Named ("edge_") = edge_out, // empty here
267 Rcpp::Named (".vx0_x") = vx0_x_out,
268 Rcpp::Named (".vx0_y") = vx0_y_out,
269 Rcpp::Named (".vx1_x") = vx1_x_out,
270 Rcpp::Named (".vx1_y") = vx1_y_out,
271 Rcpp::Named ("d") = d_out,
272 Rcpp::Named ("object_") = object_out,
273 Rcpp::Named ("highway") = highway_out,
274 Rcpp::Named ("d_weighted") = dw_out,
275 Rcpp::Named ("time") = time_out,
276 Rcpp::Named ("time_weighted") = timew_out,
277 Rcpp::Named ("old_edge_in") = old_edge_in,
278 Rcpp::Named ("old_edge_out") = old_edge_out,
279 Rcpp::_["stringsAsFactors"] = false);
280
281 return res;
282}
283
284//' rcpp_route_times
285//'
286//' @noRd
287// [[Rcpp::export]]
288Rcpp::List rcpp_route_times (const Rcpp::DataFrame graph,
289 bool left_side, int turn_penalty)
290{
291 std::unordered_map <std::string, std::pair <RTEdgeSet, RTEdgeSet> > the_edges;
292 std::unordered_set <std::string> junction_vertices;
293
294 routetimes::fill_edges (graph, the_edges, junction_vertices);
295 std::vector <OneCompoundEdge> junctions;
296 routetimes::replace_junctions (the_edges, junctions, left_side);
297
298 Rcpp::DataFrame expanded_graph = routetimes::expand_edges (graph,
299 junctions, turn_penalty);
300
301 Rcpp::CharacterVector junction_vec (junction_vertices.size ());
302 size_t i = 0;
303 for (auto j: junction_vertices)
304 junction_vec (i++) = j;
305
306 return Rcpp::List::create (
307 Rcpp::Named ("graph") = expanded_graph,
308 Rcpp::Named ("junction_vertices") = junction_vec);
309}