So, on the bike, I use a portable GPS to keep track of my speed and to track the mileage done on the bike so I know when to next put it in for service. Originally I just relied on the trip counter in the GPS, but then found that this could develop quite an error if left to tick over for a few months.
Thus, I wrote a simple CGI application in Perl and SQLite3 that would track the odometer readings. Plain, simple, and it’s worked quite well, but remembering to punch in the current odometer reading is a chore, and my stats are only as granular as I submit: if I want to see what distance I did on a particular day, I either have to have had the foresight to store readings at the start and end of that day, or I’m stuffed.
I also keep the GPX tracklogs. While the Garmin 650 is not great at handling lots of tracklogs (and for some moronic reason, they name the files “Day DD-MMM-YY HH.MM.SS.gpx”, not something sensible like “YYYY-MM-DDTHH-MM-SS.gpx”), it’s good enough that I can periodically siphon off the track logs for storage on my laptop. I then have a record of where I’ve been.
Theoretically, this also has the distance travelled, I could make a service that just consumes the GPX files, and tallies up the distances that way. Maybe even visualise heat-map style, where I go most. (No prizes for guessing “work” … but where else?)
The existing system uses SQLite, and specifically, its views, as poor man’s stored procedures. It’s hacky, inefficient, and sooner or later I’ll have performance problems. PostGIS is an extension onto PostgreSQL which supports a large number of spatial operations, including finding the length of a series of points, which is exactly the problem I’m trying to solve right now. The catch is, how do you import the data?
Enter GDAL
GDAL is a library of geographic functions for answering these kinds of questions. It ships with a utility ogr2ogr
, which can take geographic information in a variety of formats, and convert to a variety of output formats. Crucially, this tool supports consuming GPX files and writing to a PostGIS database.
Loading one file, is easy enough:
$ ogr2ogr -oo GPX_ELE_AS_25D=YES \
-dim 3 \
-gt 65536 \
-lco GEOM_TYPE=geography \
-preserve_fid \
-f PostgreSQL \
"PG:dbname=yourdb" yourfile.gpx \
tracks track_points
The arguments here were found by trial-and-error. Specifically, -oo GPX_ELE_AS_25D=YES
and -dim 3
tell ogr2ogr
to preserve the elevation in the point information (as well as keeping a copy of it in the ele
column). -lco GEOM_TYPE=geography
tells ogr2ogr
to use the geography
data type in PostGIS.
Look in the database, and you’ll see two tables, tracks
and track_points
. Sadly, you don’t get to choose the names of these (not easily anyway, there is -nln
, but it then will create one table with the given name, put the tracks in it, then blow it away and replace it with a table of the same name containing points), and there’s no foreign keys between the two.
The fun starts when you try to import a second GPX file. Run that command again, and because of -preserve_fid
, you’ll get a primary key clash. Take that away, and the track_fid
column in track_points
becomes meaningless.
If you drop -preserve_fid
, then track_fid
gets set to 0 for all points. Useless.
Importing many GPX files
Out of the box, this just wasn’t going to fly, so we needed to do things a little different. Firstly, I duplicated the schema that GDAL creates, creating my own tables which will ultimately store the data. I then used a wrapper shell script that calls psql
before and after ogr2ogr
so I can re-map the primary keys to maintain relationships.
Schema SQL
CREATE SEQUENCE public.gpx_points_ogc_fid_seq
INCREMENT 1
START 1
MINVALUE 1
MAXVALUE 2147483647
CACHE 1;
CREATE SEQUENCE public.gpx_tracks_ogc_fid_seq
INCREMENT 1
START 1
MINVALUE 1
MAXVALUE 2147483647
CACHE 1;
CREATE TABLE public.gpx_tracks
(
ogc_fid integer NOT NULL,
name character varying COLLATE pg_catalog."default",
cmt character varying COLLATE pg_catalog."default",
"desc" character varying COLLATE pg_catalog."default",
src character varying COLLATE pg_catalog."default",
link1_href character varying COLLATE pg_catalog."default",
link1_text character varying COLLATE pg_catalog."default",
link1_type character varying COLLATE pg_catalog."default",
link2_href character varying COLLATE pg_catalog."default",
link2_text character varying COLLATE pg_catalog."default",
link2_type character varying COLLATE pg_catalog."default",
"number" integer,
type character varying COLLATE pg_catalog."default",
gpxx_trackextension character varying COLLATE pg_catalog."default",
the_geog geography(MultiLineStringZ,4326),
CONSTRAINT gpx_tracks_pkey PRIMARY KEY (ogc_fid)
)
WITH (
OIDS = FALSE
)
TABLESPACE pg_default;
CREATE TABLE public.gpx_points
(
ogc_fid integer NOT NULL,
track_fid integer,
track_seg_id integer,
track_seg_point_id integer,
ele double precision,
"time" timestamp with time zone,
magvar double precision,
geoidheight double precision,
name character varying COLLATE pg_catalog."default",
cmt character varying COLLATE pg_catalog."default",
"desc" character varying COLLATE pg_catalog."default",
src character varying COLLATE pg_catalog."default",
link1_href character varying COLLATE pg_catalog."default",
link1_text character varying COLLATE pg_catalog."default",
link1_type character varying COLLATE pg_catalog."default",
link2_href character varying COLLATE pg_catalog."default",
link2_text character varying COLLATE pg_catalog."default",
link2_type character varying COLLATE pg_catalog."default",
sym character varying COLLATE pg_catalog."default",
type character varying COLLATE pg_catalog."default",
fix character varying COLLATE pg_catalog."default",
sat integer,
hdop double precision,
vdop double precision,
pdop double precision,
ageofdgpsdata double precision,
dgpsid integer,
the_geog geography(PointZ,4326),
CONSTRAINT gpx_points_pkey PRIMARY KEY (ogc_fid),
CONSTRAINT gpx_points_track_fid_fkey FOREIGN KEY (track_fid)
REFERENCES public.gpx_tracks (ogc_fid) MATCH SIMPLE
ON UPDATE RESTRICT
ON DELETE RESTRICT
)
WITH (
OIDS = FALSE
)
TABLESPACE pg_default;
The wrapper script
1
2
3 DB=tracklog
4
5 for f in "$@"; do
6 psql tracklog <<EOF
7 DROP TABLE IF EXISTS tracks;
8 DROP TABLE IF EXISTS track_points;
9 EOF
10 ogr2ogr -oo GPX_ELE_AS_25D=YES \
11 -dim 3 \
12 -gt 65536 \
13 -lco SPATIAL_INDEX=FALSE \
14 -lco GEOM_TYPE=geography \
15 -overwrite \
16 -preserve_fid \
17 -f PostgreSQL \
18 "PG:dbname=${DB}" "$f" \
19 tracks track_points
20
21
22 psql tracklog <<EOF
23 CREATE TEMPORARY TABLE track_fids AS
24 SELECT ogc_fid AS orig_fid,
25 nextval('gpx_tracks_ogc_fid_seq') AS ogc_fid
26 FROM tracks;
27
28 CREATE TEMPORARY TABLE point_fids AS
29 SELECT ogc_fid AS orig_fid,
30 nextval('gpx_points_ogc_fid_seq') AS ogc_fid
31 FROM track_points;
32
33 INSERT INTO gpx_tracks
34 SELECT track_fids.ogc_fid AS ogc_fid,
35 tracks.name as name,
36 tracks.cmt as cmt,
37 tracks."desc" as "desc",
38 tracks.src as src,
39 tracks.link1_href as link1_href,
40 tracks.link1_text as link1_text,
41 tracks.link1_type as link1_type,
42 tracks.link2_href as link2_href,
43 tracks.link2_text as link2_text,
44 tracks.link2_type as link2_type,
45 tracks."number" as "number",
46 tracks.type as type,
47 tracks.gpxx_trackextension as gpxx_trackextension,
48 tracks.the_geog as the_geog
49 FROM track_fids, tracks
50 WHERE track_fids.orig_fid=tracks.ogc_fid;
51
52 INSERT INTO gpx_points
53 SELECT point_fids.ogc_fid AS ogc_fid,
54 track_fids.ogc_fid AS track_fid,
55 track_points.track_seg_id AS track_seg_id,
56 track_points.track_seg_point_id AS track_seg_point_id,
57 track_points.ele AS ele,
58 track_points."time" AS "time",
59 track_points.magvar AS magvar,
60 track_points.geoidheight AS geoidheight,
61 track_points.name AS name,
62 track_points.cmt AS cmt,
63 track_points."desc" AS "desc",
64 track_points.src AS src,
65 track_points.link1_href AS link1_href,
66 track_points.link1_text AS link1_text,
67 track_points.link1_type AS link1_type,
68 track_points.link2_href AS link2_href,
69 track_points.link2_text AS link2_text,
70 track_points.link2_type AS link2_type,
71 track_points.sym AS sym,
72 track_points.type AS type,
73 track_points.fix AS fix,
74 track_points.sat AS sat,
75 track_points.hdop AS hdop,
76 track_points.vdop AS vdop,
77 track_points.pdop AS pdop,
78 track_points.ageofdgpsdata AS ageofdgpsdata,
79 track_points.dgpsid AS dgpsid,
80 track_points.the_geog AS the_geog
81 FROM track_points, track_fids, point_fids
82 WHERE point_fids.orig_fid=track_points.ogc_fid
83 AND track_fids.orig_fid=track_points.track_fid;
84
85 DROP TABLE tracks;
86 DROP TABLE track_points;
87 DROP TABLE track_fids;
88 DROP TABLE point_fids;
89 EOF
90 done
Getting the length of a track
Having imported all the data, we can do something like this:
SELECT ogc_fid, name,
ST_Length(the_geog, false)/1000 as dist_in_km
FROM gpx_tracks order by ogc_fid desc limit 10;
and get this:
1754 |
Day 20-JUL-18 18:09:02′ |
9.83689686312541′ |
|
1753 |
Day 15-JUL-18 09:36:16′ |
5.75919119415676′ |
|
1752 |
Day 14-JUL-18 17:12:24′ |
0.071734341651265′ |
|
1751 |
Day 14-JUL-18 17:12:23′ |
0.0729574875289383′ |
1750 |
Day 13-JUL-18 08:13:32′ |
9.88420745610283′ |
|
1749 |
Day 06-JUL-18 09:00:32′ |
9.81221316219109′ |
|
1748 |
Day 30-JUN-18 01:11:26′ |
9.77607205972035′ |
|
1747 |
Day 23-JUN-18 05:02:04′ |
19.6368592034475′ |
|
1746 |
Day 22-JUN-18 18:03:37′ |
9.91964760346248′ |
|
1745 |
Day 12-JUN-18 21:22:26′ |
0.0884092391531763′ |
Visualisation with QGIS
Turns out, this is straightforward…
- In your workspace, there’s a tree with the different layer types you can add, including PostGIS… right-click on this and select New Connection… fill in the details for your PostgreSQL database.
- Below that is XYZ Tiles…, right click again, select New Connection for OpenStreetMap, and use the URL https://a.tile.openstreetmap.org/{z}/{x}/{y}.png (also, see their policy).
- Drag the OpenStreetMap connection to your layers
- Expand the PostGIS connection you just made, and look for the
gpx_tracks
table, drag this on top of your OpenStreetMap layer.
Below is everywhere I’ve been with the GPS tracklog running. Much of what you see is the big loop a few of us did in 2012, including my trip to Ballarat for the 2012 LCA.
If I zoom in on Brisbane, unsurprisingly, some areas show up very clearly as being common haunts for me:
A bit of SQL voodoo, and I come up with this:
In orange is the territory covered on the Boulder (minus what was covered before I got the GPS), in blue the territory covered on the Talon 29 ER 0, and in red, on my current commuter (Toughroad SLR2).
Recent Comments