In [1]:
import os
import pandas as pd
import numpy as np
import sys
import shapefile
import operator

Read shapefiles


In [2]:
sf_link = os.path.join('directed', 'new_link')
sf_node = os.path.join('directed', 'new_node')

link_sf = shapefile.Reader(sf_link)
link_records = link_sf.records()
node_sf = shapefile.Reader(sf_node)
node_records = node_sf.records()

for data in link_records:
    print(data)
    break


['8290', '1', '2114', '25', '1', '0.04', '0', '8290', 1388147.792, 448212.542, 1387767.306, 449855.597]

In [3]:
link_sf.fields


Out[3]:
[('DeletionFlag', 'C', 1, 0),
 ['linkID', 'C', 80, 0],
 ['N1', 'C', 80, 0],
 ['N2', 'C', 80, 0],
 ['spd', 'C', 80, 0],
 ['lanes', 'C', 80, 0],
 ['l', 'C', 80, 0],
 ['isOD', 'C', 80, 0],
 ['note', 'C', 80, 0],
 ['X1', 'N', 11, 3],
 ['Y1', 'N', 11, 3],
 ['X2', 'N', 11, 3],
 ['Y2', 'N', 11, 3]]

In [4]:
for data in node_records:
    print(data)
    break


[1, 0]

In [5]:
node_sf.fields


Out[5]:
[('DeletionFlag', 'C', 1, 0), ['ID', 'N', 9, 0], ['isOD', 'N', 9, 0]]

In [6]:
# check if some O directly link to some D, if 0 it is good!
node_isod ={}
for data in node_records:
    node_isod[data[0]] = data[1]
_count = 0
for data in link_records:
    if node_isod[int(data[1])] and node_isod[int(data[2])]:
        print(data)
#         print("{}, {}".format(node_isod[int(data[1])], node_isod[int(data[2])]))
        _count += 1
print(_count)


0

In [7]:
# check if some link length=0, change to 0.01?
_count = 0
for data in link_records:
    if float(data[5]) == 0:
        data[5] = '0.01'
#         print(data)
        _count += 1
print(_count)


801

In [8]:
_count = 0
for data in link_records:
    if float(data[5]) == 0:
        _count += 1
print(_count)


0

In [9]:
out_path = 'input_files_MckeesRocks_SPC'

config_file_name = 'config.conf'
link_file_name = 'MNM_input_link'
node_file_name = 'MNM_input_node'
graph_file_name = 'Snap_graph'
od_file_name = 'MNM_input_od'
pathtable_file_name = 'path_table'
path_p_file_name = 'path_table_buffer'
demand_file_name = 'MNM_input_demand'

In [10]:
# ramps ID
ramp_IDs = ['2584', '2650', '2585', '2658', '2568', '2657', '2668', '2673', '2780', '2774', '2726', '2729', '2782', '2784',
            '2748', '2750', '2751', '4495', '4514', '4510', '4899', '4897', '4898', '4887', '4888', '4989', '4569', '4546',
            '4547', '4987', '4986', '4559', '4991', '4580', '4579', '4576', '4586', '4565', '5027', '5012', '5022', '5005',
            '5021', '5029', '4640', '5032', '5034', '5037', '4823', '4824', '4841', '4844', '4853', '4857', '4863', '4859',
            '4860', '4849', '4848', '4862', '4861', '4858', '4832', '4851', '4829', '5094', '5101', '5100', '4515', '4516',
            '4497', '4490', '4500', '1004490', '4496', '4512', '1005010', '5010', '5001', '1005001', '5024', '5023', '5026',
            '4587', '5031', '5030', '5033', '5035', '4641', '5038', '2322', '2229', '2294', '2304', '2301', '2298', '1002298',
            '2318', '2305', '2355', '2360', '2380', '2358', '2190', '2346', '2362', '2348', '2382', '2381', '2400', '2396',
            '2323', '2324', '1002324', '1016422', '16422', '2378', '2191', '16307', '16294', '16295', '16284', '16286', 
            '4256', '4280', '4314', '16300', '16414', '16326', '16325', '16331', '16333', '16339', '16337', '16304', '16338',
            '16345', '16343', '16409', '16371', '17777', '16365', '16373', '16355', '16383', '16377', '16387', '16382',
            '16398', '16394', '16391', '16396', '7644', '7638', '7636', '7647', '2236', '2217', '2086', '2118', '2219', 
            '1953', '1920', '1949', '2035', '1984', '2023', '2033', '2018', '1916', '2019', '2021', '2017', '2010', '2009',
            '2030', '2031', '2027', '2014', '2029', '2002', '2001', '2013', '2037', '2012', '1987', '1997', '3934', '3945',
            '16402', '4218', '3920', '3940', '3943', '3968', '4065', '3956', '3963', '3967', '16570', '4073', '4046', '4070',
            '4044', '3726', '3724', '3712', '6354', '3711', '3722', '3710', '3706', '3708', '3704', '3718', '3707', '17138',
            '17135', '6341', '6347', '6346', '6425', '6345', '6324', '6325']

