GIS & Chips

OpenLayers y Panoramio

jose — Tue, 04 May 2010 09:59:47 +0000

Tras leer el excelente post de Pepe relativo a la integración de servicios de geo-localización en un ambiente Google Maps llegué a la conclusión de que resulta muy fácil aplicarlo con dicha API, pero mi inquietud me llevó a investigar en la posibilidad de integrarlo en otra API de tipo mapping de reconocido prestigio, como es el caso de OpenLayers, que además es 100% Open Source, lo que nos permite añadir múltiples fuentes de datos (wms, wfs, gml, xml, geojson, georss, conjunto de tiles, proveedores de mapa mundiales -OpenStreetMap, GMaps, Virtual Earth, Yahoo Maps, etc.), además de esa libertad de movimientos que le confieren una cierta ventaja sobre su homónima, sin que ello suponga una crítica hacia Google Maps, que entre otras virtudes ha popularizado el uso de los servicios de mapa vía web, acercándolo a un público genérico sin apenas conocimientos.

Openlayers con Panoramio

A la hora de “copiar” el ejemplo descrito en dicho post para OpenLayer, no he planteado integrar todos los geoservicios, entre otras cosas porque algunos de ellos forman parte de la API casi desde su inicio, y no supone un gran esfuerzo incorporarlo, salvo el de YouTube, que no tuve tiempo de ver, aunque tampoco creo que haya muchos videos geoposicionados. Ejemplos de integración de KML en OL hay muchos en la red, al igual que de GML, GeoJSON o GeoRSS, y también de Flickr, que en el fondo es acceder a un documento XML. Para más información consultas los ejemplos on-line del portal Openlayers.org

De todos los geoservicios, al que más atención ha prestado ha sido al de Panoramio, entre otras cosas por la ingente cantidad de fotografías georreferenciadas que integra, lo que lo hace muy atractivo para el webmaster. Como muchos sabéis este portal dispone de un servicio REST público que nos permite ejecutar una URL con ciertos parámetros, que finalmente nos devuelve un documento estructurado en formato JSON

{
"count": 110,
"photos": [
{
"upload_date": "20 May 2007",
"owner_name": "Carlos Sieiro del Nido",
"photo_id": 2313090,
"longitude": -0.516679,
"height": 75,
"width": 100,
"photo_title": "UNIVERSIDAD DE ALICANTE-Campus de San Vicente del Raspeig.",
"latitude": 38.386611000000002,
"owner_url": "http://www.panoramio.com/user/134716",
"owner_id": 134716,
"photo_file_url": "http://mw2.google.com/mw-panoramio/photos/thumbnail/2313090.jpg",
"photo_url": "http://www.panoramio.com/photo/2313090"
},
{
"upload_date": "16 January 2007",
"owner_name": "\u00a9 www.fotoseb.es - Sebastien Pigneur Jans",
"photo_id": 453228,
"longitude": -0.51438300000000003,
"height": 66,
"width": 100,
"photo_title": "Foto Aerea Universidad Alicante 2 (Foto_Seb)",
"latitude": 38.381700000000002,
"owner_url": "http://www.panoramio.com/user/55833",
"owner_id": 55833,
"photo_file_url": "http://mw2.google.com/mw-panoramio/photos/thumbnail/453228.jpg",
"photo_url": "http://www.panoramio.com/photo/453228"
}]

}

Como podéis ver por cada foto disponemos de un array con los datos que nos interesan, y de forma particular los siguientes:

owner_name: autor
photo_title: título
latitude y longitude: coordenadas geográficas de geoposicionamiento en datum WGS 84
photo_id: identificador de foto
photo_file_url: URL de la imagen
photo_url: URL del recurso en la web de Panoramio

Integrar Panoramio en OpenLayers no resulta tan elegante como en GM, donde sólo hay que escribir 2 líneas, pero con un poco de esfuerzo y el apoyo de la comunidad seguro que alguien es capaz de escribir una clase derivada de JSON que facilite las cosas. Por ahora me limitaré a desmenuzar el código de forma secuencial para que se entienda el proceso

Lo primero que debemos de hacer es preparar nuestro proxy para que admita peticiones al host de Panoramio. En este enlace explica como hacerlo

OpenLayers.ProxyHost= "/cgi-bin/proxy.cgi?url=";

Posteriormente creamos un par de variables con la URL del servicio REST y los parámetros

url = "http://www.panoramio.com/map/get_panoramas.php";
parametros = {
   order:'popularity',
   set:'full',
   from:0,
   to:20,
   minx: minx,
   miny: miny,
   maxx: maxx,
   maxy: maxy,
   size:'thumbnail'
}

Consulta la documentación de la API de Panoramio para más información. Por ahora los parámetros que más nos interesan son:

“to”: indica el nº de fotografías que queremos cargar. Máximo 50
“minx”,”miny”,”maxx”,”maxy”: coordenadas de la BBOX donde buscará fotos. Debén de estar en coordenadas geográficas (EPSG: 4326). En mi caso los he obtenido utilizando este código:

// Obtener extensiones.
ext = map.getExtent();
var minx = ext.left;
var miny = ext.bottom;
var maxx = ext.right;
var maxy = ext.top; //alert (minx + " " + miny + " " + maxx + " " + maxy);

Por último sólo nos queda invocar a un HttpRequest de AJAX que ejecute la URL con los parámetros y los envíe a nuestra función: mostrarfotos

OpenLayers.loadURL(url, parametros, this, mostrarfotos);

Función mostrarfotos
Es la encargada de procesar el documento JSON que recibe como argumento

function mostrarfotos(response) {

El documento JSON, que todavía es una cadena de texto es parseada con OpenLayers.Format.JSON, que la convierte en un objeto que ya podemos tratar:

var json = new OpenLayers.Format.JSON();
var panoramio = json.read(response.responseText);

Creamos un array de featutes con el total de fotos (panoramio.photos.length)

var features = new Array(panoramio.photos.length);

Aplicamos un bucle “for” para recorrer cada fotografía.

for (var i = 0; i < panoramio.photos.length; i++)
{

Declaramos y almacenamos cada dato de la foto en una variable para recorrerla con comodidad.

var upload_date = panoramio.photos[i].upload_date;
var owner_name = panoramio.photos[i].owner_name;
var photo_id = panoramio.photos[i].photo_id;
var longitude =panoramio.photos[i].longitude;
var latitude = panoramio.photos[i].latitude;
var pheight = panoramio.photos[i].height;
var pwidth = panoramio.photos[i].width;
var photo_title = panoramio.photos[i].photo_title;/
var owner_url = panoramio.photos[i].owner_url;
var owner_id = panoramio.photos[i].owner_id;
var photo_file_url = panoramio.photos[i].photo_file_url;
var photo_url = panoramio.photos[i].photo_url;

Opcionalmente podemos recoger todas las fotografías en código HTML para ser ubicadas en el BODY con un DIV, en forma de “galería de fotos”

OpenLayers.Util.getElement('galeria').innerHTML +=
     " ";

Posteriormente en el BODY para representar esta galería sólo hay que escribir:

Creamos un objeto point con las coordenadas de la foto

var fpoint = new OpenLayers.Geometry.Point(longitude,latitude);

Creo un array de atributos

var atributos = {
   'upload_date' : upload_date,
   'owner_name':owner_name,
   'photo_id':photo_id,
   'longitude':longitude,
   'latitude':latitude,
   'pheight':pheight,
   'pwidth':pwidth,
   'pheight':pheight,
   'photo_title':photo_title,
   'owner_url':owner_url,
   'owner_id':owner_id,
   'photo_file_url':photo_file_url,
   'photo_url':photo_url
}

Por cada foto creo un elemento (feature) que contiene su posición y sus atributos

features[i] = new OpenLayers.Feature.Vector(fpoint,atributos);

Terminamos el bucle

Fuera del bucle, todavía en la función mostrar foto nos queda:

Definir el estilo para representar la foto. Tenemos dos opciones
a) Crear un icono único para cada fotografía

// estilo punto
var panoramio_style1 = new OpenLayers.StyleMap(OpenLayers.Util.applyDefaults({
   pointRadius: 15,
   fillColor: "red",
   fillOpacity: 1,
   strokeColor: "black",
   externalGraphic: "${photo_file_url}"
}, OpenLayers.Feature.Vector.style["default"]));

b) Por cada posición añadimos una miniatura de la foto.

