Si tiene un GeoTiff, ¿sería posible transformar un punto de Lat/Lon en una X, Y utilizando el GeoTransform?

Question

Estoy usando la biblioteca GDAL. Actualmente, puedo tomar un punto superior izquierdo y un punto superior derecho y obtener una imagen del original. Lo que me gustaría hacer ahora es tomar dos puntos WKT y convertir a coordenadas X, Y para hacer lo mismo. Me preguntaba si era posible hacer esto si conocía el GeoTransform y el sistema de coordenadas que usaba (WGS84).

Answer 1

Utilicé la transformación afín para que la imagen calcule algunas latitudes/longitudes de muestra. El único problema que tuve fue que si la imagen está orientada al norte, geoTransform [2] y geTransform [4] deben ponerse a cero al calcular Lat/Lon.

x = (int)Math.Abs(Math.Round((Latitude - geotransform[0])/geotransform[1]));
y = (int)Math.Abs(Math.Round((Longitude - geotransform[3])/geotransform[5]));

Si quería la fuerza bruta, usted podría hacer lo siguiente (lo hice y funcionó, pero esto es sólo pseudocódigo):

//Get the Pixel for the length and width, this portion is for the full image
pixelXSize = AbsoluteValue((latitudeAt(Zero)-(latitudeAt(Length)-1))/imageLength);
pixelYSize = AbsoluteValue((longitudeAt(Zero)-(LongitudeAt(Width)-1))/imageWidth);

//Calculate the x,y conversion for the points you want to calculate
x = AbsoluteValue((latitudeToConvert-latitudeAt(Zero))/pixelXSize);
y = AbsoluteValue((longitudeToConvert-longitudteAt(Zero))/pixelYSize);

Esta respuesta puede ser apagado por uno o dos píxeles Si está utilizando un geotransform, nuevamente las variables North Facing pueden estropear la respuesta que se devuelve. Hasta ahora, solo he probado esto para imágenes orientadas al norte.

Answer 2

Me encontré con esto antes también, y aquí hay una forma bastante agradable de hacer transformaciones de coordenadas.

Nota de GDAL documentation:

El sistema de coordenadas devuelto por GDALDataset :: GetProjectionRef() describe las coordenadas de referencia geográfica que implica la georreferenciación afín transforman devueltos por GDALDataset :: GetGeoTransform().

Podemos usar esto con OGRCoordinateTransformation para hacer la transformación por nosotros.

Básicamente el código se verá algo como esto:

// Load up some dataset. 
dataset = (GDALDataset *) GDALOpen(mapfile, GA_ReadOnly); 

// Define Geographic coordinate system - set it to WGS84. 
OGRSpatialReference *poSRS_Geog = new OGRSpatialReference(); 
poSRS_Geog->importFromEPSG(4326); // WGS84 

// Define Projected coordinate system - set to the GeoTransform. 
const char *sProj = dataset->GetProjectionRef(); 
OGRSpatialReference *poSRS_Proj = new OGRSpatialReference(sProj); 

// Set up the coordinate transform (geographic-to-projected). 
OGRCoordinateTransformation *poCT_Geog2Proj; 
poCT_Geog2Proj = OGRCreateCoordinateTransformation(poSRS_Geog, poSRS_Proj); 

// Now everything is set up and we set transforming coordinates! 
// Pass Lon/Lat coordinates to the Transform function: 
double x = lon; 
double y = lat; 
poCT_Geog2Proj->Transform(1, &x, &y); 

// Now x and y variables will contain the X/Y pixel coordinates. 

Así es como se puede convertir entre longitud/latitud y coordenadas de píxeles. Tenga en cuenta que puede usar matrices con Transform() y convertir varias coordenadas juntas. El primer argumento es el número de pares de coordenadas para transformar, y el segundo y tercer argumentos son punteros a las x e y. Solo transformo un par aquí.

Nota es igualmente fácil de configurar la transformada inversa:

// Set up the coordinate transform (projected-to-geographic). 
OGRCoordinateTransformation *poCT_Proj2Geog; 
poCT_Proj2Geog = OGRCreateCoordinateTransformation(poSRS_Proj, poSRS_Geog); 

Answer 3

utilizo este método:

void transformCoordinatesEPSG(OGRGeometry &geometry,int from, int to) { 
    OGRSpatialReference srcSpatialReference; 
    OGRErr error = srcSpatialReference.importFromEPSG(from); 

    #ifdef __OGRTRANSFORMDEBUG 
     qDebug() << "Import EPSG " << from << "return " << error; 
    #endif 

    OGRSpatialReference dstSpatialReference; 
    error = error | dstSpatialReference.importFromEPSG(to); 

    #ifdef __OGRTRANSFORMDEBUG 
    qDebug() << "Import EPSG " << to << "return " << error; 
    #endif 

    OGRCoordinateTransformation* coordTrans = OGRCreateCoordinateTransformation(&srcSpatialReference, &dstSpatialReference); 
    geometry.transform(coordTrans); 
} 

Para latitud/longitud que debe ser 4326.