OD_nearby = ['4064', '16571', '1004006', '4006', '4051', '1004051', '1003964', '3964', '1003962', '3962', '4067', '3959',
             '4066', '3951', '4071', '1003655', '3655', '3656', '16179', '16178', '6357', '6356', '16180', '3668', '4052',
             '1004052', '1003972', '3972', '1003909', '3909', '3908', '1003908', '3942', '3901', '3996', '1003996', '1003998',
             '3998', '1003995', '3995', '1003994', '3994', '1003991', '3991']

bridges = ['103655', '3655', '3929', '2026', '2072', '1002072', '2078', '1002078', '4216', '1004216', '16303', '2401', 
           '16810', '16423', '1002416', '2416', '2490', '1002490', '2914', '1002914', '2946', '1002946', '8126', '2999',
           '2022', '3669', '1002090', '2090', '2293', '3530', '2312', '1002312', '2564', '3425', '17084', '1017084', '5626',
           '2698', '2708', '1002708', '4922', '1004922', '4150', '7189']

link_str = '# ID Type LEN(mile) FFS(mile/h) Cap(v/hour) RHOJ(v/miles) Lane FFS_truck(mile/h) Cap_truck(v/hour) RHOJ_truck(v/miles) Convert_factor(1)\n'

_count_ctm = 0
_count_lq = 0
_count_pq = 0

veh_len_ratio = 0.7 # when jam, in average, (length occupied by 1 car) / (length occupied by 1 truck)
veh_headway_ratio = 0.8 # when max flux, in average, (length occupied by 1 car) / (length occupied by 1 truck)
convert_fc = 1.25
for data in link_records:
    _len = data[5]
    _lane = data[4]
    _cap_car = 2600
    _cap_truck = _cap_car * veh_headway_ratio
    _rhoJ_car = 280
    _rhoJ_truck = _rhoJ_car * veh_len_ratio
    _ffs_car = float(data[3])
    if (data[0] in ramp_IDs) or (data[0] in OD_nearby) or (data[0] in bridges):
        _ffs_car = 60
    if _ffs_car <= 15:
        _ffs_car = 30
        _cap_car = 2200
        _cap_truck = _cap_car * veh_headway_ratio
        _ffs_truck = _ffs_car * 0.9
    elif _ffs_car <= 25:
        _ffs_car = 40
        _cap_car = 2300
        _cap_truck = _cap_car * veh_headway_ratio
        _ffs_truck = _ffs_car * 0.85
    elif _ffs_car <= 35:
        _ffs_car = 50
        _cap_car = 2500
        _cap_truck = _cap_car * veh_headway_ratio
        _ffs_truck = _ffs_car * 0.8
    else:
        _ffs_car += 15
        _ffs_truck = _ffs_car * 0.8
    
    if data[6] == '1': # is_OD
        _type = 'PQ'
        _len = 1
        _lane = 1
        _cap_car = 99999
        _cap_truck = 99999
        _ffs_car = 99999
        _ffs_truck = 99999
        _rhoJ_car = 99999
        _rhoJ_truck = 99999
        _count_pq += 1
    elif (_ffs_car * 1600 / 3600 * 5 * 2) < (float(_len) * 1600):
        _type = 'CTM'
        _count_ctm += 1
    else:
        _type = 'LQ'
        _count_lq += 1
        _cap_car = 2600
        _cap_truck = _cap_car * veh_headway_ratio
        _rhoJ_car = 280
        _rhoJ_truck = _rhoJ_car * veh_len_ratio
        
    useful = [data[0], _type, _len, _ffs_car, _cap_car, _rhoJ_car, _lane, _ffs_truck, _cap_truck, _rhoJ_truck, convert_fc]
    newline = ' '.join([str(e) for e in useful]) + '\n'
    link_str += newline

