[gdal-dev] OGR Geometry methods

Chaitanya kumar CH chaitanya.ch at gmail.com
Tue Jul 12 02:07:45 EDT 2011


Marius,

Please do file a ticket with a small shapefile that shows this error.
http://trac.osgeo.org/gdal/newticket

On Tue, Jul 12, 2011 at 8:30 AM, Marius Jigmond
<mariusjigmond at hotmail.com>wrote:

> **
> 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:
>
> #!/usr/bin/python
>
> from osgeo import ogr
> import sys, math, time, os
>
> aquifer = '/data/romania/judeteEPSG3844.shp'
> aqDS = ogr.Open(aquifer)
> sql = 'select * from judeteEPSG3844 where DENJUD = "CLUJ"'
> aqLayer = aqDS.ExecuteSQL(sql)
> feat = aqLayer.GetFeature(0)
> print aqLayer.GetFeatureCount()
> print feat.GetField('DENJUD')
> aqDS.ReleaseResultSet(aqLayer)
>
> The result:
> marius at mobi:~/temp$ run_sql.py
> 1
> TELEORMAN
> marius at mobi:~/temp$
>
> This explains why none of my centroids were remotely close to the polygon.
>
> Is there something wrong with my query or should I file a bug?
>
> -marius
>
>
> On Mon, 2011-07-11 at 20:46 -0500, Marius Jigmond wrote:
>
> I suppose a piece of code speaks louder :):
>
> xsect = False
> for i in range(gridLayer.GetFeatureCount()):
>   feat = gridLayer.GetFeature(i)
>   geom = feat.GetGeometryRef()
>   point = geom.Centroid()
>   for j in range(aqLayer.GetFeatureCount()):
>     aqfeat = aqLayer.GetFeature(j)
>     aqgeom = aqfeat.GetGeometryRef()
>     if point.Intersect(aqgeom):
>       xsect = True
>       print 'ibound = 1'
>       break
>   if xsect:
>     feat.SetField('IBOUND', 1)
>     gridLayer.SetFeature(feat)
>     xsect = False
>
> -marius
>
> On Mon, 2011-07-11 at 20:34 -0500, Marius Jigmond wrote:
>
> Hi everyone,
>
> 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.
>
> -marius
>
> _______________________________________________
> gdal-dev mailing listgdal-dev at lists.osgeo.orghttp://lists.osgeo.org/mailman/listinfo/gdal-dev
>
>  _______________________________________________
> gdal-dev mailing listgdal-dev at lists.osgeo.orghttp://lists.osgeo.org/mailman/listinfo/gdal-dev
>
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>



-- 
Best regards,
Chaitanya kumar CH.

+91-9494447584
17.2416N 80.1426E
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.osgeo.org/pipermail/gdal-dev/attachments/20110712/d0ae8fa8/attachment.html


More information about the gdal-dev mailing list