MongoDB的地埋空间数据存储、空间索引以及空间查询

一、关于MongoDB

在众多NoSQL数据库，MongoDB是一个优秀的产品。其官方介绍如下： 
MongoDB (from "humongous") is a scalable, high-performance, open source, document-oriented database.

看起来，十分诱人！值得说明的是，MongoDB的document是以BSON（Binary JSON）格式存储的，完全支持Schema Free。这对地理空间数据是十分友好的。因为有著名的GeoJSON可供使用。另外OGR库也支持将Geometry类型导出为JSON格式。

本文将尝试使用OGR库把Shapefile导入到MongoDB存储，然后建立空间索引，进行空间查询。

著名的Foursquare使用了MongoDB数据库。

二、开发环境

MongoDB+Python+Pymongo+GDAL for Python

关于MongoDB和Python安装，本文不做介绍。
在继续本文之前，请先启动你的MongoDB服务器。本文默认采用如下服务器参数： 
Server:localhost 
Post:27017 
Database Name:gisdb

三、将shapefile导入到MongoDB

这里我直接提供代码，代码中已经有比较详尽的注释了。代码基本源于“引文1”，只是做了些改动，将MongoDB的Geometry的存储格式由wkt改成json。你可直接复制并运行下面的代码，当然需要修改一下Shapefile路径和MongoDB服务器相关参数。

import os 
import sys 
import json 
from pymongo import json_util 
from pymongo.connection import Connection 
from progressbar import ProgressBar 
from osgeo import ogr

def shp2mongodb(shape_path, mongodb_server, mongodb_port, mongodb_db, mongodb_collection, append, query_filter): 
        """Convert a shapefile to a mongodb collection""" 
        print ‘Converting a shapefile to a mongodb collection ‘ 
        driver = ogr.GetDriverByName(‘ESRI Shapefile’) 
        print ‘Opening the shapefile %s…’ % shape_path 
        ds = driver.Open(shape_path, 0) 
        if ds is None: 
                print ‘Can not open’, ds 
                sys.exit(1) 
        lyr = ds.GetLayer() 
        totfeats = lyr.GetFeatureCount() 
        lyr.SetAttributeFilter(query_filter) 
        print ‘Starting to load %s of %s features in shapefile %s to MongoDB…’ % (lyr.GetFeatureCount(), totfeats, lyr.GetName()) 
        print ‘Opening MongoDB connection to server %s:%i…’ % (mongodb_server, mongodb_port) 
        connection = Connection(mongodb_server, mongodb_port) 
        print ‘Getting database %s’ % mongodb_db 
        db = connection[mongodb_db] 
        print ‘Getting the collection %s’ % mongodb_collection 
        collection = db[mongodb_collection] 
        if append == False: 
                print ‘Removing features from the collection…’ 
                collection.remove({}) 
        print ‘Starting loading features…’ 
        # define the progressbar 
        pbar = ProgressBar(maxval=lyr.GetFeatureCount()).start() 
        k=0 
        # iterate the features and access its attributes (including geometry) to store them in MongoDb 
        feat = lyr.GetNextFeature() 
        while feat: 
                mongofeat = {} 
                geom = feat.GetGeometryRef() 
                mongogeom = geom.ExportToJson() 
                # store the geometry data with json format 
                mongofeat[‘geom‘] = json.loads(mongogeom,object_hook=json_util.object_hook)
                # iterate the feature’s  fields to get its values and store them in MongoDb 
                feat_defn = lyr.GetLayerDefn() 
                for i in range(feat_defn.GetFieldCount()): 
                        value = feat.GetField(i) 
                        if isinstance(value, str): 
                                value = unicode(value, "gb2312") 
                        field = feat.GetFieldDefnRef(i) 
                        fieldname = field.GetName() 
                        mongofeat[fieldname] = value 
                # insert the feature in the collection 
                collection.insert(mongofeat) 
                feat.Destroy() 
                feat = lyr.GetNextFeature() 
                k = k + 1 
                pbar.update(k) 
        pbar.finish() 
        print ‘%s features loaded in MongoDb from shapefile.’ % lyr.GetFeatureCount() 
        
        
input_shape = ‘/home/evan/data/map/res4_4m/XianCh_point.shp’ 
mongodb_server = ‘localhost’ 
mongodb_port = 27017 
mongodb_db = ‘gisdb’ 
mongodb_collection = ‘xqpoint’ 
filter = ”

