3

I'd like to create a little server for us, that tells us the elevation of a point in a GeoTIFF file. I already got it working basically using the following:

dataset = gdal.Open('test_file.tif', gdal.GA_ReadOnly)
geotransform = dataset.GetGeoTransform()
px, py = convertLatLngToPixels(geotransform, lat, lng)
val = self.dataset.GetRasterBand(1).ReadAsArray(px, py, 1, 1)

Now the issues is, that our real raster file, is a VRT that combines several GeoTIFFs. Using my code now on that VRT does not work any more...

Does anybody know any solution to using a VRT file?

Improve this question
4
  • Study the source code of gdal.org/programs/gdallocationinfo.html#gdallocationinfo or use that as is.
    – user30184
    Jul 5, 2019 at 8:53
  • But I'd need it for python... But i guess i could use the result of the lifonly option and then continue inside my code? But the issue will be that gdallocationinfo -lifonly is not really fast... Is there a fast way to find the file that I'd need to use inside python?
    – Georg
    Jul 5, 2019 at 9:44
  • 3
    What do you mean "does not work anymore"? Throws an error? Gives the wrong value? ... I don't see why it wouldn't work on a vrt as opposed to a tif.
    – Jon
    Jul 19, 2019 at 15:20
  • Borrow code from github.com/OSGeo/gdal/blob/master/gdal/swig/python/gdal-utils/… maybe?
    – ThomasG77
    Oct 16, 2021 at 23:11

2 Answers 2

Reset to default
0

As it was said in comments, you can find the answer in gdallocationinfo sources. It does the following thing:

In [37]: band = ds.GetRasterBand(i)

In [38]: band.GetMetadataItem('Pixel_1000_1000', 'LocationInfo')                                                                                                                                               
Out[38]: '<LocationInfo><File>/files/Files/GIS/raster_nsk/test.tif</File></LocationInfo>'

In [39]: band.GetMetadataItem('Pixel_992644_2000', 'LocationInfo')                                                                                                                                             
Out[39]: '<LocationInfo><File>/files/Files/GIS/raster_nsk/abc1.tif</File><File>/files/Files/GIS/raster_nsk/abc2.tif</File></LocationInfo>'
Improve this answer
1
  • The questioner wants pixel values, your answer seems to show something else.
    – bugmenot123
    Jul 22, 2019 at 14:50
0

Maybe you find an answer for your question here: Profile Tool has too high resource consumption: Alternative or work around?

In my solution I'm using a VRT file that combines ~5000 8MB TIF files to a virtual DEM mosaic for building a terrain profile, and I'm really impressed about the processing speed of Python GDAL.

Improve this answer

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.