Análisis de Google Maps

En noviembre de 2010 se promulgó La Ley 26.651 que establece la obligatoriedad de utilizar el mapa bicontinental de la República Argentina.

Solución al uso inminente de Google Maps API

Con argenmap.jquery pretendemos ofrecer la transición a mapas con datos argentinos sin hacer que los desarrolladores de otros sitios se vean en la necesidad de optar por Google Maps en contraste con un servicio local que pueda parecer obsoleto en contraste con el de Google.

Facilitar el contraste entre Google Maps y los datos del IGN

Uno de los principales objetivos es poner a dispocisión de muchos herramientas que permitan cotejar nuestros datos contra los de un sistema tan importante como lo es Google Maps.

Todos usamos Google Maps

Sabemos que por más aplicaciones web personalizadas que uno haga siempre en algún momento los usuarios buscan cotejar nuestros datos con Google.
En el caso del IGN, incluso a veces es una herramienta de ayuda anexa en el momento de la revisión geográfica.
No por el software en sí como interfaz gráfica, sino por la infraestructura de imágenes (satelitarias y vectoriales de calles) que ofrecen . Tanto en casa como en el trabajo siempre están disponibles.

¿Por qué GMaps API y no Openlayers por ejemplo?

Tranquilamente podríamos utilizar OpenLayers o Mapstraction como base para argenmap.jquery pero creemos que es más importante dar una solución concrecta a los problemas que la masificación de Google Maps API presenta a la representación de la geografía nacional.
Argenmap.jquery no tiene por objetivo ofrecer una alternativa a todo el servicio de Google Maps (compuesto Tecnologías WEB + APIs + infraestructura de servidores) en su totalidad.

Para dar una demostración práctica de este problema, tengamos en cuenta la siguiente ubicación expresada en coordenadas WGS84 (latitud, longitud):

(52.4827802220782, -5.625 )

Este punto, casualmente, es la coordenada WGS84 del punto que se encuenta en la esquina inferior izquierda de la tesela 031311 de Bing Maps. Si proyectás esta ubicación a la proyección oficial de Google / Bing Maps (EPSG:3857) utilizando la librería de reproyección Proj.NET, obtenés el siguiente resultado:

(-626172.13571216376, 6853955.508199729)

Podés confirmar esto usando la herramienta online en http://www.synectics-tc.com/resources/downloads/coordinate-reprojection.html que está basada en la librería Proj.Net. (Librería que sin embargo utiliza el código viejo 3785 para referirse a esta proyección).

De cualquier maner,a si transformás este punto utilizando, la última versión (esto es importante) de la herramienta CS2CS de esta manera:

echo -5.625 52.4827802220782 | cs2cs -f “%.10f” +init=epsg:4326 +to +init=epsg:3857

entonces obtenés el siguiente resultado en lugar del anterior:

(-626172.1357121646, 6887893.4928337997)

Observá que, mientras la coordenada X se mantiene igual, la coordenada Y difiere en algo así como 25km entre las dos conversiones.

El inconveniente #1 no necesita explicación realmente. Mientras más masivo se hace el uso de datos geoespaciales, es inevitable que nuevos programadores quieran hacer cosas con datos geoespaciales y mostrarlos sobre un mapa de Google Maps o Bing Maps es obviamente una de las primeras cosas que intentan hacer. Estas herramientas ofrecen interfaces de usuario grandiosas e intuitivas para navegar y analizar datos geoespaciales y siendo que WGS84 es la fuente de datos más prolífica, parece razonable esperar que la conversión WGS84 <->Bing / Google no fuese tan difícil de lograr.

En cuanto a las otras cuestiones...

Confusiones acerca de los códigos EPSG

Los sistemas de coordenadas se definen utilizando un set de parámetros, incluyendo el tipo de proyección utilizada (Mercator, Albers, etc), desplazamiento de las coordenadas (Falso Este, Falso Norte), la unidad de medida (metro, Pie, etc.), el datum subyacente (WGS84, NAD27, etc) y muchos más.

Como es bastante molesto tener que seguir citando este listado completo de parámetros, cada sistema es comunmente identificado por un código (un número entero) único. (también llamado SRID). El set de SRIDs asignados a los sistemas de referencia espaciales de uso común fue originalmente desarrollado por el Grupo de Investigaciones Petroleras Europeas (European Petroleum Survey Group) y en general, estos códigos reciben el nombre de Códigos EPSG.

"EPSG": 900913

