Category Archives: GIS

Sistemas de Información Geográfica

Geocamp ES 2017

Quedan poco más de dos semanas para la GeocampES 2017, que este año se va a celebrar el sábado 16 de Septiembre en Almería, organizada por Geoinquietos AlmeríaHacklab Almería y UNIA. Como cada año la idea es organizar un evento sencillo, divertido y muy técnico y participativo. No hay agenda, siguiendo el modelo de desconferencia habitual, al llegar allí se espera que un buen porcentaje de los asistentes se postulen para salir a contar algo: puede ser una charla convencional, una demostracíon de algo en lo que estás trabajando, una dinámica de grupo que te parezca interesante, etc. Cualquier actividad con una componente geo bien definida (o no tanto) es bienvenida.

La idea de Geocamp ES se tomó siguiendo el liderazgo de los compañeros de Portugal. A diferencia sus geocamps en pequeños pueblos más o menos alejados de todo, las anteriores Geocamp españolas se han celebrado en ciudades. En cualquier caso este año nos vamos a Almería, que está un poco más difícil de acceder (sin exagerar) pero sigue siendo una ciudad que seguro va a ofrecernos un espacio agradable y totalmente adecuado para pasarlo bien aprendiendo.

Si lo que te he contado te parece mínimamente interesante, anímate y prepara un fin de semana por el sur de España con el resto de la comunidad geoespacial.

Apúntate aquí.


Aggregating points: JSON on SQL and loops on infowindows

NOTE: I’ll use CARTO but you can apply all this to any webmapping technology backed by a modern database.

Get all the data

So we start with the typical use case where we have a one to many relationship like this:

    select e.cartodb_id,
           l.cartodb_id as locaction_id,
      from locations l
inner join employees e
        on e.location = l.location
  order by location

Easy peasy, we have a map with many stacked points. From here you can jump to this excellent post by James Milner about dense point maps. My example is not about having thousands of scattered points that at certain zoom levels overlap. Mine is a small set of locations but many points “stacking” on them. In this case you can do two things: aggregate or not. When you aggregate you pay a prize for readability: reducing all your data to those locations and maybe using visual variables to show counts or averages or any other aggregated value and finally try to use the interactivity of your map to complete the picture.

So at this point we have something like this map, no aggregation yet, but using transparency we can see where CARTO has many employees. We could also use a composite operation instead of transparency to modify the color of the stacked points.

Stacking points using transparency

Stacking points using transparency

Aggregate and count

OK, let’s do a GROUP BY the geometry and an aggregation like counting. At least now we know how many people are there but that’s all, we loose the rest of the details.

    select l.the_geom_webmercator,
           min(e.cartodb_id) as cartodb_id,
           count(1) as counts
      from locations l
inner join employees e
        on e.location = l.location
  group by l.the_geom_webmercator
Grouping by location and counting

Grouping by location and counting

Aggregate one field

But in my case, with CARTO we have PostgreSQL at hand so we can do way more than that. PostgreSQL has many many cool features, handling JSON types is one of them. Mix that with the fact that almost all template systems for front-end applications allow you to iterate over JavaScript Objects and you have a winner here.

So we can combine the json_agg function with MustacheJS iteration over objects to allow rendering the names of our employees.

    select l.the_geom_webmercator,
           min(e.cartodb_id) as cartodb_id,
           json_agg(e.firstname) as names, -- JSON aggregation
           count(1) as counts
      from locations l
inner join employees e
        on e.location = l.location
  group by l.the_geom_webmercator,l.location

And this bit of HTML and Mustache template to create a list of employees we can add to the infowindow template:

<ul style="margin:1em;list-style-type: disc;max-height:10em;">
{{#names}}<li class="CDB-infowindow-title">{{.}}</li>{{/names}}

List of employees on the infowindow

We could do this without JSON types, composing all the markup in the SQL statement but that’s generating quite a lot of content to move to the frontend and of course making the whole thing way harder to maintain.

Aggregate several fields

At this point we can repeat the same function for the rest of the fields but we need to iterate them separatedly. It’d be way better if we could create JSON objects with all the content we want to maintain in a single output field we could iterate on our infowindow. With PostgreSQL we can do this with the row_to_json function and nesting an inner double SELECT to give the properties names. We can use directly row_to_json(row(field1,field2,..)) but then our output fields would have generic names.

    select l.the_geom_webmercator,
           min(e.cartodb_id) as cartodb_id,
           count(1) as counts,
             SELECT r
               FROM (
                 SELECT photourl as photo,
                        coalesce(preferredname,firstname,'') as name
             ) r
           ),true)) as data
      from solutions.bamboo_locations l
