позиция пикселя / массива в лат длинный гдал Python

Question

Я пытаюсь преобразовать позиции в растре, представляющем .tif, в соответствующие глобальные координаты.Конвертируя весь массив в tif и загружая его в QGIS, на все ссылаются нормально, но при использовании приведенного ниже метода расчета для отдельных точек есть небольшое смещение (к востоку-северо-востоку в полученных координатах ....

raster.tif использует ETRS 89 UTM Zone 32N

У кого-нибудь есть идеи?

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))

Answer 1

Нет причин делать это вручную.Вам просто нужно установить растерио, библиотеку python, основанную на GDAL и numpy.Даже если изображение имеет растеризацию смещения / поворота, об этом позаботятся.

Вам просто нужно сделать:

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

Обе функции возвращают кортеж пикселей и координат соответственно.

Answer 2

Я думаю, что проблема может заключаться в том, что xoff и yoff содержат координаты верхнего левого угла самого верхнего левого пикселя, и вам нужно вычислить координаты центра пикселя.

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)