The following module might help, use the mapToPixel() function. #---Imports try: from osgeo import gdal except ImportError: import gdal
def mapToPixel(mx,my,gt): ''' Convert map to pixel coordinates @param mx Input map x coordinate (double) @param my Input map y coordinate (double) @param gt Input geotransform (six doubles) @return px,py Output coordinates (two doubles) ''' if gt[2]+gt[4]==0: #Simple calc, no inversion required px = (mx - gt[0]) / gt[1] py = (my - gt[3]) / gt[5] else: px,py=ApplyGeoTransform(mx,my,InvGeoTransform(gt)) return int(px+0.5),int(py+0.5) def pixelToMap(px,py,gt): ''' Convert pixel to map coordinates @param px Input pixel x coordinate (double) @param py Input pixel y coordinate (double) @param gt Input geotransform (six doubles) @return mx,my Output coordinates (two doubles) ''' mx,my=ApplyGeoTransform(px,py,gt) return mx,my def ApplyGeoTransform(inx,iny,gt): ''' Apply a geotransform @param inx Input x coordinate (double) @param iny Input y coordinate (double) @param gt Input geotransform (six doubles) @return outx,outy Output coordinates (two doubles) ''' outx = gt[0] + inx*gt[1] + iny*gt[2] outy = gt[3] + inx*gt[4] + iny*gt[5] return (outx,outy) def InvGeoTransform(gt_in): ''' ************************************************************************ * InvGeoTransform(gt_in) ************************************************************************ ** * Invert Geotransform. * * This function will invert a standard 3x2 set of GeoTransform coefficients. * * @param gt_in Input geotransform (six doubles - unaltered). * @return gt_out Output geotransform (six doubles - updated) on success, * None if the equation is uninvertable. ''' # ************************************************************************ ****** # * This code ported from GDALInvGeoTransform() in gdaltransformer.cpp # * as it isn't exposed in the python SWIG bindings until GDAL 1.7 # * copyright & permission notices included below as per conditions. # # ************************************************************************ ****** # * $Id: gdaltransformer.cpp 15024 2008-07-24 19:25:06Z rouault $ # * # * Project: Mapinfo Image Warper # * Purpose: Implementation of one or more GDALTrasformerFunc types, including # * the GenImgProj (general image reprojector) transformer. # * Author: Frank Warmerdam, warmer...@pobox.com # * # ************************************************************************ ****** # * Copyright (c) 2002, i3 - information integration and imaging # * Fort Collin, CO # * # * Permission is hereby granted, free of charge, to any person obtaining a # * copy of this software and associated documentation files (the "Software"), # * to deal in the Software without restriction, including without limitation # * the rights to use, copy, modify, merge, publish, distribute, sublicense, # * and/or sell copies of the Software, and to permit persons to whom the # * Software is furnished to do so, subject to the following conditions: # * # * The above copyright notice and this permission notice shall be included # * in all copies or substantial portions of the Software. # * # * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS # * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, # * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL # * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER # * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING # * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER # * DEALINGS IN THE SOFTWARE. # ************************************************************************ **** # we assume a 3rd row that is [1 0 0] # Compute determinate det = gt_in[1] * gt_in[5] - gt_in[2] * gt_in[4] if( abs(det) < 0.000000000000001 ): return inv_det = 1.0 / det # compute adjoint, and divide by determinate gt_out = [0,0,0,0,0,0] gt_out[1] = gt_in[5] * inv_det gt_out[4] = -gt_in[4] * inv_det gt_out[2] = -gt_in[2] * inv_det gt_out[5] = gt_in[1] * inv_det gt_out[0] = ( gt_in[2] * gt_in[3] - gt_in[0] * gt_in[5]) * inv_det gt_out[3] = (-gt_in[1] * gt_in[3] + gt_in[0] * gt_in[4]) * inv_det return gt_out ------ If you have received this transmission in error please notify us immediately by return e-mail and delete all copies. If this e-mail or any attachments have been sent to you in error, that error does not constitute waiver of any confidentiality, privilege or copyright in respect of information in the e-mail or attachments. Please consider the environment before printing this email. ------ _______________________________________________ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev