Hi all,

I think I found what might be a bug while trying to optimize some code
that uses Shapely.  In my program's innermost loop I create a Polygon
object and pass it to the 'intersects' method of another geometry.
When I profiled my code, I found that it was spending a lot of time in
Polygon.__init__, so I've been looking for a way to just set the
coordinates that make up the polygon instead of recreating the object
each time.  Which would make sense because the polygons in this inner
loop always have the same dimensionality and number of vertices- only
the coordinate values change.

Anyway, I'm using Shapely 1.2.7.  I'm wondering why it doesn't work to
do something like this:

>>> p1=Polygon(((-1, -1), (1, -1), (0, 0)))
>>> p2=Polygon(((-1, 1), (1, 1), (1, 2), (-1, 2)))
>>> p2.intersects(p1)
False
>>> list(p1.exterior.coords)
[(-1.0, -1.0), (1.0, -1.0), (0.0, 0.0), (-1.0, -1.0)]
>>> cs=lgeos.GEOSGeom_getCoordSeq(p1.exterior.__geom__)
>>> lgeos.GEOSCoordSeq_setX(cs, 2, 0)
1
>>> lgeos.GEOSCoordSeq_setY(cs, 2, 2)
1
>>> list(p1.exterior.coords)
[(-1.0, -1.0), (1.0, -1.0), (0.0, 2.0), (-1.0, -1.0)]   # okay, looks
like that worked
>>> p2.intersects(p1)
False    # but this is not the expected result
>>> p1 = Polygon(p1.exterior)   # rebuild the polygon
>>> p2.intersects(p1)    # you sure about that?
True    # ah ha!

I actually got this CoordSeq_setX/Y trick to work in an older version
of Shapely (1.2.3), but I found it would stop working the first time I
put a new object reference in p1, or if I fork()ed!  Which seems
bizarre.  In both cases the libgeos was 3.2.2.

I realize I'm messing around with interfaces I probably shouldn't be,
so forgive me for asking for help with a nasty hack!  But when I did
have it working with the older Shapely version, it resulted in a very
significant speedup for my application.  So if anyone knows what
causes this or maybe a better way to achieve the same result, I'd
really appreciate the advice.

Thanks,
Kyle Cronan
_______________________________________________
Community mailing list
[email protected]
http://lists.gispython.org/mailman/listinfo/community

Reply via email to