|
20 | 20 |
|
21 | 21 | ''' |
22 | 22 | import os |
23 | | -import re |
24 | 23 | import operator |
25 | 24 | import arcpy |
26 | | -import pandas as pd |
27 | 25 | from MHN import MasterHighwayNetwork # Custom class for MHN processing functionality |
28 | 26 |
|
29 | 27 | # ----------------------------------------------------------------------------- |
@@ -1052,247 +1050,6 @@ def get_scen_line_ids(): |
1052 | 1050 | # Node extra attribute CSVs |
1053 | 1051 | exec(open(os.path.join(MHN.src_dir, 'transit_node_extra_attributes.py')).read()) |
1054 | 1052 |
|
1055 | | - |
1056 | | -# -------------------------------------------------------------- |
1057 | | -# generate @busveq and write into highway batchin files |
1058 | | -# -------------------------------------------------------------- |
1059 | | -# write @busveq to highway batchin files |
1060 | | - |
1061 | | -arcpy.AddMessage('-- CREATING @BUSVEQ --') |
1062 | | - |
1063 | | -''' |
1064 | | -This section creates @busveq to place into the highway extra attribute batchin files. |
1065 | | -This is done by taking headway info from the bus batchin files. |
1066 | | -First, the script reads the highway and transit times of day, determines what fraction of the |
1067 | | -day each time of day occupies, then derives the fraction that represents how much of each |
1068 | | -transit TOD that a highway TOD occupies. (For example, transit TOD 1 is 6pm-6am, and |
1069 | | -highway TOD 1 is 8pm-6am, so highway TOD 1 takes up 80% of transit TOD 1, and 0.8 is stored.) |
1070 | | -
|
1071 | | -Bus volumes on links are estimated using bus itineraries and headway information. Volumes are |
1072 | | -derived for each transit TOD, then these volumes are multiplied by the fraction described above to derive the |
1073 | | -approximate volume for highway TOD. |
1074 | | -
|
1075 | | -Highway TOD bus volumes are then converted to bus vehicle equivalents (by multiplying *3). |
1076 | | -The highway extra attribute batchin files (with suffix '.l2') are then read, reconstructed, |
1077 | | -and overwritten to include @busveq as an extra attribute. |
1078 | | -''' |
1079 | | - |
1080 | | -# generate_hour_list() converts MHN.tod_periods to lists of hours instead of description+sql query |
1081 | | -# e.g.: |
1082 | | -# before: {'highway' : {'1': ('8PM-6AM', '"STARTHOUR">=20 OR "STARTHOUR"<=5'), ...}} |
1083 | | -# after: {'highway' : {'1': [20, 21, 22, 23, 0, 1, 2, 3, 4, 5], ...}} |
1084 | | - |
1085 | | -def generate_hour_list(time_range): |
1086 | | - start, end = time_range.split('-') |
1087 | | - start_hour = int(start[:-2]) |
1088 | | - end_hour = int(end[:-2]) |
1089 | | - if start.endswith('PM') and start_hour != 12: |
1090 | | - start_hour += 12 |
1091 | | - if end.endswith('PM') and end_hour != 12: |
1092 | | - end_hour += 12 |
1093 | | - if start_hour > end_hour: |
1094 | | - hours = list(range(start_hour, 24))+list(range(0, end_hour)) |
1095 | | - else: |
1096 | | - hours = list(range(start_hour, end_hour)) |
1097 | | - if end_hour == 24: |
1098 | | - hours.append(0) |
1099 | | - return hours |
1100 | | - |
1101 | | -# Generate hour lists for each time period by calling function above |
1102 | | -hour_lists = {} |
1103 | | -for mode, periods in MHN.tod_periods.items(): |
1104 | | - hour_lists[mode] = {} |
1105 | | - for period, time_range in periods.items(): |
1106 | | - hour_lists[mode][period] = generate_hour_list(time_range[0]) |
1107 | | - |
1108 | | -arcpy.AddMessage('hour_lists: \n{}'.format(hour_lists)) |
1109 | | - |
1110 | | -#using hour_lists generated above, calculate % of each hwy tod inside a given trnt tod |
1111 | | -#output is of the format: {hwy_tod: [trnt_tod: fraction]} |
1112 | | -#where hwy_tod is number 1-8, trnt_tod is number 1-4, and |
1113 | | -#fraction is (#hours in hwy_tod / #hours in trnt_tod) |
1114 | | - |
1115 | | -hwy_trnt_ratio = {} |
1116 | | -for hwy_period, hwy_hours in hour_lists['highway'].items(): |
1117 | | - trnt_periods = [] |
1118 | | - trnt_hours = [] |
1119 | | - for trnt_period_num, trnt_pd_hours in hour_lists['transit'].items(): |
1120 | | - if set(hwy_hours).intersection(trnt_pd_hours): |
1121 | | - trnt_periods.append(trnt_period_num) |
1122 | | - trnt_hours.extend(trnt_pd_hours) |
1123 | | - fraction = float(len(hwy_hours)) / float(len(trnt_hours)) |
1124 | | - hwy_trnt_ratio[hwy_period] = [trnt_periods[0], fraction] |
1125 | | - |
1126 | | -arcpy.AddMessage('\nhwy_trnt_ratio: \n{}'.format(hwy_trnt_ratio)) |
1127 | | - |
1128 | | -## -- |
1129 | | -## for each scenario year: (1) derive bus volumes by transit time-of-day, |
1130 | | -## then (2) convert volumes to highway time-of-day |
1131 | | -## -- |
1132 | | - |
1133 | | -for scen in scen_list: |
1134 | | - arcpy.AddMessage('scenario {}...'.format(scen)) |
1135 | | - #cycle through each transit time of day and derive volumes on roadways |
1136 | | - scen_dfs = [] |
1137 | | - |
1138 | | - ## GATHER AND COMBINE TRANSIT ITINERARY FILES FOR {TOD} IN {SCEN}-- |
1139 | | - for tod in MHN.tod_periods['transit'].keys(): |
1140 | | - |
1141 | | - itinerary = os.path.join(tran_path, scen, 'bus.itinerary_{}'.format(tod)) |
1142 | | - |
1143 | | - #read itinerary file for tod_period |
1144 | | - with open(itinerary, 'r') as file: |
1145 | | - lines = file.readlines() |
1146 | | - |
1147 | | - #grab only route info (starts with 'a' in itinerary) |
1148 | | - rteinfo = [] |
1149 | | - for line in lines: |
1150 | | - if line.startswith('a'): |
1151 | | - #regular expression splits the info into columns using double-space as delimiter (unless it's within quotation marks) |
1152 | | - rte = re.split(''' (?=(?:[^'"]|'[^']*'|"[^"]*")*$)''', line) |
1153 | | - |
1154 | | - #cleanup |
1155 | | - rte = [item.replace('\n','') for item in rte] #remove newline character |
1156 | | - rte = [item.replace('"','') for item in rte] #remove quotation marks |
1157 | | - rte = [item.replace("'","") for item in rte] |
1158 | | - rte.append(lines.index(line)) #stores the row in the itinerary that the route info is located |
1159 | | - rteinfo.append(rte) |
1160 | | - |
1161 | | - rteinfo = [rte[1:] for rte in rteinfo] #remove the 'a' column |
1162 | | - rteinfo.insert(0,['ROUTE_ID','MODE','VEHTYPE','HDWY','SPEED','DESC','INDEX']) #create a row of column labels |
1163 | | - |
1164 | | - #next, grab the nodes of each route itinerary |
1165 | | - indices = [[line[-1], line[0]] for line in rteinfo[1:]] #grab the row in the itinerary file that each route starts on |
1166 | | - itin = [] #will include list of all the links along with their route info |
1167 | | - #grab the nodes that each itinerary traverses |
1168 | | - for rte in indices: |
1169 | | - rte_name = rte[1] |
1170 | | - rte_index_num = rte[0] |
1171 | | - itin_start = rte_index_num+3 #the line in `indices` that corresponds to the first line of a route's itinerary |
1172 | | - |
1173 | | - #find the last line in the route's itinerary |
1174 | | - #if rte is not last in the list, the final row is the row just before the next rte |
1175 | | - if rte != indices[-1]: |
1176 | | - itin_end = indices[indices.index(rte)+1][0]-1 |
1177 | | - rte_itin = lines[itin_start: itin_end] |
1178 | | - #if rte is last in the list, read to the end |
1179 | | - else: |
1180 | | - rte_itin = lines[itin_start:] |
1181 | | - |
1182 | | - #clean up itin files; we just want the nodes per rte |
1183 | | - rte_itin = [rte.lstrip() for rte in rte_itin] |
1184 | | - rte_itin = [rte.split(' ') for rte in rte_itin] |
1185 | | - itin_nodes = [rte[0] for rte in rte_itin] |
1186 | | - |
1187 | | - #make i-j pairs (named 'ANODE' & 'BNODE' here) from the list of nodes |
1188 | | - rte_itin_rebuild = [] |
1189 | | - for node in itin_nodes[:-1]: |
1190 | | - nextnode = itin_nodes[itin_nodes.index(node)+1] |
1191 | | - link = [rte_name, node, itin_nodes[itin_nodes.index(node)+1]] # link = [ROUTE_ID, ANODE, BNODE] |
1192 | | - rte_itin_rebuild.append(link) |
1193 | | - |
1194 | | - |
1195 | | - itin.append(rte_itin_rebuild) |
1196 | | - |
1197 | | - #prep for pandas |
1198 | | - itin = [link for route in itin for link in route] #flatten to list of [ROUTE_ID, ANODE, BNODE] |
1199 | | - itin.insert(0, ['ROUTE_ID', 'ANODE', 'BNODE']) #add column labels |
1200 | | - |
1201 | | - |
1202 | | - itin_df = pd.DataFrame(data=itin[1:], columns=itin[0]) |
1203 | | - rteinfo_df = pd.DataFrame(data=rteinfo[1:], columns=rteinfo[0]) |
1204 | | - |
1205 | | - |
1206 | | - #merge into one and cleanup |
1207 | | - final = pd.merge(itin_df, rteinfo_df, how='left', on='ROUTE_ID') |
1208 | | - |
1209 | | - final.drop(labels=['VEHTYPE', 'SPEED', 'DESC','INDEX'], axis=1, inplace=True) |
1210 | | - final['HDWY'] = final['HDWY'].astype(float) |
1211 | | - for column in ['ANODE','BNODE']: |
1212 | | - final[column] = final[column].astype(int) |
1213 | | - |
1214 | | - |
1215 | | - #actual volumes analysis |
1216 | | - final['trnt_tod'] = str(tod) #transit time period |
1217 | | - final['hrs'] = len(hour_lists['transit'][str(tod)]) #hours in time period |
1218 | | - |
1219 | | - # arcpy.AddMessage('time period hours: {}'.format(len(hour_lists['transit'][str(tod)]))) |
1220 | | - |
1221 | | - final['vol'] = final['hrs'] * (60 / final['HDWY']) |
1222 | | - |
1223 | | - scen_dfs.append(final) |
1224 | | - |
1225 | | - scen_vol = pd.concat(scen_dfs) |
1226 | | - |
1227 | | - #convert transit tods to highway tods |
1228 | | - hwy_vols = [] |
1229 | | - for hwy_tod in hwy_trnt_ratio: |
1230 | | - tod_hwy_vol = scen_vol.loc[scen_vol['trnt_tod']==hwy_trnt_ratio[hwy_tod][0]].copy() |
1231 | | - tod_hwy_vol['hwyvol'] = tod_hwy_vol['vol'] * hwy_trnt_ratio[hwy_tod][1] |
1232 | | - tod_hwy_vol['hwy_tod'] = hwy_tod |
1233 | | - hwy_vols.append(tod_hwy_vol) |
1234 | | - # arcpy.AddMessage('hwy_tod={} \ntrnt_tod={}'.format(hwy_tod, hwy_trnt_ratio[hwy_tod][0])) |
1235 | | - # arcpy.AddMessage("hwyvol=tod_hwy_vol['vol']*{}".format(hwy_trnt_ratio[hwy_tod][1])) |
1236 | | - |
1237 | | - |
1238 | | - scen_hwy_vols = pd.concat(hwy_vols, ignore_index=True) |
1239 | | - |
1240 | | - arcpy.AddMessage('scen {} total volumes before/after conversion to hwy tods:'.format(scen)) |
1241 | | - arcpy.AddMessage('before -- {}'.format(scen_vol['vol'].sum())) |
1242 | | - arcpy.AddMessage('after -- {}'.format(scen_hwy_vols['hwyvol'].sum())) |
1243 | | - if round(scen_vol['vol'].sum(),1) != round(scen_hwy_vols['hwyvol'].sum(),1): |
1244 | | - raise ValueError('Problem converting transit bus volumes to highway!') |
1245 | | - |
1246 | | - scen_hwy_vols_bylink = scen_hwy_vols.groupby(['ANODE','BNODE','hwy_tod']).agg({'hwyvol':'sum'}).reset_index() |
1247 | | - scen_hwy_vols_bylink['@busveq'] = scen_hwy_vols_bylink['hwyvol'] * 3 |
1248 | | - scen_hwy_vols_bylink = scen_hwy_vols_bylink[['ANODE','BNODE','@busveq','hwy_tod']].copy() |
1249 | | - |
1250 | | - #create total daily and append to scen_hwy_vols_bylink |
1251 | | - tot_daily_vol = scen_hwy_vols_bylink.groupby(['ANODE','BNODE']).agg({'@busveq':'sum'}).reset_index() |
1252 | | - tot_daily_vol['hwy_tod'] = '0' |
1253 | | - scen_hwy_vols_bylink = pd.concat([scen_hwy_vols_bylink, tot_daily_vol]) |
1254 | | - |
1255 | | - #add total daily to the iterator used below |
1256 | | - tod_iter = sorted(MHN.tod_periods['highway'].keys()) |
1257 | | - tod_iter.append('0') |
1258 | | - |
1259 | | - # BEGIN TOD ITERATOR (for writing to hwy .l2 files) -- |
1260 | | - arcpy.AddMessage('writing busveq to .l2 files...') |
1261 | | - for tod in tod_iter: |
1262 | | - arcpy.AddMessage('tod {}...'.format(tod)) |
1263 | | - #extra attributes highway batchin file |
1264 | | - hwy_l2 = os.path.join(hwy_path, scen, '{}0{}.l2'.format(scen, tod)) |
1265 | | - |
1266 | | - busveqs = scen_hwy_vols_bylink #renamed for convenience |
1267 | | - busveqs = busveqs.loc[busveqs['hwy_tod']==tod, ['ANODE','BNODE','@busveq']].copy() |
1268 | | - busveqs['ijpair'] = busveqs['ANODE'].astype('str') + '-' + busveqs['BNODE'].astype('str') |
1269 | | - |
1270 | | - with open(hwy_l2, 'r') as file: |
1271 | | - lines = file.readlines() |
1272 | | - if lines[0].endswith('@busveq'): |
1273 | | - raise ValueError('.l2 file already contains @busveq! Re-run "Generate Highway Files."') |
1274 | | - lines[0] = lines[0].rstrip() + ',@busveq\n' #add @busveq column label |
1275 | | - |
1276 | | - #add busveq values to each ij pair |
1277 | | - for i in range(1, len(lines)): |
1278 | | - ijnode = lines[i].strip().split()[:2] #first two columns are inode and jnode |
1279 | | - ijnode = str(ijnode[0]) + '-' + str(ijnode[1]) #make id to match with busveqs |
1280 | | - if ijnode in busveqs['ijpair'].values: |
1281 | | - num = busveqs.loc[busveqs['ijpair']==ijnode, '@busveq'].iloc[0] |
1282 | | - else: |
1283 | | - num = 0 |
1284 | | - lines[i] = lines[i].rstrip() + ' {}\n'.format(num) |
1285 | | - |
1286 | | - #overwrite hwy extra attributes batchin file |
1287 | | - with open(hwy_l2, 'w') as file: #overwrite |
1288 | | - file.writelines(lines) |
1289 | | - #-- END TOD ITERATOR |
1290 | | - |
1291 | | -#-- END SCEN ITERATOR |
1292 | | - |
1293 | | - |
1294 | | -arcpy.AddMessage('Done with @busveq!') |
1295 | | - |
1296 | 1053 | # ----------------------------------------------------------------------------- |
1297 | 1054 | # Clean up script-level data. |
1298 | 1055 | # ----------------------------------------------------------------------------- |
|
0 commit comments