print ‘Importing data to mongodb…’ 
shp2mongodb(input_shape, mongodb_server, mongodb_port, mongodb_db, mongodb_collection, False, filter)

四、MongoDB中空间数据的存储格式

在MongoDB的Shell中执行： 
>db.xqpoint.findOne() 
结果如下：

{ 
    "_id" : ObjectId("4dc82e7f7de36a5ceb000000"), 
    "PERIMETER" : 0, 
    "NAME" : "漠河县", 
    "PYNAME" : "Mohe Xian", 
    "AREA" : 0, 
    "ADCODE93" : 232723, 
    "CNTYPT_ID" : 31, 
    "CNTYPT_" : 1, 
    "geom" : { 
        "type" : "Point", 
        "coordinates" : [ 
            122.53233, 
            52.968872 
        ] 
    }, 
    "ID" : 1031, 
    "PN" : 1, 
    "CLASS" : "AI" 
}

这便是一个document，使用JSON格式，一目了然。其中的"geom"即为Geometry类型的数据，即地理空间数据，也是采用JSON格式存储，这样后续的空间索引与空间查询将十分方便。

MongoDB原生地支持了空间索引与空间查询，这一点比PostgreSQL方便，不再需要使用PostGIS进行空间扩展了。至于性能，我还没测试，在此不敢妄加评论。

五、在MongoDB中建立空间索引

>db.xqpoint.ensureIndex({‘geom.coordinates’:’2d’})

是不是十分简单？其它参数及用法请自行查看MongoDB手册。

六、在MongoDB中进行空间查询

>db.xqpoint.find({"geom.coordinates":[122.53233,52.968872]})

即可查询到上述“莫河县”这个点。当然，像这种精确查询，实际应用并不多。实际应用的空间查询大多为范围查询。MongoDB支持邻域查询（$near），和范围查询（$within）。

1. 邻域查询($near)

>db.xqpoint.find({"geom.coordinates":{$near:[122,52]}}) 
上述查询语句查询点[122,52]附近的点，MongoDB默认返回附近的100个点，并按距离排序。你也可以用limit()指定返回的结果数量，如：>db.xqpoint.find({"geom.coordinates":{$near:[122,52]}}).limit(5)

另外，你也可以指定一个最大距离，只查询这个距离内的点。 
>db.xqpoint.find({"geom.coordinates":{$near:[122,52],$maxDistance:5}}).limit(5)

MongoDB的find()方法可很方便的进行查询，同时MongoDB也提供了geoNear命令，用于邻域查询。 
>db.runCommand({geoNear:"xqpoint",near:[122,56],num:2}) 
上述语句用于查询[122,56]点附近的点，并只返回2个点。结果如下：

{ 
    "ns" : "gisdb.xqpoint", 
    "near" : "1110011000111101111100010000011000111101111100010000", 
    "results" : [ 
        { 
            "dis" : 3.077515616588727, 
            "obj" : { 
                "_id" : ObjectId("4dc82e7f7de36a5ceb000000"), 
                "PERIMETER" : 0, 
                "NAME" : "漠河县", 
                "PYNAME" : "Mohe Xian", 
                "AREA" : 0, 
                "ADCODE93" : 232723, 
                "CNTYPT_ID" : 31, 
                "CNTYPT_" : 1, 
                "geom" : { 
                    "type" : "Point", 
                    "coordinates" : [ 
                        122.53233, 
                        52.968872 
                    ] 
                }, 
                "ID" : 1031, 
                "PN" : 1, 
                "CLASS" : "AI" 
            } 
        }, 
        { 
            "dis" : 4.551319677334594, 
            "obj" : { 
                "_id" : ObjectId("4dc82e7f7de36a5ceb000001"), 
                "PERIMETER" : 0, 
                "NAME" : "塔河县", 
                "PYNAME" : "Tahe Xian", 
                "AREA" : 0, 
                "ADCODE93" : 232722, 
                "CNTYPT_ID" : 66, 
                "CNTYPT_" : 2, 
                "geom" : { 
                    "type" : "Point", 
                    "coordinates" : [ 
                        124.7058, 
                        52.340332 
                    ] 
                }, 
                "ID" : 1059, 
                "PN" : 1, 
                "CLASS" : "AI" 
            } 
        } 
    ], 
    "stats" : { 
        "time" : 0, 
        "btreelocs" : 85, 
        "nscanned" : 85, 
        "objectsLoaded" : 4, 
        "avgDistance" : 3.814417646961661, 
        "maxDistance" : 4.551319677334594 
    }, 
    "ok" : 1 
}

