Build a GeoJSON file from a collection of GPX files.

gpx2geojson.py

import dateutil.parser
import fiona
import glob
import json
import math
import pytz
import re
import sys

from pathlib import Path

def guess_activity(file):
    m = re.search("\w\.(.+?)\.", file)
    if (m is None):
        return "hike"
    if (m.group(1) == "car"):
        return "drive"
    if (m.group(1) == "paddle"):
        return "paddle"
    if (m.group(1) == "ski"):
        return "ski"
    return None

assert guess_activity("123.gpx") == "hike"
assert guess_activity("123.car.gpx") == "drive"
assert guess_activity("123.paddle.gpx") == "paddle"
assert guess_activity("123.ski.gpx") == "ski"
assert guess_activity("123.foo.gpx") == None

def distance(x, y):
    return math.sqrt((x[0] - y[0])**2 + (x[1] - y[1])**2)

assert distance((-1, -1), (1, 1)) == 2.8284271247461903

def to_feature(timestamp, source, activity, coordinates):
    return { 
        "type": "Feature", 
        "properties" : { 
            "timestamp": timestamp.isoformat(),
            "source": source, 
            "activity": activity,
        }, 
        "geometry": { 
            "type": "LineString", 
            "coordinates": coordinates 
        } 
    }

def gpx_to_geojson(path, target, tz=pytz.timezone("America/Los_Angeles")):
    features = []
    for file in Path(path).rglob("*.gpx"):
        activity = guess_activity(file.name)
        print(f"{file} -> {activity}")
        if not activity:
            continue
        with fiona.open(file, layer="track_points") as points:
            timestamps = []
            coordinates = []
            for point in points:
                if not point["properties"]["time"]:
                    print(f"  missing timestamps")
                    break
                t = dateutil.parser.isoparse(point["properties"]["time"] + "Z").astimezone(tz)
                c = point["geometry"]["coordinates"]
                if c[0] == 0 or c[1] == 0:
                    print(f"  invalid coordinates: {c}")
                    continue
                if len(coordinates) > 1:
                    t_delta = t - timestamps[-1]
                    c_delta = distance(c, coordinates[-1])
                    if t_delta.total_seconds() > 3600 or c_delta > 0.05: # ~5km
                        features.append(to_feature(timestamps[0], file.name, activity, coordinates))
                        timestamps = []
                        coordinates = []
                timestamps.append(t)
                coordinates.append(c)
            if len(coordinates) > 1:
                features.append(to_feature(timestamps[0], file.name, activity, coordinates))

    with open(target, "w") as out:
        json.dump({ 
            "type": "FeatureCollection", 
            "features": features 
        }, out, indent=2)

if __name__ == "__main__":
    assert len(sys.argv) == 3
    gpx_to_geojson(sys.argv[1], sys.argv[2])