MongoDB地理空间数据存储及检索

以前写过MySQL空间索引简单使用,测试也是可用的,当时没有测试效率问题,由于存储的矢量数据都只是一个四点的多边形而已。此次使用mongoDB来作一个行政区划检索的功能,记录一下使用的过程。ios

参考资料:mongodb

一、存入地理数据

MongoDB存储的数据是bson结构,因此只要你的数据符合这个结构都是能够存储的,可是要支持空间索引,就必须按照它的规定来。
早期版本的(2.6以前)仅仅支持简单的点数据的索引,也就是filed:[x,y]这样的结构,这个适用范围太有限了。如今的版本支持GeoJSON形式的数据类型,且支持OCG的空间数据查询模型,使用上很是方便。shell

GeoJSON数据存入

MongoDB要求把GeoJSON格式的数据以子文档的形式存入,但实际上并非存入一个完整GeoJSON对象,只须要其中的typecoordinates两个字段就能够了。数据库

下面以存入一个含有地理空间数据的文档为例,把全部支持的GeoJSON对象类型都作个示例。
这里假设存储一个县的信息,数据都以json格式表示。json

{
      "xian":"潜山县",
      "sheng":"安徽"
}

一、Ponit 点数据

如今假设为这个文档添加上中心点位置,那么这个文档就变成了以下样子:api

{
      "xian":"潜山县",
      "sheng":"安徽"
      "center":{ "type":"Point","coordinates":[116.45,30.72]}
}

二、LineString 线数据(多段线)

如今加上一个到省会合肥的路径连线,那么文档就变成了以下样子:数组

{
      "xian":"潜山县",
      "sheng":"安徽"
      "center":{ "type":"Point","coordinates":[116.45,30.72]},
      "toShengHui":{
          "type":"LineString","
          coordinates":[[116.55944824218749,30.58827267102698],
          [116.87667846679689,30.791396195188927],[116.96594238281249,31.038815104128687],
          [117.18292236328124,31.264465555752835],[117.22412109375,31.819230730326613]]}
}

三、 Polygon 多边形数据

多边形是当前地理信息领域应用的比较多的数据类型。
多边形描述的是一个面对象,其由两部分组成,一个外壳shell和0或多个內洞holes
外壳和內洞的表示都是一个闭合的线环(LinearRing),表示一个闭合曲面。外壳包括的区域表示在多边形内的区域,內洞表示的区域则是从外壳表示的区域中排除的区域,以下图所示的样子,蓝色的是外壳,绿色部分的是內洞。
bash

由于行政边界的涉及到的点太多,因此这里就只添加一个外包框做为示例:
由于GeoJSON中使用bbox字段(四个double值的数组)来描述范围外包框,因此这里不能使用bbox来存储一个多边形,不然将报错error inserting documents: location object expected, location array not in correct format。但奇怪的是,这个bbox用于搜索的时候倒是无效的。
这里就不使用內洞了,有內洞的状况就是coordinates数组中后面加上线环便可。session

{
      "xian":"潜山县",
      "sheng":"安徽"
      "center":{ "type":"Point","coordinates":[116.45,30.72]},
      "toShengHui":{
          "type":"LineString","
          coordinates":[[116.55944824218749,30.58827267102698],
          [116.87667846679689,30.791396195188927],[116.96594238281249,31.038815104128687],
          [117.18292236328124,31.264465555752835],[117.22412109375,31.819230730326613]]},
      "box":{
         "type":"Polygon",
         "coordinates":[[[116.40701293945311,30.454001045389525],
         [116.77505493164062,30.454001045389525],[116.77505493164062,30.76248901825541],
         [116.40701293945311,30.76248901825541],[116.40701293945311,30.454001045389525]]]}
}

四、MultiPoint多点、MultiLineString多线、MultiPolygon多多边形

对于多点、多线和多多边形,与单个的区别,也就是将coordinates成员改成一个数组形式,存入多个单个形式的坐标数据。

五、GeometryCollection 几何集合

几何集合就是多个几何对象的集合,就是一个数组里面放多个几何对象。

六、全国区县行政区划入库示例

一、首先下载全国的性质边界矢量数据,这个能够从https://www.gadm.org/download_country_v3.html下载。由于中国的矢量数据中没有台湾和香港澳门的数据,能够下在后合成一个。这份数据还有一些其余的小问题,这里就不提了。这也是我能找到的免费数据中较好的一份。