当然，我们也可附加条件查询条件，如查询[122,56]附近的且"PYNAME"为"Tahe Xian"的点： 
>db.runCommand({geoNear:"xqpoint",near:[122,56],num:2,query:{"PYNAME":"Tahe Xian"}) 
返回结果如下：

{ 
    "ns" : "gisdb.xqpoint", 
    "near" : "1110011000111101111100010000011000111101111100010000", 
    "results" : [ 
        { 
            "dis" : 4.551319677334594, 
            "obj" : { 
                "_id" : ObjectId("4dc82e7f7de36a5ceb000001"), 
                "PERIMETER" : 0, 
                "NAME" : "塔河县", 
                "PYNAME" : "Tahe Xian", 
                "AREA" : 0, 
                "ADCODE93" : 232722, 
                "CNTYPT_ID" : 66, 
                "CNTYPT_" : 2, 
                "geom" : { 
                    "type" : "Point", 
                    "coordinates" : [ 
                        124.7058, 
                        52.340332 
                    ] 
                }, 
                "ID" : 1059, 
                "PN" : 1, 
                "CLASS" : "AI" 
            } 
        } 
    ], 
    "stats" : { 
        "time" : 45, 
        "btreelocs" : 2095, 
        "nscanned" : 2096, 
        "objectsLoaded" : 2096, 
        "avgDistance" : 4.551319677334594, 
        "maxDistance" : 4.551319677334594 
    }, 
    "ok" : 1 
}

2. 范围查询($within)

MongoDB的$within操作符支持的形状有$box（矩形）,$center（圆形）,$polygon（多边形，包括凹多边形和凸多边形）。所有的范围查询，默认是包含边界的。

查询一个矩形范围，需要指定矩形的左下角和右上角两个坐标点，如下： 
> box = [[80,40],[100,50]] 
> db.xqpoint.find({"geom.coordinates":{$within:{$box:box}}})

查询一个圆形范围，需要指定圆心坐标和半径，如下： 
> center = [80,44] 
> radius =5 
> db.xqpoint.find({"geom.coordinates":{$within:{$center:[center,radius]}}})

查询一个多边形范围，需要指定多边形的各个顶点，可以通过一个顶点数组或一系列点对象指定。其中，最后一个点是默认与第一个点连接的。如下： 
> polygon1 = [[75,35],[80,35],[80,45],[60,40]] 
> db.xqpoint.find({"geom.coordinates":{$within:{$polygon:polygon1}}}) 
或者 
> polygon2 = ｛a:{75,35},b:{80,35},c:{80,45},d:{60,40}} 
> db.xqpoint.find({"geom.coordinates":{$within:{$polygon:polygon2}}})

注意：MongoDB 1.9及以上版本才支持多边形范围查询。

P.S. MongoDB还支持复合索引，球面模型（可简单理解为投影吧），多位置文档（Multi-location Documents，即一个文档中包括多个Geometry），可参见“引文2”或MongoDB手册。

七、参考资料

引文1：http://www.paolocorti.net/2009/12/06/using-mongodb-to-store-geographic-data/ 
引文2：http://www.mongodb.org/display/DOCS/Geospatial+Indexing

时间： 2024-10-09 09:22:38

MongoDB的地埋空间数据存储、空间索引以及空间查询的相关文章

Mysql空间数据，空间索引，Spatial Data，Spatial Index

前文: 这两天因为项目原因看了一下MySQL的空间索引,发现网上的资料不多,查了一下官方文档,为了强化记忆做了一个简单的翻译.基本上理解了mysql空间索引的要点.谨以此纪. Extensions for Spatial Data Open Geospatial Consortium (OGC) 是一个由超过两百五十个公司,机构,大学组成的致力于发展管理空间数据的解决方案的组织. OGC 发布了OpenGIS® Implementation Standard for Geographic inf

MongoDB系列五（地理空间索引与查询）.