inner join solutions.bamboo_employees e
        on e.location = l.location
  group by l.the_geom_webmercator,l.location
  order by counts asc

With this query now we have a data field with an array of objects with the display name and web address for the employee picture. Easy now to compose this in a simple infowindow where you can see the faces and names of my colleagues.

<div style="column-count:3;">
<span style="display:inline-block;margin-bottom:5px;">
  <img style="height:35px;" src="{{photo}}"/> 
  <span style="font-size:0.55em;">{{name}}</span>


Adding pictures and names

That’s it. You can do even more if you retrieve all the data directly from your database and render on the frontend, for example if you use D3 you probably can do fancy symbolizations and interactions.

One final note is that if you use UTF grids (like in these maps with CARTO) you need to be conservative with the amount of content you put on your interactivity because with medium and big datasets this can make your maps slow and too heavy for the front-end. On those cases you may want to change to an interactivity that works like WMS GetFeatureInfo workflow, where you retrieve the information directly from the backend when the user clicks on the map, instead of retrieving everything when loading your tiles.

Check the map below and how the interactions show the aggregated contents. What do you think of this technique? Any other procedure to display aggregated data that you think is more effective?

How a daily digest of geospatial links is distributed

TL;DR If you are interested on getting a daily digest of geospatial links subscribe to this mailing list or this atom feed. Take «daily» with a grain of salt.

Over the last six years Raf Roset, one of my favourite geonerds out there, has been sending all the cool stuff he founds about our geospatial world to Barcelona mailing list at OSGeo mailman server. He started circa 2011 sending one link per mail, but in 2013-04-03 he started to make a daily digest. A gun burst in Spanish is called Ráfaga so the joke was really at hand when someone proposed to call those digests that way.

Time passes, September 2014 and I ask Raf to send them also to Valencia mailing list, since most people there understand Catalan and the content was too good to be enjoyed only by our loved neighbours. Finally in January 2015 I decide to start translating them into Spanish and send them also to Spanish and Seville mailing lists.

Then in May I join CARTO and @jatorre thinks is a good idea if I can send them to the whole company mailing list so after some weeks I stop translating them into Spanish. Since that day I only do it English, trying to follow Raf lead everyday translating his mails and forwarding them to CARTO internal mailing list and the rest of the OSGeo ones.

Also at June I decided to put those mails in a simple website so the Ráfagas would also be accessible on GitHub and a static jekyll website so anyone could use the Atom feed to reach them.

Final chapter, in July I also decide to create a dedicated mailing list just for those people who are only interested in receiving those digest mails, obviously thinking in a broader audience, not just my fellow friends from Spain. I think at some point I will stop sending them to the Spanish lists because normally Ráfagas don’t fire any discussion and I’m sending the same message to three lists. To be fair they sometimes provoke discussions at CARTO mailing list. By the way I’m almost certain the full team has a filter to move them to their archives and they think I’m just an annoying spammer (a couple of times I’ve changed the subject just to troll them xDDD).

To conclude I want to post here my daily Ráfagas experience:

  • Raf is an early bird and sends the digest in the morning, I copy the contents into a shared Google Doc where a group of collaborators help me on translating the content. It may seem not a lot of effort, but doing this every single day needs a team. Really.
  • I go to my favorite text editor, put the translated content into a new file and start a local server to check the website renders properly.
  • If everything is OK I copy the rendered content and send it to CARTO and OSGeo mailing lists
  • I commit and Push to the GitHub repo so the website is updated along with the feed.
  • I archive Raf’s mail from my inbox.

Creating a Ráfaga

That’s it. Raf you are a formidable example of perseverance and I hope you’ll have the energy to keep giving us all those contents for many years. Thanks mate!

Mapping Party en Poio, Pontevedra

El sábado 11 de abril en el Centro Xaime Illa de Raxó (Poio) llega la fiesta de las fiestas, la Mapping Party Poio’15 con el objetivo de pasarlo bien mientras aumentamos los datos geográficos de la costa de Poio con una licencia libre.

