Skip to content
Permalink
Browse files

use native GDAL functions

  • Loading branch information
alexbruy committed Sep 19, 2014
1 parent 12dc59f commit a809dcce141fb89cf04a44c729a678bb43afd28a
Showing with 3 additions and 70 deletions.
  1. +3 −70 python/plugins/processing/tools/raster.py
@@ -52,80 +52,13 @@ def scanraster(layer, progress):


def mapToPixel(mX, mY, geoTransform):
"""Convert map coordinates to pixel coordinates.
@param mX Input map X coordinate (double)
@param mY Input map Y coordinate (double)
@param geoTransform Input geotransform (six doubles)
@return pX, pY Output coordinates (two doubles)
"""

if geoTransform[2] + geoTransform[4] == 0:
pX = (mX - geoTransform[0]) / geoTransform[1]
pY = (mY - geoTransform[3]) / geoTransform[5]
else:
(pX, pY) = applyGeoTransform(mX, mY, invertGeoTransform(geoTransform))
(pX, pY) = gdal.ApplyGeoTransform(
gdal.InvGeoTransform(geoTransform)[1], mX, mY)
return (int(pX), int(pY))


def pixelToMap(pX, pY, geoTransform):
"""Convert pixel coordinates to map coordinates.
@param pX Input pixel X coordinate (double)
@param pY Input pixel Y coordinate (double)
@param geoTransform Input geotransform (six doubles)
@return mX, mY Output coordinates (two doubles)
"""

(mX, mY) = applyGeoTransform(pX + 0.5, pY + 0.5, geoTransform)
return (mX, mY)


def applyGeoTransform(inX, inY, geoTransform):
"""Apply a geotransform to coordinates.
@param inX Input coordinate (double)
@param inY Input coordinate (double)
@param geoTransform Input geotransform (six doubles)
@return outX, outY Output coordinates (two doubles)
"""

outX = geoTransform[0] + inX * geoTransform[1] + inY * geoTransform[2]
outY = geoTransform[3] + inX * geoTransform[4] + inY * geoTransform[5]
return (outX, outY)


def invertGeoTransform(geoTransform):
"""Invert standard 3x2 set of geotransform coefficients.
@param geoTransform Input GeoTransform (six doubles - unaltered)
@return outGeoTransform Output GeoTransform (six doubles - updated)
on success, None if the equation is uninvertable
"""

# We assume a 3rd row that is [1 0 0]
# Compute determinate
det = geoTransform[1] * geoTransform[5] - geoTransform[2] * geoTransform[4]

if abs(det) < 0.000000000000001:
return

invDet = 1.0 / det

# Compute adjoint and divide by determinate
outGeoTransform = [0, 0, 0, 0, 0, 0]
outGeoTransform[1] = geoTransform[5] * invDet
outGeoTransform[4] = -geoTransform[4] * invDet

outGeoTransform[2] = -geoTransform[2] * invDet
outGeoTransform[5] = geoTransform[1] * invDet

outGeoTransform[0] = (geoTransform[2] * geoTransform[3] - geoTransform[0]
* geoTransform[5]) * invDet
outGeoTransform[3] = (-geoTransform[1] * geoTransform[3] + geoTransform[0]
* geoTransform[4]) * invDet

return outGeoTransform
return gdal.ApplyGeoTransform(geoTransform, pX + 0.5, pY + 0.5)


class RasterWriter:

0 comments on commit a809dcc

Please sign in to comment.
You can’t perform that action at this time.