pixel/array position to lat long gdal Python

Question

I am trying to convert positions in a raster representing a .tif to the corresponding global coordinates. Converting the whole array to a tif and load it to QGIS everything is referenced fine, but using the below calculation method for single points there is a slight offset (to east-north-east in the resulting coordinates....

raster.tif uses ETRS 89 UTM Zone 32N

Does anyone has an idea?

from osgeo import ogr, gdal, osr
import numpy as np


raster = gdal.Open("rasters/raster.tif")
raster_array = np.array(raster.ReadAsArray())


def pixel2coord(x, y):
     xoff, a, b, yoff, d, e = raster.GetGeoTransform()
     xp = a * x + b * y + xoff
     yp = d * x + e * y + yoff
     return(xp, yp)


print(pixel2cood(500,598))

Anesthetic

Answer 1

I think the issue might be that the xoff and yoff contain coordinates of the upper left corner of the most upper left pixel, and you need to calculate coordinates of the pixel's center.

def pixel2coord(x, y):
    xoff, a, b, yoff, d, e = raster.GetGeoTransform()

    xp = a * x + b * y + a * 0.5 + b * 0.5 + xoff
    yp = d * x + e * y + d * 0.5 + e * 0.5 + yoff
    return(xp, yp)

Scuff

Answer 2

There is no reason to do this manually. You just need to install rasterio, a python library based on GDAL and numpy. Even if the image has an offset/rotation rasterio will take care of it.

You just have to do:

import rasterio

with rasterio.open('rasters/raster.tif') as map_layer:
    coords2pixels = map_layer.index(235059.32,810006.31) #input lon,lat
    pixels2coords = map_layer.xy(500,598)  #input px, py

Both functions return a tuple of the pixels and coordinates respectively.

Scarab

They're looking for global coordinates though, the .xy method returns the coordinates in the particular projection being used. – Shiftless

Also this requires an additional package rasterio, which will require extra steps to retrieve information from a pure GDAL context (rasterio is after all an abstraction layer on top of GDAL and not just another Python binding for it). – Etter