Los sistemas de referencia de coordenadas (CRS) o sistemas de referencia espacial (SRS) son complicados, al menos para mí que no soy geógrafo y tengo un conocimiento superficial de GIS: existen más de 13.000 sistemas diferentes y siempre me cuesta saber en qué sistema están las ubicaciones de un conjunto de datos georreferenciados. Siempre he tenido la intuición de que es algo complejísimo, sin saber exactamente a qué nivel. Solo cuando leí este maravilloso post de 2011 sobre el tema (que por supuesto no hace falta leer para conseguir cambiar de CRS un conjunto de datos), empecé a tomar consciencia de la magnitud de la complejidad: el post explica la incapacidad de las distintas herramientas de transformación para conseguir precisión en las conversiones, dando resultados que pueden diferir en torno a 20km. Apasionante lectura también para personas interesadas en saber curiosidades como de dónde vienen los códigos EPSG, hoy uno de los estándares para sistemas de referencia de coordenadas.
Pero voy a tema. A mí que trabajo casi siempre desarrollando para web me interesa en general que los datos geográficos usen latitud y longitud. Es el sistema de coordenadas que se suele usar en web y en aplicaciones GPS: Google Maps, Bing Maps, OpenStreetMap usan latitud y longitud. Si los datos que recibo usan un sistema cartesiano, tengo que transformarlos.
Detectar el sistema de referencia de coordenadas de un conjunto de datos
Lo primero que hay que hacer es saber en qué sistema de referencia están. Esto puede no ser fácil. Para ello hay que mirar los meta-datos del conjunto de datos y no todos los formatos guardan información sobre el CRS. Cuando no sé en qué sistema están los datos, uso QGIS: abro el archivo con los datos y luego lo exporto a geojson. Este formato guarda meta-datos sobre el CRS. Por ejemplo:
{ "type": "FeatureCollection", "name": "projects", "crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:EPSG::3857" } }, "features": []} |
Otra forma de detectar el CRS de un conjunto de datos es usar el programa de línea de comandos projinfo, que forma parte del paquete de software PROJ. Sin embargo projinfo únicamente detecta el CRS de archivos .prj (información geoespacial de los archivos shape) o .geotiff, hasta donde he podido leer en la documentación. Así que es QGIS el que me suele ayudar en esta tarea.
Cambiar el sistema de referencia de coordenadas de un conjunto de datos
Hay multitud de maneras de cambiar el CRS de un conjunto de datos. Yo busco herramientas que permitan automatizar la tarea dentro de los procesos de limpieza de datos que hago. Para este tipo de procesos, los programitas que integran el paquete de software PROJ me funcionan de maravilla.
PROJ integra librerías para varios lenguajes de programación. Yo he usado dos: cs2cs en la línea de comandos, y pyproj en scripts de Python.
Transformar EPSG:3857 (Web Mercator) a EPSG:4326 (WGS 84) con cs2cs
EPSG:4386 es el sistema que usan la mayoría de cartografías web como OpenStreetMap o Google Maps. Para transformar datos en otro CRS, por ejemplo EPSG:3857, a EPSG:4326 con cs2cs en la línea de comandos:
echo "20037508.34 20037508.34 100" | cs2cs -f "%.6f" EPSG:3857 EPSG:4326 |
El modificador -f permite especificar el formato de salida tras la conversión. Si no se usa -f, cs2cs devuelve por omisión grados, minutos y segundos. Para obtener grados decimales, que es lo que yo necesito para web, hay que añadir -f "%.6f" donde 6 es el número de cifras decimales.
cs2cs devuelve tres coordenadas: x y z. En proyecciones planas la tercera es siempre 0, por supuesto:
85.051129 180.000000 0.000000 |
cs2cs puede transformar un conjunto de coordenadas si se le suministra un archivo con pares por líneas. Por ejemplo, el siguiente listado en EPSG:3857 alojado en el archivo coords.list:
echo coords.list231362.441622 5082306.625994-97756.310090 5108535.370440-674547.189064 4487249.462584-220554.193010 5360852.515295-36260.920048 5255629.363755245097.478426 5071621.614563-455214.864321 4964279.036953-53978.135458 4628117.767552 |
Se puede transformar masivamente de la siguiente manera:
cs2cs -f "%.6f" EPSG:3857 EPSG:4326 coords.list41.472631 2.078364 0.00000041.648932 -0.878160 0.00000037.344603 -6.059560 0.00000043.320300 -1.981272 0.00000042.628723 -0.325737 0.00000041.400673 2.201748 0.00000040.673346 -4.089265 0.00000038.343871 -0.484894 0.000000 |
Transformar EPSG:3857 (Web Mercator) a EPSG:4326 (WGS 84) con pyproj en Python
Para transformaciones de CRS dentro de un proceso de limpieza o estructuración de datos más amplio, quizás tenga más sentido usar una librería como pyproj en Python. La transformación en Python se podría hacer de la siguiente manera, suponiendo que el conjunto de datos esté en formato CSV:
import pandas as pdfrom pyproj import Transformer# define la conversión de EPSG:3857 a EPSG:4326transformer = Transformer.from_crs("EPSG:3857", "EPSG:4326", always_xy=True)# variables para el dataset de origen y el transformadoinput_csv = "input.csv"output_csv = "output.csv"# abre el archivo de origendf = pd.read_csv(input_csv)# variables para las columnas con la información geográficax_col = "longitude"y_col = "latitude"# hace la conversióndf[[x_col, y_col]] = df.apply(lambda row: transformer.transform(row[x_col], row[y_col]), axis=1, result_type="expand")# guarda los datos transformadosdf.to_csv(output_csv, index=False) |