<HTML>
Hi Vincent,<br>
<br>
Thanks for the repliy; unfortunately <br>
<br>
gdal_array.BandWriteArray(ds.GetRasterBand(1),a) <br>
(BandWriteAsArray doesn't exist)<br>
<br>
has the same effect.<br>
<br>
Tom<br>
<br>
<br>
<span style="font-weight: bold;">On Fri 28/10/11 10:57 , Vincent Schut schut@sarvision.nl sent:<br>
</span><blockquote style="BORDER-LEFT: #F5F5F5 2px solid; MARGIN-LEFT: 5px; MARGIN-RIGHT:0px; PADDING-LEFT: 5px; PADDING-RIGHT: 0px">
<defanged_meta content="text/html; charset=UTF-8" http-equiv="Content-Type">
<defanged_body text="#000000" bgcolor="#FFFFFF">
Tom van der Putte wrote:
<blockquote cite="mid:4335.1319790184@vdputte.nl" type="cite">
Hi
List,<br>
<br>
I'm trying to save a numpy array to a GTiff, which should be quite
straightforward, considering all the examples available on the
web. However, I get an file that is not readable by Qgis, is
regarded by ArcMap to have only 0 and gdalinfo displays nothing on
the values.<br>
<br>
I start with a ndarray, dtype = int32, shape = ( 10,10). It has
integer values between 1 and 10. The code is as follows:<br>
<br>
import numpy as np<br>
from osgeo import gdal<br>
from osgeo import osr<br>
<br>
#np array<br>
a = np.ones((10,10))<br>
<br>
#coord-system<br>
sr_str = 'LOCAL_CS["arbitrary"]'<br>
sr = osr.SpatialReference( sr_str )<br>
<br>
#driver<br>
driver = gdal.GetDriverByName("GTiff")<br>
<br>
#dataset<br>
ds = driver.Create("D:\\temp\\output.tif", 10, 10, 1,
gdal.GDT_Byte)<br>
ds.SetProjection( sr_str )<br>
ds.SetGeoTransform((0, 1, 0,0, 0,1))<br>
<br>
#write and close<br>
ds.GetRasterBand(1).WriteArray(a)<br>
ds = None<br>
</blockquote>
<br>
Tom,<br>
<br>
don't know if you can use WriteArray this way... what I always do,
is:<br>
<br>
from osgeo import gdal_array<br>
# to write numpy array:<br>
gdal_array.BandWriteAsArray(ds.GetRasterBand(1), a)<br>
<br>
Best,<br>
Vincent.<br>
<blockquote cite="mid:4335.1319790184@vdputte.nl" type="cite">
<br>
The corresponding gdalinfo is:<br>
*********<br>
Raster dataset parameters:<br>
Projection:
LOCAL_CS["arbitrary",UNIT["metre",1,AUTHORITY["EPSG","9001"]]]<br>
RasterCount: 1<br>
RasterSize (10,10)<br>
Using driver GeoTIFF<br>
Metadata:<br>
0: AREA_OR_POINT=Area<br>
<br>
Image Structure Metadata:<br>
0: INTERLEAVE=BAND<br>
<br>
Corner Coordinates:<br>
Upper Left (0, 0)<br>
Lower Left (0, 10)<br>
Upper Right (10, 0)<br>
Lower Right (10, 10)<br>
Center (5, 5)<br>
<br>
Coordinate System is:<br>
LOCAL_CS["arbitrary",<br>
UNIT["metre",1,<br>
AUTHORITY["EPSG","9001"]]]<br>
Band 1 :<br>
DataType: Byte<br>
ColorInterpretation: Gray<br>
Description:<br>
Size (10,10)<br>
BlockSize (10,10)<br>
*********<br>
<br>
I have tried both version 1.6.3 (I'm trying to make a Qgis plugin,
and thus restricted to 1.6.3) as well as 1.8.1. I have also tried
ds = driver.Create("D:\\temp\\output.tif", 10, 10, 1,
gdal.GDT_Int32), but that didn't make much difference.<br>
<br>
I tried to make a Surfer grid (GSBG), and this somewhat worked: it
gave a nice picture, but the values were off the charts.<br>
<br>
Any clues as to what I'm missing here? Cheers,<br>
<br>
Tom van der Putte<br>
<br>
<br>
<fieldset class="mimeAttachmentHeader"></fieldset>
<br>
<pre wrap="">_______________________________________________
gdal-dev mailing list
<a class="moz-txt-link-abbreviated" href="javascript:top.opencompose('gdal-dev@lists.osgeo.org','','','')">gdal-dev@lists.osgeo.org</a>
<a class="moz-txt-link-freetext" href="http://lists.osgeo.org/mailman/listinfo/gdal-dev">http://lists.osgeo.org/mailman/listinfo/gdal-dev</a></pre>
</blockquote>
</defanged_body></defanged_meta></blockquote></HTML>