Creación de las imágenes georreferenciadas

El resultado final de la selección de puntos en NightCitiesISS y su correspondencia con las coordenadas, permite crear un fichero en formato GeoTIFF o un KMZ para que pueda ser utilizado con el software de gestión GIS. Cada imagen puede ser tratada como una capa que se superpone a la cartografía y sobre la que se pueden llevar a cabo análisis espaciales

El siguiente script permite leer la imagen original sin georrerenciar, conectarse a NightCitiesISS para descargarse las mediciones de los usuarios, y finalmente, generar un fichero que permite la georreferenciación y su visualización usando QGIS o GlobalMapper. Estos son dos programas ampliamente utilizados por la comunidad GIS, el primero de ellos Open Source. Por último, se genera un shell script que utiliza la librería GDAL para crear un KMZ y una imagen GeoTIFF


In [ ]:
import urllib2
import json
import asciitable
import time
import Image


idISS = 'ISS030-E-211378'

dirImagenesISS = 'images/'
dirGeoTIFF = 'geotiff/'
dirPuntosQGIS = 'puntosQGIS/'
dirPuntosGlobal = 'puntosGlobal/'
dirScriptsGDAL = 'scriptsGdal/'

def getKey(item):
	return item[0]

hayProxy = False

# cargo las tareas acabadas de NightCitiesISS para usar todos los puntos disponibles
if hayProxy == True:
	proxy = urllib2.ProxyHandler({'http': 'http://usuario:clave@proxy.empresa.es:8080'})
	auth = urllib2.HTTPBasicAuthHandler()
	opener = urllib2.build_opener(proxy, auth, urllib2.HTTPHandler)
	urllib2.install_opener(opener)
else:
	opener = urllib2.build_opener()	

req = urllib2.Request("http://crowdcrafting.org/api/taskrun?app_id=1712&limit=10000")
f = opener.open(req)
json = json.loads(f.read())

# lista con el conjunto de imágenes a tratar
lista = asciitable.read('../completa.csv') 

puntos = []

#Fichero de puntos para QGIS
gcp = open(dirPuntosQGIS + idISS + ".points", "w")
linea = "mapX,mapY,pixelX,pixelY,enable\n"
gcp.writelines(linea)

#Fichero de puntos para Global Mapper
fglobal = open(dirPuntosGlobal + idISS + ".gcp", "w")

# Fichero para script de GDAL 
scriptGdal = open(dirScriptsGDAL + idISS + ".sh", "w")

#imagenes que si se han reconocido
imgSi = []
punto = 0
for i in range(len(json)):
	if json[i]['info']['LONLAT'] != '':
		
		link = json[i]['info']['img_big'].split('/')
		idiss = link[8].split('.')[0]

		if idiss == idISS:
			
			im=Image.open(dirImagenesISS + idiss + '.jpg')
			dimensiones = im.size # (width,height) tuple
			xMax = dimensiones[0]
			yMax = dimensiones[1]
			#print 'Dimensiones ' + str(dimensiones[0]) + ',' + str(dimensiones[1])

			linea1GDAL = 'gdal_translate -of GTiff  '
			posicionesGDAL = ''

			for k in range(len(json[i]['info']['LONLAT'])):
				
				punto = punto + 1

				x = json[i]['info']['XY'][k].split(';')[0].split(',')[0]
				y = json[i]['info']['XY'][k].split(';')[0].split(',')[1]

				#x = xMax - int(x)
				y = yMax - int(y)

				lineaGCP = json[i]['info']['LONLAT'][k].split(';')[0] + ',' + str(x) + ',' + str(y) + ',' + '1\n'
				#print lineaGCP
				
				gcp.writelines(lineaGCP)

				lineaGlobal = str(x) + ',' + str(y) + ',' + json[i]['info']['LONLAT'][k].split(';')[0] + ',"punto' + str(punto) + '",0\n'  
				fglobal.writelines(lineaGlobal)

				
				# genero el script GDAL

				posicionesGDAL = ' -gcp ' +  str(x) + ' ' + str(y) + ' ' + json[i]['info']['LONLAT'][k].split(';')[0].replace(',', ' ') + ' ' + "origen" "tmpDestino\n"

			linea1GDAL = linea1GDAL + posicionesGDAL + '"' + idiss + 'jpg"' '"tmp/' + idiss + '.jpg"\n'
			linea2GDAL = 'gdalwarp -r near -tps -co COMPRESS=NONE "/tmp/' + idiss + '.jpg" "' + dirGeoTIFF + idiss + '.tif"\n'
			
			scriptGdal.writelines(linea1GDAL)
			scriptGdal.writelines(linea2GDAL)

			
