4

Given a GeoTIFF height map image, how can an altitude value be read at a given latitude/longitude using C++ and GDAL?

Ralf Ebert
  • 181
  • 6

1 Answers1

4

At the command line this can be achieved using:

gdallocationinfo -valonly [image] -geoloc lon lat

The same in C++:

GDALRasterImage image(src_image);
int value = image.valueAt(lat, lon);

GDALRasterImage.h:

#pragma once

#include <gdal.h>


class GDALRasterImage {
public:
  GDALRasterImage(const char* filename);
  ~GDALRasterImage();
  int valueAt(double lat, double lon);


private:
  GDALDatasetH dataset;
  GDALRasterBandH band;
  double inverseTransform[6];
};

GDALRasterImage.cpp:

#import "GDALRasterImage.h"

GDALRasterImage::GDALRasterImage(const char* filename) {
  GDALAllRegister();
  this->dataset = GDALOpenEx(filename, GDAL_OF_RASTER, nullptr, nullptr, nullptr);
  assert(dataset != nullptr);


// Assume there is only one band in the raster source and use that
  assert(GDALGetRasterCount(dataset) == 1);
  this->band = GDALGetRasterBand(dataset, 1);


// Get the inverse geo transform to map from geo location -> pixel location
  double datasetTransform[6] = {};
  assert(GDALGetGeoTransform(dataset, datasetTransform) == CE_None);
  assert(GDALInvGeoTransform(datasetTransform, inverseTransform));
}


GDALRasterImage::~GDALRasterImage() {
  GDALClose(dataset);
}


int GDALRasterImage::valueAt(double lat, double lon) {
  int x = static_cast<int>(floor(inverseTransform[0] + inverseTransform[1] * lon + inverseTransform[2] * lat));
  int y = static_cast<int>(floor(inverseTransform[3] + inverseTransform[4] * lon + inverseTransform[5] * lat));


int32_t pixelValue;
  assert(GDALRasterIO(band, GF_Read, x, y, 1, 1, &pixelValue, 1, 1, GDT_Int32, 0, 0) == CE_None);
  return pixelValue;
}

See also:

Ralf Ebert
  • 181
  • 6