[Gdal-dev] -spat switch in ogr2ogr

[Frank Warmerdam]

Brian Timoney wrote:
> Frank Warmerdam wrote:
>> Brian Timoney wrote:
>>>
>>> All:
>>>
>>> At the risk of disproving the adage that the only dumb question is 
>>> the one left unasked, I'm having some trouble applying the -spat 
>>> switch using ogr2ogr from the command line.
>>>
>>> I have a table of US cities in PostGIS called..."cities" in lat/long 
>>> NAD 83 (EPSG:4269)
>>>
>>> I want to pull cities from the middle of the country from a BBOX of  
>>> (-114,30,-100,45) and throw them into a new shapefile.
>>>
>>> So I hazarded...
>>>
>>> ogr2ogr mid_city.shp PG:'dbname=mydb user=me' cities -spat '-114 30 
>>> -100 45'
>>>
>>> No dice (empty shapefile).  There are 75+ features that fall inside 
>>> the bounding box.
>>>
>>> I've also tried double quotes and no quotes around the xmin ymin xmax 
>>> ymax values.
>>
>> Brian,
>>
>> It should work with:
>>
>>  ogr2ogr mid_city.shp PG:'dbname=mydb user=me' cities -spat -114 30 
>> -100 45
>>
>> Would it be fair for me to guess that the following works fine?
>>
>>  ogr2ogr mid_city.shp PG:'dbname=mydb user=me' cities
>>
>> So if my suggested approach (no quotes on -spat argument) doesn't 
>> work, then
>> I'm at a bit of a loss.  I guess I'd have to actually test things at 
>> my end!
>>
>> By the way, you might want to set the CPL_DEBUG environment variable 
>> to ON.
>> This will produce various debug output of the actual SQL being sent to
>> postgis.  This might give clues.
>>
>> Best regards,
> 
> Frank:
> 
> OK, switching the debugging on--
> 
> $ CPL_DEBUG=ON ogr2ogr mid_city.shp PG:'dbname=mydb user=me' cities 
> -spat -114 30 -100 45
> 
> Gave me what ogr2ogr is passing to PostGIS.  The key part is the WHERE 
> statement--
> 
> WHERE the_geom && GeometryFromText('BOX3D(-114.000000 
> 30.000000,-100.00000 45.000000)'::box3d,-2))
> 
> I fooled around replicating the WHERE statement in PgAdmin's SQL window 
> and ran into the following issues--
> 
> 1)  The GeometryFromText function no longer honors 'BOX3D'
> 
> http://postgis.refractions.net/pipermail/postgis-users/2005-July/008631.html 
> 
> 
> 2)  Following the suggestion above, I tried
> 
> WHERE the_geom && SetSRID('BOX3D(-114 30,-100 45)::box3d,-1)
> 
> and received an error for "Operation on two geometries with different 
> SRIDs"
> 
> 3)  Since my PostGIS table has an SRID of 4269, I tried
> 
> WHERE the_geom && SetSRID('BOX3D(-114 30,-100 45)::box3d,4269)
> 
> and everything worked swimmingly....

Brian,

OK, there are two issues here.  One is that the SRID was wrong.  The value
-2 is a special marker kept around till the SRID can be looked up.  The
looking up was broken (recently?) and this may have been the cause of that
problem.  the other issue is use of SetSRID(...) and that seems to have been
fixed in the development version of the code too, though I can't see exactly
when from the log messages.

An ogrinfo with -spat for me does:
OGR_PG: PQexec(DECLARE OGRPGLayerReader CURSOR for SELECT "ogc_fid", 
AsEWKT("wkb_geometry"), "city", "id", "mapid", "milltype", "mslink", "owner", 
"state", "yrly_capacity" FROM "Mills" WHERE wkb_geometry && 
SetSRID('BOX3D(831000.000000000000 3525000.000000000000, 832000.000000000000 
3526000.000000000000)'::box3d,32767) )

and works right.

Would you possibly be willing to try a GDAL night cvs snapshot and see if it
works properly for you?  If not I'd appreciate your filing a bug and I can
look into it more deeply.

Best regards,
-- 
---------------------------------------+--------------------------------------
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | President OSGeo, http://osgeo.org