Distancia mínima desde los centroides de una cuadrícula a un vectorial de puntos

En un foro que trata temas relacionados con SIGs, un usuario pretende determinar la distancia entre las celdas de una cuadrícula y un shapefile de líneas (hidrografía) y que la distancia mínima encontrada quede incorporada en la tabla atributiva de la cuadrícula. Para resolver este problema tomé, en primera instancia, el uso del comando de GRASS r.distance para el cálculo de las distancias mínimas pero existe un equivalente para su contra parte vectorial: v.distance (para el cual verifiqué su existencia una vez que terminé el código que ocupa este post). La desventaja de usar r.distance es que alarga innecesariamente el código al tener que rasterizar los puntos individuales que corresponden a los centroides de la cuadrícula al barrer el área total.

Para usar el código hay que extraer los nodos del vectorial de línea que corresponde a la hidrografía y rasterizarlo por un valor constante antes de incorporarlo al ambiente de GRASS. Lo mismo hay que hacer con el primer centroide de la cuadrícula (en este caso de 1000×1000). La ejecución del programa permite barrer la cuadrícula de oeste a este y de norte a sur, partiendo del primer centroide, y determina en cada caso la distancia mínima y la almacena en una lista vacía antes de incorporar los valores en la columna distancia de la tabla atributiva de la cuadrícula; mediante el comando db.execute. Esta columna es creada automáticamente si el código detecta que no existe. A continuación el programa (py_distance.py):

#!/usr/bin/env python
#-*using coding: utf-8

from os import system
import grass.script as grass

system("clear")

#Número de filas y columnas de la cuadrícula
rows = 22
cols = 18

#Tamaño de celda
celda = 1000

#Define el arreglo de valores que van a ser introducidos a la tabla atributiva
#En este caso rows*cols

lista = []

#Define el archivo temporal para vectorial punto
tempfile = grass.tempfile()
 
#Loop para el cálculo de las distancias

print "Espere..."

#Coordenadas del centroide de la primera celda de la cuadrícula
x0 = 658201.833231
y0 = 1061243.17669

x = x0
y = y0

#Formato ascii del punto de partida
puntos=grass.read_command("v.out.ascii", 
				                    input="punto", 
				                    format="standard", 
				                    quiet="True")

puntos=puntos.split("\n")

for i in range(rows):
	for j in range(cols):

		#Imprime archivo temporal para el vectorial tipo punto

		pfile=open(tempfile, 'w')  #Abre el archivo temporal		 

		for i in range(len(puntos)):
			if i==11:
				puntos[i] = str(x) + "   " + str(y)
			print >> pfile, str(puntos[i])

		pfile.close() #Cierra el archivo temporal
		
		#Importa el punto creado al formato nativo de GRASS
 
		grass.run_command("v.in.ascii", 
				            format="standard", 
				            input=tempfile, 
				            out="punto2", 
				            quiet="True", 
				            overwrite="True")

		#Convierte el vectorial punto a ráster 

		grass.run_command("v.to.rast", 
							input="punto2", 
							value="1", 
							output="punto_raster", 
							use="val", 
							overwrite="True", 
							quiet="True")
			
		#Remueve el vectorial punto
		grass.run_command("g.remove",
							vect="punto2",
							quiet="True")

		#Determina la distancia mínima al punto creado (desplazado 1000 m con relación al anterior)
		salida = grass.read_command("r.distance",
									maps="punto_raster,annual_precipitation3_raster",
									fs=" ",
									quiet="True")

		salida = list(salida.split(" "))

		#El segundo término de esta sub lista temporal es la distancia mínima
		tmp = float(salida[2])

		lista.append(tmp)
		
		x += celda

	x = x0
	y -= celda

#Se obtiene información sobre las columnas de la cuadrícula
info= grass.read_command("v.info", 
					flags="c", 
					map="cuadricula",
					quiet="True")

info= info.split("\n")

#Crear columna distancia si no existe
val = 0

for strings in info:
	if 'distancia' in strings:
		val = 1
		print "La columna distancia existe"
		break

if val == 0:
	grass.run_command("v.db.addcol", 
					    map="cuadricula", 
					    columns="distancia double precision")

	print "Se creó la columna distancia"

print "Se está escribiendo el arreglo a la tabla..."

#Introduce los valores del arreglo lista en la columna distancia
for i in range(rows*cols):
	input="UPDATE cuadricula SET distancia=" + str(lista[i]) + " WHERE cat = " + str(i+1)
	grass.write_command("db.execute", stdin = input)

print "Listo!"

#Remueve el archivo temporal
system ("rm " + tempfile)

La imagen siguiente refleja las cuatro menores distancias que fueron incorporadas a la tabla atributiva, resaltadas en amarillo, donde los puntos de la “hidrografía” prácticamente coinciden con el centroide de las celdas y permiten corroborar que se ejecutó de manera adecuada.

distancia2

En la imagen que se tiene a continuación se observan los valores de distancia menores que 1000 m (espaciamiento de la cuadrícula), también resaltados en amarillo, y su distribución es una comprobación adicional de la buena ejecución del código.

distancia1

Esta entrada fue publicada en Código Python, GRASS, QGIS, Software Libre. Guarda el enlace permanente.

Una respuesta a Distancia mínima desde los centroides de una cuadrícula a un vectorial de puntos

  1. Pingback: Distancia mínima desde los centroides de una cuadrícula a un vectorial de puntos con v.distance |

Responder

Por favor, inicia sesión con uno de estos métodos para publicar tu comentario:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Cerrar sesión / Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Cerrar sesión / Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Cerrar sesión / Cambiar )

Google+ photo

Estás comentando usando tu cuenta de Google+. Cerrar sesión / Cambiar )

Conectando a %s