[gdal-dev] (no subject)

[Frank Warmerdam]

Discourse Maps wrote:
> I am simply trying to open a raster image, manipulate in NumPy and then spit
> it back out in GDAL.
> 
> I am having trouble setting null data...
> I tried: ds.GetRasterBand(1).SetNoDataValue( -9999 )
> but nothing seems to work.
...

Dear Discourse Maps,

Presumably you should be setting the nodata value on the
first band of the "dst_ds", not the source dataset.  But
I don't see any SetNoDataValue() call at all in your script.

I would suggest adding it right after the Create() call.

I would also note that the HFA driver does not have a COMPRESS=LZW
creation option, though you can use [ 'COMPRESSED=YES' ] to use
the grid compression scheme supported by the format.

> dst_format = 'HFA'
> dst_datatype = gdal.GDT_Int16
> dst_options = ['COMPRESS=LZW']
> dst_file='sequest_final.img'
> dst_xsize = nx
> dst_ysize = ny
> dst_nbands = 1
> 
> geoTransform = ds.GetGeoTransform()
> proj = ds.GetProjection()
> 
> 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(ny):
> nxarray = data[j,:]
> nxarray.shape = (1,nx)
> print nxarray
> dst_ds.GetRasterBand(1).WriteArray(nxarray,0,j)
> dst_ds.SetGeoTransform(geoTransform)
> dst_ds.SetProjection(proj)
> dst_ds = None

Your email client has reformatted the script so it is hard
to see the indentation.  But I would advise strongly against
doing SetGeoTransform() and SetProjection() many times on the
same dataset.

Best regards,
-- 
---------------------------------------+--------------------------------------
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | Geospatial Programmer for Rent