1# Eclipse SUMO, Simulation of Urban MObility; see https://eclipse.org/sumo 2# Copyright (C) 2016-2019 German Aerospace Center (DLR) and others. 3# SUMOPy module 4# Copyright (C) 2012-2017 University of Bologna - DICAM 5# This program and the accompanying materials 6# are made available under the terms of the Eclipse Public License v2.0 7# which accompanies this distribution, and is available at 8# http://www.eclipse.org/legal/epl-v20.html 9# SPDX-License-Identifier: EPL-2.0 10 11# @file mapmatching-11-shortestproblem.py 12# @author Joerg Schweizer 13# @date 14# @version $Id$ 15 16""" 17This plugin provides methods to match GPS traces to a SUMO road network and to generate routes. 18 19The general structure is as follows: 20GPSPoints 21routes 22 23Traces 24 25""" 26import os 27import sys 28from xml.sax import saxutils, parse, handler 29 30import time 31import string 32import numpy as np 33from xml.sax import saxutils, parse, handler 34 35if __name__ == '__main__': 36 try: 37 APPDIR = os.path.dirname(os.path.abspath(__file__)) 38 except: 39 APPDIR = os.path.dirname(os.path.abspath(sys.argv[0])) 40 SUMOPYDIR = os.path.join(APPDIR, '..', '..') 41 sys.path.append(SUMOPYDIR) 42 43 44from coremodules.modules_common import * 45import agilepy.lib_base.classman as cm 46import agilepy.lib_base.arrayman as am 47import agilepy.lib_base.xmlman as xm 48#from agilepy.lib_base.misc import get_inversemap 49from agilepy.lib_base.geometry import get_dist_point_to_segs 50from agilepy.lib_base.processes import Process # ,CmlMixin,ff,call 51from coremodules.network.network import SumoIdsConf 52from coremodules.network import routing 53from coremodules.demand.demand import Trips, Routes, DemandobjMixin 54 55try: 56 try: 57 import pyproj 58 except: 59 from mpl_toolkits.basemap import pyproj 60 try: 61 from shapely.geometry import Polygon, MultiPolygon, MultiLineString, Point, LineString, MultiPoint, asLineString, asMultiPoint 62 from shapely.ops import cascaded_union 63 except: 64 print 'Import error: No shapely module available.' 65except: 66 print 'Import error: in order to run the traces plugin please install the following modules:' 67 print ' mpl_toolkits.basemap and shapely' 68 print 'Please install these modules if you want to use it.' 69 print __doc__ 70 raise 71 72TRIPPUROPSES = {'unknown': -1, 73 'HomeToWork': 1, 74 'WorkToHome': 2, 75 'HomeToSchool': 3, 76 'SchoolToHome': 4, 77 'SchoolToHome': 5, 78 'Leisure': 6, 79 'Other': 7, 80 } 81 82OCCUPATIONS = {'unknown': -1, 83 'Worker': 1, 84 'Student': 2, 85 'Employee': 3, 86 'Public employee': 4, 87 'Selfemployed': 5, 88 'Pensioneer': 6, 89 'Other': 7 90 } 91 92DEVICES = {'unknown': -1, 93 'cy-android': 1, 94 'cy-windows': 2, 95 'cy-ios': 3, 96 'cy-web-gpx': 4, 97 'cy-web-manual': 5, 98 'Other': 7 99 } 100 101WEEKDAYCHOICES = {'Monday-Friday': [0, 1, 2, 3, 4], 102 'Monday-Saturday': [0, 1, 2, 3, 4, 5], 103 'Saturday, Sunday': [5, 6], 104 'all': [0, 1, 2, 3, 4, 5, 6], 105 } 106 107COLOR_MATCHED_ROUTE = np.array([1.0, 0.4, 0.0, 0.6], np.float32) 108COLOR_SHORTEST_ROUTE = np.array([0.23529412, 1.0, 0.0, 0.6], np.float32) 109# Structure GPS traces 110 111 112# 2016 113 114# UserID TripID TimeStamp Start DT Distance ECC AvgSpeed TrackType Sex Year Profession Frequent User ZIP Source TypeOfBike TipeOfTrip Max Spd 115# 57249bcd88c537874f9fa1ae 57515edc88c537576ca3e16f 1464945480 2016-06-03T09:18:00.000Z 4.75 4.75 10.16 urban bicycle F 1999 Studente yes cy-web-gpx MyBike HomeToSchool 25.45 116# 5550800888c53765217661aa 574ec13788c537476aa3e122 1464797040 2016-06-01T16:04:00.000Z 2.93 2.93 5.87 urban bicycle F 1972 Worker yes 40136 cy-web-gpx -- -- 0 117# 5535586e88c53786637b23c6 574f3c3688c537c877a3e139 1464796080 2016-06-01T15:48:00.000Z 8.95 8.95 7.84 urban bicycle M 1973 Lavoratore yes 40133 cy-web-gpx MyBike HomeToWork 34.08 118# 5550800888c53765217661aa 574ec12b88c537812fa3e118 1464781800 2016-06-01T11:50:00.000Z 3.75 3.75 14.99 urban bicycle F 1972 Worker yes 40136 cy-web-gpx -- -- 0 119# 5550800888c53765217661aa 574ec11988c5370944a3e107 1464778980 2016-06-01T11:03:00.000Z 2.77 2.77 11.09 urban bicycle F 1972 Worker yes 40136 cy-web-gpx -- -- 0 120# 55e5edfc88c5374b4974fdbc 574e98c988c5378163a3e11f 1464765060 2016-06-01T07:11:00.000Z 11.09 11.14 15.63 urban bicycle 1982 SI cy-web-gpx MyBike HomeToWork 27.9 121 122 123# csv 2016 124# TripID, TimeStamp,Latitude, Longitude, Altitude, Distance, Speed, Type 125# 574e98c988c5378163a3e11f,1462347278,44.52606,11.27617,78,0.027255420500625783,5,<start|mid|end>, 126# 5741cdd388c537f10192ee97, 1463926725,44.50842,11.3604,101.0417,0.01615021486623964,3.483146,mid 127# timestamp: 1464160924 after 1970 128 129 130# 2015 csv 131 132# workouts 133# UserID TripID TimeStamp Start DT Distance AvgSpeed TrackType Sex Year Profession Frequent User ZIP 134# 54eb068de71f393530a9a74d 54eb0737e71f394c2fa9a74d 1424692504 Mon, 23 Feb 2015 11:55:04 GMT 0 0 urban bicycle M 1987 Developer no 135# 54eb9374e71f39f02fa9a750 5505cb04e71f39542e25e2d4 1426442994 Sun, 15 Mar 2015 18:09:54 GMT 0 0.7 urban bicycle M 1974 Worker yes 40128 136 137 138#TripID, TimeStamp,Date, Latitude, Longitude, Altitude, Distance, Speed, Type 139# 54eb0737e71f394c2fa9a74d,1424692509,"Mon, 23 Feb 2015 11:55:09 GMT",44.499096,11.361185,49.419395,0,0.000815,start 140 141 142# https://docs.python.org/2/library/time.html 143# >>> t = time.mktime((2017,03,07,13,46,0,0,0,0)) 144# >>> lt = time.localtime(t) 145# >>> time.strftime("%a, %d %b %Y %H:%M:%S +0000", lt) 146#'Tue, 07 Mar 2017 13:46:00 +0000' 147# >>> time.strftime("%a, %d %b %Y %H:%M:%S", lt) 148#'Tue, 07 Mar 2017 13:46:00' 149 150def load_results(filepath, parent=None, logger=None): 151 # typically parent is the scenario 152 results = cm.load_obj(filepath, parent=parent) 153 if logger is not None: 154 results.set_logger(logger) 155 return results 156 157 158def calc_seconds(t_data, 159 sep_date_clock=' ', sep_date='-', sep_clock=':', 160 is_float=False): 161 """ 162 Returns time in seconds after 1/1/1970. 163 Time format for time data string used: 164 2012-05-02 12:57:08.0 165 """ 166 167 if len(t_data.split(sep_date_clock)) != 2: 168 return -1 169 (date, clock) = t_data.split(sep_date_clock) 170 171 if (len(clock.split(sep_clock)) == 3) & (len(date.split(sep_date)) == 3): 172 173 (year_str, month_str, day_str) = date.split(sep_date) 174 # print ' year_str,month_str,day_str',year_str,month_str,day_str 175 (hours_str, minutes_str, seconds_str) = clock.split(sep_clock) 176 # print ' hours_str,minutes_str,seconds_str',hours_str,minutes_str,seconds_str 177 178 t = time.mktime((int(year_str), int(month_str), int(day_str), 179 int(hours_str), int(minutes_str), int(float(seconds_str)), -1, -1, -1)) 180 181 # print 'calc_seconds',t 182 # print ' t_data' 183 # print ' tupel',int(year_str),int(month_str),int(day_str), int(hours_str),int(minutes_str),int(float(seconds_str)),0,0,0 184 if is_float: 185 return t 186 else: 187 return int(t) 188 else: 189 return -1 190 191 192def get_routelinestring(route, edges): 193 # print 'get_routelinestring' 194 # TODO: we could make the np.sum of shapes if shapes where lists 195 shape = [] 196 for id_edge in route: 197 # print ' ',edges.shapes[id_edge ],type(edges.shapes[id_edge ]) 198 shape += list(edges.shapes[id_edge]) 199 # print ' ',shape 200 201 return LineString(shape) 202 203 204def get_boundary(coords): 205 # print 'get_boundary',len(coords),coords.shape 206 # print ' coords',coords 207 if len(coords) == 0: 208 return [0, 0, 0, 0] 209 else: 210 x_min, y_min = coords[:, :2].min(0) 211 x_max, y_max = coords[:, :2].max(0) 212 213 return [x_min, y_min, x_max, y_max] 214 215 216class BirgilMatcher(Process): 217 def __init__(self, ident, mapmatching, logger=None, **kwargs): 218 print 'BirgilMatcher.__init__' 219 220 # TODO: let this be independent, link to it or child?? 221 222 self._init_common(ident, 223 parent=mapmatching, 224 name='Birgillito Map matching', 225 logger=logger, 226 info='Birgillito Map matching.', 227 ) 228 229 attrsman = self.set_attrsman(cm.Attrsman(self)) 230 231 self.width_buffer_max = attrsman.add(cm.AttrConf('width_buffer_max', kwargs.get('width_buffer_max', 40.0), 232 groupnames=['options'], 233 perm='rw', 234 name='Max. buffer width', 235 unit='m', 236 info='GPS points are only valid if they are within the given maximum buffer width distant from a network edge.', 237 )) 238 239 self.n_points_min = attrsman.add(cm.AttrConf('n_points_min', kwargs.get('n_points_min', 5), 240 groupnames=['options'], 241 perm='rw', 242 name='Min. point number', 243 info='Minimum number of valid GPS points. Only if this minimum number is reached, a map-matching attempt is performed.', 244 )) 245 246 # self.dist_min = attrsman.add(cm.AttrConf( 'dist_min',kwargs.get('dist_min',3.0), 247 # groupnames = ['options'], 248 # name = 'Min. dist', 249 # unit = 'm', 250 # info = 'Minimum distance used in the calculation of the edge cost function.', 251 # )) 252 253 self.weight_cumlength = attrsman.add(cm.AttrConf('weight_cumlength', kwargs.get('weight_cumlength', 0.01), 254 groupnames=['options'], 255 name='Cumlength weight', 256 info='Cost function weight of cumulative length.', 257 )) 258 259 self.weight_angle = attrsman.add(cm.AttrConf('weight_angle', kwargs.get('weight_angle', 20.0), 260 groupnames=['options'], 261 name='Angle weight', 262 info='Cost function weight of angular coincidence.', 263 )) 264 265 self.weight_access = attrsman.add(cm.AttrConf('weight_access', kwargs.get('weight_access', 5.0), 266 groupnames=['options'], 267 name='Access weight', 268 info='Cost function weight of access level for matched mode.', 269 )) 270 271 self.n_routes_follow = attrsman.add(cm.AttrConf('n_routes_follow', kwargs.get('n_routes_follow', 20), 272 groupnames=['options'], 273 name='Followed routes', 274 info='Number of routes which are followed in parallel.', 275 )) 276 277 def do(self): 278 logger = self.get_logger() 279 timeconst = 0.95 280 281 trips = self.parent.trips 282 #routes = trips.get_routes() 283 edges = self.get_edges() 284 vtypes = self.parent.get_scenario().demand.vtypes 285 286 fstar = edges.get_fstar(is_return_arrays=True, is_ignor_connections=True) 287 times = edges.get_times(id_mode=2, is_check_lanes=False, speed_max=5.555) 288 ids_trip = trips.get_ids_selected() 289 290 ids_mode = vtypes.ids_mode[trips.ids_vtype[ids_trip]] 291 292 n_trips = len(ids_trip) 293 logger.w('Start mapmatching of %d GPS traces with Birgillito method' % n_trips) 294 295 distancesmap = self.parent.get_distancesmap() 296 accesslevelsmap = self.parent.get_accesslevelsmap() 297 298 n_trip_matched = 0 299 time_match_trace_av = 3.0 # initial gues of match time 300 301 # ,vtypes.speed_max[ids_vtype] 302 for id_trip, ids_point, id_mode in zip(ids_trip, trips.ids_points[ids_trip], ids_mode): 303 tick_before = time.time() 304 route, length_route, length_route_mixed, length_route_exclusive, duration_gps, lengthindex, err_dist, t_match, ids_point_edgeend, is_connected = \ 305 self.match_trip_birgil(id_trip, ids_point, fstar, distancesmap[id_mode], accesslevelsmap[id_mode]) 306 307 trips.set_matched_route(id_trip, route, 308 length_route, 309 length_route_mixed, 310 length_route_exclusive, 311 duration_gps, 312 lengthindex, 313 err_dist, 314 t_match, 315 is_connected=is_connected, 316 ids_point_edgeend=ids_point_edgeend) 317 n_trip_matched += 1 318 319 time_match_trace = time.time()-tick_before 320 321 time_match_trace_av = timeconst*time_match_trace_av + (1.0-timeconst)*time_match_trace 322 time_end_est = (n_trips-n_trip_matched)*time_match_trace_av 323 324 progress = 100.0*n_trip_matched/float(n_trips) 325 logger.w("Matched %d/%d traces, %.2f%% comleted. Avg match time = %.2fs, Terminated in %.1fh" % 326 (n_trip_matched, n_trips, progress, time_match_trace_av, float(time_end_est)/3600.0), key='message') 327 logger.w(progress, key='progress') 328 329 def get_edges(self): 330 return self.parent.get_scenario().net.edges 331 332 def match_trip_birgil(self, id_trip, ids_point, fstar, costs, accesslevels): 333 # TODO: record time intervals 334 # calculate initial and final position on edge to get more precise 335 # triplength 336 print 79*'=' 337 print 'match_trip_birgil', id_trip 338 tick = time.time() 339 routes = [] 340 route = None 341 route_mindist = None 342 intervals_route = [] 343 length_route = -1.0 344 length_bikeway = -1.0 345 length_mindist = -1.0 346 matchindex = -1.0 347 lengthindex = -1.0 348 err_dist = -1.0 349 err_dist_alg = -1.0 350 t_match = -1.0 351 duration_gps = -1.0 352 duration_est = -1.0 353 354 points = self.parent.points 355 #trips = self.parent.trips 356 #ids_point = trips.ids_points[id_trip] 357 #pointset = self.pointsets.get(id_trip) 358 359 # create a multi point object with points of traces 360 # should no longer be necessary 361 #tracepoints = MultiPoint(pointset.get_coords().tolist()) 362 363 # make a array with time stamps of all points 364 pointtimes = points.timestamps[ids_point] 365 # segind_to_id_edge = self._segind_to_id_edge #<<<<<<<<<<< 366 coords = points.coords[ids_point][:, :2] 367 368 x1, y1, x2, y2 = self.parent.get_segvertices_xy() 369 edges = self.get_edges() 370 get_ids_edge_from_inds_seg = edges.get_ids_edge_from_inds_seg 371 get_dist_point_to_edge = edges.get_dist_point_to_edge 372 #edgehit_counts = np.zeros(len(self._net.getEdges()),int) 373 # segind_to_id_edge = self._segind_to_id_edge 374 # for p in pointset.get_coords(): 375 # print ' ',p 376 # hitinds = get_dist_point_to_segs(p,self._x1,self._y1,self._x2,self._y2, is_ending=True)< self.width_buffer**2 377 # edgehit_counts[segind_to_id_edge[hitinds]]+= 1 378 379 # --------------------------------------------------------------------- 380 # 1. initialization 381 382 # Set of initial edges 383 384 # create a list with lists for each point 385 # Each list contans the edges, where the point is in the edge buffer 386 #edges_containing_point = np.zeros(len(ids_points), object) 387 #edgedists_to_point = np.zeros(len(ids_points), object) 388 389 ids_edge_initial = [] 390 ids_edge_final = [] 391 dists_initial = [] 392 n_points = len(ids_point) 393 ind_point_initial = -1 394 ind_point_final = -1 395 t_point_initial = -1 396 t_point_final = -1 397 398 # print ' search initial edges' 399 ind_point = 0 400 for p in coords: 401 dists2 = get_dist_point_to_segs(p, x1, y1, x2, y2) 402 inds_hit = np.flatnonzero(dists2 < self.width_buffer_max**2) 403 # print ' id_point,inds_hit',ids_point[ind_point],inds_hit 404 #ids_edge = get_ids_edge_from_inds_seg(inds_hit) 405 406 if len(inds_hit) > 0: 407 if ind_point_initial < 0: 408 # print ' inds_hit',ind_point,inds_hit 409 ids_edge_initial = set(get_ids_edge_from_inds_seg(inds_hit)) 410 # print ' ids_edge_initial',ids_edge_initial 411 #dists_initial = np.sqrt(dists2[inds_hit]) 412 ind_point_initial = ind_point 413 t_point_initial = pointtimes[ind_point] 414 # print ' ind_point, inds_hit, ind_point_initial',ind_point,inds_hit,ind_point_initial,ids_edge_initial 415 else: 416 ids_edge_final = set(get_ids_edge_from_inds_seg(inds_hit)) 417 # print ' ids_edge_final',ids_edge_final 418 ind_point_final = ind_point 419 t_point_final = pointtimes[ind_point] 420 # print ' ind_point, inds_hit, ind_point_final',ind_point,inds_hit,ind_point_initial,ids_edge_final 421 422 ind_point += 1 423 424 n_points_eff = ind_point_final - ind_point_initial 425 426 duration_gps = pointtimes[ind_point_final] - pointtimes[ind_point_initial] 427 428 # if self._logger: 429 # self._logger.w( '>>match_trip_birgil : n_points_eff=%d, len(ids_edge_initial)=%d,len(ids_edge_final)=%d'%(n_points_eff, len(ids_edge_initial),len(ids_edge_final)) ) 430 # print ' completed init:' 431 print ' ids_edge_initial', ids_edge_initial 432 print ' ids_edge_final', ids_edge_final 433 434 # print ' ind_point_initial,ind_point_final,n_points_eff',ind_point_initial,ind_point_final,n_points_eff 435 print ' id_point_initial=%d,id_point_final=%d,n_points_eff=%d' % ( 436 ids_point[ind_point_initial], ids_point[ind_point_final], n_points_eff) 437 if (ind_point_initial < 0) | (ind_point_final < 0) | (n_points_eff < self.n_points_min): 438 print 'ABOARD: insufficient valid points' 439 return [], 0.0, 0.0, -1.0, -1.0, -1.0, -1.0, 0.0, [], False 440 441 # print ' coords',coords 442 # create initial routeset 443 routelist = [] 444 445 #id_edge_unique = [] 446 for id_edge in ids_edge_initial: 447 448 cost0 = get_dist_point_to_edge(coords[ind_point_initial], id_edge) 449 450 # print ' coords of initial point:',coords[ind_point_initial] 451 # print ' cost0, id_edge',cost0, id_edge 452 453 #get_ind_seg_from_id_edge(self, id_edge) 454 455 length_edge = edges.lengths[id_edge] 456 accesslevel = accesslevels[id_edge] 457 cost_tot = cost0+self.weight_cumlength * length_edge - self.weight_access * accesslevel 458 routelist.append([cost_tot, cost0, 0.0, length_edge, accesslevel, length_edge, [id_edge], []]) 459 460 # print ' initial routelist',routelist 461 462 # --------------------------------------------------------------------- 463 # 2. main loop through rest of the points 464 465 # for i in xrange(ind_point+1, n_points):#coords[ind_point+1:]: 466 #is_connected = True 467 ind_point = ind_point_initial+1 468 length_gps = 0.0 469 while (ind_point < ind_point_final): # & is_connected: 470 point = coords[ind_point] 471 delta_point = point-coords[ind_point-1] 472 dist_interpoint = np.sqrt(np.sum(delta_point**2)) 473 phi_point = np.arctan2(delta_point[1], delta_point[0]) 474 length_gps += dist_interpoint 475 476 # print 79*'_' 477 # print ' check ID point %d, dist_interpoint=%.2fm, length_gps=%.2fm, phi_point %d deg'%(ids_point[ind_point],dist_interpoint,length_gps,phi_point/np.pi*180) 478 479 routelist_new = [] 480 ind_route = 0 481 for routeinfo in routelist: 482 costs_tot, cost, length_partial, length_cum, accesslevel, length_edge, ids_edge, ids_point_edgeend = routeinfo 483 484 if len(ids_edge) == 0: 485 break 486 487 length_partial += dist_interpoint 488 id_edge_last = ids_edge[-1] 489 # minimum distance point-edge 490 dist_point_edge = get_dist_point_to_edge(point, id_edge_last, 491 is_detect_initial=False, 492 is_detect_final=True,) 493 #dist_point_edge = self.get_dist_point_edge(coords[ind_point], id_edge_last, is_detect_final = True) 494 # if dist_point_edge is not a number (nan) the point is outside projection of edge 495 496 # print 79*'-' 497 # print ' route ',ind_route,'costs_tot',costs_tot,'dist_point_edge',dist_point_edge, len(ids_edge),id_edge_last 498 # print ' cost= %.3f, length_cum=%.3f'%(cost,length_cum,) 499 # print ' Xost= %.3f, Xength_cum=%.3f'%(cost,self.weight_cumlength *length_cum,) 500 # print ' ids_edge',ids_edge 501 # if length_partial > self.alpha * length_edge: 502 if (np.isnan(dist_point_edge)) | (dist_point_edge > self.width_buffer_max): 503 # point is beyond edge 504 505 # has reached end of edge 506 length_partial = length_partial-length_edge 507 508 # get FSTAR 509 # print ' fstar',id_edge_last,fstar[id_edge_last] 510 ids_edge_next_all = fstar[id_edge_last] 511 512 # filter allowed edges, by verifying positivness of 513 # travel costs 514 ids_edge_next = ids_edge_next_all[costs[ids_edge_next_all] >= 0] 515 # print ' parse ids_edge_next',ids_edge_next 516 517 if len(ids_edge_next) == 0: 518 # route not connected 519 cost_edge = np.Inf # punish route with infinite costs 520 length_edge = np.Inf 521 length_cum = np.Inf 522 accesslevel = -1 523 ids_edge_new = [] 524 ids_point_edgeend_new = [] 525 else: 526 527 id_edge = ids_edge_next[0] 528 accesslevel = accesslevels[id_edge] 529 length_edge = edges.lengths[id_edge] 530 length_cum += length_edge 531 cost_edge = self._get_cost_birgil(point, id_edge, phi_point, 532 accesslevel, get_dist_point_to_edge) 533 #self.get_dist_point_edge(coords[ind_point], ind_edge) 534 ids_edge_new = ids_edge + [id_edge] 535 ids_point_edgeend_new = ids_point_edgeend + [ids_point[ind_point]] 536 # print ' add to route id_edge',id_edge,ids_edge_new 537 538 # save updated data of current route 539 cost0 = cost+cost_edge 540 cost_tot = cost0+self.weight_cumlength * length_cum 541 routeinfo[:] = cost_tot, cost0, length_partial, length_cum, accesslevel, length_edge, ids_edge_new, ids_point_edgeend_new 542 543 # add new children if forward star greater 1 544 if len(ids_edge_next) > 1: 545 for id_edge in ids_edge_next[1:]: 546 #ind_edge = self._edge_to_id_edge[edge] 547 length_edge = edges.lengths[id_edge] 548 length_cum += length_edge 549 accesslevel = accesslevels[id_edge] 550 cost_edge = self._get_cost_birgil( 551 point, id_edge, phi_point, accesslevel, get_dist_point_to_edge) 552 ids_edge_new = ids_edge + [id_edge] 553 ids_point_edgeend_new = ids_point_edgeend + [ids_point[ind_point]] 554 # print ' new route id_edge',id_edge,ids_edge_new 555 cost0 = cost+cost_edge 556 cost_tot = cost0+self.weight_cumlength * length_cum 557 routelist_new.append([cost_tot, cost0, length_partial, length_cum, 558 accesslevel, length_edge, ids_edge_new, ids_point_edgeend_new]) 559 560 else: 561 # has not reached end of current edge 562 # save updated data of current route 563 if len(ids_edge) > 0: 564 #dist += self.get_dist_point_edge(coords[ind_point], ids_edge[-1]) 565 cost += self._get_cost_birgil(point, id_edge_last, phi_point, 566 accesslevel, get_dist_point_to_edge) 567 cost_tot = cost+self.weight_cumlength * length_cum 568 routeinfo[:] = cost_tot, cost, length_partial, length_cum, accesslevel, length_edge, ids_edge, ids_point_edgeend 569 570 ind_route += 1 571 572 routelist += routelist_new 573 # print ' routelist',routelist 574 routelist.sort() 575 576 # cut list to self.n_routes_follow 577 if len(routelist) > self.n_routes_follow+1: 578 routelist = routelist[:self.n_routes_follow+1] 579 580 ind_point += 1 581 582 # recover final edge objects of winner route! 583 costs_tot, cost, length_partial, length_cum, accesslevel, length_edge, ids_edge, ids_point_edgeend = routelist[ 584 0] 585 586 # print ' ids_edge[-1]',ids_edge[-1],ids_edge[-1] in ids_edge_final 587 # print ' ids_edge_final',ids_edge_final 588 t_match = time.time() - tick 589 590 # -------------------------------------------------------------------- 591 # post matching analisis 592 print '\n'+79*'-' 593 print ' matched route:', ids_edge 594 # print ' ids_point_edgeend',ids_point_edgeend 595 route = ids_edge 596 if len(route) == 0: 597 print 'ABOARD: route contains no edges ids_edge=', ids_edge 598 return [], 0.0, 0.0, -1.0, -1.0, -1.0, -1.0, 0.0, [], False 599 600 length_route = np.sum(costs[ids_edge]) 601 length_route_mixed = 0.0 602 length_route_exclusive = 0.0 603 for id_edge, length in zip(ids_edge, costs[ids_edge]): 604 accesslevel = accesslevels[id_edge] 605 if accesslevel == 1: 606 length_route_mixed += length 607 elif accesslevel == 2: 608 length_route_exclusive += length 609 610 if length_gps > 0: 611 lengthindex = length_route/length_gps 612 else: 613 lengthindex = -1.0 614 615 if ids_edge[-1] in ids_edge_final: # dist!=np.Inf: DID WE ARRIVE? 616 print 'SUCCESS: target', ids_edge[-1], 'reached' 617 is_connected = True 618 else: 619 print 'DISCONNECTED: last matched edge', ids_edge[-1], ' did not reach final edges', ids_edge_final 620 is_connected = False 621 622 # check dist error 623 inds_point = xrange(ind_point_initial, ind_point_final) 624 n_points = len(inds_point) 625 626 if (n_points >= 2) & (len(route) > 0): 627 # print ' initialize distance error measurement',len(route) 628 dist_points_tot = 0.0 629 routestring = get_routelinestring(route, edges) 630 # print ' routestring =',routestring 631 632 # print ' coords',coords[inds_point].tolist() 633 tracepoints = MultiPoint(coords[inds_point].tolist()) 634 635 # measure precision 636 #routeset = set(route) 637 #n_points_boundary = 0 638 ind_tracepoint = 0 639 for ind_point in inds_point: # xrange(n_points): 640 d = routestring.distance(tracepoints[ind_tracepoint]) 641 dist_points_tot += d 642 # print ' v distance measurement',d,tracepoints[ind_tracepoint],coords[ind_point]#,n_points#,n_points_eff 643 644 ind_tracepoint += 1 645 646 err_dist = dist_points_tot/float(n_points) 647 648 # print ' check',(lengthindex<self.lengthindex_max),(lengthindex >= self.lengthindex_min),(duration_gps>self.duration_min),(length_route>self.dist_min),(len(route)>0) 649 # print ' twice',(lengthindex<self.lengthindex_max)&(lengthindex >= self.lengthindex_min)&(duration_gps>self.duration_min)&(length_route>self.dist_min)&(len(route)>0) 650 651 return route, length_route, length_route_mixed, length_route_exclusive, duration_gps, lengthindex, err_dist, t_match, ids_point_edgeend, is_connected 652 653 def _get_cost_birgil(self, p, id_edge, phi_point, accesslevel, get_dist_point_to_edge): 654 # print '_get_cost_birgil p, id_edge',p, id_edge 655 656 dist, segment = get_dist_point_to_edge(p, id_edge, 657 is_return_segment=True) 658 x1, y1, x2, y2 = segment 659 660 # print ' x1,y1',x1,y1 661 # print ' p',p[0],p[1] 662 # print ' x2,y2',x2,y2 663 # print ' dist',dist 664 #seginds = self._id_edge_to_seginds[ind_edge] 665 #dists2 = get_dist_point_to_segs( p, self._x1[seginds], self._y1[seginds], self._x2[seginds], self._y2[seginds]) 666 667 #ind_segind = np.argmin(dists2) 668 #segind = seginds[ind_segind] 669 #x1,y1,x2,y2 = (self._x1[segind], self._y1[segind], self._x2[segind], self._y2[segind]) 670 phi_seg = np.arctan2(y2-y1, x2-x1) 671 #dist_ps = max(dist, self.dist_min) 672 # print ' x1,y1,x2,y2',x1,y1,x2,y2 673 # print ' dist =%.1f phi_point=%.2f,phi_seg=%.2f,phi_point-phi_seg=%.2f, delta_phi=%.2f'%(dist,phi_point/np.pi*180,phi_seg/np.pi*180,(phi_point-phi_seg)/np.pi*180,np.clip(phi_point-phi_seg,-np.pi/2,np.pi/2)/np.pi*180) 674 #cost = dist_ps*np.abs(np.sin(np.abs(phi_point-phi_seg))) 675 cost_angle = np.abs(np.sin(np.clip(phi_point-phi_seg, -np.pi/2, np.pi/2))) 676 cost = dist + self.weight_angle * cost_angle - self.weight_access * accesslevel 677 678 # print ' cost_birgil=%.3f, dist=%.3f, cost_angle=%.3f, accesslevel=%d'%(cost,dist,cost_angle,accesslevel) 679 # print ' cost_birgil=%.3f, Xist=%.3f, Xost_angle=%.3f, Xccesslevel=%.3f'%(cost,dist,self.weight_angle *cost_angle,self.weight_access*accesslevel) 680 # print ' cost_birgil=',cost 681 return cost 682 683 def get_dist_point_edge(self, p, ind_edge, is_detect_final=False): 684 seginds = self._id_edge_to_seginds[ind_edge] 685 dists2 = get_dist_point_to_segs(p, self._x1[seginds], self._y1[seginds], self._x2[seginds], 686 self._y2[seginds], is_ending=True, is_detect_final=is_detect_final) 687 if is_detect_final: 688 is_final = np.isnan(dists2) 689 if np.all(is_final): # point outside final segction of edge 690 return np.nan 691 else: 692 return np.sqrt(np.min(dists2[~is_final])) 693 694 return np.sqrt(np.min(dists2)) 695 696 def get_edge_info(self, id_trip, log=None): 697 # print 'get_edge_info for trace %s with c_length=%.6f'%(id_trip,self.c_length) 698 pointset = self.pointsets.get(id_trip) 699 700 # create a multi point object with points of traces 701 # should no longer be necessary 702 #tracepoints = MultiPoint(pointset.get_coords().tolist()) 703 704 # make a array with time stamps of all points 705 pointtimes = pointset.get_times() 706 ids_points = pointset.get_ids() 707 segind_to_id_edge = self._segind_to_id_edge 708 #edgehit_counts = np.zeros(len(self._net.getEdges()),int) 709 # segind_to_id_edge = self._segind_to_id_edge 710 # for p in pointset.get_coords(): 711 # print ' ',p 712 # hitinds = get_dist_point_to_segs(p,self._x1,self._y1,self._x2,self._y2, is_ending=True)< self.width_buffer**2 713 # edgehit_counts[segind_to_id_edge[hitinds]]+= 1 714 715 # create a list with lists for each point 716 # Each list contans the edges, where the point is in the edge buffer 717 edges_containing_point = np.zeros(len(ids_points), object) 718 edgedists_to_point = np.zeros(len(ids_points), object) 719 ind_point = 0 720 for p in pointset.get_coords(): 721 dists2 = get_dist_point_to_segs(p, self._x1, self._y1, self._x2, self._y2, is_ending=True) 722 inds_hit = np.nonzero(dists2 < self.width_buffer**2) 723 edges_containing_point[ind_point] = self._edges[segind_to_id_edge[inds_hit]] 724 edgedists_to_point[ind_point] = np.sqrt(dists2[inds_hit]) 725 726 # print ' ind_point,p,inds_hit',ind_point,p,len(np.nonzero(inds_hit)[0])#,edges_containing_point[ind_point] 727 ind_point += 1 728 729 # this is a dictionary with edge instance as key and 730 # the weight as value. 731 occurrences = {} 732 for edge in self._edges: 733 occurrences[edge] = 0.0 734 735 id_mode = self.id_mode 736 737 # 2. determine part of weight inversely proportional to the number 738 # of GPS points contained in the buffer 739 ind = 0 740 time_start = +10**10 741 time_end = -10**10 742 743 ind_point_start = -1 744 ind_point_end = -1 745 intervals = {} 746 for pointedges in edges_containing_point: 747 n_edges = len(pointedges) 748 if n_edges > 0: 749 t_point = pointtimes[ind] 750 if t_point < time_start: 751 time_start = t_point 752 ind_point_start = ind 753 754 elif t_point > time_end: 755 time_end = t_point 756 ind_point_end = ind 757 758 # compute partial weight for containing GPS points in edge buffer 759 wp = 1.0/n_edges 760 761 for edge in pointedges: 762 # assign weight component to edges 763 occurrences[edge] += wp 764 765 # determine intervals per edge 766 if intervals.has_key(edge): 767 t_edge_start, t_edge_end = intervals[edge] 768 if t_point < t_edge_start: 769 t_edge_start = t_point 770 elif t_point > t_edge_end: 771 t_edge_end = t_point 772 intervals[edge] = (t_edge_start, t_edge_end) 773 774 else: 775 intervals[edge] = (t_point, t_point) 776 777 ind += 1 778 779 edges_start = set(edges_containing_point[ind_point_start]) 780 edges_end = set(edges_containing_point[ind_point_end]) 781 782 # 1. calculate length proportional edge weight 783 c_length = float(self.c_length) 784 c_bike = float(self.c_bike) 785 weights = {} 786 787 for edge in self._edges: 788 len_edge = edge.getLength() 789 p_empty = c_length*len_edge 790 791 l = len_edge 792 allowed = edge.getLanes()[0].getAllowed() 793 if len(allowed) > 0: 794 if id_mode in allowed: 795 l *= c_bike # weight for dedicated lanes 796 797 if occurrences[edge] > 0.0: 798 # compute edge length proportional weight 799 weights[edge] = l/occurrences[edge] 800 801 else: 802 weights[edge] = l/p_empty 803 804 # if self._logger: 805 # self._logger.w( ' New: duration in network %d s from %ds to %ds'%(time_end-time_start,time_start,time_end)) 806 807 #log.w( ' Check dimensions len(weights)=%d,len(intervals)=%d'%(len(weights), len(intervals))) 808 # print ' intervals',intervals 809 #log.w( ' Possible start edges') 810 # for edge in edges_start: 811 # log.w( ' %s'%edge.getID()) 812 #log.w( ' Possible end edges') 813 # for edge in edges_end: 814 # log.w( ' %s'%edge.getID()) 815 # 816 #log.w('get_edge_info done.') 817 818 #self.intervals.set(id_trip, [ time_start,time_end]) 819 return weights, edges_start, edges_end, intervals, edges_containing_point, ind_point_start, ind_point_end 820 821 822def get_colvalue(val, default=0.0): 823 824 if (len(val) > 0) & (val != 'NULL'): 825 return float(val) 826 else: 827 return default 828 829 830class FilterMixin(Process): 831 def _init_filter_preview(self): 832 attrsman = self.get_attrsman() 833 trips = self.parent.trips 834 self.filterresults = attrsman.add(cm.FuncConf('filterresults', 835 'filterpreview', # function attribute of object 836 '%d/%d' % (len(trips.get_ids_selected()), len(trips)), 837 name='Preview selected trips', 838 groupnames=['options', 'results'], 839 )) 840 841 def filterpreview(self): 842 """ 843 Previews selected trips after filtering. 844 """ 845 trips = self.parent.trips 846 n_trips = len(trips) 847 n_sel_current = len(trips.get_ids_selected()) 848 n_sel_eliminate = len(self.filter_ids()) 849 n_sel_after = n_sel_current-n_sel_eliminate 850 return '%d/%d (currently %d/%d)' % (n_sel_after, n_trips, n_sel_current, n_trips) 851 852 def _init_traceoptions(self, **kwargs): 853 attrsman = self.get_attrsman() 854 self.dist_trip_min = attrsman.add(cm.AttrConf('dist_trip_min', kwargs.get('dist_trip_min', 100.0), 855 groupnames=['options'], 856 perm='rw', 857 name='Min. trip distance', 858 unit='m', 859 info='Minimum distance of one trip. Shorter trips will not be selected.', 860 )) 861 862 self.dist_trip_max = attrsman.add(cm.AttrConf('dist_trip_max', kwargs.get('dist_trip_max', 25000.0), 863 groupnames=['options'], 864 perm='rw', 865 name='Max. trip distance', 866 unit='m', 867 info='Maximum distance of one trip. Shorter trips will not be selected.', 868 )) 869 870 self.duration_trip_min = attrsman.add(cm.AttrConf('duration_trip_min', kwargs.get('duration_trip_min', 30.0), 871 groupnames=['options'], 872 perm='rw', 873 name='Min. trip duration', 874 unit='s', 875 info='Minimum duration of one trip. Trips with shorter duration will not be selected.', 876 )) 877 878 self.duration_trip_max = attrsman.add(cm.AttrConf('duration_trip_max', kwargs.get('duration_trip_max', 999999.0), 879 groupnames=['options'], 880 perm='rw', 881 name='Max. trip duration', 882 unit='s', 883 info='Maximum duration of one trip. Trips with longer duration will not be selected.', 884 )) 885 886 self.speed_trip_min = attrsman.add(cm.AttrConf('speed_trip_min', kwargs.get('speed_trip_min', 1.0), 887 groupnames=['options'], 888 perm='rw', 889 name='Min. av. trip speed', 890 unit='m/s', 891 info='Minimum average trip speed. Trips with lower average speed will not be selected.', 892 )) 893 894 self.speed_trip_max = attrsman.add(cm.AttrConf('speed_trip_max', kwargs.get('speed_trip_max', 50.0), 895 groupnames=['options'], 896 perm='rw', 897 name='Max. av. trip speed', 898 unit='m/s', 899 info='Maximum average trip speed. Trips with higher average speed will not be selected.', 900 )) 901 902 def _init_filter_time(self, **kwargs): 903 attrsman = self.get_attrsman() 904 905 self.hour_from_morning = attrsman.add(cm.AttrConf('hour_from_morning', kwargs.get('hour_from_morning', 5), 906 groupnames=['options'], 907 perm='rw', 908 name='From morning hour', 909 unit='h', 910 info='Keep only morning trips which start after this hour.', 911 )) 912 913 self.hour_to_morning = attrsman.add(cm.AttrConf('hour_to_morning', kwargs.get('hour_to_morning', 9), 914 groupnames=['options'], 915 perm='rw', 916 name='To morning hour', 917 unit='h', 918 info='Keep only morning trips which start before this hour.', 919 )) 920 921 self.hour_from_evening = attrsman.add(cm.AttrConf('hour_from_evening', kwargs.get('hour_from_evening', 15), 922 groupnames=['options'], 923 perm='rw', 924 name='From evening hour', 925 unit='h', 926 info='Keep only evening trips which start after this hour.', 927 )) 928 929 self.hour_to_evening = attrsman.add(cm.AttrConf('hour_to_evening', kwargs.get('hour_to_evening', 19), 930 groupnames=['options'], 931 perm='rw', 932 name='To evening hour', 933 unit='h', 934 info='Keep only evening trips which start before this hour.', 935 )) 936 937 self.weekdays = attrsman.add(cm.AttrConf('weekdays', kwargs.get('weekdays', WEEKDAYCHOICES['all']), 938 groupnames=['options'], 939 name='Weekdays', 940 choices=WEEKDAYCHOICES, 941 info='Keep only trips at the given weekdays.', 942 )) 943 944 def filter_time(self, timestamps): 945 """ 946 timestamps is an array with timestamps. 947 Function returns a binary array, with a True value for each 948 time stamp that does NOT SATISFY the specified time constrants. 949 """ 950 # print 'filter_time' 951 # print ' self.hour_from_morning,self.hour_to_morning',self.hour_from_morning,self.hour_to_morning 952 localtime = time.localtime 953 inds_elim = np.zeros(len(timestamps), dtype=np.bool) 954 i = 0 955 for timestamp in timestamps: 956 957 dt = localtime(timestamp) 958 is_keep = dt.tm_wday in self.weekdays 959 # print ' dt.tm_wday,self.weekdays',dt.tm_wday,self.weekdays 960 h = dt.tm_hour 961 is_keep &= (h > self.hour_from_morning) & (h < self.hour_to_morning)\ 962 | (h > self.hour_from_evening) & (h < self.hour_to_evening) 963 964 # print ' is_keep,w,h=',is_keep,dt.tm_wday,h 965 inds_elim[i] = not is_keep 966 i += 1 967 968 return inds_elim 969 970 def is_timestamp_ok(self, timestamp): 971 dt = time.localtime(timestamp) 972 is_ok = dt.tm_wday in self.weekdays 973 h = dt.tm_hour 974 is_ok &= (h > self.hour_from_morning) & (h < self.hour_to_morning)\ 975 | (h > self.hour_from_evening) & (h < self.hour_to_evening) 976 977 return is_ok 978 979 def filter_ids(self): 980 """ 981 Returns an array of ids to be eliminated or deselected. 982 To be overridden 983 """ 984 return [] 985 986 987class PostMatchfilter(FilterMixin): 988 def __init__(self, mapmatching, logger=None, **kwargs): 989 print 'PostMatchfilter.__init__' 990 self._init_common('postmatchfilter', 991 parent=mapmatching, 992 name='Post matchfilter', 993 logger=logger, 994 info='Removes matched tripe with defined characteristics from current selection.', 995 ) 996 997 attrsman = self.set_attrsman(cm.Attrsman(self)) 998 999 info_lengthindex = "Length index is length of matched route divided by length of line interpolated GPS points in percent." 1000 self.lengthindex_min = attrsman.add(cm.AttrConf('lengthindex_min', kwargs.get('lengthindex_min', 80.0), 1001 groupnames=['options'], 1002 name='Min. length index', 1003 unit='%', 1004 info='Minimum allowed length index. '+info_lengthindex, 1005 )) 1006 self.lengthindex_max = attrsman.add(cm.AttrConf('lengthindex_max', kwargs.get('lengthindex_max', 110.0), 1007 groupnames=['options'], 1008 name='Min. length index', 1009 unit='%', 1010 info='Maximum allowed length index '+info_lengthindex, 1011 )) 1012 1013 info_error_dist = 'The distance error is the average distance between the GPS points and the matched route.' 1014 self.error_dist_max = attrsman.add(cm.AttrConf('error_dist_max', kwargs.get('error_dist_max', 10000.0), 1015 groupnames=['options'], 1016 name='Max. distance err.', 1017 unit='mm', 1018 info='Maximum allowed distance error. '+info_error_dist, 1019 )) 1020 1021 self._init_filter_time(**kwargs) 1022 self._init_filter_preview() 1023 1024 def filter_ids(self): 1025 """ 1026 Returns an array of ids to be eliminated or deselected. 1027 """ 1028 # print 'filter_ids' 1029 trips = self.parent.trips 1030 ids_selected = trips.get_ids_selected() 1031 inds_eliminate = np.logical_or( 1032 trips.lengthindexes[ids_selected] < self.lengthindex_min, 1033 trips.lengthindexes[ids_selected] > self.lengthindex_max, 1034 trips.errors_dist[ids_selected] > self.error_dist_max, 1035 # (trips.lengths_route_matched[ids_selected]<1.0), # too many args?? 1036 ) 1037 inds_eliminate = np.logical_or(inds_eliminate, trips.lengths_route_matched[ids_selected] < 1.0) 1038 # print ' lengths_route_matched',trips.lengths_route_matched[ids_selected] 1039 # print ' lengthindexes',self.lengthindex_min,self.lengthindex_max,trips.lengthindexes[ids_selected] 1040 # print ' errors_dist',self.error_dist_max,trips.errors_dist[ids_selected] 1041 # print ' lengthindex_min',trips.lengthindexes[ids_selected]<self.lengthindex_min 1042 # print ' lengthindex_max',trips.lengthindexes[ids_selected]>self.lengthindex_max 1043 # print ' lengths_route_matched',trips.errors_dist[ids_selected]>self.error_dist_max 1044 # print ' inds_eliminate',inds_eliminate 1045 1046 inds_eliminate = np.logical_or(inds_eliminate, self.filter_time(trips.timestamps[ids_selected])) 1047 # print ' inds_eliminate',inds_eliminate 1048 return ids_selected[inds_eliminate] 1049 1050 def do(self): 1051 1052 # execute filtering 1053 self.parent.trips.are_selected[self.filter_ids()] = False 1054 1055 1056class TripGeomfilter(FilterMixin): 1057 def __init__(self, mapmatching, logger=None, **kwargs): 1058 print 'TripGeomfilter.__init__' 1059 self._init_common('tripGeomfilter', 1060 parent=mapmatching, 1061 name='Geometry trip filter', 1062 logger=logger, 1063 info='Removes trips with defined geometric characteristics from current selection.', 1064 ) 1065 1066 attrsman = self.set_attrsman(cm.Attrsman(self)) 1067 1068 self.dist_point_max = attrsman.add(cm.AttrConf('dist_point_max', kwargs.get('dist_point_max', 1000.0), 1069 groupnames=['options'], 1070 perm='rw', 1071 name='Max. point dist.', 1072 unit='m', 1073 info='Keep only traces where the distance between all successive points is below this maximal distance.', 1074 )) 1075 1076 self.duration_point_max = attrsman.add(cm.AttrConf('duration_point_max', kwargs.get('duration_point_max', 300.0), 1077 groupnames=['options'], 1078 perm='rw', 1079 name='Max. point duration.', 1080 unit='s', 1081 info='Keep only traces where the duration between all successive points is below this maximal duration.', 1082 )) 1083 1084 self.const_return_max = attrsman.add(cm.AttrConf('const_return_max', kwargs.get('const_return_max', 0.3), 1085 groupnames=['options'], 1086 perm='rw', 1087 name='Max. return const', 1088 info='Keep only traces where the user is not returning more than a share of this constant from the maximum (line of site) distance from the origin.', 1089 )) 1090 1091 self._init_filter_preview() 1092 1093 def do(self): 1094 1095 # execute filtering 1096 self.parent.trips.are_selected[self.filter_ids()] = False 1097 1098 def filter_ids(self): 1099 """ 1100 Returns an array of ids to be eliminated or deselected. 1101 """ 1102 c_cutoff = 1.0 - self.const_return_max 1103 print 'TripGeomfilter.do c_cutoff', c_cutoff 1104 dist_point_max = self.dist_point_max 1105 duration_point_max = self.duration_point_max 1106 trips = self.parent.trips 1107 points = self.parent.points 1108 1109 ids_trips = trips.get_ids_selected() 1110 1111 ids_points = trips.ids_points[ids_trips] 1112 # print ' ids_trips',ids_trips 1113 # print ' ids_points',ids_points 1114 1115 intersects_boundaries = self.parent.get_scenario().net.intersects_boundaries 1116 1117 inds_elim = np.zeros(len(ids_trips), dtype=np.bool) 1118 j = 0 1119 for id_trip, ids_point in zip(ids_trips, ids_points): 1120 1121 coords = points.coords[ids_point] 1122 times = points.timestamps[ids_point] 1123 # print 79*'-' 1124 # print ' check split ',id_trip 1125 # print ' ids_point', ids_point 1126 # print ' times',times 1127 # print ' coords', coords 1128 # print ' duration_point_max',duration_point_max 1129 #dist = pointset.get_distance() 1130 1131 if ids_point is None: 1132 # this happens if the points of a trip in the workout file 1133 # have not been imported for some reason 1134 is_eliminate = True 1135 1136 elif not intersects_boundaries(get_boundary(coords)): 1137 is_eliminate = True 1138 1139 elif np.any(times < 0): 1140 is_eliminate = True 1141 else: 1142 dist_max = 0.0 1143 is_eliminate = False 1144 i = 0 1145 n = len(times) 1146 x, y, z = coords[0] 1147 while (not is_eliminate) & (i < (n-1)): 1148 i += 1 1149 dist_point = np.sqrt((coords[i, 0]-coords[i-1, 0])**2 1150 + (coords[i, 1]-coords[i-1, 1])**2) 1151 1152 if dist_point > dist_point_max: 1153 is_eliminate = True 1154 1155 # print ' times[i]-times[i-1]',times[i]-times[i-1] 1156 1157 if (times[i]-times[i-1]) > duration_point_max: 1158 is_eliminate = True 1159 1160 d = np.sqrt((x-coords[i, 0])**2 + (y-coords[i, 1])**2) 1161 if d > dist_max: 1162 dist_max = d 1163 1164 elif d < c_cutoff*dist_max: 1165 is_eliminate = True 1166 1167 inds_elim[j] = is_eliminate 1168 j += 1 1169 1170 return ids_trips[inds_elim] 1171 1172 1173class EccTracesImporter(FilterMixin): 1174 def __init__(self, mapmatching, logger=None, **kwargs): 1175 print 'EccTracesImporter.__init__', mapmatching.get_ident() 1176 self._init_common('traceimporter', 1177 parent=mapmatching, 1178 name='ECC Trace Importer', 1179 logger=logger, 1180 info='Import workouts and GPS points of a European cycling challange.', 1181 ) 1182 1183 attrsman = self.set_attrsman(cm.Attrsman(self)) 1184 1185 scenario = mapmatching.get_scenario() 1186 rootfilepath = scenario.get_rootfilepath() 1187 1188 # here we ged classes not vehicle type 1189 # specific vehicle type within a class will be generated later 1190 modechoices = scenario.net.modes.names.get_indexmap() 1191 1192 # print ' modechoices',modechoices 1193 self.id_mode = attrsman.add(am.AttrConf('id_mode', modechoices['bicycle'], 1194 groupnames=['options'], 1195 choices=modechoices, 1196 name='Mode', 1197 info='Transport mode to be matched.', 1198 )) 1199 1200 self.workoutsfilepath = attrsman.add( 1201 cm.AttrConf('workoutsfilepath', kwargs.get('workoutsfilepath', rootfilepath+'.workouts.csv'), 1202 groupnames=['options'], 1203 perm='rw', 1204 name='Workout file', 1205 wildcards='CSV file (*.csv)|*.csv', 1206 metatype='filepath', 1207 info="""CSV text file with workout database.""", 1208 )) 1209 1210 self.pointsfilepath = attrsman.add(cm.AttrConf('pointsfilepath', kwargs.get('pointsfilepath', rootfilepath+'.points.csv'), 1211 groupnames=['options'], 1212 perm='rw', 1213 name='Points file', 1214 wildcards='CSV file (*.csv)|*.csv', 1215 metatype='filepath', 1216 info="CSV text file with GPS point database.", 1217 )) 1218 1219 self.year = attrsman.add(cm.AttrConf('year', kwargs.get('year', 2014), 1220 groupnames=['options'], 1221 choices={'2014': 2014, 'from 2015': 2015}, 1222 perm='rw', 1223 name='Year of challange', 1224 info='Year of challange is used to identify the correct database format.', 1225 )) 1226 1227 self._init_traceoptions(**kwargs) 1228 1229 self.sep_column_workout = attrsman.add(cm.AttrConf('sep_column_workout', kwargs.get('sep_column_workout', ','), 1230 groupnames=['options'], 1231 perm='rw', 1232 name='Workoutdata seperator', 1233 info='Workout column seperator of CSV file', 1234 )) 1235 1236 self.sep_column_points = attrsman.add(cm.AttrConf('sep_column_points', kwargs.get('sep_column_points', ','), 1237 groupnames=['options'], 1238 perm='rw', 1239 name='Point data seperator', 1240 info='Pointdata column seperator of CSV file', 1241 )) 1242 1243 self._init_filter_time(**kwargs) 1244 1245 def do(self): 1246 print 'TraceImporter.do' 1247 if self.year == 2014: 1248 self.import_workouts_2014() 1249 self.import_points_2014() 1250 1251 self.parent.points.project() 1252 1253 def import_workouts_2014(self): 1254 print 'import_workouts_2014' 1255 # 2014 ecomondo workouts 1256 # id pointDBNode pointPathId startTime distance duration sport calories maxSpeed altitudeMin altitudeMax metersAscent metersDescent 1257 # 329308466 7 37073516 2014-05-01 19:00:00 26 15600 1 1182.64 NULL NULL NULL NULL NULL 1258 # 0 1 2 3 4 5 6 7 8 9 10 11 12 1259 1260 j_id, j_node, j_id_trip, j_time, j_dist, j_duration = range(6) 1261 n_cols = 13 1262 null = 'NULL' 1263 trips = self.parent.trips 1264 # persons = self.parent.persons # no person data in 2014 :( 1265 ids_person_sumo = {} 1266 1267 ids_trips = [] 1268 scenario = self.parent.get_scenario() 1269 1270 get_vtype_for_mode = scenario.demand.vtypes.get_vtype_for_mode 1271 1272 f = open(self.workoutsfilepath, 'r') 1273 #if self._logger: self._logger.w('import_workouts_2014 %s'%os.path.basename(self.workoutsfilepath)) 1274 i_line = 0 1275 sep = self.sep_column_workout 1276 1277 #self.get_logger().w(100.0*self.simtime/self.duration, key ='progress') 1278 for line in f.readlines()[1:]: 1279 # if True:#i_line>1: 1280 cols = line.split(sep) 1281 # print ' len(cols)',len(cols),n_cols 1282 if len(cols) == n_cols: 1283 # row is complete 1284 1285 # if cols[j_dist] != null: 1286 dist = get_colvalue(cols[j_dist])*1000.0 1287 duration = get_colvalue(cols[j_duration]) 1288 if duration > 0: 1289 speed_av = dist/duration 1290 else: 1291 speed_av = 0.0 1292 1293 if (dist > self.dist_trip_min)\ 1294 & (dist < self.dist_trip_max)\ 1295 & (duration > self.duration_trip_min)\ 1296 & (speed_av > self.speed_trip_min)\ 1297 & (speed_av < self.speed_trip_max): 1298 # print 'parametric conditions verified' 1299 timestamp = calc_seconds(cols[j_time]) 1300 if timestamp is not None: 1301 # print ' valid time stamp',timestamp 1302 if self.is_timestamp_ok(timestamp): 1303 id_trip = trips.make(id_sumo=cols[j_id_trip], 1304 id_vtype=get_vtype_for_mode(self.id_mode), 1305 timestamp=timestamp, 1306 distance_gps=dist, 1307 duration_gps=duration, 1308 speed_average=speed_av, 1309 ) 1310 ids_trips.append(id_trip) 1311 1312 def import_points_2014(self): 1313 print 'import_points_2014' 1314 # csv2014 1315 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1316 # 4,61565791,23648171762,2013-05-01 06:33:58,44.501085,11.372906,NULL,0,NULL,2,NULL 1317 1318 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1319 # 7,37073516,15138524460,NULL,44.51579,11.36257,NULL,NULL,NULL,NULL,NULL 1320 1321 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1322 # 7,37073516,15138524460,NULL,44.51579,11.36257,NULL,NULL,NULL,NULL,NULL 1323 1324 # Endomondo export format 1325 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1326 # 0 1 2 3 4 5 6 7 8 9 10 1327 ind_id_path = 1 1328 ind_id_point = 2 1329 ind_time = 3 1330 ind_lat = 4 1331 ind_lon = 5 1332 ind_alt = 6 1333 ind_dist = 7 1334 ind_speed = 10 1335 1336 n_cols = 11 1337 #TripID, TimeStamp,Date, Latitude, Longitude, Altitude, Distance, Speed, Type 1338 1339 trips = self.parent.trips 1340 points = self.parent.points 1341 1342 exist_id_trip_sumo = trips.ids_sumo.has_index 1343 get_id_trip = trips.ids_sumo.get_id_from_index 1344 1345 sep = self.sep_column_points 1346 1347 f = open(self.pointsfilepath, 'r') 1348 if self._logger: 1349 self._logger.w('import_points_2014 %s' % os.path.basename(self.pointsfilepath)) 1350 i_line = 0 1351 id_trip_sumo = None 1352 id_trip = -1 1353 is_valid_trip = False 1354 1355 n_points_imported = 0 1356 for line in f.readlines(): 1357 1358 cols = line.split(sep) 1359 # print ' len(cols)',len(cols),n_cols 1360 if len(cols) == n_cols: 1361 id_trip_sumo_current = cols[ind_id_path] 1362 # print ' id_trip_sumo_current,id_trip_sumo',id_trip_sumo_current,id_trip_sumo,is_valid_trip 1363 if id_trip_sumo_current != id_trip_sumo: 1364 # this point is part of new trip 1365 1366 if is_valid_trip: 1367 # print ' store past points for valid trip',id_trip 1368 ids_point = points.add_rows( 1369 timestamps=timestamps, 1370 ids_trip=id_trip*np.ones(len(timestamps), dtype=np.int32), 1371 longitudes=longitudes, 1372 latitudes=latitudes, 1373 altitudes=altitudes, 1374 ) 1375 1376 trips.set_points(id_trip, ids_point) 1377 1378 # check if new trip is valid 1379 if exist_id_trip_sumo(id_trip_sumo_current): 1380 1381 is_valid_trip = True # start recording 1382 id_trip = get_id_trip(id_trip_sumo_current) 1383 # print ' found trip',id_trip,id_trip_sumo_current,' exisits-> record' 1384 1385 # ids_point_sumo = [] # useless? 1386 timestamps = [] 1387 ids_trip = [] 1388 longitudes = [] 1389 latitudes = [] 1390 altitudes = [] 1391 1392 else: 1393 # print ' trip',id_trip_sumo_current,'does not exisit' 1394 is_valid_trip = False 1395 1396 id_trip_sumo = id_trip_sumo_current 1397 1398 if is_valid_trip: 1399 # print ' store point timestamp',cols[ind_time] 1400 # current point belongs to a valid trip 1401 # ids_point_sumo.append(cols[ind_id_point]) 1402 timestamps.append(calc_seconds(cols[ind_time])) 1403 # ids_trip.append(id_trip) 1404 longitudes.append(get_colvalue(cols[ind_lon])) 1405 latitudes.append(get_colvalue(cols[ind_lat])) 1406 altitudes.append(get_colvalue(cols[ind_alt])) 1407 1408 else: 1409 print 'WARNING: inconsistent columns in line %d, file %s' % (i_line, filepath) 1410 1411 i_line += 1 1412 1413 # register points of last trip 1414 if is_valid_trip: 1415 # store past points for valid trip 1416 ids_point = points.add_rows( 1417 timestamps=timestamps, 1418 ids_trip=id_trip*np.ones(len(timestamps), dtype=np.int32), 1419 longitudes=longitudes, 1420 latitudes=latitudes, 1421 altitudes=altitudes, 1422 ) 1423 1424 trips.set_points(id_trip, ids_point) 1425 1426 # self.odtab.print_rows() 1427 1428 f.close() 1429 1430 1431class GpxImporter(FilterMixin): 1432 def __init__(self, mapmatching, logger=None, **kwargs): 1433 print 'GpxImporter.__init__', mapmatching.get_ident() 1434 self._init_common('gpximporter', 1435 parent=mapmatching, 1436 name='GPX Importer', 1437 logger=logger, 1438 info='Import GPS traces from GPX files.', 1439 1440 ) 1441 1442 attrsman = self.set_attrsman(cm.Attrsman(self)) 1443 1444 scenario = mapmatching.get_scenario() 1445 rootfilepath = scenario.get_rootfilepath() 1446 1447 # here we ged classes not vehicle type 1448 # specific vehicle type within a class will be generated later 1449 modechoices = scenario.net.modes.names.get_indexmap() 1450 1451 # print ' modechoices',modechoices 1452 self.id_mode = attrsman.add(am.AttrConf('id_mode', modechoices['bicycle'], 1453 groupnames=['options'], 1454 choices=modechoices, 1455 name='Mode', 1456 info='Transport mode to be matched.', 1457 )) 1458 1459 self.filepath = attrsman.add( 1460 cm.AttrConf('filepath', kwargs.get('filepath', rootfilepath+'.gpx'), 1461 groupnames=['options'], 1462 perm='rw', 1463 name='filename', 1464 wildcards='GPX file (*.gpx)|*.gpx|XML file (*.xml)|*.xml', 1465 metatype='filepath', 1466 info="""Path and file name of GPX file with GPS traces in XML format.""", 1467 )) 1468 self._init_traceoptions(**kwargs) 1469 1470 def do(self): 1471 """ 1472 Reads endomondo gpx xml file and stores point data in traces table. 1473 If there is no traces table, one will be initialized and returned 1474 """ 1475 mapmatching = self.parent 1476 #scenario = mapmatching.get_scenario() 1477 logger = self.get_logger() 1478 1479 get_vtype_for_mode = mapmatching.get_scenario().demand.vtypes.get_vtype_for_mode 1480 1481 parser = GpxParser(mapmatching.trips, mapmatching.points, 1482 logger=logger, 1483 id_vtype=get_vtype_for_mode(self.id_mode), 1484 dist_trip_min=self.dist_trip_min, 1485 dist_trip_max=self.dist_trip_max, 1486 duration_trip_min=self.speed_trip_min, 1487 duration_trip_max=self.duration_trip_max, 1488 speed_trip_min=self.speed_trip_min, 1489 speed_trip_max=self.speed_trip_max, 1490 ) 1491 1492 parse(self.filepath, parser) 1493 ids_trips = parser.get_ids_trip() 1494 if logger: 1495 logger.w('imported %d traces, project coordinates...' % len(ids_trips)) 1496 1497 # recalculate projection for network in scenario 1498 # for id_trace in ids_traces: 1499 # traces.pointsets.get(id_trace).project( traces._proj, traces.offset) 1500 1501 mapmatching.points.project() 1502 1503 if logger: 1504 logger.w('imported %d trips done.' % len(ids_trips)) 1505 1506 1507class GpxParser(handler.ContentHandler): 1508 """Reads endomondo gpx xml file and parses lat,lon, ele and time. 1509 """ 1510 1511 def __init__(self, trips, points, logger=None, 1512 id_vtype=0, 1513 dist_trip_min=10.0, 1514 dist_trip_max=50000.0, 1515 duration_trip_min=0.5, 1516 duration_trip_max=999999.0, 1517 speed_trip_min=0.1, 1518 speed_trip_max=100.0, 1519 ): 1520 1521 self._logger = logger 1522 self._points = points 1523 self._trips = trips 1524 self._id_vtype = id_vtype 1525 self._dist_trip_min = dist_trip_min 1526 self._dist_trip_max = dist_trip_max 1527 self._duration_trip_min = speed_trip_min 1528 self._duration_trip_max = duration_trip_max 1529 self._speed_trip_min = speed_trip_min 1530 self._speed_trip_max = speed_trip_max 1531 1532 self.reset_trip() 1533 1534 def reset_trip(self): 1535 self._lons = [] 1536 self._lats = [] 1537 self._eles = [] 1538 self._times = [] 1539 self._is_record_segment = False 1540 self._is_record_point = False 1541 self._is_time = False 1542 self._is_ele = False 1543 1544 self._ids_trip = [] 1545 self._ids_point = [] 1546 1547 def get_ids_trip(self): 1548 return self._ids_trip 1549 1550 def get_ids_point(self): 1551 return self._ids_point 1552 1553 def startElement(self, name, attrs): 1554 # if self._pointset is None: 1555 1556 if name == 'trkpt': 1557 self._is_record_point = True 1558 self._lons.append(float(attrs['lon'])) 1559 self._lats.append(float(attrs['lat'])) 1560 1561 if name == 'time': 1562 self._is_time = True 1563 1564 if name == 'ele': 1565 self._is_ele = True 1566 1567 if name == 'trkseg': 1568 self._is_record_segment = True 1569 1570 if name == 'trk': 1571 self._is_record_trip = True 1572 1573 def characters(self, content): 1574 1575 if self._is_time: 1576 # print ' got time',content[:-1] 1577 1578 if self._is_record_point: 1579 # 2013-05-30T16:42:33Z 2014-01-26T12:32:21Z 1580 self._times.append(calc_seconds(content[:-1], 1581 sep_date_clock='T', 1582 sep_date='-', 1583 sep_clock=':', 1584 is_float=False)) 1585 1586 if self._is_ele: 1587 # print 'characters content',content 1588 self._eles.append(float(content)) 1589 1590 def endElement(self, name): 1591 if name == 'ele': 1592 self._is_ele = False 1593 1594 if name == 'time': 1595 self._is_time = False 1596 1597 if name == 'trkpt': 1598 # print 'endElement: point',self._lon,self._lat,self._ele,self._time 1599 #self._is_ele = False 1600 #self._is_time = False 1601 self._is_record_point = False 1602 1603 if name == 'trk': 1604 # here we could stop recording a segment 1605 # but currently we join all segments together 1606 pass 1607 1608 if name == 'trk': 1609 # print 'endElement',len(self._lons) 1610 1611 # print ' self._lons',self._lons 1612 # print ' self._lats',self._lats 1613 # print ' self._times',self._times 1614 1615 # trip recording ends 1616 n_points = len(self._lons) 1617 1618 if (n_points > 0): 1619 timestamp = self._times[0] 1620 duration = self._times[-1]-timestamp 1621 # print ' timestamp',timestamp 1622 # print ' duration',duration,self._duration_trip_min,self._duration_trip_max 1623 # TODO: make this filter a functin of filtermixin 1624 if (duration > self._duration_trip_min)\ 1625 & (duration < self._duration_trip_max): 1626 1627 id_trip = self._trips.make(id_vtype=self._id_vtype, 1628 timestamp=timestamp, 1629 #distance_gps = dist, 1630 duration_gps=duration 1631 #speed_average = speed_av, 1632 ) 1633 self._ids_trip.append(id_trip) 1634 1635 if len(self._eles) == n_points: 1636 altitudes = self._eles 1637 else: 1638 altitudes = np.zeros(n_points, dtype=np.float32) 1639 1640 ids_point = self._points.add_rows( 1641 timestamps=self._times, 1642 ids_trip=id_trip*np.ones(n_points, dtype=np.int32), 1643 longitudes=self._lons, 1644 latitudes=self._lats, 1645 altitudes=altitudes, 1646 ) 1647 1648 self._trips.set_points(id_trip, ids_point) 1649 self._trips.ids_sumo[id_trip] = str(id_trip) 1650 1651 self.reset_trip() 1652 1653 1654class GpsTrips(Trips): 1655 def __init__(self, ident, mapmatching, **kwargs): 1656 # print 'Trips.__init__' 1657 self._init_objman(ident=ident, 1658 parent=mapmatching, 1659 name='Trips', 1660 info='Table with GPS trips, matched routes and alternative routes.', 1661 xmltag=('trips', 'trip', 'ids_sumo'), 1662 version=0.0, 1663 **kwargs) 1664 1665 self._init_attributes() 1666 self._init_constants() 1667 1668 def _init_attributes(self): 1669 1670 self.add_col(SumoIdsConf('GPS Trip', xmltag='id', info='GPS trip data.')) 1671 1672 self.add_col(am.ArrayConf('are_selected', default=True, 1673 dtype=np.bool, 1674 groupnames=['parameters', ], 1675 name='selected', 1676 symbol='Sel.', 1677 info='Selected for being processed (example mapmatching, export, etc).', 1678 )) 1679 1680 self.add_col(am.ArrayConf('timestamps', default=0, 1681 dtype=np.int, 1682 perm='r', 1683 groupnames=['parameters', 'gps'], 1684 name='timestamp', 1685 unit='s', 1686 metatype='datetime', 1687 info='Timestamp when trip started in seconds after 01 January 1970.', 1688 )) 1689 1690 self.timestamps.metatype = 'datetime' 1691 self.timestamps.set_perm('r') 1692 1693 self.add_col(am.ArrayConf('durations_gps', default=0.0, 1694 dtype=np.float32, 1695 groupnames=['parameters', 'gps'], 1696 name='GPS duration', 1697 unit='s', 1698 info='Time duration measure with GPS points.', 1699 )) 1700 1701 self.add_col(am.ArrayConf('distances_gps', default=0.0, 1702 dtype=np.float32, 1703 groupnames=['parameters', 'gps'], 1704 name='GPS distance', 1705 unit='m', 1706 info='Distance measure with GPS points.', 1707 )) 1708 1709 self.add_col(am.ArrayConf('speeds_average', default=0.0, 1710 dtype=np.float32, 1711 groupnames=['parameters', 'gps'], 1712 name='Av. speed', 1713 unit='m/s', 1714 info='Average speed based on GPS info.', 1715 )) 1716 self.add_col(am.ArrayConf('speeds_max', default=0.0, 1717 dtype=np.float32, 1718 groupnames=['parameters', 'gps'], 1719 name='Max. speed', 1720 unit='m/s', 1721 info='Maximum speed based on GPS info.', 1722 )) 1723 1724 # Trips._init_attributes(self) 1725 self.add_col(am.IdsArrayConf('ids_vtype', self.get_obj_vtypes(), 1726 groupnames=['state'], 1727 name='Type', 1728 info='Vehicle type.', 1729 xmltag='type', 1730 )) 1731 1732 self.add_col(am.ArrayConf('ids_purpose', default=TRIPPUROPSES['unknown'], 1733 dtype=np.int32, 1734 groupnames=['parameters', 'gps'], 1735 choices=TRIPPUROPSES, 1736 name='Purpose', 1737 info='Trip purpose ID', 1738 )) 1739 1740 self.add_col(am.ArrayConf('ids_device', default=DEVICES['unknown'], 1741 dtype=np.int32, 1742 groupnames=['parameters', 'gps'], 1743 choices=DEVICES, 1744 name='Devices', 1745 info='Device ID', 1746 )) 1747 self.add(cm.ObjConf(Routes('routes', self, self.parent.get_scenario().net))) 1748 1749 self.add_col(am.IdsArrayConf('ids_route_matched', self.get_routes(), 1750 groupnames=['results'], 1751 name='ID matched route', 1752 info='Route ID of mached route.', 1753 )) 1754 1755 self.add_col(am.IdsArrayConf('ids_route_shortest', self.get_routes(), 1756 groupnames=['results'], 1757 name='ID shortest route', 1758 info='Route ID of shortest route.', 1759 )) 1760 1761 self.add_col(am.ArrayConf('lengths_gpsroute_matched', default=-1.0, 1762 dtype=np.float32, 1763 groupnames=['results'], 1764 name='Matched GPS length', 1765 symbol='L_{MGPS}', 1766 unit='m', 1767 info='Length of the matched part of the GPS trace, measured by linear interpolation of GPS points. Note the only a fraction of the GPS trace ma be within the given network.', 1768 )) 1769 1770 self.add_col(am.ArrayConf('lengths_route_matched', default=-1.0, 1771 dtype=np.float32, 1772 groupnames=['results'], 1773 name='Matched length', 1774 symbol='L_{M}', 1775 unit='m', 1776 info='Length of the matched part of the GPS trace, measured by summing the length of edges of the matched route. Note the only a fraction of the GPS trace ma be within the given network.', 1777 )) 1778 1779 self.add_col(am.ArrayConf('durations_route_matched', default=-1.0, 1780 dtype=np.float32, 1781 groupnames=['results'], 1782 name='Matched duration', 1783 unit='s', 1784 info='Duration of the matched part of the GPS trace. This is the difference in timestamps between last and first GPS point of the matched route. Note the only a fraction of the GPS trace ma be within the given network.', 1785 )) 1786 1787 self.add_col(am.ArrayConf('lengths_route_shortest', default=-1.0, 1788 dtype=np.float32, 1789 groupnames=['results'], 1790 name='Shortest length', 1791 symbol='L_{S}', 1792 unit='m', 1793 info='Length of the shortest route. Shortest route is connecting the first matched edge and the final matched edge.', 1794 )) 1795 1796 self.add_col(am.ArrayConf('lengths_route_matched_mixed', default=-1.0, 1797 dtype=np.float32, 1798 groupnames=['results'], 1799 name='Matched length mixed access', 1800 symbol='L_{MM}', 1801 unit='m', 1802 info='Length of the matched part of the GPS trace. Note the only a fraction of the GPS trace ma be within the given network.', 1803 )) 1804 1805 self.add_col(am.ArrayConf('lengths_route_matched_exclusive', default=-1.0, 1806 dtype=np.float32, 1807 groupnames=['results'], 1808 name='Matched length exclusive access', 1809 symbol='L_{ME}', 1810 unit='m', 1811 info='Length of the matched part of the GPS trace. Note the only a fraction of the GPS trace ma be within the given network.', 1812 )) 1813 1814 self.add_col(am.ArrayConf('lengthindexes', default=-1.0, 1815 dtype=np.float32, 1816 groupnames=['results'], 1817 name='Length index', 1818 unit='%', 1819 info='Length index is the length of the matched route divided by length of line-interpolated GPS points.', 1820 )) 1821 1822 self.add_col(am.ArrayConf('errors_dist', default=-1.0, 1823 dtype=np.float32, 1824 groupnames=['results'], 1825 name='Distance error', 1826 unit='mm', 1827 info='The distance error is the average distance between the GPS points and the matched route.', 1828 )) 1829 1830 self.add_col(am.ArrayConf('times_computation', default=-1.0, 1831 dtype=np.float32, 1832 groupnames=['results'], 1833 name='Computation time', 1834 unit='ms', 1835 info='Computation time of the match algorithm.', 1836 )) 1837 1838 self.add_col(am.ArrayConf('are_match_connected', default=False, 1839 dtype=np.bool, 1840 groupnames=['results', ], 1841 name='Match connected', 1842 #symbol = 'Match conn.', 1843 info='The matched route connects first and last network edge where GPS points have been detected.', 1844 )) 1845 1846 self.add_col(am.IdlistsArrayConf('ids_points', self.parent.points, 1847 #groupnames = ['_private'], 1848 name='Point IDs', 1849 info="GPS point IDs.", 1850 )) 1851 1852 self.add_col(am.IdlistsArrayConf('ids_points_edgeend', self.parent.points, 1853 groupnames=['results', ], 1854 name='Edge endpoint IDs', 1855 info="This is a list of GPS point IDs which represent the last point associated with each matched edge.", 1856 )) 1857 1858 def get_obj_vtypes(self): 1859 return self.parent.get_scenario().demand.vtypes 1860 1861 def get_routes(self): 1862 return self.routes.get_value() 1863 1864 def clear_routes(self): 1865 self.get_routes().clear() 1866 for attrconf in self.get_group('results'): 1867 attrconf.reset() 1868 1869 def select_all(self): 1870 self.are_selected.get_value()[:] = True 1871 1872 def unselect_all(self): 1873 self.are_selected.get_value()[:] = False 1874 1875 def invert_selection(self): 1876 self.are_selected.get_value()[:] = np.logical_not(self.are_selected.get_value()[:]) 1877 1878 def delete_unselected(self): 1879 ids_del = self.select_ids(np.logical_not(self.are_selected.get_value())) 1880 points = self.parent.points 1881 for ids_point in self.ids_points[ids_del]: 1882 if ids_point is not None: 1883 points.del_rows(ids_point) 1884 self.del_rows(ids_del) 1885 1886 def set_matched_route(self, id_trip, route, 1887 length_matched=0.0, 1888 length_route_mixed=0.0, 1889 length_route_exclusive=0.0, 1890 duration_matched=0.0, 1891 lengthindex=-1.0, 1892 error_dist=-1.0, 1893 comptime=0.0, 1894 is_connected=False, 1895 ids_point_edgeend=[], 1896 ): 1897 1898 if len(route) > 0: 1899 id_route = self.ids_route_matched[id_trip] 1900 if id_route >= 0: 1901 # already a matched route existant 1902 self.get_routes().ids_edges[id_route] = route 1903 self.get_routes().colors[id_route] = COLOR_MATCHED_ROUTE 1904 # self.get_routes().set_row( id_route, 1905 # ids_edges = route, 1906 # colors = COLOR_MATCHED_ROUTE, 1907 # ) 1908 1909 # self.get_routes().set_row( id_route, 1910 # ids_trip = id_trip, 1911 # ids_edges = route, 1912 # #costs = duration_matched, 1913 # #probabilities = 1.0, 1914 # ) 1915 else: 1916 id_route = self.get_routes().add_row(ids_trip=id_trip, 1917 ids_edges=route, 1918 #costs = duration_matched, 1919 #probabilities = 1.0, 1920 colors=COLOR_MATCHED_ROUTE, 1921 ) 1922 else: 1923 id_route = -1 1924 1925 # print 'set_matched_route id_trip', id_trip,'id_route', id_route 1926 # print ' ids_point_edgeend',ids_point_edgeend 1927 #self.ids_route_matched[id_trip] = id_route 1928 1929 self.set_row(id_trip, 1930 ids_route_matched=id_route, 1931 lengths_gpsroute_matched=length_matched/lengthindex, 1932 durations_route_matched=duration_matched, 1933 lengths_route_matched=length_matched, 1934 lengths_route_matched_mixed=length_route_mixed, 1935 lengths_route_matched_exclusive=length_route_exclusive, 1936 lengthindexes=100*lengthindex, 1937 errors_dist=1000 * error_dist, 1938 times_computation=1000*comptime, 1939 are_match_connected=is_connected 1940 ) 1941 # TODO: do this extra! this is a bug!! 1942 # if included in set_row, only first value in list is taken!!! 1943 self.ids_points_edgeend[id_trip] = ids_point_edgeend 1944 1945 # print ' ids_points_edgeend[id_trip]',self.ids_points_edgeend[id_trip] 1946 1947 return id_route 1948 1949 def make(self, **kwargs): 1950 # if self.ids_sumo.has_index(id_sumo): 1951 # id_trip = self.ids_sumo.get_id_from_index(id_sumo) 1952 # #self.set_row(id_sumo, **kwargs) 1953 # return id_trip 1954 # else: 1955 id_trip = self.add_row(ids_sumo=kwargs.get('id_sumo', None), 1956 ids_vtype=kwargs.get('id_vtype', None), 1957 timestamps=kwargs.get('timestamp', None), 1958 distances_gps=kwargs.get('distance_gps', None), 1959 durations_gps=kwargs.get('duration_gps', None), 1960 speeds_average=kwargs.get('speed_average', None), 1961 speeds_max=kwargs.get('speed_max', None), 1962 ids_purpose=kwargs.get('id_purpose', None), 1963 ids_device=kwargs.get('id_device', None), 1964 ids_points=kwargs.get('ids_point', None), 1965 ) 1966 return id_trip 1967 1968 def set_points(self, id_trip, ids_point): 1969 self.ids_points[id_trip] = ids_point 1970 1971 def get_ids_selected(self): 1972 return self.select_ids(self.are_selected.get_value()) 1973 1974 def route_shortest(self, dist_modespecific=5.0, c_modespecific=0.9, 1975 is_ignor_connections=False): 1976 """ 1977 Shortest path routint. 1978 """ 1979 print 'route_shortest', dist_modespecific, c_modespecific 1980 # TODO: if too mant vtypes, better go through id_modes 1981 exectime_start = time.clock() 1982 scenario = self.parent.get_scenario() 1983 net = scenario.net 1984 edges = net.edges 1985 vtypes = scenario.demand.vtypes 1986 routes = self.get_routes() 1987 ids_edges = [] 1988 ids_trip = [] 1989 costs = [] 1990 1991 distancesmap = self.parent.get_distancesmap() 1992 accesslevelsmap = self.parent.get_accesslevelsmap() 1993 1994 # delete current 1995 ids_with_shortes = self.select_ids(np.logical_and( 1996 self.are_selected.get_value(), self.ids_route_shortest.get_value() >= 0)) 1997 routes.del_rows(self.ids_route_shortest[ids_with_shortes]) 1998 self.ids_route_shortest[ids_with_shortes] = -1 1999 2000 fstar = edges.get_fstar(is_ignor_connections=is_ignor_connections) 2001 for id_vtype in self.get_vtypes(): 2002 id_mode = vtypes.ids_mode[id_vtype] 2003 2004 # no routing for pedestrians 2005 if id_mode != net.modes.get_id_mode('pedestrian'): 2006 weights = distancesmap[id_mode] 2007 2008 # this will subtract some meters dependent on 2009 # access-level of the edge 2010 accesslevels = accesslevelsmap[id_mode] 2011 ids_edge = edges.get_ids() 2012 are_valid = weights > 3*dist_modespecific 2013 weights[ids_edge] -= dist_modespecific * are_valid[ids_edge] * accesslevels[ids_edge] 2014 weights[ids_edge] *= 1 - (1-c_modespecific) * are_valid[ids_edge] * (accesslevels[ids_edge] == 2) 2015 2016 ids_trip_vtype = self.select_ids(np.logical_and(self.ids_vtype.get_value( 2017 ) == id_vtype, self.are_selected.get_value(), self.ids_route_matched.get_value() >= 0)) 2018 #ids_trip_vtype = self.get_trips_for_vtype(id_vtype) 2019 # print ' id_vtype,id_mode',id_vtype,id_mode#,ids_trip_vtype 2020 # print ' weights',weights 2021 2022 #ids_edge_depart = self.ids_edge_depart[ids_trip_vtype] 2023 #ids_edge_arrival = self.ids_edge_arrival[ids_trip_vtype] 2024 2025 for id_trip, id_route in zip(ids_trip_vtype, self.ids_route_matched[ids_trip_vtype]): 2026 route_matched = routes.ids_edges[id_route] 2027 2028 cost, route = routing.get_mincostroute_edge2edge(route_matched[0], 2029 route_matched[-1], 2030 weights=weights, 2031 fstar=fstar) 2032 if len(route) > 0: 2033 ids_edges.append(route) 2034 ids_trip.append(id_trip) 2035 costs.append(cost) 2036 2037 ids_route = routes.add_rows(ids_trip=ids_trip, 2038 ids_edges=ids_edges, 2039 costs=costs, 2040 colors=COLOR_SHORTEST_ROUTE*np.ones((len(ids_trip), 4), dtype=np.float32), 2041 ) 2042 2043 self.ids_route_shortest[ids_trip] = ids_route 2044 self.lengths_route_shortest[ids_trip] = costs 2045 #self.add_routes(ids_trip, ids_route) 2046 print ' exectime', time.clock()-exectime_start 2047 return ids_trip, ids_route 2048 2049 def get_speedprofile(self, id_trip): 2050 2051 points = self.parent.points 2052 ids_point = self.ids_points[id_trip] 2053 ids_point_edgeend = self.ids_points_edgeend[id_trip] 2054 ids_edge = self.get_routes().ids_edges[self.ids_route_matched[id_trip]] 2055 id_edge_current = -1 2056 n_edges = len(ids_edge) 2057 2058 positions_gps = [] 2059 times_gps = [] 2060 ids_edges_profile = [] 2061 ids_nodes_profile = [] 2062 2063 for id_point, coord, timestamp in zip(ids_point, points.coords[ids_point], points.timestamps[ids_point]): 2064 #get_pos_from_coord(id_edge_current, coord) 2065 2066 if id_point in ids_point_edgeend: 2067 ind = ids_point_edgeend.index(id_point) 2068 if ind == n_edges-1: 2069 id_edge_current = -1 2070 else: 2071 id_edge_current = ids_edge[ind+1] 2072 2073 def get_flows(self, is_shortest_path=False): 2074 """ 2075 Determine the total number of vehicles for each edge. 2076 returns ids_edge and flows 2077 """ 2078 ids_edges = self.get_routes.ids_edges 2079 counts = np.zeros(np.max(self.get_net().edges.get_ids())+1, int) 2080 2081 ids_trip = self.get_ids_selected() 2082 if not is_shortest_path: 2083 ids_route = self.ids_route_matched[ids_trip] 2084 else: 2085 ids_route = self.ids_route_shortest[ids_trip] 2086 inds_valid = np.flatnonzero(ids_route > 0) 2087 for id_trip, id_route in zip(ids_trip[inds_valid], ids_route[inds_valid]): 2088 counts[ids_edges[id_route][:]] += 1 2089 2090 ids_edge = np.flatnonzero(counts) 2091 2092 return ids_edge, counts[ids_edge].copy() 2093 2094 def get_ids_route_selected(self): 2095 # TODO: here we could append direct routes 2096 # print 'get_ids_route_selected' 2097 ids_route_matched = self.ids_route_matched[self.get_ids_selected()] 2098 ids_route_shortest = self.ids_route_shortest[self.get_ids_selected()] 2099 ids_route = np.concatenate([ids_route_matched[ids_route_matched >= 0], 2100 ids_route_shortest[ids_route_shortest >= 0]]) 2101 # print ' ids_route',ids_route 2102 # print ' ids_route_matched[ids_route_matched >= 0] ',ids_route_matched[ids_route_matched >= 0] 2103 return np.concatenate([ids_route_matched[ids_route_matched >= 0], ids_route_shortest[ids_route_shortest >= 0]]) 2104 # return ids_route_matched[ids_route_matched >= 0] 2105 2106 2107class GpsPoints(am.ArrayObjman): 2108 """ 2109 Contains data of points of a single trace. 2110 """ 2111 2112 def __init__(self, ident, mapmatching, **kwargs): 2113 # print 'PrtVehicles vtype id_default',vtypes.ids_sumo.get_id_from_index('passenger1') 2114 self._init_objman(ident=ident, 2115 parent=mapmatching, 2116 name='GPS Points', 2117 info='GPS points database.', 2118 version=0.0, 2119 **kwargs) 2120 2121 self._init_attributes() 2122 self._init_constants() 2123 2124 def _init_attributes(self): 2125 scenario = self.get_scenario() 2126 2127 # ident is the path id of the trace 2128 2129 # the actual numpy arrays are stored in .cols 2130 self.add_col(am.ArrayConf('longitudes', default=0.0, 2131 dtype=np.float32, 2132 groupnames=['parameters', ], 2133 perm='rw', 2134 name='Longitude', 2135 symbol='Lon', 2136 unit='deg', 2137 info='Longitude of point', 2138 )) 2139 2140 self.add_col(am.ArrayConf('latitudes', default=0.0, 2141 groupnames=['parameters', ], 2142 dtype=np.float32, 2143 perm='rw', 2144 name='Latitude', 2145 symbol='Lat', 2146 unit='deg', 2147 info='Latitude of point', 2148 )) 2149 2150 self.add_col(am.ArrayConf('altitudes', default=-100000.0, 2151 groupnames=['parameters', ], 2152 dtype=np.float32, 2153 perm='rw', 2154 name='Altitude', 2155 symbol='Alt', 2156 unit='m', 2157 info='Altitude of point', 2158 )) 2159 2160 self.add_col(am.ArrayConf('radii', 10.0, 2161 dtype=np.float32, 2162 groupnames=['parameters', ], 2163 perm='rw', 2164 name='Radius', 2165 unit='m', 2166 info='Point radius, representing the imprecision of the point, which depends on the recording device ane the environment.', 2167 )) 2168 2169 self.add_col(am.ArrayConf('coords', default=[0.0, 0.0, 0.0], 2170 groupnames=['parameters', ], 2171 dtype=np.float32, 2172 perm='rw', 2173 name='Coordinate', 2174 symbol='x,y,z', 2175 unit='m', 2176 info='Local 3D coordinate of point', 2177 )) 2178 2179 self.add_col(am.ArrayConf('timestamps', default=0.0, 2180 dtype=np.float, 2181 groupnames=['parameters', ], 2182 perm='r', 2183 name='timestamp', 2184 symbol='t', 2185 metatype='datetime', 2186 unit='s', 2187 digits_fraction=2, 2188 info='Time stamp of point in seconds after 01 January 1970.', 2189 )) 2190 # test: 2191 self.timestamps.metatype = 'datetime' 2192 self.timestamps.set_perm('r') 2193 2194 def set_trips(self, trips): 2195 self.add_col(am.IdsArrayConf('ids_trip', trips, 2196 groupnames=['state'], 2197 name='Trip ID', 2198 info='ID of trips to which this point belongs to.', 2199 )) 2200 # put in geometry filter 2201 # self.add_col(am.ArrayConf( 'are_inside_boundary', default =False, 2202 # groupnames = ['parameters',], 2203 # perm='r', 2204 # name = 'in boundary', 2205 # info = 'True if this the data point is within the boundaries of the road network.', 2206 # )) 2207 2208 def get_ids_selected(self): 2209 """ 2210 Returns point ids of selected traces 2211 """ 2212 # print 'GpsPoints.get_ids_selected' 2213 # print ' ??ids_points = ',self.select_ids(self.parent.trips.are_selected[self.ids_trip.get_value()] ) 2214 # TODO: why is this working??? do we need trips.ids_points???? 2215 return self.select_ids(self.parent.trips.are_selected[self.ids_trip.get_value()]) 2216 # return self.get_ids()#self.select_ids(self.get_ids() 2217 # return self.select_ids( 2218 2219 def get_scenario(self): 2220 return self.parent.get_scenario() 2221 2222 def get_coords(self): 2223 """ 2224 Returns an array of x,y coordinates of all points. 2225 """ 2226 return self.coords.get_value() 2227 2228 def get_times(self): 2229 """ 2230 Returns an array of time stamps of all points. 2231 """ 2232 return self.timestamp.get_value() 2233 2234 def get_interval(self): 2235 if len(self) == 0: 2236 return [0.0, 0.0] 2237 else: 2238 timestamps = self.timestamps.get_value() 2239 return [timestamps.min(), timestamps.max()] 2240 2241 def get_duration(self): 2242 ts, te = self.get_interval() 2243 return te-ts 2244 2245 def get_distance(self): 2246 v = self.get_coords() 2247 # return np.linalg.norm( self.cols.coords[0:-1] - self.cols.coords[1:] ) 2248 return np.sum(np.sqrt((v[1:, 0]-v[0:-1, 0])**2 + (v[1:, 1]-v[0:-1, 1])**2)) 2249 2250 def get_boundary(self): 2251 if len(self) == 0: 2252 return [0, 0, 0, 0] 2253 else: 2254 x_min, y_min = self.get_coords().min(0) 2255 x_max, y_max = self.get_coords().max(0) 2256 2257 return [x_min, y_min, x_max, y_max] 2258 2259 def project(self, proj=None, offset=None): 2260 if proj is None: 2261 proj, offset = self.parent.get_proj_and_offset() 2262 x, y = proj(self.longitudes.get_value(), self.latitudes.get_value()) 2263 self.get_coords()[:] = np.transpose(np.concatenate( 2264 ([x+offset[0]], [y+offset[1]], [self.altitudes.get_value()]), axis=0)) 2265 2266 2267class GpsPersons(am.ArrayObjman): 2268 2269 def __init__(self, ident, mapmatching, **kwargs): 2270 # print 'PrtVehicles vtype id_default',vtypes.ids_sumo.get_id_from_index('passenger1') 2271 self._init_objman(ident=ident, 2272 parent=mapmatching, 2273 name='GPS Persons', 2274 info='GPS person database.', 2275 **kwargs) 2276 2277 self._init_attributes() 2278 2279 def _init_attributes(self): 2280 trips = self.parent.trips 2281 2282 # TODO: add/update vtypes here 2283 self.add_col(SumoIdsConf('User', xmltag='id')) 2284 2285 self.add_col(am.IdlistsArrayConf('ids_trips', trips, 2286 #groupnames = ['_private'], 2287 name='Trip IDs', 2288 info="IDs of trips made by this vehicle. This is a collection of recorded trips associated with this person.", 2289 )) 2290 2291 self.add_col(am.ArrayConf('ids_genders', default=0, 2292 dtype=np.int, 2293 groupnames=['parameters'], 2294 name='gender', 2295 info='Gender of person.', 2296 )) 2297 2298 self.add_col(am.ArrayConf('years_birth', default=-1, 2299 dtype=np.int, 2300 groupnames=['parameters'], 2301 name='birth year', 2302 info='Year when person has been born.', 2303 )) 2304 2305 self.add_col(am.ArrayConf('ids_occupation', default=OCCUPATIONS['unknown'], 2306 dtype=np.int32, 2307 choices=OCCUPATIONS, 2308 groupnames=['parameters'], 2309 name='occupation', 2310 info='Tupe of occupation', 2311 )) 2312 2313 self.add_col(am.ArrayConf('are_frequent_user', False, 2314 dtype=np.int32, 2315 groupnames=['parameters'], 2316 name='frequent user', 2317 info='If true, this person is a frequent user of the recorded transport mode.', 2318 )) 2319 2320 self.add_col(am.ArrayConf('zips', -1, 2321 dtype=np.int32, 2322 groupnames=['parameters'], 2323 name='ZIP', 2324 info='ZIP code of persons home.', 2325 )) 2326 2327 def make(self, id_sumo, **kwargs): 2328 if self.ids_sumo.has_index(id_sumo): 2329 id_pers = self.ids_sumo.get_id_from_index(id_sumo) 2330 #self.set_row(id_pers, **kwargs) 2331 return id_pers 2332 else: 2333 id_pers = self.add_row(ids_sumo=id_sumo, 2334 ids_trips=kwargs.get('ids_trip', None), 2335 ids_genders=kwargs.get('id_gender', None), 2336 years_birth=kwargs.get('year_birth', None), 2337 ids_occupation=kwargs.get('id_occupation', None), 2338 are_frequent_user=kwargs.get('is_frequent_user', None), 2339 zips=kwargs.get('zip', None), 2340 ) 2341 2342 return id_pers 2343 2344 2345class Mapmatching(DemandobjMixin, cm.BaseObjman): 2346 def __init__(self, ident, demand=None, 2347 name='Mapmatching', info='Mapmatching functionality.', 2348 **kwargs): 2349 2350 self._init_objman(ident=ident, parent=demand, 2351 name=name, info=info, **kwargs) 2352 2353 attrsman = self.set_attrsman(cm.Attrsman(self)) 2354 2355 self._init_attributes() 2356 self._init_constants() 2357 2358 def get_scenario(self): 2359 return self.parent.parent 2360 2361 def clear_all(self): 2362 self.trips.clear() 2363 self.points.clear() 2364 self.persons.clear() 2365 self._init_constants() 2366 2367 def clear_routes(self): 2368 self.trips.clear_routes() 2369 2370 def _init_attributes(self): 2371 print 'Mapmatching._init_attributes' 2372 attrsman = self.get_attrsman() 2373 2374 self.points = attrsman.add(cm.ObjConf(GpsPoints('points', self))) 2375 self.trips = attrsman.add(cm.ObjConf(GpsTrips('trips', self))) 2376 self.points.set_trips(self.trips) 2377 self.persons = attrsman.add(cm.ObjConf(GpsPersons('persons', self))) 2378 2379 def _init_constants(self): 2380 self._proj = None 2381 self._segvertices_xy = None 2382 self._distancesmap = None 2383 self._accesslevelsmap = None 2384 2385 attrsman = self.get_attrsman() 2386 attrsman.do_not_save_attrs(['_segvertices_xy', '_proj', '_distancesmap', '_accesslevelsmap']) 2387 2388 def get_proj_and_offset(self): 2389 if self._proj is None: 2390 net = self.get_scenario().net 2391 proj_params = str(net.get_projparams()) 2392 # try: 2393 self._proj = pyproj.Proj(proj_params) 2394 # except: 2395 # proj_params ="+proj=utm +zone=32 +ellps=WGS84 +datum=WGS84 +units=m +no_defs" 2396 # self._proj = pyproj.Proj(self.proj_params) 2397 2398 self._offset = net.get_offset() 2399 2400 return self._proj, self._offset 2401 2402 def get_segvertices_xy(self): 2403 if self._segvertices_xy is None: 2404 self._segvertices_xy = self.get_scenario().net.edges.get_segvertices_xy() 2405 2406 return self._segvertices_xy 2407 2408 def get_distancesmap(self, is_check_lanes=False): 2409 """ 2410 Returns a dictionary where key is id_mode and 2411 value is a distance-lookup table, mapping id_edge to edge distance 2412 """ 2413 # print 'get_distancesmap',self._distancesmap is None 2414 2415 if self._distancesmap is None: 2416 2417 vtypes = self.get_scenario().demand.vtypes 2418 edges = self.get_scenario().net.edges 2419 ids_vtype = self.trips.ids_vtype[self.trips.get_ids_selected()] 2420 ids_mode = vtypes.ids_mode[ids_vtype] 2421 # print ' ids_mode',ids_mode 2422 2423 self._distancesmap = {} 2424 for id_mode in set(ids_mode): 2425 #ids_vtype_mode = vtypes.select_by_mode(id_mode) 2426 2427 self._distancesmap[id_mode] = edges.get_distances(id_mode=id_mode, 2428 is_check_lanes=is_check_lanes, 2429 ) 2430 2431 # print ' len(self._distancesmap)',len(self._distancesmap) 2432 return self._distancesmap 2433 2434 def get_accesslevelsmap(self): 2435 """ 2436 Returns a dictionary where key is id_mode and 2437 value is a distance-lookup table, mapping id_edge to edge distance 2438 """ 2439 2440 if self._accesslevelsmap is None: 2441 vtypes = self.get_scenario().demand.vtypes 2442 edges = self.get_scenario().net.edges 2443 ids_vtype = self.trips.ids_vtype[self.trips.get_ids_selected()] 2444 ids_mode = vtypes.ids_mode[ids_vtype] 2445 2446 self._accesslevelsmap = {} 2447 for id_mode in set(ids_mode): 2448 #ids_vtype_mode = vtypes.select_by_mode(id_mode) 2449 self._accesslevelsmap[id_mode] = edges.get_accesslevels(id_mode) 2450 return self._accesslevelsmap 2451 2452 2453class Matchresults(cm.BaseObjman): 2454 def __init__(self, ident, mapmatching, 2455 name='Mapmatching results', 2456 info='Results of mapmatching analysis.', 2457 **kwargs): 2458 2459 # make results a child of process or of wxgui 2460 # use these objects to access matched trips 2461 2462 self._init_objman(ident, parent=mapmatching, name=name, 2463 info=info, **kwargs) 2464 attrsman = self.set_attrsman(cm.Attrsman(self)) 2465 2466 self._init_attributes() 2467 2468 def _init_attributes(self): 2469 attrsman = self.get_attrsman() 2470 mapmatching = self.parent 2471 # self.routesresults = attrsman.add(cm.ObjConf( Routesresults('routesresults', 2472 # self, mapmatching.trips.routes), 2473 # groupnames = ['Route results'], 2474 # )) 2475 2476 # add trip results from all demand objects 2477 print 'Matchresults._init_attributes' 2478 2479 def config(self, resultobj, **kwargs): 2480 # attention: need to check whether already set 2481 # because setattr is set explicitely after add 2482 if not hasattr(self, resultobj.get_ident()): 2483 if kwargs.has_key('groupnames'): 2484 kwargs['groupnames'].append('Results') 2485 else: 2486 kwargs['groupnames'] = ['Results'] 2487 attrsman = self.get_attrsman() 2488 attrsman.add(cm.ObjConf(resultobj, **kwargs)) 2489 setattr(self, resultobj.get_ident(), resultobj) 2490 2491 def get_scenario(self): 2492 return self.parent.get_scenario() 2493 2494 def clear_all(self): 2495 self.clear() 2496 2497 def save(self, filepath=None, is_not_save_parent=True): 2498 if filepath is None: 2499 self.get_scenario().get_rootfilepath()+'.mmatch.obj' 2500 # parent will not be saved because in no_save set 2501 cm.save_obj(self, filepath, is_not_save_parent=is_not_save_parent) 2502 2503 2504class Nodeinfolists(am.ListArrayConf): 2505 def format_value(self, _id, show_unit=False, show_parentesis=False): 2506 if show_unit: 2507 unit = ' '+self.format_unit(show_parentesis) 2508 else: 2509 unit = '' 2510 # return repr(self[_id])+unit 2511 2512 #self.min = minval 2513 #self.max = maxval 2514 #self.digits_integer = digits_integer 2515 #self.digits_fraction = digits_fraction 2516 val = self[_id] 2517 tt = type(val) 2518 2519 if tt in (np.int, np.int32, np.float64): 2520 return str(val)+unit 2521 2522 elif tt in (np.float, np.float32, np.float64): 2523 if hasattr(self, 'digits_fraction'): 2524 digits_fraction = self.digits_fraction 2525 else: 2526 digits_fraction = 3 2527 s = "%."+str(digits_fraction)+"f" 2528 return s % (val)+unit 2529 2530 else: 2531 return str(val)+unit 2532 2533 2534class Accesslevellists(am.ListArrayConf): 2535 def format_value(self, _id, show_unit=False, show_parentesis=False): 2536 if show_unit: 2537 unit = ' '+self.format_unit(show_parentesis) 2538 else: 2539 unit = '' 2540 # return repr(self[_id])+unit 2541 2542 #self.min = minval 2543 #self.max = maxval 2544 #self.digits_integer = digits_integer 2545 #self.digits_fraction = digits_fraction 2546 val = self[_id] 2547 tt = type(val) 2548 2549 if tt in (np.int, np.int32, np.float64): 2550 return str(val)+unit 2551 2552 elif tt in (np.float, np.float32, np.float64): 2553 if hasattr(self, 'digits_fraction'): 2554 digits_fraction = self.digits_fraction 2555 else: 2556 digits_fraction = 3 2557 s = "%."+str(digits_fraction)+"f" 2558 return s % (val)+unit 2559 2560 else: 2561 return str(val)+unit 2562 2563 2564class Edgesresults(am.ArrayObjman): 2565 def __init__(self, ident, parent, 2566 name='Edges results', 2567 info='Table with results from edges that are part of matched routes or alternative routes.', 2568 **kwargs): 2569 2570 self._init_objman(ident=ident, 2571 parent=parent, # main results object 2572 info=info, 2573 name=name, 2574 **kwargs) 2575 2576 # self.add(cm.AttrConf( 'datapathkey',datapathkey, 2577 # groupnames = ['_private'], 2578 # name = 'data pathkey', 2579 # info = "key of data path", 2580 # )) 2581 2582 self.add_col(am.IdsArrayConf('ids_edge', parent.get_scenario().net.edges, 2583 groupnames=['state'], 2584 is_index=True, 2585 name='Edge ID', 2586 info='ID of network edge.', 2587 )) 2588 self._init_attributes() 2589 2590 def _init_attributes(self): 2591 2592 self.add_col(am.ArrayConf('speeds_average', default=0.0, 2593 dtype=np.float32, 2594 groupnames=['results'], 2595 name='Avg. speed', 2596 unit='m/s', 2597 info='Average speed on this edge.', 2598 )) 2599 2600 self.add_col(am.ArrayConf('durations_tot_matched', default=0.0, 2601 dtype=np.float32, 2602 groupnames=['results'], 2603 name='Tot. duration', 2604 unit='s', 2605 info='Total time of matched routes spent on this edge.', 2606 )) 2607 2608 self.add_col(am.ArrayConf('numbers_tot_matched', default=0, 2609 dtype=np.int32, 2610 groupnames=['results'], 2611 name='number matched', 2612 info='Total number of matched routes crossing this edge.', 2613 )) 2614 2615 self.add_col(am.ArrayConf('numbers_tot_shortest', default=0, 2616 dtype=np.int32, 2617 groupnames=['results'], 2618 name='number shortest', 2619 info='Total number of shortest routes crossing this edge.', 2620 )) 2621 2622 self.add_col(am.ArrayConf('differences_dist_tot_shortest', default=0.0, 2623 dtype=np.float32, 2624 groupnames=['results'], 2625 name='Tot. length difference', 2626 unit='m', 2627 info='Sum of length differences between matched route length and the corrisponding shortest route using this edge.', 2628 )) 2629 2630 def init_for_routes(self, routes): 2631 """ 2632 Initializes a row for each edge ID in routes 2633 """ 2634 ids_edge_init = set() 2635 for ids_edge in routes.ids_edges.get_value(): 2636 if ids_edge is not None: 2637 ids_edge_init.update(ids_edge) 2638 ids_edge_init = list(ids_edge_init) 2639 # return result ids and respective edge ids 2640 ids_edgeresmap = np.zeros(np.max(ids_edge_init)+1, dtype=np.int32) 2641 ids_edgeresmap[ids_edge_init] = self.add_rows(n=len(ids_edge_init), ids_edge=ids_edge_init) 2642 2643 return ids_edgeresmap 2644 2645 2646class Routesresults(am.ArrayObjman): 2647 def __init__(self, ident, parent, 2648 name='Route results', 2649 info='Table with results from analysis of each matched route.', 2650 **kwargs): 2651 2652 self._init_objman(ident=ident, 2653 parent=parent, # main results object 2654 info=info, 2655 name=name, 2656 **kwargs) 2657 2658 # self.add(cm.AttrConf( 'datapathkey',datapathkey, 2659 # groupnames = ['_private'], 2660 # name = 'data pathkey', 2661 # info = "key of data path", 2662 # )) 2663 2664 self.add_col(am.IdsArrayConf('ids_route', parent.parent.trips.routes, 2665 groupnames=['state'], 2666 is_index=True, 2667 name='ID route', 2668 info='ID of route.', 2669 )) 2670 self._init_attributes() 2671 2672 def _init_attributes(self): 2673 2674 self.add_col(am.ArrayConf('distances', default=-1.0, 2675 dtype=np.float32, 2676 groupnames=['results'], 2677 name='distance', 2678 unit='m', 2679 info='Length of the route.', 2680 )) 2681 2682 self.add_col(am.ArrayConf('durations', default=0.0, 2683 dtype=np.float32, 2684 groupnames=['results'], 2685 name='duration', 2686 unit='s', 2687 info='Time duration of the route.', 2688 )) 2689 2690 self.add_col(am.ArrayConf('lengths_mixed', default=-1.0, 2691 dtype=np.float32, 2692 groupnames=['results'], 2693 name='Length mixed access', 2694 #symbol = 'L_{MIX}', 2695 unit='m', 2696 info='Length of roads with mixed access.', 2697 )) 2698 2699 self.add_col(am.ArrayConf('lengths_exclusive', default=-1.0, 2700 dtype=np.float32, 2701 groupnames=['results'], 2702 name='Length exclusive access', 2703 #symbol = 'L_{EXC}', 2704 unit='m', 2705 info='Length of roads with exclusive access.', 2706 )) 2707 2708 self.add_col(am.ArrayConf('lengths_low_priority', default=-1.0, 2709 dtype=np.float32, 2710 groupnames=['results'], 2711 name='Length low priority', 2712 #symbol = 'L_{EXC}', 2713 unit='m', 2714 info='Length of low priority roads. These are roads with either speed limit to 30 or below, or classified as residential, or exclusive bike or pedestrian ways. This correspond to SUMO priorities 1-6.', 2715 )) 2716 2717 self.add_col(am.ArrayConf('lengths_overlap_matched', default=-1.0, 2718 dtype=np.float32, 2719 groupnames=['results'], 2720 name='Length overlap', 2721 #symbol = 'L_{EXC}', 2722 unit='m', 2723 info='Length of overlap with matched route.', 2724 )) 2725 2726 self.add_col(am.ArrayConf('numbers_nodes', default=-1, 2727 dtype=np.int32, 2728 groupnames=['results'], 2729 name='number of nodes', 2730 #symbol = 'L_{EXC}', 2731 info='Total number of nodes.', 2732 )) 2733 2734 self.add_col(am.ArrayConf('numbers_nodes_tls', default=-1, 2735 dtype=np.int32, 2736 groupnames=['results'], 2737 name='number of TL nodes', 2738 #symbol = 'L_{EXC}', 2739 info='Total number of traffic light controlled nodes.', 2740 )) 2741 2742 # self.add_col(am.ListArrayConf( 'intervallists', 2743 # groupnames = ['_private'], 2744 # perm='rw', 2745 # name = 'Intervals', 2746 # unit = 's', 2747 # info = 'Time interval as tuple for each edge of the route. These time intervals are obtained from timestamps of GPS points which can be perpendicularly projected to the respective route edge.', 2748 # )) 2749 2750 # self.add_col(am.ListArrayConf( 'edgedurationlists', 2751 # groupnames = ['_private'], 2752 # perm='rw', 2753 # name = 'Edge durations', 2754 # unit = 's', 2755 # info = 'Duration on each edge of the route. These times are obtained from timestamps of GPS points which can be perpendicularly projected to the respective route edge.', 2756 # )) 2757 2758 # self.add_col(am.ListArrayConf( 'nodedurationlists', 2759 # groupnames = ['_private'], 2760 # perm='rw', 2761 # name = 'Node durations', 2762 # unit = 's', 2763 # info = 'Duration on each intersection of the route. These times are obtained from timestamps of GPS points on the edges before and after the respective node.', 2764 # )) 2765 2766 # self.add_col(Nodeinfolists( 'nodeinfolists', 2767 # groupnames = ['_private'], 2768 # perm='rw', 2769 # name = 'Node info', 2770 # info = 'Node info for each node on the route. Node info is defined as (node type, number of legs, accesslevel) and for traffic light nodes (node type, total cycle time, green time).', 2771 # )) 2772 2773 # self.nodetypemap = {\ 2774 # "priority":0, 2775 # "traffic_light":1, 2776 # "right_before_left":2, 2777 # "unregulated":3, 2778 # "priority_stop":4, 2779 # "traffic_light_unregulated":5, 2780 # "allway_stop":6, 2781 # "rail_signal":7, 2782 # "zipper":8, 2783 # "traffic_light_right_on_red":9, 2784 # "rail_crossing":10, 2785 # "dead_end":11, 2786 # } 2787 2788 # self.add_col(Accesslevellists( 'accesslevellists', 2789 # groupnames = ['_private'], 2790 # perm='rw', 2791 # name = 'Accesslevels', 2792 # info = 'Accesslevels for each edge of the route. -1 = no access, 0 = all have access, 1 = reserved access together with other modes, 2 = exclusive access.', 2793 # )) 2794 2795 2796class Shortestrouter(Process): 2797 def __init__(self, ident, mapmatching, logger=None, **kwargs): 2798 print 'Shortestrouter.__init__' 2799 2800 # TODO: let this be independent, link to it or child?? 2801 2802 self._init_common(ident, 2803 parent=mapmatching, 2804 name='Shortest Router', 2805 logger=logger, 2806 info='Shortest path router.', 2807 ) 2808 2809 attrsman = self.set_attrsman(cm.Attrsman(self)) 2810 2811 self.c_modespecific = attrsman.add(cm.AttrConf('c_modespecific', kwargs.get('c_modespecific', 0.9), 2812 groupnames=['options'], 2813 perm='rw', 2814 name='Mode constant', 2815 info='The Mode constant is multiplied with the the distance of each edge if vehicle has exclusive access.', 2816 )) 2817 2818 self.dist_modespecific = attrsman.add(cm.AttrConf('dist_modespecific', kwargs.get('dist_modespecific', 0.0), 2819 groupnames=['options'], 2820 perm='rw', 2821 name='Mode dist', 2822 unit='m', 2823 info='The mode distance is multiplied with the accesslevel and added to the distance of each edge.', 2824 )) 2825 2826 self.is_ignor_connections = attrsman.add(cm.AttrConf('is_ignor_connections', kwargs.get('is_ignor_connections', True), 2827 groupnames=['options'], 2828 perm='rw', 2829 name='Ignore connections', 2830 info='Ignore connections means that the vehicle can always make turns into all possible edges, whether they are allowed or not.', 2831 )) 2832 2833 def do(self): 2834 print 'Shortestrouter.do' 2835 # links 2836 mapmatching = self.parent 2837 trips = mapmatching.trips 2838 trips.route_shortest(dist_modespecific=self.dist_modespecific, 2839 c_modespecific=self.c_modespecific, 2840 is_ignor_connections=self.is_ignor_connections, 2841 ) 2842 2843 2844class Routesanalyzer(Process): 2845 def __init__(self, ident, mapmatching, results=None, logger=None, **kwargs): 2846 print 'Routesanalyzer.__init__' 2847 2848 # TODO: let this be independent, link to it or child?? 2849 2850 if results is None: 2851 self._results = Matchresults('matchresults', mapmatching) 2852 else: 2853 self._results = results 2854 2855 self._init_common(ident, 2856 parent=mapmatching, 2857 name='Routes Analyzer', 2858 logger=logger, 2859 info='Determines attributes of matched route and alternative routes.', 2860 ) 2861 2862 attrsman = self.set_attrsman(cm.Attrsman(self)) 2863 2864 self.priority_max_low = attrsman.add(cm.AttrConf('priority_max_low', kwargs.get('priority_max_low', 6), 2865 groupnames=['options'], 2866 perm='rw', 2867 name='Max. low priority', 2868 info='Maximum priority of an edge, such that it is still considered low priority.', 2869 )) 2870 2871 self.timeloss_intersection = attrsman.add(cm.AttrConf('timeloss_intersection', kwargs.get('timeloss_intersection', 5.0), 2872 groupnames=['options'], 2873 perm='rw', 2874 name='Timeloss intersection', 2875 info='Estimated timeloss at intersections due to waiting or reduction of speed. This value is used to estimate travel time of route alternatives.', 2876 )) 2877 2878 self.timeloss_tl = attrsman.add(cm.AttrConf('timeloss_tl', kwargs.get('timeloss_tl', 15.0), 2879 groupnames=['options'], 2880 perm='rw', 2881 name='Timeloss TL', 2882 info='Additional estimated timeloss at traffic light intersections due to waiting at the red light. This value is used to estimate travel time of route alternatives.', 2883 )) 2884 2885 self.edgedurations_tot_min = attrsman.add(cm.AttrConf('edgedurations_tot_min', kwargs.get('edgedurations_tot_min', 10.0), 2886 groupnames=['options'], 2887 perm='rw', 2888 name='min tot. edge duration', 2889 info='Minimum total duration of all users in an edge in order to calculate average speed. A too small value distors the results.', 2890 )) 2891 2892 results = self.get_results() 2893 results.config(Routesresults('routesresults_matched', results, 2894 name='Results of matched routes'), name='Results of matched routes') 2895 results.config(Routesresults('routesresults_shortest', results, 2896 name='Results of shortest routes'), name='Results of shortest routes') 2897 results.config(Routesresults('routesresults_matched_nonoverlap', results, name='Results of matched, non-overlapping routes'), 2898 name='Results of matched routes using only edges which do not overlap with the shortest route.') 2899 results.config(Routesresults('routesresults_shortest_nonoverlap', results, name='Results of shortest, non-overlapping routes'), 2900 name='Results of shortest routes using only edges which do not overlap with the matched route.') 2901 2902 results.config(Edgesresults('edgesresults', results)) 2903 2904 def get_results(self): 2905 return self._results 2906 2907 def do(self): 2908 print 'Routesanalyzer.do' 2909 # results 2910 results = self.get_results() 2911 routesresults_shortest = results.routesresults_shortest 2912 routesresults_matched = results.routesresults_matched 2913 2914 routesresults_shortest_nonoverlap = results.routesresults_shortest_nonoverlap 2915 routesresults_matched_nonoverlap = results.routesresults_matched_nonoverlap 2916 2917 edgesresults = results.edgesresults 2918 2919 # links 2920 mapmatching = self.parent 2921 trips = mapmatching.trips 2922 2923 routes = trips.get_routes() 2924 scenario = mapmatching.get_scenario() 2925 edges = scenario.net.edges 2926 vtypes = scenario.demand.vtypes 2927 points = mapmatching.points 2928 timestamps = points.timestamps 2929 2930 priority_max_low = self.priority_max_low 2931 timeloss_intersection = self.timeloss_intersection 2932 timeloss_tl = self.timeloss_tl 2933 2934 distancesmap = mapmatching.get_distancesmap() 2935 accesslevelsmap = mapmatching.get_accesslevelsmap() 2936 priorities = edges.priorities 2937 ids_tonode = edges.ids_tonode 2938 nodetypes = scenario.net.nodes.types 2939 tlstype = nodetypes.choices['traffic_light'] 2940 ids_tls = scenario.net.nodes.ids_tls 2941 2942 ids_trip_sel = trips.get_ids_selected() 2943 2944 ids_route_sel = trips.ids_route_matched[ids_trip_sel] 2945 inds_valid = np.flatnonzero(ids_route_sel > 0) 2946 2947 ids_route = ids_route_sel[inds_valid] 2948 ids_trip = ids_trip_sel[inds_valid] 2949 2950 ids_route_shortest = trips.ids_route_shortest[ids_trip] 2951 2952 ids_vtype = trips.ids_vtype[ids_trip] 2953 ids_mode = vtypes.ids_mode[ids_vtype] 2954 2955 #ind = 0 2956 2957 print ' analyzing %d trips' % len(ids_trip) 2958 2959 routesresults_shortest.clear() 2960 routesresults_matched.clear() 2961 routesresults_shortest_nonoverlap.clear() 2962 routesresults_matched_nonoverlap.clear() 2963 2964 edgesresults.clear() 2965 ids_res = routesresults_matched.add_rows(n=len(ids_trip)) 2966 routesresults_shortest.add_rows(ids=ids_res) 2967 routesresults_matched_nonoverlap.add_rows(ids=ids_res) 2968 routesresults_shortest_nonoverlap.add_rows(ids=ids_res) 2969 2970 ids_edgeresmap = edgesresults.init_for_routes(routes) 2971 2972 for id_trip, id_route_matched, id_route_shortest, id_mode, id_res, dist_gps, duration_gps\ 2973 in zip(ids_trip, 2974 ids_route, 2975 ids_route_shortest, 2976 ids_mode, 2977 ids_res, 2978 trips.lengths_gpsroute_matched[ids_trip], 2979 trips.durations_route_matched[ids_trip]): 2980 #ids_point = trips.ids_points[id_trip] 2981 distances = distancesmap[id_mode] 2982 2983 ids_edge_matched = np.array(routes.ids_edges[id_route_matched][1:], dtype=np.int32) 2984 ids_points_edgeend = np.array(trips.ids_points_edgeend[id_trip], dtype=np.int32) # [1:] 2985 accesslevels_matched = accesslevelsmap[id_mode][ids_edge_matched] 2986 priorities_matched = priorities[ids_edge_matched] 2987 nodetypes_matched = nodetypes[ids_tonode[ids_edge_matched]] 2988 # print ' len(ids_edge_matched)',len(ids_edge_matched) 2989 # print ' len(ids_points_edgeend)',len(ids_points_edgeend) 2990 # print ' ids_edge_matched',ids_edge_matched 2991 # print ' accesslevels_matched',accesslevels_matched 2992 # print ' tlstype',tlstype 2993 # print ' accesslevels_mixed',accesslevels_matched==1 2994 # print ' accesslevels_mixed',np.flatnonzero(accesslevels_matched==1) 2995 # print ' distances', 2996 # print ' nodetypes_matched',nodetypes_matched 2997 dist_matched = np.sum(distances[ids_edge_matched]) 2998 duration_matched = timestamps[ids_points_edgeend[-1]] - timestamps[ids_points_edgeend[0]] 2999 3000 n_nodes_matched = len(nodetypes_matched) 3001 n_tls_matched = np.sum(nodetypes_matched == tlstype) 3002 lineduration_matched = duration_matched-n_nodes_matched*timeloss_intersection-n_tls_matched*timeloss_tl 3003 linespeed_matched = dist_matched/lineduration_matched 3004 3005 # do edge by edge analyses if at least more than 1 edge 3006 if len(ids_edge_matched) > 1: 3007 # cut off last edge 3008 3009 # print ' ids_edge_matched',ids_edge_matched.shape,ids_edge_matched 3010 # print ' ids_points_edgeend',ids_points_edgeend.shape,ids_points_edgeend 3011 ids_edge = ids_edge_matched[:-1] 3012 ids_points_ee = ids_points_edgeend # [:-1] 3013 3014 # print ' ids_edge',ids_edge.shape 3015 # print ' ids_points_ee',ids_points_ee[1:].shape,ids_points_ee[1:] 3016 durations_edge = timestamps[ids_points_ee[1:]] - timestamps[ids_points_ee[:-1]] 3017 # print ' durations_edge',durations_edge.shape,durations_edge 3018 edgesresults.numbers_tot_matched[ids_edgeresmap[ids_edge]] += 1 3019 edgesresults.durations_tot_matched[ids_edgeresmap[ids_edge]] += durations_edge 3020 #dists_edge = distances[ids_edge] 3021 # for id_edge, duration_edge, dists_edge, id_point_from, id_point_to in zip(ids_edge,durations_edge,ids_points_ee[:1],ids_points_ee[1:]): 3022 # edgesresults.numbers_tot_matched[id_edge] += 1 3023 # edgesresults.durations_tot_matched[id_edge] += duration_edge 3024 # identify gps points with id_point_from, id_point_to 3025 # and compute line speed, waiting time.... 3026 3027 routesresults_matched.set_row(id_res, 3028 ids_route=id_route_matched, 3029 distances=dist_matched, 3030 durations=duration_matched, 3031 lengths_mixed=np.sum(distances[ids_edge_matched[accesslevels_matched == 1]]), 3032 lengths_exclusive=np.sum( 3033 distances[ids_edge_matched[accesslevels_matched == 2]]), 3034 lengths_low_priority=np.sum( 3035 distances[ids_edge_matched[priorities_matched <= priority_max_low]]), 3036 lengths_overlap_matched=np.sum(distances[ids_edge_matched]), 3037 numbers_nodes=n_nodes_matched, 3038 numbers_nodes_tls=n_tls_matched, 3039 ) 3040 3041 if id_route_shortest >= 0: 3042 ids_edge_shortest = np.array(routes.ids_edges[id_route_shortest], dtype=np.int32)[1:] 3043 accesslevels_shortest = accesslevelsmap[id_mode][ids_edge_shortest] 3044 priorities_shortest = priorities[ids_edge_shortest] 3045 nodetypes_shortest = nodetypes[ids_tonode[ids_edge_shortest]] 3046 3047 ids_edge_overlap = list(set(ids_edge_matched).intersection(ids_edge_shortest)) 3048 ids_edge_shortest_nonoverlap = np.array( 3049 list(set(ids_edge_shortest).difference(ids_edge_matched)), dtype=np.int32) 3050 ids_edge_match_nonoverlap = np.array( 3051 list(set(ids_edge_matched).difference(ids_edge_shortest)), dtype=np.int32) 3052 dist_shortest = np.sum(distances[ids_edge_shortest]) 3053 n_nodes_shortest = len(nodetypes_shortest) 3054 n_tls_shortest = np.sum(nodetypes_shortest == tlstype) 3055 3056 edgesresults.numbers_tot_shortest[ids_edgeresmap[ids_edge_shortest]] += 1 3057 edgesresults.differences_dist_tot_shortest[ids_edgeresmap[ids_edge_shortest_nonoverlap] 3058 ] += dist_matched-dist_shortest 3059 3060 routesresults_shortest.set_row(id_res, 3061 ids_route=id_route_shortest, 3062 distances=dist_shortest, 3063 durations=dist_shortest/linespeed_matched+timeloss_intersection * n_nodes_shortest + timeloss_tl * n_tls_shortest, 3064 lengths_mixed=np.sum( 3065 distances[ids_edge_shortest[accesslevels_shortest == 1]]), 3066 lengths_exclusive=np.sum( 3067 distances[ids_edge_shortest[accesslevels_shortest == 2]]), 3068 lengths_low_priority=np.sum( 3069 distances[ids_edge_shortest[priorities_shortest <= priority_max_low]]), 3070 numbers_nodes=n_nodes_shortest, 3071 numbers_nodes_tls=n_tls_shortest, 3072 lengths_overlap_matched=np.sum(distances[ids_edge_overlap]), 3073 ) 3074 3075 accesslevels_nonoverlap = accesslevelsmap[id_mode][ids_edge_match_nonoverlap] 3076 # print ' ids_edge_match_nonoverlap',ids_edge_match_nonoverlap 3077 # print ' accesslevels_nonoverlap==1',np.flatnonzero(accesslevels_nonoverlap==1),np.flatnonzero(accesslevels_nonoverlap==1).dtype 3078 # print ' ids_edge_match_nonoverlap[accesslevels_nonoverlap==1]',ids_edge_match_nonoverlap[np.flatnonzero(accesslevels_nonoverlap==1)] 3079 priorities_nonoverlap = priorities[ids_edge_match_nonoverlap] 3080 nodetypes_nonoverlap = nodetypes[ids_tonode[ids_edge_match_nonoverlap]] 3081 3082 dist_nonoverlap = np.sum(distances[ids_edge_match_nonoverlap]) 3083 n_nodes_nonoverlap = len(nodetypes_nonoverlap) 3084 n_tls_nonoverlap = np.sum(nodetypes_nonoverlap == tlstype) 3085 3086 routesresults_matched_nonoverlap.set_row(id_res, 3087 ids_route=id_route_matched, 3088 distances=dist_nonoverlap, 3089 durations=dist_nonoverlap/linespeed_matched+timeloss_intersection * n_nodes_nonoverlap + timeloss_tl * n_tls_nonoverlap, 3090 lengths_mixed=np.sum( 3091 distances[ids_edge_match_nonoverlap[accesslevels_nonoverlap == 1]]), 3092 lengths_exclusive=np.sum( 3093 distances[ids_edge_match_nonoverlap[accesslevels_nonoverlap == 2]]), 3094 lengths_low_priority=np.sum( 3095 distances[ids_edge_match_nonoverlap[priorities_nonoverlap <= priority_max_low]]), 3096 lengths_overlap_matched=np.sum( 3097 distances[ids_edge_match_nonoverlap]), 3098 numbers_nodes=n_nodes_nonoverlap, 3099 numbers_nodes_tls=n_tls_nonoverlap, 3100 ) 3101 3102 accesslevels_nonoverlap = accesslevelsmap[id_mode][ids_edge_shortest_nonoverlap] 3103 priorities_nonoverlap = priorities[ids_edge_shortest_nonoverlap] 3104 nodetypes_nonoverlap = nodetypes[ids_tonode[ids_edge_shortest_nonoverlap]] 3105 3106 dist_nonoverlap = np.sum(distances[ids_edge_shortest_nonoverlap]) 3107 n_nodes_nonoverlap = len(nodetypes_nonoverlap) 3108 n_tls_nonoverlap = np.sum(nodetypes_nonoverlap == tlstype) 3109 3110 routesresults_shortest_nonoverlap.set_row(id_res, 3111 ids_route=id_route_matched, 3112 distances=dist_nonoverlap, 3113 durations=dist_nonoverlap/linespeed_matched+timeloss_intersection * n_nodes_nonoverlap + timeloss_tl * n_tls_nonoverlap, 3114 lengths_mixed=np.sum( 3115 distances[ids_edge_shortest_nonoverlap[accesslevels_nonoverlap == 1]]), 3116 lengths_exclusive=np.sum( 3117 distances[ids_edge_shortest_nonoverlap[accesslevels_nonoverlap == 2]]), 3118 lengths_low_priority=np.sum( 3119 distances[ids_edge_shortest_nonoverlap[priorities_nonoverlap <= priority_max_low]]), 3120 lengths_overlap_matched=np.sum( 3121 distances[ids_edge_shortest_nonoverlap]), 3122 numbers_nodes=n_nodes_nonoverlap, 3123 numbers_nodes_tls=n_tls_nonoverlap, 3124 ) 3125 3126 ids_valid = edgesresults.select_ids(edgesresults.durations_tot_matched.get_value() > self.edgedurations_tot_min) 3127 edgesresults.speeds_average[ids_valid] = edges.lengths[edgesresults.ids_edge[ids_valid]] / \ 3128 (edgesresults.durations_tot_matched[ids_valid]/edgesresults.numbers_tot_matched[ids_valid]) 3129