Merge olusanya to mal23-car-ownership

Scripts/assignment/assignment_period.py

@@ -200,7 +200,7 @@ def calc_transit_cost(self, fares, peripheral_cost, mapping):
                 "trip_part_transit_work", "board_cost")
             # If the number of visits is less than 1, there seems to
             # be an easy way to avoid visiting this transit zone
-            has_visited[transit_zone] = (nr_visits >= 1)
+            has_visited[transit_zone] = (nr_visits > 0.99)
         for centroid in network.centroids():
             # Add transit zone of destination to visited
             has_visited[centroid.label][:, mapping[centroid.number]] = True

Scripts/datatypes/histogram.py

@@ -13,7 +13,7 @@ def __init__(self):
         self.histogram = pandas.Series(0, index)
 
     def add(self, dist):
-        self.histogram.iloc[numpy.searchsorted(self._u, dist, "right")] += 1
+        self.histogram.iat[numpy.searchsorted(self._u, dist, "right")] += 1
 
     def count_tour_dists(self, tours, dists):
         self.histogram[:], _ = numpy.histogram(dists, intervals, weights=tours)

Scripts/datatypes/purpose.py

@@ -1,6 +1,7 @@
 import numpy
 import pandas
 
+import parameters.zone as param
 from parameters.destination_choice import secondary_destination_threshold
 import models.logit as logit
 import models.generation as generation
@@ -35,25 +36,13 @@ def __init__(self, specification, zone_data, resultdata=None):
         self.dest = specification["dest"]
         self.area = specification["area"]
         self.sources = []
-        if self.area == "metropolitan":
-            l = 0
-            m = zone_data.first_surrounding_zone
-            u = zone_data.first_peripheral_zone
-        if self.area == "peripheral":
-            l = zone_data.first_peripheral_zone
-            m = None
-            u = zone_data.nr_zones
-        if self.area == "all":
-            l = 0
-            m = zone_data.first_surrounding_zone
-            u = zone_data.nr_zones
-        if self.area == "external":
-            l = zone_data.first_external_zone
-            m = None
-            u = None
-        self.bounds = slice(l, u)
-        self.lbounds = slice(l, m)
-        self.ubounds = slice(m, u)
+        zone_numbers = zone_data.zone_numbers
+        zone_intervals = param.purpose_areas[self.area]
+        self.bounds = slice(*zone_numbers.searchsorted(
+            [zone_intervals[0], zone_intervals[-1]]))
+        sub_intervals = zone_numbers[self.bounds].searchsorted(zone_intervals)
+        self.sub_bounds = [slice(sub_intervals[i-1], sub_intervals[i])
+            for i in range(1, len(sub_intervals))]
         self.zone_data = zone_data
         self.resultdata = resultdata
         self.model = None

Scripts/datatypes/tour.py

@@ -4,7 +4,7 @@
 import parameters.car as param
 import parameters.zone as zone_param
 from parameters.assignment import assignment_classes, vot_inv
-from parameters.impedance_transformation import divided_classes, trips_month
+from parameters.impedance_transformation import divided_classes, transit_trips_per_month
 
 
 class Tour:
@@ -200,14 +200,12 @@ def calc_cost(self, impedance):
 
     def _get_cost(self, impedance, mtx_type):
         """Get cost and time components from tour dest choice."""
-        if self.mode in divided_classes:
-            ass_class = "{}_{}".format(
-                self.mode, assignment_classes[self.purpose.name])
-        else:
-            ass_class = self.mode
+        demand_type = assignment_classes[self.purpose.name]
+        ass_class = ("{}_{}".format(self.mode, demand_type)
+            if self.mode in divided_classes else self.mode)
         cost = 0
         try:
-            if assignment_classes[self.purpose_name] == "work":
+            if demand_type == "work":
                 departure_imp = impedance["aht"][mtx_type][ass_class]
                 sec_dest_imp = impedance["iht"][mtx_type][ass_class]
                 return_imp = impedance["iht"][mtx_type][ass_class]