二、下载的数据可使用GDAL或QGIS工具将其转换为geojson格式文档,也能够转换,直接写程序来读取。我把转换后的程序再通过了一次简化,由于全部的边界线都是MutilPolygon,而大多数边界是仅仅一个Polygon的,因此我把能转换的都转换成了Polygon

三、由于MongoDB中存储字段bbox为GeoJSON中的数组形式在查询的时候会有问题,因此我把它改成了多边形。但bbox字段又不支持多边形,因此改成以box字段来存储。由于县级行政区划的边界都比较复杂,点比较多,在查询的时候会比较慢,因此使用$and来先查询box在查询geometry会比较快。

处理好的数据能够在这里下载 连接: https://pan.baidu.com/s/1f2c9FEQhkfDzC1dZHn6XmA 密码:5a2u

四、处理完成以后把全部的json对象按行写入了到了(县级边界.json.txt)文件中,由于每行都是一个json对象,可使用mongoimport将其导入数据库中。可是直接导入会有一个问题,就是数据量太大,没法一次写入。因此我先拆分红了几个小文件。

split -l 512 县级边界.json.txt -d -a 1 中国县级行政边界_

执行上面命令后生成了5个小文件,而后逐个导入到mongodb便可。

mongoimport --host 192.168.0.28 --db us --collection xzbj --file 中国县级行政边界_0
2018-10-12T14:14:45.165+0800    connected to: 192.168.0.28
2018-10-12T14:14:47.386+0800    imported 512 documents
mongoimport --host 192.168.0.28 --db us --collection xzbj --file 中国县级行政边界_1
mongoimport --host 192.168.0.28 --db us --collection xzbj --file 中国县级行政边界_2
mongoimport --host 192.168.0.28 --db us --collection xzbj --file 中国县级行政边界_3
mongoimport --host 192.168.0.28 --db us --collection xzbj --file 中国县级行政边界_4

bbox只能以数组的形式存在,如图所示,我在处理的时候已经改用box多边形表示了。
bbox存在的形式

二、建立地理索引

MongoDB的空间索引有三种,2dsphere2dgeoHaystack

对于某些地理空间查询操做,必须有相应的索引才行

2.一、2dsphere索引

官网文档:https://docs.mongodb.com/manual/core/2dsphere/

2dsphere索引支持查询球面几何实体对象。2dsphere是MongoDB地理空间索引支持全部查询:用于查询、交点和接近。地理空间查询的更多信息,参见地理空间查询。
2dsphere索引支持点数据(包括传统点数据方式和GeoJSON的Point方式)。2dsphere索引因为是球面索引,因此仅仅支持经纬度数据,坐标系为WGS84
版本2之后的MongoDB支持附加GeoJSON对象包括MultiPoint、MultiLineString、MultiPolygon和 GeometryCollection,具体的参考官方文档。

建立一个2dsphere索引的语句以下:

db.collection.createIndex( { <location field> : "2dsphere" } )

这里有一个问题,就是建立的时候,有MultiPolygon等不支持的几何类型的时候会出现错误"errmsg" : "Can't extract geo keys...,但看官网文档,这个是能够附加的类型,但目前没有找到相关的文档。

2.二、2d索引

官网文档:https://docs.mongodb.com/manual/core/2d/

2d索引对存储为二维平面上的点的数据使用。 2d索引适用于MongoDB 2.2及更早版本中使用的旧坐标对。

应对仅在下列状况下使用2d索引
    您的数据库具备MongoDB 2.2或更早版本的遗留遗留坐标对,以及
    您不打算将任何位置数据存储为GeoJSON对象。

建立一个2d索引的语句以下:

db.collection.createIndex( { <location field> : "2d" } )

要在具备非简单排序规则的集合上建立2d索引,必须在建立索引时显式指定{collation:{locale:“simple”}}

2.三、geoHaystacks索引

官网文档:https://docs.mongodb.com/manual/core/geohaystack/

geoHaystack索引是一种特殊索引,优化小面积内的返回结果。geoHaystack索引可提升在平面进行几何(geometry)查询的性能。
对于使用球面的几何查询,2dsphere索引是一个比geoHaystack索引更好的选择。2dsphere索引容许字段从新排序。
geoHaystack索引要求第一个字段为location字段。此外,geoHaystack索引仅可经过命令使用,所以始终一次返回全部结果。

三、检索地理数据

检索这个也是看官方的文档比较快,这里只是一个简单的介绍。
官方文档:Geospatial Query Operators

3.1地理空间模型

