CODEWITHSUNDEEP
c++gdalogr
René Orozco
René Orozco

Reputation: 11

Reproject XY points of GeoTiff to Lat/Long using GDAL/OGR C++

I'm using the GDAL/OGR API (C++) to reproject XY points into Lat/Long coordinates from a given Tiff.

I'm able to achieve it right now, but, the coordenates that i am getting do not correspond to the actual coordinates i am trying to get.

I realized that the method is not reading the input dataset, so, there is not an "starting point".

In fact, using the gdaltranslate tool from the command line, gives the right corrdinates, because i am sending the Tiff as a parameter. Unfortunately, i cannot afford the use of this command.

Below is my code (so far) and the example of gdaltranslate, as the outputs:

C++ Method:

void reproject_coords(string map_biomass) {
        GDALDataset *dataset;
        GDALAllRegister();
        string ds = map_biomass;
        dataset = (GDALDataset *) GDALOpen(ds.c_str(), GA_ReadOnly);

        OGRSpatialReference source(dataset->GetMetadata()), target;
        OGRCoordinateTransformation *poCT;

        target.importFromEPSG(4326);
        //source.AutoIdentifyEPSG();
        //source.SetWellKnownGeogCS("32643");
        double x = 1794, y = 6644;

        //source.dumpReadable();
        poCT = OGRCreateCoordinateTransformation(&source, &target);

        if( poCT == NULL || !poCT->Transform( 1, &x, &y ) )
            printf( "Transformation failed.\n" );
        else
            printf( "(%f,%f)\n",
                    x, y );
}

Input: x = 1794 y = 6644
Output: (70.527326,0.059926)

Command Line:

gdaltransform Karnataka_biomass.tif -t_srs EPSG:4326 -s_srs EPSG:32643

Input: x = 1794 y = 6644 Output: 75.7461474090405 12.4826172522718

The question is if is there any way to "read" the Tiff in the C++ Method in order to get the real Lat/Long coords. I need?

Upvotes: 1

Views: 2429

Answers (1)

Sergio Monteleone
Sergio Monteleone

Reputation: 2886

You creating the source spatial reference the wrong way. You should use the method GetProjectionRef() not GetMetaData():

OGRSpatialReference source(dataset->GetProjectionRef()), target;

This should initialize the source spatial reference with the same SRS in your GeoTIFF.

If you just need to transform points, there is no need to use the GeoTIFF though. You can just initialize two OGRSpatialReference objects and transform an OGRPoint.

Update #1

Here is an example:

GDALAllRegister();

OGRSpatialReference source, target;

source.importFromEPSG(32643);
target.importFromEPSG(4326);

OGRPoint p;
p.setX(1794);
p.setY(6644);
p.assignSpatialReference(&source);

p.transformTo(&target);

// transformed coordinates
cout << p.getX() << " | " << p.getY();

However, please note that your results are correct, and you will get the same with the code above. The reason you see different values when you look at the transformed GeoTIFF is because gdaltransform will probably need to rotate / rescale it, and since an image is always rectangular no matter what, it may choose a different upper left corner to accomplish this.

Update #2

The OP clarified he doesn't want to project a single point but rather the whole image. This task can be accomplished using the GDALWarp API (a tutorial here).

As an alternative, it is still possible to use the gdaltransform utility, invoking it via a system() call.

Upvotes: 3

Related Questions