[gdal-dev] Re: Intersection of polygon layers

[Simon Lyngby Kokkendorff]

Hi Yosuke,

  If you use MultiPolygon geometries (with your polygons as the parts), geos
(via ogr or shapely) should be able to do the job in one go, which should
speed things up considerably...

  Cheers,
   Simon


On Sat, Nov 6, 2010 at 5:04 AM, yosuke kimura <yosukesabai at gmail.com> wrote:

> All,
>
> I think I found what I should have been using.  pick a feature from
> one layer, use its geometry as an argument of SetSpatialFilter method
> of the other layer.  Seems to be yielding same answer quicker...
>
> Francis,
>
> Thank you for suggestion.  I looked at the Shapely documentation
> http://gispython.org/shapely/docs/1.0/manual.html . It wasn't obvious
> to me advantage of the library over ogr, except that it looks more
> pythonic.  I even didn't find mention of layer vs. feature, and to
> solve my problem it is a step back to me?
>
>
>
> On Nov 5, 4:54 pm, Francis Markham <u2546... at anu.edu.au> wrote:
> > If you're using Python, I'd recommend checking out the shapely package (
> http://pypi.python.org/pypi/Shapely):
> >
> > You can load your shapefile using ogr, and then import the geometry to
> > shapely where you can do the geometry predicate calculations:
> >
> > >>> import ogr
> > >>> from shapely.wkb import loads
> > >>> source = ogr.Open("/tmp/world_borders.shp")
> > >>> borders = source.GetLayerByName("world_borders")
> > >>> feature = borders.GetNextFeature()
> > >>> loads(feature.GetGeometryRef().ExportToWkb())
> >
> > Good luck,
> >
> > Francis
> >
> > On 6 November 2010 06:08, yosuke kimura <yosukesa... at gmail.com> wrote:
> >
> > > All,
> >
> > > I am trying to calculate intersection of two polygon layers with
> > > polygons.  I searched threads, found following thread, and what I got
> > > was that OGR/GEOS provide me to intersection of two features, but I
> > > have to write something in order to perform this operation for two
> > > layers with multiple features.
> >
> > > Open source vector geoprocessing libraries?
> > >http://groups.google.com/group/gdal/browse_thread/thread/7b5a4f461a4d.
> ..
> >
> > > Intersect of two shapefiles
> > >http://groups.google.com/group/gdal/browse_thread/thread/4b178c776aec.
> ..
> >
> > > My questions are:
> >
> > > 1. Did I miss something?  Is there development of library which takes
> > > two polygon layers to calculate intersection or other basic geometry
> > > operations that i find in OGR?
> >
> > > 2. Assuming no to above, I got started writing code using python
> > > binding.  Not having too much experience, I codes to go through all
> > > possible pair (complete bipartite), test with intersect, then run
> > > intersection if it returns true.  I am wondering there is faster way
> > > to do this.  I am thinking something like get envelope polygon of each
> > > feature first, and then test intersection of them and then perform
> > > intersect operation only when first test passes (i assume that it is
> > > easier to test intersect if there are only handful of vertices as
> > > oppose to hundreds of vertices in my case).
> >
> > > Right now I am trying with layer of hundreds of thousands of polygons
> > > vs. thousands of polygon, and it seems to be taking 10s of hours to
> > > finish in my environment.
> >
> > > Thank you,
> >
> > > Yosuke Kimura
> > > _______________________________________________
> > > gdal-dev mailing list
> > > gdal-... at lists.osgeo.org
> > >http://lists.osgeo.org/mailman/listinfo/gdal-dev
> >
> >
> >
> > _______________________________________________
> > gdal-dev mailing list
> > gdal-... 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
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.osgeo.org/pipermail/gdal-dev/attachments/20101109/1bdf0fa9/attachment-0001.html