[gdal-dev] Python RasterizeLayer and gdal_rasterize

Kai Muehlbauer kai.muehlbauer at uni-bonn.de
Wed Jul 24 06:46:15 PDT 2019


Hi Ivan,

I'm doing something similar in my application.

I remember that there were some issues with all black images, but I do
not recall the details. One thing, though, from comparing your snipped
with my code, I'm explicitely flushing the rasterband of the target dataset:

band = mem_ds.GetRasterBand(1)
band.FlushCache()

before doing the rasterization. There are other minor differences. If
you want to have a look into the code I can provide the github links.

Cheers,
Kai

Am 24.07.19 um 15:18 schrieb Ivan Lucena:
> Hi there,
> 
> I am having a problem with gdal.RasterizeLayer() in Python. I look at
> the documentation and code example in forums and tutorial but I can't
> find a solution to my problem.
> 
> My code is supposed to loops thru a shape files, with countries for
> example, and create a boolean image of each country but the images
> always come as all black pixels. Itshould use GDAL and OGR memory
> drivers but in order to debug the issue, I am saving one GEOJSON and one
> GTIFF file for each country.
> 
> What is most intriguing is that I can run gdal_rasterize with any of the
> countries pairs of file and it works as expected. I open it on QGIS and
> it all looks perfectly georeferenced too. The resolution is coming from
> another remote sensing image that I am using as reference. And again,
> everything works fine with the gdal_rasterize command, so it should be
> right.
> 
> Here is my code:
> 
>         geom_ds =
> ogr.GetDriverByName('GEOJSON').CreateDataSource('E:\TMP\MEMORY_%s.JSON'
> % name)
>         geom_lyr = geom_ds.CreateLayer('MEMORY_LAYER',
> geom_type=ogr.wkbPolygon, srs=ref_srs)​
>         geom_feat = ogr.Feature(geom_lyr.GetLayerDefn())​
>         geom_feat.SetGeometry(geom)​
>        geom_lyr.CreateFeature(geom_feat)​
> 
>         mem_ds =
> gdal.GetDriverByName('GTIFF').Create('E:\TMP\MEMORY_%s.TIF' % name,
> xcount, ycount, 1, gdal.GDT_Byte)​
>         mem_ds.SetGeoTransform((​xmin, xResolut, 0,​ ymax, 0, yResolut,​))​
>         mem_ds.SetProjection(ref_ds.GetProjectionRef())​
>         mem_ds.FlushCache()​
>    ​
>         err = gdal.RasterizeLayer(mem_ds, [1], geom_lyr, burn_values=[1])​
>         if err != 0:​
>             print(err)​
> 
> Here is the command line:
> 
>     (base) E:\TMP>gdal_rasterize MEMORY_Peru.JSON" "MEMORY_Peru.TIF" -burn 1
> 
> 
> I am using GDAL 2.3.3 and Python 3.7.3 installed by Anaconda3 (64-bit).
> The function gdal.UseExceptions() was called but there is nothing been
> reported.
> 
> I might be asking too much, but if someone had the same experience and
> have found a solution and want to share, that would be awesome. Or maybe
> there is an obvious error on my code that you can point out. I tried
> several things like flushing and closing and re-opening the files,  to
> make sure that the data is ready for the functions, but nothing really
> helped.
> 
> One other thing is, once I got it to work, I am thinking that I don't
> need to prepare the little geometry file if I could call
> gdal.rasterizelayer() the original shape with a *query* for each country
> at the time. That should be more efficient, I guess.
> 
> Thank you for your attention and have a nice day everybody,
> 
> Ivan
> 
> 
> 
> 
> 
> 
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev
> 

-- 
Kai Muehlbauer
Institute for Geosciences and Meteorology University of Bonn
Auf dem Huegel 20       | +49 228 739083
D-53121 Bonn            | kai.muehlbauer at uni-bonn.de


More information about the gdal-dev mailing list