# if not os.path.isdir(out_path):
#     os.makedirs(out_path)
# f = open(os.path.join(out_path, link_file_name), 'w')
# f.write(link_str)
# f.close()

print(_count_pq)
print(_count_ctm)
print(_count_lq)


1382
3788
10940

MNM_input_node


In [11]:
node_str = '# ID Type Convert_factor\n'
node_to_DMOND = {}
node_to_DMDND = {}
_cur_ind_DMOND = 100001
_cur_ind_DMDND = 150001
node_count = 0
# In-out node
for data in node_records:
    _id = data[0]
    if data[1] == 1: # is_OD
        node_to_DMOND[_id] = _cur_ind_DMOND
        node_to_DMDND[_id] = _cur_ind_DMDND
        _cur_ind_DMOND += 1
        _cur_ind_DMDND += 1
        node_count += 2
    else:
        useful = [_id, 'FWJ', convert_fc]
        newline = ' '.join([str(e) for e in useful]) + '\n'
        node_str += newline
        node_count += 1

# DMOND
for v in sorted(node_to_DMOND.items(), key=operator.itemgetter(1)):
    useful = [v[1], 'DMOND', convert_fc]
    newline = ' '.join([str(e) for e in useful]) + '\n'
    node_str += newline
# DMOND
for v in sorted(node_to_DMDND.items(), key=operator.itemgetter(1)):
    useful = [v[1], 'DMDND', convert_fc]
    newline = ' '.join([str(e) for e in useful]) + '\n'
    node_str += newline
    
# if not os.path.isdir(out_path):
#     os.makedirs(out_path)
# f = open(os.path.join(out_path, node_file_name), 'w')
# f.write(node_str)
# f.close()

node_count

for _id in [3577, 7425, 1963, 7393, 223]:
    print(node_to_DMOND[_id])
    print(node_to_DMDND[_id])


100123
150123
100251
150251
100063
150063
100247
150247
100010
150010

Snap_graph


In [23]:
graph_str = '# EdgeId FromNodeId ToNodeId\n'
for data in link_records:
    _id = data[0]
    _from = int(data[1])
    _to = int(data[2])
    if _from in node_to_DMOND:
        _from = node_to_DMOND[_from]
    if _to in node_to_DMDND:
        _to = node_to_DMDND[_to]
        
    useful = [_id, _from, _to]
    newline = ' '.join([str(e) for e in useful]) + '\n'
    graph_str += newline

if not os.path.isdir(out_path):
    os.makedirs(out_path)
f = open(os.path.join(out_path, graph_file_name), 'w')
f.write(graph_str)
f.close()

MNM_input_od


In [24]:
od_str = '# Origin_ID <-> node_ID\n'
all_origins = []
for v in sorted(node_to_DMOND.items(), key=operator.itemgetter(1)):
    od_str += ' '.join([str(e) for e in [v[1] + 100000, v[1]]]) + '\n'
    all_origins.append(v[1] + 100000)

