Smathermather's Weblog

Remote Sensing, GIS, Ecology, and Oddball Techniques

Posts Tagged ‘Voronoi’

Voronoi in PostGIS

Posted by smathermather on December 21, 2013

PostGIS has for a time had ST_DelaunayTriangles for Delaunay Triangulation, and since 2.1 apparently can even natively create a 2.5D TIN using said function, which is pretty cool. I think with SFCGAL, we will eventually have true 3D TINs as well.

Overlay of Delaunay and Voronoi

We’re still waiting on native Vororoi support for PostGIS though. According to http://trac.osgeo.org/postgis/wiki/UsersWikiPostgreSQLPostGIS

“GEOS 3.5 is still in development but will have Voronoi which will be a feature in 2.2 only if you have GEOS 3.5.”

Vishal Tiwari through a Google of Summer of Code under the mentorship of Sandro Santilli completed the port of Voronoi to GEOS from JTS Topology Suite. Now we just need to wait for PostGIS 2.2….

In the mean time,

http://geogeek.garnix.org/2012/04/faster-voronoi-diagrams-in-postgis.html

One caveat– this python function doesn’t always provide just Voronoi but some artifact polygons.

Voronoi polygons with some artifacts

Once we have a table with Voronoi, we can filter out just the true Voronoi cells by counting the number of original points we find within them, and only return the polygons which contain a single point:

CREATE TABLE true_voronoi AS
WITH tvoronoi AS (
	SELECT COUNT(*), v.geom
		FROM voronoi v, input_points p
		WHERE ST_Intersects(v.geom, p.geom)
		GROUP BY v.geom
		)
SELECT the_geom FROM tvoronoi WHERE count = 1;

voronoi without extra polygons

But it’s still not a perfect solution. I can’t wait for PostGIS 2.2….

Posted in Analysis, Database, PostGIS, PostgreSQL, SQL | Tagged: , , | 4 Comments »

Way beyond red-dot fever– bees hives and overlapping home ranges

Posted by smathermather on April 6, 2013

Last week, post Boston Code Sprint, I spent a couple of hours playing with bee data, specifically bee keeper data for Norfolk County Massachusetts.

My friend Eric (a bee keeper of 4 hives in said county) says that worker bees can fly as far as 3 miles for nectar, but after that approximate limit, they hit diminishing returns relative to the calories they burn.

Proximity of bee hives is important for reasons of competition (although bees are quite friendly with each other on the flower) but also disease transfer, which in this age of not-so-easy bee keeping, is an important factor.  So we set out to map bees.

With a list of addresses, QGIS, and the mmqgis extension (see Plugins:Fetch Python Plugins)  which uses Google services for address lookup, we get a nice little set of points for beekeepers in Norfolk County.

01 points

Having a good context map for any such set of points is important.  This is the part of the project I worried most about.  Fortunately, we also have the OpenLayers Plugin, which allows us to add a number of web maps as background (BTW, when we want print quality maps down the road, this will cause us suffering, but for now it is a nice stand-in for a good base layer).  We decided the colorless Stamen Toner map was perfect for our purposes.

02 points toner

Now, let’s buffer these points– give those bees 1.5 miles of free reign, and see where the overlaps in range are (we could make these 3 mile radii, but you get the point):

03 buffer

The map tells us something of overlapping ranges, etc., but it seems like we could use something like we could do something more sophisticated.  Anyone who has followed my blog for a while will not be surprised by the employment of Voronoi polygons… .

04 voronoi

Eric, having suffered an explanation of different ways to calculate Voronoi earlier in the week, remembered the relationship between buffers and Voronoi, and asked if we could combine them.   In this case, picture soap bubbles of your youth.  Blow a soap bubble for each bee point, and the walls of bubbles in between form Voronoi polygons.  In our case, we use cones instead of soap bubbles, but after a couple hours of playing in PovRay, we made some code like this:

#version 3.6;

// Pov Includes
#include “colors.inc”
#include “transforms.inc”
#include “cone_coords.inc”

#declare Camera_Size = 500000; //freekin huge

background {color <0, 0, 0>}

//Set the center of image
#declare scene_center_x=759496;
#declare scene_center_y=2861454;

#declare Camera_Location = <scene_center_x,5,scene_center_y> ;
#declare Camera_Lookat = <scene_center_x,0,scene_center_y> ;

// Use orthographic camera for true plan view
camera {
orthographic
location Camera_Location
look_at Camera_Lookat
right Camera_Size*x
up Camera_Size*y
}

// Union all the cones together into one object

#declare coner = union {
cone { <0,0,0>,15840,<0,5,0>,0

} };

union {

#declare LastIndex = dimension_size(cone_coords, 1)-2;
#declare Index = 0;
#while(Index <= LastIndex)
object {
coner
translate cone_coords[Index]
}
#declare Index = Index + 1;
#end

// Pigment cones according to distance from camera
pigment {
gradient x color_map {
[0 color Blue ] // Blue to orange will help us employ blue as the alpha values later
[1 color Orange ]
}
scale <vlength(Camera_Location-Camera_Lookat),1,1>
Reorient_Trans(x, Camera_Lookat-Camera_Location)
translate Camera_Location
}
finish {ambient 1}
}

We get a funky result (yes, this is the result I hoped for):

05 raw pov blue 06 dialog 1

A little transparency, alternate QGIS blending mode, and using the blue band as a transparency band give us a nice output:

07 final 1

We’ll add in those Voronoi boundaries to emphasize the edge of these areas:

07 final

And since this is really only meant for Norfolk County bee keepers, and doesn’t represent folks outside that area, we’ll use a mask with a blending mode of “Subtract” to really offset this map:

08 mask norfolk 09 norfolk masked

I will be interested to see if this helps them look at disease control and spread.  In case you are wondering, Eric is luck enough to have an almost dedicated area for his bees, which includes a huge wetland complex (think: lots of flowers) that is (in his words) all his own.  Well, all his bees own, anyway.

Posted in Analysis, Optics, Other, POV-Ray, QGIS | Tagged: , , , | 2 Comments »