You are here:GeoTux»Geo-Blogs»Librerías para Geomática»Transformaciones afines y de semejanza en R

Estadísticas

Invitados: 20
Usuarios registrados: 2491
Usuarios en línea:
-
Registrados hoy:
-

Registro

Redifusión (RSS)

Blogs y Noticias:
Recibe las actualizaciones en Geo-Noticias y Geo-Blogs

Recibir por e-mail
Recibir Geo-Noticias y Geo-Blogs por e-mail

¿Qué es esto?

Lunes 07 de Mayo de 2012 05:00

Transformaciones afines y de semejanza en R

Written by  German Carrillo
Rate this item
(1 Vote)

En este post doy un ejemplo de cómo usar mi primer paquete de R, llamado 'vec2dtransf', para realizar transformaciones afínes y de semejanza sobre datos vectoriales.

For English click here.

Casos de uso

Las transformaciones afínes y de semejanza son útiles cuando se integran datos espaciales de distintas fuentes. Suele ocurrir que los vectores de un conjunto de datos (llamémoslo 'A') no coinciden con un conjunto de datos base ('B'), el cual puede ser raster o vector. En dicho escenario uno buscaría reposicionar 'A' tomando como referencia 'B'.

Hay varios casos en los que esta situación ocurre, por ejemplo:

  • 'A' no tiene información de su proyección y no hay ninguna pista que la indique.
  • 'A' ya fue proyectado pero simplemente no coincide con 'B' (piensen en datos con datum local cuando se migran a WGS84).
  • 'A' tiene un sistema de coordenadas arbitrario o falso y va a ser integrado con datos de carácter nacional. Este puede ser el caso de proyectos arqueológicos o de ingeniería.
  • 'A' fue digitalizado a partir de una imagen distorsionada, como un mapa antiguo escaneado.
  • 'A' fue digitalizado a partir de una imagen de satélite mal georreferenciada (recientemente ocurrió en OpenStreetMap).

 

El paquete 'vec2dtransf'

El paquete de R 'vec2dtransf' provee clases para definir y realizar transformaciones afínes y de semejanza sobre datos vectoriales, concretamente en objetos 'sp' (objetos creados por el paquete de R 'sp'). Las transformaciones se pueden definir a partir de puntos de control o directamente desde parámetros de transformación. Si se proveen más puntos de control que los requeridos, se ejecutan mínimos cuadrados, permitiendo obtener residuales y RMSE (error medio cuadrático).

Las transformaciones de semejanza permiten rotar, desplazar y escalar geometrías, mientras que las afines permiten rotar, desplazar, escalar (incluso aplicando diferentes factores de escala en cada eje) y torcer geometrías. Por lo menos dos (2) puntos de control son requeridos para la de semejanza y tres (3) para la afín. Ver (Knippers, 2009) para una introducción a las fórmulas de cada transformación.

 

Ejemplo

Como ejemplo, podemos tomar estos dos conjuntos de datos: Departamento de La Guajira (conjunto de datos 'A', fuente: SIGOT) y límite de Colombia (conjunto de datos 'B'). Los datos se incluyen solo para propósitos de demostración. Los datos comparten el sistema de referencia espacial especificado por el código EPSG 3116.

Importémoslos en R a través del paquete rgdal (primero hay que extraer los datos y ajustar la ruta a ellos cuando se requiera):

 

1
2
3
require("rgdal")
l = readOGR(dsn="/path/to/the/data/", layer="La_Guajira")
r = readOGR(dsn="/path/to/the/data/", layer="Colombia")

 

Así es como lucen (La Guajira: línea roja, Límite de Colombia: Polígono verde con borde negro):

 

4
5
plot(l,col="red")
plot(r,add=TRUE,col="lightGreen")

 

 

Podemos también mirarlos de cerca:

 

6
7
plot(r,col="lightGreen",xlim=c(1206309,1264449),ylim=c(1832839,1864827))
plot(l,col="red",add=TRUE)

 

Puede notarse que los datos simplemente no coinciden, a pesar de compartir el mismo sistema de referencia. En este ejemplo llevaremos el límite de La Guajira a coincidir con el límite de Colombia.

Puntos de control

Para reposicionar el conjunto de datos 'A' tomando el 'B' como base, seleccionamos puntos presentes en ambos conjuntos de datos y disponemos las coordenadas de esos puntos en cuatro columnas, así:

X SOURCE        Y SOURCE        X TARGET        Y TARGET

