Thanks for responding Anton, I have been playing with the datatype but so far it hasn't fixed the problem. Also I do have access to ArcGIS to check the output.<br>I will keep trying, this is where I am at currently:<br>
<br>from osgeo import gdal<br>from osgeo import gdal_array<br>from osgeo import osr<br><br>myarray=myarray.astype(N.float32)<br>xmin,ymin,xmax,ymax=[139.8,-39.2,150.0,-33.6]<br>ncols,nrows=[193,106]<br>xres=(xmax-xmin)/float(ncols)<br>
yres=(ymax-ymin)/float(nrows)<br>geotransform=(xmin,xres,0,ymax,0, -yres)<br>dst_ds = gdal.GetDriverByName('GTiff').Create('E:/test/rasterise/myraster.tif',ncols, nrows, 1 ,gdal.GDT_Float32)<br>dst_ds.SetGeoTransform(geotransform)<br>
srs = osr.SpatialReference()<br>srs.ImportFromEPSG(4326) #WGS84 lat long.<br>dst_ds.SetProjection( srs.ExportToWkt() )<br>dst_ds.GetRasterBand(1).WriteArray(myarray)<br><br>dst_ds=None<br><br><br><div class="gmail_quote">
On Wed, Feb 1, 2012 at 7:12 PM, Anton Korosov <span dir="ltr"><<a href="mailto:anton.korosov@nersc.no">anton.korosov@nersc.no</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0pt 0pt 0pt 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
Probably the problem is with data type. You obviously have data of float type and you give GDT_Byte to the gdal.GetDriverByName('GTiff').<u></u>Create().<br>
Either try to create array with Byte data or give GDT_Float32. In the latter case the produced geotiff will have 32 bits fo each pixel. It is not supported by simple viewers but should be read with ArcGIS.<div class="im">
<br>
<br>
On 02/01/2012 03:31 AM, questions anon wrote:<br>
</div><blockquote class="gmail_quote" style="margin:0pt 0pt 0pt 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class="im">
no I quickly checked that!<br>
myarray is:<br>
[[ 16.15035553 16.14380074 16.15581551 ..., 18.06388149 18.08930645<br>
18.08825245]<br>
[ 16.2154911 16.21180592 16.23977184 ..., 18.1085537 18.12040272<br>
18.12342682]<br>
[ 16.32851467 16.29202938 16.28964043 ..., 18.16753635 18.14905453<br>
18.16632977]<br>
...,<br>
[ 30.50812759 30.58384018 30.66707535 ..., 26.31020527 26.47789459<br>
27.11495361]<br>
[ 30.76499577 30.76497536 30.79138317 ..., 26.45928288 26.64059887<br>
<a href="tel:27.03641129" value="+12703641129" target="_blank">27.03641129</a>]<br>
[ 31.01263275 30.96269417 30.9933857 ..., <a href="tel:26.78247185" value="+12678247185" target="_blank">26.78247185</a> <a href="tel:26.77845631" value="+12677845631" target="_blank">26.77845631</a><br>
<a href="tel:26.97975636" value="+12697975636" target="_blank">26.97975636</a>]]<br>
<br>
On Wed, Feb 1, 2012 at 1:25 PM, Kyle Shannon <<a href="mailto:ksshannon@gmail.com" target="_blank">ksshannon@gmail.com</a><br></div><div class="im">
<mailto:<a href="mailto:ksshannon@gmail.com" target="_blank">ksshannon@gmail.com</a>>> wrote:<br>
<br>
Is 'myarray' full of zeros?<br>
<br>
/**<br>
*<br>
* Kyle Shannon<br></div>
* <a href="mailto:ksshannon@gmail.com" target="_blank">ksshannon@gmail.com</a> <mailto:<a href="mailto:ksshannon@gmail.com" target="_blank">ksshannon@gmail.com</a>><div class="im"><br>
*<br>
*/<br>
<br>
<br>
<br>
On Tue, Jan 31, 2012 at 19:20, questions anon<br></div><div class="im">
<<a href="mailto:questions.anon@gmail.com" target="_blank">questions.anon@gmail.com</a> <mailto:<a href="mailto:questions.anon@gmail.com" target="_blank">questions.anon@gmail.<u></u>com</a>>> wrote:<br>
<br>
thanks Frank, following your instructions with:<br>
<br>
src_ds=gdal_array.OpenArray(<u></u>myarray)<br>
dst_ds =<br>
gdal.GetDriverByName('GTiff').<u></u>Create('E:/test/rasterise/<u></u>mynewraster.tif',ncols,<br>
nrows, 1 ,gdal.GDT_Byte)<br>
dst_ds.SetGeoTransform(<u></u>geotransform)<br>
dst_ds.GetRasterBand(1).<u></u>WriteArray(myarray)<br>
<br>
I do not receive any error messages but the tiff produced are<br>
all just zeros. Is there a step I am missing?<br>
thanks<br>
<br>
<br>
On Wed, Feb 1, 2012 at 12:07 PM, Frank Warmerdam<br></div><div class="im">
<<a href="mailto:warmerdam@pobox.com" target="_blank">warmerdam@pobox.com</a> <mailto:<a href="mailto:warmerdam@pobox.com" target="_blank">warmerdam@pobox.com</a>>> wrote:<br>
<br>
On Tue, Jan 31, 2012 at 4:38 PM, questions anon<br></div>
<<a href="mailto:questions.anon@gmail.com" target="_blank">questions.anon@gmail.com</a> <mailto:<a href="mailto:questions.anon@gmail.com" target="_blank">questions.anon@gmail.<u></u>com</a>>><div>
<div class="h5"><br>
wrote:<br>
> I need to output my numpy array as a raster so that<br>
someone else can access<br>
> the data in ArcGIS. So basically the steps I need are:<br>
> read numpy array into gdal<br>
> convert to raster<br>
> use latitude and longitude and array size to set projection<br>
><br>
> I am really struggling with gdal because I can't seem to<br>
find enough<br>
> documentation about each step to understand what it is doing.<br>
> Here are some of the steps I think I need:<br>
><br>
> myarray=myarray<br>
> #the extent and shape of my array<br>
> xmin,ymin,xmax,ymax=[139.8,-<u></u>39.2,150.0,-33.6]<br>
> ncols,nrows=[193,106]<br>
> xres=(xmax-xmin)/float(ncols)<br>
> yres=(ymax-ymin)/float(nrows)<br>
> geotransform=(xmin,xres,0,<u></u>ymax,0, -yres)<br>
><br>
> from osgeo import gdal<br>
> from osgeo import gdal_array<br>
><br>
> src_ds=gdal_array.OpenArray(<u></u>myarray)<br>
><br>
> dst_ds =<br>
><br>
gdal.GetDriverByName('GTiff').<u></u>Create('E:/test/rasterise/<u></u>mynewraster.tif',ncols,<br>
> nrows, 1 ,gdal.GDT_Byte)<br>
> dst_rb = dst_ds.GetRasterBand(0)<br>
> dst_ds.SetGeoTransform(<u></u>geotransform)<br>
> output = gdal.RasterizeLayer(dst_ds)<br>
<br>
Dear "Questions Anon",<br>
<br>
There are a variety of Python related information at:<br>
<br>
<a href="http://trac.osgeo.org/gdal/wiki/GdalOgrInPython" target="_blank">http://trac.osgeo.org/gdal/<u></u>wiki/GdalOgrInPython</a><br>
<br>
One serious issue with the above is that<br>
"gdal.RasterizeLayer()" is used to<br>
turn vector data into raster data - for instance to<br>
rasterize polygon features<br>
into an existing raster file. I think you want to write<br>
your array into a<br>
raster file. I think you can just call<br>
"dst_ds.WriteArray(myarray)".<br>
<br>
Alternatively if that method does not exist on the dataset, you<br>
can write the one band like:<br>
<br>
dst_ds.GetRasterBand(1).<u></u>WriteArray(myarray)<br>
<br>
Some of the samples referenced from the GdalOgrInPython<br>
should be<br>
helpful though I understand it can be hard to know where to<br>
look.<br>
<br>
Best regards,<br>
--<br>
------------------------------<u></u>---------+--------------------<u></u>------------------<br>
I set the clouds in motion - turn up | Frank Warmerdam,<br></div></div>
<a href="mailto:warmerdam@pobox.com" target="_blank">warmerdam@pobox.com</a> <mailto:<a href="mailto:warmerdam@pobox.com" target="_blank">warmerdam@pobox.com</a>><div class="im"><br>
light and sound - activate the windows |<br></div>
<a href="http://pobox.com/%7Ewarmerdam" target="_blank">http://pobox.com/~warmerdam</a> <<a href="http://pobox.com/%7Ewarmerdam" target="_blank">http://pobox.com/%7Ewarmerdam</a><u></u>><div class="im">
<br>
and watch the world go round - Rush | Geospatial Software<br>
Developer<br>
<br>
<br>
<br>
______________________________<u></u>_________________<br>
gdal-dev mailing list<br></div>
<a href="mailto:gdal-dev@lists.osgeo.org" target="_blank">gdal-dev@lists.osgeo.org</a> <mailto:<a href="mailto:gdal-dev@lists.osgeo.org" target="_blank">gdal-dev@lists.osgeo.<u></u>org</a>><br>
<a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/<u></u>mailman/listinfo/gdal-dev</a><div class="im"><br>
<br>
<br>
<br>
<br>
<br>
______________________________<u></u>_________________<br>
gdal-dev mailing list<br>
<a href="mailto:gdal-dev@lists.osgeo.org" target="_blank">gdal-dev@lists.osgeo.org</a><br>
<a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/<u></u>mailman/listinfo/gdal-dev</a><br>
</div></blockquote><div class="HOEnZb"><div class="h5">
______________________________<u></u>_________________<br>
gdal-dev mailing list<br>
<a href="mailto:gdal-dev@lists.osgeo.org" target="_blank">gdal-dev@lists.osgeo.org</a><br>
<a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/<u></u>mailman/listinfo/gdal-dev</a><br>
</div></div></blockquote></div><br>