Python gdal stopped working when using ogr Within, Contains or other

Question

Update: By further ivestigation I found out it must be caused by somehow corrupted geometry. But even if I run:

if drillhole[1].IsValid():

It is causing the crash. So not sure how to check it.

I am trying to test if lines are within polygons. Using shapely it works for me very well except of speed - I have ten thousands of lines and the same for polygons. I just wanted to test if ogr can do this faster but no luck to make it work.

import ogr
#Load cells as polygons with slight buffer
data_source = ogr.Open(file_path_cells)
source_layer = data_source.GetLayer()
source_layer.ResetReading()
cells = []
for source_feature in source_layer:
  feature_id = source_feature.GetFieldAsInteger(0)
  feature_geometry = source_feature.geometry()
  feature_geometry = feature_geometry.Buffer(0.1, quadsecs = 3)
  cells.append((feature_id,feature_geometry))

#Load drillholes as lines wihtin cells
data_source = ogr.Open(file_path_drillholes)
source_layer = data_source.GetLayer()
source_layer.ResetReading()
drillholes = []
for source_feature in source_layer:
  feature_elev = source_feature.GetFieldAsInteger("elev")
  feature_geometry = source_feature.geometry()
  drillholes.append((feature_elev,feature_geometry))  
for cell in cells:
  for drillhole in drillholes:
    if cell[1].Contains(drillhole[1]):
      print("yes")
    else:
      print("no")

Any idea what can be wrong with the line?:

if cell[1].Contains(drillhole[1]):

On both my Windows 7 machines I allways get python.exe stopped working... Problem signature:

  Problem Event Name:   APPCRASH
  Application Name: python.exe
  Application Version:  0.0.0.0
  Application Timestamp:    5193f3af
  Fault Module Name:    gdal111.dll
  Fault Module Version: 1.11.2.0
  Fault Module Timestamp:   54e65215
  Exception Code:   c0000005
  Exception Offset: 00000000005e5fb3
  OS Version:   6.1.7601.2.1.0.768.3

Or is there some other faster way or way how to improve the speed of this using shapely?

Answer 1

The crash is most likely a well documented gotcha.

Similar work can be done with Shapely, which doesn't have that gotcha. Also, the query in the for loop can be done faster by using an Rtree index, there are a few Q&A's related to shapely + Rtree out there.