Las coordenadas fuente (source) corresponden al conjunto de datos 'A', es decir, el conjunto de datos a ser transformado, mientras que las coordenadas destino (target) corresponden al conjunto de datos base ('B').

Los puntos de control deben estar bien distribuidos sobre la región a ser reposicionada (Iliffe, pp.135-137). Para este ejemplo, he seleccionado 16 puntos que están incluidos en el paquete 'vec2dtransf' como datos de ejemplo, así que podemos acceder a ellos ejecutando:

 

 8
 9
10
11
require("vec2dtransf")
data(control.points)
options(digits=10) # To see the decimals
control.points

 

Estos puntos de control están bien distribuidos sobre la región a ser transformada:

 

12
13
plot(control.points[2:3],pch=19)
plot(l,add=TRUE,col="red")

 

 

Transformando coordenadas

Primero definimos el objeto transformación afín a partir de los puntos de control (solo pasamos las columnas de las coordenadas):

14
aft = AffineTransformation(control.points[2:5])


Ahora calculamos los parámetros de la transformación:

15
calculateParameters(aft)


Podemos ver los parámetros calculados de esta forma:

16
getParameters(aft)


Y pueso que usamos más puntos de control que los requeridos para una transformación afín, podemos acceder a los residuales y al RMSE (en metros), los cuales resultan de los mínimos cuadrados.

17
18
getResiduals(aft)
getRMSE(aft)

 

Los residuales son las diferencias entre las coordenadas fuente transformadas y las coordenadas destino  de los puntos de control, mientras que el RMSE mide la desviación general de las coordenadas  fuente transformadas con respecto a las coordenadas destino de los puntos de control. Más información sobre residuales y RMSE puede encontrarse en (Iliffe, 2008).

Para ver el efecto global de la transformación, podríamos plotear la transformación de una grilla regular, con base en (UC Davis Soil Resource Laboratory, 2007):

 

19
plotGridTransformation(aft,bbox(l),100)

 

Y agregar el conjunto de datos 'A':

 

20
plot(l,add=TRUE,col="red")

 

 

Finalmente, realizamos la transformación afín sobre el conjunto de datos 'A':

 

21
newLines = applyTransformation(aft,l)

 

Así es como lucen los datos ahora (el conjunto de datos transformado está en color azul):

 

22
23
plot(newLines,col="blue")
plot(r,add=TRUE,col="lightGreen")

 

Y en la vista detallada:

 

24
25
plot(r,col="lightGreen",xlim=c(1206309,1264449),ylim=c(1832839,1864827))
plot(newLines,col="blue",add=TRUE)

 

 

Mucho mejor, ¿cierto? Comparémoslo con el conjunto de datos original (el mal posicionado):

 

26
plot(l,col="red",add=TRUE)

 

 

Usando el paquete rgdal podemos exportar el resultado como, por ejemplo, un Shapefile, en este caso a un directorio temporal:

 

27
writeOGR(newLines, tempdir(), "transformed_Guajira", driver="ESRI Shapefile")

 

Referencias

 

Last modified on Lunes 09 de Julio de 2012 08:24

Comentarios  

 
0 # Traducción correctatuxman 09-07-2012 08:27
Gracias a Geowarrior por recordarme que 'similarity transformation' se traduce al español como 'transformación de semejanza' y no como 'transformación de similaridad'. Ya he actualzado el post con dicha corrección.
Responder | Responder con una citación | Citar
 
 
0 # transform xy into geographic coordinatesnuri 02-05-2013 11:07
thanks for the detailed information. I am trying to use your package. I have 9 xy coordinates and 4 of it has geographical coordinate. I transform those 4 points, but now is it possible to identify the rest of the xy coordinates (5 remaining) into geographic coordinates?
Responder | Responder con una citación | Citar
 
 
0 # Re:tuxman 02-05-2013 15:23
Hi nurl,

please take into account that the approach here exposed deals with cartesian coordinates (i.e., x,y). If you need to convert x,y coordinates into geographical ones, you need to consider a different approach, which, by the way, is nowadays automated in most GIS software. Please have a look at [1] for theoretical aspects on coordinate transformations and to software such as QGIS [2] for the practical part.

Regards,

Tuxman
-----------
[1] tinyurl.com/d6vtsaa
[2] qgis.org/
Responder | Responder con una citación | Citar
 

Escribir un comentario


Código de seguridad
Refescar

 

¿Dónde nos leen?