MongoDB的地理空间查询能够实现球面或平面几何实体对象的查询。

  • 2dsphere索引仅仅支持球面查询(即把点坐标数据当作球面经纬度处理)。
  • 2d索引支持平面查询(即将点坐标数据当作平面直角坐标系点坐标处理)和一些球面查询,虽然2d索引支持一些球面查询,可是对这些球面查询使用2d索引可能会致使错误,这样的数据尽可能优先使用2dsphere索引。

下面列出每一个地理空间操做使用的地理空间查询运算符,支持的查询和相关说明:

操做符 参数类型 索引 支持查询
$near
邻近查询
GeoJSON质心点在这个line和下一个line, 2dsphere 球面
$near legacy coordinates 2d 平面
$nearSphere GeoJSON点 2dsphere 球面
$nearSphere legacy coordinates 2d 球面
$geoWithin
内部查询
GeoJSON几何对象
{ $geometry: … }
球面
$geoWithin { $box: … } 2d 平面
$geoWithin { $polygon: … } 2d 平面
$geoWithin { $center: … } 2d 平面
$geoWithin { $centerSphere: … } 2d/2dsphere 球面
$geoIntersects
相交查询
{ $geometry: … }
多边形或多多边形
球面

还有$geoNear这个操做符,这里就不摘录了,直接去官网看好了。
由于这些操做符都比较简单,这里只单独介绍一下最有用的$geoIntersects操做符。

$geoIntersects操做符使用
选择地理空间数据与指定GeoJSON对象相交的文档; 即数据和指定对象的交集是非空的。
$geoIntersects运算符使用$geometry运算符来指定一个GeoJSON对象做为参数,使用默认坐标系(CRS)指定GeoJSON多边形或多边形的使用语法以下:

{
     空间数据字段名: {
     $geoIntersects: {
        $geometry: {
           type: "<GeoJSON对象类型>" ,
           coordinates: [ <coordinates> ]
        }
     }
  }
}

对于$geoIntersects查询,当指定的GeoJSON的几何对象大于半个球面时,使用默认的坐标系(CRS)会致使互补的几何对象在查询结果中。
3.0版中的新功能:要指定单环的GeoJSON多边形使用自定义MongoDB CRS,使用如下语法在$geometry表达式中指定自定义MongoDB CRS:

{
  <location field>: {
     $geoIntersects: {
        $geometry: {
           type: "Polygon" ,
           coordinates: [ <coordinates> ],
           crs: {
              type: "name",
              properties: { name: "urn:x-mongodb:crs:strictwinding:EPSG:4326" }
           }
        }
     }
  }
}

自定义MongoDB CRS使用逆时针顺序包覆,并容许$geoIntersects支持具备单环GeoJSON多边形的查询,该多边形的面积大于或等于单个半球。若是指定的多边形小于单个半球,则带有MongoDB CRS的$ geoIntersects的行为与默认的CRS相同。“Big” Polygons参考

若是指定纬度和经度坐标,请先列出经度而后列出纬度:
    有效经度值介于-180和180之间(包括二者)。
    有效纬度值介于-90和90之间(包括二者)。

3.二、查询示例(使用全国县级行政边界数据)

上面的数据导入以后,写几个查询的例子来测试一下。

3.2.一、使用$geoIntersects查询相交的区域

$geoIntersects用于查询与给定参数有相交区域的记录,只能查询geojson形式表示的位置字段,有没有索引均可以查,速度与几何对象复杂程度有关。

查询代码以下

db.getCollection('xzbj').find({
    "geometry": {
        "$geoIntersects": {
            "$geometry": {
                "type": "Polygon",
                "coordinates": [[[116.24633789062499, 40.168380093142446], [116.17492675781251, 40.15998434802335], [116.1199951171875, 40.057052221322], [116.09527587890624, 40.002371935876475], [116.1474609375, 39.890772566959534], [116.10626220703124, 39.70929962338767], [116.3177490234375, 39.69662085337441], [116.57592773437499, 39.7642140375156], [116.6912841796875, 39.86969567045658], [116.69677734375, 39.99605985169435], [116.62261962890624, 40.094882122321145], [116.6143798828125, 40.13899044275822], [116.43310546875, 40.15788524950653], [116.28753662109375, 40.18097176388719], [116.24633789062499, 40.168380093142446]]]
            }
        }
    }
})

这些点是沿着北京六环线画的一个多边形,查询的速度比较慢,耗时达到7.12秒。

