GPX to line overlaid atop USGS Urban Area photo

The other day I wrote about a script to simplify GPS plots for plotting polylines using Google Maps. The momentum from that finally propelled me to learn enough about this stuff to be able to have the script output this:

image

(The inlined image is resized to fit into my Wordpress layout, it links to the full image. I also rendered using the topographic map tile set.)

This is a 4 meters/pixel resolution render of the same track used in the previous post to TerraServer-served USGS Urban Area photos. To make it, the script:

  1. Parse the GPX file into lat/long (same as it did before
  2. Convert all the lat/long points to Universal Transverse Mercator (UTM) points — this is a convenient coordinate system, particularly since it is what TerraServer uses..
  3. Compute the TerraServer image tiles corrosponding to the UTM coordinates covered by the track. This is straightforward division. Since I wanted 4mpp resolution, I used Scale 12 when hitting the TerraServer tile URL. This meant converting to TerraServer tiles is just “divide by 800”
  4. Fetch all the tiles and assemble them into an image (the script caches the tiles it fetches so repeated runs covering the same geographic area are fast)
  5. Plot the track over the map

pygps.org was a big help here, for both the LatLongUTMconversion library and the source code to the mapviewer, which clued me into the easy way to fetch TerraServer tiles. I had spent quality time browsing the TerraServer SOAP documentation, but never quite got a SOAP implementation for Python to work.

I modified the LatLongUTMconversion library to work with Numeric Python arrays.

UTM is also convenient because it dispenses with lots of messy distance-on-sphere distance calculations to figure out the distances between points — UTM coordinates are in meters and use a square grid. So I can use the usual Cartesian distance calculation for distances. My code is probably little dodgy on UTM zone boundaries but fortunately I have no tracks that cross a zone boundary so I can put off dealing with that for the time being.

The track shown there is also interpolated using a spline. I dug the code for splines out of ancient C++ code I wrote back in the last millenium for a college independent study. After I converted it to Python, it fit into my script and worked. The track is smoother and most of it looks better but it behaves badly in some cases.

I need to do some (more) research about how this works so I can improve it. Topofusion does what I want, I just need to learn how it works. TopoFusion is a great app — it’s what I actually use to manage to look at most of my GPS tracks. Working on this script is about learning about mapping and visualizing rather than addressing some problem with TopoFusion.

Making this all work has been educational. Most of what I’ve done I could use USGS DOQ and DRG files rather than TerraServer as a map/photo image source. TerraServer is very convenient, though.