Este taller está organizado por la asociación de software libre Xeopesca y cuenta con la colaboración del Concello de Poio y las asociaciones SCD Raxó y ACDeM Armadiña.


  • 10:00-11:00 Presentación de OSM, organización de equipos para la zona a cartografiar.
  • 11:00-14:00 Trabajo de campo con explicaciones de como emplear osmAndroid.
  • 14:00 -16:00 Comida.
  • 16:00-20:00 Trabajar con las computadoras para el volcado de los datos OSM
  • 20:00-20:30 Clausura del curso.


El número de asistentes será de 25. La selección de los candidatos se realizará por orden de inscripción. Se recomienda la disposición de cualquiera de los siguientes dispositivos:  GPS, teléfono con GPS y cámara digital.

Formulario de Inscripción

Para inscribirse a Mapping Party Poio’15 cubre el formulario. (ver aquí) ( .

Material fungible

Se hará entrega a cada uno de los asistentes de un bloc de notas, un bolígrafo y un lápiz.

Redes SociaLes

Establecemos el hashtag #mappingpartypoio para seguir el evento a través de las redes sociales. Además también puedes seguir la  Mapping Party Poio’15 a través de twitter mediante el hashtag #mappingpartypoio o en la  página de facebook de XeoPesca.




The Null Island Algorithm

We geomaticians like to gather around a mythical place called Null Island. This island has everything: airports, train stations, hotels, postcodes, all kinds of shops, a huge lot of geocoded addresses, and whatever geographical feature ends up with null coordinates due to whatever buggy geoprocessing pipeline and ends up in the (0,0) coordinates.

But earlier this year, some geonerds such as @mizmay and @schuyler realized that there is no one Null Island, but one Null Island per datum / coordinate system (depending on who you ask). And @smathermather had the spare time to find out how the “Null Archipielago” looks like:

(Null archipielago image by @smathermather, containing Map tiles by Stamen Design, under CC BY 3.0. Data by OpenStreetMap, under CC BY SA)

Fast forward a few months. I received an e-mail from one of my university peers, asking for help with a puzzle:

A friend of mine received a puzzle with some coordinates. He has to find a place on earth represented by 861126.41, 941711.64.

It’s supposed to be a populated place. Any ideas?

Well, off the top of my head, those looked slightly like UTM coordinates – two digits after the decimal point, suggesting centimeter precision… but the easting is way off the valid range for UTM coordinates.

And I realized this is the Null Archipielago problem, all over again; but instead of plotting (0,0) on a map, let’s plot all points having (861126.41, 941711.64) coordinates in any reference system.

Cue PostGIS. We can create a point in every CRS like so:

select srid, ST_GeomFromText('POINT(861126.41 941711.64)',srid) as geom
 from spatial_ref_sys;

Note the complete absence of PL/SQL in there.

But it will be much easier to work with the data if all the points are in our beloved EPSG:4326 latitude-longitude coordinate system. And while we’re at it, let’s materialize that data into a table:

select srid, ST_Transform(ST_GeomFromText('POINT(861126.41 941711.64)',srid),4326) as geom
 from spatial_ref_sys;

But there is a problem with this – the PostGIS query will crash due to some CRSs having an empty Proj4 string. This took me a while to trace and fix:

select srid, ST_Transform(ST_GeomFromText('POINT(861126.41 941711.64)',srid),4326) as geom
 from spatial_ref_sys where proj4text!='';

And now we can take this data out into a file… but once again, there’s a catch: some of the coordinates are out of bounds and represented by (∞, ∞) coordinate pair. Even though file formats can handle ∞/-∞ values (good thing we know how IEEE floating point format works, right folks?), some mapping software can not accommodate for these values. And I’m looking at you, CartoDB upload page.

In this particular case, there are only points for (∞, ∞) so the data can be cleaned up in just one pass:

delete from archipielago where ST_X(geom)>180;

Then just add a tiny bit of CartoDB magic, and publish a map:


I still don’t know if the original puzzle has anything to do with any obscure used-in-the-real-world CRS, but at least it’s worth a try.

Cómo instalar la toolbox de LAStools en QGIS

El conjunto de herramientas LAStools para procesar datos LiDAR están disponibles en el software propiertario ArcGIS desde abril de 2012, pero no ha sido hasta el FOSS4G13 de Nottingham que no han estado disponibles para QGIS. Las herramientas han sido testadas en QGIS 1.8.0-Lisboa y 2.0.1-Dufour. En esta entrada os enseñamos el modo de installar LAStools en QGIS y, básicamente, es la traducción y adaptación del manual de instalación que se encuentra en el blog de Martin Isenburg aka rapidlasso.

Nos centraremos en instalar las herramientas en la versión de QGIS 2.0.1-Dufour, por lo que lo primero que necesitamos es descargar e installar dicha versión. A partir de entonces, seguiremos las siguientes instrucciones:

  1. Si tienes abierto QGIS, ciérralo

  2. Borra la carpeta donde está contenido el plugin de LiDAR en sextante. La ruta para los usuarios de Windows es:

    "C:\Program Files\QGIS Dufour\apps\qgis\python\plugins\processing\lidar"

    La ruta para los usuarios de GNU/Linux es

  3. Copia la carpeta lidar que está en este archivo comprimido en el lugar de la carpeta borrada.

    Nota: Los usuarios de QGIS 1.8.0-Lisboa tienen que sustituir la carpeta lidar dentro de:

    "C:\Program Files\QGIS Dufour\apps\qgis\python\plugins\sextante\lidar"

    o si utilizas GNU/Linux, dentro de:


    por la que se encuentra en este otro archivo comprimido.

  4. Descarga la versión más reciente de las LAStools que se encuentra en

  5. Extrae la carpeta lastools. Si estás usando Windows, procura no estraerlo en ninguna ruta con espacios, si estás en GNU/Linux tendrás que compilar las herramientas.

  6. Inicia QGIS, y si encuentras algún error cuando se carge el script de Python, repite los pasos de 1 a 3 con más cuidado 😉

  7. Habilita el toolbox que se encuentra bajo el menú procesing como se muestra en la siguiente figura:


  8. Cambia el toolbox de Simplified Interface a Advanced Interface, imágenes de la izquierda y derecha respectivamente:

  9. Abre el submenú Options and configurations de la pestaña Processing como se muestra a continuación:

  10. Activa la casilla Activate que se encuentra en Providers -> Tools for LiDAR data como se muestra en la siguiente figura, e introduce la ruta de la carpeta donde estén alojadas las LAStools en LAStools folder:

  11. Ahora deberías ver el conjunto de herramientas Tools for LiDAR data en el toolbox y todas las LAStools como en esta figura:

  12. Inicia cualquier comando mediante un doble click y rellena las opciones. En la siguiente figura se muestra la interfaz de lasinfo.

Tened en cuenta que, al distrubuirse el código de manera libre sólo de unas pocas herramientas, los usuarios de GNU/Linux tendrán menos funcionalidades disponibles. No así los usuarios de Windows, cuyas herramientas se distribuyen como shareware.

Martin Isenburg agradece a Victor Olaya por crear todo el entorno de sextante para crear nuevos plugins y por los ejemplos de cómo crear módulos. Al igual que él, yo también agradezo a Victor todo el trabajo porque también estoy trasteando con la creación de módulos dentro de sextante 😛

Herramientas Libres para trabajar con datos LiDAR

En los últimos años ha proliferado el uso del LiDAR como técnica topográfica. Básicamente, consiste en un telémetro láser que mide el tiempo que tarda una pulso láser en ir y volver después de haber rebotado en un objeto. De este modo consigue hallar la distancia entre el instrumento y el objeto. Es decir, es sencillamente un distanciómetro, pero con la particularidad de que puede llegar a medir unos 100000 puntos por segundo (100 MHz). Si, además, incorporamos a los equipos de medición un GPS que nos dé la posición y un sistema inercial que nos de la orientación si estamos en movimiento, podemos dar coordenadas globales, normalmente en el sistema WGS84, a todos los puntos medidos. Por tanto, tendremos lo que se denominan nubes de puntos de los cuales conoceremos su posición en un sistema de referencia global, además de otras características relativas al objeto, como la intensidad o diferentes ecos de retornos, o referentes a la medición como ángulo de emisión del pulso, tiempo o distancia relativa al sensor.

Las compañías que desarrollan instrumentación LiDAR, ya sean telémetros aerotransportados, dando lugar a lo que se denomina en inglés Airborne Laser Scanning (ALS), o telémetros terrestes, Terrestrial Laser Scanning (TLS), también desarrollan su propio software destinado a:

  1. Extraer los datos del equipo de medida y ofrecer los datos en algún
    formato propio o, al menos, conocido.
  2. En algunos casos, hacer post-proceso de dichos datos y obtener productos
    cartográficos finales.

Sin embargo, como en el ecosistema de los GIS, existe un software privativo que monopoliza casi todo el mercado del software para realizar los trabajos de post-proceso. No voy a dar el nombre de ninguno de estos paquetes porque el objetivo de esta entrada es justamente la contraria, exponer las posibilidades existentes para utilizar software libre a la hora de trabajar con datos LiDAR.

Lamentablemente, para extraer los datos de la mayoría de los equipos actuales no hay alternativas y necesitamos utilizar obligatoriamente software privativo. Comprar la licencia no es el verdadero problema, porque si tenemos dinero para comprar un equipo que cuesta varias decenas de miles de euros es que también podemos comprar al menos una licencia por un par de miles. Lo peor es lo que verdaderamente implica el software privativo, es decir, que no eres libre de hacer lo que quieras con producto adquirido. Sin embargo, si alguien quiere hacerse un telémetro láser, puede utilizar el módulo para la captura de datos de las librerías PCL, ya que soporta algunos dispositivos conocidos. Pero sobre la librería PCL y otras más hablaremos más datelladamente en otra entrada. Además, en la mayoría de los casos, cuando podamos acceder a datos LiDAR, estos estarán ya en un fomato conocido. Lo lógico es que obtengamos los datos en el formato LAS que es el formato estándar que define la ASPRS (American Society for Photogrammetry and Remote Sensing).

Hay varias librerías libres para la lectura y escritura de archivos en este formato LAS. La decana de ellas es LASlib. Está desarrollada y mantenida por Martin Isenburg y está escrita en C++. Está licenciada bajo LGPL, por lo que se puede utilizar en otros paquetes, aunque sean privativos. Al descargar estas librerías y compilarlas genera unas herramientas llamadas LAStools que sirven para la gestión de archivos LAS (las2las, lasmerge), para creación de LAS a partir de archivos de texto (txt2las) o archivos de texto a partir de LAS (las2txt), para dar información sobre archivos (lasinfo, lasprecision, lasdiff) o para crear un índice espacial de los puntos dentro de los archivos (lasindex). En las últimas versiones, también se crea la herramienta laszip que sirve para comprimir archivo LAS. El formato de salida es LAZ y el archivo de comprimido ocupa SOLO entre el 7% y 20% del tamaño del archivo original. A todo esto hay que añadir que con la librería LASlib también se distribuye otra librería para leer y escribir formatos LAZ, también en C++ y también con licencia LGPL. Para los usuarios de windows, además, están disponibles otras herramientas precompiladas, de las cuales no voy a dar detalles porque no son libres, sino que son shareware.

De la librería LASlib se hizo un fork y nació la librería libLAS, que está bajo el auspicio de OSGeo. También escrita en C++ pero incorporan bindings para una gran cantidad de lenguajes de programación. También incorpora herramientas para la gestión de archivos en formato LAS y texto, que se llaman de igual manera que las LAStools, aunque la utilización de los comandos pueda variar. Las diferencias que existen actualmente entre ambas librerías se pueden encontrar aquí.

La librería SPDLib es bastante reciente. Tanto es así que empezó a desarrollarse en el verano de 2011. El formato estándar que utiliza para trabajar con datos LiDAR se denomina precisamente SPD (Sorted Pulse Data) y está basado en el formato HDF. Es un formato de datos ordenados e indexado que está optimizado para el acceso rápido a los datos y en el que es posible trabajar con toda la señal del pulso de retorno, en inglés full waveform, no sólo con ecos discretos, como hace las librerías anteriores. Y precisamente ésta es una de sus grandes virtudes. Por lo demás, está en consonancia con las librerías ya vistas. Está desarrollada en C++, con bindings para python e IDL, y con licencia GPL. Al compilar se crean una serie de herramientas para la gestión, manipulación e información de archivos SPD, así como una herramienta imprescindible para transformar entre formatos LAS y SPD. Otra de las ventajas que presenta es que incorpora utilidades para generar modelos digitales en cualquier formato soportado por gdal, y para aplicaciones forestales. ¡Todas ellas libres! Sin embargo, al ser desarrollada principalmente por una sola persona y ser tan reciente, uno de los problemas que presenta es la escasa y, en algunos casos, inexacta documentación.

Sin duda, existen más librerías capaces de trabajar con datos LiDAR pero con otros propósitos distintos a los que hemos cubierto aquí. Las intentaremos tratar en otras entradas. Manteneos atentos.

PD: Mientras escribo estas líneas me llega un tuit en el que hablan de otra libreria
en python para leer y escribir datos LiDAR en formatos LAS. Se llaman laspy 🙂