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-10-preresults.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 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 150 151def calc_seconds(t_data, 152 sep_date_clock=' ', sep_date='-', sep_clock=':', 153 is_float=False): 154 """ 155 Returns time in seconds after 1/1/1970. 156 Time format for time data string used: 157 2012-05-02 12:57:08.0 158 """ 159 160 if len(t_data.split(sep_date_clock)) != 2: 161 return -1 162 (date, clock) = t_data.split(sep_date_clock) 163 164 if (len(clock.split(sep_clock)) == 3) & (len(date.split(sep_date)) == 3): 165 (year_str, month_str, day_str) = date.split(sep_date) 166 (hours_str, minutes_str, seconds_str) = clock.split(sep_clock) 167 t = time.mktime((int(year_str), int(month_str), int(day_str), 168 int(hours_str), int(minutes_str), int(float(seconds_str)), -1, -1, -1)) 169 170 # print 'calc_seconds',t 171 # print ' t_data' 172 # print ' tupel',int(year_str),int(month_str),int(day_str), int(hours_str),int(minutes_str),int(float(seconds_str)),0,0,0 173 if is_float: 174 return t 175 else: 176 return int(t) 177 else: 178 return -1 179 180 181def get_routelinestring(route, edges): 182 # print 'get_routelinestring' 183 # TODO: we could make the np.sum of shapes if shapes where lists 184 shape = [] 185 for id_edge in route: 186 # print ' ',edges.shapes[id_edge ],type(edges.shapes[id_edge ]) 187 shape += list(edges.shapes[id_edge]) 188 # print ' ',shape 189 190 return LineString(shape) 191 192 193def get_boundary(coords): 194 # print 'get_boundary',len(coords),coords.shape 195 # print ' coords',coords 196 if len(coords) == 0: 197 return [0, 0, 0, 0] 198 else: 199 x_min, y_min = coords[:, :2].min(0) 200 x_max, y_max = coords[:, :2].max(0) 201 202 return [x_min, y_min, x_max, y_max] 203 204 205class BirgilMatcher(Process): 206 def __init__(self, ident, mapmatching, logger=None, **kwargs): 207 print 'BirgilMatcher.__init__' 208 209 # TODO: let this be independent, link to it or child?? 210 211 self._init_common(ident, 212 parent=mapmatching, 213 name='Birgillito Map matching', 214 logger=logger, 215 info='Birgillito Map matching.', 216 ) 217 218 attrsman = self.set_attrsman(cm.Attrsman(self)) 219 220 self.width_buffer_max = attrsman.add(cm.AttrConf('width_buffer_max', kwargs.get('width_buffer_max', 40.0), 221 groupnames=['options'], 222 perm='rw', 223 name='Max. buffer width', 224 unit='m', 225 info='GPS points are only valid if they are within the given maximum buffer width distant from a network edge.', 226 )) 227 228 self.n_points_min = attrsman.add(cm.AttrConf('n_points_min', kwargs.get('n_points_min', 5), 229 groupnames=['options'], 230 perm='rw', 231 name='Min. point number', 232 info='Minimum number of valid GPS points. Only if this minimum number is reached, a map-matching attempt is performed.', 233 )) 234 235 # self.dist_min = attrsman.add(cm.AttrConf( 'dist_min',kwargs.get('dist_min',3.0), 236 # groupnames = ['options'], 237 # name = 'Min. dist', 238 # unit = 'm', 239 # info = 'Minimum distance used in the calculation of the edge cost function.', 240 # )) 241 242 self.weight_cumlength = attrsman.add(cm.AttrConf('weight_cumlength', kwargs.get('weight_cumlength', 0.01), 243 groupnames=['options'], 244 name='Cumlength weight', 245 info='Cost function weight of cumulative length.', 246 )) 247 248 self.weight_angle = attrsman.add(cm.AttrConf('weight_angle', kwargs.get('weight_angle', 20.0), 249 groupnames=['options'], 250 name='Angle weight', 251 info='Cost function weight of angular coincidence.', 252 )) 253 254 self.weight_access = attrsman.add(cm.AttrConf('weight_access', kwargs.get('weight_access', 5.0), 255 groupnames=['options'], 256 name='Access weight', 257 info='Cost function weight of access level for matched mode.', 258 )) 259 260 self.n_routes_follow = attrsman.add(cm.AttrConf('n_routes_follow', kwargs.get('n_routes_follow', 20), 261 groupnames=['options'], 262 name='Followed routes', 263 info='Number of routes which are followed in parallel.', 264 )) 265 266 def do(self): 267 logger = self.get_logger() 268 timeconst = 0.95 269 270 trips = self.parent.trips 271 #routes = trips.get_routes() 272 edges = self.get_edges() 273 vtypes = self.parent.get_scenario().demand.vtypes 274 275 fstar = edges.get_fstar(is_return_arrays=True, is_ignor_connections=True) 276 times = edges.get_times(id_mode=2, is_check_lanes=False, speed_max=5.555) 277 ids_trip = trips.get_ids_selected() 278 279 ids_mode = vtypes.ids_mode[trips.ids_vtype[ids_trip]] 280 281 n_trips = len(ids_trip) 282 logger.w('Start mapmatching of %d GPS traces with Birgillito method' % n_trips) 283 284 distancesmap = self.parent.get_distancesmap() 285 accesslevelsmap = self.parent.get_accesslevelsmap() 286 287 n_trip_matched = 0 288 time_match_trace_av = 3.0 # initial gues of match time 289 290 # ,vtypes.speed_max[ids_vtype] 291 for id_trip, ids_point, id_mode in zip(ids_trip, trips.ids_points[ids_trip], ids_mode): 292 tick_before = time.time() 293 route, length_route, length_route_mixed, length_route_exclusive, duration_gps, lengthindex, err_dist, t_match, ids_point_edgeend, is_connected = \ 294 self.match_trip_birgil(id_trip, ids_point, fstar, distancesmap[id_mode], accesslevelsmap[id_mode]) 295 296 trips.set_matched_route(id_trip, route, 297 length_route, 298 length_route_mixed, 299 length_route_exclusive, 300 duration_gps, 301 lengthindex, 302 err_dist, 303 t_match, 304 is_connected=is_connected, 305 ids_point_edgeend=ids_point_edgeend) 306 n_trip_matched += 1 307 308 time_match_trace = time.time()-tick_before 309 310 time_match_trace_av = timeconst*time_match_trace_av + (1.0-timeconst)*time_match_trace 311 time_end_est = (n_trips-n_trip_matched)*time_match_trace_av 312 313 progress = 100.0*n_trip_matched/float(n_trips) 314 logger.w("Matched %d/%d traces, %.2f%% comleted. Avg match time = %.2fs, Terminated in %.1fh" % 315 (n_trip_matched, n_trips, progress, time_match_trace_av, float(time_end_est)/3600.0), key='message') 316 logger.w(progress, key='progress') 317 318 def get_edges(self): 319 return self.parent.get_scenario().net.edges 320 321 def match_trip_birgil(self, id_trip, ids_point, fstar, costs, accesslevels): 322 # TODO: record time intervals 323 # calculate initial and final position on edge to get more precise 324 # triplength 325 print 79*'=' 326 print 'match_trip_birgil', id_trip 327 tick = time.time() 328 routes = [] 329 route = None 330 route_mindist = None 331 intervals_route = [] 332 length_route = -1.0 333 length_bikeway = -1.0 334 length_mindist = -1.0 335 matchindex = -1.0 336 lengthindex = -1.0 337 err_dist = -1.0 338 err_dist_alg = -1.0 339 t_match = -1.0 340 duration_gps = -1.0 341 duration_est = -1.0 342 343 points = self.parent.points 344 #trips = self.parent.trips 345 #ids_point = trips.ids_points[id_trip] 346 #pointset = self.pointsets.get(id_trip) 347 348 # create a multi point object with points of traces 349 # should no longer be necessary 350 #tracepoints = MultiPoint(pointset.get_coords().tolist()) 351 352 # make a array with time stamps of all points 353 pointtimes = points.timestamps[ids_point] 354 # segind_to_id_edge = self._segind_to_id_edge #<<<<<<<<<<< 355 coords = points.coords[ids_point][:, :2] 356 357 x1, y1, x2, y2 = self.parent.get_segvertices_xy() 358 edges = self.get_edges() 359 get_ids_edge_from_inds_seg = edges.get_ids_edge_from_inds_seg 360 get_dist_point_to_edge = edges.get_dist_point_to_edge 361 #edgehit_counts = np.zeros(len(self._net.getEdges()),int) 362 # segind_to_id_edge = self._segind_to_id_edge 363 # for p in pointset.get_coords(): 364 # print ' ',p 365 # hitinds = get_dist_point_to_segs(p,self._x1,self._y1,self._x2,self._y2, is_ending=True)< self.width_buffer**2 366 # edgehit_counts[segind_to_id_edge[hitinds]]+= 1 367 368 # --------------------------------------------------------------------- 369 # 1. initialization 370 371 # Set of initial edges 372 373 # create a list with lists for each point 374 # Each list contans the edges, where the point is in the edge buffer 375 #edges_containing_point = np.zeros(len(ids_points), object) 376 #edgedists_to_point = np.zeros(len(ids_points), object) 377 378 ids_edge_initial = [] 379 ids_edge_final = [] 380 dists_initial = [] 381 n_points = len(ids_point) 382 ind_point_initial = -1 383 ind_point_final = -1 384 t_point_initial = -1 385 t_point_final = -1 386 387 # print ' search initial edges' 388 ind_point = 0 389 for p in coords: 390 dists2 = get_dist_point_to_segs(p, x1, y1, x2, y2) 391 inds_hit = np.flatnonzero(dists2 < self.width_buffer_max**2) 392 # print ' id_point,inds_hit',ids_point[ind_point],inds_hit 393 #ids_edge = get_ids_edge_from_inds_seg(inds_hit) 394 395 if len(inds_hit) > 0: 396 if ind_point_initial < 0: 397 # print ' inds_hit',ind_point,inds_hit 398 ids_edge_initial = set(get_ids_edge_from_inds_seg(inds_hit)) 399 # print ' ids_edge_initial',ids_edge_initial 400 #dists_initial = np.sqrt(dists2[inds_hit]) 401 ind_point_initial = ind_point 402 t_point_initial = pointtimes[ind_point] 403 # print ' ind_point, inds_hit, ind_point_initial',ind_point,inds_hit,ind_point_initial,ids_edge_initial 404 else: 405 ids_edge_final = set(get_ids_edge_from_inds_seg(inds_hit)) 406 # print ' ids_edge_final',ids_edge_final 407 ind_point_final = ind_point 408 t_point_final = pointtimes[ind_point] 409 # print ' ind_point, inds_hit, ind_point_final',ind_point,inds_hit,ind_point_initial,ids_edge_final 410 411 ind_point += 1 412 413 n_points_eff = ind_point_final - ind_point_initial 414 415 duration_gps = pointtimes[ind_point_final] - pointtimes[ind_point_initial] 416 417 # if self._logger: 418 # 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)) ) 419 # print ' completed init:' 420 print ' ids_edge_initial', ids_edge_initial 421 print ' ids_edge_final', ids_edge_final 422 423 # print ' ind_point_initial,ind_point_final,n_points_eff',ind_point_initial,ind_point_final,n_points_eff 424 print ' id_point_initial=%d,id_point_final=%d,n_points_eff=%d' % ( 425 ids_point[ind_point_initial], ids_point[ind_point_final], n_points_eff) 426 if (ind_point_initial < 0) | (ind_point_final < 0) | (n_points_eff < self.n_points_min): 427 print 'ABOARD: insufficient valid points' 428 return [], 0.0, 0.0, -1.0, -1.0, -1.0, -1.0, 0.0 429 430 # print ' coords',coords 431 # create initial routeset 432 routelist = [] 433 434 #id_edge_unique = [] 435 for id_edge in ids_edge_initial: 436 437 cost0 = get_dist_point_to_edge(coords[ind_point_initial], id_edge) 438 439 # print ' coords of initial point:',coords[ind_point_initial] 440 # print ' cost0, id_edge',cost0, id_edge 441 442 #get_ind_seg_from_id_edge(self, id_edge) 443 444 length_edge = edges.lengths[id_edge] 445 accesslevel = accesslevels[id_edge] 446 cost_tot = cost0+self.weight_cumlength * length_edge - self.weight_access * accesslevel 447 routelist.append([cost_tot, cost0, 0.0, length_edge, accesslevel, length_edge, [id_edge], []]) 448 449 # print ' initial routelist',routelist 450 451 # --------------------------------------------------------------------- 452 # 2. main loop through rest of the points 453 454 # for i in xrange(ind_point+1, n_points):#coords[ind_point+1:]: 455 #is_connected = True 456 ind_point = ind_point_initial+1 457 length_gps = 0.0 458 while (ind_point < ind_point_final): # & is_connected: 459 point = coords[ind_point] 460 delta_point = point-coords[ind_point-1] 461 dist_interpoint = np.sqrt(np.sum(delta_point**2)) 462 phi_point = np.arctan2(delta_point[1], delta_point[0]) 463 length_gps += dist_interpoint 464 465 # print 79*'_' 466 # 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) 467 468 routelist_new = [] 469 ind_route = 0 470 for routeinfo in routelist: 471 costs_tot, cost, length_partial, length_cum, accesslevel, length_edge, ids_edge, ids_point_edgeend = routeinfo 472 473 if len(ids_edge) == 0: 474 break 475 476 length_partial += dist_interpoint 477 id_edge_last = ids_edge[-1] 478 # minimum distance point-edge 479 dist_point_edge = get_dist_point_to_edge(point, id_edge_last, 480 is_detect_initial=False, 481 is_detect_final=True,) 482 #dist_point_edge = self.get_dist_point_edge(coords[ind_point], id_edge_last, is_detect_final = True) 483 # if dist_point_edge is not a number (nan) the point is outside projection of edge 484 485 # print 79*'-' 486 # print ' route ',ind_route,'costs_tot',costs_tot,'dist_point_edge',dist_point_edge, len(ids_edge),id_edge_last 487 # print ' cost= %.3f, length_cum=%.3f'%(cost,length_cum,) 488 # print ' Xost= %.3f, Xength_cum=%.3f'%(cost,self.weight_cumlength *length_cum,) 489 # print ' ids_edge',ids_edge 490 # if length_partial > self.alpha * length_edge: 491 if (np.isnan(dist_point_edge)) | (dist_point_edge > self.width_buffer_max): 492 # point is beyond edge 493 494 # has reached end of edge 495 length_partial = length_partial-length_edge 496 497 # get FSTAR 498 # print ' fstar',id_edge_last,fstar[id_edge_last] 499 ids_edge_next_all = fstar[id_edge_last] 500 501 # filter allowed edges, by verifying positivness of 502 # travel costs 503 ids_edge_next = ids_edge_next_all[costs[ids_edge_next_all] >= 0] 504 # print ' parse ids_edge_next',ids_edge_next 505 506 if len(ids_edge_next) == 0: 507 # route not connected 508 cost_edge = np.Inf # punish route with infinite costs 509 length_edge = np.Inf 510 length_cum = np.Inf 511 accesslevel = -1 512 ids_edge_new = [] 513 ids_point_edgeend_new = [] 514 else: 515 516 id_edge = ids_edge_next[0] 517 accesslevel = accesslevels[id_edge] 518 length_edge = edges.lengths[id_edge] 519 length_cum += length_edge 520 cost_edge = self._get_cost_birgil(point, id_edge, phi_point, 521 accesslevel, get_dist_point_to_edge) 522 #self.get_dist_point_edge(coords[ind_point], ind_edge) 523 ids_edge_new = ids_edge + [id_edge] 524 ids_point_edgeend_new = ids_point_edgeend + [ids_point[ind_point]] 525 # print ' add to route id_edge',id_edge,ids_edge_new 526 527 # save updated data of current route 528 cost0 = cost+cost_edge 529 cost_tot = cost0+self.weight_cumlength * length_cum 530 routeinfo[:] = cost_tot, cost0, length_partial, length_cum, accesslevel, length_edge, ids_edge_new, ids_point_edgeend_new 531 532 # add new children if forward star greater 1 533 if len(ids_edge_next) > 1: 534 for id_edge in ids_edge_next[1:]: 535 #ind_edge = self._edge_to_id_edge[edge] 536 length_edge = edges.lengths[id_edge] 537 length_cum += length_edge 538 accesslevel = accesslevels[id_edge] 539 cost_edge = self._get_cost_birgil( 540 point, id_edge, phi_point, accesslevel, get_dist_point_to_edge) 541 ids_edge_new = ids_edge + [id_edge] 542 ids_point_edgeend_new = ids_point_edgeend + [ids_point[ind_point]] 543 # print ' new route id_edge',id_edge,ids_edge_new 544 cost0 = cost+cost_edge 545 cost_tot = cost0+self.weight_cumlength * length_cum 546 routelist_new.append([cost_tot, cost0, length_partial, length_cum, 547 accesslevel, length_edge, ids_edge_new, ids_point_edgeend_new]) 548 549 else: 550 # has not reached end of current edge 551 # save updated data of current route 552 if len(ids_edge) > 0: 553 #dist += self.get_dist_point_edge(coords[ind_point], ids_edge[-1]) 554 cost += self._get_cost_birgil(point, id_edge_last, phi_point, 555 accesslevel, get_dist_point_to_edge) 556 cost_tot = cost+self.weight_cumlength * length_cum 557 routeinfo[:] = cost_tot, cost, length_partial, length_cum, accesslevel, length_edge, ids_edge, ids_point_edgeend 558 559 ind_route += 1 560 561 routelist += routelist_new 562 # print ' routelist',routelist 563 routelist.sort() 564 565 # cut list to self.n_routes_follow 566 if len(routelist) > self.n_routes_follow+1: 567 routelist = routelist[:self.n_routes_follow+1] 568 569 ind_point += 1 570 571 # recover final edge objects of winner route! 572 costs_tot, cost, length_partial, length_cum, accesslevel, length_edge, ids_edge, ids_point_edgeend = routelist[ 573 0] 574 575 # print ' ids_edge[-1]',ids_edge[-1],ids_edge[-1] in ids_edge_final 576 # print ' ids_edge_final',ids_edge_final 577 t_match = time.time() - tick 578 579 # -------------------------------------------------------------------- 580 # post matching analisis 581 print '\n'+79*'-' 582 print ' matched route:', ids_edge 583 # print ' ids_point_edgeend',ids_point_edgeend 584 route = ids_edge 585 if len(route) == 0: 586 print 'ABOARD: route contains no edges ids_edge=', ids_edge 587 return [], 0.0, 0.0, -1.0, -1.0, -1.0, -1.0, 0.0, [], False 588 589 length_route = np.sum(costs[ids_edge]) 590 length_route_mixed = 0.0 591 length_route_exclusive = 0.0 592 for id_edge, length in zip(ids_edge, costs[ids_edge]): 593 accesslevel = accesslevels[id_edge] 594 if accesslevel == 1: 595 length_route_mixed += length 596 elif accesslevel == 2: 597 length_route_exclusive += length 598 599 if length_gps > 0: 600 lengthindex = length_route/length_gps 601 else: 602 lengthindex = -1.0 603 604 if ids_edge[-1] in ids_edge_final: # dist!=np.Inf: DID WE ARRIVE? 605 print 'SUCCESS: target', ids_edge[-1], 'reached' 606 is_connected = True 607 else: 608 print 'DISCONNECTED: last matched edge', ids_edge[-1], ' did not reach final edges', ids_edge_final 609 is_connected = False 610 611 # check dist error 612 inds_point = xrange(ind_point_initial, ind_point_final) 613 n_points = len(inds_point) 614 615 if (n_points >= 2) & (len(route) > 0): 616 # print ' initialize distance error measurement',len(route) 617 dist_points_tot = 0.0 618 routestring = get_routelinestring(route, edges) 619 # print ' routestring =',routestring 620 621 # print ' coords',coords[inds_point].tolist() 622 tracepoints = MultiPoint(coords[inds_point].tolist()) 623 624 # measure precision 625 #routeset = set(route) 626 #n_points_boundary = 0 627 ind_tracepoint = 0 628 for ind_point in inds_point: # xrange(n_points): 629 d = routestring.distance(tracepoints[ind_tracepoint]) 630 dist_points_tot += d 631 # print ' v distance measurement',d,tracepoints[ind_tracepoint],coords[ind_point]#,n_points#,n_points_eff 632 633 ind_tracepoint += 1 634 635 err_dist = dist_points_tot/float(n_points) 636 637 # print ' check',(lengthindex<self.lengthindex_max),(lengthindex >= self.lengthindex_min),(duration_gps>self.duration_min),(length_route>self.dist_min),(len(route)>0) 638 # print ' twice',(lengthindex<self.lengthindex_max)&(lengthindex >= self.lengthindex_min)&(duration_gps>self.duration_min)&(length_route>self.dist_min)&(len(route)>0) 639 640 return route, length_route, length_route_mixed, length_route_exclusive, duration_gps, lengthindex, err_dist, t_match, ids_point_edgeend, is_connected 641 642 def _get_cost_birgil(self, p, id_edge, phi_point, accesslevel, get_dist_point_to_edge): 643 # print '_get_cost_birgil p, id_edge',p, id_edge 644 645 dist, segment = get_dist_point_to_edge(p, id_edge, 646 is_return_segment=True) 647 x1, y1, x2, y2 = segment 648 649 # print ' x1,y1',x1,y1 650 # print ' p',p[0],p[1] 651 # print ' x2,y2',x2,y2 652 # print ' dist',dist 653 #seginds = self._id_edge_to_seginds[ind_edge] 654 #dists2 = get_dist_point_to_segs( p, self._x1[seginds], self._y1[seginds], self._x2[seginds], self._y2[seginds]) 655 656 #ind_segind = np.argmin(dists2) 657 #segind = seginds[ind_segind] 658 #x1,y1,x2,y2 = (self._x1[segind], self._y1[segind], self._x2[segind], self._y2[segind]) 659 phi_seg = np.arctan2(y2-y1, x2-x1) 660 #dist_ps = max(dist, self.dist_min) 661 # print ' x1,y1,x2,y2',x1,y1,x2,y2 662 # 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) 663 #cost = dist_ps*np.abs(np.sin(np.abs(phi_point-phi_seg))) 664 cost_angle = np.abs(np.sin(np.clip(phi_point-phi_seg, -np.pi/2, np.pi/2))) 665 cost = dist + self.weight_angle * cost_angle - self.weight_access * accesslevel 666 667 # print ' cost_birgil=%.3f, dist=%.3f, cost_angle=%.3f, accesslevel=%d'%(cost,dist,cost_angle,accesslevel) 668 # print ' cost_birgil=%.3f, Xist=%.3f, Xost_angle=%.3f, Xccesslevel=%.3f'%(cost,dist,self.weight_angle *cost_angle,self.weight_access*accesslevel) 669 # print ' cost_birgil=',cost 670 return cost 671 672 def get_dist_point_edge(self, p, ind_edge, is_detect_final=False): 673 seginds = self._id_edge_to_seginds[ind_edge] 674 dists2 = get_dist_point_to_segs(p, self._x1[seginds], self._y1[seginds], self._x2[seginds], 675 self._y2[seginds], is_ending=True, is_detect_final=is_detect_final) 676 if is_detect_final: 677 is_final = np.isnan(dists2) 678 if np.all(is_final): # point outside final segction of edge 679 return np.nan 680 else: 681 return np.sqrt(np.min(dists2[~is_final])) 682 683 return np.sqrt(np.min(dists2)) 684 685 def get_edge_info(self, id_trip, log=None): 686 # print 'get_edge_info for trace %s with c_length=%.6f'%(id_trip,self.c_length) 687 pointset = self.pointsets.get(id_trip) 688 689 # create a multi point object with points of traces 690 # should no longer be necessary 691 #tracepoints = MultiPoint(pointset.get_coords().tolist()) 692 693 # make a array with time stamps of all points 694 pointtimes = pointset.get_times() 695 ids_points = pointset.get_ids() 696 segind_to_id_edge = self._segind_to_id_edge 697 #edgehit_counts = np.zeros(len(self._net.getEdges()),int) 698 # segind_to_id_edge = self._segind_to_id_edge 699 # for p in pointset.get_coords(): 700 # print ' ',p 701 # hitinds = get_dist_point_to_segs(p,self._x1,self._y1,self._x2,self._y2, is_ending=True)< self.width_buffer**2 702 # edgehit_counts[segind_to_id_edge[hitinds]]+= 1 703 704 # create a list with lists for each point 705 # Each list contans the edges, where the point is in the edge buffer 706 edges_containing_point = np.zeros(len(ids_points), object) 707 edgedists_to_point = np.zeros(len(ids_points), object) 708 ind_point = 0 709 for p in pointset.get_coords(): 710 dists2 = get_dist_point_to_segs(p, self._x1, self._y1, self._x2, self._y2, is_ending=True) 711 inds_hit = np.nonzero(dists2 < self.width_buffer**2) 712 edges_containing_point[ind_point] = self._edges[segind_to_id_edge[inds_hit]] 713 edgedists_to_point[ind_point] = np.sqrt(dists2[inds_hit]) 714 715 # print ' ind_point,p,inds_hit',ind_point,p,len(np.nonzero(inds_hit)[0])#,edges_containing_point[ind_point] 716 ind_point += 1 717 718 # this is a dictionary with edge instance as key and 719 # the weight as value. 720 occurrences = {} 721 for edge in self._edges: 722 occurrences[edge] = 0.0 723 724 id_mode = self.id_mode 725 726 # 2. determine part of weight inversely proportional to the number 727 # of GPS points contained in the buffer 728 ind = 0 729 time_start = +10**10 730 time_end = -10**10 731 732 ind_point_start = -1 733 ind_point_end = -1 734 intervals = {} 735 for pointedges in edges_containing_point: 736 n_edges = len(pointedges) 737 if n_edges > 0: 738 t_point = pointtimes[ind] 739 if t_point < time_start: 740 time_start = t_point 741 ind_point_start = ind 742 743 elif t_point > time_end: 744 time_end = t_point 745 ind_point_end = ind 746 747 # compute partial weight for containing GPS points in edge buffer 748 wp = 1.0/n_edges 749 750 for edge in pointedges: 751 # assign weight component to edges 752 occurrences[edge] += wp 753 754 # determine intervals per edge 755 if intervals.has_key(edge): 756 t_edge_start, t_edge_end = intervals[edge] 757 if t_point < t_edge_start: 758 t_edge_start = t_point 759 elif t_point > t_edge_end: 760 t_edge_end = t_point 761 intervals[edge] = (t_edge_start, t_edge_end) 762 763 else: 764 intervals[edge] = (t_point, t_point) 765 766 ind += 1 767 768 edges_start = set(edges_containing_point[ind_point_start]) 769 edges_end = set(edges_containing_point[ind_point_end]) 770 771 # 1. calculate length proportional edge weight 772 c_length = float(self.c_length) 773 c_bike = float(self.c_bike) 774 weights = {} 775 776 for edge in self._edges: 777 len_edge = edge.getLength() 778 p_empty = c_length*len_edge 779 780 l = len_edge 781 allowed = edge.getLanes()[0].getAllowed() 782 if len(allowed) > 0: 783 if id_mode in allowed: 784 l *= c_bike # weight for dedicated lanes 785 786 if occurrences[edge] > 0.0: 787 # compute edge length proportional weight 788 weights[edge] = l/occurrences[edge] 789 790 else: 791 weights[edge] = l/p_empty 792 793 # if self._logger: 794 # self._logger.w( ' New: duration in network %d s from %ds to %ds'%(time_end-time_start,time_start,time_end)) 795 796 #log.w( ' Check dimensions len(weights)=%d,len(intervals)=%d'%(len(weights), len(intervals))) 797 # print ' intervals',intervals 798 #log.w( ' Possible start edges') 799 # for edge in edges_start: 800 # log.w( ' %s'%edge.getID()) 801 #log.w( ' Possible end edges') 802 # for edge in edges_end: 803 # log.w( ' %s'%edge.getID()) 804 # 805 #log.w('get_edge_info done.') 806 807 #self.intervals.set(id_trip, [ time_start,time_end]) 808 return weights, edges_start, edges_end, intervals, edges_containing_point, ind_point_start, ind_point_end 809 810 811def get_colvalue(val, default=0.0): 812 813 if (len(val) > 0) & (val != 'NULL'): 814 return float(val) 815 else: 816 return default 817 818 819class FilterMixin(Process): 820 def _init_filter_preview(self): 821 attrsman = self.get_attrsman() 822 trips = self.parent.trips 823 self.filterresults = attrsman.add(cm.FuncConf('filterresults', 824 'filterpreview', # function attribute of object 825 '%d/%d' % (len(trips.get_ids_selected()), len(trips)), 826 name='Preview selected trips', 827 groupnames=['options', 'results'], 828 )) 829 830 def filterpreview(self): 831 """ 832 Previews selected trips after filtering. 833 """ 834 trips = self.parent.trips 835 n_trips = len(trips) 836 n_sel_current = len(trips.get_ids_selected()) 837 n_sel_eliminate = len(self.filter_ids()) 838 n_sel_after = n_sel_current-n_sel_eliminate 839 return '%d/%d (currently %d/%d)' % (n_sel_after, n_trips, n_sel_current, n_trips) 840 841 def _init_filter_time(self, **kwargs): 842 attrsman = self.get_attrsman() 843 844 self.hour_from_morning = attrsman.add(cm.AttrConf('hour_from_morning', kwargs.get('hour_from_morning', 5), 845 groupnames=['options'], 846 perm='rw', 847 name='From morning hour', 848 unit='h', 849 info='Keep only morning trips which start after this hour.', 850 )) 851 852 self.hour_to_morning = attrsman.add(cm.AttrConf('hour_to_morning', kwargs.get('hour_to_morning', 9), 853 groupnames=['options'], 854 perm='rw', 855 name='To morning hour', 856 unit='h', 857 info='Keep only morning trips which start before this hour.', 858 )) 859 860 self.hour_from_evening = attrsman.add(cm.AttrConf('hour_from_evening', kwargs.get('hour_from_evening', 15), 861 groupnames=['options'], 862 perm='rw', 863 name='From evening hour', 864 unit='h', 865 info='Keep only evening trips which start after this hour.', 866 )) 867 868 self.hour_to_evening = attrsman.add(cm.AttrConf('hour_to_evening', kwargs.get('hour_to_evening', 19), 869 groupnames=['options'], 870 perm='rw', 871 name='To evening hour', 872 unit='h', 873 info='Keep only evening trips which start before this hour.', 874 )) 875 876 self.weekdays = attrsman.add(cm.AttrConf('weekdays', kwargs.get('weekdays', WEEKDAYCHOICES['all']), 877 groupnames=['options'], 878 name='Weekdays', 879 choices=WEEKDAYCHOICES, 880 info='Keep only trips at the given weekdays.', 881 )) 882 883 def filter_time(self, timestamps): 884 """ 885 timestamps is an array with timestamps. 886 Function returns a binary array, with a True value for each 887 time stamp that does NOT SATISFY the specified time constrants. 888 """ 889 # print 'filter_time' 890 localtime = time.localtime 891 inds_elim = np.zeros(len(timestamps), dtype=np.bool) 892 i = 0 893 for timestamp in timestamps: 894 895 dt = localtime(timestamp) 896 is_keep = dt.tm_wday in self.weekdays 897 h = dt.tm_hour 898 is_keep &= (h > self.hour_from_morning) & (h < self.hour_to_morning)\ 899 | (h > self.hour_from_evening) & (h < self.hour_to_evening) 900 901 # print ' is_keep,w,h=',is_keep,dt.tm_wday,h 902 inds_elim[i] = not is_keep 903 i += 1 904 905 return inds_elim 906 907 def is_timestamp_ok(self, timestamp): 908 dt = time.localtime(timestamp) 909 is_ok = dt.tm_wday in self.weekdays 910 h = dt.tm_hour 911 is_ok &= (h > self.hour_from_morning) & (h < self.hour_to_morning)\ 912 | (h > self.hour_from_evening) & (h < self.hour_to_evening) 913 914 return is_ok 915 916 def filter_ids(self): 917 """ 918 Returns an array of ids to be eliminated or deselected. 919 To be overridden 920 """ 921 return [] 922 923 924class PostMatchfilter(FilterMixin): 925 def __init__(self, mapmatching, logger=None, **kwargs): 926 print 'PostMatchfilter.__init__' 927 self._init_common('postmatchfilter', 928 parent=mapmatching, 929 name='Post matchfilter', 930 logger=logger, 931 info='Removes matched tripe with defined characteristics from current selection.', 932 ) 933 934 attrsman = self.set_attrsman(cm.Attrsman(self)) 935 936 info_lengthindex = "Length index is length of matched route divided by length of line interpolated GPS points in percent." 937 self.lengthindex_min = attrsman.add(cm.AttrConf('lengthindex_min', kwargs.get('lengthindex_min', 80.0), 938 groupnames=['options'], 939 name='Min. length index', 940 unit='%', 941 info='Minimum allowed length index. '+info_lengthindex, 942 )) 943 self.lengthindex_max = attrsman.add(cm.AttrConf('lengthindex_max', kwargs.get('lengthindex_max', 110.0), 944 groupnames=['options'], 945 name='Min. length index', 946 unit='%', 947 info='Maximum allowed length index '+info_lengthindex, 948 )) 949 950 info_error_dist = 'The distance error is the average distance between the GPS points and the matched route.' 951 self.error_dist_max = attrsman.add(cm.AttrConf('error_dist_max', kwargs.get('error_dist_max', 10000.0), 952 groupnames=['options'], 953 name='Max. distance err.', 954 unit='mm', 955 info='Maximum allowed distance error. '+info_error_dist, 956 )) 957 958 self._init_filter_time() 959 self._init_filter_preview() 960 961 def filter_ids(self): 962 """ 963 Returns an array of ids to be eliminated or deselected. 964 """ 965 # print 'filter_ids' 966 trips = self.parent.trips 967 ids_selected = trips.get_ids_selected() 968 inds_eliminate = np.logical_or( 969 trips.lengthindexes[ids_selected] < self.lengthindex_min, 970 trips.lengthindexes[ids_selected] > self.lengthindex_max, 971 trips.errors_dist[ids_selected] > self.error_dist_max, 972 ) 973 # print ' inds_eliminate',inds_eliminate 974 inds_eliminate = np.logical_or(inds_eliminate, self.filter_time(trips.timestamps[ids_selected])) 975 # print ' inds_eliminate',inds_eliminate 976 return ids_selected[inds_eliminate] 977 978 def do(self): 979 980 # execute filtering 981 self.parent.trips.are_selected[self.filter_ids()] = False 982 983 984class TripGeomfilter(FilterMixin): 985 def __init__(self, mapmatching, logger=None, **kwargs): 986 print 'TripGeomfilter.__init__' 987 self._init_common('tripGeomfilter', 988 parent=mapmatching, 989 name='Geometry trip filter', 990 logger=logger, 991 info='Removes trips with defined geometric characteristics from current selection.', 992 ) 993 994 attrsman = self.set_attrsman(cm.Attrsman(self)) 995 996 self.dist_point_max = attrsman.add(cm.AttrConf('dist_point_max', kwargs.get('dist_point_max', 1000.0), 997 groupnames=['options'], 998 perm='rw', 999 name='Max. point dist.', 1000 unit='m', 1001 info='Keep only traces where the distance between all successive points is below this maximal distance.', 1002 )) 1003 1004 self.duration_point_max = attrsman.add(cm.AttrConf('duration_point_max', kwargs.get('duration_point_max', 300.0), 1005 groupnames=['options'], 1006 perm='rw', 1007 name='Max. point duration.', 1008 unit='s', 1009 info='Keep only traces where the duration between all successive points is below this maximal duration.', 1010 )) 1011 1012 self.const_return_max = attrsman.add(cm.AttrConf('const_return_max', kwargs.get('const_return_max', 0.3), 1013 groupnames=['options'], 1014 perm='rw', 1015 name='Max. return const', 1016 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.', 1017 )) 1018 1019 self._init_filter_preview() 1020 1021 def do(self): 1022 1023 # execute filtering 1024 self.parent.trips.are_selected[self.filter_ids()] = False 1025 1026 def filter_ids(self): 1027 """ 1028 Returns an array of ids to be eliminated or deselected. 1029 """ 1030 c_cutoff = 1.0 - self.const_return_max 1031 print 'TripGeomfilter.do c_cutoff', c_cutoff 1032 dist_point_max = self.dist_point_max 1033 duration_point_max = self.duration_point_max 1034 trips = self.parent.trips 1035 points = self.parent.points 1036 1037 ids_trips = trips.get_ids_selected() 1038 1039 ids_points = trips.ids_points[ids_trips] 1040 # print ' ids_trips',ids_trips 1041 # print ' ids_points',ids_points 1042 1043 intersects_boundaries = self.parent.get_scenario().net.intersects_boundaries 1044 1045 inds_elim = np.zeros(len(ids_trips), dtype=np.bool) 1046 j = 0 1047 for id_trip, ids_point in zip(ids_trips, ids_points): 1048 1049 coords = points.coords[ids_point] 1050 times = points.timestamps[ids_point] 1051 # print 79*'-' 1052 # print ' check split ',id_trip 1053 # print ' ids_point', ids_point 1054 # print ' times',times 1055 # print ' coords', coords 1056 # print ' duration_point_max',duration_point_max 1057 #dist = pointset.get_distance() 1058 1059 if ids_point is None: 1060 # this happens if the points of a trip in the workout file 1061 # have not been imported for some reason 1062 is_eliminate = True 1063 1064 elif not intersects_boundaries(get_boundary(coords)): 1065 is_eliminate = True 1066 1067 elif np.any(times < 0): 1068 is_eliminate = True 1069 else: 1070 dist_max = 0.0 1071 is_eliminate = False 1072 i = 0 1073 n = len(times) 1074 x, y, z = coords[0] 1075 while (not is_eliminate) & (i < (n-1)): 1076 i += 1 1077 dist_point = np.sqrt((coords[i, 0]-coords[i-1, 0])**2 1078 + (coords[i, 1]-coords[i-1, 1])**2) 1079 1080 if dist_point > dist_point_max: 1081 is_eliminate = True 1082 1083 # print ' times[i]-times[i-1]',times[i]-times[i-1] 1084 1085 if (times[i]-times[i-1]) > duration_point_max: 1086 is_eliminate = True 1087 1088 d = np.sqrt((x-coords[i, 0])**2 + (y-coords[i, 1])**2) 1089 if d > dist_max: 1090 dist_max = d 1091 1092 elif d < c_cutoff*dist_max: 1093 is_eliminate = True 1094 1095 inds_elim[j] = is_eliminate 1096 j += 1 1097 1098 return ids_trips[inds_elim] 1099 1100 1101class TripTimefilter(Process): 1102 def __init__(self, mapmatching, logger=None, **kwargs): 1103 print 'Timefilter.__init__' 1104 self._init_common('triptimefilter', 1105 parent=mapmatching, 1106 name='Trip time filter', 1107 logger=logger, 1108 info='Keeps only trips with specified time charactaristict in selection.', 1109 ) 1110 1111 def do(self): 1112 c_cutoff = 1.0 - self.const_return_max 1113 print 'TripGeomfilter.do c_cutoff', c_cutoff 1114 dist_point_max = self.dist_point_max 1115 duration_point_max = self.duration_point_max 1116 trips = self.parent.trips 1117 #points = self.parent.points 1118 ids_trips = trips.get_ids_selected() 1119 1120 #ids_points = trips.ids_points[ids_trips] 1121 # print ' ids_trips',ids_trips 1122 # print ' ids_points',ids_points 1123 1124 for id_trip, ids_point in zip(ids_trips, ids_points): 1125 1126 coords = points.coords[ids_point] 1127 times = points.timestamps[ids_point] 1128 print 79*'-' 1129 print ' check split ', id_trip 1130 print ' ids_point', ids_point 1131 print ' times', times 1132 print ' duration_point_max', duration_point_max 1133 #dist = pointset.get_distance() 1134 1135 if np.any(times < 0): 1136 is_eliminate = True 1137 else: 1138 dist_max = 0.0 1139 is_eliminate = False 1140 i = 0 1141 n = len(times) 1142 x, y, z = coords[0] 1143 while (not is_eliminate) & (i < (n-1)): 1144 i += 1 1145 dist_point = np.sqrt((coords[i, 0]-coords[i-1, 0])**2 1146 + (coords[i, 1]-coords[i-1, 1])**2) 1147 1148 if dist_point > dist_point_max: 1149 is_eliminate = True 1150 1151 print ' times[i]-times[i-1]', times[i]-times[i-1] 1152 1153 if (times[i]-times[i-1]) > duration_point_max: 1154 is_eliminate = True 1155 1156 d = np.sqrt((x-coords[i, 0])**2 + (y-coords[i, 1])**2) 1157 if d > dist_max: 1158 dist_max = d 1159 1160 elif d < c_cutoff*dist_max: 1161 is_eliminate = True 1162 1163 trips.are_selected[id_trip] = not is_eliminate 1164 1165 1166class EccTracesImporter(FilterMixin): 1167 def __init__(self, mapmatching, logger=None, **kwargs): 1168 print 'EccTracesImporter.__init__', mapmatching.get_ident() 1169 self._init_common('traceimporter', 1170 parent=mapmatching, 1171 name='ECC Trace Importer', 1172 logger=logger, 1173 info='Import workouts and GPS points of a European cycling challange.', 1174 ) 1175 1176 attrsman = self.set_attrsman(cm.Attrsman(self)) 1177 1178 scenario = mapmatching.get_scenario() 1179 rootfilepath = scenario.get_rootfilepath() 1180 1181 # here we ged classes not vehicle type 1182 # specific vehicle type within a class will be generated later 1183 modechoices = scenario.net.modes.names.get_indexmap() 1184 1185 # print ' modechoices',modechoices 1186 self.id_mode = attrsman.add(am.AttrConf('id_mode', modechoices['bicycle'], 1187 groupnames=['options'], 1188 choices=modechoices, 1189 name='Mode', 1190 info='Transport mode to be matched.', 1191 )) 1192 1193 self.workoutsfilepath = attrsman.add( 1194 cm.AttrConf('workoutsfilepath', kwargs.get('workoutsfilepath', rootfilepath+'.workouts.csv'), 1195 groupnames=['options'], 1196 perm='rw', 1197 name='Workout file', 1198 wildcards='CSV file (*.csv)|*.csv', 1199 metatype='filepath', 1200 info="""CSV text file with workout database.""", 1201 )) 1202 1203 self.pointsfilepath = attrsman.add(cm.AttrConf('pointsfilepath', kwargs.get('pointsfilepath', rootfilepath+'.points.csv'), 1204 groupnames=['options'], 1205 perm='rw', 1206 name='Points file', 1207 wildcards='CSV file (*.csv)|*.csv', 1208 metatype='filepath', 1209 info="CSV text file with GPS point database.", 1210 )) 1211 1212 self.year = attrsman.add(cm.AttrConf('year', kwargs.get('year', 2014), 1213 groupnames=['options'], 1214 choices={'2014': 2014, 'from 2015': 2015}, 1215 perm='rw', 1216 name='Year of challange', 1217 info='Year of challange is used to identify the correct database format.', 1218 )) 1219 1220 self.dist_trip_min = attrsman.add(cm.AttrConf('dist_trip_min', kwargs.get('dist_trip_min', 100.0), 1221 groupnames=['options'], 1222 perm='rw', 1223 name='Min. trip distance', 1224 unit='m', 1225 info='Minimum distance of one trip. Shorter trips will not be imported.', 1226 )) 1227 1228 self.dist_trip_max = attrsman.add(cm.AttrConf('dist_trip_max', kwargs.get('dist_trip_max', 25000.0), 1229 groupnames=['options'], 1230 perm='rw', 1231 name='Max. trip distance', 1232 unit='m', 1233 info='Maximum distance of one trip. Shorter trips will not be imported.', 1234 )) 1235 1236 self.duration_trip_min = attrsman.add(cm.AttrConf('duration_trip_min', kwargs.get('duration_trip_min', 30.0), 1237 groupnames=['options'], 1238 perm='rw', 1239 name='Min. trip duration', 1240 unit='s', 1241 info='Minimum duration of one trip. Trips with shorter duration will not be imported.', 1242 )) 1243 1244 self.speed_trip_min = attrsman.add(cm.AttrConf('speed_trip_min', kwargs.get('speed_trip_min', 1.0), 1245 groupnames=['options'], 1246 perm='rw', 1247 name='Min. av. trip speed', 1248 unit='m/s', 1249 info='Minimum average trip speed. Trips with lower average speed will not be imported.', 1250 )) 1251 1252 self.speed_trip_max = attrsman.add(cm.AttrConf('speed_trip_max', kwargs.get('speed_trip_max', 50.0), 1253 groupnames=['options'], 1254 perm='rw', 1255 name='Max. av. trip speed', 1256 unit='m/s', 1257 info='Maximum average trip speed. Trips with higher average speed will not be imported.', 1258 )) 1259 1260 self.sep_column_workout = attrsman.add(cm.AttrConf('sep_column_workout', kwargs.get('sep_column_workout', ','), 1261 groupnames=['options'], 1262 perm='rw', 1263 name='Workoutdata seperator', 1264 info='Workout column seperator of CSV file', 1265 )) 1266 1267 self.sep_column_points = attrsman.add(cm.AttrConf('sep_column_points', kwargs.get('sep_column_points', ','), 1268 groupnames=['options'], 1269 perm='rw', 1270 name='Point data seperator', 1271 info='Pointdata column seperator of CSV file', 1272 )) 1273 1274 self._init_filter_time(**kwargs) 1275 1276 def do(self): 1277 print 'TraceImporter.do' 1278 if self.year == 2014: 1279 self.import_workouts_2014() 1280 self.import_points_2014() 1281 1282 self.parent.points.project() 1283 1284 def import_workouts_2014(self): 1285 print 'import_workouts_2014' 1286 # 2014 ecomondo workouts 1287 # id pointDBNode pointPathId startTime distance duration sport calories maxSpeed altitudeMin altitudeMax metersAscent metersDescent 1288 # 329308466 7 37073516 2014-05-01 19:00:00 26 15600 1 1182.64 NULL NULL NULL NULL NULL 1289 # 0 1 2 3 4 5 6 7 8 9 10 11 12 1290 1291 j_id, j_node, j_id_trip, j_time, j_dist, j_duration = range(6) 1292 n_cols = 13 1293 null = 'NULL' 1294 trips = self.parent.trips 1295 # persons = self.parent.persons # no person data in 2014 :( 1296 ids_person_sumo = {} 1297 1298 ids_trips = [] 1299 scenario = self.parent.get_scenario() 1300 1301 get_vtype_for_mode = scenario.demand.vtypes.get_vtype_for_mode 1302 1303 f = open(self.workoutsfilepath, 'r') 1304 #if self._logger: self._logger.w('import_workouts_2014 %s'%os.path.basename(self.workoutsfilepath)) 1305 i_line = 0 1306 sep = self.sep_column_workout 1307 1308 #self.get_logger().w(100.0*self.simtime/self.duration, key ='progress') 1309 for line in f.readlines()[1:]: 1310 # if True:#i_line>1: 1311 cols = line.split(sep) 1312 # print ' len(cols)',len(cols),n_cols 1313 if len(cols) == n_cols: 1314 # row is complete 1315 1316 # if cols[j_dist] != null: 1317 dist = get_colvalue(cols[j_dist])*1000.0 1318 duration = get_colvalue(cols[j_duration]) 1319 if duration > 0: 1320 speed_av = dist/duration 1321 else: 1322 speed_av = 0.0 1323 1324 if (dist > self.dist_trip_min)\ 1325 & (dist < self.dist_trip_max)\ 1326 & (duration > self.duration_trip_min)\ 1327 & (speed_av > self.speed_trip_min)\ 1328 & (speed_av < self.speed_trip_max): 1329 # print 'parametric conditions verified' 1330 timestamp = calc_seconds(cols[j_time]) 1331 if timestamp is not None: 1332 # print ' valid time stamp',timestamp 1333 if self.is_timestamp_ok(timestamp): 1334 id_trip = trips.make(id_sumo=cols[j_id_trip], 1335 id_vtype=get_vtype_for_mode(self.id_mode), 1336 timestamp=timestamp, 1337 distance_gps=dist, 1338 duration_gps=duration, 1339 speed_average=speed_av, 1340 ) 1341 ids_trips.append(id_trip) 1342 1343 def import_points_2014(self): 1344 print 'import_points_2014' 1345 # csv2014 1346 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1347 # 4,61565791,23648171762,2013-05-01 06:33:58,44.501085,11.372906,NULL,0,NULL,2,NULL 1348 1349 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1350 # 7,37073516,15138524460,NULL,44.51579,11.36257,NULL,NULL,NULL,NULL,NULL 1351 1352 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1353 # 7,37073516,15138524460,NULL,44.51579,11.36257,NULL,NULL,NULL,NULL,NULL 1354 1355 # Endomondo export format 1356 # pointDBNode,pointPathId,id,timestamp,latitude,longitude,altitude,distance,heartRate,instruction,speed 1357 # 0 1 2 3 4 5 6 7 8 9 10 1358 ind_id_path = 1 1359 ind_id_point = 2 1360 ind_time = 3 1361 ind_lat = 4 1362 ind_lon = 5 1363 ind_alt = 6 1364 ind_dist = 7 1365 ind_speed = 10 1366 1367 n_cols = 11 1368 #TripID, TimeStamp,Date, Latitude, Longitude, Altitude, Distance, Speed, Type 1369 1370 trips = self.parent.trips 1371 points = self.parent.points 1372 1373 exist_id_trip_sumo = trips.ids_sumo.has_index 1374 get_id_trip = trips.ids_sumo.get_id_from_index 1375 1376 sep = self.sep_column_points 1377 1378 f = open(self.pointsfilepath, 'r') 1379 if self._logger: 1380 self._logger.w('import_points_2014 %s' % os.path.basename(self.pointsfilepath)) 1381 i_line = 0 1382 id_trip_sumo = None 1383 id_trip = -1 1384 is_valid_trip = False 1385 1386 n_points_imported = 0 1387 for line in f.readlines(): 1388 1389 cols = line.split(sep) 1390 # print ' len(cols)',len(cols),n_cols 1391 if len(cols) == n_cols: 1392 id_trip_sumo_current = cols[ind_id_path] 1393 # print ' id_trip_sumo_current,id_trip_sumo',id_trip_sumo_current,id_trip_sumo,is_valid_trip 1394 if id_trip_sumo_current != id_trip_sumo: 1395 # this point is part of new trip 1396 1397 if is_valid_trip: 1398 # print ' store past points for valid trip',id_trip 1399 ids_point = points.add_rows( 1400 timestamps=timestamps, 1401 ids_trip=id_trip*np.ones(len(timestamps), dtype=np.int32), 1402 longitudes=longitudes, 1403 latitudes=latitudes, 1404 altitudes=altitudes, 1405 ) 1406 1407 trips.set_points(id_trip, ids_point) 1408 1409 # check if new trip is valid 1410 if exist_id_trip_sumo(id_trip_sumo_current): 1411 1412 is_valid_trip = True # start recording 1413 id_trip = get_id_trip(id_trip_sumo_current) 1414 # print ' found trip',id_trip,id_trip_sumo_current,' exisits-> record' 1415 1416 # ids_point_sumo = [] # useless? 1417 timestamps = [] 1418 ids_trip = [] 1419 longitudes = [] 1420 latitudes = [] 1421 altitudes = [] 1422 1423 else: 1424 # print ' trip',id_trip_sumo_current,'does not exisit' 1425 is_valid_trip = False 1426 1427 id_trip_sumo = id_trip_sumo_current 1428 1429 if is_valid_trip: 1430 # print ' store point timestamp',cols[ind_time] 1431 # current point belongs to a valid trip 1432 # ids_point_sumo.append(cols[ind_id_point]) 1433 timestamps.append(calc_seconds(cols[ind_time])) 1434 # ids_trip.append(id_trip) 1435 longitudes.append(get_colvalue(cols[ind_lon])) 1436 latitudes.append(get_colvalue(cols[ind_lat])) 1437 altitudes.append(get_colvalue(cols[ind_alt])) 1438 1439 else: 1440 print 'WARNING: inconsistent columns in line %d, file %s' % (i_line, filepath) 1441 1442 i_line += 1 1443 1444 # register points of last trip 1445 if is_valid_trip: 1446 # store past points for valid trip 1447 ids_point = points.add_rows( 1448 timestamps=timestamps, 1449 ids_trip=id_trip*np.ones(len(timestamps), dtype=np.int32), 1450 longitudes=longitudes, 1451 latitudes=latitudes, 1452 altitudes=altitudes, 1453 ) 1454 1455 trips.set_points(id_trip, ids_point) 1456 1457 # self.odtab.print_rows() 1458 1459 f.close() 1460 1461 1462class GpsTrips(Trips): 1463 def __init__(self, ident, mapmatching, **kwargs): 1464 # print 'Trips.__init__' 1465 self._init_objman(ident=ident, 1466 parent=mapmatching, 1467 name='Trips', 1468 info='Table with GPS trips, matched routes and alternative routes.', 1469 xmltag=('trips', 'trip', 'ids_sumo'), 1470 version=0.0, 1471 **kwargs) 1472 1473 self._init_attributes() 1474 self._init_constants() 1475 1476 def _init_attributes(self): 1477 1478 self.add_col(SumoIdsConf('GPS Trip', xmltag='id', info='GPS trip data.')) 1479 1480 self.add_col(am.ArrayConf('are_selected', default=True, 1481 dtype=np.bool, 1482 groupnames=['parameters', ], 1483 name='selected', 1484 symbol='Sel.', 1485 info='Selected for being processed (example mapmatching, export, etc).', 1486 )) 1487 1488 self.add_col(am.ArrayConf('timestamps', default=0, 1489 dtype=np.int, 1490 perm='r', 1491 groupnames=['parameters', 'gps'], 1492 name='timestamp', 1493 unit='s', 1494 metatype='datetime', 1495 info='Timestamp when trip started in seconds after 01 January 1970.', 1496 )) 1497 1498 self.timestamps.metatype = 'datetime' 1499 self.timestamps.set_perm('r') 1500 1501 self.add_col(am.ArrayConf('durations_gps', default=0.0, 1502 dtype=np.float32, 1503 groupnames=['parameters', 'gps'], 1504 name='GPS duration', 1505 unit='s', 1506 info='Time duration measure with GPS points.', 1507 )) 1508 1509 self.add_col(am.ArrayConf('distances_gps', default=0.0, 1510 dtype=np.float32, 1511 groupnames=['parameters', 'gps'], 1512 name='GPS distance', 1513 unit='m', 1514 info='Distance measure with GPS points.', 1515 )) 1516 1517 self.add_col(am.ArrayConf('speeds_average', default=0.0, 1518 dtype=np.float32, 1519 groupnames=['parameters', 'gps'], 1520 name='Av. speed', 1521 unit='m/s', 1522 info='Average speed based on GPS info.', 1523 )) 1524 self.add_col(am.ArrayConf('speeds_max', default=0.0, 1525 dtype=np.float32, 1526 groupnames=['parameters', 'gps'], 1527 name='Max. speed', 1528 unit='m/s', 1529 info='Maximum speed based on GPS info.', 1530 )) 1531 1532 # Trips._init_attributes(self) 1533 self.add_col(am.IdsArrayConf('ids_vtype', self.get_obj_vtypes(), 1534 groupnames=['state'], 1535 name='Type', 1536 info='Vehicle type.', 1537 xmltag='type', 1538 )) 1539 1540 self.add_col(am.ArrayConf('ids_purpose', default=TRIPPUROPSES['unknown'], 1541 dtype=np.int32, 1542 groupnames=['parameters', 'gps'], 1543 choices=TRIPPUROPSES, 1544 name='Purpose', 1545 info='Trip purpose ID', 1546 )) 1547 1548 self.add_col(am.ArrayConf('ids_device', default=DEVICES['unknown'], 1549 dtype=np.int32, 1550 groupnames=['parameters', 'gps'], 1551 choices=DEVICES, 1552 name='Devices', 1553 info='Device ID', 1554 )) 1555 self.add(cm.ObjConf(Routes('routes', self, self.parent.get_scenario().net))) 1556 1557 self.add_col(am.IdsArrayConf('ids_route_matched', self.get_routes(), 1558 groupnames=['results'], 1559 name='ID matched route', 1560 info='Route ID of mached route.', 1561 )) 1562 1563 self.add_col(am.IdsArrayConf('ids_route_shortest', self.get_routes(), 1564 groupnames=['results'], 1565 name='ID shortest route', 1566 info='Route ID of shortest route.', 1567 )) 1568 1569 self.add_col(am.ArrayConf('lengths_route_matched', default=-1.0, 1570 dtype=np.float32, 1571 groupnames=['results'], 1572 name='Matched length', 1573 symbol='L_{M}', 1574 unit='m', 1575 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.', 1576 )) 1577 1578 self.add_col(am.ArrayConf('lengths_route_shortest', default=-1.0, 1579 dtype=np.float32, 1580 groupnames=['results'], 1581 name='Shortest length', 1582 symbol='L_{S}', 1583 unit='m', 1584 info='Length of the shortest route. Shortest route is connecting the first matched edge and the final matched edge.', 1585 )) 1586 1587 self.add_col(am.ArrayConf('lengths_route_matched_mixed', default=-1.0, 1588 dtype=np.float32, 1589 groupnames=['results'], 1590 name='Matched length mixed access', 1591 symbol='L_{MM}', 1592 unit='m', 1593 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.', 1594 )) 1595 1596 self.add_col(am.ArrayConf('lengths_route_matched_exclusive', default=-1.0, 1597 dtype=np.float32, 1598 groupnames=['results'], 1599 name='Matched length exclusive access', 1600 symbol='L_{ME}', 1601 unit='m', 1602 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.', 1603 )) 1604 1605 self.add_col(am.ArrayConf('durations_route_matched', default=-1.0, 1606 dtype=np.float32, 1607 groupnames=['results'], 1608 name='Matched duration', 1609 unit='s', 1610 info='Duration of the matched part of the GPS trace. Note the only a fraction of the GPS trace ma be within the given network.', 1611 )) 1612 1613 self.add_col(am.ArrayConf('lengthindexes', default=-1.0, 1614 dtype=np.float32, 1615 groupnames=['results'], 1616 name='Length index', 1617 unit='%', 1618 info='Length index is the length of the matched route divided by length of line-interpolated GPS points.', 1619 )) 1620 1621 self.add_col(am.ArrayConf('errors_dist', default=-1.0, 1622 dtype=np.float32, 1623 groupnames=['results'], 1624 name='Distance error', 1625 unit='mm', 1626 info='The distance error is the average distance between the GPS points and the matched route.', 1627 )) 1628 1629 self.add_col(am.ArrayConf('times_computation', default=-1.0, 1630 dtype=np.float32, 1631 groupnames=['results'], 1632 name='Computation time', 1633 unit='ms', 1634 info='Computation time of the match algorithm.', 1635 )) 1636 1637 self.add_col(am.ArrayConf('are_match_connected', default=False, 1638 dtype=np.bool, 1639 groupnames=['results', ], 1640 name='Match connected', 1641 #symbol = 'Match conn.', 1642 info='The matched route connects first and last network edge where GPS points have been detected.', 1643 )) 1644 1645 self.add_col(am.IdlistsArrayConf('ids_points', self.parent.points, 1646 #groupnames = ['_private'], 1647 name='Point IDs', 1648 info="GPS point IDs.", 1649 )) 1650 1651 self.add_col(am.IdlistsArrayConf('ids_points_edgeend', self.parent.points, 1652 groupnames=['results', ], 1653 name='Edge endpoint IDs', 1654 info="This is a list of GPS point IDs which represent the last point associated with each matched edge.", 1655 )) 1656 1657 def get_obj_vtypes(self): 1658 return self.parent.get_scenario().demand.vtypes 1659 1660 def get_routes(self): 1661 return self.routes.get_value() 1662 1663 def clear_routes(self): 1664 self.get_routes().clear() 1665 for attrconf in self.get_group('results'): 1666 attrconf.reset() 1667 1668 def select_all(self): 1669 self.are_selected.get_value()[:] = True 1670 1671 def unselect_all(self): 1672 self.are_selected.get_value()[:] = False 1673 1674 def invert_selection(self): 1675 self.are_selected.get_value()[:] = np.logical_not(self.are_selected.get_value()[:]) 1676 1677 def delete_unselected(self): 1678 ids_del = self.select_ids(np.logical_not(self.are_selected.get_value())) 1679 points = self.parent.points 1680 for ids_point in self.ids_points[ids_del]: 1681 if ids_point is not None: 1682 points.del_rows(ids_point) 1683 self.del_rows(ids_del) 1684 1685 def set_matched_route(self, id_trip, route, 1686 length_matched=0.0, 1687 length_route_mixed=0.0, 1688 length_route_exclusive=0.0, 1689 duration_matched=0.0, 1690 lengthindex=-1.0, 1691 error_dist=-1.0, 1692 comptime=0.0, 1693 is_connected=False, 1694 ids_point_edgeend=[], 1695 ): 1696 1697 if len(route) > 0: 1698 id_route = self.ids_route_matched[id_trip] 1699 if id_route >= 0: 1700 # already a matched route existant 1701 self.get_routes().ids_edges[id_route] = route 1702 self.get_routes().colors[id_route] = COLOR_MATCHED_ROUTE 1703 # self.get_routes().set_row( id_route, 1704 # ids_edges = route, 1705 # colors = COLOR_MATCHED_ROUTE, 1706 # ) 1707 1708 # self.get_routes().set_row( id_route, 1709 # ids_trip = id_trip, 1710 # ids_edges = route, 1711 # #costs = duration_matched, 1712 # #probabilities = 1.0, 1713 # ) 1714 else: 1715 id_route = self.get_routes().add_row(ids_trip=id_trip, 1716 ids_edges=route, 1717 #costs = duration_matched, 1718 #probabilities = 1.0, 1719 colors=COLOR_MATCHED_ROUTE, 1720 ) 1721 else: 1722 id_route = -1 1723 1724 # print 'set_matched_route id_trip', id_trip,'id_route', id_route 1725 # print ' ids_point_edgeend',ids_point_edgeend 1726 #self.ids_route_matched[id_trip] = id_route 1727 1728 self.set_row(id_trip, 1729 ids_route_matched=id_route, 1730 lengths_route_matched=length_matched, 1731 lengths_route_matched_mixed=length_route_mixed, 1732 lengths_route_matched_exclusive=length_route_exclusive, 1733 durations_route_matched=duration_matched, 1734 lengthindexes=100*lengthindex, 1735 errors_dist=1000 * error_dist, 1736 times_computation=1000*comptime, 1737 are_match_connected=is_connected 1738 ) 1739 # TODO: do this extra! this is a bug!! 1740 # if included in set_row, only first value in list is taken!!! 1741 self.ids_points_edgeend[id_trip] = ids_point_edgeend 1742 1743 # print ' ids_points_edgeend[id_trip]',self.ids_points_edgeend[id_trip] 1744 1745 return id_route 1746 1747 def make(self, **kwargs): 1748 # if self.ids_sumo.has_index(id_sumo): 1749 # id_trip = self.ids_sumo.get_id_from_index(id_sumo) 1750 # #self.set_row(id_sumo, **kwargs) 1751 # return id_trip 1752 # else: 1753 id_trip = self.add_row(ids_sumo=kwargs.get('id_sumo', None), 1754 ids_vtype=kwargs.get('id_vtype', None), 1755 timestamps=kwargs.get('timestamp', None), 1756 distances_gps=kwargs.get('distance_gps', None), 1757 durations_gps=kwargs.get('duration_gps', None), 1758 speeds_average=kwargs.get('speed_average', None), 1759 speeds_max=kwargs.get('speed_max', None), 1760 ids_purpose=kwargs.get('id_purpose', None), 1761 ids_device=kwargs.get('id_device', None), 1762 ids_points=kwargs.get('ids_point', None), 1763 ) 1764 return id_trip 1765 1766 def set_points(self, id_trip, ids_point): 1767 self.ids_points[id_trip] = ids_point 1768 1769 def get_ids_selected(self): 1770 return self.select_ids(self.are_selected.get_value()) 1771 1772 def route_shortest(self): 1773 """ 1774 Shortest path routint. 1775 """ 1776 print 'route_shortest' 1777 # TODO: if too mant vtypes, better go through id_modes 1778 exectime_start = time.clock() 1779 scenario = self.parent.get_scenario() 1780 net = scenario.net 1781 edges = net.edges 1782 vtypes = scenario.demand.vtypes 1783 routes = self.get_routes() 1784 ids_edges = [] 1785 ids_trip = [] 1786 costs = [] 1787 1788 distancesmap = self.parent.get_distancesmap() 1789 1790 fstar = edges.get_fstar() 1791 for id_vtype in self.get_vtypes(): 1792 id_mode = vtypes.ids_mode[id_vtype] 1793 1794 # no routing for pedestrians 1795 if id_mode != net.modes.get_id_mode('pedestrian'): 1796 weights = distancesmap[id_mode] 1797 # weights = edges.get_distances( id_mode = id_mode, 1798 # is_check_lanes = is_check_lanes) 1799 1800 ids_trip_vtype = self.select_ids(np.logical_and(self.ids_vtype.get_value( 1801 ) == id_vtype, self.are_selected.get_value(), self.ids_route_matched.get_value() >= 0)) 1802 #ids_trip_vtype = self.get_trips_for_vtype(id_vtype) 1803 # print ' id_vtype,id_mode',id_vtype,id_mode#,ids_trip_vtype 1804 # print ' weights',weights 1805 1806 #ids_edge_depart = self.ids_edge_depart[ids_trip_vtype] 1807 #ids_edge_arrival = self.ids_edge_arrival[ids_trip_vtype] 1808 1809 for id_trip, id_route in zip(ids_trip_vtype, self.ids_route_matched[ids_trip_vtype]): 1810 route_matched = routes.ids_edges[id_route] 1811 1812 cost, route = routing.get_mincostroute_edge2edge(route_matched[0], 1813 route_matched[-1], 1814 weights=weights, 1815 fstar=fstar) 1816 if len(route) > 0: 1817 ids_edges.append(route) 1818 ids_trip.append(id_trip) 1819 costs.append(cost) 1820 1821 ids_route = self.get_routes().add_rows(ids_trip=ids_trip, 1822 ids_edges=ids_edges, 1823 costs=costs, 1824 colors=COLOR_SHORTEST_ROUTE * 1825 np.ones((len(ids_trip), 4), dtype=np.float32), 1826 ) 1827 1828 self.ids_route_shortest[ids_trip] = ids_route 1829 self.lengths_route_shortest[ids_trip] = costs 1830 #self.add_routes(ids_trip, ids_route) 1831 print ' exectime', time.clock()-exectime_start 1832 return ids_trip, ids_route 1833 1834 def get_speedprofile(self, id_trip): 1835 1836 points = self.parent.points 1837 ids_point = self.ids_points[id_trip] 1838 ids_point_edgeend = self.ids_points_edgeend[id_trip] 1839 ids_edge = self.get_routes().ids_edges[self.ids_route_matched[id_trip]] 1840 id_edge_current = -1 1841 n_edges = len(ids_edge) 1842 1843 positions_gps = [] 1844 times_gps = [] 1845 ids_edges_profile = [] 1846 ids_nodes_profile = [] 1847 1848 for id_point, coord, timestamp in zip(ids_point, points.coords[ids_point], points.timestamps[ids_point]): 1849 #get_pos_from_coord(id_edge_current, coord) 1850 1851 if id_point in ids_point_edgeend: 1852 ind = ids_point_edgeend.index(id_point) 1853 if ind == n_edges-1: 1854 id_edge_current = -1 1855 else: 1856 id_edge_current = ids_edge[ind+1] 1857 1858 def get_flows(self, is_shortest_path=False): 1859 """ 1860 Determine the total number of vehicles for each edge. 1861 returns ids_edge and flows 1862 """ 1863 ids_edges = self.get_routes.ids_edges 1864 counts = np.zeros(np.max(self.get_net().edges.get_ids())+1, int) 1865 1866 ids_trip = self.get_ids_selected() 1867 if not is_shortest_path: 1868 ids_route = self.ids_route_matched[ids_trip] 1869 else: 1870 ids_route = self.ids_route_shortest[ids_trip] 1871 inds_valid = np.flatnonzero(ids_route > 0) 1872 for id_trip, id_route in zip(ids_trip[inds_valid], ids_route[inds_valid]): 1873 counts[ids_edges[id_route][:]] += 1 1874 1875 ids_edge = np.flatnonzero(counts) 1876 1877 return ids_edge, counts[ids_edge].copy() 1878 1879 def get_ids_route_selected(self): 1880 # TODO: here we could append direct routes 1881 # print 'get_ids_route_selected' 1882 ids_route_matched = self.ids_route_matched[self.get_ids_selected()] 1883 ids_route_shortest = self.ids_route_shortest[self.get_ids_selected()] 1884 ids_route = np.concatenate([ids_route_matched[ids_route_matched >= 0], 1885 ids_route_shortest[ids_route_shortest >= 0]]) 1886 # print ' ids_route',ids_route 1887 # print ' ids_route_matched[ids_route_matched >= 0] ',ids_route_matched[ids_route_matched >= 0] 1888 return np.concatenate([ids_route_matched[ids_route_matched >= 0], ids_route_shortest[ids_route_shortest >= 0]]) 1889 # return ids_route_matched[ids_route_matched >= 0] 1890 1891 1892class GpsPoints(am.ArrayObjman): 1893 """ 1894 Contains data of points of a single trace. 1895 """ 1896 1897 def __init__(self, ident, mapmatching, **kwargs): 1898 # print 'PrtVehicles vtype id_default',vtypes.ids_sumo.get_id_from_index('passenger1') 1899 self._init_objman(ident=ident, 1900 parent=mapmatching, 1901 name='GPS Points', 1902 info='GPS points database.', 1903 version=0.0, 1904 **kwargs) 1905 1906 self._init_attributes() 1907 self._init_constants() 1908 1909 def _init_attributes(self): 1910 scenario = self.get_scenario() 1911 1912 # ident is the path id of the trace 1913 1914 # the actual numpy arrays are stored in .cols 1915 self.add_col(am.ArrayConf('longitudes', default=0.0, 1916 dtype=np.float32, 1917 groupnames=['parameters', ], 1918 perm='rw', 1919 name='Longitude', 1920 symbol='Lon', 1921 unit='deg', 1922 info='Longitude of point', 1923 )) 1924 1925 self.add_col(am.ArrayConf('latitudes', default=0.0, 1926 groupnames=['parameters', ], 1927 dtype=np.float32, 1928 perm='rw', 1929 name='Latitude', 1930 symbol='Lat', 1931 unit='deg', 1932 info='Latitude of point', 1933 )) 1934 1935 self.add_col(am.ArrayConf('altitudes', default=-100000.0, 1936 groupnames=['parameters', ], 1937 dtype=np.float32, 1938 perm='rw', 1939 name='Altitude', 1940 symbol='Alt', 1941 unit='m', 1942 info='Altitude of point', 1943 )) 1944 1945 self.add_col(am.ArrayConf('radii', 10.0, 1946 dtype=np.float32, 1947 groupnames=['parameters', ], 1948 perm='rw', 1949 name='Radius', 1950 unit='m', 1951 info='Point radius, representing the imprecision of the point, which depends on the recording device ane the environment.', 1952 )) 1953 1954 self.add_col(am.ArrayConf('coords', default=[0.0, 0.0, 0.0], 1955 groupnames=['parameters', ], 1956 dtype=np.float32, 1957 perm='rw', 1958 name='Coordinate', 1959 symbol='x,y,z', 1960 unit='m', 1961 info='Local 3D coordinate of point', 1962 )) 1963 1964 self.add_col(am.ArrayConf('timestamps', default=0.0, 1965 dtype=np.float, 1966 groupnames=['parameters', ], 1967 perm='r', 1968 name='timestamp', 1969 symbol='t', 1970 metatype='datetime', 1971 unit='s', 1972 digits_fraction=2, 1973 info='Time stamp of point in seconds after 01 January 1970.', 1974 )) 1975 # test: 1976 self.timestamps.metatype = 'datetime' 1977 self.timestamps.set_perm('r') 1978 1979 def set_trips(self, trips): 1980 self.add_col(am.IdsArrayConf('ids_trip', trips, 1981 groupnames=['state'], 1982 name='Trip ID', 1983 info='ID of trips to which this point belongs to.', 1984 )) 1985 # put in geometry filter 1986 # self.add_col(am.ArrayConf( 'are_inside_boundary', default =False, 1987 # groupnames = ['parameters',], 1988 # perm='r', 1989 # name = 'in boundary', 1990 # info = 'True if this the data point is within the boundaries of the road network.', 1991 # )) 1992 1993 def get_ids_selected(self): 1994 """ 1995 Returns point ids of selected traces 1996 """ 1997 # print 'GpsPoints.get_ids_selected' 1998 # print ' ??ids_points = ',self.select_ids(self.parent.trips.are_selected[self.ids_trip.get_value()] ) 1999 # TODO: why is this working??? do we need trips.ids_points???? 2000 return self.select_ids(self.parent.trips.are_selected[self.ids_trip.get_value()]) 2001 # return self.get_ids()#self.select_ids(self.get_ids() 2002 # return self.select_ids( 2003 2004 def get_scenario(self): 2005 return self.parent.get_scenario() 2006 2007 def get_coords(self): 2008 """ 2009 Returns an array of x,y coordinates of all points. 2010 """ 2011 return self.coords.get_value() 2012 2013 def get_times(self): 2014 """ 2015 Returns an array of time stamps of all points. 2016 """ 2017 return self.timestamp.get_value() 2018 2019 def get_interval(self): 2020 if len(self) == 0: 2021 return [0.0, 0.0] 2022 else: 2023 timestamps = self.timestamps.get_value() 2024 return [timestamps.min(), timestamps.max()] 2025 2026 def get_duration(self): 2027 ts, te = self.get_interval() 2028 return te-ts 2029 2030 def get_distance(self): 2031 v = self.get_coords() 2032 # return np.linalg.norm( self.cols.coords[0:-1] - self.cols.coords[1:] ) 2033 return np.sum(np.sqrt((v[1:, 0]-v[0:-1, 0])**2 + (v[1:, 1]-v[0:-1, 1])**2)) 2034 2035 def get_boundary(self): 2036 if len(self) == 0: 2037 return [0, 0, 0, 0] 2038 else: 2039 x_min, y_min = self.get_coords().min(0) 2040 x_max, y_max = self.get_coords().max(0) 2041 2042 return [x_min, y_min, x_max, y_max] 2043 2044 def project(self, proj=None, offset=None): 2045 if proj is None: 2046 proj, offset = self.parent.get_proj_and_offset() 2047 x, y = proj(self.longitudes.get_value(), self.latitudes.get_value()) 2048 self.get_coords()[:] = np.transpose(np.concatenate( 2049 ([x+offset[0]], [y+offset[1]], [self.altitudes.get_value()]), axis=0)) 2050 2051 2052class GpsPersons(am.ArrayObjman): 2053 2054 def __init__(self, ident, mapmatching, **kwargs): 2055 # print 'PrtVehicles vtype id_default',vtypes.ids_sumo.get_id_from_index('passenger1') 2056 self._init_objman(ident=ident, 2057 parent=mapmatching, 2058 name='GPS Persons', 2059 info='GPS person database.', 2060 **kwargs) 2061 2062 self._init_attributes() 2063 2064 def _init_attributes(self): 2065 trips = self.parent.trips 2066 2067 # TODO: add/update vtypes here 2068 self.add_col(SumoIdsConf('User', xmltag='id')) 2069 2070 self.add_col(am.IdlistsArrayConf('ids_trips', trips, 2071 #groupnames = ['_private'], 2072 name='Trip IDs', 2073 info="IDs of trips made by this vehicle. This is a collection of recorded trips associated with this person.", 2074 )) 2075 2076 self.add_col(am.ArrayConf('ids_genders', default=0, 2077 dtype=np.int, 2078 groupnames=['parameters'], 2079 name='gender', 2080 info='Gender of person.', 2081 )) 2082 2083 self.add_col(am.ArrayConf('years_birth', default=-1, 2084 dtype=np.int, 2085 groupnames=['parameters'], 2086 name='birth year', 2087 info='Year when person has been born.', 2088 )) 2089 2090 self.add_col(am.ArrayConf('ids_occupation', default=OCCUPATIONS['unknown'], 2091 dtype=np.int32, 2092 choices=OCCUPATIONS, 2093 groupnames=['parameters'], 2094 name='occupation', 2095 info='Tupe of occupation', 2096 )) 2097 2098 self.add_col(am.ArrayConf('are_frequent_user', False, 2099 dtype=np.int32, 2100 groupnames=['parameters'], 2101 name='frequent user', 2102 info='If true, this person is a frequent user of the recorded transport mode.', 2103 )) 2104 2105 self.add_col(am.ArrayConf('zips', -1, 2106 dtype=np.int32, 2107 groupnames=['parameters'], 2108 name='ZIP', 2109 info='ZIP code of persons home.', 2110 )) 2111 2112 def make(self, id_sumo, **kwargs): 2113 if self.ids_sumo.has_index(id_sumo): 2114 id_pers = self.ids_sumo.get_id_from_index(id_sumo) 2115 #self.set_row(id_pers, **kwargs) 2116 return id_pers 2117 else: 2118 id_pers = self.add_row(ids_sumo=id_sumo, 2119 ids_trips=kwargs.get('ids_trip', None), 2120 ids_genders=kwargs.get('id_gender', None), 2121 years_birth=kwargs.get('year_birth', None), 2122 ids_occupation=kwargs.get('id_occupation', None), 2123 are_frequent_user=kwargs.get('is_frequent_user', None), 2124 zips=kwargs.get('zip', None), 2125 ) 2126 2127 return id_pers 2128 2129 2130class Mapmatching(cm.BaseObjman): 2131 def __init__(self, ident, demand=None, 2132 name='Mapmatching', info='Mapmatching functionality.', 2133 **kwargs): 2134 2135 self._init_objman(ident=ident, parent=demand, 2136 name=name, info=info, **kwargs) 2137 2138 attrsman = self.set_attrsman(cm.Attrsman(self)) 2139 2140 self._init_attributes() 2141 self._init_constants() 2142 2143 def get_scenario(self): 2144 return self.parent.parent 2145 2146 def clear_all(self): 2147 self.trips.clear() 2148 self.points.clear() 2149 self.persons.clear() 2150 self._init_constants() 2151 2152 def clear_routes(self): 2153 self.trips.clear_routes() 2154 2155 def _init_attributes(self): 2156 print 'Mapmatching._init_attributes' 2157 attrsman = self.get_attrsman() 2158 2159 self.points = attrsman.add(cm.ObjConf(GpsPoints('points', self))) 2160 self.trips = attrsman.add(cm.ObjConf(GpsTrips('trips', self))) 2161 self.points.set_trips(self.trips) 2162 self.persons = attrsman.add(cm.ObjConf(GpsPersons('persons', self))) 2163 2164 def _init_constants(self): 2165 self._proj = None 2166 self._segvertices_xy = None 2167 self._distancesmap = None 2168 self._accesslevelsmap = None 2169 2170 attrsman = self.get_attrsman() 2171 attrsman.do_not_save_attrs(['_segvertices_xy', '_proj', '_distancesmap', '_accesslevelsmap']) 2172 2173 def get_proj_and_offset(self): 2174 if self._proj is None: 2175 net = self.get_scenario().net 2176 proj_params = str(net.get_projparams()) 2177 # try: 2178 self._proj = pyproj.Proj(proj_params) 2179 # except: 2180 # proj_params ="+proj=utm +zone=32 +ellps=WGS84 +datum=WGS84 +units=m +no_defs" 2181 # self._proj = pyproj.Proj(self.proj_params) 2182 2183 self._offset = net.get_offset() 2184 2185 return self._proj, self._offset 2186 2187 def get_segvertices_xy(self): 2188 if self._segvertices_xy is None: 2189 self._segvertices_xy = self.get_scenario().net.edges.get_segvertices_xy() 2190 2191 return self._segvertices_xy 2192 2193 def get_distancesmap(self, is_check_lanes=False): 2194 """ 2195 Returns a dictionary where key is id_mode and 2196 value is a distance-lookup table, mapping id_edge to edge distance 2197 """ 2198 # print 'get_distancesmap',self._distancesmap is None 2199 if self._distancesmap is None: 2200 2201 vtypes = self.get_scenario().demand.vtypes 2202 edges = self.get_scenario().net.edges 2203 ids_vtype = self.trips.ids_vtype[self.trips.get_ids_selected()] 2204 ids_mode = vtypes.ids_mode[ids_vtype] 2205 # print ' ids_mode',ids_mode 2206 2207 self._distancesmap = {} 2208 for id_mode in set(ids_mode): 2209 #ids_vtype_mode = vtypes.select_by_mode(id_mode) 2210 self._distancesmap[id_mode] = edges.get_distances(id_mode=id_mode, 2211 is_check_lanes=is_check_lanes, 2212 ) 2213 # print ' len(self._distancesmap)',len(self._distancesmap) 2214 return self._distancesmap 2215 2216 def get_accesslevelsmap(self): 2217 """ 2218 Returns a dictionary where key is id_mode and 2219 value is a distance-lookup table, mapping id_edge to edge distance 2220 """ 2221 2222 if self._accesslevelsmap is None: 2223 vtypes = self.get_scenario().demand.vtypes 2224 edges = self.get_scenario().net.edges 2225 ids_vtype = self.trips.ids_vtype[self.trips.get_ids_selected()] 2226 ids_mode = vtypes.ids_mode[ids_vtype] 2227 2228 self._accesslevelsmap = {} 2229 for id_mode in set(ids_mode): 2230 #ids_vtype_mode = vtypes.select_by_mode(id_mode) 2231 self._accesslevelsmap[id_mode] = edges.get_accesslevels(id_mode) 2232 return self._accesslevelsmap 2233 2234 2235class Matchresults(cm.BaseObjman): 2236 def __init__(self, ident, mapmatching, 2237 name='Mapmatching results', 2238 info='Results of mapmatching and successive analyses.', 2239 **kwargs): 2240 2241 # make results a child of process or of wxgui 2242 # use these objects to access matched trips 2243 2244 self._init_objman(ident, parent=mapmatching, name=name, 2245 info=info, **kwargs) 2246 attrsman = self.set_attrsman(cm.Attrsman(self)) 2247 2248 self._init_attributes() 2249 2250 def _init_attributes(self): 2251 attrsman = self.get_attrsman() 2252 mapmatching = self.parent 2253 self.routesresults = attrsman.add(cm.ObjConf(Routesresults('routesresults', 2254 self, mapmatching.trips.routes), 2255 groupnames=['Route results'], 2256 )) 2257 2258 # add trip results from all demand objects 2259 print 'Matchresults._init_attributes' 2260 2261 def get_scenario(self): 2262 return self.parent.get_scenario() 2263 2264 def clear_all(self): 2265 self.clear() 2266 2267 def save(self, filepath=None, is_not_save_parent=True): 2268 if filepath is None: 2269 self.get_scenario().get_rootfilepath()+'.mmatch.obj' 2270 # parent will not be saved because in no_save set 2271 cm.save_obj(self, filepath, is_not_save_parent=is_not_save_parent) 2272 2273 2274class Nodeinfolists(am.ListArrayConf): 2275 def format_value(self, _id, show_unit=False, show_parentesis=False): 2276 if show_unit: 2277 unit = ' '+self.format_unit(show_parentesis) 2278 else: 2279 unit = '' 2280 # return repr(self[_id])+unit 2281 2282 #self.min = minval 2283 #self.max = maxval 2284 #self.digits_integer = digits_integer 2285 #self.digits_fraction = digits_fraction 2286 val = self[_id] 2287 tt = type(val) 2288 2289 if tt in (np.int, np.int32, np.float64): 2290 return str(val)+unit 2291 2292 elif tt in (np.float, np.float32, np.float64): 2293 if hasattr(self, 'digits_fraction'): 2294 digits_fraction = self.digits_fraction 2295 else: 2296 digits_fraction = 3 2297 s = "%."+str(digits_fraction)+"f" 2298 return s % (val)+unit 2299 2300 else: 2301 return str(val)+unit 2302 2303 2304class Accesslevellists(am.ListArrayConf): 2305 def format_value(self, _id, show_unit=False, show_parentesis=False): 2306 if show_unit: 2307 unit = ' '+self.format_unit(show_parentesis) 2308 else: 2309 unit = '' 2310 # return repr(self[_id])+unit 2311 2312 #self.min = minval 2313 #self.max = maxval 2314 #self.digits_integer = digits_integer 2315 #self.digits_fraction = digits_fraction 2316 val = self[_id] 2317 tt = type(val) 2318 2319 if tt in (np.int, np.int32, np.float64): 2320 return str(val)+unit 2321 2322 elif tt in (np.float, np.float32, np.float64): 2323 if hasattr(self, 'digits_fraction'): 2324 digits_fraction = self.digits_fraction 2325 else: 2326 digits_fraction = 3 2327 s = "%."+str(digits_fraction)+"f" 2328 return s % (val)+unit 2329 2330 else: 2331 return str(val)+unit 2332 2333 2334class Routesresults(am.ArrayObjman): 2335 def __init__(self, ident, parent, routes, 2336 is_add_default=True, 2337 name='Route results', 2338 info='Table with results from analysis of each matched route.', 2339 **kwargs): 2340 2341 self._init_objman(ident=ident, 2342 parent=parent, # main results object 2343 info=info, 2344 name=name, 2345 **kwargs) 2346 2347 # self.add(cm.AttrConf( 'datapathkey',datapathkey, 2348 # groupnames = ['_private'], 2349 # name = 'data pathkey', 2350 # info = "key of data path", 2351 # )) 2352 2353 self.add_col(am.IdsArrayConf('ids_route', routes, 2354 groupnames=['state'], 2355 is_index=True, 2356 name='ID route', 2357 info='ID of matched route.', 2358 )) 2359 self._init_attributes() 2360 2361 def _init_attributes(self): 2362 2363 self.add_col(am.ListArrayConf('intervallists', 2364 groupnames=['_private'], 2365 perm='rw', 2366 name='Intervals', 2367 unit='s', 2368 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.', 2369 )) 2370 2371 self.add_col(am.ListArrayConf('edgedurationlists', 2372 groupnames=['_private'], 2373 perm='rw', 2374 name='Edge durations', 2375 unit='s', 2376 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.', 2377 )) 2378 2379 self.add_col(am.ListArrayConf('nodedurationlists', 2380 groupnames=['_private'], 2381 perm='rw', 2382 name='Node durations', 2383 unit='s', 2384 info='Duration on each intersection of the route. These times are obtained from timestamps of GPS points which can be perpendicularly projected to the respective route edge.', 2385 )) 2386 2387 self.add_col(Nodeinfolists('nodeinfolists', 2388 groupnames=['_private'], 2389 perm='rw', 2390 name='Node info', 2391 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).', 2392 )) 2393 2394 # self.nodetypemap = {\ 2395 # "priority":0, 2396 # "traffic_light":1, 2397 # "right_before_left":2, 2398 # "unregulated":3, 2399 # "priority_stop":4, 2400 # "traffic_light_unregulated":5, 2401 # "allway_stop":6, 2402 # "rail_signal":7, 2403 # "zipper":8, 2404 # "traffic_light_right_on_red":9, 2405 # "rail_crossing":10, 2406 # "dead_end":11, 2407 # } 2408 2409 self.add_col(Accesslevellists('accesslevellists', 2410 groupnames=['_private'], 2411 perm='rw', 2412 name='Accesslevels', 2413 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.', 2414 )) 2415 2416 # self.add_col(am.IdsArrayConf('ids_route_shortest', self.get_routes(), 2417 # groupnames = ['results'], 2418 # name = 'ID shortest route', 2419 # info = 'Route ID of shortest route. Shortest route is connecting the first matched edge and the final matched edge.', 2420 # )) 2421 2422 # self.add_col(am.ArrayConf('lengths_route_shortest', default = -1.0, 2423 # dtype = np.float32, 2424 # groupnames = ['results'], 2425 # name = 'Shortest length', 2426 # symbol = 'L_{S}', 2427 # unit = 'm', 2428 # info = 'Length of the shortest route. Shortest route is connecting the first matched edge and the final matched edge.', 2429 # )) 2430 2431 2432class Routesanalysis(Process): 2433 def __init__(self, ident, mapmatching, results=None, logger=None, **kwargs): 2434 print 'Routesanalysis.__init__' 2435 2436 # TODO: let this be independent, link to it or child?? 2437 2438 if results is None: 2439 self._results = mapmatching.Matchresults('matchresults', 2440 mapmatching) 2441 else: 2442 self._results = results 2443 2444 self._init_common(ident, 2445 parent=mapmatching, 2446 name='Birgillito Map matching', 2447 logger=logger, 2448 info='Birgillito Map matching.', 2449 ) 2450 2451 attrsman = self.set_attrsman(cm.Attrsman(self)) 2452 2453 self.width_buffer_max = attrsman.add(cm.AttrConf('width_buffer_max', kwargs.get('width_buffer_max', 40.0), 2454 groupnames=['options'], 2455 perm='rw', 2456 name='Max. buffer width', 2457 unit='m', 2458 info='GPS points are only valid if they are within the given maximum buffer width distant from a network edge.', 2459 )) 2460 2461 self.n_points_min = attrsman.add(cm.AttrConf('n_points_min', kwargs.get('n_points_min', 5), 2462 groupnames=['options'], 2463 perm='rw', 2464 name='Min. point number', 2465 info='Minimum number of valid GPS points. Only if this minimum number is reached, a map-matching attempt is performed.', 2466 )) 2467 2468 def get_results(self): 2469 return self._results 2470 2471 def do(self): 2472 mapmatching = self.parent 2473 trips = mapmatching.trips 2474 routes = trips.routes 2475 edges = mapmatching.get_scenario().net.edges 2476 points = mapmatching.points 2477 2478 ids_trips_sel = trips.get_ids_selected() 2479 ids_route_sel = trips.ids_route_matched[ids_trips_sel] 2480 inds_valid = np.flatnonzero(ids_route_sel > 0) 2481 2482 ids_route = ids_route_sel[inds_valid] 2483 ids_trips = ids_trips_sel[inds_valid] 2484 2485 ids_vtype = trips.ids_vtype[ids_trip] 2486 ids_mode = vtypes.ids_mode[ids_vtype] 2487 2488 for id_trip, id_route in zip(ids_trip[inds_valid], ids_route[inds_valid]): 2489 ids_point = trips.ids_points[id_trip] 2490 2491 2492class GpxParser(handler.ContentHandler): 2493 """Reads endomondo gpx xml file and parses lat,lon, ele and time. 2494 """ 2495 2496 def __init__(self, scenario, traces, log=None): 2497 2498 self._traces = traces 2499 2500 self._pointset = None 2501 self._id_point = None 2502 2503 self._lon = None 2504 self._lat = None 2505 self._ele = None 2506 self._time = None 2507 self._time_flag = False 2508 self._ele_flag = False 2509 self._ids_traces = [] 2510 2511 def get_ids_traces(self): 2512 return self._ids_traces 2513 2514 def get_traces(self): 2515 return self._traces 2516 2517 def startElement(self, name, attrs): 2518 if self._pointset is None: 2519 if name == 'trkseg': 2520 # start new trace with new pointse 2521 if len(self._traces) == 0: 2522 id_trip = '0' 2523 else: 2524 id_trip = str(np.max(self._traces.get_ids())+1) 2525 print 'startElement id_trip', id_trip, self._traces.get_ids(), self._traces 2526 self._pointset = PointSet(id_trip, self._traces) 2527 self._traces.set_row(id_trip, pointsets=self._pointset) 2528 self._ids_traces.append(id_trip) 2529 2530 if self._id_point is None: 2531 if name == 'trkpt': 2532 2533 self._id_point = self._pointset.suggest_id() 2534 # print 'startElement trkpt',self._id_point 2535 self._lon = float(attrs['lon']) 2536 self._lat = float(attrs['lat']) 2537 2538 if name == 'time': 2539 self._time_flag = True 2540 2541 if name == 'ele': 2542 self._ele_flag = True 2543 2544 def characters(self, content): 2545 if self._time_flag: 2546 # print ' got time',content[:-1] 2547 2548 # 2013-05-30T16:42:33Z 2549 self._time = calc_seconds(content[:-1], 2550 sep_date_clock='T', sep_date='-', sep_clock=':', 2551 is_float=False) 2552 2553 if self._ele_flag: 2554 self._ele = float(content) 2555 2556 def endElement(self, name): 2557 if name == 'ele': 2558 self._ele_flag = False 2559 2560 if name == 'time': 2561 self._time_flag = False 2562 2563 if name == 'trkpt': 2564 # print 'endElement:add point to pointset',self._id_point,self._lon,self._lat,self._ele,self._time 2565 2566 if (self._time is not None)\ 2567 & (self._id_point is not None)\ 2568 & (self._lon is not None) & (self._lat is not None)\ 2569 & (self._pointset is not None): 2570 if (self._ele is not None): 2571 self._pointset.set_row(self._id_point, 2572 lon=self._lon, 2573 lat=self._lat, 2574 altitude=self._ele, 2575 timestamp=self._time, 2576 ids_point=self._id_point 2577 ) 2578 else: 2579 # no altutude provided, may be corrected later 2580 self._pointset.set_row(self._id_point, 2581 lon=self._lon, 2582 lat=self._lat, 2583 timestamp=self._time, 2584 ids_point=self._id_point 2585 ) 2586 2587 self._id_point = None 2588 self._lon = None 2589 self._lat = None 2590 self._lat = None 2591 self._ele = None 2592 2593 self._ele_flag = False 2594 self._time_flag = False 2595