"""Replicate every number behind bartanalysis.com/stationmap from BART's raw data. Standard library only. Downloads BART's public hourly origin-destination file and writes three CSVs into the current directory: table1_raw_trips.csv station_in, station_out, date, hour, people (every row BART publishes: exact faregate counts, one row per station-pair per exit-hour per date) table2_hourly_od.csv year, daytype, hour, station_in, station_out, avg_people_per_day (averaged across the year's weekdays / weekend+holiday days) table3_link_flows.csv year, daytype, hour, link_from, link_to, avg_people_per_day (every trip routed along its shortest path; riders added to each track link crossed) Usage: python3 replicate.py [year] (default 2025) Data source: https://afcweb.bart.gov/ridership/origin-destination/ (one gzipped CSV per year, no login, CC-BY. "Hour" is the EXIT hour.) Holidays: to avoid dependencies this script treats only fixed-date federal holidays plus scanned floating ones via a small built-in table for 2019/2025; results match the site to within rounding. The full pipeline (with the python-holidays package, input hashing, and cross-checks) is the analysis notebook this site links to. """ import collections import csv import datetime import gzip import sys import urllib.request YEAR = sys.argv[1] if len(sys.argv) > 1 else "2025" URL = f"https://afcweb.bart.gov/ridership/origin-destination/date-hour-soo-dest-{YEAR}.csv.gz" # US federal holidays for the two years the site covers (extend for other years) FEDERAL_HOLIDAYS = { "2019": ["2019-01-01", "2019-01-21", "2019-02-18", "2019-05-27", "2019-07-04", "2019-09-02", "2019-10-14", "2019-11-11", "2019-11-28", "2019-12-25"], "2025": ["2025-01-01", "2025-01-20", "2025-02-17", "2025-05-26", "2025-06-19", "2025-07-04", "2025-09-01", "2025-10-13", "2025-11-11", "2025-11-27", "2025-12-25"], } holidays = set(FEDERAL_HOLIDAYS.get(YEAR, [])) # the BART track network: every pair of directly-connected stations (2025 system) SEGMENTS = [ ("EMBR","MONT"),("MONT","POWL"),("POWL","CIVC"),("CIVC","16TH"),("16TH","24TH"), ("24TH","GLEN"),("GLEN","BALB"),("BALB","DALY"),("DALY","COLM"),("COLM","SSAN"), ("SSAN","SBRN"),("SBRN","SFIA"),("SBRN","MLBR"),("SFIA","MLBR"),("WOAK","EMBR"), ("WOAK","12TH"),("WOAK","LAKE"),("12TH","LAKE"),("12TH","19TH"),("19TH","MCAR"), ("MCAR","ASHB"),("ASHB","DBRK"),("DBRK","NBRK"),("NBRK","PLZA"),("PLZA","DELN"), ("DELN","RICH"),("MCAR","ROCK"),("ROCK","ORIN"),("ORIN","LAFY"),("LAFY","WCRK"), ("WCRK","PHIL"),("PHIL","CONC"),("CONC","NCON"),("NCON","PITT"),("PITT","PCTR"), ("PCTR","ANTC"),("LAKE","FTVL"),("FTVL","COLS"),("COLS","SANL"),("SANL","BAYF"), ("COLS","OAKL"),("BAYF","CAST"),("CAST","WDUB"),("WDUB","DUBL"),("BAYF","HAYW"), ("HAYW","SHAY"),("SHAY","UCTY"),("UCTY","FRMT"),("FRMT","WARM"),("WARM","MLPT"), ("MLPT","BERY"), ] neighbors = collections.defaultdict(list) for a, b in SEGMENTS: neighbors[a].append(b) neighbors[b].append(a) def shortest_path(origin, destination): """Fewest-stations path (breadth-first search).""" came_from = {origin: origin} frontier = [origin] while frontier: nxt = [] for station in frontier: for nb in neighbors[station]: if nb not in came_from: came_from[nb] = station nxt.append(nb) frontier = nxt path = [destination] while path[-1] != origin: path.append(came_from[path[-1]]) return list(reversed(path)) def daytype_of(date_text): day = datetime.date.fromisoformat(date_text) return "weekend+holiday" if (day.weekday() >= 5 or date_text in holidays) else "weekday" print("downloading", URL) local = f"date-hour-soo-dest-{YEAR}.csv.gz" urllib.request.urlretrieve(URL, local) # ---- table 1: the raw rows, plus year totals for tables 2 and 3 ---- year_totals = collections.Counter() # (daytype, hour, in, out) -> people, whole year days_in_bucket = collections.defaultdict(set) with gzip.open(local, "rt") as raw, open("table1_raw_trips.csv", "w", newline="") as out1: writer = csv.writer(out1) writer.writerow(["station_in", "station_out", "date", "hour", "people"]) for date_text, hour_text, origin, destination, people in csv.reader(raw): writer.writerow([origin, destination, date_text, hour_text, people]) daytype = daytype_of(date_text) days_in_bucket[daytype].add(date_text) if origin != destination: year_totals[(daytype, int(hour_text), origin, destination)] += int(people) print("wrote table1_raw_trips.csv") # ---- table 2: average people per day for each (station in, station out, hour) ---- with open("table2_hourly_od.csv", "w", newline="") as out2: writer = csv.writer(out2) writer.writerow(["year", "daytype", "hour", "station_in", "station_out", "avg_people_per_day"]) for (daytype, hour, origin, destination), total in sorted(year_totals.items()): writer.writerow([YEAR, daytype, hour, origin, destination, round(total / len(days_in_bucket[daytype]), 2)]) print("wrote table2_hourly_od.csv") # ---- table 3: route every pair along its path; sum people per track link ---- link_flow = collections.Counter() # (daytype, hour, from, to) -> people, whole year path_cache = {} for (daytype, hour, origin, destination), total in year_totals.items(): if origin not in neighbors or destination not in neighbors: continue # station outside the 2025 network map if (origin, destination) not in path_cache: path_cache[(origin, destination)] = shortest_path(origin, destination) path = path_cache[(origin, destination)] for a, b in zip(path, path[1:]): link_flow[(daytype, hour, a, b)] += total with open("table3_link_flows.csv", "w", newline="") as out3: writer = csv.writer(out3) writer.writerow(["year", "daytype", "hour", "link_from", "link_to", "avg_people_per_day"]) for (daytype, hour, a, b), total in sorted(link_flow.items()): writer.writerow([YEAR, daytype, hour, a, b, round(total / len(days_in_bucket[daytype]), 2)]) print("wrote table3_link_flows.csv") print("days:", {k: len(v) for k, v in days_in_bucket.items()})