Simple Rest API for storing “point” observations

Database stuff

First I’ll describe the database that backs the service.

PostGIS database backend

Here is how to make a simple table in PostgreSQL, that can store geo-tagged “observations”. It uses a hstore type for key-value pairs and a geography point for the GPS dot. It’s very versatile, and could store anything from bird observations to endomondo like GPS tracks.

CREATE TABLE observations(
    id SERIAL PRIMARY KEY, 
    utc_timestamp TIMESTAMP, 
    geog GEOGRAPHY(Point, 4326), 
    kvp HSTORE
);

Continue reading “Simple Rest API for storing “point” observations”

Examples of querying a OSM PostgreSQL table with the hstore tags column

I have found some tutorials (1 and 2) for using the hstore column type in PostgreSQL. This blog post is specific to querying OpenStreetMap data using the tags column in relations ways and nodes (the tags column is a hstore).

Continue reading “Examples of querying a OSM PostgreSQL table with the hstore tags column”

Import OSM data into PostGIS using Osmosis

This blog post is made as a “note-to-self” so that I can remember the procedure. You are of course welcome to read along. It’s does nothing fancy, simply imports the planet.osm file into PostGIS using Osmosis with the Snapshot Schema.

Step by step

Assuming Osmosis is installed (if not download osmosis), and a planet.osm file has been download.

Continue reading “Import OSM data into PostGIS using Osmosis”

Plotting data on maps with matplotlib

I’m learning about matplotlib, and actually just bought the book Matplotlib for Python Developers.

Geographical plots

Browsing stackoverflow, the matplotlib homepage, and other resources, I eventually came by this stackoverflow post, which mentions BaseMap. Since the data that I’m plotting is inherently geographical, it makes sense to show the data on a map.

There are several nice examples on the basemap Github page.

Heatmaps

Often I want to create heatmaps of the data, using matplotlib.

On stackoverflow there are several posts on this topic:

There are different colormaps available for matplotlib, if you want to try different colorschemes.

Pyramidal tile cache cheat sheet

This table lists the number of total tiles for increasing zoom-level in a tile cache. The tile cache is assumed to be pyramidal: \left|z_{i}\right| = n \Rightarrow \left|z_{i+1}\right| = n \dot 4.

level 1: 1 tile
level 2: 5 tiles
level 3: 21 tiles
level 4: 85 tiles
level 5: 341 tiles
level 6: 1,365 tiles
level 7: 5,461 tiles
level 8: 21,845 tiles
level 9: 87,381 tiles
level 10: 349,525 tiles
level 11: 1,398,101 tiles
level 12: 5,592,405 tiles
level 13: 22,369,621 tiles
level 14: 89,478,485 tiles
level 15: 357,913,941 tiles
level 16: 1,431,655,765 tiles
level 17: 5,726,623,061 tiles
level 18: 22,906,492,245 tiles
level 19: 91,625,968,981 tiles
level 20: 366,503,875,925 tiles
level 21: 1,466,015,503,701 tiles

Hello world of raster creation with GDAL and Python

Mostly so I myself can remember how to do it, here is how to create a random geotiff with GDAL in Python

Note: the width and height are given in opposite order in the GDAL raster and numpy arrays!

import osr
import numpy
import gdal
import math
 
width = 4000
height = 3000
 
format = "GTiff"
driver = gdal.GetDriverByName( format )
 
dst_ds = driver.Create( "test.tiff", width, height, 1, gdal.GDT_Byte )
 
dst_ds.SetGeoTransform( [ 444720, 30, 0, 3751320, 0, -30 ] )
 
srs = osr.SpatialReference()
srs.ImportFromEPSG(25832)
dst_ds.SetProjection( srs.ExportToWkt() )
 
raster = numpy.zeros( (height, width), dtype=numpy.uint32 )
color_range = 2**8
seed = math.pi**10
for i in range(height):
	for j in range(width):
		color = (seed*i*j) % color_range
		raster[i][j] = color
 
dst_ds.GetRasterBand(1).WriteArray( raster )

It’s kind of slow, so perhaps the operation can be speeded up somehow? The result looks kind of nice though (image created with width and height both 4000):

So, not completely random.