[gdal-dev] Python error (missing?) when performing RasterizeLayer on specific polygon with two very close vertices

Tue Sep 9 09:08:16 PDT 2014

I have tried with GDAL 1.11 and Python 3.2 (again the builds are from GISInternals), and it gives me the same result. Python just terminated without any warning or error.

Both versions of GDAL are x64.

Kind regards, Casper

From: gdal-dev-bounces at lists.osgeo.org [mailto:gdal-dev-bounces at lists.osgeo.org] On Behalf Of Casper Børgesen (CABO)
Sent: 9. september 2014 17:16
To: gdal-dev at lists.osgeo.org
Subject: [gdal-dev] Python error (missing?) when performing RasterizeLayer on specific polygon with two very close vertices

I have the following polygon:

POLYGON ((539500.00004634191 6314429.5954758488,539499.999 6314429.5964747742,539500.0 6314500.001,539563.49086194986 6314500.001,539574.70556674246 6314491.5008581309,539591.23494657199 6314478.9664945751,539600.001 6314472.1168122729,539600.00014411798 6314460.5601839637,539543.07124903298 6314452.2691510962,539546.30058158806 6314438.619163854,539507.3702219734 6314429.939004397,539505.68004823104 6314429.8589810915,539500.00004634191 6314429.5954758488))

I have created a memory vector layer and added this polygon as the only feature.

    projection = ...
    geometry = ...
    ogr_mem_drv = ogr.GetDriverByName('Memory')

    mem_vector_ds = ogr_mem_drv.CreateDataSource('')
    mem_vector_layer = mem_vector_ds.CreateLayer('', srs=projection)
    feature = ogr.Feature(mem_vector_layer.GetLayerDefn())
    feature.SetGeometry(geometry.Clone())
    mem_vector_layer.CreateFeature(feature)

    gdal_mem_drv = gdal.GetDriverByName('MEM')
    mem_raster_ds = gdal_mem_drv.Create('',5,5,1,gdal.GDT_Byte)
    mem_raster_ds.SetGeoTransform([...])
    mem_raster_ds.SetProjection(mem_vector_layer.GetSpatialRef().ExportToWkt())

    err = gdal.RasterizeLayer(mem_raster_ds, [1], mem_vector_layer, None, None, [1], ['ALL_TOUCHED=TRUE'])

When running the last line with gdal.RasterizeLayer(...) my program waits a couple of seconds and then just terminates without any exception or error.

If  I remove vertex number two from the above WKT it works. Vertex 1 and 2 are very close (located in the lower left corner of the polygon).

I am running GDAL 1.10.1 and Python 3.2, builds are stable releases from http://www.gisinternals.com/sdk/.

What could be the reason for this strange error?

Kind regards, Casper
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20140909/621090f3/attachment.html>