all_destinations = []
od_str += '# Dest_ID <-> node_ID\n'
for v in sorted(node_to_DMDND.items(), key=operator.itemgetter(1)):
    od_str += ' '.join([str(e) for e in [v[1] + 100000, v[1]]]) + '\n'
    all_destinations.append(v[1] + 100000)

if not os.path.isdir(out_path):
    os.makedirs(out_path)
f = open(os.path.join(out_path, od_file_name), 'w')
f.write(od_str)
f.close()

MNM_input_demand


In [25]:
demand_str = '# Origin_ID Destination_ID <demand by interval> <truck demand by interval>\n'

# Randomly generate some demand for debug
# 20 intervals * 15 min/interval = 5 hours
np.random.seed(0)

for o in all_origins:
    for d in all_destinations:
        s_car = np.random.uniform(0, 0.1, 20)
        vari = np.random.uniform(-0.05, 0.05, 20)
        s_truck = s_car * (1 + vari) / 20
        demand_str += ' '.join([str(e) for e in ([o, d] + list(s_car) + list(s_truck))]) + '\n'

if not os.path.isdir(out_path):
    os.makedirs(out_path)
f = open(os.path.join(out_path, demand_file_name), 'w')
f.write(demand_str)
f.close()

path_table & path_table_buffer


In [26]:
# To be completed...
path_table_str = ''

# if not os.path.isdir(out_path):
#     os.makedirs(out_path)
# f = open(os.path.join(out_path, pathtable_file_name), 'w')
# f.write(path_table_str)
# f.close()


path_p_str = ''

# if not os.path.isdir(out_path):
#     os.makedirs(out_path)
# f = open(os.path.join(out_path, path_p_file_name), 'w')
# f.write(path_p_str)
# f.close()

config.conf


In [29]:
pathtable_loc = os.path.join('..', '..', 'data', 'input_files_MckeesRocks_SPC', 'path_table')
pathtable_file = open(pathtable_loc, 'r')
num_path = len(pathtable_file.readlines())

DTA_FIELDS = [('network_name', 'Snap_graph'), ('unit_time', 5), ('total_interval', -1), ('assign_frq', 180),
              ('start_assign_interval', 0), ('max_interval', 4), ('flow_scalar', 2), ('num_of_link', len(link_records)),
              ('num_of_node', node_count), ('num_of_O', len(all_origins)), ('num_of_D', len(all_destinations)),
              ('OD_pair', len(all_origins) * len(all_destinations)), ('routing_type', 'Fixed'), ('adaptive_ratio', 0.50)]

STAT_FIELDS = [('rec_mode', 'LRn'), ('rec_mode_para', 2), ('rec_folder', 'record'), ('rec_volume', 0), 
               ('volume_load_automatic_rec', 0), ('volume_record_automatic_rec', 0), ('rec_tt', 0), 
               ('tt_load_automatic_rec', 0), ('tt_record_automatic_rec', 0)]

ADAPTIVE_FIELDS = [('route_frq', 180)]

FIXED_FIELDS = [('path_file_name', 'path_table'), ('num_path', num_path), ('choice_portion', 'Buffer'), ('route_frq', 180)]

config_str = '[DTA]\n'
for v in DTA_FIELDS:
    config_str += "{} = {}\n".format(str(v[0]), str(v[1]))

config_str += '\n[STAT]\n'
for v in STAT_FIELDS:
    config_str += "{} = {}\n".format(str(v[0]), str(v[1]))

config_str += '\n[ADAPTIVE]\n'
for v in ADAPTIVE_FIELDS:
    config_str += "{} = {}\n".format(str(v[0]), str(v[1]))
    
config_str += '\n[FIXED]\n'
for v in FIXED_FIELDS:
    config_str += "{} = {}\n".format(str(v[0]), str(v[1]))
    
if not os.path.isdir(out_path):
    os.makedirs(out_path)
f = open(os.path.join(out_path, config_file_name), 'w')
f.write(config_str)
f.close()

In [ ]: