Map Challenge – Day 2 (lines): Straightening the A2

Day 2 of the 30daymapchallenge focused on cleaning and preparing the A2 motorway data before creating a straightened schematic line. Raw geospatial data is rarely perfect, and getting the order, connectivity, and attributes right is essential for meaningful visualization. So I decided to pick up with the data from day 1 focusing on toll data in Germany.

Getting the data

Data was provided by Toll Collect on the Verkehrsportal.

open data site on the “Verkehrsportal”

First I converted the csv to a geojson and filtered all non A2 segments from the set:

import pandas as pd
import geojson
from shapely import wkt

# --- Settings ---
input_file = "Lkw_Befahrungen_1m/Daten/abschnitte_2025-10-02_prc_2025-10-11_exp_2025-10-10.csv"
output_file = "a2_lines.geojson"

# --- Read the semicolon-separated CSV ---
df = pd.read_csv(input_file, sep=';')

# --- Filter only rows where bundesfernstrasse == 'A2' ---
df = df[df['bundesfernstrasse'] == 'A2']

# --- Create GeoJSON features ---
features = []
for _, row in df.iterrows():
    try:
        geometry = wkt.loads(row['wkt'])
        feature = geojson.Feature(
            geometry=geojson.loads(geojson.dumps(geometry.__geo_interface__)),
            properties={
                "datum": row["datum"],
                "abschnitt_id": row["abschnitt_id"],
                "mautknoten_von": row["mautknoten_name_von"],
                "mautknoten_nach": row["mautknoten_name_nach"],
                "bundesfernstrasse": row["bundesfernstrasse"],
                "anzahl_befahrungen": row["anzahl_befahrungen"],
                "anzahl_einfahrten": row["anzahl_einfahrten"],
                "anzahl_ausfahrten": row["anzahl_ausfahrten"],
                "laenge_km": row["laenge_km"],
                "fahrleistung_km": row["fahrleistung_km"]
            }
        )
        features.append(feature)
    except Exception as e:
        print(f"Skipping row {row['abschnitt_id']}: {e}")

# --- Create GeoJSON FeatureCollection ---
geojson_data = geojson.FeatureCollection(features)

# --- Write to file ---
with open(output_file, "w", encoding="utf-8") as f:
    geojson.dump(geojson_data, f, ensure_ascii=False, indent=2)

print(f"Exported {len(features)} features to {output_file}")

Cleaning the data in QGIS

Before plotting, I needed to address several issues:

  • Join attributes by location: To create a single line feature for both directions I used this in QGIS
  • Remove duplicate geometries: The result from above still has both directions inside the data, so I removed duplicate geometries
  • Remove individual ramps mislabelled as A2: Some small connecting ramps were still labeled as “A2”; I filtered these out to avoid misleading spikes in the schematic.
  • Correct geometry directions: A few segments were digitized backwards by the data provider or my processes above; I reversed their direction ensures consistent east → west ordering.
  • Check non-overlapping segments and manually sum values: In some cases, a single motorway segment was split into two in opposite directions. I merged them manually, and traffic counts summed to maintain accuracy.

This preprocessing work is crucial — “garbage in, garbage out” is especially true when visualizing linear transport data.

From messy geometry to a straight line

With the cleaned dataset, I created a straightened schematic:

  • Segments are laid out in sequence from east to west.
  • Segment width corresponds to length; color indicates traffic intensity (anzahl_befahrungen_sum).
  • Labels are shown for the first node, the last node, and every 4th segment, providing context without clutter. This also indicate “density” as more densely populated areas result in denser segments

This approach is inspired by metro maps: clear connectivity while still showing relative magnitude.

Adding geographic context

To connect the schematic back to reality, I added an overview map below the straightened line:

  • Shows actual A2 geometries on a clean basemap.
  • Start and end nodes are labeled, so viewers can relate schematic and geographic views.
  • The map is aligned in width with the schematic and automatically scaled to show all segments.

The map creation script

We all use it, so did I: I’ve used and co-worked with chatGPT to come up with this. And I don’t feel guilty as it is the result of a conversation (which I would have with others as well). The whole process from Idea to “visual” took me around 2 hours.

import geopandas as gpd
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
from shapely.geometry import LineString
import contextily as ctx

# === Load data ===
gdf = gpd.read_file("a2_sum.geojson")

# Reproject to Web Mercator for both basemap and length calc
if not gdf.crs or gdf.crs.to_epsg() != 25832:
    gdf = gdf.to_crs(25832)

gdf = gdf.copy()
gdf["start_x"] = gdf.geometry.apply(lambda geom: geom.coords[0][0])
gdf = gdf.sort_values("start_x", ascending=True).reset_index(drop=True)

# Now recompute linear positions for straightened plot
gdf["length_km"] = gdf.geometry.length / 1000.0
gdf["cum_length"] = gdf["length_km"].cumsum()
gdf["start_pos"] = gdf["cum_length"] - gdf["length_km"]

# Reorder GeoDataFrame
#gdf = gdf.iloc[ordered_indices].reset_index(drop=True)