var panoramio_style2 = new OpenLayers.StyleMap(OpenLayers.Util.applyDefaults({
   pointRadius: 7,
   fillColor: "red",
   fillOpacity: 1,
   strokeColor: "black",
   externalGraphic: "panoramio-marker.png"
}, OpenLayers.Feature.Vector.style["default"]));

Imagen: Ejemplos de estilos

Estilos para Panoramio

Ahora creamos la vector layer asignándole el estilo escogido

var vectorPano = new OpenLayers.Layer.Vector("Panoramio fotos", {
   styleMap: panoramio_style2
});

Añadimos los elementos:

vectorPano.addFeatures(features);

Añadimos la layer al mapa

map.addLayer(vectorPano);

Definimos para terminar la función el comportamiento del evento que se desencadenará cuando seleccionemos ( click con el ratón) cada fotografía

selectControl = new OpenLayers.Control.SelectFeature(vectorPano,
   {onSelect: onFeatureSelect, onUnselect: onFeatureUnselect});
map.addControl(selectControl);
selectControl.activate();
}

Como ves, se ejecutará la función “onFeatureSelect” cada vez que se selecciona un feature, y “onFeatureUnselect” cuando no hay nada seleccionado.

Con este código tenemos toda la lógica para obtener un vector layer con su posición y atributos. Ahora sólo nos queda definir el comportamiento de los eventos. En este caso vamos a crear un “popup” por cada elemento.

// popups
function onPopupClose(evt) {
   selectControl.unselect(selectedFeature);
}
function onFeatureSelect(feature) {
   selectedFeature = feature;

   // HTML del PopUp
      mensaje = "
"+
	      ""+ feature.attributes.photo_title + "
" +
	      "
" +
	      "autor: "+feature.attributes.owner_name +"";

   popup = new OpenLayers.Popup.FramedCloud("chicken",
      feature.geometry.getBounds().getCenterLonLat(),
      null,
      mensaje,
      null, true, onPopupClose);
   feature.popup = popup;
   map.addPopup(popup);
}
function onFeatureUnselect(feature) {
   map.removePopup(feature.popup);
   feature.popup.destroy();
   feature.popup = null;
}

Ya lo tenemos todo listo. Sólo nos queda verlo en un navegador. ¿A que esperas?

Ver demo en proyección Mercator con fondo de OpenStreetMap
Ver demo en WGS84

Nota de interés:

Por último nos gustaría hacer una consideración de especial relevancia. Este código tiene sentido cuando todas las layers (WMS y vector layer de Panoramio) se encuentran en la misma proyección, en este caso se utiliza el datum WGS84 en coordenadas geográficas (EPSG: 4326). En el caso de que utilices en el mapa otra proyección debes de recordar que has de pasarle a la URL de Panoramio los parámetros de caja en WGS84, y lo mismo ocurre cuando el servicio te devuelve los puntos de cada foto, que ahora debes de hacer lo contrario, es decir, convertirla de EPSG:4326 a la proyección que estés utilizando.

API de Google + OpenStreetMap y otros servicios georreferenciados

pepe — Wed, 24 Mar 2010 11:22:00 +0000

Voy a tratar de diseñar un visor de cartografía utilizando el API de Google Maps ( si quereis saber más sobre este tema, podeis leer artículos en esta misma web que hacen referencia: Google Web Toolkit & Google Maps, Integración de un servicio WMS en GoogleMaps usando GWT ) con la cartografía de OpenStreetMap y además voy a incluir datos de servicios de tipo georeferenciado (panoramio, youtube, flickr, kml, wikipedia,…) en dicho mapa.

Ya aviso que, una vez que muestro todos los servicios en el mapa, queda un mapa demasiado lleno de elementos, pero eso es justo lo que pretendo con este visor, tener todos los elementos en el mismo mapa y que se visualicen tanto fotos, videos, artículos de la wikipedia,kml,… a la vez, permitiéndo incluso añadir muchos más elementos.

En primer lugar, para mostrar el mapa de Google Maps, necesito registrar la clave previamente para el sitio web en cuestion que voy a tratar.

Y una vez registrada la clave ya se puede realizar el visor, la programación del visor sería la siguiente:

//Hojas de estilo de google maps

El visor completo con todos los datos estaría accesible desde la siguiente URL:
http://www.gisandchips.org/demos/pepe/osm/test.html

Como digo, este mapa incluye tantos elementos que a veces es un poco extraño, pero es justo lo que pretendía, incluir muchos servicios en un mismo mapa.

Cartografía temática II: Cálculo de intervalos de clase con PL/R

jose — Thu, 04 Feb 2010 14:53:22 +0000

Población por municipios de Madrid (2008)

Como ya vimos en la primera parte de esta mini-serie dedicada a la cartografía temática, hasta ahora sólo tenemos resuelta una parte, la asignación de paletas de colores para nuestros mapas o semiología cartográfica, pero queda otra parte por resolver, la lógica estadística, es decir, calcular cuales serán los valores de cada uno de los intervalos de clase que intervienen en una distribución numérica.

Como ya es habitual en www.gisandchips.org tenemos una especial predilección por almacenar nuestra fuente de datos en tablas espaciales contenidas en una base de datos PostgreSQL-PostGIS, y esta vez no va a ser menos.

Generalmente, el cálculo de los puntos de ruptura de una distribución ha sido tratado habitualmente desde el punto de vista del cliente (gvSIG, Quantum GIS, ArcGIS, etc), descargándose a la aplicación los datos oportunos (shapefiles, PostGIS, etc) para escoger posteriormente algún método de agrupamiento o clustering. Esta solución resulta apropiada cuando los datos a representar proceden de diferentes fuentes, cambios de geometría constantes realizamos muchos temáticos distintos, o cualquier otra variable a considerar. Pero sin embargo, cuando necesitamos representar la misma realidad espacial (municipios, provincias, comunidades autónomas, continentes, etc) desde diferentes aspectos (demografía, economía, agriculturas, etc) no resulta tán cómodo delegar esta función a la aplicación cliente.

