1 /*!
2 \file vector/neta/timetables.c
3
4 \brief Network Analysis library - timetables
5
6 Shortest path using timetables.
7
8 (C) 2009-2010 by Daniel Bundala, and the GRASS Development Team
9
10 This program is free software under the GNU General Public License
11 (>=v2). Read the file COPYING that comes with GRASS for details.
12
13 \author Daniel Bundala (Google Summer of Code 2009)
14 */
15
16 #include <stdio.h>
17 #include <stdlib.h>
18 #include <grass/gis.h>
19 #include <grass/vector.h>
20 #include <grass/dbmi.h>
21 #include <grass/glocale.h>
22 #include <grass/dgl/graph.h>
23 #include <grass/neta.h>
24
25 /*!
26 \brief Get number of distinct elements
27
28 \param driver DB driver
29 \param sql SQl string
30 \param[out] list of lengths
31 \param[out] list of ids
32
33 \return number of distinct elements
34 \return -1 on failure
35 */
NetA_init_distinct(dbDriver * driver,dbString * sql,int ** lengths,int ** ids)36 int NetA_init_distinct(dbDriver * driver, dbString * sql, int **lengths,
37 int **ids)
38 {
39 int count, last, cur, result, index, more;
40 dbCursor cursor;
41 dbTable *table;
42 dbColumn *column;
43 dbValue *value;
44
45 if (db_open_select_cursor(driver, sql, &cursor, DB_SEQUENTIAL) != DB_OK) {
46 G_warning(_("Unable to open select cursor: %s"), db_get_string(sql));
47 return -1;
48 }
49 /*TODO: check column types */
50
51 count = last = 0;
52 /*count number of distinct routes */
53 table = db_get_cursor_table(&cursor);
54 column = db_get_table_column(table, 0);
55 while (db_fetch(&cursor, DB_NEXT, &more) == DB_OK && more) {
56 value = db_get_column_value(column);
57 cur = db_get_value_int(value);
58 if (count == 0 || cur != last) {
59 last = cur;
60 count++;
61 }
62 }
63 result = count;
64 db_close_cursor(&cursor);
65
66 *lengths = (int *)G_calloc(count, sizeof(int));
67 *ids = (int *)G_calloc(count, sizeof(int));
68 if (!*lengths || !*ids) {
69 G_warning(_("Out of memory"));
70 return -1;
71 }
72 db_open_select_cursor(driver, sql, &cursor, DB_SEQUENTIAL);
73 count = index = 0;
74 /*calculate the lengths of the routes */
75 table = db_get_cursor_table(&cursor);
76 column = db_get_table_column(table, 0);
77 while (db_fetch(&cursor, DB_NEXT, &more) == DB_OK && more) {
78 value = db_get_column_value(column);
79 cur = db_get_value_int(value);
80 if (count != 0 && cur != last)
81 index++;
82 if (count == 0 || cur != last)
83 (*ids)[index] = cur;
84 (*lengths)[index]++;
85 last = cur;
86 count++;
87 }
88 return result;
89 }
90
cmp_int(const void * a,const void * b)91 static int cmp_int(const void *a, const void *b)
92 {
93 return *(int *)a - *(int *)b;
94 }
95
96 /*!
97 \brief Initialises timetable from a database
98
99 \param In pointer to Map_info structure
100 \param route_layer layer number of routes
101 \param walk_layer layer number of walkers
102 \param route_id id of route
103 \param times list of timestamps
104 \param to_stop ?
105 \param walk_length walk length as string
106 \param timetable pointer to neta_timetable
107 \param route_ids list of route ids
108 \param stop_ids lits of stop ids
109
110 \return 0 on success
111 \return non-zero value on failure
112 */
NetA_init_timetable_from_db(struct Map_info * In,int route_layer,int walk_layer,char * route_id,char * times,char * to_stop,char * walk_length,neta_timetable * timetable,int ** route_ids,int ** stop_ids)113 int NetA_init_timetable_from_db(struct Map_info *In, int route_layer,
114 int walk_layer, char *route_id, char *times,
115 char *to_stop, char *walk_length,
116 neta_timetable * timetable, int **route_ids,
117 int **stop_ids)
118 {
119 int more, i, stop, route, time, *stop_pnt, stop1, stop2;
120 dbString sql;
121 dbCursor cursor;
122 dbTable *table;
123 dbColumn *column1, *column2, *column3;
124 dbValue *value;
125 char buf[2000];
126
127 dbDriver *driver;
128 struct field_info *Fi;
129
130 Fi = Vect_get_field(In, route_layer);
131 driver = db_start_driver_open_database(Fi->driver, Fi->database);
132 if (driver == NULL)
133 G_fatal_error(_("Unable to open database <%s> by driver <%s>"),
134 Fi->database, Fi->driver);
135
136
137 db_init_string(&sql);
138 sprintf(buf, "select %s from %s order by %s", route_id, Fi->table,
139 route_id);
140 db_set_string(&sql, buf);
141 timetable->routes =
142 NetA_init_distinct(driver, &sql, &(timetable->route_length),
143 route_ids);
144 if (timetable->routes < 0)
145 return 1;
146
147 sprintf(buf, "select %s from %s order by %s", Fi->key, Fi->table,
148 Fi->key);
149 db_set_string(&sql, buf);
150 timetable->stops =
151 NetA_init_distinct(driver, &sql, &(timetable->stop_length), stop_ids);
152 if (timetable->stops < 0)
153 return 1;
154
155 timetable->route_stops =
156 (int **)G_calloc(timetable->routes, sizeof(int *));
157 timetable->route_times =
158 (int **)G_calloc(timetable->routes, sizeof(int *));
159 timetable->stop_routes =
160 (int **)G_calloc(timetable->stops, sizeof(int *));
161 timetable->stop_times = (int **)G_calloc(timetable->stops, sizeof(int *));
162 timetable->walk_length = (int *)G_calloc(timetable->stops, sizeof(int));
163 timetable->walk_stops = (int **)G_calloc(timetable->stops, sizeof(int *));
164 timetable->walk_times = (int **)G_calloc(timetable->stops, sizeof(int *));
165 if (!timetable->route_stops || !timetable->route_times ||
166 !timetable->stop_routes || !timetable->stop_times ||
167 !timetable->walk_length) {
168 G_warning(_("Out of memory"));
169 return 2;
170 }
171
172 for (i = 0; i < timetable->routes; i++) {
173 timetable->route_stops[i] =
174 (int *)G_calloc(timetable->route_length[i], sizeof(int));
175 timetable->route_times[i] =
176 (int *)G_calloc(timetable->route_length[i], sizeof(int));
177 if (!timetable->route_stops[i] || !timetable->route_times[i]) {
178 G_warning(_("Out of memory"));
179 return 2;
180 }
181
182 timetable->route_length[i] = 0;
183 }
184
185 for (i = 0; i < timetable->stops; i++) {
186 timetable->stop_routes[i] =
187 (int *)G_calloc(timetable->stop_length[i], sizeof(int));
188 timetable->stop_times[i] =
189 (int *)G_calloc(timetable->stop_length[i], sizeof(int));
190 if (!timetable->stop_routes[i] || !timetable->stop_times[i]) {
191 G_warning(_("Out of memory"));
192 return 2;
193 }
194 timetable->walk_length[i] = 0;
195 timetable->stop_length[i] = 0;
196 }
197
198 sprintf(buf, "select %s, %s, %s from %s order by %s", Fi->key, route_id,
199 times, Fi->table, times);
200 db_set_string(&sql, buf);
201
202 if (db_open_select_cursor(driver, &sql, &cursor, DB_SEQUENTIAL) != DB_OK) {
203 G_warning(_("Unable to open select cursor: %s"), db_get_string(&sql));
204 return 1;
205 }
206
207
208 table = db_get_cursor_table(&cursor);
209 column1 = db_get_table_column(table, 0);
210 column2 = db_get_table_column(table, 1);
211 column3 = db_get_table_column(table, 2);
212 while (db_fetch(&cursor, DB_NEXT, &more) == DB_OK && more) {
213 value = db_get_column_value(column1);
214 stop = db_get_value_int(value);
215 value = db_get_column_value(column2);
216 route = db_get_value_int(value);
217 value = db_get_column_value(column3);
218 time = db_get_value_int(value);
219 stop =
220 (int *)bsearch(&stop, *stop_ids, timetable->stops, sizeof(int),
221 cmp_int) - (*stop_ids);
222 route =
223 (int *)bsearch(&route, *route_ids, timetable->routes, sizeof(int),
224 cmp_int) - (*route_ids);
225
226 timetable->stop_routes[stop][timetable->stop_length[stop]] = route;
227 timetable->stop_times[stop][timetable->stop_length[stop]++] = time;
228
229 timetable->route_stops[route][timetable->route_length[route]] = stop;
230 timetable->route_times[route][timetable->route_length[route]++] =
231 time;
232 }
233 db_close_cursor(&cursor);
234
235 if (walk_layer != -1) {
236
237 Fi = Vect_get_field(In, walk_layer);
238 sprintf(buf, "select %s, %s, %s from %s", Fi->key, to_stop,
239 walk_length, Fi->table);
240 db_set_string(&sql, buf);
241
242 if (db_open_select_cursor(driver, &sql, &cursor, DB_SEQUENTIAL) !=
243 DB_OK) {
244 G_warning(_("Unable to open select cursor: %s"),
245 db_get_string(&sql));
246 return 1;
247 }
248
249 table = db_get_cursor_table(&cursor);
250 column1 = db_get_table_column(table, 0);
251 column2 = db_get_table_column(table, 1);
252 while (db_fetch(&cursor, DB_NEXT, &more) == DB_OK && more) {
253 value = db_get_column_value(column2);
254 stop = db_get_value_int(value);
255 stop_pnt =
256 (int *)bsearch(&stop, *stop_ids, timetable->stops,
257 sizeof(int), cmp_int);
258 if (stop_pnt) {
259 value = db_get_column_value(column1);
260 stop = db_get_value_int(value);
261 stop_pnt =
262 (int *)bsearch(&stop, *stop_ids, timetable->stops,
263 sizeof(int), cmp_int);
264 if (stop_pnt) {
265 stop = stop_pnt - (*stop_ids);
266 timetable->walk_length[stop]++;
267 }
268 }
269 }
270 db_close_cursor(&cursor);
271
272 for (i = 0; i < timetable->stops; i++) {
273 timetable->walk_stops[i] =
274 (int *)G_calloc(timetable->walk_length[i], sizeof(int));
275 timetable->walk_times[i] =
276 (int *)G_calloc(timetable->walk_length[i], sizeof(int));
277 if (!timetable->walk_stops[i] || !timetable->walk_times[i]) {
278 G_warning(_("Out of memory"));
279 return 2;
280 }
281 timetable->walk_length[i] = 0;
282 }
283
284 if (db_open_select_cursor(driver, &sql, &cursor, DB_SEQUENTIAL) !=
285 DB_OK) {
286 G_warning(_("Unable to open select cursor: %s"),
287 db_get_string(&sql));
288 return 1;
289 }
290
291 table = db_get_cursor_table(&cursor);
292 column1 = db_get_table_column(table, 0);
293 column2 = db_get_table_column(table, 1);
294 column3 = db_get_table_column(table, 2);
295 while (db_fetch(&cursor, DB_NEXT, &more) == DB_OK && more) {
296 value = db_get_column_value(column2);
297 stop = db_get_value_int(value);
298 stop_pnt =
299 (int *)bsearch(&stop, *stop_ids, timetable->stops,
300 sizeof(int), cmp_int);
301 if (stop_pnt) {
302 stop2 = stop_pnt - (*stop_ids);
303 value = db_get_column_value(column1);
304 stop = db_get_value_int(value);
305 stop_pnt =
306 (int *)bsearch(&stop, *stop_ids, timetable->stops,
307 sizeof(int), cmp_int);
308 if (stop_pnt) {
309 stop1 = stop_pnt - (*stop_ids);
310 value = db_get_column_value(column3);
311 time = db_get_value_int(value);
312 timetable->walk_stops[stop1][timetable->
313 walk_length[stop1]] = stop2;
314 timetable->walk_times[stop1][timetable->
315 walk_length[stop1]++] = time;
316 }
317 }
318 }
319 db_close_cursor(&cursor);
320 }
321 db_close_database_shutdown_driver(driver);
322
323 return 0;
324 }
325
326 typedef struct
327 {
328 int v;
329 int conns;
330 } neta_heap_data;
331
new_heap_data(int conns,int v)332 static neta_heap_data *new_heap_data(int conns, int v)
333 {
334 neta_heap_data *d =
335 (neta_heap_data *) G_calloc(1, sizeof(neta_heap_data));
336 d->v = v;
337 d->conns = conns;
338 return d;
339 }
340
341 /*!
342 \brief Update Dijkstra structures
343
344 \param olc_conns old connection
345 \param new_conns new connection
346 \param to old 'to' node
347 \param new_dst new 'to' node
348 \param v ?
349 \param route id of route
350 \param rows ?
351 \param update ?
352 \param[out] result pointer to neta_timetable_result structure
353 \param heap ?
354 */
NetA_update_dijkstra(int old_conns,int new_conns,int to,int new_dst,int v,int route,int rows,int update,neta_timetable_result * result,dglHeap_s * heap)355 void NetA_update_dijkstra(int old_conns, int new_conns, int to, int new_dst,
356 int v, int route, int rows, int update,
357 neta_timetable_result * result, dglHeap_s * heap)
358 {
359 if (result->dst[new_conns][to] == -1 ||
360 result->dst[new_conns][to] > new_dst) {
361 result->dst[new_conns][to] = new_dst;
362 result->prev_stop[new_conns][to] = v;
363 result->prev_route[new_conns][to] = route;
364 result->prev_conn[new_conns][to] = old_conns;
365 if (update) {
366 dglHeapData_u heap_data;
367
368 heap_data.pv = (void *)new_heap_data(new_conns, to);
369 dglHeapInsertMin(heap, new_dst, ' ', heap_data);
370 }
371 }
372 }
373
374 /*!
375 \brief Computes the earliest arrival time.
376
377 Computes the earliest arrival time to to_stop from from_stop
378 starting at start_time, or -1 if no path exists
379
380 \param timetable pointer to neta_timetable structure
381 \param from_stop 'from' node
382 \param to_stop 'to' stop
383 \param start_time start timestamp
384 \param min_change ?
385 \param max_changes ?
386 \param walking_change ?
387 \param[out] result pointer to neta_timetable_result
388
389 \return ?
390 \return -1 on error
391 */
NetA_timetable_shortest_path(neta_timetable * timetable,int from_stop,int to_stop,int start_time,int min_change,int max_changes,int walking_change,neta_timetable_result * result)392 int NetA_timetable_shortest_path(neta_timetable * timetable, int from_stop,
393 int to_stop, int start_time, int min_change,
394 int max_changes, int walking_change,
395 neta_timetable_result * result)
396 {
397 int i, j;
398 dglHeap_s heap;
399
400 int opt_conns, rows = 1;
401
402 if (max_changes != -1)
403 rows = max_changes + 2;
404
405 result->rows = rows;
406 result->dst = (int **)G_calloc(rows, sizeof(int *));
407 result->prev_stop = (int **)G_calloc(rows, sizeof(int *));
408 result->prev_route = (int **)G_calloc(rows, sizeof(int *));
409 result->prev_conn = (int **)G_calloc(rows, sizeof(int *));
410
411 if (!result->dst || !result->prev_stop || !result->prev_route ||
412 !result->prev_conn) {
413 G_warning(_("Out of memory"));
414 return -1;
415 }
416
417 for (i = 0; i < rows; i++) {
418 result->dst[i] = (int *)G_calloc(timetable->stops, sizeof(int));
419 result->prev_stop[i] = (int *)G_calloc(timetable->stops, sizeof(int));
420 result->prev_route[i] =
421 (int *)G_calloc(timetable->stops, sizeof(int));
422 result->prev_conn[i] = (int *)G_calloc(timetable->stops, sizeof(int));
423 if (!result->dst[i] || !result->prev_stop[i] || !result->prev_route[i]
424 || !result->prev_conn[i]) {
425 G_warning(_("Out of memory"));
426 return -1;
427 }
428 }
429
430 if (from_stop == to_stop) {
431 result->dst[0][to_stop] = start_time;
432 result->prev_route[0][to_stop] = result->prev_stop[0][to_stop] = -1;
433 result->routes = 0;
434 return start_time;
435 }
436
437 result->routes = -1;
438 if (walking_change > 1)
439 walking_change = 1;
440 if (walking_change < 0 || max_changes == -1)
441 walking_change = 0;
442 dglHeapInit(&heap);
443
444 for (i = 0; i < rows; i++)
445 for (j = 0; j < timetable->stops; j++)
446 result->dst[i][j] = result->prev_stop[i][j] =
447 result->prev_route[i][j] = -1;
448
449 result->dst[0][from_stop] = start_time - min_change;
450 result->prev_stop[0][from_stop] = result->prev_route[0][from_stop] = -1;
451 dglHeapData_u heap_data;
452
453 heap_data.pv = (void *)new_heap_data(0, from_stop);
454 dglHeapInsertMin(&heap, start_time - min_change, ' ', heap_data);
455
456 while (1) {
457 dglInt32_t v, dist, conns;
458 dglHeapNode_s heap_node;
459 int new_conns, walk_conns, update;
460
461 if (!dglHeapExtractMin(&heap, &heap_node))
462 break;
463 v = ((neta_heap_data *) (heap_node.value.pv))->v;
464 conns = ((neta_heap_data *) (heap_node.value.pv))->conns;
465 dist = heap_node.key;
466
467 if (dist > result->dst[conns][v])
468 continue;
469 if (v == to_stop)
470 break;
471 new_conns = (max_changes == -1) ? 0 : (conns + 1);
472 walk_conns = conns + walking_change;
473
474 /*walking */
475 if (walk_conns < rows) {
476 /* update = (max_changes == -1 || walk_conns <= max_changes); */
477 update = 1;
478 for (i = 0; i < timetable->walk_length[v]; i++) {
479 int to = timetable->walk_stops[v][i];
480 int new_dst = dist + timetable->walk_times[v][i];
481
482 NetA_update_dijkstra(conns, walk_conns, to, new_dst, v, -2,
483 rows, update, result, &heap);
484 }
485 }
486
487 if (new_conns >= rows)
488 continue;
489 /*process all routes arriving after dist+min_change */
490 for (i = 0; i < timetable->stop_length[v]; i++)
491 if (timetable->stop_times[v][i] >= dist + min_change) {
492 int route = timetable->stop_routes[v][i];
493
494 /*find the index of v on the route */
495 for (j = 0; j < timetable->route_length[route]; j++)
496 if (timetable->route_stops[route][j] == v)
497 break;
498 j++;
499 for (; j < timetable->route_length[route]; j++) {
500 int to = timetable->route_stops[route][j];
501
502 NetA_update_dijkstra(conns, new_conns, to,
503 timetable->route_times[route][j], v,
504 route, rows, 1, result, &heap);
505 }
506 }
507 }
508 dglHeapFree(&heap, NULL);
509 opt_conns = -1;
510 for (i = 0; i < rows; i++)
511 if (result->dst[i][to_stop] != -1 &&
512 (opt_conns == -1 ||
513 result->dst[opt_conns][to_stop] > result->dst[i][to_stop]))
514 opt_conns = i;
515 result->routes = opt_conns;
516
517 if (opt_conns != -1)
518 return result->dst[opt_conns][to_stop];
519 return -1;
520 }
521
522 /*!
523 \brief Get time
524
525 Get time when route "route" arrives at stop "stop" or -1.
526
527 \param timetable pointer to neta_timetable structure
528 \param stop 'stop' node id
529 \param route route id
530
531 \return time
532 \return -1 if not found
533 */
NetA_timetable_get_route_time(neta_timetable * timetable,int stop,int route)534 int NetA_timetable_get_route_time(neta_timetable * timetable, int stop,
535 int route)
536 {
537 int i;
538
539 for (i = 0; i < timetable->route_length[route]; i++)
540 if (timetable->route_stops[route][i] == stop)
541 return timetable->route_times[route][i];
542 return -1;
543 }
544
545 /*!
546 \brief Free neta_timetable_result structure
547
548 \param result pointer to neta_timetable_result structure
549 */
NetA_timetable_result_release(neta_timetable_result * result)550 void NetA_timetable_result_release(neta_timetable_result * result)
551 {
552 int i;
553
554 for (i = 0; i < result->rows; i++) {
555 G_free(result->dst[i]);
556 G_free(result->prev_stop[i]);
557 G_free(result->prev_route[i]);
558 }
559 G_free(result->dst);
560 G_free(result->prev_stop);
561 G_free(result->prev_route);
562 }
563