下面添加外包框过滤,在进行相交比较,加快查询速度。这里要注意我是要的box字段内存的也是一个Polygon而不是GeoJSON中的bbox。

db.getCollection('xzbj').find({
    "$and": [{
        "box": {
            "$geoIntersects": {
                "$geometry": {
                    "type": "Polygon",
                    "coordinates": [[[116.0870361328125, 39.69873414348139], [116.71600341796874, 39.69873414348139], [116.71600341796874, 40.17257757632168], [116.0870361328125, 40.17257757632168], [116.0870361328125, 39.69873414348139]]]
                }
            }
        }
    },
    {
        "geometry": {
            "$geoIntersects": {
                "$geometry": {
                    "type": "Polygon",
                    "coordinates": [[[116.24633789062499, 40.168380093142446], [116.17492675781251, 40.15998434802335], [116.1199951171875, 40.057052221322], [116.09527587890624, 40.002371935876475], [116.1474609375, 39.890772566959534], [116.10626220703124, 39.70929962338767], [116.3177490234375, 39.69662085337441], [116.57592773437499, 39.7642140375156], [116.6912841796875, 39.86969567045658], [116.69677734375, 39.99605985169435], [116.62261962890624, 40.094882122321145], [116.6143798828125, 40.13899044275822], [116.43310546875, 40.15788524950653], [116.28753662109375, 40.18097176388719], [116.24633789062499, 40.168380093142446]]]
                }
            }
        }
    }]
})

这样查询的速度就大大提高了,仅耗时41毫秒就完成了检索。通过测试,$and数组中的元素顺序对检索速度没有影响,不知道是否是与字段名的排序有关系。

3.2.3 使用$geoWithin查询

使用上和$geoIntersects差很少,查询速度上也差很少。

db.getCollection('xzbj').find({
    "geometry": {
        "$geoWithin": {
            "$geometry": {
                "type": "Polygon",
                "coordinates": [[[93.69140625, 28.76765910569123], [113.37890625, 28.76765910569123], [113.37890625, 39.30029918615029], [93.69140625, 39.30029918615029], [93.69140625, 28.76765910569123]]]
            }
        }
    }
})

能够指定参数几何对象的坐标系。

db.getCollection('xzbj').find({
    "box": {
        "$geoWithin": {
            "$geometry": {
                "type": "Polygon",
                "coordinates": [[[93.69140625, 28.76765910569123], [113.37890625, 28.76765910569123], [113.37890625, 39.30029918615029], [93.69140625, 39.30029918615029], [93.69140625, 28.76765910569123]]],
                "crs": {
                    "type": "name",
                    "properties": {
                        "name": "urn:x-mongodb:crs:strictwinding:EPSG:4326"
                    }
                }
            }
        }
    }
})

上面使用Polygon做为查询参数,没有索引也能够查。若是使用$center(中心点加半径)和$box查询,则仅有存在2d索引的状况才能查询。
使用$centerSphere做为查询参数的时候,有没有索引均可以查,并且速度很快。能够查询GeoJSON表示的字段,没有类型限制。
$centerShpere的参数也是一个(先经度后纬度)加一个弧度,示例以下:

db.getCollection('xzbj').find({
    "box": {
        "$geoWithin": { "$centerSphere": [ [116.3623809,39.9013085] ,0.008 ] }
    }
})

上的很快就搜索到了,由于我对box字段建了2dSphere索引。把搜索的字段换为geometry就比较慢了,由于没有建索引。

3.2.2 其余的

其余的这里就不记录了,须要的时候查询便可。


下面是我对QGIS转出的geojson进行提取的代码,放在这里作个存档。由于这是简单的使用一次,没有处理各类错误异常等。

#include <iostream>
#include <rapidjson/filereadstream.h>
#include <rapidjson/rapidjson.h>
#include <rapidjson/document.h>
#include <rapidjson/stringbuffer.h>
#include <rapidjson/writer.h>
#include <Poco/Net/HTTPClientSession.h>
#include <Poco/Net/HTTPRequest.h>
#include <Poco/Net/HTTPResponse.h>
#include <Poco/URI.h>
#include <Poco/StreamCopier.h>

