Spatial analysis with SQL: 5- aggregation function (ArcGis Union)

In Article previous we have discussed some SQL functions for spatial processing . In this article we will see a function a little more complicated to put at work: the aggregation function for two layers of polygons. This operation is the one called Union in ArcGis. The purpose of the operation is to obtain a layer with all the intersections of the entities of the two original layers, and, also, all non-intersecting entities.

Let’s use the following example. We have a first layer “rectangles”

 

A second layer “circles”

 

That overlaps spatially:

 

What we want to obtain is a layer containing all the polygons of both layers but creating polygons to overlay areas:

 

We will proceed by following the same method used in the previous article, i.e. by developing an SQL query including several sub queries.

In the first place, we will develop a query that will build the polygons of the overlapping areas (the intersections) of both layers. These polygons appear in yellow in the next picture.

 

We will resume the request used in the previous article:

select circle.id as id1, rectangles.id as id2, st_intersection (circle.geom, rectangles.geom ) as the_geom
from circle, rectangles
where st_intersects ( circle.geom, rectangles.geom )
group by circle.id, rectangles.id

You will notice that the select is built not just to retrieve intersected geometries, but also the identifiers of the original entities (circle.id and rectangles.id). Of course, you can recover other attributes, or all by changing this part of the request.

To these polygons, we will add the rest of the areas occupied by the rectangle layer (in yellow in the following picture)

 

 

select 0 as id1, rectangles.id as id2, st_difference ( rectangles.geom,

(select multi ST_ (ST_Union (the_geom)) as the_geom from (
select st_intersection (circle.geom, rectangles.geom ) as the_geom
from circle, rectangles
where st_intersects ( circle.geom, rectangles.geom ))
as the_geom)) as the_geom
from rectangles

What we get to do with this request is to extract all the intersection areas (as in the previous query) and convert them in a single multi-polygon. The request removes the difference of this multi- polygon with the layer of rectangles.

Finally, these polygons, we are going to add the rest of the areas occupied by the layer circles to these polygons (yellow in the following image)

 

We will proceed likewise the rectangles, by using a multi- polygon of intersection areas and calculating the difference with the layer circles:

 

select circle.id as id1, 0 as id2 , st_difference ( circle.geom , (
select st_multi ( st_union ( the_geom )) as the_geom from (
select st_intersection ( circle.geom, rectangles.geom ) as the_geom
from circle, rectangles
where st_intersects ( circle.geom, rectangles.geom ))
as the_geom )) as the_geom
from circle

The union of these three requests produces the outcome layer we expected :

 

select *
from
(select circle.id as id1, rectangles.id as id2, st_intersection ( circle.geom, rectangles.geom ) as the_geom
from circle, rectangles
where st_intersects ( circle.geom, rectangles.geom )
group by circle.id, rectangles.id
union all
select 0 as id1, rectangles.id as id2, st_difference ( rectangles.geom , (
select st_multi ( st_union ( the_geom )) as the_geom from (
select st_intersection ( circle.geom, rectangles.geom ) as the_geom
from circle, rectangles
where st_intersects ( circle.geom, rectangles.geom ))
as the_geom )) as the_geom
from rectangles
union all
select circle.id as id1, 0 as id2, st_difference ( circle.geom , (
select st_multi ( st_union ( the_geom )) as the_geom from (
select st_intersection ( circle.geom, rectangles.geom ) as the_geom
from circle, rectangles
where st_intersects ( circle.geom, rectangles.geom ))
as the_geom )) as the_geom
from circle
) as T1  

 

Displaying marine charts with S57 QGis

The sources of information for the coastal regions are not limited to classical terrestrial sources. The marine charts can provide information very helpful and they allow creating rich and impactful documents.

The most adapted sources to quickly integrate a great deal of information are the marine electronic charts (ENC). In a single data set you find all information contained in a chart: coastline, bathymetry, navigation aids, etc…

Unlike the USA and other countries where the charts are free and available to download, in France one has to pay in order to gain access to them, which slows very much their dissemination and use.

