¿Cómo realizo una geotransformación y proyección de 3 parámetros en Sql Server 2012?

10

Tengo una tabla con columnas de latitud longitud (NAD27). Calculo otras dos columnas, X e Y, que representan la ubicación de Web Mercator (WGS84).

Actualmente estoy usando un Arcmap para hacer esto, aplicando la geotransformación recomendada para el área de estudio, la geotransformación de 3 parámetros (geocéntrica) , para pasar de NAD27 a WGS84.

Me gustaría hacer esto completamente dentro de Sql Server 2012. Por lo que puedo decir, Sql Server no admite transformaciones de datos de fábrica. ¿Alguien sabe de una biblioteca SQL que admite esta geotransformación? Me gustaría simplemente usar los mismos coeficientes en SQL que estoy usando actualmente en Arcmap.

También necesito proyectar desde WGS84 lat / long en web mercator. Veo esta fórmula implementada en javascript , pero si alguien tiene un procedimiento almacenado SQL que hace esto, sería genial.

Kirk Kuykendall
fuente
Que yo sepa, no hay una solución OO que funcione en este momento para las transformaciones de datos. La forma más fácil de compilarlo en la base de datos sería usar sharpmap.codeplex.com lib- O tomar el código existente y convertirlo a T-SQL que probé ...
simplexio
@simplexio Gracias, ¿tuvo suerte con la conversión T-SQL?
Kirk Kuykendall
¿Qué tan precisa quieres que sean tus coordenadas convertidas? ¿O la precisión importa tanto?
Mintx
@Mintx Me gustaría reproducir los mismos resultados que obtengo actualmente con Arcmap.
Kirk Kuykendall el
1
Por supuesto. Si puede cambiar db a PostGIS, tiene soporte de re-transformación. El servidor MS SQL puede ser un buen DB y tiene un buen soporte, pero pierdo postgresq cuando hablamos de herramientas
prefabricadas

Respuestas:

5

Con respecto al javascript para SQL, esta es probablemente la forma en que manejarías eso:

SELECT  FromX, 
        FromY, 
        CASE WHEN FromX > 180 THEN NULL ELSE FromX * 0.017453292519943295 * 6378137.0 END AS mercatorX_lon2,
        CASE WHEN FromY > 90 THEN NULL ELSE 3189068.5 * LOG((1.0 + SIN(FromY * 0.017453292519943295)) / (1.0 - SIN(FromY * 0.017453292519943295))) END AS mercatorY_lat2
FROM TABLENAME

Creo que lo siguiente responderá a su primera pregunta. Requerirá bastante verificación de errores. Para ayudar, puede encontrar la ecuación original aquí: http://www.colorado.edu/geography/gcraft/notes/datum/gif/molodens.gif

--fromTheta :column --radians
--fromLamda :column --radians
--fromH     :column --meters

DECLARE @fromA float = 6378206.4        --radius of earth, meters
DECLARE @fromF float =1.0/294.9786982   --Flattening
DECLARE @toA float =6378137.0           --radius of earth, meters
DECLARE @toF float = 1.0/298.257223563  --Flattening
DECLARE @dA float = @toA - @fromA       --change in equatorial radius
DECLARE @dX float = -8.0                --change in X, meters
DECLARE @dY float = 160.0               --change in Y, meters
DECLARE @dZ float = 176.0               --change in Z, meters
DECLARE @dF float = @toF-@fromF         --change in flattening
DECLARE @fromES float = 2.0*@fromF - @fromF*@fromF --first eccentricity squared
DECLARE @bda float = 1.0-@fromF         --polar radius divided by equatorial radius

--RM = (@fromA*(1-@fromES)/POWER(1-@fromES*sin(fromTheta)*sin(fromTheta), 1.5))

--RN = (@fromA/SQRT(1.00-@fromES*sin(fromTheta)*sin(fromTheta)))

SELECT 

((((-@dX*sin(fromTheta)*cos(fromLamda)-@dY*sin(fromTheta)*sin(fromLamda))+@dZ*cos(fromTheta))+@dA*(@fromA/SQRT(1.00-@fromES*sin(fromTheta)*sin(fromTheta)))*@fromES*sin(fromTheta)*cos(fromTheta)/@fromA)+@df*((@fromA*(1-@fromES)/POWER(1-@fromES*sin(fromTheta)*sin(fromTheta), 1.5))/@bda+(@fromA/SQRT(1.00-@fromES*sin(fromTheta)*sin(fromTheta)))*@bda)*sin(fromTheta)*cos(fromTheta))/((@fromA*(1-@fromES)/POWER(1-@fromES*sin(fromTheta)*sin(fromTheta), 1.5)) + fromH) AS deltaTheta,
(-@dX*sin(fromLamda)+@dY*cos(fromLamda))/((((@fromA/SQRT(1.00-@fromES*sin(fromTheta)*sin(fromTheta))) +fromH) * cos(fromTheta)) AS deltaLamda,
@dX*cos(fromTheta)*cos(fromLamda)+@dY*cos(fromTheta)*sin(fromLamda)+@dZ*sin(fromTheta)-@da*@fromA/(@fromA/SQRT(1.00-@fromES*sin(fromTheta)*sin(fromTheta)))+@dF*@bda*(@fromA/SQRT(1.00-@fromES*sin(fromTheta)*sin(fromTheta)))*sin(fromTheta)*sin(fromTheta) AS deltaH

FROM TABLENAME

Editar: un par de variables que deberían haber sido nombres de columna, y una coma y paréntesis faltantes.

Editar: un paréntesis más.

He probado esta fórmula y funciona utilizando puntos aleatorios contra la transformación de ArcGIS. Recuerde que sus unidades pueden estar en pies / grados. Recuerde también que estos resultados son deltas, por lo que deberá agregarlos a sus valores para obtener sus resultados finales.

ike
fuente
1
Gracias, creo que los deltas XYZ deben aplicarse después de la conversión de lat, largo en espacio XYZ donde el origen de los ejes XY y Z está en el centro de la tierra.
Kirk Kuykendall
Voy a imprimir ese gif y enmarcarlo en la pared frente a mi escritorio.
nickves
@KirkKuykendall Este método es el Molodensky abreviado, donde los deltas que obtienes están en segundos de arco y se pueden aplicar a tus lat / long iniciales para obtener la traducción a tu dato objetivo. No conozco su AOI, pero geocéntrico es generalmente la forma menos precisa (¡pero más fácil!) De obtener de NAD27-> WGS84.
Mintx
También tenga en cuenta los @dX @dY @dZvalores de ike que pueden ser diferentes dependiendo del NAD_1927_To_WGS_1984método geocéntrico que haya elegido.
Mintx