diff options
| -rw-r--r-- | fietsboek/geo.py | 132 |
1 files changed, 132 insertions, 0 deletions
diff --git a/fietsboek/geo.py b/fietsboek/geo.py new file mode 100644 index 0000000..c665705 --- /dev/null +++ b/fietsboek/geo.py @@ -0,0 +1,132 @@ +from dataclasses import dataclass +from itertools import islice +from math import sqrt, sin, cos, radians + +import gpxpy + + +# WGS-84 equatorial radius, also called the semi-major axis. +# https://en.wikipedia.org/wiki/Earth_radius +EARTH_RADIUS = 6378137.0 +"""Radius of the earth, in meters.""" + +# https://en.wikipedia.org/wiki/Preferred_walking_speed +MOVING_THRESHOLD = 1.1 +"""Speed which is considered to be the moving threshold, in m/s.""" + + +@dataclass +class MovementData: + duration: float = 0.0 + """Duration of the path, in seconds.""" + + moving_duration: float = 0.0 + """Duration spent moving, in seconds.""" + + stopped_duration: float = 0.0 + """Duration spent stopped, in seconds.""" + + length: float = 0.0 + """Length of the path, in meters.""" + + average_speed: float = 0.0 + """Average speed, in m/s.""" + + maximum_speed: float = 0.0 + """Maximum speed, in m/s.""" + + uphill: float = 0.0 + """Uphill elevation, in meters.""" + + downhill: float = 0.0 + """Downhill elevation, in meters.""" + + +@dataclass(slots=True) +class Point: + longitude: float + latitude: float + elevation: float + time_offset: float + + def distance(self, other: "Point") -> float: + """Returns the distance between this point and the given other point in + meters. + """ + r_1 = EARTH_RADIUS + self.elevation + r_2 = EARTH_RADIUS + other.elevation + # The formula assumes that 0° is straight upward, but 0° in geo + # coordinates is actually on the equator plane. + t_1 = radians(90 - self.latitude) + t_2 = radians(90 - other.latitude) + p_1 = radians(self.longitude) + p_2 = radians(other.longitude) + # See + # https://en.wikipedia.org/wiki/Spherical_coordinate_system#Distance_in_spherical_coordinates + # While this is not the Haversine formula for distances along the + # circle curvature, it allows us to take the elevation into account, + # and for most GPS point differences that we encounter it should be + # enough. + radicand = ( + r_1**2 + + r_2**2 - + 2 * r_1 * r_2 * ( + sin(t_1) * sin(t_2) * cos(p_1 - p_2) + + cos(t_1) * cos(t_2) + ) + ) + if radicand < 0.0: + return 0.0 + return sqrt(radicand) + + +class Path: + @classmethod + def from_gpx(cls, gpx: gpxpy.gpx.GPX) -> "Path": + points = [] + start_time = None + + for track in gpx.tracks: + for segment in track.segments: + for point in segment.points: + if start_time is None: + start_time = point.time + + time_offset = (point.time - start_time).total_seconds() + points.append(Point( + longitude=point.longitude, + latitude=point.latitude, + elevation=point.elevation, + time_offset=time_offset, + )) + + return cls(points) + + def __init__(self, points: list[Point]): + self.points = points + + def _point_pairs(self): + return zip(self.points, islice(self.points, 1, None)) + + def movement_data(self) -> MovementData: + """Returns the movement data.""" + movement_data = MovementData() + for a, b in self._point_pairs(): + distance = a.distance(b) + time = b.time_offset - a.time_offset + speed = distance / time + elevation = b.elevation - a.elevation + + movement_data.length += distance + if speed >= MOVING_THRESHOLD: + movement_data.moving_duration += time + else: + movement_data.stopped_duration += time + movement_data.maximum_speed = max(movement_data.maximum_speed, speed) + if elevation > 0.0: + movement_data.uphill += elevation + else: + movement_data.downhill += -elevation + movement_data.duration = b.time_offset + movement_data.average_speed = movement_data.length / movement_data.moving_duration + return movement_data |