# Ensure numeric color attribute
color_col = "anzahl_befahrungen_sum"
if color_col not in gdf.columns:
    raise ValueError(f"Column '{color_col}' missing in dataset!")

# === Compute segment lengths and cumulative positions ===
gdf["length_km"] = gdf.geometry.length / 1000.0
gdf["cum_length"] = gdf["length_km"].cumsum()
gdf["start_pos"] = gdf["cum_length"] - gdf["length_km"]

# === Create a straightened version ===
straight_lines = []
for _, row in gdf.iterrows():
    start_x = row["start_pos"]
    end_x = row["cum_length"]
    straight_lines.append(LineString([(start_x, 0), (end_x, 0)]))
gdf["straight_geom"] = straight_lines

# === Prepare color scale ===
norm = mcolors.Normalize(vmin=gdf[color_col].min(), vmax=gdf[color_col].max())
cmap = plt.cm.viridis

# === Create figure with two subplots ===
# === Create figure with 2 vertical subplots ===
fig, (ax_main, ax_map) = plt.subplots(
    2, 1, figsize=(18, 6),
    gridspec_kw={"height_ratios": [1, 0.7]}
)
# === Main straightened plot ===
for idx, row in gdf.iterrows():
    color = cmap(norm(row[color_col]))
    ax_main.plot(
        [row["start_pos"], row["cum_length"]],
        [0, 0],
        color=color,
        linewidth=6,
        solid_capstyle="butt",
    )

# Label every 20th feature
for idx, row in gdf.iloc[::4].iterrows():
    name_from = str(row.get("mautknoten_nach", ""))
    start_x = row["start_pos"]

    if name_from:
        ax_main.text(
            start_x,
            0.2,
            name_from,
            rotation=45,
            ha="left",
            va="bottom",
            fontsize=7,
            clip_on=False,
        )

last_row = gdf.iloc[-1]
end_x = last_row["cum_length"]
end_name = str(last_row.get("mautknoten_nach", last_row.get("mautknoten_nach", "End")))

ax_main.text(
    end_x,
    -0.2,                    # slightly above the line
    end_name,
    fontsize=8,
    fontweight="bold",
    ha="right",
    va="bottom",
    color="darkblue",
    clip_on=False
)
first_row = gdf.iloc[0]
first_x = first_row["cum_length"]
first_name = str(first_row.get("mautknoten_von", last_row.get("mautknoten_von", "Start")))

ax_main.text(
    first_x,
    -0.2,                    # slightly above the line
    first_name,
    fontsize=8,
    fontweight="bold",
    ha="left",
    va="bottom",
    color="darkblue",
    clip_on=False
)


ax_main.set_ylim(-0.5, 0.6)
ax_main.set_yticks([])
ax_main.set_xlabel("Kilometres along A2")
ax_main.set_title("A2 – Linear Visualization (Straightened)")

sm = plt.cm.ScalarMappable(cmap=cmap, norm=norm)
sm._A = []
plt.colorbar(sm, ax=ax_main, label="Number of trips (sum of directions)")

# === Overview Map ===
gdf.plot(ax=ax_map, color="red", linewidth=2)
ctx.add_basemap(ax_map, crs=gdf.crs.to_string(), source=ctx.providers.CartoDB.Positron)
ax_map.set_axis_off()

# Adjust bottom subplot limits to show the full A2 geometry
x0, y0, x1, y1 = gdf.total_bounds
ax_map.set_xlim(x0-10000, x1+10000)
ax_map.set_ylim(y0-10000, y1+10000)
ax_map.set_aspect('equal', adjustable='box')  # preserves horizontal stretching
ax_map.set_axis_off()
ax_map.set_title("A2 Overview Map")
# --- Add labels for first and last segment ---
first_geom = gdf.geometry.iloc[0]
last_geom  = gdf.geometry.iloc[-1]

first_name = str(gdf["mautknoten_von"].iloc[0])
last_name  = str(gdf["mautknoten_nach"].iloc[-1])

# Use the centroid of each segment for label placement
ax_map.text(
    first_geom.centroid.x, first_geom.centroid.y + (y1 - y0)*0.01,
    first_name, fontsize=8, fontweight='bold', ha='center', va='bottom', color='black'
)
ax_map.text(
    last_geom.centroid.x, last_geom.centroid.y + (y1 - y0)*0.01,
    last_name, fontsize=8, fontweight='bold', ha='center', va='bottom', color='black'
)

plt.tight_layout()
plt.show()

plt.tight_layout()
plt.show()
linear view of the A2
The resulting map

Key takeaways

  • Data cleaning is critical: joins, filtering, duplicate removal, direction correction, and merging ensure your schematic reflects reality.
  • Straightening the line highlights traffic patterns immediately.
  • Combining schematic and geographic maps provides clarity and context.

Tomorrow, I work on polygons…

5 1 vote
Article Rating
Subscribe
Notify of
guest

This site uses Akismet to reduce spam. Learn how your comment data is processed.

0 Comments
Inline Feedbacks
View all comments