How to execute spatial queries


This document was published with and applies to ArcGIS 9.3.
A 10 version also exists. A 9.2 version also exists.
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 (finding all streets within an envelope) or searching for features based on their relation to other features (finding all farmland intersected by a particular stream). The topic, How to query geodatabase tables and feature classes, introduces spatial querying—this topic expands upon that information with illustrations, code examples, 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. 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

The following are two different APIs available to developers who need 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 previously known. Also, 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 before executing the spatial query.

Prerequisites for code examples in this topic

The following considerations apply to all code examples in this topic:

Examples of spatial queries

This section shows several spatial query illustrations and code exmples that use the ISpatialFilter and IRelationalOperator interfaces. For an introduction to using the ISpatialFilter interface and the IQueryFilter interface, which the ISpatialFilter interface extends, see How to query geodatabase tables and feature classes. The examples shown include the following scenarios:
 
Finding features within a polygon
This scenario shows how to execute a spatial query to find polyline features that intersect a polygon—specifically, major highways that pass through Iowa. See the following illustration:
 
 
Retrieve the geometry for Iowa, assuming that the ObjectID of the feature is known. When the feature is found, its geometry can be used as the query geometry in a spatial filter. Since the query locates all highways that pass through it, and not necessarily those contained entirely within it, set the spatial filter's relationship to "intersects."
 
The following code example 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 ObjectID.
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 shown.
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);
}

// Discard the cursors as they are no longer needed.
Marshal.ReleaseComObject(highwayCursor);

[VB.NET]
' Get the feature and its geometry given an ObjectID.
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 shown.
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

' Discard the cursors as they are no longer needed.
Marshal.ReleaseComObject(highwayCursor)
Querying for features with relation to multiple geometries
This scenario 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 scenario shows how to find parcel features within a known set of block features. The block features are not adjacent, but their ObjectIDs are known.
 
The following illustration shows the blocks and parcels for the search area—limit the spatial query to find the parcels within the blue highlighted blocks:
 
 
Create a new geometry bag and fill it with the geometries of the three block features (in this scenario, it is assumed the ObjectIDs 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.
// A non-recycling cursor is used, as the features' geometries are 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);
}

// Discard the cursors as they are no longer needed.
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.
' A non-recycling cursor is used, as the features' geometries are 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

' Discard the cursors as they are no longer needed.
Marshal.ReleaseComObject(blocksCursor)
When using a geometry bag as the query geometry, create a spatial index on the geometry bag to allow rapid access to the contained geometries during the spatial query. See the following code example:
 

[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 following code example 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 the three blocks:
 

[C#]
// Create the spatial filter. 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. 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 that can be accomplished by 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 be used as the query geometry of a spatial filter to find all features within the specified distance of the feature. This scenario shows how to find cities that have populations of 500,000 or more within 500 kilometers of Osaka, Japan.
 
The following illustration shows city data with a 500 kilometer buffer to search within (it's assumed no buffer exists):
 
 
Retrieve Osaka's geometry (it's assumed the feature's ObjectID is known) and apply 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. It is assumed that the cities' feature class is using a metric spatial reference and the units are meters. See the following code example:
 

[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)
Create a spatial feature by using the buffer as the query geometry with a "contains" spatial relationship. A WhereClause 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 less than 4 have at least 500,000 people. Iterating through the cursor's results return the names of the cities within the buffer. See the following code example:
 

[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 shows the relationships represented by each character:
 
 
Requested geometry
Interior
Boundary
Exterior
Query geometry
Interior
1
2
3
Boundary
4
5
6
Exterior
7
8
9
 
Some examples of SpatialRelDescription strings include the strings in the following table (more can be found in the previously mentioned ISpatialFilter.SpatialRelDescription). The accompanying illustrations show the query geometry of each relationship with dashed-blue borders and the polygons that satisfy the relationship are highlighted in red. How spatial relationships are evaluated can vary depending on the geometry types of 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 disjointed.
 
One scenario where building a string (such as those in the preceding table) 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 ensure the interiors and boundaries of both geometries only intersect the interiors and boundaries of the other, and that neither intersect the other's exterior. The Exterior-Exterior relationship is never tested, as the exteriors of two geometries always intersect.
 
This scenario shows how to find identical polygons in a temporal dataset containing several features for each state, where the features of each have identical geometries but differ in date and population attributes. Although this scenario is contrived (the state name attribute can be used to achieve the same result), the process is 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 ObjectID):
 

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

// Create a spatial filter with a shape comparison language 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 ObjectID (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
Use IRelationalOperator to verify the existence of specific spatial relationships between two high-level geometries. IRelationalOperator is used by casting a geometry to the interface, and providing one of the methods with a comparison geometry. The interface's methods have Boolean return types that indicate if the specific relationship exists.
 
This scenario shows how to determine if a planned highway comes into contact with a freshwater marsh, given a feature class of roads being considered for construction and a feature class of vegetation communities. See the following illustration:
 
 
On the preceding illustration, it is easy to determine the highway (red line) intersects the marsh (selected polygon), but assuming the process requires automation, the following code example shows how to determine if 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 if
// 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 if
' 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)
IRelationalOperator3D 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 database management system (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 (with ISpatialCacheManager2 and ISpatialCacheManager3) provides methods and properties to fill the spatial cache given an extent, checks the extent and status of the cache, and empties it when it is not needed.
 
The following illustration shows a situation in which spatial caching should be used; four spatial queries are executed within a common extent (indicated by the red-dashed line):
 
 
Do the following requirements to use a spatial cache with the queries:
 
The following code example shows how this is done:
 

[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 x,y tolerance of a feature class when evaluating spatial relationships. Different x,y tolerance values can produce different results for relational and topological operations. For example, two geometries can be classified as disjointed (features physically not touching) with a small x,y tolerance, but a larger x,y tolerance value may cause them to be classified as intersecting. This is because the x,y tolerance is taken into consideration when evaluating spatial relationships between objects in the geodatabase. See the following illustration:
 
 

Summary

The choice to use ISpatialFilter or IRelationalOperator for spatial queries depends on what is known prior to executing the query and the type of results needed. 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), use a spatial filter. When trying to verify that a spatial relationship exists (or doesn't exist) between two specific geometries, a relational operator is a better choice.
 
Consider spatial caching and tolerance 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. A feature class x,y tolerance can have an effect on the results of a spatial query and should be considered when executing spatial queries, particularly with feature classes that have unusually large x,y tolerances.


See Also:

How to query geodatabase tables and feature classes
Shape comparison language