TLDR;
gdal_rasterize has begun silently producing wrong output at least since GDAL
1.11.0.
Whenever >500 polygons are present in the database query result (and perhaps
under other circumstances), only 500 will be drawn.
Hi,
I am rasterizing many maps and areas of maps in SRID 25833 using gdal_rasterize
with a postgis database.
When I use GDAL 1.11, the output is missing much of the data, despite
gdal_rasterize reporting completing successfully.
With GDAL 1.09 and GDAL 1.10, the output is correct and complete and identical.
I'm posting here to gather more information about the bug and also to raise an
alarm for anyone using postgis/gdal_rasterize (1.11) in their data workflow
over the last year.
Background of bug
--------
Versions:
POSTGIS="2.1.5 r13152" GEOS="3.4.2-CAPI-1.8.2 r3921" PROJ="Rel. 4.8.0, 6 March
2012" GDAL="GDAL 1.11.0, released 2014/04/16" (from EPEL RPM)
A. The source table in postgis is registered in geometry_columns correctly.
B. I run this command in Centos linux.
gdal_rasterize -te 1060000 7865000 1120000 7940000 -a jord -init 255 -tr 100
100 -a_nodata 255 -a_srs 'EPSG:25833' -co 'TILED=YES' -co 'BIGTIFF=IF_SAFER'
-co 'COMPRESS=PACKBITS' -co 'PREDICTOR=1' -co 'ZLEVEL=1' -ot Byte
PG:"host=XXXX'' dbname='XXXX' user='XXXX' password='XXXX' tables='XXXX'" -sql
"SELECT xxxx FROM xxxx WHERE st_intersects(geo,
ST_GeometryFromText('POLYGON((1060000 7865000,1060000 7940000,1120000 7940000,
1120000 7865000, 1060000 7865000))',25833)) " xxxx.tif
C. The SQL query (-sql), when run directly on postgres with COUNT(), produces
this output.
count
-------
3959
(1 row)
D. GDAL 1.9.2 and GDAL 1.10.0 produce a TIF file of 12238 bytes which is
visually complete.
E. GDAL 1.11.0 produces a TIF file of 10920 bytes which is missing many
polygons.
F. Console output from gdal_rasterize in all 3 cases is as follows;
0...10...20...30...40...50...60...70...80...90...100 - done.
G. PG server log output for GDAL 1.9.2 and GDAL 1.10.0
(nothing...)
PG server log output for GDAL 1.11.0
< 2015-02-09 13:52:11.720 CET >STATEMENT: FETCH 500 in OGRPGLayerReader0xd5ce50
< 2015-02-09 13:52:11.720 CET >ERROR: cursor "ogrpglayerreader0xd5ce50" does
not exist
< 2015-02-09 13:52:11.720 CET >STATEMENT: CLOSE OGRPGLayerReader0xd5ce50
< 2015-02-09 13:52:18.657 CET >ERROR: cursor "ogrpglayerreader0x1fd8e50" does
not exist
H. Debug mode output for GDAL 1.11.0
OGR: OGROpen(...) succeeded as PostgreSQL.
GDAL: QuietDelete(xxxx.tif) invoking Delete()
GDAL: GDALOpen(xxxx.tif, this=0x259a8f0) succeeds as GTiff.
GDAL: GDALDefaultOverviews::OverviewScan()
GDAL: GDALClose(xxxx.tif, this=0x259a8f0)
GDAL: GDALDriver::Create(GTiff,xxxx.tif,600,750,1,Byte,0x25744e0)
PG: Primary key name (FID): gid
PG: Using column 'gid' as FID for table 'xxxx_25833'
GDAL: Rasterizer operating on 1 swaths of 750 scanlines.
0...10...20...30...40...50...60...70...80...90...100 - done.
PG: 500 features read on layer 'sql_statement'.
GDAL: GDALClose(xxxx.tif, this=0x2598930)
I. Debug mode output for GDAL 1.10.0
OGR: OGROpen(xxxx) succeeded as PostgreSQL.
GDAL: QuietDelete(xxxx) invoking Delete()
GDAL: GDALOpen(xxxx) succeeds as GTiff.
GDAL: GDALDefaultOverviews::OverviewScan()
GDAL: GDALClose(xxxx.tif, this=0xeb5030)
GDAL: GDALDriver::Create(GTiff,xxxx.tif,600,750,1,Byte,0xe8ede0)
PG: Primary key name (FID): gid
PG: Using column 'gid' as FID for table 'dmk_25833'
PG: Layer 'xxxx_25833' geometry type: GEOMETRY:Unknown (any), Dim=2
0...10...20...30...40...50...60...70...80...90...100 - done.
PG: 3959 features read on layer 'sql_statement'.
GDAL: GDALClose(xxxx.tif, this=0xeb2fd0)
As well as local unix sockets I've tested TCP connections to a remote PG
database from GDAL1.10 and GDAL1.11 clients, and I've tested with other data
layers. Every time, I'm getting the same problem e.g. 500 features out 8000.
Has anyone else experienced anything similar?
I've spent some time looking through the source but I wasn't able to determine
exactly where the problem is occuring.
It seems to somehow be related to the select page size (500 rows) being used in
place of the total feature count; or perhaps some logic mistake causing the
reads to finish after 1 page is complete.
Graeme Bell
_______________________________________________
gdal-dev mailing list
[email protected]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