Pongamos un ejemplo para clarificar el tema: Imaginemos que disponemos de datos de población agregados a nivel de municipios, donde existen tantas columnas como censos de población se hayan hecho, y deseamos generar en un portal de Internet mapas temáticos a requerimiento del usuario (año del censo, nº de clases, método de agrupamiento o clustering y paleta de colores). Crear un mapa coroplético de la población lleva implícito un conocimiento de la distribución de los datos (valores mínimos y máximos de la serie, media), y sobre todo las medidas de dispersión o de variabilidad: desviación típica, varianza, co-varianza, coeficiente de correlación de Pearson. Con estos parámetros estamos en condiciones de definir cual será el método de cálculo de intervalo de clase más apropiado para nuestra distribución, para posteriormente asignar la paleta de colores acorde al fenómeno a representar (este último aspecto ya lo tenemos solucionado con nuestra clase PHPcolorBrewer). Con una aplicación de GIS de escritorio sólo podemos realizar un mapa para visualizarlo en pantalla o imprimirlo, pero en un ambiente web(GIS) este forma de actuar carece de lógica, pues estamos limitando al usuario la posibilidad de elegir estos parámetros (método de clustering, número de clases y paleta de colores), ya que es muy complicado obtener datos de un GIS y realizar los cálculos de los intervalos de clases en la própia página web (JavaScript, PHP).

Flujo cartografía temática

En este artículo vamos a solucionar esta cuestión proponiendo un método que nos permita desde el servidor, es decir, desde el gestor de base de datos, calcular los intervalos de clase. Para ello no nos queda más remedio que recurrir a la creación de funciones o procedimientos almacenados que contemplen dicha lógica. Realizarlo con PLPg/SQL es una solución, pero programar estas funciones requiere de mucho esfuerzo al no formar parte de las funciones nativas del servidor (por ejemplo, no existe una función que calcule la desviación típica). ¿Por qué no utilizamos un lenguaje especializado en el manejo estadístico? En este sentido, PL/R reúne todas las condiciones para convertirse en el lenguaje de programación elegido. Una prueba de la potencia de este lenguaje, en este caso R lo podemos ver con este ejemplo.

Implementación en R
Tenemos la población activa (valores en miles) de todas las provincias españolas, obtenidas del INE, en el año 1996, disponibles en un objeto vector:

c(126.4,137.9,570.9,192.1,409.7,60.2,251.3,325.9,2109.1,145.8,159.4,425.5,207.8,194.3,
172.5,301.5,443,69,252.4,300,55.2,296.4,158.2,84.5,244.3,200.3,141.7,174.2,2180.7,
475.8,438,224.6,161.4,71.2,320.8,390,104.7,144.3,328.3,56.5,646.5,36.3,242,52.8,196.5,
919.6,211.6,473.9,70.2,356.7,50)

En una sesión de R se crearía de esta manera

poblacion.1996 <- c(126.4,137.9,570.9,192.1,409.7,60.2,251.3,325.9,2109.1,145.8,159.4,
425.5,207.8,194.3,172.5,301.5,443,69,252.4,300,55.2,296.4,158.2,84.5,244.3,200.3,141.7,
174.2,2180.7,475.8,438,224.6,161.4,71.2,320.8,390,104.7,144.3,328.3,56.5,646.5,36.3,
242,52.8,196.5,919.6,211.6,473.9,70.2,356.7,50)

Y se podría ver simplemente escribiendo:

poblacion.1996

Resultado

[1] 126.4 137.9 570.9 192.1 409.7 60.2 251.3 325.9 2109.1 145.8
[11] 159.4 425.5 207.8 194.3 172.5 301.5 443.0 69.0 252.4 300.0
[21] 55.2 296.4 158.2 84.5 244.3 200.3 141.7 174.2 2180.7 475.8
[31] 438.0 224.6 161.4 71.2 320.8 390.0 104.7 144.3 328.3 56.5
[41] 646.5 36.3 242.0 52.8 196.5 919.6 211.6 473.9 70.2 356.7
[51] 50.0

Para ver un sumario estadístico básico de dicho objeto escribimos:

summary(poblacion.1996)

Con este resultado:

Min. 1st Qu. Median Mean 3rd Qu. Max.
36.3 139.8 207.8 320.8 342.5 2181.0

Entrando ya en materia, vamos a calcular 5 clases distribuidas por el método de intervalos iguales, es decir, (valor máximo- valor mínimo) / nº intervalos:

seq(min(poblacion.1996),max(poblacion.1996),length.out=(5+1))

Resultado

[1] 36.30 465.18 894.06 1322.94 1751.82 2180.70

Es decir, las clases quedarían agrupadas de esta forma:

clase 1: 36.30 – 465.18
clase 2: 465.18 – 894.06
clase 3: 894.06 – 1322.94
clase 4: 1322.94 – 1751.82
clase 5: 1751.82 – 2180.70

Este método, aunque sencillo de implementar (incluso desde la parte cliente) ofrece una clasificación muy distorsionada cuando los valores son muy dispares, por lo que necesitamos recurrir a otros métodos más eficientes, como son los cuantiles. En R este cálculo queda reducido a la mínima expressión:

quantile(poblacion.1996)

Resultado:

0% 25% 50% 75% 100%
36.3 139.8 207.8 342.5 2180.7

Esta sería nuestra clasificación

1º quantil: 36.3 – 139.8
2º quantil: 139.8 – 207.8
3º quantil: 207.8 – 342.5
4º quantil: 342.5 – 2180.7

¡Más sencillo imposible!

El problema de trabajar directamente en R estriba en la dificultad para obtener datos desde otras fuentes, aunque siempre podemos recurrir a paquetes externos para traer los datos de shapefiles o tablas de PostgreSQL a R

Implementación en PostgreSQL

Como ya hemos visto en anteriores artículos sobre PL/R, esta lenguaje procedural es muy sencillo de utilizar. Sólo hay que copiar el código de R y realizar pequeños cambios para adaptarlo a la sintaxis de PL/R.

Un aspecto que debemos tener muy claro es que casi todos los cálculos en R de una distribución deben de estar dispuestos en forma de un vector. Como es lógico, PostgreSQL no soporta este tipo de dato, pero el lenguaje PL/R si que nos permite una conversión hacia el tipo de datos “array”, ofreciéndonos además una función que realiza todo el trabajo (array_accum)

Por ejemplo, si queremos obtener un array con los datos de la poblacion activa en 1996 (columna t1996t1) de la tabla ine.poblacion_activa, esta es la forma de conseguirlo con SQL:

SELECT array_accum(t1996t1) FROM ine.poblacion_activa;

Y este el resultado:

{126.4,137.9,570.9,192.1,409.7,60.2,251.3,325.9,2109.1,145.8,
159.4,425.5,207.8,194.3,172.5,301.5,443,69,252.4,300,55.2,296.4,
158.2,84.5,244.3,200.3,141.7,174.2,2180.7,475.8,438,224.6,161.4,
71.2,320.8,390,104.7,144.3,328.3,56.5,646.5,36.3,242,52.8,196.5,
919.6,211.6,473.9,70.2,356.7,50}

Con esta pequeña aclaración ya podemos utilizar cualquier columna con datos numéricos para realizar funciones de clustering. Las siguientes funciones de PL/R reproducen los cálculos realizados anteriormente con R:

Función para el método de intervalos iguales:

CREATE OR REPLACE FUNCTION r_equal(double precision[], integer)
RETURNS double precision[] AS
$BODY$
sort(seq(min(arg1),max(arg1),length.out=(arg2+1)))

$BODY$
LANGUAGE 'plr' VOLATILE STRICT
COST 100;

Vamos a analizar esta función.

Tiene dos parámetros de entrada:

double precision[]: El array con los datos de entrada que queremos utilizar para calcular los intervalos, utilizando la función array_accum
integer: el nº de intervalos de clase que deseamos obtener