If you are an ArcGis user, you count with an extension developed by NOAA to integrate ENC data in ArcMap. For QGis, there is not such extension. Continue reading “Displaying marine charts with S57 QGis”

ArcGis Model Builder: FOR and WHILE iterators

In the previous article ( ArcGIS Model Builder: expand a model to a set of input data ) we have addressed the most common iterators used with Model Builder. These iterators make it possible to apply a processing model to a set of input layers.

But there are two other types of iterators available:

  • the iterators FOR and WHILE
  • the iterators that browse entities of a single layer

In this article we will discuss the iterators FOR and WHILE Continue reading “ArcGis Model Builder: FOR and WHILE iterators”

1-n joins in QGis: virtual layers

One of the most common techniques in GIS in order to treat non-spatial data is that of the join. We use a geographical layer having geometries, as support for another layer having no location data. To this end, we will use a common field to both tables, that allows “the joining” the records from the geographic table to those of the non- geographic table. Therefore we created a new virtual table where the attributes of the non- geographic table can be used to be mapped with the geometry of the first table. For this join to work it is necessary that for every registration of geographical table corresponds a record in the non- geographic table, and only one. Because if more than one is found, the join cannot work. In that case we have a “relationship” between the tables, never a “juncture”.

For a long time, GIS software has, only, allowed to work with joins. The relationships do not stand most spatial operators applicable to joins. Continue reading “1-n joins in QGis: virtual layers”

How to create Tanaka’s illuminated contours with QGis

The Tanaka’s contours are a method of displaying terrain developed by Professor Tanaka Kitiro in 1950. Tanaka called his technique the “Relief Contours Method”. However, it is usually called the “Illuminated Contour Method” or the “Tanaka Method”.

This method applies a light source north-west at a contour map. The result is a 3D style representation of the terrain. The method involves changing the width and colour of the contour lines based on their position relative to the light source: The contour lines that face the light source are drawn white while those found in the shadow are traced in black; the lines perpendicular to the light source are drawn finer than the parallels. In summary, this method gives a realistic shadowy effect to the terrain maps and allows the viewer to better understand the features of the landscape. Continue reading “How to create Tanaka’s illuminated contours with QGis”

Model Builder: how to convert a batch model to an interactive model

The principle of Model Builder is to create geoprocessing tools, which are by definition “batch” tools. This means that there is no   interaction between the user and the processing, once launched. You define the input data and you retrieve the output data.
In some occasions this aspect is binding. Imagine that you have to run a model in using only a portion of the data set and if you use the Clip tool you will be duplicating the data, or else you want to select an area on the map and apply your geoprocessing model only to entities completely contained in this area. There is another possibility, the one we are going to follow in this example, we will use a template with a specific tool that requires an entry point, such as “Field of View” and you want to enter manually the, or the, locations of the observer. Continue reading “Model Builder: how to convert a batch model to an interactive model”

How to calculate local averages in the Voronoï polygons with Qgis

The Voronoi maps are built based on a series of polygons created around the location of each sampling point.

The Voronoï polygons are created in such a way that each location in a polygon is closer to the sampling point present in this polygon than any other sampling point .

In the article regarding the Exploratory data analysis , we discussed how to use the ArcGIS Geostatistical Analyst to build and analyse the Voronoi polygons.

In this article we will discuss how to get the same results with QGis. Continue reading “How to calculate local averages in the Voronoï polygons with Qgis”

Share your resources (symbols, images, scripts) for Qgis with the new Resource Sharing plugin

The QGis resource sharing plugin is now available. This plugin makes it possible to share all type of resources for QGis: libraries of symbols , styles, images and scripts. These items can be stored   remotely (for now the plugin supports Github and Bitbucket as public repositories), the local file system, or an http server with zipped collections.

The installation is performed in a classic way with the Manager extensions ( Extensions menu- > Install / Manage extensions ).

Continue reading “Share your resources (symbols, images, scripts) for Qgis with the new Resource Sharing plugin”