How to execute spatial queries

Summary A spatial query is a query that returns features based on their spatial relationship with a query geometry. This includes searching for features given an extent (i.e. finding all of the streets within an envelope) or searching for features based on their relation to other features (i.e. finding all farmland intersected by a particular stream). The article "How to query geodatabase tables and feature classes" introduces spatial querying – this article expands upon that information with illustrated examples, code samples, and detailed information concerning optimizing queries with spatial caching.

Development licensing Deployment licensing
ArcView ArcView
ArcEditor ArcEditor
ArcInfo ArcInfo
Engine Developer Kit Engine Runtime

To use the code in this topic, reference the following namespaces via the using (C#) or Imports (VB .NET) statements. It is also necessary to add the corresponding references to the project to gain access to these application programming interfaces (APIs).

In this topic

Approaches to executing spatial queries on features and feature classes

There are two different APIs available to a developer who needs to query the spatial attributes of a feature class and the spatial relationships between features in a feature class.
 
 
 
 
 
 
Both types of spatial queries offer advantages in different applications. In general, relational operators are ideal for discrete geometry-on-geometry comparisons where the features being compared are known up front. On the other hand, spatial filters are good when working within the larger scope of a feature class and there is little to no information about the input features prior to executing the spatial query.
 

Prerequisites for code examples in this article
The following considerations apply to all code examples this article:
 

Examples of spatial queries

This section provides several illustrated examples of spatial queries along with code samples, using the ISpatialFilter and IRelationalOperator interfaces. An introduction to using the ISpatialFilter interface – as well as the IQueryFilter interface, which the ISpatialFilter interface extends – can be found in the article, How to query geodatabase tables and feature classes. The examples shown involve the following scenarios:
 
Find features within a polygon
This example shows how to execute a spatial query that will find polyline features that intersect a polygon – specifically, major highways that pass through Indiana.
 
The illustration below shows several major highways running through Iowa:
 
 
 
The first step is to retrieve the geometry for the state of Iowa, assuming that the Object ID of the feature is known. Once the feature has been found, its geometry can be used as the query geometryin a spatial filter. Since the query should locate all highways that pass through it, and not necessarily those contained entirely within, the spatial filter’s relationship should be set to “intersects”.
 
The following code shows how to use the state’s geometry as a query geometry, construct a spatial filter using it, execute the query with the spatial filter, then iterate through the results, displaying the name of each highway:
 

[C#] 
// Get the feature and its geometry given an Object ID.
IFeature stateFeature = stateFeatureClass.GetFeature(14);
IGeometry queryGeometry = stateFeature.Shape;

// Create the spatial filter. "highwayFeatureClass" is the feature class containing
// the highway data. Set the SubFields property to "FULL_NAME", as only that field is
// going to be displayed.
ISpatialFilter spatialFilter = new SpatialFilterClass();
spatialFilter.Geometry = queryGeometry;
spatialFilter.GeometryField = highwayFeatureClass.ShapeFieldName;
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelIntersects;
spatialFilter.SubFields = "FULL_NAME";

// Find the position of the "FULL_NAME" field in the highway feature class.
int nameFieldPosition = highwayFeatureClass.FindField("FULL_NAME");

// Execute the query and iterate through the cursor's results.
IFeatureCursor highwayCursor = highwayFeatureClass.Search(spatialFilter, false);
IFeature highwayFeature = null;
while ((highwayFeature = highwayCursor.NextFeature()) != null)
{
  String name = Convert.ToString(highwayFeature.get_Value(nameFieldPosition));
  Console.WriteLine("Highway found: {0}", name);
}

// The cursors is no longer needed, so dispose of it.
Marshal.ReleaseComObject(highwayCursor);

[VB.NET] 
' Get the feature and its geometry given an Object ID.
Dim stateFeature As IFeature = stateFeatureClass.GetFeature(14)
Dim queryGeometry As IGeometry = stateFeature.Shape

' Create the spatial filter. "highwayFeatureClass" is the feature class containing
' the highway data. Set the SubFields property to "FULL_NAME", as only that field is
' going to be displayed.
Dim spatialFilter As ISpatialFilter = New SpatialFilterClass()
spatialFilter.Geometry = queryGeometry
spatialFilter.GeometryField = highwayFeatureClass.ShapeFieldName
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelIntersects
spatialFilter.SubFields = "FULL_NAME"

' Find the position of the "FULL_NAME" field in the highway feature class.
Dim nameFieldPosition As Integer = highwayFeatureClass.FindField("FULL_NAME")

' Execute the query and iterate through the cursor's results.
Dim highwayCursor As IFeatureCursor = highwayFeatureClass.Search(spatialFilter, False)
Dim highwayFeature As IFeature = highwayCursor.NextFeature()
While Not highwayFeature Is Nothing
    Dim Name As String = Convert.ToString(highwayFeature.Value(nameFieldPosition))
    Console.WriteLine("Highway found: {0}", Name)
    highwayFeature = highwayCursor.NextFeature()
End While

' The cursors is no longer needed, so dispose of it.
Marshal.ReleaseComObject(highwayCursor)
 
Query for features with relation to multiple geometries
This example shows how to use a geometry bag as a query geometry. A geometry bag is a single high-level geometry that stores a collection of geometries. The use of a query bag as the query geometry in a single query is an alternative to performing multiple queries for each geometry. This example shows how to find parcel features within a known set of block features. The block features are not adjacent, but their Object IDs are known.
 
The image below shows the blocks and the parcels for our search area – the spatial query should be limited to finding the parcels within the blue highlighted blocks.
 
 
 
The first step in the process is to create a new geometry bag and fill it with the geometries of the three block features. For this example, it’s assumed the Object IDs of the blocks are known. See the following code example:
 

[C#] 
// Create a new geometry bag and give it the same spatial reference as the
// blocks feature class.
IGeometryBag geometryBag = new GeometryBagClass();
IGeometryCollection geometryCollection = (IGeometryCollection)geometryBag;
IGeoDataset geoDataset = (IGeoDataset)blocksFeatureClass;
ISpatialReference spatialReference = geoDataset.SpatialReference;
geometryBag.SpatialReference = spatialReference;

// Get a feature cursor for the three blocks and put their geometries into the geometry bag.
// Note that a non-recycling cursor is used, as the features' geometries are being stored
// for later use.
int[] blockObjectIDs = 
{
  11043, 11049, 11057
};
IFeatureCursor blocksCursor = blocksFeatureClass.GetFeatures(blockObjectIDs,
  false);
IFeature blockFeature = null;
object missingType = Type.Missing;
while ((blockFeature = blocksCursor.NextFeature()) != null)
{
  geometryCollection.AddGeometry(blockFeature.Shape, ref missingType, ref
    missingType);
}

// The cursor is no longer needed, so dispose of it.
Marshal.ReleaseComObject(blocksCursor);

[VB.NET] 
' Create a new geometry bag and give it the same spatial reference as the
' blocks feature class.
Dim geometryBag As IGeometryBag = New GeometryBagClass()
Dim geometryCollection As IGeometryCollection = CType(geometryBag, IGeometryCollection)
Dim geoDataset As IGeoDataset = CType(blocksFeatureClass, IGeoDataset)
Dim spatialReference As ISpatialReference = geoDataset.SpatialReference
geometryBag.SpatialReference = spatialReference

' Get a feature cursor for the three blocks and put their geometries into the geometry bag.
' Note that a non-recycling cursor is used, as the features' geometries are being stored
' for later use.
Dim blockObjectIDs As Integer() = {11043, 11049, 11057}
Dim blocksCursor As IFeatureCursor = blocksFeatureClass.GetFeatures(blockObjectIDs, False)
Dim blockFeature As IFeature = blocksCursor.NextFeature()
Dim missingType As Object = Type.Missing
While Not blockFeature Is Nothing
    geometryCollection.AddGeometry(blockFeature.Shape, missingType, missingType)
    blockFeature = blocksCursor.NextFeature()
End While

' The cursor is no longer needed, so dispose of it.
Marshal.ReleaseComObject(blocksCursor)
 
When using a geometry bag as the query geometry, a spatial index should be created on the geometry bag to allow rapid access to the contained geometries during the spatial query:
 

[C#] 
// Cast the geometry bag to the ISpatialIndex interface and call the Invalidate method
// to generate a new spatial index.
ISpatialIndex spatialIndex = (ISpatialIndex)geometryBag;
spatialIndex.AllowIndexing = true;
spatialIndex.Invalidate();

[VB.NET] 
' Cast the geometry bag to the ISpatialIndex interface and call the Invalidate method
' to generate a new spatial index.
Dim spatialIndex As ISpatialIndex = CType(geometryBag, ISpatialIndex)
spatialIndex.AllowIndexing = True
spatialIndex.Invalidate()
Now the geometry bag can be used as the query geometry in a new spatial filter. The code below shows how to use the spatial filter with a spatial relationship of “contains” (for parcels completely contained by the blocks) to find the number of parcels inside of the three blocks:
 

[C#] 
// Create the spatial filter. Note that the SubFields property specifies that only
// the Shape field is retrieved, since the features' attributes aren't being inspected.
ISpatialFilter spatialFilter = new SpatialFilterClass();
spatialFilter.Geometry = geometryBag;
spatialFilter.GeometryField = parcelsFeatureClass.ShapeFieldName;
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelContains;
spatialFilter.SubFields = "Shape";

// Use IFeatureClass.FeatureCount to get a parcel count.
int parcelCount = parcelsFeatureClass.FeatureCount(spatialFilter);
Console.WriteLine("Parcels in the three blocks: {0}", parcelCount);

[VB.NET] 
' Create the spatial filter. Note that the SubFields property specifies that only
' the Shape field is retrieved, since the features' attributes aren't being inspected.
Dim spatialFilter As ISpatialFilter = New SpatialFilterClass()
spatialFilter.Geometry = geometryBag
spatialFilter.GeometryField = parcelsFeatureClass.ShapeFieldName
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelContains
spatialFilter.SubFields = "Shape"

' Use IFeatureClass.FeatureCount to get a parcel count.
Dim parcelCount As Integer = parcelsFeatureClass.FeatureCount(spatialFilter)
Console.WriteLine("Parcels in the three blocks: {0}", parcelCount)
 
Buffering and querying
Finding features within a certain distance of other features is a common task, and can be accomplished using a buffer. A buffer is a polygon that encloses a point, line or polygon at a specified distance. After buffering a feature, the buffer can then be used as the query geometry of a spatial filter to find all of the features within the specified distance of the feature. This example shows how to find the cities that have populations of 500,000 people or more within 500 kilometers of Osaka, Japan.
 
The image below shows the city data, along with the 500-kilometer buffer to be searched within (of course, the example assumes no such buffer already exists):
 
 
 
The first step is retrieving Osaka’s geometry (for this example, it’s assumed the feature’s Object ID is known) and applying a buffer to it. The buffer method takes an inbound argument in the same units as the spatial reference of the feature class being buffered. For the sake of simplicity it is assumed that the cities feature class is using a metric spatial reference and the units are meters.
 

[C#] 
// Find the feature for Osaka and get its geometry.
IFeature osakaFeature = citiesFeatureClass.GetFeature(2263);
IGeometry osakaGeometry = osakaFeature.Shape;

// Use the ITopologicalOperator interface to create a buffer.
ITopologicalOperator topoOperator = (ITopologicalOperator)osakaGeometry;
IGeometry buffer = topoOperator.Buffer(500000);

[VB.NET] 
' Find the feature for Osaka and get its geometry.
Dim osakaFeature As IFeature = citiesFeatureClass.GetFeature(2263)
Dim osakaGeometry As IGeometry = osakaFeature.Shape

' Use the ITopologicalOperator interface to create a buffer.
Dim topoOperator As ITopologicalOperator = CType(osakaGeometry, ITopologicalOperator)
Dim buffer As IGeometry = topoOperator.Buffer(500000)
A spatial filter can now be created, using the buffer as the query geometry, with a “Contains” spatial relationship. A where clause can also be applied, to remove smaller cities from the search – the cities dataset contains a “POP_RANK” integer field, with a lower number indicating higher population. Cities with a population rank of 3 or less have at least 500,000 people. Iterating through the cursor’s results will return the names of the cities within the buffer.
 

[C#] 
// Create the spatial filter.
ISpatialFilter spatialFilter = new SpatialFilterClass();
spatialFilter.Geometry = buffer;
spatialFilter.GeometryField = citiesFeatureClass.ShapeFieldName;
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelContains;
spatialFilter.SubFields = "CITY_NAME, POP_CLASS";
spatialFilter.WhereClause = "POP_RANK < 4";

// Find the position of the CITY_NAME and POP_CLASS fields in the feature class.
int cityNamePosition = citiesFeatureClass.FindField("CITY_NAME");
int popClassPosition = citiesFeatureClass.FindField("POP_CLASS");

// Execute the query.
IFeatureCursor featureCursor = citiesFeatureClass.Search(spatialFilter, true);
IFeature feature = null;
while ((feature = featureCursor.NextFeature()) != null)
{
  String cityName = Convert.ToString(feature.get_Value(cityNamePosition));
  String popClass = Convert.ToString(feature.get_Value(popClassPosition));
  Console.WriteLine("{0} --- Population: {1}", cityName, popClass);
}

[VB.NET] 
' Create the spatial filter.
Dim spatialFilter As ISpatialFilter = New SpatialFilterClass()
spatialFilter.Geometry = buffer
spatialFilter.GeometryField = citiesFeatureClass.ShapeFieldName
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelContains
spatialFilter.SubFields = "CITY_NAME, POP_CLASS"
spatialFilter.WhereClause = "POP_RANK < 4"

' Find the position of the CITY_NAME and POP_CLASS fields in the feature class.
Dim cityNamePosition As Integer = citiesFeatureClass.FindField("CITY_NAME")
Dim popClassPosition As Integer = citiesFeatureClass.FindField("POP_CLASS")

' Execute the query.
Dim featureCursor As IFeatureCursor = citiesFeatureClass.Search(spatialFilter, True)
Dim feature As IFeature = featureCursor.NextFeature()
While Not feature Is Nothing
    Dim cityName As String = Convert.ToString(feature.Value(cityNamePosition))
    Dim popClass As String = Convert.ToString(feature.Value(popClassPosition))
    Console.WriteLine("{0} --- Population: {1}", cityName, popClass)
End While
 
Defining a spatial relationship with a Shape Comparison Language string
In most cases, values from the esriSpatialRelEnum enumeration – such as esriSpatialRelTouches and esriSpatialRelWithin – can be used to define the appropriate spatial relationship for a query. If a query requires a spatial relationship that isn’t defined by the enumeration, a special value from the enumeration (esriSpatialRelRelation) and a Shape Comparison Language string can be used to define any spatial relationship. A filter’s shape comparison string can be set with the ISpatialFilter.SpatialRelDescription property. Each character in the string represents a relationship between the query geometry and the geometry being tested, and can have a value of ‘T’ (true), ‘F’ (false), or ‘*’ (not tested). The following table lists the relationships represented by each character:
 
 
Requested Geometry
Interior
Boundary
Exterior
Query Geometry
Interior
1
2
3
Boundary
4
5
6
Exterior
9
8
7
 
Some examples of SpatialRelDescription strings include the following strings (more can be found in the Component Help for ISpatialFilter.SpatialRelDescription). The accompanying illustrations show the query geometry of each relationship with dashed-blue borders and the polygons that satisfy the relationship highlighted in red. Note that how spatial relationships are evaluated can vary depending on the geometry types of both the query geometry and the feature class being queried.
 
Illustration
String
Description
T********
 
The interiors of the geometries must intersect.
T***T****
The boundaries of the geometries must intersect and their interiors must intersect.
F***T****
The boundaries of the geometries must intersect and their interiors must not intersect.
FF*FF****
The geometries must be completely disjoint.
 
One scenario where building a string like those listed above can be useful is when a spatial query is needed to find identical geometries. The string used depends on the type of geometry being compared. When trying to find points with identical geometries, T******** can be used, because if two points share the same interior, they must be coincident. With polylines and polygons, TFFFTFFF* should be used to make sure that the interiors and boundaries of both geometries only intersect the interiors and boundaries of the other, and that neither intersect the other’s exterior. Note that the Exterior-Exterior relationship is never tested, as the exteriors of two geometries will always intersect.
 
This example will show how to find identical polygons in a temporal dataset containing several features for each US state, where the features of each have identical geometries but differ in date and population attributes. Although this example is contrived (the state name attribute could be used to achieve the same result), the process would be the same for datasets with coincident geometries and no common attributes to search by.
 
The following code example shows how to create a selection set containing features with the same geometry as a known feature (given its Object ID):
 

[C#] 
// Get the feature with the known Object ID (California).
IFeature caFeature = featureClass.GetFeature(4);
IGeometry caGeometry = caFeature.Shape;

// Create a spatial filter with a SCL spatial relationship string.
ISpatialFilter spatialFilter = new SpatialFilterClass();
spatialFilter.Geometry = caGeometry;
spatialFilter.GeometryField = featureClass.ShapeFieldName;
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelRelation;
spatialFilter.SpatialRelDescription = "TFFFTFFF*";
spatialFilter.SubFields = "Start_Date, Population";

// Find the positions of the Start_Date and Population fields in the class.
int startDatePosition = featureClass.FindField("Start_Date");
int populationPosition = featureClass.FindField("Population");

// Execute the query.
ISelectionSet selectionSet = featureClass.Select(spatialFilter,
  esriSelectionType.esriSelectionTypeSnapshot,
  esriSelectionOption.esriSelectionOptionNormal, null);

[VB.NET] 
' Get the feature with the known Object ID (California).
Dim caFeature As IFeature = featureClass.GetFeature(4)
Dim caGeometry As IGeometry = caFeature.Shape

' Create a spatial filter with a 9IM spatial relationship.
Dim spatialFilter As ISpatialFilter = New SpatialFilterClass()
spatialFilter.Geometry = caGeometry
spatialFilter.GeometryField = featureClass.ShapeFieldName
spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelRelation
spatialFilter.SpatialRelDescription = "TFFFTFFF*"
spatialFilter.SubFields = "Start_Date, Population"

' Find the positions of the Start_Date and Population fields in the class.
Dim startDatePosition As Integer = featureClass.FindField("Start_Date")
Dim populationPosition As Integer = featureClass.FindField("Population")

' Execute the query.
Dim selectionSet As ISelectionSet = featureClass.Select(spatialFilter, esriSelectionType.esriSelectionTypeSnapshot, _
                                    esriSelectionOption.esriSelectionOptionNormal, Nothing)
 
Evaluating a specific spatial relationship using IRelationalOperator
The IRelationalOperator interface can be used to verify the existence of specific spatial relationships between two high-level geometries. The interface is used by casting a geometry to the interface, and providing one of the methods with a comparison geometry. The interface’s methods all have boolean return types that indicate whether or not the specific relationship actually exists.
 
The following scenario shows how to determine whether a planned highway will come into contact with a freshwater marsh, given a feature class of roads being considered for construction and a feature class of vegetation communities.
 
 
 
From the image above, it is easy to tell that the highway (the red line) intersects the marsh (the selected polygon), but assuming the process requires automation, the code example below shows how to determine whether the two intersect, assuming the ObjectIDs of both features are known:
 

[C#] 
// Get the marsh and highway features from their respective classes.
IFeature marshFeature = vegFeatureClass.GetFeature(518);
IFeature highwayFeature = roadsFeatureClass.GetFeature(39);

// Get the geometries of the two features.
IGeometry marshGeometry = marshFeature.Shape;
IGeometry highwayGeometry = highwayFeature.Shape;

// Cast the highway's geometry to IRelationalOperator and determine whether or
// not it crosses the marsh's geometry.
IRelationalOperator relationalOperator = (IRelationalOperator)highwayGeometry;
Boolean crosses = relationalOperator.Crosses(marshGeometry);
Console.WriteLine("Highway crosses marsh: {0}", crosses);

[VB.NET] 
' Get the marsh and highway features from their respective classes.
Dim marshFeature As IFeature = vegFeatureClass.GetFeature(518)
Dim highwayFeature As IFeature = roadsFeatureClass.GetFeature(39)

' Get the geometries of the two features.
Dim marshGeometry As IGeometry = marshFeature.Shape
Dim highwayGeometry As IGeometry = highwayFeature.Shape

' Cast the highway's geometry to IRelationalOperator and determine whether or
' not it crosses the marsh's geometry.
Dim relationalOperator As IRelationalOperator = CType(highwayGeometry, IRelationalOperator)
Dim crosses As Boolean = relationalOperator.Crosses(marshGeometry)
Console.WriteLine("Highway crosses marsh: {0}", crosses)
The IRelationalOperator3D interface provides similar functionality for evaluating whether Z-aware geometries are coincident, but with a single method - Disjoint3D.
 

Using spatial caching to optimize spatial queries

Client-side caching of feature values within a certain extent is known as spatial caching. Spatial caching can significantly improve performance when multiple spatial queries are performed in a common extent, as it reduces the number of DBMS round trips required. An example of where the spatial cache functionality is used in ArcGIS Desktop is ArcMap’s Map Cache.
 
The ISpatialCacheManager interface (along with ISpatialCacheManager2 and ISpatialCacheManager3) provides methods and properties to fill the spatial cache given an extent, check the extent and status of the cache, and to empty it when no longer needed. The image below illustrates an example of a situation where spatial caching should be used; four spatial queries are executed within a common extent (indicated with the red-dashed line):
 
 
 
To use a spatial cache with these queries, the required steps should be taken:
 
The following code example shows how to do this:

[C#] 
// Open the feature classes used by the queries.
IFeatureClass blocksFeatureClass = featureWorkspace.OpenFeatureClass("Blocks");
IFeatureClass parcelsFeatureClass = featureWorkspace.OpenFeatureClass("Parcels")
  ;

// Fill the spatial cache.
ISpatialCacheManager spatialCacheManager = (ISpatialCacheManager)
  featureWorkspace;

// Check if the cache has been filled.
if (!spatialCacheManager.CacheIsFull)
{
  // If not full, fill the cache.
  spatialCacheManager.FillCache(cacheExtent);
}

// Execute spatial queries.

// Empty the cache.
spatialCacheManager.EmptyCache();

[VB.NET] 
' Open the feature classes used by the queries.
Dim blocksFeatureClass As IFeatureClass = featureWorkspace.OpenFeatureClass("Blocks")
Dim parcelsFeatureClass As IFeatureClass = featureWorkspace.OpenFeatureClass("Parcels")

' Fill the spatial cache.
Dim spatialCacheManager As ISpatialCacheManager = CType(featureWorkspace, ISpatialCacheManager)

' Check if the cache has been filled.
If spatialCacheManager.CacheIsFull Then
    ' If not full, fill the cache.
    spatialCacheManager.FillCache(cacheExtent)
End If

' Execute spatial queries.

' Empty the cache.
spatialCacheManager.EmptyCache()
 

Spatial relationships and feature class tolerance
It is important to consider the XY tolerance of a feature class when evaluating the spatial relationships. Different XY tolerances values can produce different results for relational and topological operations. For example, two geometries might be classified as disjoint (features physically not touching) with a small XY tolerance, but a larger XY tolerance value might cause them to be classified as intersecting. This is because the XY tolerance is taken into consideration when evaluating spatial relationships between objects in the Geodatabase.
 
 

Summary
The choice to use ISpatialFilter or IRelationalOperator for spatial queries depends on what is known prior to executing the query, as well as what type of results are desired. If the goal of the query is to find features that satisfy a spatial relationship with a single geometry (or a collection, if a geometry bag is used), a spatial filter should be used. When trying to verify that a spatial relationship exists (or doesn’t exist) between two specific geometries, a relational operator would be a better choice.
 
Spatial caching and tolerance should be considered when executing spatial queries. If multiple spatial queries are going to be executed within a common extent, the use of spatial caching can significantly increase performance by reducing round trips to the data source. XY Tolerance of a feature class can have an effect on the results of a spatial query, and should be kept in mind when executing spatial queries, particularly with feature classes that have unusually large XY tolerances.


See Also:

How to query geodatabase tables and feature classes