Indexof

IndexofVector Skeletonization: How to Generate Centerlines Using OGR and Python › Last update: Mar 18, 2026@poothangAbout › #VectorSkeletonization

Linear Approximation: Generating Centerlines Directly via OGR

In transport and utility mapping, we often encounter "dual-line" or "casing" geometries—two parallel lines representing the edges of a road or a wide pipe. Converting these into a single, routable centerline is a common spatial challenge. While many GIS workflows involve a "Line-to-Polygon-to-Centerline" conversion, this approach is computationally expensive and prone to topological errors during the conversion steps. This tutorial focuses on a direct vector approach using OGR (GDAL) and Python, utilizing medial axis approximation and interpolation techniques to derive centerlines while bypassing the creation of intermediate polygon features.

Table of Content

Purpose

Generating centerlines directly from existing line strings serves several technical needs:

  • Network Topology: Creating a routable graph where only a single line represents a multi-lane corridor.
  • Database Optimization: Reducing the storage footprint of linear assets by 50% without losing the central trajectory.
  • Annotation Placement: Providing a central anchor point for dynamic labeling in web maps.

Geometric Concepts: Medial Axis vs. Interpolation

Without polygons, we rely on Voronoi-based Medial Axis approximation. By densifying the vertices of the two parallel lines and calculating the Voronoi vertices between them, the points that fall roughly equidistant from both lines form the skeleton. Alternatively, for simple parallel lines, Linear Interpolation at specific intervals can be used to calculate a mid-point sequence.

Step-by-Step: Extracting Centerlines with Python/OGR

1. Prepare the Input Geometries

Ensure your "casing" lines are part of the same feature or can be easily grouped. In OGR, we load these as wkbLineString.

2. Densify the Vertices

A centerline is only as accurate as its vertex density. Use the Segmentize method to add vertices every few meters.

from osgeo import ogr

# Open the source
ds = ogr.Open("road_edges.shp", 1)
layer = ds.GetLayer()

for feature in layer:
    geom = feature.GetGeometryRef()
    # Add vertices every 1 meter
    geom.Segmentize(1.0)

3. Generate the Mid-Point Skeleton

Using the densified vertices of Line A and Line B, we find the nearest neighbor for each vertex and calculate the midpoint. In a more robust script, this involves the Distance and Interpolate functions.

# Conceptual logic for midpoint calculation
def get_midpoint(p1, p2):
    return ogr.Geometry(ogr.wkbPoint)
    x = (p1.GetX() + p2.GetX()) / 2
    y = (p1.GetY() + p2.GetY()) / 2
    # Create new line point...

4. Reconstruct the LineString

Assemble the resulting mid-points into a new ogr.wkbLineString and apply a smoothing filter (Laplacian or moving average) to remove any "zig-zag" artifacts caused by uneven vertex distribution.

Use Case: Hydrographic Channel Mapping

A maritime agency has high-resolution shoreline vectors for a narrow canal and needs a navigation track for vessels.

  • The Challenge: The canal is 50km long; converting this to a polygon for skeletonization crashes the local memory.
  • The Action: Using the OGR Segmentize and a midpoint interpolation script, they process the lines in 1km chunks.
  • The Result: A smooth navigation centerline is produced that remains perfectly centered regardless of the canal's varying width.

Best Results

Technique Pros Cons
Linear Interpolation Fastest; works well for strictly parallel lines. Fails at sharp bends or varying widths.
Voronoi Skeleton High accuracy; handles complex shapes. Requires external libraries (like SciPy) to filter noise.
OGR Segmentize Native to GDAL; preserves projection data. Adds significant vertex overhead to the dataset.

FAQ

Can I use OGR's Buffer tool to find a centerline?

Not directly. Buffering creates polygons. To stay in the vector domain, you would need to use a "negative buffer" (erosion), but this eventually collapses the geometry rather than providing a clean line string.

What happens at 'dead ends'?

In a non-polygon workflow, dead ends often result in "dangling" vertices. You must implement a threshold to prune lines that are shorter than the average width of the dual-line casing.

Does this work on MultiLineStrings?

Yes, but you must ensure the segments are ordered sequentially. If your lines are digitized in opposite directions, the midpoint calculation will "cross over" itself, creating a knot. Always verify line directionality using ST_Line_Locate_Point logic.

Disclaimer

Centerline extraction without a constrained polygon area is a mathematical approximation. It assumes that your two input lines are logically paired. If a third line enters the proximity, the interpolation logic will fail unless sophisticated spatial partitioning is used. March 2026.

Tags: OGR_Geometry, Centerline_Extraction, Spatial_Python, Vector_Algorithms



What’s new

Close [x]
Loading special offers...