<html>
<head>
<style><!--
.hmmessage P
{
margin:0px;
padding:0px
}
body.hmmessage
{
font-size: 10pt;
font-family:Tahoma
}
--></style>
</head>
<body class='hmmessage'><div dir='ltr'>
Strangely enough, ExecuteSQL works (it returns an osgeo.ogr.layer object with one feature) with<br>
<br>
'select FID from judeteEPSG3844 where DENJUD = "CLUJ"'<br>
<br>
but later a GetField('FID') request returns "Illegal field requested in GetField()"<br>
<br>
-marius<br>
<br><br><div><hr id="stopSpelling">Date: Tue, 12 Jul 2011 11:36:42 -0400<br>Subject: Re: [gdal-dev] OGR Geometry methods<br>From: luke.peterson@gmail.com<br>To: mariusjigmond@hotmail.com<br><br>does executing<div><br></div><div>sql = 'select FID from judeteEPSG3844 where DENJUD = "CLUJ"'</div><div><br></div><div>return anything?<br clear="all">-----<br>Luke Peterson<br>312-533-1051(m)<br>
<br><br><div class="ecxgmail_quote">On Tue, Jul 12, 2011 at 11:30 AM, Marius Jigmond <span dir="ltr"><<a href="mailto:mariusjigmond@hotmail.com">mariusjigmond@hotmail.com</a>></span> wrote:<br><blockquote class="ecxgmail_quote" style="border-left:1px #ccc solid;padding-left:1ex">
<div><div dir="ltr">
Luke,<br>
<br>
GetFIDColumn returns an empty string (not supported according to the docs). I tried the following SQL statement:<br>
sql = 'select DENJUD from judeteEPSG3844 where DENJUD = "CLUJ"'<br>
but it's still returning that same wrong feature.<br><font color="#888888">
<br>
-marius<br>
</font><div class="ecxhm"><br><br></div><div><div class="ecxhm"><hr>Date: Tue, 12 Jul 2011 09:13:00 -0400<br>Subject: Re: [gdal-dev] OGR Geometry methods<br>From: <a href="mailto:luke.peterson@gmail.com">luke.peterson@gmail.com</a><br>
To: <a href="mailto:mariusjigmond@hotmail.com">mariusjigmond@hotmail.com</a></div><div><div></div><div class="h5"><br><br><br>
-----<br>
Luke Peterson<br>
- sent via mobile device -<br>
On Jul 11, 2011 11:00 PM, "Marius Jigmond" <<a href="mailto:mariusjigmond@hotmail.com">mariusjigmond@hotmail.com</a>> wrote:<br>
><br>
> After some more investigation that is likely NOT the issue. I have an ExecuteSQL statement which selects a certain polygon based on an attribute value. Unfortunately it seems to return the wrong feature. The feature I query for is unique so a duplicate is out of the question. Here's the code:<br>
><br>
> #!/usr/bin/python<br>
><br>
> from osgeo import ogr<br>
> import sys, math, time, os<br>
><br>
> aquifer = '/data/romania/judeteEPSG3844.shp'<br>
> aqDS = ogr.Open(aquifer)<br>
> sql = 'select * from judeteEPSG3844 where DENJUD = "CLUJ"'<br>
> aqLayer = aqDS.ExecuteSQL(sql)<br>
> feat = aqLayer.GetFeature(0)<br>
> print aqLayer.GetFeatureCount()<br>
> print feat.GetField('DENJUD')<br>
> aqDS.ReleaseResultSet(aqLayer)<br>
><br>
> The result:<br>
> marius@mobi:~/temp$ run_sql.py <br>
> 1<br>
> TELEORMAN<br>
> marius@mobi:~/temp$ <br>
><br>
> This explains why none of my centroids were remotely close to the polygon.<br>
><br>
> Is there something wrong with my query or should I file a bug?<br>
><br>
Can you select the feature by its FID instead? (apologies, I'm banging this out on my phone ... probably won't run but you should get the idea):<br>
aqLayer = aqDS.GetLayerByIndex(0) #gets the whole shapefile<br>
aqFID = aqLayer.GetFIDColumn() # I forget the method name here... hope this is right but may have to pull in the layer def<br>
sql = 'select %s from judeteEPSG3844 where DENJUD = "CLUJ"' % aqFID # specifically asking for just the FID field in the resultset<br>
sqlResults = aqDS.ExecuteSQL(sql) # sql results same as before<br>
featFID = sqlResults.GetFeature(0).GetFieldAsInt(aqFID) #this should be the actual FID of the results (problems would occur if your sql returned more than one hit)<br>
feat = aqLayer.GetFeature(featFID) # using the fid to grab the feature from the 0th layer in the shapefile<br>
This may be a long way around ... <br><br>
> -marius<br>
><br>
><br>
> On Mon, 2011-07-11 at 20:46 -0500, Marius Jigmond wrote:<br>
>><br>
>> I suppose a piece of code speaks louder :):<br>
>><br>
>> xsect = False<br>
>> for i in range(gridLayer.GetFeatureCount()):<br>
>> feat = gridLayer.GetFeature(i)<br>
>> geom = feat.GetGeometryRef()<br>
>> point = geom.Centroid()<br>
>> for j in range(aqLayer.GetFeatureCount()):<br>
>> aqfeat = aqLayer.GetFeature(j)<br>
>> aqgeom = aqfeat.GetGeometryRef()<br>
>> if point.Intersect(aqgeom):<br>
>> xsect = True<br>
>> print 'ibound = 1'<br>
>> break<br>
>> if xsect:<br>
>> feat.SetField('IBOUND', 1)<br>
>> gridLayer.SetFeature(feat)<br>
>> xsect = False<br>
>><br>
>> -marius<br>
>><br>
>> On Mon, 2011-07-11 at 20:34 -0500, Marius Jigmond wrote:<br>
>>><br>
>>> Hi everyone,<br>
>>><br>
>>> I am trying to test whether centroids of polygons lie/intersect within another polygon. I have tried Intersect, Within, and Contains but they always return false. Should these methods work for my intended purpose or do I need to implement a point in polygon function? Thanks.<br>
>>><br>
>>> -marius<br>
>>><br>
>>> _______________________________________________<br>
>>> gdal-dev mailing list<br>
>>> <a href="mailto:gdal-dev@lists.osgeo.org">gdal-dev@lists.osgeo.org</a><br>
>>> <a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/mailman/listinfo/gdal-dev</a><br>
>><br>
>> _______________________________________________<br>
>> gdal-dev mailing list<br>
>> <a href="mailto:gdal-dev@lists.osgeo.org">gdal-dev@lists.osgeo.org</a><br>
>> <a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/mailman/listinfo/gdal-dev</a><br>
><br>
><br>
> _______________________________________________<br>
> gdal-dev mailing list<br>
> <a href="mailto:gdal-dev@lists.osgeo.org">gdal-dev@lists.osgeo.org</a><br>
> <a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/mailman/listinfo/gdal-dev</a><br>
<br></div></div></div> </div></div>
</blockquote></div><br></div></div> </div></body>
</html>