Wiki page
[tesselations-4.0] by
sandro
2012-08-29 22:23:14.
D 2012-08-29T22:23:14.463
L tesselations-4.0
P 5233d7fdaca430825f1a3cd83bc2c489d37b6890
U sandro
W 15065
x<h2>Tesselation-related SQL functions supported in version 4.0.0</h2>
Back to <a href="https://www.gaia-gis.it/fossil/libspatialite/wiki?name=switching-to-4.0#tessellation">main page</a>
<h3>Generalities about tessellations</h3>
<i>Tessellation is the process of creating a two-dimensional plane using the repetition of a geometric shape with no overlaps and no gaps.</i>
<a href="http://en.wikipedia.org/wiki/Tessellation">read more</a><br><br>
A tessellation could be eventually based on indentical <i>cells</i>, all of exactly the same shape and size.
In this case we'll have a <i>regular tessellation</i>.<br>
Just a quick recall of elementary geometry; there are simply <b>3</b> regular polygonal shapes we can use in order to get a regular tessellation: the <b>equilateral triangle</b>, the <b>square</b> and the <b>regular hexagon</b>.<br><br>
On the other way many tessellations aren't regular at all, because each single cell has an individual size and shape of its own.<br>
It's really interesting to note that many <i>natural shapes</i> closely resemble a tessellation: going from biology to cristallography since geology and landscapes it's not at all difficult to identify many natural tessellation examples <i>on the wild</i>.
So it's not at all surprising to discover that tessellations are often useful in geography as well.
<h3>Setting up a testbed DB</h3>
In this short tutorial we'll use a very simple SpatiaLite DB, just containing the following Geometry tables:<ul>
<li>Administrative boundaries for Italy regions: you can <a href="http://www.istat.it/it/files/2011/04/reg2011.zip">download</a> this dataset from ISTAT (released on CC-BY license terms).<br></li>
<li>Main Italian populated places (namely, Local Councils): you can <a href="http://download.geonames.org/export/dump/allCountries.zip">download</a> this dataset from GeoNames (released on CC-BY license terms).<br>
Please note: this one is a worldwide dataset; italian populated places have then been extracted imposing the SQL clause<br>
<b>WHERE county_code = 'IT' AND population > 0</b></li>
</ul><br>
Just to keep any example as simple as possible, both datasets have been referenced into the <b>SRID=23032 - ED50 / UMT32 zone N</b>
<br><br><hr>
<h3>Creating regular tessellations in Spatial SQL</h3>
<table cellspacing="4" cellpadding="4">
<tr><td colspan="2">
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<i>generic prototypes</i>:<hr>
ST_xxxxGrid( input <i>Geometry</i>, size <i>double precision</i> ) : grid <i>Geometry</i><br>
ST_xxxxGrid( input <i>Geometry</i>, size <i>double precision</i>, edges_only <i>boolean</i> ) : grid <i>Geometry</i><br>
ST_xxxxGrid( input <i>Geometry</i>, size <i>double precision</i>, edges_only <i>boolean</i>, origin <i>Geometry</i> ) : grid <i>Geometry</i>
</td></tr></table>
<ul>
<li>the <b>input</b> Geometry is always expected to be a Polygon or a MultiPolygon, and will be exactly covered by the return grid.</li>
<li>the <b>size</b> argument identifies the edge length of the grid cell.</li>
<li>the facultative <b>edges_only</b> argument will be interpreted as follows:<ul>
<li>if <b>FALSE</b> (<i>default value</i>) a MultiPolygon will be returned.</li>
<li>if <b>TRUE</b> a MultiLinestring will be returned (simply representing the cells edges).</li>
</ul></li>
<li>the facultative <b>origin</b> Geometry is always assumed to be a Point, and will identify the grid's origin.<br>
By default a <b>(0, 0)</b> origin will be assumed.</li>
</ul></td></tr>
<tr><td>
<h4>using Square cells</h4>
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<b>SELECT ST_SquareGrid(geometry, 10000)<br>
FROM regions<br>WHERE cod_reg = 9;</b>
</td></tr>
</table><br>
This SQL query will return a regular grid (square cells) covering Tuscany (<b>cod_reg=9</b>).<br>
Each grid's cell will have an edge length of exactly <b>10 Km</b>
</td>
<td><img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/square-grid.png" alt="square grid"></td></tr>
<tr><td>
<h4>using Triangular cells</h4>
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<b>SELECT ST_TriangularGrid(geometry, 10000)<br>
FROM regions<br>WHERE cod_reg = 9;</b>
</td></tr>
</table></td>
<td><img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/tri-grid.png" alt="triangular grid"></td></tr>
<tr><td>
<h4>using Hexagonal cells</h4>
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<b>SELECT ST_HexagonalGrid(geometry, 10000)<br>
FROM regions<br>WHERE cod_reg = 9;</b>
</td></tr>
</table></td>
<td><img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/hex-grid.png" alt="hexagonal grid"></td></tr>
</table>
<br><hr>
<h3>Delaunay Triangulations</h3>
Triangulations simply represent a special case of tessellations: in this case all cells are represented by generic <b>triangles</b>, not necessarily of the equilateral kind.<br>
A <b>Delaunay Triangulation</b> is very peculiar, imposing a specific constraint.<br>
All triangles in a Delaunay triangulation must satisfy the <u><i>empty circle</i></u> property: i.e. for each edge we can find a circle containing the edge's endpoints but not containing any other points.
<a href="http://en.wikipedia.org/wiki/Delaunay_triangulation">read more</a><br><br>
Computing a Delaunay Triangulation representing many points is a very complex opeation, and one possibily imposing a huge computational load and may be requiring a long time.
Happily enough many highly efficient alghorithms have been already developed for the Dalaunay problem.<br>
The next-to-come <a href="http://trac.osgeo.org/geos/">GEOS 3.4.0</a> will support Delaunay Triangulations; this GEOS version is currently still under active development, but the <i>experimental</i> base-code (<i>trunk</i>) seems to be stable enough to be safely tested.
SpatiaLite already relies on GEOS for many tasks, so integrating a smooth support for Delaunay Triangulation as well wasn't at all difficult.<br><br>
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<i>prototype</i>:<hr>
ST_DelaunayTriangulation( input <i>Geometry</i> ) : delaunay <i>Geometry</i><br>
ST_DelaunayTriangulation( input <i>Geometry</i>, edges_only <i>boolean</i> ) : delaunay <i>Geometry</i><br>
ST_DelaunayTriangulation( input <i>Geometry</i>, edges_only <i>boolean</i>, tolerance <i>double precision</i> ) : delaunay <i>Geometry</i>
</td></tr></table>
<ul>
<li>the <b>input</b> Geometry can be of absolutely arbitrary type; all Linestrings and / or Polygowns will eventually be dissolved into Points corresponding to vertices.<br>
So after all ST_DelaunayTriangulation() will always process a MultiPoint.</li>
<li>the facultative <b>edges_only</b> argument will be interpreted as follows:<ul>
<li>if <b>FALSE</b> (<i>default value</i>) a MultiPolygon will be returned.</li>
<li>if <b>TRUE</b> a MultiLinestring will be returned (simply representing the triangles edges).</li>
</ul></li>
<li>the facultative <b>tolerance</b> argument is intended to normalize the input point-set, suppressing repeated (or too much close) points.<br>
By default a <b>0.0</b> tolerance will be assumed.</li>
</ul></td></tr></table>
<table cellspacing="4" cellpadding="4">
<tr><td>
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<b>SELECT ST_DelaunayTriangulation(ST_Collect(geometry))<br>
FROM italy_populated_places;</b>
</td></tr>
</table><br>
This SQL query will return a Delaunay Triangulation based on Italy's populated places (about 8,000+ Points).<br>
The visual example simply covers Tuscany, so the ensure an easy readibility.
</td>
<td><img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/delaunay.png" alt="delaunay triangulation"></td></tr>
</table>
<br><hr>
<h3>Voronoj Diagrams</h3>
<table cellspacing="4" cellpadding="4">
<tr><td><img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/voronoj-1.png" alt="voronoj-delaunay relationship">
</td><td>
The <b>Voronoj Diagram</b> simply is the <i><u>dual graph</u></i> of the <b>Delaunay Triangulation</b>.
A Voronoj Diagram still is a tessellation, and cells in a Voronoj can have any arbitrary polygonal (irregular) shape.<br><br>
This figure clearly shows the relation joining the Delaunay Triangulation and the corresponding Voronoj Diagram.<br><br>
Each single Voronoj's cell is obtained by connecting all the <i>circumcenters</i> of adjacent Delaunay's triangles; consequently each Voronoj cell surely contains one (and only one) Delaunay's node.<br><br>
The cell in the Voronoj Diagram presents an interesting property: all points falling within the same cell are ensured to be nearest to the Delaunay's node placed on the cell itself than to any other Delaunay's node placed in a different cell.<br>
So the Voronoj Diagram is a very effective conceptual tool allowing to divide an arbitrary space region in many <i>rational</i> cells, and is thus widely usued on many applicative fields.<br>
This including Geography, obviously.<br><br>
<a href="http://en.wikipedia.org/wiki/Voronoi_tessellation">read more</a>
<br></td></tr></table>
<a href="http://trac.osgeo.org/geos/">GEOS 3.4.0</a> will not support Voronoj; so the SpatiaLite own support is based on an original implementation (based in turn on the Delaunay support made available by GEOS).<br><br>
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<i>prototype</i>:<hr>
ST_VoronojDiagram( input <i>Geometry</i> ) : voronoj <i>Geometry</i><br>
ST_VoronojDiagram( input <i>Geometry</i>, edges_only <i>boolean</i> ) : voronoj <i>Geometry</i><br>
ST_VoronojDiagramn( input <i>Geometry</i>, edges_only <i>boolean</i>, extra_frame_size <i>double precision</i> ) : voronoj <i>Geometry</i><br>
ST_VoronojDiagramn( input <i>Geometry</i>, edges_only <i>boolean</i>, extra_frame_size <i>double precision</i>, tolerance <i>double precision</i> ) : voronoj <i>Geometry</i>
</td></tr></table>
<ul>
<li>the <b>input</b> Geometry can be of absolutely arbitrary type; all Linestrings and / or Polygowns will eventually be dissolved into Points corresponding to vertices.<br>
So after all ST_VoronojDiagram() (exactly as ST_DelaunayTriangulation) will always process a MultiPoint.</li>
<li>the facultative <b>edges_only</b> argument will be interpreted as follows:<ul>
<li>if <b>FALSE</b> (<i>default value</i>) a MultiPolygon will be returned.</li>
<li>if <b>TRUE</b> a MultiLinestring will be returned (simply representing the triangles edges).</li>
</ul></li>
<li>the facultative <b>extra_frame_size</b> xxxxxx
<li>the facultative <b>tolerance</b> argument is intended to normalize the input point-set, suppressing repeated points (simply used when internally computing the Delaunay Triangulation).<br>
By default a <b>0.0</b> tolerance will be assumed.</li>
</ul></td></tr></table>
<table cellspacing="4" cellpadding="4">
<tr><td>
<table bgcolor="#f0f0f0" cellspacing="4" cellpadding="8">
<tr><td>
<b>SELECT ST_VoronojDiagram(ST_Collect(geometry))<br>
FROM italy_populated_places;</b>
</td></tr>
</table><br>
This SQL query will return a Voronoj Diagram based on Italy's populated places (about 8,000+ Points).<br>
The visual example simply covers Tuscany, so the ensure an easy readibility.<br>
All <i>populated places</i> (aka cell <i>seeds</i>) are explicitly represented.
</td>
<td><img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/voronoj-2.png" alt="voronoj diagram"></td></tr>
</table>
<br><hr>
<h3>Delaunay Triangulation, Convex Hull and Concave Hull</h3>
<table cellspacing="4" cellpadding="4">
<tr><td><img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/convex.png" alt="delaunay-convexhull relationship">
</td><td>
An intresting point absolutely worth to be explicitly noted: the boundary of any <b>Delaunay Triangulation</b> exactly corresponds to the <b>Convex Hull</b> for the same input Geometry.
<a href="http://en.wikipedia.org/wiki/Convex_hull">read more</a><br><br>
And this in turn opens the way to a further consideration: we could purposely simplify a Delaunay Triangulation so to get a <b>concave hull</b>.<br>
You can get more extensive informations about this approach from <a href="http://grass.osgeo.org/wiki/Create_concave_hull">here</a><br><br>
<u>Please note well</u>: the <b>convex hull</b> concept corresponds to a robust and formal mathematical definition. On the other side the <b>concave hull</b> is a much more vague and indetermined notion.<br>
There is one and only one ConvexHull for a given Geometry; but many ConcaveHulls are possible.
Choosing the one or the other is much more a matter of personal taste than a mathematical operation formally defined.<br>
Computing a ConcaveHull always is an inherently arbitrary and heuristic process.
</td></tr>
<tr><td>
<h4>The SpatiaLite's own approach to ConcaveHull</h4>
<ol>
<li>the Delaunay Triangulation corresponding to the input Geometry will be computed first.</li>
<li>then the statistical distribution of all triangle edge's lengths will be evaluated, so to determine <b>σ</b> (aka the <i>standard deviation</i>)</li>
<li>a second pass will now examine yet again all Delaunay's triangles; any triangle presenting at least one edge longer than <b>σ * factor</b> will be discarded.</li>
<li>and finally all filtered triangles will be dissolved so to form the ConcaveHull to be returned.</li>
</ol>
<u>Please note</u>: by setting an appropriate value to <b>factor</b> you can freely influence how much <i>aggressive</i> the filtering pass will be.<br>
Adopting a very high <b>factor</b> value practically means applying a very bland filtering (very few triangles will be discarded, and you'll consenquently get a rather <i>convex</i> shape).
On the other side adopting a very low <b>factor</b> values will apply a very strong filtering (many triangles will be now discarded, and you'll get a very <i>concave</i> shape).<br><br>
<u>Useful constants</u>: assuming a perfectly normal distribution of edge lengths (a by far unrealistic assumption for real world datasets), <b>3σ</b> will imply suppressing about <b>0.1%</b> triangles (only the few ones presenting abnormally lengthy edges), <b>2σ</b> corresponds to about <b>2.1%</b>, and <b>1σ</b> roughly corresponds to <b>15.8%</b>.<br>
Usually values ranging between <b>3σ</b> and <b>2σ</b> are the most appropriate to be used.<br><br>
The figure shows what you can actually get by applying a <b>3σ</b> filtering to Italian Populated Places; for clarity all Italian Regions are represented in yellow, and the ConvexHull is represented in blue.
</td><td>
<img border="1" src="http://www.gaia-gis.it/gaia-sins/write-view-pics/concave3.png" alt="concavehull-1">
</td></tr>
</table>
<br><hr>
Back to <a href="https://www.gaia-gis.it/fossil/libspatialite/wiki?name=switching-to-4.0#tesselation">main page</a>
Z 8c93429b7bb0b79e8541dc97ba670ddc