[gdal-dev] Creating Geotiff with WriteArray in Python

Lucena, Ivan ivan.lucena at pmldnet.com
Thu Mar 20 13:38:24 EDT 2008


32/8...

Lucena, Ivan wrote:
> Also, remember to check if your GDAL distribution supports Bigtiff if 
> the file size is greater than 4GB (13700x13300x32=5,830,720,000 ?).
> 
> gaffigan at sfos.uaf.edu wrote:
>> Hello, Matt.  I think the array argument to the WriteArray function must
>> be 2D.  In your example, nxarray should have shape (1,nx), instead of
>> (nx,).  So, on the relevant line, try the following instead:
>>   dst_ds.GetRasterBand(1).WriteArray(reshape(nxarray,(1,nx)),0,j)
>>
>> Hope this gets you closer.
>>
>> Steve
>>
>>> Dear GDAL/Python Users,
>>>
>>> I am struggling to make the following python/GDAL script work. I am
>> reading/converting large Palsar images to GeoTiff, and therefore I have
>> to be concerned with memory. I believe I have the following script very
>> close to working. The problem is related to the line
>> "dst_ds.GetRasterBand(1).WriteArray(nxarray,0,j)" - I am trying to write
>> an imported integer array(y,x) , 13700x13300, to a tif file row by row
>> using a for loop structure as a converted float. I constructed this code
>> from found examples, ie, 4bandto3band.py. However,
>>> "WriteArray" is giving me trouble.  I believe the problem is something
>> very simple; the arguments in WriteArray are not quite right. I am quite
>> new to GDAL and Python. Some help would be very much
>>> appreciated. Below are the errors I receive and my script.
>>> thanks,
>>> Matt
>>>
>>>
>>> Traceback (most recent call last):
>>>   File "./test.py", line 43, in <module>
>>>     dst_ds.GetRasterBand(1).WriteArray(nxarray,0,j)
>>>   File "/usr/lib/python2.5/site-packages/gdal.py", line 876, in
>> WriteArray
>>>     return gdalnumeric.BandWriteArray( self, array, xoff, yoff )
>>>   File "/usr/lib/python2.5/site-packages/gdalnumeric.py", line 187, in
>>> BandWriteArray
>>>     xsize = array.shape[1]
>>> IndexError: tuple index out of range
>>>
>>>
>>>
>>>
>>> #! /usr/bin/env python
>>>
>>> import Numeric
>>> from Numeric import *
>>> import pickle, array, sys, os, shutil
>>> import gdal, gdalnumeric
>>> from gdalconst import *
>>> import numpy
>>> from math import *
>>>
>>>
>>> file =
>>> "china/hangzhou/processed/ALPSRP021230580/IMG-HH-ALPSRP021230580-H1.5GUA" 
>>>
>>>
>>> gd=gdal.Open(file)
>>> array=gd.ReadAsArray()
>>> dims = array.shape
>>> nx=dims[1]
>>> ny=dims[0]
>>>
>>> nxarray=zeros(nx,Float)
>>> cal_factor=-83.0
>>>
>>> dst_format = 'GTiff'
>>> dst_datatype = gdal.GDT_Float32
>>> dst_options = ['COMPRESS=LZW']
>>> dst_file='foo.tif'
>>> dst_xsize = nx
>>> dst_ysize = ny
>>> dst_nbands = 1
>>>
>>> driver = gdal.GetDriverByName( dst_format )
>>> dst_ds = driver.Create(dst_file, dst_xsize, dst_ysize, dst_nbands,
>> dst_datatype, dst_options)
>>>
>>> for j in range(0,ny-1):
>>>   for i in range(0,nx-1):
>>>
>>>     nxarray[i]=0.0
>>>     tmp = array[j,i] * 1.0
>>>     if (tmp > 0):
>>>       nxarray[i] = 10.0 * log10(pow(tmp,2)) + cal_factor
>>>
>>>   dst_ds.GetRasterBand(1).WriteArray(nxarray,0,j)
>>>
>>> dst_ds = None
>>> _______________________________________________
>>> gdal-dev mailing list
>>> gdal-dev at lists.osgeo.org
>>> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>>>
>>
>>
>>
>>
>> _______________________________________________
>> gdal-dev mailing list
>> gdal-dev at lists.osgeo.org
>> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>>
>>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
> 
> 


More information about the gdal-dev mailing list