scriptGdal.close()
gcp.close()
fglobal.close()

Como resultado obtenemos el fichero de puntos que relaciona la posición de cada pixel de la imagen con su posición geográfica, listo para usar en QGIS


In [ ]:
mapX,mapY,pixelX,pixelY,enable
-6.29923900000000003,36.53782000000000352,1451,2331,1
-6.27469199999999994,36.53312999999999988,1408,2206,1
-6.30627699999999969,36.52857900000000058,1513,2324,1
-6.26748200000000022,36.48946300000000065,1609,2033,1
-6.16431299999999993,36.5218190000000007,1173,1700,1
-6.22362199999999977,36.56829499999999911,1093,2098,1
-6.11581900000000012,36.66419299999999737,297,1941,1
-6.35837699999999995,36.6163969999999992,1201,2837,1
-6.42721300000000006,36.74621700000000146,706,3584,1
-6.44231900000000035,36.73810100000000034,789,3622,1
-6.08933000000000035,36.2724410000000006,2226,577,1
-6.06687400000000032,36.28812099999999674,2101,522,1
-6.20306599999999975,36.38529900000000339,1965,1422,1
-6.14375699999999991,36.42949999999999733,1590,1323,1
-6.12581799999999976,36.69798999999999722,147,2104,1
-6.13118300000000005,36.70036400000000043,153,2129,1
-6.28443299999999994,36.53378599999999921,1419,2242,1

También obtenemos el mismo resultado para usar con GlobalMapper


In [ ]:
1451,2331,-6.299239,36.537820,"punto1",0
1408,2206,-6.274692,36.533130,"punto2",0
1513,2324,-6.306277,36.528579,"punto3",0
1609,2033,-6.267482,36.489463,"punto4",0
1173,1700,-6.164313,36.521819,"punto5",0
1093,2098,-6.223622,36.568295,"punto6",0
297,1941,-6.115819,36.664193,"punto7",0
1201,2837,-6.358377,36.616397,"punto8",0
706,3584,-6.427213,36.746217,"punto9",0
789,3622,-6.442319,36.738101,"punto10",0
2226,577,-6.089330,36.272441,"punto11",0
2101,522,-6.066874,36.288121,"punto12",0
1965,1422,-6.203066,36.385299,"punto13",0
1590,1323,-6.143757,36.429500,"punto14",0
147,2104,-6.125818,36.697990,"punto15",0
153,2129,-6.131183,36.700364,"punto16",0
1419,2242,-6.284433,36.533786,"punto17",0
1954,2952,10.235532,36.804439,"punto18",0

El shell script quedaría de la siguiente forma


In [ ]:
gdal_translate -of GTiff  "ISS030-E-209446jpg""tmp/ISS030-E-209446.jpg"
gdalwarp -r near -tps -co COMPRESS=NONE "/tmp/ISS030-E-209446.jpg" "geotiff/ISS030-E-209446.tif"
gdal_translate -of GTiff   -gcp 1590 1323 -6.143757 36.429500 origentmpDestino
"ISS030-E-209446jpg""tmp/ISS030-E-209446.jpg"
gdalwarp -r near -tps -co COMPRESS=NONE "/tmp/ISS030-E-209446.jpg" "geotiff/ISS030-E-209446.tif"
gdal_translate -of GTiff   -gcp 1419 2242 -6.284433 36.533786 origentmpDestino
"ISS030-E-209446jpg""tmp/ISS030-E-209446.jpg"
gdalwarp -r near -tps -co COMPRESS=NONE "/tmp/ISS030-E-209446.jpg" "geotiff/ISS030-E-209446.tif"
gdal_translate -of GTiff  "ISS030-E-209446jpg""tmp/ISS030-E-209446.jpg"
gdalwarp -r near -tps -co COMPRESS=NONE "/tmp/ISS030-E-209446.jpg" "geotiff/ISS030-E-209446.tif"
gdal_translate -of GTiff   -gcp 1954 2952 10.235532 36.804439 origentmpDestino
"ISS030-E-209446jpg""tmp/ISS030-E-209446.jpg"
gdalwarp -r near -tps -co COMPRESS=NONE "/tmp/ISS030-E-209446.jpg" "geotiff/ISS030-E-209446.tif"
gdal_translate -of GTiff  "ISS030-E-209446jpg""tmp/ISS030-E-209446.jpg"
gdalwarp -r near -tps -co COMPRESS=NONE "/tmp/ISS030-E-209446.jpg" "geotiff/ISS030-E-209446.tif"