Los programadores comenzaron a proyectar sus datos geoespaciales para superponerlos sobre Google o Bing Maps desde que su API lo permite, que fue aproximadamente en el año 2005. Desafortundamente, la gente de EPSG se comportaron despectivamente con respecto al sistema utilizado por Microsoft y Google (más adelante, podés encontrar más información) y se negaron a asignarle un código EPSG oficial. En un comunicado (reportado por Morten Nielsen) se podía leer:

"Hemos revisado el sistema de referencia de coordenadas utilizado por Microsoft, Google, etc y creemos que está fallado, técnicamente hablando. No vamos a devaluar el set de datos EPSG incluyendo semejante geodesia y cartografía inapropiada"

Como resultado de esto, durante los primeros años, los desarrolladores necesitaron una manera alternativa de referirse al set de parámetros de proyección utilizados en los mapas de Google, y el código 900913 (es decir, Google) se hizo popular. Vale la pena notar que este código nunca fue desarrollado, ni soportado, por el EPSG (por eso las comillas en el título de esta sección).

El problema que surgió fue, a falta de un organismo oficial que supervise el estándar, que aparecieron muchas variantes diferentes —todas usando en general los mismos parámetros, pero con diferentes convenciones en su denominación, etc... y estas variantes tenían que ser agregadas en las bases de datos que, de otra manera, solo tenían las listas de los sistemas oficialmente reconocidos por el EPSG.

EPSG: 3785

Al EPSG le llevó un tiempo masticar el sistema utilizado por Microsoft y Google, pero eventualmente, a mediados del 2008 se decidieron a agregar el sistema a su registros de códigos, asignándole el código EPSG:3785. Hubo unas cuatnas alegrías entre quienes participan de la comunidad geoespacial, y la gente comenzó a recodificar el viejo 900913 para utilizar el nuevo código "oficial" en su lugar.

El set de parámetros de proyección correspondiente al EPSG:3785 era el siguiente:

PROJCS["Popular Visualisation CRS / Mercator",
  GEOGCS["Popular Visualisation CRS",
    DATUM["Popular Visualisation Datum",
      SPHEROID["Popular Visualisation Sphere", 6378137, 0, AUTHORITY["EPSG",7059]],
      TOWGS84[0, 0, 0, 0, 0, 0, 0],
      AUTHORITY["EPSG",6055]],
    PRIMEM["Greenwich", 0, AUTHORITY["EPSG", "8901"]],
    UNIT["degree", 0.0174532925199433, AUTHORITY["EPSG", "9102"]],
    AXIS["E", EAST],
    AXIS["N", NORTH],
    AUTHORITY["EPSG",4055]],
  PROJECTION["Mercator"],
  PARAMETER["False_Easting", 0],
  PARAMETER["False_Northing", 0],
  PARAMETER["Central_Meridian", 0],
  PARAMETER["Latitude_of_origin", 0],
  UNIT["metre", 1, AUTHORITY["EPSG", "9001"]],
  AXIS["East", EAST],
  AXIS["North", NORTH],
  AUTHORITY["EPSG",3785]]					
					

(Como dato extra, la proyección era listada en el registro EPSG con una nota que decía "No es un sistema geodésico reconocido").

EPSG: 3857

Sin embargo, menos de un año luego de que el código oficial EPS:3785 fuera presentado, el EPSG realizó un nuevo pedido de cambio (EPSG::2008.114, de Diciembre de 2008), dando de baja el código 3785 y presentando un nuevo código, 3857, con los siguiente parámetros correspondientes:

