[gdal-dev] Intersection of two Layers

[Kai Muehlbauer]

Am 22.01.2016 um 12:23 schrieb Ari Jolma:
> 22.01.2016, 12:48, Kai Muehlbauer kirjoitti:
>>
>>
>> Am 22.01.2016 um 10:52 schrieb Ari Jolma:
>>> 22.01.2016, 10:06, Kai Muehlbauer kirjoitti:
>>>>
>>>>
>>>> To quickly check if there is an improvement I used the python shapely
>>>> package. Converting the outer loop geometry to an shapely prepared
>>>> geometry and calling 'contains()' with the inner loop geometries to
>>>> find those fully contained geometries I get an improvement of nearly
>>>> 40% to the former ogr only code.
>>>>
>>>> As Even pointed out, there is only OGRPreparedGeometryIntersects()
>>>> (GEOSPreparedIntersects) exposed within ogrgeometry.cpp. In my use
>>>> case I would need GEOSPreparedContains to mimic the shapely behaviour.
>>>
>>> I made a small test with the attached diff, which prepares the
>>> geometries of the input layer.
>>>
>>> The input layer had two polygons (boundaries of two Finnish
>>> municipalities) and the method layer had 3091 points, all in either one
>>> of the polygons.
>>>
>>> The improvement is ~15%.
>>>
>>> I think the prepared geometries could be used by default easily.
>>> Prepared contains test would need addition of the respective GEOS
>>> method.
>>
>> Thanks Ari, I very much appreciate your efforts in testing layer
>> intersection. If the prepared contains test could be done within the
>> layer intersection then those features can "just be copied" without
>> the need of costly intersection calculation. This would give a real
>> speedup.
>
> It seems so, I added the contains test and in this case the total
> speedup is ~80 %.
>
>>
>> So, it would need the addition of the OGRPreparedGeometryContains
>> method (GEOSPreparedContains) and then enhancing your patch.
>> Would it work the way as in the attached diff? Or did I miss something?
>
> Attached is my diff, which basically copies to add-to-result code, which
> might be done better.
>
>>
>> What are the next steps to get this feature fully implemented?
>
> The thing is that most of the other layer methods could be enhanced this
> way, which means some more work. I believe the C++ API can be enhanced
> for version 2.1.

This is very good news! If I understand correctly, I could use your 
patch and self-compile gdal to have access to the improved 
layer-intersection function via the python bindings. I'll give this a try.

I also created a ticket here: https://trac.osgeo.org/gdal/ticket/6323 to 
keep track of this enhancement.

Cheers,
Kai