1.环境

-
VS2015

-
C#

-
DotSpatial

2.范围索引

2.1.简介

范围索引即在记录每个空间实体的坐标时，同时记录每个空间实体的最大和最小坐标（即是指该空间实体的最小外包矩形）。在通过一个查询范围查询包含在其中的空间实体时，根据空间实体的MBR（最小外包矩形），先筛选出落入查询范围内的空间实体，再对筛选的结果进行进一步的精细查询，最终得到真正与查询范围相交的空间实体。

图1 范围索引示例图

2.2.建立索引

使用C#的字典对象，其key为索引，value为空间实体。

1	`//接口ISpatialIndex using System.Collections.Generic; namespace SpatialIndex { interface ISpatialIndex { Dictionary<T1, T2> GetSpatialIndexDictionary<T1,T2>(string shpPath); } }`

`//MBR索引的实现 using System.Collections.Generic; using System.Linq; using DotSpatial.Data; using DotSpatial.Topology; namespace SpatialIndex { class MbrIndex:ISpatialIndex { public Dictionary<T1, T2> GetSpatialIndexDictionary<T1, T2>(string shpPath) { // 使用FeatureSet打开文件 IFeatureSet featureSet = FeatureSet.Open(shpPath); // 循环遍历各个Feature,获取其MBR，并建立MBR与要素间的对应关系 Dictionary<IEnvelope, IFeature> extentIndexDict = featureSet.Features.ToDictionary(feature => feature.Envelope); object obj = (object)extentIndexDict; Dictionary<T1, T2> dict = (Dictionary<T1, T2>) obj; return dict; } } }`

3.格网索引

3.1.简介

格网索引的基本思想是将区域划分成相同大小的网格，记录每个网格内所包含的空间实体在数据库中的地址。格网索引编码可以按照行列顺序进行编码，也可使用莫顿码等其他方式来进行编码。

图2 格网索引示例图

3.2.建立索引

同样使用C#的字典对象，其key为索引，但value为空间实体集合。

`//格网索引实体类 namespace SpatialIndex { class GridIndexInfo { /// <summary> /// 格网索引的行列数，须为2的幂，默认值为16 /// </summary> public int RowColNum { get; set; } = 16; /// <summary> /// 格网索引中行列间距， /// [0]列间距 /// [1]行间距 /// </summary> public double[] RowColInterval { get; set; } /// <summary> /// 格网索引的范围 /// [0]最上方，maxY /// [1]最右侧，maxX /// [2]最下方，minY /// [3]最左侧，minX /// </summary> public double[] GridExtent { get; set; } } }`

`//莫顿码生成 using System; using System.Collections.Generic; namespace SpatialIndex { class Morton { /// <summary> /// 获取指定行列数的莫顿码 /// </summary> /// <param name="rowIndex">行数</param> /// <param name="colIndex">列数</param> /// <returns></returns> public static string GetMortonCode(int rowIndex, int colIndex) { string morton = ""; List<char> mc = new List<char>(); string rowBit = Convert.ToString(rowIndex, 2); string colBit = Convert.ToString(colIndex, 2); while (rowBit.Length > colBit.Length) { colBit = 0 + colBit; } while (colBit.Length > rowBit.Length) { rowBit = 0 + rowBit; } for (int i = 0; i < rowBit.Length; i++) { mc.Add(rowBit[i]); mc.Add(colBit[i]); } for (int i = 0; i < mc.Count; i++) { morton += mc[i].ToString(); } int mottonInt = Convert.ToInt32(morton, 2); return mottonInt.ToString(); } /// <summary> /// 位于一定行列数内的所有莫顿码 /// </summary> /// <param name="minColIndex">最小列数</param> /// <param name="maxColIndex">最大列数</param> /// <param name="minRowIndex">最小行数</param> /// <param name="maxRowIndex">最大行数</param> /// <returns></returns> public static List<string> GetMortonList(int minColIndex, int maxColIndex, int minRowIndex, int maxRowIndex) { List<string> mortonList = new List<string>(); for (int i = minRowIndex; i <= maxRowIndex; i++) { for (int j = minColIndex; j <= maxColIndex; j++) { mortonList.Add(GetMortonCode(i, j)); } } return mortonList; } } }`

`//格网索引的实现 using System.Collections.Generic; using System.Linq; using DotSpatial.Data; using DotSpatial.Topology; namespace SpatialIndex { class GridIndex: GridIndexInfo,ISpatialIndex { public Dictionary<T1, T2> GetSpatialIndexDictionary<T1, T2>(string shpPath) { IFeatureSet featureSet = FeatureSet.Open(shpPath); //格网索引 Dictionary<string, List<IFeature>> gridIndexDict = new Dictionary<string, List<IFeature>>(); //shapefile要素类的范围 double top = featureSet.Extent.MaxY; double right = featureSet.Extent.MaxX; double bottom = featureSet.Extent.MinY; double left = featureSet.Extent.MinX; double xInterval = (right - left) / RowColNum; double yInterval = (top - bottom) / RowColNum; double bottomCache = bottom; for (int i = 0; i < RowColNum; i++) { double leftcache = left; for (int j = 0; j < RowColNum; j++) { //获取格网对应的莫顿码 string morton = Morton.GetMortonCode(i, j); //创建格网 Envelope enve = new Envelope(leftcache, leftcache + xInterval, bottomCache, bottomCache + yInterval); //如果多边形和格网相交，将其添加到索引集合中 List<IFeature> list = featureSet.Features.Where(feature => feature.Intersects(enve)).ToList(); gridIndexDict.Add(morton, list); leftcache = leftcache + xInterval; } bottomCache = bottomCache + yInterval; } base.RowColInterval = new double[] {xInterval, yInterval}; base.GridExtent = new double[] {top, right, bottom, left}; object obj = (object) gridIndexDict; Dictionary<T1, T2> dict = (Dictionary<T1, T2>) obj; return dict; } } }`