PROJCS["WGS 84 / Pseudo-Mercator",
  GEOGCS["WGS 84",
    DATUM["WGS_1984",
      SPHEROID["WGS 84", 6378137, 298.257223563, AUTHORITY["EPSG","7030"]],
      AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich", 0, AUTHORITY["EPSG","8901"]],
    UNIT["degree", 0.01745329251994328, AUTHORITY["EPSG","9122"]],
    AUTHORITY["EPSG","4326"]],
  PROJECTION["Mercator_1SP"],
  PARAMETER["central_meridian",0],
  PARAMETER["scale_factor",1],
  PARAMETER["false_easting",0],
  PARAMETER["false_northing",0],
  UNIT["metre", 1, AUTHORITY["EPSG","9001"]],
  AXIS["X",EAST],
  AXIS["Y",NORTH],
  AUTHORITY["EPSG","3857"]]					
					

No es poco usual que los códigos EPSG se actualicen de vez en cuando, pero este cambio fue particularmente raro — primero porque ocurrió casi al instante después de que el primer código fuera presentado, y segundo porque el código de "reemplazo" 3857 era una permutación cercana del código anterior 3785 (¿Coincidencia?)

Aunque uno siempre recuerda los códigos EPSG que usa frecuentamente, generalmente se hace dificultoso recordar si 3857 o 3785 es el nuevo y muchos sitios en la web listan los parámetros de uno bajo el título de otro de manera incorrecta. Además, a pesar del hecho de que el registro EPSG normalmente ofrece detalles de los códigos viejos, el 3785 desapareció misteriosamente de la base de datos por completo, con la única salvedad de una nota que dice "OGP revisó su enfoque a la descripción de CRS de Visualización Popular."...

Una vez resuelto el misterio de cuál código EPSG es cuál, nos vemos en la problemática de entender las discrepancias entre los diferentes resultados de proyección que los códigos anteriores representan.

Los parámetros de proyección definidos para el EPSG:3785 describen una proyección de Mercator de coordenadas geográficas definidas sobre un modelo esférico de la Tierra, de radio de 6.378.137m, como se especifica en la siguiente línea:

SPHEROID["Popular Visualisation Sphere", 6378137, 0, AUTHORITY["EPSG","7059"]]						
						

La proyección EPSG:3857, en contraste, sugiere una proyección Mercator de coordenadas definidas sobre el elipsoide WGS84, con un achatamiento (f) de 1/298,257223563:

SPHEROID["WGS 84", 6378137, 298.257223563, AUTHORITY["EPSG","7030"]]
					

El problema es que ninguna de las dos describe realmente la proyección utilizada por Google y Bing Maps, en la cual las coordenadas geográficas subyacentes están definidas utilizando el elipsoide WGS84 (como en la EPSG:3857) pero proyectadas como si estuvieran definidas sobre una esfera (como en el EPSG:3785).

El problema surge porque no hay una manera estandarizada de representar este sistema de "doble"-elipsoide, en el cual se utiliza un elipsoide para interpretar las coordenadas geográficas diferente al utilizado para proyectarlas.

Diferentes librerías de software manejan esta anomalía de varias maneras, lo que explica la diferencia entre las coordenadas obtenidas al utilizar esta proyección en varias herramientas

Los sistemas que tratan las coordenadas ingresadas como si estuviesen medidas sobre una esfera son incorrectos, pero también los son los que proyectan el elipsoide WGS84 directamente utilizando la proyección de Mercator. De hecho, esta segunda clase de errores explica el aparente "desplazamiento de 20km" en la coordenada Y, que solo se aprecia en los links que nombramos anteriormente. Al ser el elipsoide WGS84 oblongo (como una esfera sobre la cual alguien se sentó) y menos "alto" que una esfera, los objetos geográficos en un mapa Mercator proyectados directamente sobre el elipsoide WGS84 parecerán estar más cercanos al ecuador que aquellos proyectados desde una esfera (es decir, ubicados un tanto a sur en el hemisferio norte o un tanto al norte en el hemisferio sur).

PROJCS["Popular Visualisation CRS / Mercator",
 GEOGCS["Popular Visualisation CRS",
  DATUM["WGS84",
    SPHEROID["WGS84", 6378137.0, 298.257223563, AUTHORITY["EPSG","7059"]],
  AUTHORITY["EPSG","6055"]],
 PRIMEM["Greenwich", 0, AUTHORITY["EPSG", "8901"]],
 UNIT["degree", 0.0174532925199433, AUTHORITY["EPSG", "9102"]],
 AXIS["E", EAST], AXIS["N", NORTH], AUTHORITY["EPSG","4055"]],
PROJECTION["Mercator"],
PARAMETER["semi_minor",6378137],
PARAMETER["False_Easting", 0],
PARAMETER["False_Northing", 0],
PARAMETER["Central_Meridian", 0],
PARAMETER["Latitude_of_origin", 0],
UNIT["metre", 1, AUTHORITY["EPSG", "9001"]],
AXIS["East", EAST], AXIS["North", NORTH],
AUTHORITY["EPSG","3785"]]					
					

double[] WGS84_a_GoogleBing(double lon, double lat) {
	double x = lon * 20037508.34 / 180;
	double y = Math.Log(Math.Tan((90 + lat) * Math.PI / 360)) / (Math.PI / 180);
	y = y * 20037508.34 / 180;
	return new double[] {x, y};
}
double[] GoogleBing_a_WGS84Mercator (double x, double y) {
	double lon = (x / 20037508.34) * 180;
	double lat = (y / 20037508.34) * 180;

	lat = 180/Math.PI * (2 * Math.Atan(Math.Exp(lat * Math.PI / 180)) - Math.PI / 2);
	return new double[] {lon, lat};
}