int main()
{

    FILE* fp = fopen("./xianjie.geojson","rb");
    char buffer[65536];
    rapidjson::FileReadStream is(fp,buffer,sizeof(buffer));
    rapidjson::Document doc;
    doc.ParseStream(is);
    if(doc.HasParseError()){
        std::cout<<"parse error:"<< doc.GetParseError()<<std::endl;
        return 0;
    }

    rapidjson::Value::ConstMemberIterator iter_features = doc.FindMember("features");
    if(iter_features == doc.MemberEnd()){
        std::cerr<<"没有找到 iter_features"<<std::endl;
        return 0;
    }
    const auto& features = iter_features->value.GetArray();
    // std::cout<<"features size = "<< features.Size()<<std::endl;

    for(size_t i =0;i<features.Size();++i){
        auto& obj = features[i];
        rapidjson::Value::ConstMemberIterator iter_properties = obj.FindMember("properties");
        if(iter_properties == obj.MemberEnd()){
            std::cout<<"features["<<i<<"] format error"<<std::endl;
            continue;
        }
        std::string sheng = iter_properties->value["NL_NAME_1"].GetString();
        size_t pos = sheng.rfind('|');
        if(pos != std::string::npos){
            sheng = sheng.substr(pos+1);
        }
        std::string di = iter_properties->value["NL_NAME_2"].GetString();
        pos = di.find('|');
        if(pos != std::string::npos){
            di = di.substr(0,pos);
        }
        
        std::string xian = di;

        if(iter_properties->value.HasMember("NL_NAME_3") && iter_properties->value["NL_NAME_3"].IsString()){
            xian = iter_properties->value["NL_NAME_3"].GetString();
            pos = xian.find('|');
            if(pos != std::string::npos){
                xian = xian.substr(0,pos);
            }
        }
        std::string fullname = sheng + "-" + di + "-" + x2;
        rapidjson::Value::ConstMemberIterator iter_geometry = obj.FindMember("geometry");
        if(iter_geometry == obj.MemberEnd()){
            std::cout<<"Failed "<<fullname<<std::endl;
            continue;
        }
        std::string geometry;
        {
            rapidjson::StringBuffer buffer;
            rapidjson::Writer<rapidjson::StringBuffer> w(buffer);
            if(iter_geometry->value["coordinates"].GetArray().Size() == 1){
                iter_geometry->value["coordinates"].GetArray()[0].Accept(w);
                geometry.append("{\"type\":\"Polygon\",\"coordinates\":");
                geometry.append(buffer.GetString(),buffer.GetSize());
                geometry.append("}");
            }
            else{
                obj.Accept(w);
                geometry.assign(buffer.GetString(),buffer.GetSize());
            }
        }
        double bbox[4];
        {
            bbox[0] = obj["bbox"].GetArray()[0].GetDouble();
            bbox[1] = obj["bbox"].GetArray()[1].GetDouble();
            bbox[2] = obj["bbox"].GetArray()[2].GetDouble();
            bbox[3] = obj["bbox"].GetArray()[3].GetDouble();
        }
        std::cout<<"{\"name\":\""<<fullname<<"\",\"type\":\"Feature\","
            <<"\"properties\":{\"sheng\":\""<<sheng<<"\",\"di\":\""<<di<<"\",\"xian\":\""<<xian
            <<"\"},\"box\":{\"type\":\"Polygon\",\"coordinates\":[["
            <<"["<<bbox[0]<<","<<bbox[1]<<"],"
            <<"["<<bbox[2]<<","<<bbox[1]<<"],"
            <<"["<<bbox[2]<<","<<bbox[3]<<"],"
            <<"["<<bbox[0]<<","<<bbox[3]<<"],"
            <<"["<<bbox[0]<<","<<bbox[1]<<"]]]}"
            /*<<"\"},\"bbox\":["<<bbox[0]<<","<<bbox[1]<<","<<bbox[2]<<","<<bbox[3]<<"]"*/
            <<",\"geometry\":"<<geometry<<"}\n";

#ifdef CouchDB_Insert
        Poco::Net::HTTPClientSession session("192.168.0.28",5984);
        Poco::Net::HTTPRequest request(Poco::Net::HTTPRequest::HTTP_PUT ,"/xzbj/" + fullname);
        request.setContentType("application/json");
        request.setContentLength((int)geometry.size());
        request.set("Authorization","Basic 用户名密码base64");
        std::ostream& ss = session.sendRequest(request);
        ss.write(geometry.data(),geometry.size());
        Poco::Net::HTTPResponse response;
        std::istream& rs = session.receiveResponse(response);
        std::cout<<"\n\n"<<request.getURI()<<"\t\t"<<response.getStatus()<<std::endl;
        Poco::StreamCopier copier;
        copier.copyStream(rs,std::cout);
#endif

    }
    return 0;
}
相关文章
相关标签/搜索