[gdal-dev] numpy array to raster

Anton Korosov anton.korosov at nersc.no
Wed Feb 1 03:12:17 EST 2012


Probably the problem is with data type. You obviously have data of float 
type and you give GDT_Byte to the gdal.GetDriverByName('GTiff').Create().
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.

On 02/01/2012 03:31 AM, questions anon wrote:
> no I quickly checked that!
> myarray is:
> [[ 16.15035553  16.14380074  16.15581551 ...,  18.06388149  18.08930645
>     18.08825245]
>   [ 16.2154911   16.21180592  16.23977184 ...,  18.1085537   18.12040272
>     18.12342682]
>   [ 16.32851467  16.29202938  16.28964043 ...,  18.16753635  18.14905453
>     18.16632977]
>   ...,
>   [ 30.50812759  30.58384018  30.66707535 ...,  26.31020527  26.47789459
>     27.11495361]
>   [ 30.76499577  30.76497536  30.79138317 ...,  26.45928288  26.64059887
>     27.03641129]
>   [ 31.01263275  30.96269417  30.9933857  ...,  26.78247185  26.77845631
>     26.97975636]]
>
> On Wed, Feb 1, 2012 at 1:25 PM, Kyle Shannon <ksshannon at gmail.com
> <mailto:ksshannon at gmail.com>> wrote:
>
>     Is 'myarray' full of zeros?
>
>     /**
>       *
>       * Kyle Shannon
>       * ksshannon at gmail.com <mailto:ksshannon at gmail.com>
>       *
>       */
>
>
>
>     On Tue, Jan 31, 2012 at 19:20, questions anon
>     <questions.anon at gmail.com <mailto:questions.anon at gmail.com>> wrote:
>
>         thanks Frank, following your instructions with:
>
>         src_ds=gdal_array.OpenArray(myarray)
>         dst_ds =
>         gdal.GetDriverByName('GTiff').Create('E:/test/rasterise/mynewraster.tif',ncols,
>         nrows, 1 ,gdal.GDT_Byte)
>         dst_ds.SetGeoTransform(geotransform)
>         dst_ds.GetRasterBand(1).WriteArray(myarray)
>
>         I do not receive any error messages but the tiff produced are
>         all just zeros. Is there a step I am missing?
>         thanks
>
>
>         On Wed, Feb 1, 2012 at 12:07 PM, Frank Warmerdam
>         <warmerdam at pobox.com <mailto:warmerdam at pobox.com>> wrote:
>
>             On Tue, Jan 31, 2012 at 4:38 PM, questions anon
>             <questions.anon at gmail.com <mailto:questions.anon at gmail.com>>
>             wrote:
>              > I need to output my numpy array as a raster so that
>             someone else can access
>              > the data in ArcGIS. So basically the steps I need are:
>              > read numpy array into gdal
>              > convert to raster
>              > use latitude and longitude and array size to set projection
>              >
>              > I am really struggling with gdal because I can't seem to
>             find enough
>              > documentation about each step to understand what it is doing.
>              > Here are some of the steps I think I need:
>              >
>              > myarray=myarray
>              > #the extent and shape of my array
>              > xmin,ymin,xmax,ymax=[139.8,-39.2,150.0,-33.6]
>              > ncols,nrows=[193,106]
>              > xres=(xmax-xmin)/float(ncols)
>              > yres=(ymax-ymin)/float(nrows)
>              > geotransform=(xmin,xres,0,ymax,0, -yres)
>              >
>              > from osgeo import gdal
>              > from osgeo import gdal_array
>              >
>              > src_ds=gdal_array.OpenArray(myarray)
>              >
>              > dst_ds =
>              >
>             gdal.GetDriverByName('GTiff').Create('E:/test/rasterise/mynewraster.tif',ncols,
>              > nrows, 1 ,gdal.GDT_Byte)
>              > dst_rb = dst_ds.GetRasterBand(0)
>              > dst_ds.SetGeoTransform(geotransform)
>              > output = gdal.RasterizeLayer(dst_ds)
>
>             Dear "Questions Anon",
>
>             There are a variety of Python related information at:
>
>             http://trac.osgeo.org/gdal/wiki/GdalOgrInPython
>
>             One serious issue with the above is that
>             "gdal.RasterizeLayer()" is used to
>             turn vector data into raster data - for instance to
>             rasterize polygon features
>             into an existing raster file.  I think you want to write
>             your array into a
>             raster file.  I think you can just call
>             "dst_ds.WriteArray(myarray)".
>
>             Alternatively if that method does not exist on the dataset, you
>             can write the one band like:
>
>               dst_ds.GetRasterBand(1).WriteArray(myarray)
>
>             Some of the samples referenced from the GdalOgrInPython
>             should be
>             helpful though I understand it can be hard to know where to
>             look.
>
>             Best regards,
>             --
>             ---------------------------------------+--------------------------------------
>             I set the clouds in motion - turn up   | Frank Warmerdam,
>             warmerdam at pobox.com <mailto:warmerdam at pobox.com>
>             light and sound - activate the windows |
>             http://pobox.com/~warmerdam <http://pobox.com/%7Ewarmerdam>
>             and watch the world go round - Rush    | Geospatial Software
>             Developer
>
>
>
>         _______________________________________________
>         gdal-dev mailing list
>         gdal-dev at lists.osgeo.org <mailto: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