--- a/python/traffic_engineering.py	Thu Apr 18 15:29:33 2013 -0400
+++ b/python/traffic_engineering.py	Fri Dec 05 12:13:53 2014 -0500
@@ -12,6 +12,20 @@
 # Simulation
 #########################
 
+def generateTimeHeadways(meanTimeHeadway, simulationTime):
+    '''Generates the time headways between arrivals 
+    given the meanTimeHeadway and the negative exponential distribution
+    over a time interval of length simulationTime (assumed to be in same time unit as headway'''
+    from random import expovariate
+    headways = []
+    totalTime = 0
+    flow = 1/meanTimeHeadway
+    while totalTime < simulationTime:
+        h = expovariate(flow)
+        headways.append(h)
+        totalTime += h
+    return headways
+
 class Vehicle:
     '''Generic vehicle class
     1D coordinates for now'''
@@ -159,11 +173,16 @@
             self.types = types
             self.proportions = proportions
             self.equivalents = equivalents
-            self.nLanes = nLanes # unused
+            self.nLanes = nLanes
         else:
             print('Proportions do not sum to 1')
             pass
 
+    def checkProtected(self, opposedThroughMvt):
+        '''Checks if this left movement should be protected,
+        ie if one of the main two conditions on left turn is verified'''
+        return self.volume >= 200 or self.volume*opposedThroughMvt.volume/opposedThroughMvt.nLanes > 50000
+
     def getPCUVolume(self):
         '''Returns the passenger-car equivalent for the input volume'''
         v = 0
@@ -182,15 +201,6 @@
     def getTVUVolume(self):
         return self.mvtEquivalent*self.volume.getPCUVolume()    
 
-class IntersectionApproach: # should probably not be used
-    def __init__(self, leftTurnVolume, throughVolume, rightTurnVolume):
-        self.leftTurnVolume = leftTurnVolume
-        self.throughVolume = throughVolume
-        self.rightTurnVolume = rightTurnVolume
-
-    def getTVUVolume(self, leftTurnEquivalent = 1, throughEquivalent = 1, rightTurnEquivalent = 1):
-        return self.leftTurnVolume.getPCUVolume()*leftTurnEquivalent+self.throughVolume.getPCUVolume()*throughEquivalent+self.rightTurnVolume.getPCUVolume()*rightTurnEquivalent
-
 class LaneGroup:
     '''Class that represents a group of mouvements'''
 
@@ -204,11 +214,6 @@
     def getCharge(self, saturationVolume):
         return self.getTVUVolume()/(self.nLanes*saturationVolume)
 
-def checkProtectedLeftTurn(leftMvt, opposedThroughMvt):
-    '''Checks if one of the main two conditions on left turn is verified
-    The lane groups should contain left and through movement'''
-    return leftMvt.volume >= 200 or leftMvt.volume*opposedThroughMvt.volume/opposedThroughMvt.nLanes > 50000
-
 def optimalCycle(lostTime, criticalCharge):
     return (1.5*lostTime+5)/(1-criticalCharge)
 
@@ -226,9 +231,7 @@
         self.saturationVolume = saturationVolume
 
     def computeCriticalCharges(self):
-        self.criticalCharges = []
-        for phase in self.phases:
-            self.criticalCharges.append(max([lg.getCharge(self.saturationVolume) for lg in phase]))
+        self.criticalCharges = [max([lg.getCharge(self.saturationVolume) for lg in phase]) for phase in self.phases]
         self.criticalCharge = sum(self.criticalCharges)
         
     def computeOptimalCycle(self):
@@ -242,7 +245,7 @@
         return self.C
 
     def computeEffectiveGreen(self):
-        from numpy import round
+        #from numpy import round
         #self.computeCycle() # in case it was not done before
         effectiveGreenTime = self.C-self.lostTime
         self.effectiveGreens = [round(c*effectiveGreenTime/self.criticalCharge,1) for c in self.criticalCharges]
@@ -263,6 +266,15 @@
     '''Computes the uniform delay'''
     return 0.5*cycleLength*(1-float(effectiveGreen)/cycleLength)/(1-float(effectiveGreen*saturationDegree)/cycleLength)
 
+def overflowDelay(T, X, c, k=0.5, I=1):
+    '''Computes the overflow delay (HCM)
+    T in hours
+    c capacity of the lane group
+    k default for fixed time signal
+    I=1 for isolated intersection (Poisson arrival)'''
+    from math import sqrt
+    return 900*T*(X - 1 + sqrt((X - 1)**2 + 8*k*I*X/(c*T)))
+
 #########################
 # misc
 #########################