[gdal-dev] OGR Spatialite driver and spatial index
Hi, Does the OGR SQLite/Spatialite driver really utilise the spatial index if such exists in the Spatialite database? The documentation suggests so http://gdal.org/ogr/drv_sqlite.html. However, I haven't been able to see any difference in my speed tests with ogr yet. I made a set of tests with ogrinfo against a Spatialite database with and without spatial index. These two tests took the same 30 seconds to run both with and without spatial index. ogrinfo OSM_Finland.sqlite -sql select geometry, osm_id ,highway,ref, name, tunnel from osm_line where highway is not null -spat 389116 6677305 389579 6677661 ogrinfo OSM_Finland.sqlite -sql select geometry, osm_id ,highway,ref, name, tunnel from osm_line where highway is not null AND MBRIntersects(geometry, BuildMBR(389116,6677305,389579,6677661)) This one takes only 2 seconds if spatial index table is available. ogrinfo OSM_Finland.sqlite -sql SELECT geometry, osm_id, highway,ref, name, tunnel FROM osm_line WHERE highway IS NOT NULL AND ROWID IN (SELECT pkid FROM idx_osm_line_GEOMETRY WHERE xmax 389116 AND xmin 389579 AND ymax 6677305 AND ymin 6677661) The Spatialite man Alessandro Furieri wrote this information: SpatiaLite isn't PostGIS: you *must* explicitly write your SQL Queries in such a way to access the corresponding Spatial Index table as appropriate. Please note well: in SQLite/SpatiaLite the R*Tree Spatial Index simply is another table between many others. The SQL engine has absolutely no idea that a strict correlation exists between the geometry table and the corresponding R*Tree. So you are explicitly required to define an explicit sub-query in order to inquiry the R*Tree. Is the OGR Spatialite driver clever enough for making a conclusion If exists spatial index table 'idx_table_GEOMETRY' for layer 'table' do use if when building bounding box filter? I found a blog posting about this thing and the author tells how he could make SharpMap to render blisteringly fast from Spatialite http://epsg27700.blogspot.com/2009/08/adventures-with-spatialite.html I would not call my OpenStreetMap rendering from Spatialite with Mapserver 6.0 (MS4W 3.0.3) blistering at the moment. The highway layer takes more than 5 minutes to render even when zoomed very close. The corresponding SQL query which is utilising the idx_osm_line_GEOMETRY takes less than a second to run directly through Spatialite-GUI so there is absolutely something to improve. -Jukka Rahkonen- ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] OGR Spatialite driver and spatial index
Selon Jukka Rahkonen jukka.rahko...@mmmtike.fi: Hi, ogrinfo OSM_Finland.sqlite -sql select geometry, osm_id ,highway,ref, name, tunnel from osm_line where highway is not null AND MBRIntersects(geometry, BuildMBR(389116,6677305,389579,6677661)) The MBRIntersects() stuff is exactly what the OGR driver does, when you specify -spat, so it is not surprising that they have the same runtime. This one takes only 2 seconds if spatial index table is available. ogrinfo OSM_Finland.sqlite -sql SELECT geometry, osm_id, highway,ref, name, tunnel FROM osm_line WHERE highway IS NOT NULL AND ROWID IN (SELECT pkid FROM idx_osm_line_GEOMETRY WHERE xmax 389116 AND xmin 389579 AND ymax 6677305 AND ymin 6677661) The Spatialite man Alessandro Furieri wrote this information: SpatiaLite isn't PostGIS: you *must* explicitly write your SQL Queries in such a way to access the corresponding Spatial Index table as appropriate. Please note well: in SQLite/SpatiaLite the R*Tree Spatial Index simply is another table between many others. The SQL engine has absolutely no idea that a strict correlation exists between the geometry table and the corresponding R*Tree. So you are explicitly required to define an explicit sub-query in order to inquiry the R*Tree. Is the OGR Spatialite driver clever enough for making a conclusion If exists spatial index table 'idx_table_GEOMETRY' for layer 'table' do use if when building bounding box filter? Hum, the OGR driver was indeed written a bit too naively, assuming that MBRIntersects() would be similar to the operator in PostGIS. Apparently not. Strange that this hasn't been reported before, but I assume you can only see the difference with big amount of data, like your OSM database (and the MBRIntersects() test must be faster than evaluating the spatial filter on OGR side, so that's perhaps why it wasn't seen at implementation time). Anyway, would you mind opening a Trac ticket about this ? Changing the request to use the spatial index table should not be too difficult. Thanks, Even ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] reading postgis raster in mode 2 error
2011/8/26 Jorge Arévalo jorge.arev...@deimos-space.com: -- Forwarded message -- From: Ricardo Filipe Soares Garcia da ricardo.garcia.si...@gmail.com Date: Wed, Aug 24, 2011 at 6:07 PM Subject: reading postgis raster in mode 2 error To: gdal-dev@lists.osgeo.org Hi list I am trying out the postgis raster driver. I'm running Ubuntu 11.04 with gdal 1.8 from the ubuntugis repository. As for postgis I've compiled a snapshot of postgis 2.0 from svn, as offered on the postgis website. Following the directions in the driver info page[1], I've been trying to get a gdalinfo on the katrina.tif file. When I run it using mode='1' everything goes as expected and I get info on each of the tiles. but when using mode='2' I am getting this error: gdalinfo -mm -stats -checksum PG:host='localhost' dbname='gis_testes' user='gisuser' password='resusig' table='katrina' mode='2' ERROR 1: Error, the ONE_RASTER_PER_TABLE mode can't be applied if the raster rows don't have the same metadata for band 1 gdalinfo failed - unable to open 'PG:host='localhost' dbname='gis_testes' user='test_user' password='test_pass' table='katrina' mode='2''. When importing the file into postgis I ran the following command: /usr/lib/postgresql/8.4/bin/raster2pgsql.py -r ~/Downloads/katrina.tif -t katrina -k 64x64 -o katrina.sql -s 4326 -I Thanks for your help [1] - http://trac.osgeo.org/gdal/wiki/frmts_wtkraster.html -- ___ ___ __ Ricardo Garcia Silva Hi Ricardo, That error is caused because this query select (foo.md).* from (select distinct st_bandmetadata(rast, n) as md from katrina) as foo (being 'n' the band number). returns more than one result, and it shouldn't. Could you please execute that query in a client and send me back the results? Best regards, and sorry for the inconvenience -- Jorge Arévalo Internet Mobility Division, DEIMOS jorge.arev...@deimos-space.com http://es.linkedin.com/in/jorgearevalo80 http://mobility.grupodeimos.com/ http://gis4free.wordpress.com http://geohash.org/ezjqgrgzz0g Hi Jorge, list Running the query you suggest I get this output: gis_testes= select (foo.md).* from (select distinct st_bandmetadata(rast, 1) as md from katrina) as foo; pixeltype | hasnodatavalue | nodatavalue | isoutdb | path ---++--+-+-- 8BUI | f | -1.82132e-05 | f | 8BUI | f | -1.81834e-05 | f | 8BUI | f | -1.80478e-05 | f | 8BUI | f | -1.80471e-05 | f | (4 rows) Can it be that the nodatavalue is messing things up, since it is not the same for all the blocks? Thanks for your help -- ___ ___ __ Ricardo Garcia Silva ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Re: OGR Spatialite driver and spatial index
Even Rouault even.rouault at mines-paris.org writes: Hum, the OGR driver was indeed written a bit too naively, assuming that MBRIntersects() would be similar to the operator in PostGIS. Apparently not. Strange that this hasn't been reported before, but I assume you can only see the difference with big amount of data, like your OSM database (and the MBRIntersects() test must be faster than evaluating the spatial filter on OGR side, so that's perhaps why it wasn't seen at implementation time). Anyway, would you mind opening a Trac ticket about this ? Changing the request to use the spatial index table should not be too difficult. Ticket [GDAL] #4212 done with sample requests which should show the speed difference and a download link to a test database which is created with ogr2ogr v. 1.8.1. -Jukka Rahkonen- ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] reading postgis raster in mode 2 error
On Mon, Aug 29, 2011 at 12:15 PM, Ricardo Filipe Soares Garcia da ricardo.garcia.si...@gmail.com wrote: 2011/8/26 Jorge Arévalo jorge.arev...@deimos-space.com: -- Forwarded message -- From: Ricardo Filipe Soares Garcia da ricardo.garcia.si...@gmail.com Date: Wed, Aug 24, 2011 at 6:07 PM Subject: reading postgis raster in mode 2 error To: gdal-dev@lists.osgeo.org Hi list I am trying out the postgis raster driver. I'm running Ubuntu 11.04 with gdal 1.8 from the ubuntugis repository. As for postgis I've compiled a snapshot of postgis 2.0 from svn, as offered on the postgis website. Following the directions in the driver info page[1], I've been trying to get a gdalinfo on the katrina.tif file. When I run it using mode='1' everything goes as expected and I get info on each of the tiles. but when using mode='2' I am getting this error: gdalinfo -mm -stats -checksum PG:host='localhost' dbname='gis_testes' user='gisuser' password='resusig' table='katrina' mode='2' ERROR 1: Error, the ONE_RASTER_PER_TABLE mode can't be applied if the raster rows don't have the same metadata for band 1 gdalinfo failed - unable to open 'PG:host='localhost' dbname='gis_testes' user='test_user' password='test_pass' table='katrina' mode='2''. When importing the file into postgis I ran the following command: /usr/lib/postgresql/8.4/bin/raster2pgsql.py -r ~/Downloads/katrina.tif -t katrina -k 64x64 -o katrina.sql -s 4326 -I Thanks for your help [1] - http://trac.osgeo.org/gdal/wiki/frmts_wtkraster.html -- ___ ___ __ Ricardo Garcia Silva Hi Ricardo, That error is caused because this query select (foo.md).* from (select distinct st_bandmetadata(rast, n) as md from katrina) as foo (being 'n' the band number). returns more than one result, and it shouldn't. Could you please execute that query in a client and send me back the results? Best regards, and sorry for the inconvenience -- Jorge Arévalo Internet Mobility Division, DEIMOS jorge.arev...@deimos-space.com http://es.linkedin.com/in/jorgearevalo80 http://mobility.grupodeimos.com/ http://gis4free.wordpress.com http://geohash.org/ezjqgrgzz0g Hi Jorge, list Running the query you suggest I get this output: gis_testes= select (foo.md).* from (select distinct st_bandmetadata(rast, 1) as md from katrina) as foo; pixeltype | hasnodatavalue | nodatavalue | isoutdb | path ---++--+-+-- 8BUI | f | -1.82132e-05 | f | 8BUI | f | -1.81834e-05 | f | 8BUI | f | -1.80478e-05 | f | 8BUI | f | -1.80471e-05 | f | (4 rows) Can it be that the nodatavalue is messing things up, since it is not the same for all the blocks? Thanks for your help -- ___ ___ __ Ricardo Garcia Silva ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev Hi, Yes, it's a problem with metadata precission. I found this error in a different part of the code. I'm not sure if is caused because changes in PostGIS Raster code, but I'm aware of it. I applied a temporal patch in r23006. Could you please update your GDAL copy? Best regards, -- Jorge Arévalo Internet Mobility Division, DEIMOS jorge.arev...@deimos-space.com http://es.linkedin.com/in/jorgearevalo80 http://mobility.grupodeimos.com/ http://gis4free.wordpress.com http://geohash.org/ezjqgrgzz0g ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Projection definition issue with Mapinfo files
Hello, One of QGIS user report an issue with mapInfo file (TAB file). If you understand french you can see the thread here: http://georezo.net/forum/viewtopic.php?pid=197022 I try to sum up the issue in english in this email. The process is basicly as follow: 1. check projection information with ogrinfo DATUM[GRS_80, SPHEROID[GRS 80,6378137,298.257222101], 2. update projection information file with ogr2ogr: ogr2ogr -f MapInfo file -s_srs EPSG:2154 -t_srs EPSG:2154 new.tab origin.tab 3. check projection information with ogrinfo DATUM[WGS_1984, SPHEROID[WGS 84,6378137,298.257223563], 1 and 3 show differents projection information! If I try to change file format to shp, it seems that projection information are correct. I can share tab files if needed. Any idea? Thanks. Y. -- Yves Jacolin http://yjacolin.gloobe.org ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Re: OGR Spatialite driver and spatial index
Le lundi 29 août 2011 12:19:26, Jukka Rahkonen a écrit : Even Rouault even.rouault at mines-paris.org writes: Hum, the OGR driver was indeed written a bit too naively, assuming that MBRIntersects() would be similar to the operator in PostGIS. Apparently not. Strange that this hasn't been reported before, but I assume you can only see the difference with big amount of data, like your OSM database (and the MBRIntersects() test must be faster than evaluating the spatial filter on OGR side, so that's perhaps why it wasn't seen at implementation time). Anyway, would you mind opening a Trac ticket about this ? Changing the request to use the spatial index table should not be too difficult. Ticket [GDAL] #4212 done with sample requests which should show the speed difference and a download link to a test database which is created with ogr2ogr v. 1.8.1. Improvement commited in r23008. Before : $ time ogrinfo berlin.sqlite osm_polygon -spat 1489000 6899000 149 690 -al -so real0m0.619s user0m0.580s sys 0m0.040s After : $ time ogrinfo berlin.sqlite osm_polygon -spat 1489000 6899000 149 690 -al -so real0m0.053s user0m0.020s sys 0m0.040s And for comparison, the same for a postgis instance : $ time ogrinfo pg:dbname=autotest osm_polygon -spat 1489000 6899000 149 690 -al -so real0m0.168s user0m0.050s sys 0m0.000s However, there is no gain to expect from the above change for a request such as : ogrinfo OSM_Finland.sqlite -sql select geometry, osm_id ,highway,ref, name, tunnel from osm_line where highway is not null -spat 389116 6677305 389579 6677661 Indeed, when you specify -sql, the driver (and to my knowlegde, this is true for all other OGR drivers) makes no attempt to merge the spatial filter, so the spatial filter is evaluated on OGR side, which will be rather slow. So, you have 2 possibilities, either use the -where highway is not null on a layer so that OGR can merge the attribute and spatial filters, or incorporate at hand the spatial filter inside the SQL query specified with -sql. -Jukka Rahkonen- ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Re: OGR Spatialite driver and spatial index
Even Rouault even.rouault at mines-paris.org writes: Improvement commited in r23008. Before : $ time ogrinfo berlin.sqlite osm_polygon -spat 1489000 6899000 149 690 -al -so real 0m0.619s user 0m0.580s sys 0m0.040s After : $ time ogrinfo berlin.sqlite osm_polygon -spat 1489000 6899000 149 690 -al -so real 0m0.053s user 0m0.020s sys 0m0.040s So eleven times faster than it used to be and now three times faster than PostGIS? I would also call it improvement. -Jukka Rahkonen- ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Re: gdal_rasterize output off by one pixel/cell
Thanks Even for making the fix! Is there documentation anywhere for how to apply changes like this? I apologize but I am just getting into this type of development. I thought I was able to install the latest nightly build today, but I'm not sure. Running Mac OS X. If anyone has any pointers I would really appreciate them off-list. Thanks, Matt On 28 Aug, 2011, at 11:20 PM, Eli Adam wrote: Even, This works for me and the test case included with the initial report. Thank you very much, Eli Even Rouault even.roua...@mines-paris.org 08/28/11 5:25 AM I've applied a fix in trunk in http://trac.osgeo.org/gdal/changeset/22998 . Please test and confirm if it works for you. - Matthew Burton-Kelly, M.S. Graduate Student Department of Geology and Geological Engineering University of North Dakota Mobile: (802) 922-3696 matthew.burton.ke...@my.und.edu matthew.burtonke...@gmail.com http://www.protichnoctem.com http://und.academia.edu/MatthewBurtonKelly About thirty years ago there was much talk that geologists ought only to observe and not theorize; and I well remember someone saying that at this rate a man might as well go into a gravelpit and count the pebbles and describe the colors. How odd it is that anyone should not see that all observation must be for or against some view if it is to be of any service! -Charles Darwin, in an 1861 letter to Henry Fawcett. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Re: gdal_rasterize output off by one pixel/cell
On 29-08-2011 22:10, Burton-Kelly, Matthew wrote: Thanks Even for making the fix! Is there documentation anywhere for how to apply changes like this? I apologize but I am just getting into this type of development. I thought I was able to install the latest nightly build today, but I'm not sure. Running Mac OS X. If anyone has any pointers I would really appreciate them off-list. Thanks, Matt Matt, I never tried it as I always (re)build from SVN but that implies having previously built the dependencies, but homebrew (no Fink, no MacPorts, HOMEBREW) has a GDAL formula that contemplates svn build. Maybe it's worth to give it a try. Joaquim ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Set color table in NITF
Hello, I am trying to create a NITF that uses a color lookup table. I found the option to add an RGB/LUT to the image, which didnt complain, however, when I had my LUT data to the image, it all seems to be kind of black and white. I know when I was displaying the images in my viewer (which I wrote in Qt), I had a function to create a 256 element color table, that changed some of the values. How would I modify the color table in the NITF? Thanks ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Set color table in NITF
Le mardi 30 août 2011 00:14:12, Cole, Derek a écrit :
Hello,
I am trying to create a NITF that uses a color lookup table.
I found the option to add an RGB/LUT to the image, which didnt complain,
however, when I had my LUT data to the image, it all seems to be kind of
black and white.
I know when I was displaying the images in my viewer (which I wrote in Qt),
I had a function to create a 256 element color table, that changed some of
the values. How would I modify the color table in the NITF?
Use SetColorTable() on the raster band object.
Python example :
from osgeo import gdal
ds = gdal.GetDriverByName('NITF').Create('test.ntf', 1, 1, options =
['IREP=RGB/LUT', 'LUT_SIZE=3'])
ct = gdal.ColorTable()
ct.SetColorEntry( 0, (0,0,0) )
ct.SetColorEntry( 1, (10,10,10) )
ct.SetColorEntry( 2, (20, 20, 20))
ds.GetRasterBand(1).SetColorTable(ct)
ds = None
$ gdalinfo test.ntf
Driver: NITF/National Imagery Transmission Format
Files: test.ntf
Size is 1, 1
[...]
Band 1 Block=1x1 Type=Byte, ColorInterp=Palette
NoData Value=3
Color Table (RGB with 4 entries)
0: 0,0,0,255
1: 10,10,10,255
2: 20,20,20,255
3: 0,0,0,0
Thanks
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Projection definition issue with Mapinfo files
Le lundi 29 août 2011 13:34:22, Yves Jacolin a écrit : Hello, One of QGIS user report an issue with mapInfo file (TAB file). If you understand french you can see the thread here: http://georezo.net/forum/viewtopic.php?pid=197022 This looks really similar to the issue raised in ticket http://trac.osgeo.org/gdal/ticket/3292 that you reported almost 2 years ago ;-) but that did not see any progress unfortunately :-( Although not a specialist of the mapinfo format, I think the fix for EPSG:2154 should be straightforward (adding an entry in a table) once we know the mapinfo code that matches to the Reseau_Geodesique_Francais_1993 datum. One way to get it would be to attach the COMMUNEL93.TAB (and the other related files that make dataset) mentionned in http://georezo.net/forum/viewtopic.php?pid=197022 to the ticket. I try to sum up the issue in english in this email. The process is basicly as follow: 1. check projection information with ogrinfo DATUM[GRS_80, SPHEROID[GRS 80,6378137,298.257222101], 2. update projection information file with ogr2ogr: ogr2ogr -f MapInfo file -s_srs EPSG:2154 -t_srs EPSG:2154 new.tab origin.tab 3. check projection information with ogrinfo DATUM[WGS_1984, SPHEROID[WGS 84,6378137,298.257223563], 1 and 3 show differents projection information! If I try to change file format to shp, it seems that projection information are correct. I can share tab files if needed. Any idea? Thanks. Y. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Exporting functions under Windows
Hi, I came across this problem before, and I don't remember id I solved. The problem is then gdal18.dll is not exporting many functions. I use default build using VC 2010 express edition. For instance instead of function GDALAllRegister library exports _GDALAllRegister@0 GetProcAddress() fails in first case, but works for second. I am using C interface. Am I missing something very basic? Thank you Mikhail ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Exporting functions under Windows
I think I figured it out. I added -DCPL_DISABLE_STDCALL to the OPTFLAGS and commented out STDCALL=YES After this it all worked. Thank you Mikhail On 8/29/2011 5:29 PM, Mikhail Tchernychev wrote: Hi, I came across this problem before, and I don't remember id I solved. The problem is then gdal18.dll is not exporting many functions. I use default build using VC 2010 express edition. For instance instead of function GDALAllRegister library exports _GDALAllRegister@0 GetProcAddress() fails in first case, but works for second. I am using C interface. Am I missing something very basic? Thank you Mikhail ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