El resultado (RETURNS) será siempre una estructura de tipo array.
El cuerpo de la función se reduciría a una sola línea que recibe los dos argumentos de la función (arg1 y arg2):
sort(seq(min(arg1),max(arg1),length.out=(arg2+1)))

El uso de esta función es muy sencillo:

SELECT r_equal( array_accum(t1996t1), 5 ) FROM ine.poblacion_activa;

Nos devolverá el siguiente array:

{36.3,465.18,894.06,1322.94,1751.82,2180.7}

Para aquellos más versados en PLPg/SQL habrán advertido que realizar la función en este lenguaje requiere de algunas líneas más de código, que pueden ser bastantes más si queremos utilizar otros métodos, como el de cuantiles, que reproducimos a continuación.

CREATE OR REPLACE FUNCTION r_quantile(double precision[], integer)
RETURNS double precision[] AS
$BODY$
quantile(arg1,probs = seq(0, 1, 1/arg2))
$BODY$
LANGUAGE 'plr' VOLATILE STRICT
COST 100;

Ejemplo de uso para 5 intervalos de clase:

SELECT r_quantile(array_accum(t1996t1),5) FROM ine.poblacion_activa;

Devuelve:

{36.3,104.7,174.2,251.3,409.7,2180.7}

Ahora que tenemos implementado un método para calcular los intervalos de clase ya estamos en condiciones de desarrollar un auténtico motor de cartografía temática (corológica) en un ambiente web que combine la capacidad de cálculo de PL/R con la semiología cartográfica implementada en PHPcolorBrewer que sea capaz de “pintar” geometrias de una tabla espacial de PostGIS. El resultado de esta coctelera explosiva lo podeis ver en las siguientes demostraciones on-line:

EJEMPLOS

Test 1: Distribución de la población activa por provincias (valores en miles de personas)
Este sencillo ejemplo nos permite ver las posibilidades de integración de varias piezas de software (PostGIS, PL/R. PHP MapScripts y PHPcolorBrewer). Los cálculos de los intervalos de clase utilizan el método “Kmeans”, y el usuario tiene la posibilidad de cambiar la paleta de colores y el nº de clases, así como el año que desea consultar.

Test 2: Distribución de la tasa de ocupación/desempleo por Comunidades Autónomas(valores en miles de personas)
Similar al anterior, pero aplicado a CC.AA. Ofrece la posibilidad de múltiples consultas (número de ocupados/desempleados diferenciando hombre, mujeres y total)

Test 3: Distribución de la población activa por provincias (valores en miles de personas) según diferentes métodos de clustering
Este ejemplo es similar al primero, pero hemos añadido la posibilidad de que el usuario selecciones el método de cálculo de intervalos de clase desea. Es sin duda, un buen ejemplo de la potencia que tiene PL/R al permitirnos elegir hasta 8 métodos diferentes:

Test 4: Varios mapas temáticos en una sóla página con todos los métodos anteriores
El objeto de este ejemplo es demostrar de una sóla vez como se representa un mismo fenómeno desde diferentes métodos de clustering, para ver cual es el que mejor se adapta a nuestras necesidades.

Test 5: Mapa temático municipal sensible con HTML dinámico
Se trata del ejemplo más completo de todos, que nos permite seleccionar una provincia y ver la distribución de la población (total, hombres o mujeres) por municipio. También permite seleccionar todos los métodos indicados, así como la paleta de colores y nº de intervalos.

Test 6: Todas las paletas PHPcolorBrewer aplicadas a un mapa temático
El objetivo de este ejemplo es demostrar en una sóla página como un mismo fenómeno es representado con todas las paleta de PHPcolorBrewer, para que el usuario seleccione la que más le convenga.

Cartografía temática I: PHPcolorBrewer

jose — Thu, 04 Feb 2010 14:51:30 +0000

La creación de cartografía coroplética es quizás uno de los aspectos más recurrentes que podemos encontrar en un ambiente webgis. Todos hemos visto o realizado alguna vez con nuestro programa de GIS favorito un mapa de distribución de la población o cualquier otra variable por municipios, provincias, autonomías o estados. En su preparación siempre intervienen tres elementos: una base cartográfica con datos asociados para representar, un método estadístico para agrupar los datos en intervalos de clase o conjunto de datos, y finalmente una simbología cartográfica aplicada a dichos grupos que represente de una manera clara el fenómeno que deseamos destacar. Los datos a tratar siempre serán de tipo numérico y harán referencia a datos cuantitativos.

PHPcolorBrewer palettes

NOTA: La segunda parte de este artículo lo puede consultar en: Cartografía temática II: Cálculo de intervalos de clase con Pl/R

Sobre este tipo de cartografía siempre han circulado muchos tópicos del tipo: “la mejor manera de representar los datos es con un mapa”, con una clara alusión al uso de mapas en detrimento de largos listados o tablas cuya lectura resulta menos intuitiva que un mapa con colores degradados. También encontramos declaraciones con una percepción negativa sobre los mapas temáticos, “la mejor forma de ocultar datos es con un mapa”, aludiendo esta vez al hecho incuestionable de que un fenómeno dado puede dar diferentes lecturas en función del mayor o menor acierto en la elección de los intervalos de clase, en el sentido de que muchos intervalos dificultan la lectura y pocos generalizan demasiado los datos.

Cuando trabajamos con aplicaciones de escritorio GIS nuestra labor se reduce a indicar el fenómeno a representar, el número de intervalos y como mucho elegir una de las gamas de paletas, para finalmente obtener una distribución más o menos aceptable donde es raro que el usuario intervenga en el proceso estadístico que implica decidir las rupturas de clase. De igual forma ocurre con la gama de colores, donde solemos aceptar los predefinidos, o como mucho definimos una degradación de colores.

Quizás sea por deformación profesional, pero suelo advertir una cierta dejadez del usuario respecto al tratamiento de los colores que intervienen en una distribución normalizada en mucha de la cartografía coroplética existente, utilizando, en la mayoría de los casos, los valores que por defecto nos proporciona la aplicación. Para cualquiera que haya estudiado cartografía, sabe que la simbolización del color es uno de los temás más complejos, que muchas veces delegamos en la aplicación de escritorio, en la que intervienen factores tales como la representatividad de los colores, legibilidad, idoneidad en función de su uso (por ejemplo en un mapa pensado para su impresión no todos los colores son imprimibles), adecuación del número de intervalos al fenómeno estudiado, capacidad del usuario para distinguir colores, y por supuesto la elección de un sistema de agrupamiento (clustering) adecuado.