@@ -229,8 +227,10 @@ def _get_cost(self, impedance, mtx_type):
             pass
         # scale transit costs from month to day
         if self.mode == "transit" and mtx_type == "cost":
-            idx = int(self.position[0] > self.purpose.zone_data.first_surrounding_zone)
-            cost /= trips_month[ass_class][idx]
+            idx = numpy.searchsorted(
+                [bounds.stop for bounds in self.purpose.sub_bounds],
+                self.position[0], side="right")
+            cost /= transit_trips_per_month[self.purpose.area][demand_type][idx]
         return cost
 
     def __str__(self):

Scripts/dev-config.json

@@ -1,5 +1,5 @@
 {
-    "HELMET_VERSION": "v4.1.0-beta.3",
+    "HELMET_VERSION": "v4.1.1",
     "LOG_LEVEL": "INFO",
     "LOG_FORMAT": "TEXT",
     "SCENARIO_NAME": "test",

Scripts/models/logit.py

@@ -29,8 +29,7 @@ def __init__(self, zone_data, purpose, resultdata, is_agent_model):
         self.resultdata = resultdata
         self.purpose = purpose
         self.bounds = purpose.bounds
-        self.lbounds = purpose.lbounds
-        self.ubounds = purpose.ubounds
+        self.sub_bounds = purpose.sub_bounds
         self.zone_data = zone_data
         self.dest_exps = {}
         self.mode_exps = {}
@@ -112,13 +111,12 @@ def _add_constant(self, utility, b):
         """
         try: # If only one parameter
             utility += b
-        except ValueError: # Separate params for cap region and surrounding
-            if utility.ndim == 1: # 1-d array calculation
-                utility[self.lbounds] += b[0]
-                utility[self.ubounds] += b[1]
-            else: # 2-d matrix calculation
-                utility[self.lbounds, :] += b[0]
-                utility[self.ubounds, :] += b[1]
+        except ValueError: # Separate sub-region parameters
+            for i, bounds in enumerate(self.sub_bounds):
+                if utility.ndim == 1: # 1-d array calculation
+                    utility[bounds] += b[i]
+                else: # 2-d matrix calculation
+                    utility[bounds, :] += b[i]
 
     def _add_impedance(self, utility, impedance, b):
         """Adds simple linear impedances to utility.
@@ -139,9 +137,9 @@ def _add_impedance(self, utility, impedance, b):
         for i in b:
             try: # If only one parameter
                 utility += b[i] * impedance[i]
-            except ValueError: # Separate params for cap region and surrounding
-                utility[self.lbounds, :] += b[i][0] * impedance[i][self.lbounds, :]
-                utility[self.ubounds, :] += b[i][1] * impedance[i][self.ubounds, :]
+            except ValueError: # Separate sub-region parameters
+                for j, bounds in enumerate(self.sub_bounds):
+                    utility[bounds, :] += b[i][j] * impedance[i][bounds, :]
         return utility
 
     def _add_log_impedance(self, exps, impedance, b):
@@ -168,11 +166,10 @@ def _add_log_impedance(self, exps, impedance, b):
         for i in b:
             try: # If only one parameter
                 exps *= numpy.power(impedance[i] + 1, b[i])
-            except ValueError: # Separate params for cap region and surrounding
-                exps[self.lbounds, :] *= numpy.power(
-                    impedance[i][self.lbounds, :] + 1, b[i][0])
-                exps[self.ubounds, :] *= numpy.power(
-                    impedance[i][self.ubounds, :] + 1, b[i][1])
+            except ValueError: # Separate sub-region parameters
+                for j, bounds in enumerate(self.sub_bounds):
+                    exps[bounds, :] *= numpy.power(
+                        impedance[i][bounds, :] + 1, b[i][j])
         return exps
 
     def _add_zone_util(self, utility, b, generation=False):
@@ -196,15 +193,13 @@ def _add_zone_util(self, utility, b, generation=False):
         for i in b:
             try: # If only one parameter
                 utility += b[i] * zdata.get_data(i, self.bounds, generation)
-            except ValueError: # Separate params for cap region and surrounding
-                data_cap_region = zdata.get_data(i, self.lbounds, generation)
-                data_surrounding = zdata.get_data(i, self.ubounds, generation)
-                if utility.ndim == 1: # 1-d array calculation
-                    utility[self.lbounds] += b[i][0] * data_cap_region
-                    utility[self.ubounds] += b[i][1] * data_surrounding
-                else: # 2-d matrix calculation
-                    utility[self.lbounds, :] += b[i][0] * data_cap_region
-                    utility[self.ubounds, :] += b[i][1] * data_surrounding
+            except ValueError: # Separate sub-region parameters
+                for j, bounds in enumerate(self.sub_bounds):
+                    data = zdata.get_data(i, bounds, generation)
+                    if utility.ndim == 1: # 1-d array calculation
+                        utility[bounds] += b[i][j] * data
+                    else: # 2-d matrix calculation
+                        utility[bounds, :] += b[i][j] * data
         return utility
 
     def _add_sec_zone_util(self, utility, b, orig=None, dest=None):
@@ -344,8 +339,8 @@ def calc_individual_prob(self, mod_mode, dummy):
         try:
             self.mode_exps[mod_mode] *= numpy.exp(b)
         except ValueError:
-            self.mode_exps[mod_mode][self.lbounds] *= numpy.exp(b[0])
-            self.mode_exps[mod_mode][self.ubounds] *= numpy.exp(b[1])
+            for i, bounds in enumerate(self.sub_bounds):
+                self.mode_exps[mod_mode][bounds] *= numpy.exp(b[i])
         mode_expsum = numpy.zeros_like(self.mode_exps[mod_mode])
         for mode in self.mode_choice_param:
             mode_expsum += self.mode_exps[mode]
@@ -399,8 +394,9 @@ def calc_individual_mode_prob(self, is_car_user, zone):
             # Convert utility into euros
             money_utility = 1 / b
         except TypeError:
-            # Separate params for cap region and surrounding
-            money_utility = 1 / b[0] if self.lbounds.stop < zone else 1 / b[1]
+            # Separate sub-region parameters
+            money_utility = (1 / b[0] if self.sub_bounds[0].stop < zone
+                else 1 / b[1])
         money_utility /= self.mode_choice_param["car"]["log"]["logsum"]
         accessibility = -money_utility * logsum
         return probs, accessibility

Scripts/parameters/assignment.py

@@ -362,7 +362,7 @@
 years_average_day_factor = 0.85
 # Factor for converting day traffic into 7:00-22:00 traffic
 share_7_22_of_day = 0.9
-# Noice zone width as function of start noise
+# Noise zone width as function of start noise
 noise_zone_width = {
     (0, 55): lambda x: 5,
     (55, 65): lambda x: 10 + 31./10*x,

Scripts/parameters/impedance_transformation.py

@@ -1,8 +1,18 @@
 ### IMPEDANCE TRANSFORMATION PARAMETERS ###
 
-trips_month = {
-    "transit_work": (60.0, 44.0),
-    "transit_leisure": (30.0, 30.0),
+transit_trips_per_month = {
+    "metropolitan": {
+        "work": (60.0, 44.0),
+        "leisure": (30.0, 30.0),
+    },
+    "peripheral": {
+        "work": (44.0,),
+        "leisure": (30.0,),
+    },
+    "all": {
+        "work": (60.0, 44.0),
+        "leisure": (30.0, 30.0),
+    },
 }
 
 impedance_share = {

Scripts/parameters/zone.py

@@ -93,6 +93,15 @@
         "area": "all",
     },
 )
+# Tour purpose zone intervals
+# Some demand models have separate sub-region parameters,
+# hence need sub-intervals defined.
+purpose_areas = {
+    "metropolitan": (0, 6000, 16000),
+    "peripheral": (16000, 31000),
+    "all": (0, 6000, 31000),
+    "external": (31031, 40000),
+}
 areas = {
     "helsinki_cbd": (0, 999),
     "helsinki_other": (1000, 1999),

Scripts/tests/unit/test_logit.py

@@ -45,8 +45,7 @@ class Purpose:
             },
         }
         pur.bounds = slice(0, 9)
-        pur.lbounds = slice(0, 7)
-        pur.ubounds = slice(7, 9)
+        pur.sub_bounds = [slice(0, 7), slice(7, 9)]
         pur.zone_numbers = METROPOLITAN_ZONES
         for i in ("hw", "hc", "hu", "hs", "ho"):
             pur.name = i

Scripts/transform/impedance_transformer.py

@@ -1,7 +1,7 @@
 from collections import defaultdict
 import numpy
 
-from parameters.impedance_transformation import impedance_share, divided_classes, trips_month
+import parameters.impedance_transformation as param
 from parameters.assignment import assignment_classes
 
 
@@ -32,9 +32,10 @@ def transform(self, purpose, impedance):
         cols = (purpose.bounds if purpose.name == "hoo"
             else slice(0, purpose.zone_data.nr_zones))
         day_imp = {}
+        impedance_share = param.impedance_share
         for mode in impedance_share[purpose.name]:
             day_imp[mode] = defaultdict(float)
-            if mode in divided_classes:
+            if mode in param.divided_classes:
                 ass_class = "{}_{}".format(
                     mode, assignment_classes[purpose.name])
             else:
@@ -47,9 +48,11 @@ def transform(self, purpose, impedance):
                         day_imp[mode][mtx_type] += share[0] * imp[rows, cols]
                         day_imp[mode][mtx_type] += share[1] * imp[cols, rows].T
         # transit cost to eur per day
-        transit_class = "{}_{}".format("transit", assignment_classes[purpose.name])
+        trips_month = (param.transit_trips_per_month
+            [purpose.area][assignment_classes[purpose.name]])
         trips_per_month = numpy.full_like(
-            day_imp["transit"]["cost"], trips_month[transit_class][0])
-        trips_per_month[purpose.ubounds, :] = trips_month[transit_class][1]
+            day_imp["transit"]["cost"], trips_month[0])
+        for i in range(1, len(purpose.sub_bounds)):
+            trips_per_month[purpose.sub_bounds[i], :] = trips_month[i]
         day_imp["transit"]["cost"] /= trips_per_month
         return day_imp

