Marey / string-line diagram for Amtrak Lincoln Service
Rail — CHI-STL
Example from the compendium of canonical charts
Python Code
"""Rail — Marey / string-line diagram for Amtrak Lincoln Service (CHI-STL)."""
import numpy as np
import pandas as pd
import plotly.graph_objects as go
import io, zipfile, requests, warnings, math
GTFS_URL = "https://content.amtrak.com/content/gtfs/GTFS.zip"
STATION_NAMES = ["Chicago Union Station", "Joliet", "Pontiac",
"Bloomington-Normal", "Springfield", "Carlinville", "Alton",
"St. Louis"]
MILE_POSTS = [0, 37, 97, 136, 188, 228, 270, 284]
DWELL_MIN = 3 # minutes at each intermediate stop
def haversine(lat1, lon1, lat2, lon2):
R = 3958.8
phi1, phi2 = math.radians(lat1), math.radians(lat2)
dphi = math.radians(lat2 - lat1)
dlam = math.radians(lon2 - lon1)
a = math.sin(dphi/2)**2 + math.cos(phi1)*math.cos(phi2)*math.sin(dlam/2)**2
return 2 * R * math.asin(math.sqrt(a))
def parse_time(s):
try:
parts = str(s).strip().split(":")
return int(parts[0]) * 60 + int(parts[1]) + int(parts[2]) / 60
except Exception:
return None
def train_times(dep_min, mile_posts, speed_mph, dwell_min):
"""Return list of (arrival_min, departure_min) pairs per station.
Works for both southbound (ascending mileposts) and northbound (descending)."""
origin_mp = mile_posts[0]
dwell_elapsed = 0
pairs = []
for j, mp in enumerate(mile_posts):
travel = abs(mp - origin_mp) / speed_mph * 60
arrival = dep_min + travel + dwell_elapsed
if j > 0 and j < len(mile_posts) - 1:
departure = arrival + dwell_min
dwell_elapsed += dwell_min
else:
departure = arrival
pairs.append((arrival, departure))
return pairs
def add_train(fig, dep_min, mile_posts, speed_mph, dwell_min, color, name):
"""Add a single train trace with visible station dwell segments."""
pairs = train_times(dep_min, mile_posts, speed_mph, dwell_min)
x_pts, y_pts = [], []
for j, ((arr, dep), mp) in enumerate(zip(pairs, mile_posts)):
x_pts.append(arr / 60)
y_pts.append(mp)
if arr != dep: # show dwell as horizontal segment
x_pts.append(dep / 60)
y_pts.append(mp)
fig.add_trace(go.Scatter(
x=x_pts, y=y_pts,
mode="lines+markers",
line=dict(color=color, width=4),
marker=dict(size=9, color=color),
name=name,
hovertemplate="<b>%{fullData.name}</b><br>Time: %{x:.2f}h<br>Mile: %{y:.0f}<extra></extra>",
))
def build_synthetic_fig():
speed_mph = 60
# Southbound departures from Chicago (hours from midnight)
sb_deps = [7*60, 11*60 + 30, 16*60]
# Northbound departures from St. Louis
nb_deps = [8*60, 13*60, 17*60 + 30]
sb_colors = [VIOLET, GREEN, TEAL]
nb_colors = [TEAL, PINK, MUTED]
fig = go.Figure()
trains = sorted(
[(dep, MILE_POSTS, sb_colors[i], f"SB-{i+1}") for i, dep in enumerate(sb_deps)] +
[(dep, list(reversed(MILE_POSTS)), nb_colors[i], f"NB-{i+1}") for i, dep in enumerate(nb_deps)],
key=lambda t: t[0],
)
for dep, miles, color, name in trains:
add_train(fig, dep, miles, speed_mph, DWELL_MIN, color, name)
for mp in MILE_POSTS:
fig.add_hline(y=mp, line=dict(color="#e0e0e0", width=0.8))
fig.update_layout(
xaxis=dict(
title="Time of day",
range=[5.8, 23.5],
tickvals=list(range(6, 24, 2)),
ticktext=[f"{h}:00" for h in range(6, 24, 2)],
showgrid=False,
),
yaxis=dict(
title="Distance from Chicago (miles)",
range=[-8, 295],
showgrid=False,
showticklabels=True,
tickvals=MILE_POSTS,
ticktext=STATION_NAMES,
tickfont=dict(size=12, color="#1a1a1a"),
),
legend=dict(orientation="h", x=0, y=1.0, xanchor="left", yanchor="bottom", font=dict(size=11)),
margin=dict(t=50, b=50, l=160, r=40),
)
return fig
def try_gtfs():
"""Attempt to build figure from Amtrak GTFS."""
warnings.filterwarnings("ignore")
requests.packages.urllib3.disable_warnings()
r = requests.get(GTFS_URL, verify=False, timeout=60)
r.raise_for_status()
zf = zipfile.ZipFile(io.BytesIO(r.content))
routes_df = pd.read_csv(io.BytesIO(zf.read("routes.txt")))
trips_df = pd.read_csv(io.BytesIO(zf.read("trips.txt")))
stop_times = pd.read_csv(io.BytesIO(zf.read("stop_times.txt")))
stops = pd.read_csv(io.BytesIO(zf.read("stops.txt")))
routes_df["_name_str"] = routes_df["route_long_name"].fillna("").astype(str)
lincoln_routes = routes_df[routes_df["_name_str"].str.contains("Lincoln", case=False)]
if lincoln_routes.empty:
raise ValueError("No Lincoln routes found")
route_ids = lincoln_routes["route_id"].tolist()
lincoln_trips = trips_df[trips_df["route_id"].isin(route_ids)]
trip_ids = lincoln_trips["trip_id"].tolist()
lincoln_st = stop_times[stop_times["trip_id"].isin(trip_ids)].copy()
lincoln_st = lincoln_st.merge(
stops[["stop_id", "stop_lat", "stop_lon", "stop_name"]], on="stop_id", how="left")
lincoln_st["dep_min"] = lincoln_st["departure_time"].apply(parse_time)
lincoln_st = lincoln_st.dropna(subset=["dep_min"]).sort_values(
["trip_id", "stop_sequence"])
if lincoln_st.empty:
raise ValueError("No valid stop times")
# Build cumulative distance from reference southbound trip
ref_trip = lincoln_trips[lincoln_trips["direction_id"] == 0]["trip_id"].iloc[0]
ref_stops = lincoln_st[lincoln_st["trip_id"] == ref_trip].sort_values("stop_sequence")
lats, lons = ref_stops["stop_lat"].values, ref_stops["stop_lon"].values
cum_dist = [0.0]
for i in range(1, len(lats)):
if pd.notna(lats[i]) and pd.notna(lats[i-1]):
cum_dist.append(cum_dist[-1] + haversine(lats[i-1], lons[i-1], lats[i], lons[i]))
else:
cum_dist.append(cum_dist[-1])
seq_to_dist = dict(zip(ref_stops["stop_sequence"].values, cum_dist))
total_dist = cum_dist[-1]
station_labels = list(zip(ref_stops["stop_name"].tolist(), cum_dist))
color_cycle = [VIOLET, GREEN, TEAL, TEAL, PINK, MUTED, "#e63030", "#a0a0a0"]
fig = go.Figure()
trip_ids = sorted(trip_ids, key=lambda tid: lincoln_st[lincoln_st["trip_id"] == tid]["dep_min"].min())
for i, tid in enumerate(trip_ids):
trip = lincoln_st[lincoln_st["trip_id"] == tid].sort_values("stop_sequence")
seqs = trip["stop_sequence"].values
dists_mapped = [seq_to_dist.get(s) for s in seqs]
trip = trip.copy()
trip["dist_mi"] = dists_mapped
trip = trip.dropna(subset=["dist_mi"])
if len(trip) < 2:
continue
dir_id = lincoln_trips[lincoln_trips["trip_id"] == tid]["direction_id"].values
if len(dir_id) and dir_id[0] == 1:
trip["dist_mi"] = total_dist - trip["dist_mi"]
direction = "NB"
else:
direction = "SB"
fig.add_trace(go.Scatter(
x=trip["dep_min"] / 60,
y=trip["dist_mi"],
mode="lines+markers",
line=dict(color=color_cycle[i % len(color_cycle)], width=4),
marker=dict(size=9),
name=f"{direction}-{tid}",
hovertemplate="<b>%{fullData.name}</b><br>Time: %{x:.2f}h<br>Mile: %{y:.1f}<extra></extra>",
))
for _, dist in station_labels:
fig.add_hline(y=dist, line=dict(color="#e0e0e0", width=0.8))
fig.update_layout(
xaxis=dict(
title="Time of day",
tickvals=list(range(4, 25, 2)),
ticktext=[f"{h}:00" for h in range(4, 25, 2)],
range=[3.8, 24],
showgrid=False,
),
yaxis=dict(
title="Distance from Chicago (miles)",
range=[-5, total_dist + 10],
showgrid=False,
showticklabels=True,
tickvals=[dist for _, dist in station_labels],
ticktext=[name for name, _ in station_labels],
tickfont=dict(size=12, color="#1a1a1a"),
),
legend=dict(orientation="h", x=0, y=1.0, xanchor="left", yanchor="bottom", font=dict(size=11)),
margin=dict(t=50, b=50, l=160, r=40),
)
return fig
def generate():
print("building Amtrak Lincoln Service Marey diagram …")
try:
fig = try_gtfs()
print(" using real GTFS data")
except Exception as e:
print(f" GTFS failed ({e}), using synthetic timetable")
fig = build_synthetic_fig()
return fig
if __name__ == "__main__":
generate()Made with Plotly