21

I am VERY new to this whole GIS thing.

I am trying to convert some projected geoTiff images to WGS84 using gdal in python. I have found a post that outlines the process to transform points within projected GeoTiffs using something similar to the following:

from osgeo import osr, gdal

# get the existing coordinate system
ds = gdal.Open('path/to/file')
old_cs= osr.SpatialReference()
old_cs.ImportFromWkt(ds.GetProjectionRef())

# create the new coordinate system
wgs84_wkt = """
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.01745329251994328,
        AUTHORITY["EPSG","9122"]],
    AUTHORITY["EPSG","4326"]]"""
new_cs = osr.SpatialReference()
new_cs .ImportFromWkt(wgs84_wkt)

# create a transform object to convert between coordinate systems
transform = osr.CoordinateTransformation(old_cs,new_cs) 

#get the point to transform, pixel (0,0) in this case
width = ds.RasterXSize
height = ds.RasterYSize
gt = ds.GetGeoTransform()
minx = gt[0]
miny = gt[3] + width*gt[4] + height*gt[5] 

#get the coordinates in lat long
latlong = transform.TransformPoint(x,y) 

If I want to convert these points and create a new WGS84 GeoTiff file, is this the best way to go about it?

Is there a function that exists that will do such as task in 1 step?

2 Answers 2

26

One simpler way would be to use the GDAL command line tools:

gdalwarp infile.tif outfile.tif -t_srs "+proj=longlat +ellps=WGS84"

That can be invoked easily enough via scripting for batch jobs. This will warp the raster to a new grid, and there are other options to control the details.

Link

The target (t_srs) coordinate system can be PROJ.4 as shown, or an actual file with WKT, amongst other options. The source (-s_srs) coordinate system is assumed known, but can be set explicitly like the target.

That might be easier to get the job done if you don't have to use the GDAL API via Python.

There is a tutorial for warping with the API here, it says the support for Python is incomplete (I don't know the details): Link

0
17

As mdsumner said, it's much easier to use command line than the python bindings, unless you want to execute very complex tasks.

So, if you like python, like I do, you can run the command line tool with:

import os  
os.sys('gdalwarp infile.tif outfile.tif -t_srs "+proj=longlat +ellps=WGS84"')

or iterate through a list of files:

listoffiles=['file1, file2, file3, filen']  
for file in listoffiles:
    os.system('gdalwarp %s %s -t_srs "+proj=longlat +ellps=WGS84"' %(file,'new_'+file))  

And even use multiprocessing tools use the full power of your machine to execute big tasks.

1
  • 1
    My answer is outdated, subprocess.call is a better approach.
    – Pablo
    Jan 28, 2013 at 10:49

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.