4.查询的实现

`//查询接口 using System.Collections.Generic; using DotSpatial.Data; namespace SpatialIndex { interface IQuery { List<IFeature> GetPointQueryResult(double pointX,double pointY, double r); List<IFeature> GetRectQueryResult(double minX,double minY,double maxX,double maxY); } }`

`//范围索引查询 using System.Collections.Generic; using DotSpatial.Data; using DotSpatial.Topology; namespace SpatialIndex { class MbrIndexQuery : IQuery { public Dictionary<IEnvelope, IFeature> MbrIndexDict { get; set; } public List<IFeature> GetPointQueryResult(double pointX, double pointY, double r) { Envelope queryRect = new Envelope(pointX - r, pointX + r, pointY - r, pointY + r); // 待进一步验证的要素表 List<IFeature> features = new List<IFeature>(); foreach (KeyValuePair<IEnvelope, IFeature> mbrPair in MbrIndexDict) { // 比较查询矩形与MBR的相交关系 if (queryRect.Intersects(mbrPair.Key)) { features.Add(mbrPair.Value); } } Point queryPoint = new Point(pointX, pointY); List<IFeature> resultfeatures = new List<IFeature>(); foreach (IFeature feature in features) { // 比较查询点与空间实体的精确相交关系 if (feature.Intersects(queryPoint)) { resultfeatures.Add(feature); } } return resultfeatures; } public List<IFeature> GetRectQueryResult(double minX, double minY, double maxX, double maxY) { Envelope queryRect = new Envelope(minX, maxX, minY, maxY); List<IFeature> features = new List<IFeature>(); foreach (KeyValuePair<IEnvelope, IFeature> mbrPair in MbrIndexDict) { // 比较查询矩形与MBR的相交关系 if (queryRect.Intersects(mbrPair.Key)) { features.Add(mbrPair.Value); } } List<IFeature> resultfeatures = new List<IFeature>(); foreach (IFeature feature in features) { // 比较查询矩形与空间实体的精确相交关系 if (feature.Intersects(queryRect)) { resultfeatures.Add(feature); } } return resultfeatures; } } }`

`//格网索引查询 using System.Collections.Generic; using System.Linq; using DotSpatial.Data; using DotSpatial.Topology; namespace SpatialIndex { class GridIndexQuery : GridIndexInfo, IQuery { public Dictionary<string, List<IFeature>> GridIndexDict { get; set; } public List<IFeature> GetPointQueryResult(double pointX, double pointY, double r) { //格网索引可以根据点坐标直接计算出其所在格网的莫顿码，因此实际上可以不需要使用模糊距离r Point queryPoint = new Point(pointX, pointY); int row = (int)((pointY - base.GridExtent[2]) / base.RowColInterval[1]); int col = (int)((pointX - base.GridExtent[3]) / base.RowColInterval[0]); string gridPointMotton = Morton.GetMortonCode(row, col); List<IFeature> roughlyFeatures = GridIndexDict[gridPointMotton];//粗略结果集 List<IFeature> reFeatures = roughlyFeatures.Where(feature => feature.Contains(queryPoint)).ToList(); //有重结果集 return reFeatures; } public List<IFeature> GetRectQueryResult(double minX, double minY, double maxX, double maxY) { Envelope queryRect = new Envelope(minX, maxX, minY, maxY); int mincol = (int)((queryRect.Left() - base.GridExtent[3]) / base.RowColInterval[0]); int minrow = (int)((queryRect.Bottom() - base.GridExtent[2]) / base.RowColInterval[1]); int maxcol = (int)((queryRect.Right() - base.GridExtent[3]) / base.RowColInterval[0]); int maxrow = (int)((queryRect.Top() - base.GridExtent[2]) / base.RowColInterval[1]); List<string> gridindex = Morton.GetMortonList(mincol, maxcol, minrow, maxrow); List<IFeature> roughlyFeatures = new List<IFeature>();//粗略结果集 foreach (string key in gridindex) { roughlyFeatures.AddRange(GridIndexDict[key]); } List<IFeature> reFeatures = new List<IFeature>(); foreach (IFeature feature in roughlyFeatures) { if (queryRect.Intersects(feature.Envelope)) { reFeatures.Add(feature); } } //去重，因为01号格网包含空间实体A，02号也可能包含A HashSet<IFeature> hsResultIfeature = new HashSet<IFeature>(reFeatures); return hsResultIfeature.ToList(); } } }`