Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I have two feature layers. The first one is just a polygon and the other one consists of a series of polygons which are enclosed by the first one. Thus, some polygons of the second feature layer touch on the boundaries of the first feauture layer. I wrotre a code (the basic part shown below to identify which polygons of the second layer touches the boundary of the polygon in the first layer but while the code runs without a problem it does not identify (e.g. by giving FID with Msgbox) any such polygon whilsty it should. Any help please?

Dim pBlockBoundaryFC  As IFeatureClass
Set pBlockBoundaryFC = pBlockBoundary.FeatureClass


Dim pBlockBoundaryFields  As IFields
Set pBlockBoundaryFields = pBlockBoundaryFC.Fields

Dim pGeom As IGeometry


Dim pSpatialFilter As ISpatialFilter
Set pSpatialFilter = New SpatialFilter

Dim pBlockBoundaryCursor As IFeatureCursor
Set pBlockBoundaryCursor = pBlockBoundaryFC.Search(pSpatialFilter, False)


Set pBlockBoundaryFeature = pBlockBoundaryCursor.NextFeature
Do Until pBlockBoundaryFeature Is Nothing

Set pGeom = pBlockBoundaryFeature.Shape
With pSpatialFilter

'************** 'Conversion Dim pPolygonBoundary As ITopologicalOperator Set pPolygonBoundary = pFeature.Shape

Set pGeom = pPolygonBoundary.Boundary '*************

Set .Geometry = pGeom
    .GeometryField = "Shape"
    .SpatialRel = esriSpatialRelTouches
End With

Dim pRelOp As IRelationalOperator
Set pRelOp = pGeom

  Set pNewParcelsFeature = pNewParcelsCursor.NextFeature
  Do Until pNewParcelsFeature Is Nothing

    FID = pNewParcelsFeature.Value(intposFID)
    If pRelOp.Touches(pNewParcelsFeature.Shape) Then

    ''Report the FID of each parcel that 'touches' the exetrnal boundary
    MsgBox FID
    End If

Set pNewParcelsFeature = pNewParcelsCursor.NextFeature

Loop


Set pBlockBoundaryFeature = pBlockBoundaryCursor.NextFeature

Loop
share|improve this question

2 Answers 2

up vote 5 down vote accepted

There are a few concepts that I see you are missing. I will try to explain them.

1- The first one has to do with spatial filter usage. You should pass a filter when you are using it, otherwise, just pass nothing.

You have this code:

Dim pSpatialFilter As ISpatialFilter
Set pSpatialFilter = New SpatialFilter

Dim pBlockBoundaryCursor As IFeatureCursor
Set pBlockBoundaryCursor = pBlockBoundaryFC.Search(pSpatialFilter, False)

Why do this? It seems you want to pass Nothing to loop through every single feature (and by the way, you should make that False a True) You want to use a recycling cursor. But that is a different discussion.

2.- You are modifying the spatial filter's geometry while you are looping. That is the reason your results are all out of whack. Do not ever ever ever ever do that, unless the geometry ownership has been transfered using the FilterOwnsGeometry property.

There are ways to make this much much faster, but a simpler (non-optimal) version of what you want is (in pseudo code):

cursor = boundaryFC.search(Nothing, True)
for every feature returned by cursor
   create spatial filter
   spatialfilter.geometry =  feature.shape
   spatialfilter.SpatialRel = esriSpatialRelTouches

   secondCursor = parcelFC.search(spatialfilter)
   for every parcelFeature returned by secondCursor
     debug.print parcelFeature.objectid
   next parcelFeature

next feature

That should do the trick. There are ways to make this faster, but it will do for now.

share|improve this answer
    
Thank you very much for the remarks and the pseudocoe. –  Demetris Jan 31 '12 at 5:56
    
What Does FilterOwnsGeometry even do? –  matrixugly Oct 15 at 4:38
    
havk for internal memory management. If the filter owns the geometry, then in each loop it can assume it won't change so you can save a deep clone of it in the internals. A classic case where this happens is when you use the geometry of a feature that came from a recycling cursor. On every iteration that geometry will change. Honestly, IMHO it should just always deep clone and forget about it. –  Ragi Yaser Burhum Nov 1 at 18:30

The borders have to match exactly in order to get results from IRelationalOperator::Touches. From the ArcGIS documentation:

Two geometries are said to touch when the intersection of the geometries is non-empty, but the intersection of their interiors is empty. Touches is a Clementini operator. For example, a point touches a polyline only if the point is coincident with one of the polyline end points.

I would convert the geometry from the first layer into a line (ITopologicalOperator::Boundary) and then use the crosses operator (IRelationalOperator::Crosses) to check if they intersect.

share|improve this answer
    
I tried it but it still does not work. However, I am not sure if I did Ok. Please, have a look in the code I added between the asterisks. Thanks –  Demetris Jan 30 '12 at 12:25
    
I think you should replace "If pRelOp.Touches(pNewParcelsFeature.Shape) Then" with "If pRelOp.Crosses(pNewParcelsFeature.Shape) Then" –  Samuel Jan 30 '12 at 16:16
    
Eventually, the 'touches' worked! –  Demetris Jan 31 '12 at 6:44

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.