¿Geotransformación para la estereografía polar?

Question

Actualmente estoy trabajando para importar los datos climáticos de CANGRID (proporcionados como archivos Surfer Grid ascii, ".grd") a ArcGIS. La cuadrícula tiene un tamaño de 95 filas por 125 columnas. Metadatos proporciona la latitud del origen (esquina inferior izquierda), el tamaño de la celda (50 km), así como la proyección de las notas como estereográfica polar con el meridiano central (110 grados W) y la latitud del origen (60 grados N).

Después de intentar primero convertir el .grd a rásteres como .ascii y .flt sin éxito, he conseguido utilizar GDAL para establecer la extensión y la proyección, sin embargo el conjunto de datos no se alinea correctamente con los límites de la zona prevista. Véase la imagen de abajo.

¿Existe alguna geotransformación aceptada para la estereografía polar que pueda explicar esta falta de alineación?

Por ejemplo, ¿hay algún factor de conversión o rotación específico que deba utilizar?

Un archivo de ejemplo del conjunto de datos está aquí: "t201113.grd"

Este es el código que estoy utilizando actualmente en GDAL

ds = gdal.Open("t201113.grd")
array = ds.ReadAsArray()

x_rotation = 0
y_rotation = 0
xres = 1
yres = -1

llx = -129.8530
lly = 40.0451
ulx = -175.144
uly = 71.385

input_osr = osr.SpatialReference()
input_osr.ImportFromWkt(ds.GetProjection())

wgs84_osr = osr.SpatialReference()
wgs84_osr.ImportFromEPSG(4326)

wgs_to_nps_trans = osr.CoordinateTransformation(wgs84_osr, input_osr)
x, y, z = wgs_to_nps_trans.TransformPoint(llx,lly)

geo_transform = [ x, xres, x_rotation, y, y_rotation, yres ]

ncol = ds.RasterXSize
nrow = ds.RasterYSize

out_driver = gdal.GetDriverByName("HFA")
out_ds = out_driver.Create(t201113.img", ncol, nrow, 1, gdal.GDT_Float32)

out_ds.SetGeoTransform(geo_transform)

out_prj = 'PROJCS["North_Pole_Stereographic",GEOGCS["GCS_WGS_1984",DATUM["WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Stereographic"],PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-110.0],PARAMETER["Scale_Factor",1.0],PARAMETER["Latitude_Of_Origin",60.0],UNIT["50_Kilometers",50000.0]]'

out_ds.SetProjection(out_prj)

out_ds.GetRasterBand(1).WriteArray(array)
out_ds.GetRasterBand(1).SetNoDataValue(1.70141e+038)
out_ds.FlushCache()
out_ds = None
`

Además, aquí está la información de la proyección, tal y como se define en la entrada, es decir, desde "GetProjection()":

'PROJCS["North_Pole_Stereographic",GEOGCS["GCS_WGS_1984",DATUM["WGS_1984",SPHEROID["WGS_1984",6378137.0,298. 257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Stereographic"],PARAMETER["False_Easting",0. 0],PARÁMETRO["Falso_noroeste",0.0],PARÁMETRO["Meridiano_central",0.0],PARÁMETRO["Factor_de_escala",1.0],PARÁMETRO["Latitud_de_origen",90.0],UNIDAD["50_Kilómetros",50000.0]]'

Y la GeoTransformación de entrada:

(-0.5, 1.0, 0.0, 94.5, 0.0, -1.0)

También se proporcionan las coordenadas lat y long de la cuadrícula, y cuando se ven en el sistema de coordenadas proyectadas tienen el aspecto siguiente. Cuando la geotransformación está definida por las coordenadas de la parte inferior izquierda (amarillo) o superior derecha (rosa), puedo establecer efectivamente la extensión, pero sigue habiendo problemas de alineación en todo el raster.

Answer 1

Demasiado largo para un comentario, esto es para acompañar Respuesta de @Matej .

1. Añada los datos ".grd" a ArcGIS.

2. Utilice la función "Raster to Other Format" y convierta su archivo .grd en un formato ESRII GRID. Esto es importante porque la mayoría de las funciones de rasterización en ArcGIS son accesibles sólo para este formato, o eso, o suele ser demasiado lento cuando se utiliza en otros formatos.

3. Como ya tiene el archivo de proyección asociado al expediente. En lugar de proyectar los nuevos datos convertidos, defina su proyección. ArcToolbox > Herramientas de gestión de datos > Proyecciones y transformaciones > Definir proyección. Puedes navegar hasta la proyección estereográfica polar ESRII predefinida y ver si sus parámetros coinciden con los dados en los metadatos (no es así), por lo que puedes modificarla según @Matej. Sólo que aquí - en lugar de modificar, crea una nueva basada en la proyección NPS con el meridiano central y la latitud de origen cambiados y guárdala en el disco, luego navega a la nueva proyección y usa esa cuando definas tu proyección. Esto se debe a que su modificación sobre la marcha no estará disponible más tarde cuando quiera usarla para establecer el sistema de coordenadas para su marco de datos, lo cual es necesario para alinear sus otros archivos de límites con el mismo sistema NPS modificado.

Answer 2

La información de proyección adecuada para estos archivos se puede encontrar aquí: Definir la cadena Proj4 para el conjunto de datos CANGRD

Answer 3

No creo que sea necesario reproyectar la imagen. Sólo hacer lo siguiente:

1. definir(modificar) la proyección del mapa base
2. georeferenciar (desplazar) la imagen a la ubicación especificada

Tenga en cuenta que la imagen (grd) ya está en la proyección estereográfica del Polo Norte, que sólo le da la indicación de cómo ajustar el mapa base que se alineará con la imagen.

Paso 1 :

Modificar la proyección estereográfica original del Polo Norte (WKID: 102018) para ajustar la Latitud de Origen y el Meridiano Central:

Paso 2:

Georreferenciar el archivo grd estableciendo la esquina inferior izquierda en la coordenada especificada (lat,lon). Cuando se actualiza la georreferenciación, el archivo .gdwx se crea en la misma carpeta. Al asignar la esquina SW a (40.0451, -129.853) el contenido del archivo tiene el siguiente aspecto:

50000
0
0
-50000
-1730620.4315
2744092.9724

editar: el archivo world de arriba ha sido modificado manualmente basándose en el tamaño de la celda y la ubicación proporcionada de la esquina SW - la 5ª y 6ª línea representan la ubicación calculada del píxel superior izquierdo de la imagen. La posición de la imagen ha cambiado ligeramente.

Los valores anteriores colocan (desplazan) la imagen en el lugar especificado y definen la escala.

Y este es el resultado:

Si esto no parece estar alineado correctamente, yo cuestionaría las coordenadas proporcionadas para la esquina SW de la imagen. Si tiene acceso a las coordenadas de la esquina NE de la imagen, puede recalcular los parámetros de transformación que escalarían y rotarían la imagen entre dos (o más) puntos.