Scripts/utils/calc_noise.py

@@ -33,28 +33,27 @@ def calc_noise(self, link):
         float
             Noise zone width (m)
         """
-        traffic = sum([link[mode] for mode in self.light_modes])
+        traffic = max(sum([link[mode] for mode in self.light_modes]), 0.01)
         rlink = link.reverse_link
         if rlink is None:
             reverse_traffic = 0
         else:
             reverse_traffic = sum([rlink[mode] for mode in self.light_modes])
         cross_traffic = (param.years_average_day_factor
-                            * param.share_7_22_of_day
-                            * (traffic+reverse_traffic))
+                         * param.share_7_22_of_day
+                         * (traffic+reverse_traffic))
         heavy = sum([link[mode] for mode in self.heavy_modes])
-        traffic = max(traffic, 0.01)
         heavy_share = heavy / (traffic+heavy)
 
         # Calculate speed
         link = self.morning_network.link(link.i_node, link.j_node)
         rlink = link.reverse_link
         if reverse_traffic > 0:
             speed = (60 * 2 * link.length
-                        / (link[self.car_morning]+rlink[self.car_morning]))
+                     / (link[self.car_morning]+rlink[self.car_morning]))
         else:
             speed = (0.3*(60*link.length/link[self.car_morning])
-                        + 0.7*link.data2)
+                     + 0.7*link.data2)
         speed = max(speed, 50.0)
 
         # Calculate start noise
@@ -75,4 +74,6 @@ def calc_noise(self, link):
             if interval[0] <= start_noise < interval[1]:
                 zone_width = func[interval](start_noise - interval[0])
                 break
+        else:
+            raise ValueError("No noise interval found for link {}".format(link.id))
         return zone_width