一.经纬度表示方式 MongoDB 中对经纬度的存储有着自己的一套规范(主要是为了可以在该字段上建立地理空间索引).包括两种方式,分别是 Legacy Coordinate Pairs (这个词实在不知道怎么翻译...) 和  GeoJSON . Legacy Coordinate Pairs Legacy Coordinate Pairs 又有两种方式可以存储经纬度,可以使用数组(首选)或嵌入式文档. 数组: <field>: [<longitude>, <latitude

浅析MongoDB数据库的海量数据存储应用

[摘要]当今已进入大数据时代,特别是大规模互联网web2.0应用不断发展及云计算所需要的海量存储和海量计算发展,传统的关系型数据库已无法满足这方面的需求.随着NoSQL数据库的不断发展和成熟,可以较好地解决海量存储和海量计算方面的应用需求.本文重点描述作为NoSQL之一MongoDB数据库在海量数据存储方面的应用. 1 引言NoSQL,全称是“Not Only Sql”,指的是非关系型的数据库.这类数据库主要有这些特点:非关系型的.分布式.开源的.水平可扩展的.原始目的是为了大规模web应用,这

SQL Server 2008空间数据应用系列八：基于Bing Maps(Silverlight)的空间数据存储

原文:SQL Server 2008空间数据应用系列八:基于Bing Maps(Silverlight)的空间数据存储 友情提示,您阅读本篇博文的先决条件如下: 1.本文示例基于Microsoft SQL Server 2008 R2调测. 2.具备 Transact-SQL 编程经验和使用 SQL Server Management Studio 的经验. 3.具有使用 Microsoft Visual Studio 进行 Microsoft .NET Framework开发的经验. 4.具有

Oracle获取干净的建表DDL语句，不含其它存储、表空间、段属性

早上一个同事资讯怎么获取到建表语句而且是不带存储那种SQL.Oracle自己提供了一个函数DBMS_METADATA.GET_DDL,但是获取到的建表语句含有存储.表空间.以及一些其他段的属性.如图: 看到这个获取到的ddl语句,想通过利用Oracle函数来截取的方式获取建表语句. 思路为:1.通过get_ddl获取建表语句 abc2.将abc中的pctfree'替换成';'3.计算';'的位置4.用substr来截取abc,从开头到';'的长度 SQL如下: SELECT SUBSTR(REP

当存储无可用空间时无法启动虚拟机

今天一个网友问我一下问题: 老师这个是什么问题呀开机时提示这个 答:你的存储是不是没空间了 网友:我看一下 是存储空间不够了,如何删掉一些不要的东西呀 . 这个问题的原因很简单,是虚拟机所在存储空间不够了,导致虚拟机不能启动. 虚拟机在启动的时候,会创建一个交换文件,此交换文件与虚拟机分配 的内存同样大小,当虚拟机所在存储空间不足时就会出现这个报警. 解决方法: (1)当前这个虚拟机分配了6144MB的内存.如果存储空间不足6144MB,修改的虚拟机内存,小于存储剩余空间,可以再次启动.但这只是

如何使用sqlserver 2012 空间查询（geometry及 geography）

---恢复内容开始--- 介绍 SQL Server 2008为大地测量空间数据提供了geography数据类型为平面空间数据提供了geometry数据类型.这两个都是Microsoft .NET Framework通用语言运行时(CLR)类型,并且可以用来存储不同种类的地理元素,例如点.线和多边形.这两个数据类型都提供了你可以用来执行空间操作的属性和方法,例如计算位置间的距离和找出两者间交叉的地理特性(例如一条河流经一个城镇.) geography 数据类型 geography数据类型为空间数

MYSQL的空间查询(转帖)

SELECT x(location),y(location) FROM frddata.points; 本文将向各位介绍如何使用MySql5.x中的空间数据库,并展示一下它高效的性能(前提是正确使用). 本文适合于对SQL和MYSQL熟悉的人员. 步骤1:创建支持空间查询的表 首先来说一下如何创建一个包含空间数据的名为Points的表. CREATE TABLE `points` ( `name` varchar(20) NOT NULL DEFAULT '', `location` point

Oracle spatial空间查询的选择度分析

在上一篇中,我用一个案例演示了对于数值或字符串类型的字段,选择度的计算方法.并证明了当字段值的选择度不同时,将会影响CBO选择最终的执行计划.对于可排序的字段类型,选择度计算模型已经有很多人写博客介绍了,但空间查询的选择度怎么计算却少有人研究?Oracle会根据不同的检索范围,产生不同的选择度吗? 接下来,我们来研究一下这个问题. 创建表,并使用SDO_GEOMETRY数据类型存储矢量数据. 查看表中记录数: 创建空间索引: CREATE INDEX "TDDCSDE"."A