La doctora Cynthia A. Brewer, profesora del Departamento de Geografía en la universidad estatal de Pennsylvania (http://www.personal.psu.edu/cab38/) ha dedicado su faceta profesional al estudio del color en la cartografía coroplética, dejando como legado una serie de paletas de colores adecuadas para la representación cartográfica que han sido bautizadas con su nombre y utilizadas en muchas publicaciones de renombre, como el Census Atlas of the United States (http://www.census.gov/population/www/cen2000/censusatlas/) e implementadas en muchas otras aplicaciones de GIS (ArcGIS)

Basándonos en su trabajo hemos elaborado una sencilla clase en PHP que nos permite hacer uso de dichas paletas de una manera sencilla para el usuario, susceptibles de ser utilizadas en la elaboración de cartografía coroplética.

Descargar PHPcolorBrewer

Esta utilidad está sujeta a la Apache License 2.0

Esta clase consta de los siguientes métodos:

listColors: Listar colores de una paleta y un nº de intervalos o clases
getPalettes: Obtiene un array de paletas
getPaletteIntervals: Obtiene un array de intervalos de una paleta. Los valores siempre estarán entre de 3 a 9
getPalIntColors: Obtiene un array de intervalos de una paletas. Siempre será de 3 a 9
getColor: Obtiene un array con el color de una paleta, un intervalo y un número.

Dada la vocación web de esta clase se han creado estos métodos:

comboPalettes: Utilidad HTML para listar paletas en forma de combo
comboNumInt: Utilidad HTML para listado del nº de intervalos de una paleta en forma de combo
rgb2html: Convierte un color en formato RGB a hexadecimal para su uso en web.
tableColor: Crea una tabla HTML dada una paleta y un intervalo (de 3 a 9)

Con el objeto de facilitar al usuario su implementación se han añadido dos scripts con ejemplos de uso:
1. Listado de todas las paletas de colores: examples.palettes.colorBrewer.php
2. Ejemplos de uso: examples.colorBrewer.php

Finalmente, y como ejemplo de cartografía temática les mostramos algunas demostraciones que utilizan esta utilidad, incluidos en el siguiente post: Cartografía temática II: Cálculos de intervalos de clase con PL/R

Extracción de características estructurales de una imagen ( II ) (Semi-variograma de una imagen usando R).

benizar — Thu, 28 Jan 2010 12:45:08 +0000

En este artículo usaremos una librería de R (“fields”), que contiene métodos para calcular el semi-variograma empírico de una imagen, lo cual quiere decir que calcularemos el semi-variograma de una gran cantidad de puntos. Y a partir de ahí, analizaremos la existencia de patrones espaciales en tres parcelas agrícolas arbóreas.

Figura 1: Parcelas estudiadas en este artículo.

Espero que este sea el primero de varios posts donde aprendamos a estudiar las ortoimagenes, extrayendo características de distintos tipos como ya hicimos en un post anterior (http://www.gisandchips.org/2009/11/11/procesamiento-de-imagenes-digitales-con-c-y-una-aplicacion-para-el-analisis-de-parcelas-agricolas/ ). En aquel caso, se extrajeron algunas sencillas frecuencias a partir del cálculo de la Transformada de Hough que nos proporcionaba AForge.NET, lo cual nos servía para automatizar la decisión de si una parcela era una plantación arbórea y realizar un conteo automático del número de árboles que contenía. En dicho post, se ofrece una aplicación llamada RAPID, que incorpora una galería de imágenes de parcelas agrícolas. Hemos extraído tres de ellas para realizar los siguientes análisis. Son las siguientes:

Al estudiar las 3 parcelas que aparecen en la Figura 1 con el RAPID: La primera es identificada totalmente como agrícola y nos permite contar con bastante certidumbre el número de árboles de un modo automático. En el caso de la segunda, su estructura en hileras no permite el conteo de árboles pues están muy juntos. Y por último, en los olivos de la tercera no se cumplía la regla de una estructura normal, aunque sí que podíamos contar los árboles. Este último caso era el único donde RAPID no podía “acertar” pues identificaba la parcela como no agrícola debido a que la diferencia angular entre las direcciones principales de la Transformada de Hough no se aproximaba a 90º.

¿Qué es un variograma experimental? (una breve explicación para entender los resultados)

Los variogramas se utilizan para caracterizar la posible estructura espacial de un conjunto de datos. Podemos distinguir dos tipos de variogramas, el experimental y el modelizado. De estos nos interesa más el experimental. Este último se usa para describir los datos de una muestra, habitualmente una nube de puntos (xyz)

Matemáticamente el semi-variograma, o variograma/2, se puede definir como la mitad del promedio de las diferencias al cuadrado.

Figura 2: Fórmula para calcular el semi-variograma

Donde Np(h) es el número de pares a la distancia h, h es el incremento.

Z(xi) son los valores experimentales.

xi localizaciones donde son medidos los valores z(xi).

Aplicado a imágenes, el semi-variograma (o variograma/2), mide el grado de correlación espacial entre los píxeles de una imagen. Podemos comentar varios conceptos (ver Figura 3):

Sill o meseta: representa la varianza máxima.
Range, rango o alcance: muestra la distancia (en nuestro caso en píxeles) a la que el semi-variograma alcanza la meseta.
Nugget o “efecto pepita”, es la discontinuidad en el origen. Ésta es debida a que el semivariograma, en la práctica, no es nulo en el origen.

Figura 3: Interpretación del semivariograma

Ejemplos de cálculo en R

Os aviso de que las imágenes referidas en el código estan en formato *.ppm aunque R tiene muchas librerías capaces de convertir a este formato. También podéis utilizar OpenOffice.

En este post, aplicamos el siguiente código para calcular el semi-variograma empírico de cada una de estas imágenes:

# Instalo los packages necesarios:

#-------------------------------------------------

install.packages("pixmap", dependencies= T)

install.packages("fields", dependencies= T)

# Cargo las librerias:
#-------------------------------------------------
library(pixmap)

library(fields)

# Hago una lista con las imagenes de las parcelas y las cargo todas:
#-------------------------------------------------
imag_dir <- list.files("C:\\ ... \\Articulos Gisandchips\\ parcelas_ppm\\", full.names=T)

parcela1 <- read.pnm(imag_dir[1])

parcela2 <- read.pnm(imag_dir[2])

parcela3 <- read.pnm(imag_dir[3])

# Aislamos una banda...
#-------------------------------------------------
matriz1<-parcela1@green

matriz2<-parcela2@green

matriz3<-parcela3@green

# Calculamos el semivariograma (si queremos podemos acabar antes cambiando el alcance de 40  a 10, por ejemplo, pero los resultados pueden cambiar y quedará menos bonito (-: )
#-------------------------------------------------
vgram1_40<-vgram.matrix( matriz1, R=40) # esto llevará un poco de tiempo

vgram2_40<-vgram.matrix( matriz2, R=40) # esto llevará un poco de tiempo

vgram3_40<-vgram.matrix( matriz3, R=40) # esto llevará un poco de tiempo

# Añadimos al Plot las matrices que acabamos de calcular
#-------------------------------------------------
plot.vgram.matrix(vgram1_40)# La matriz del variograma

plot.vgram.matrix(vgram2_40)# La matriz del variograma

plot.vgram.matrix(vgram3_40)# La matriz del variograma

# Creamos una curva que ajuste bien sobre la muestra y la añadimos al variograma
#-------------------------------------------------
polyfit1_40_20 <- lm(vgram1_40$vgram ~ poly(vgram1_40$d, 20));

plot(vgram1_40$d, vgram1_40$vgram)

lines(sort(vgram1_40$d), polyfit1_40_20$fit[order(vgram1_40$d)], col=2, lwd=4)

polyfit2_40_20 <- lm(vgram2_40$vgram ~ poly(vgram2_40$d, 20));

plot(vgram2_40$d, vgram2_40$vgram)

lines(sort(vgram2_40$d), polyfit2_40_20$fit[order(vgram2_40$d)], col=2, lwd=4)

polyfit3_40_20 <- lm(vgram3_40$vgram ~ poly(vgram3_40$d, 20));

plot(vgram3_40$d, vgram3_40$vgram)

lines(sort(vgram3_40$d), polyfit3_40_20$fit[order(vgram3_40$d)], col=2, lwd=4)

El plot que obtenemos es el siguiente, o similar si es que hemos preferido cambiar algún parámetro:
Figura 4: Panel donde mostramos los resultados



En la imagen vemos las parcelas, la matriz de su variograma y el variograma, con una curva que ajustamos mucho al los datos. Es fácil interpretar que en la estructura de las parcelas se llega a apreciar cierta “ciclicidad”, ya que encontramos mesetas a distintas distancias. Cada meseta se  corresponde aproximadamente a las hileras de árboles. El alcance es el que hemos definido (40, en este caso 40 píxeles; puede que unos 20 m. en la realidad).
Viendo estas imágenes nos podemos dar cuenta de que es más fácil determinar una regla de clasificación que cuando lo hacíamos en el caso del RAPID. En este caso, nos fijamos en el número de máximos relativos de la función polinómica de ajuste del variograma. A simple vista, la primera parcela tiene unos 4, la segunda apenas 1 y la última 3 o 4. En este caso, podríamos citar una nueva regla para nuestro programa según la que a partir de 2 o 3 máximos relativos una parcela puede ser considerada una plantación arbórea.
Como nos interesa automatizar la tarea creamos una función que nos cuente los máximos relativos. Aquí viene el ejemplo aplicado a la primera parcela:

# Encontrar máximos relativos en la función 1

#--------------------------------------------------

maxRelativos1=0

hMaxRelativos1=0

for(i in 1:(length(polyfit1_40_20$fitted.values)-1))

{

if (polyfit1_40_20$fitted.values[i] > polyfit1_40_20$fitted.values[i+1] && polyfit1_40_20$fitted.values[i] > polyfit1_40_20$fitted.values[i-1])

{maxRelativos1[i]<- polyfit1_40_20$fitted.values[i]

hMaxRelativos1[i]<- vgram1_40$d[i]}

}

MaxRelativos1<-data.frame(hMaxRelativos1,maxRelativos1)

MaxRelativos1<- na.omit(MaxRelativos1)

NumMaxRelat1_40<-length(rownames(MaxRelativos1))

Consultamos los vértices de la curva que hemos dibujado mediante un bucle, de modo que si el vértice i tiene un valor superior al vértice anterior y también al vértice posterior entonces es considerado un máximo relativo y queda almacenado en la lista.
Algunos comentarios sobre todo esto:
(1) El semivariograma implica un gran esfuerzo de cálculo por parte del ordenador. Esto hace que aún habiendo varias librerías en R que obtienen el semivariograma de una nube de puntos, sean significativamente más lentas que Fields. Además, Fields hace directamente lo que necesitábamos. No obstante, en un futuro no muy lejano, intentaré desarrollar una librería en C# para reproducir un análisis de este tipo y otras cosillas relacionadas.
(2) Por cuestiones de tiempo no he tratado de crear funciones para las distintas etapas de la demostración. Puede ser un buen ejercicio tratar de implementar todo esto en funciones que permitan entre sus parámetros especificar la banda con la que trabajar, el alcance del semivariograma…
(3) Existen más posibilidades a la hora de establecer reglas, por ejemplo podríamos explorar la distancia que separa los máximos relativos para distinguir los cultivos más intensivos de los más tradicionales. Se me ocurren muchas más posibilidades.
(4) Por último, aunque no le doy mucha importancia, deciros que el número de máximos relativos no es del todo correcto pues el cero siempre aparece en el conteo (sale 5, 2, 5 y no 4, 1, 4 como he dicho más arriba), esto es porqué no he trabajado bastante el bucle, pero es fácil restar 1 a la lista final.  Mi idea principal es la de explorar el concepto y dar unos ejemplos de código a modo de ideas.
Referencias:
Para entender bien todo esto creo que es interesante ver mi artículo anterior…
Solamente os remito al siguiente tutorial de Surfer. Es bastante didáctico. http://www.goldensoftware.com/variogramTutorial.pdf
De todos modos la red está llena de materiales, apuntes de clase, libros en PDF… No tendréis problemas en encontrar información sobre los variogramas.
——————————————————————-
Si quereis contactar podeis enviarme un email (asunto: gisandchips):
Benito M. Zaragozí
benito.zaragozi@ua.es

Integración de un servicio WMS en GoogleMaps usando GWT

jpiera — Mon, 18 Jan 2010 22:26:00 +0000

Este artículo es una continuación del artículo “Google Web Toolkit & Google Maps” que se publicó en este mismo blog hace un par de meses. La idea de escribirlo surgió debido a una pregunta que nos hizo una lectora en la que preguntaba cómo se podía integrar un servidor WMS propio en Google Maps utilizando GWT.

Este artículo explica paso a paso cómo hacerlo y para ello se ha tomado como servidor de ejemplo un WMS de la NASA, pero se ha escrito el código de forma que haciendo unas pocas modificaciones se puede integrar cualquier otro servidor WMS.

Se parte de un entorno en el que se ha creado un proyecto de GWT llamado GoogleMaps-WMS utilizando para ello el entorno de desarrollo Eclipse. Además se ha creado una página de inicio y se ha añadido un mapa de GoogleMaps en el centro. Cómo crear el proyecto y cómo crear la página de inicio es algo que ya se explicó en el artículo anterior.

Antes de empezar a ver código es necesario entender cómo funciona el sistema de tileado de GoogleMaps basado en un QuadTree, que no es más que una estructura de datos muy utilizada para crear índices espaciales. El QuadTree se aplica sobre un área concreta y en el caso de GoogleMaps se aplica sobre toda la superficie terrestre que varía entre la longitud -180º a 180º y la latitud -90º a 90º (estas cifras no son del todo exactas, ya que la latitud se “acorta” unos grados en los polos)

El QuadTree divide esta superficie en 4 áreas en forma rectangular de tamaño fijo, que en su conjunto pueden ser vistas como una matriz de 2×2 y se les pueden dar nombres según la fila y la columna en la que están. De aquí en adelante a esas áreas rectangulares las vamos a llamar tiles. La nomenclatura de Google es numérica y crece de Oeste a Este y de Norte a sur de forma que el tile de la parte superior izquierda es el (0,0), el tile de la parte superior derecha es el (1,0), el de la parte inferior izquierda es el (0,1) y el de la parte inferior derecha es el (1,1). Esta matriz de 2×2 es la matriz de tiles para un nivel de zoom 1.

Cada una de esos tiles se puede descomponer a su vez en 4 partes iguales de área constante. Si dividimos los 4 tiles que hay inicialmente en 4 partes iguales obtendremos una matriz de 4×4 tiles, que nombraremos matriz de tiles de nivel 2. Podemos utilizar la misma nomenclatura que hemos utilizado para los tiles de nivel de zoom 1, siendo el tile de la esquina superior izquierda el (0,0) y el tile de la esquina inferior derecha el (3,3).

Este proceso de descomposición de tiles en 4 tiles de un nivel inferior lo podemos realizar N veces de forma que al final tendremos una estructura de datos de N niveles en los que en cada nivel habrá 4 veces más tiles que en el nivel predecesor. En el caso de GoogleMaps cada uno de estos tiles es una imagen de 256×256 píxels independientemente del nivel en el que nos encontremos, por lo que a mayor nivel de zoom más detalle tendrá la imagen. Para poder entender mejor la nomenclatura de los tiles de Google podemos consultar ésta página web.

Una vez entendido el sistema de tileado y la nomenclatura que usa GoogleMaps, vamos a empezar a escribir código. Existe en el API que proporciona Google una clase llamada TileLayer que representa un proveedor de tiles del que se nutrirá el mapa para poder visualizarse. Tenemos por tanto que crear una clase llamada WMSTileLayer que herede de TileLayer, lo que implica que hay que implementar algunos métodos:

public class WMSTileLayer extends TileLayer{
public WMSTileLayer(String baseServerUrl, String layers, String styles, Projection projection){
super();
}
@Override
public String getTileURL(Point tile, int zoom) {
return null;
}
@Override
public double getOpacity() {
return 1.0;
}
@Override
public boolean isPng() {
return true;
}

En el método getOpacity hay que devolver la transparencia de la imagen que varía entre 0 y 1. Si queremos que la imagen se visualice en el mapa tendremos que poner el valor 1.0.

En el método isPng hay que indicar si la imagen que se devuelve es png. En nuestro caso hemos devuelto true ya que las imágenes devueltas las pedimos en png.

Pero en el tercer método de la clase es dónde está la complejidad de esta clase. El método getTileURL acepta como parámetros un objeto de tipo Point y un entero. El entero hace referencia el nivel de zoom mientras que el Point es una pareja de enteros que hacen referencia al nombre del tile que hay que devolver. Si hemos entendido la estructura de QuadrTree y la nomenclatura que usa Google en los tiles podremos entender el significado de estos parámetros. A medida que el usuario se desplaza por el mapa, el componente encargado de recuperar los tiles invoca al método getTileURL para recuperar la URL de los tiles del área en la que se encuentra. Este método se invocará normalmente unas cuántas veces por cada vez que el usuario modifique el área del mapa que está visualizando.

En el ejemplo de código hemos visto cómo en la implementación del método getTileURL se devolvía el valor null y ahora lo vamos a modificar para que devuelva la URL del tile que se está pidiendo en los argumentos del método y que tendrá que ser descargado desde un servidor WMS. Podríamos haber creado un componente específico para un servidor WMS concreto, pero en lugar de eso hemos implementado un componente algo más genérico para que se pueda utilizar con otros servidores de mapas.

Para ello hemos añadido a la clase WMSTileLayer algunos parámetros que se establecen en el constructor. Estos parámetros son la dirección del servidor WMS, la lista de capas a cargar y los estilos con los que se desea visualizar la capa. Cómo conocer los valores de la capa(s) que soporta un servidor o los valores del estilo(s) soportados está fuera del alcance de este artículo, pero bastaría con hacer la operación GetCapabilities sobre el servidor.

El código resultante quedaría del siguiente modo:

public class WMSTileLayer extends TileLayer{
private String baseServerUrl;
private String layers;
private String styles;
private Projection projection;
private String getMapUrl = null;
public WMSTileLayer(String baseServerUrl, String layers, String styles, Projection projection){
super(null, 0, 20);
this.baseServerUrl = baseServerUrl;
this.layers = layers;
this.styles = styles;
//construct the base URL
String urlBase = baseServerUrl + "?" +
"FORMAT=image/png&SERVICE=WMS&REQUEST=GetMap&WIDTH" +
"=256&HEIGHT=256&VERSION=1.1.0&" +
"STYLES=" + styles + "&" +
"LAYERS="+ layers + "&";
getMapUrl = urlBase;
}
@Override
public String getTileURL(Point tile, int zoom) {
return null;
}
@Override
public double getOpacity() {
return 1.0;
}
@Override
public boolean isPng() {
return true;
}

Además de los parámetros comentados anteriormente, podemos ver que se ha incluido un cuarto parámetro llamado projection, que se utiliza para identificar la proyección y deberá ser el mismo que el utilizado para crear el componente contenedor de mapas de GoogleMaps.

En el constructor lo que se está haciendo es básicamente crear la URL que servirá de base para lanzar la petición GetMap a un servidor WMS.Todavía faltan algunos parámetros para que la petición sea completa, pero estos parámetros los vamos a añadir en el método getTileURL que es el último método que nos queda por ver para terminar con la explicación.

Para simplificar el problema, vamos a suponer que el servidor WMS que estamos utilizando soporta el sistema de coordenadas en latitud/longitud (EPSG:4326). El problema será entonces pasar un tile de GoogleMaps definido mediante un nombre y un nivel de zoom a un par de coordenadas que definan el mismo tile en latitud/longitud. Para ello utilizamos el API de GoogleMaps y una vez tengamos ese tile, creamos la Url que lo contiene, de modo que el método getTileURL quedará de la siguiente forma:

@Override
public String getTileURL(Point tile, int zoom) {
//Transforma las coordenadas a coordenadas Lat/Lon
Point tileIndexLLPoint = Point.newInstance(tile.getX() * 256, (tile.getY() + 1) * 256);
Point tileIndexURPoint = Point.newInstance(((tile.getX() + 1) * 256), (tile.getY()) * 256);
LatLng llPoint = projection.fromPixelToLatLng(tileIndexLLPoint, zoom, true);
LatLng urPoint = projection.fromPixelToLatLng(tileIndexURPoint, zoom, true);
String url =  getMapUrl + "SRS=EPSG:4326&BBOX=" +
	llPoint.getLongitude() + "," +
	llPoint.getLatitude() + "," +
	urPoint.getLongitude() + "," +
	urPoint.getLatitude();
return url;
}

En las variables llPoint y urPoint se calculan las coordenadas de la esquina inferior izquierda y de la esquina superior derecha del tile que en latitud/longitud equivale al tile de GoogleMaps. Añadiendo estas coordenadas a la URL mediante el parámetro BBOX y añadiendo el parámetro SRS se forma la URL que se devuelve en el método y que será usada por el componente gráfico para recuperar las imágenes y mostrarlas en pantalla. Ahora ya sólo falta crear un objeto WMSTileLayer y asociarlo al componente gráfico de GoogleMaps. La clase que hace esto es GoogleMaps_WMS, y el código es el sigueinte:

public class GoogleMaps_WMS implements EntryPoint {
 private MapWidget map;
 private String wmsURL = "http://wms.jpl.nasa.gov/wms.cgi";
 private String wmsLayers = "global_mosaic_base";
 private String wmsStyles = "pseudo";
 private String mapName = "Nasa";
 public void onModuleLoad() {
//Coordenadas centrales para posicionar el mapa
LatLng latLonCenter = LatLng.newInstance(0,0);
//Creamos el mapa con un tamaÃ±o fijo
map = new MapWidget(latLonCenter, 2);
map.setSize("520px", "520px");
Projection projection = new MercatorProjection(20);
//Creamos el servicio WMS que nutrirá a GoogleMaps
MapType myWMSMap = new MapType(new TileLayer[] {
new WMSTileLayer(wmsURL, wmsLayers, wmsStyles, projection)},
projection,
mapName);
//Añadimos el servicio WMS al mapa
map.addMapType(myWMSMap);
//Añadimos algunas opciones al mapa
MapUIOptions options = MapUIOptions.newInstance(Size.newInstance(400,400));
options.setMapTypeControl(true);
map.setUI(options);
// Add the map to the HTML host page
RootPanel.get("mapContainer").add(map);
}
}

En el ejemplo, la clase tiene algunos atributos:

wmsURL: para definir la URL del servicio WMS.
wmsLayers: lista separada por comas de las capas a mostrar.
wmsStyle: estilo con el que se desea mostrar la capa
mapName: nombre del mapa que aparecerá en el interfaz de usuario de GoogleMaps

Cambiando estos atributos se puede acceder a cualquier otro servidor WMS.

Cómo evitar la inserción de geometrías duplicadas en PostgreSQL

josetomas — Tue, 22 Dec 2009 11:51:02 +0000

Una de las tareas habituales en el tratamiento de inconsistencias cartográficas es la detección y corrección de geometrías duplicadas. En muchos casos eliminarlas es necesario para evitar incoherencias en nuestra base de datos geográfica, sobre todo en cuanto al cómputo de frecuencias, longitudes y superficies. Por otra parte la eliminación a posteriori, también plantea problemas de decisión que pueden complicar el código a desarrollar y desembocar en tiempos dilatados de proceso. De modo que, en ocasiones merece la pena plantearse la detección de geometrías duplicadas a priori. Implementar semejante mecanismo de control del lado del cliente puede resultar relativamente sencillo, pero lo que aquí os propongo es centralizar la lógica y evitar la inserción de geometrías duplicadas desde el lado del servidor PostgreSQL/PostGIS.

Planteado el problema, la primera solución que viene a la mente pasa por un “trigger” que haga la llamada a la correspondiente función de evaluación y lance un mensaje de error en caso necesario. Sin embargo, si lanzar el mensaje de error no es absolutamente necesario y lo que perseguimos es que una transacción de, por ejemplo, 1000 geometrías donde 50 son duplicadas, se complete con los 950 registros válidos insertados, entonces podemos usar una regla o RULE. Cada dos por tres estoy revisitando el sistema de reglas de PostgreSQL y no deja de sorprenderme por su versatilidad y elegancia. Echad un vistazo a la solución propuesta y decidid por vosotros mismos:

CREATE OR REPLACE RULE skip_duplicate AS
 ON INSERT TO table_name
 WHERE 0 <> (SELECT count(*)
 FROM table_name
 WHERE geom_field ~= new.geom_field) DO INSTEAD NOTHING;

Como podéis comprobar, hay 2 claves para interpretar esta regla:

1) Es una regla condicional (CREATE RULE … WHERE … DO …), y en la condición integramos una subquery que emplea el operador de igualdad geométrica (~=) para evaluar si la geometría de entrada (new.geom_field) es un duplicado. Recordad que, para acelerar el proceso, los operadores de PostGIS explotan, si está presente, el índice GIST sobre el campo de geometría.

2) Con INSTEAD NOTHING indicamos al motor de reescritura de expresiones subyacente al sistema de reglas que, cuando se cumpla la condición anterior, simplemente vacíe la expresión entrante y no haga nada. Por tanto, la sentencia INSERT se omite y deja de formar parte de la transacción actual.

¡Eso es todo!

Uso de “dissolve” con POSTGIS

Miguel — Thu, 17 Dec 2009 10:36:51 +0000

Uno de los problemas que tuve al manejar los datos geométricos del INE fue que la capa de municipios constaba de geometría de tipo POLYGON. Lo que necesitaba conseguir era agrupar esa geometría simple y convertirla en MULTIPOLYGON. Eso me dio la idea para el siguiente consejo práctico.

Para ello una simple ojeada a la documentación de POSTGIS y encontré la función que necesitaba, st_union. Esta función nos une la geometría en función de un atributo seleccionado.

Voy a mostrar un ejemplo sencillo, pero muy pedagógico. En este ejemplo partimos de tres polígonos simples, su geometría es de tipo POLYGON y quiero unirlos por un atributo común, en su caso el código.

Vemos en la imagen anterior los 3 polígonos, dos de los cuales tienen el mismo atributo (campo código).

La SQL quedaría de la siguiente manera:

select st_union(geometria) from poligonos group by codigo;

El problema viene con el resultado de esa consulta, y es la heterogeneidad del tipo de geometría:

Para que eso no ocurra, añadimos la función st_multi, y conseguimos una homogeneidad del tipo de geometría, convirtiendo todas las geometrías en MULTIPOLYGON:

select st_multi(st_union(geometria)) from poligonos group by codigo;

En caso de que los polígonos sean adyacentes, pasaría a ser un único polígono, es decir, elimina la frontera entre ambos.

Ahora comprobamos que al usar st_union, simplemente nos ha creado un polígono a partir de un atributo común, manteniendo la geometría de tipo POLYGON y en este caso ya no hace falta usar st_multi.

Google Web Toolkit & Google Maps

jpiera — Fri, 20 Nov 2009 14:40:09 +0000

Google Web Toolkit (GWT) es un framework creado por Google que permite al programador hacer aplicaciones web utilizando el lenguaje de programación Java. Existe un plugin de Eclipse que se utiliza para facilitar la generación de aplicaciones GWT donde se puede depurar código utilizando el depurador de código del entorno. Una vez que la aplicación esta depurada, se tiene que pasar por el compilador GWT que genera código HTML y Javascript. La idea es poder generar páginas que utilicen la tecnología Ajax (HTML + Javascript asíncrono), sin ser necesario escribir ni una sola línea de Javascript.

Google Maps es un servicio de mapas que ofrece Google que tiene un API que permite integrarlo en cualquier página web utilizando Javascript. Pero si estamos utilizando GWT, existe además un jar que contiene un API para Java lo que permite a un programador de GWT poder integrar la tecnología de Google Maps en sus páginas web. Este artículo explica ese API y mediante ejemplos.

En este tutorial vamos a crear un proyecto de GWT y le vamos a añadir un mapa de Google Maps. Partiremos de un entorno en el que se ha instalado un Eclipse y se ha instalado el plugin de GWT. A continuación hay que crear un proyecto del tipo “Google -> Web Application Project” que tiene que tener un nombre (GoogleMapsDemo) y un paquete raíz(org.gisandchips.gwt.maps). No es necesario utilizar los mismos nombres que se usan en el ejemplo, pero si se hace así se podrá seguir este tutorial sin la necesidad modificar nada.

El wizard de creación de un proyecto GWT crea un proyecto que se puede ejecutar: para ello seleccionamos la opción “Run -> Debug Configurations -> Web Application” y seleccionaremos GoogleMapsDemo que es el lanzador que ha creado el wizard por nosotros. Pulsamos en “Debug” y si todo va bien se abrirá un navegador que muestra un formulario con un único campo de texto y un botón para enviar.

Podemos rellenar el campo de texto con cualquier valor y al pulsar el botón se abrirá una ventana dónde el servidor nos da la bienvenida. ¿Y qué tiene esto de especial? Lo que realmente ha ocurrido es que el cliente ha enviado una petición asíncrona (¡el famoso Ajax!) y ha esperado a que el servidor le devuelva una respuesta para mostrarla.

Ahora que ya sabemos crear y ejecutar una aplicación en GWT vamos a crear una aplicación utilizando Google Maps. Primero hay que descargarse la librería de Google Maps para GWT desde la página web de Google (gwt-maps-XXX), hay que descomprimirla, y hay que copiar el jar gwt-maps a la carpeta “war -> WEB-INF -> lib” de nuestro proyecto. A continuación hay que modificar el path del proyecto (“Project -> Properties -> Java Build Path -> Libraries -> Add Jars”) y añadir la librería que acabamos de bajar.

Con esto ya tenemos acceso desde el Eclipse a las clases de Google Maps, pero para tener acceso en tiempo de ejecución tenemos que editar el fichero “org.gisandchips.gwt.maps.GoogleMapsDemo.gwt.xml” y añadir una línea para incluir la dependencia en el proyecto. Al final el fichero quedará de esta forma: