GEOS库学习之五：与GDAL/OGR结合使用

要学习GEOS库，肯定绕不开地理方面的东西。如果需要判断的两个多边形或几何图形，不是自己创建的，而是来自shapefile文件，那就得将GEOS库和GDAL/OGR库结合使用了。实际上只需要OGR就行了，但OGR和GDAL是放在一起的。OGR库用来读取和输出shapefile（shp)文件，geos库用来判断空间关系。

需要注意的地方：OGR里面的几何图形类OGRGeometry和GOES里面的类Geometry基本上是一样的，函数也大体相同。OGRGeometry类里面也有类似于Disjoint(),Touches(),Overlaps()这样的函数。但是这些函数是花架子，没有GEOS库的支持是用不成的。如果只有OGR库而没有GEOS库，运行这些函数就会提示GEOS support not enabled.这样的错误，也得不到正确的结果。

在OGR的官方文档中，对Overlaps()函数有这样一句话描述：This method is built on the GEOS library, check it for the definition of the geometry operation. If OGR is built without the GEOS library, this method will always fail, issuing a CPLE_NotSupported error.

大致意思就是这个函数是建立在GEOS库的基础上，没有GEOS库的支持，这个函数运行会出错。

1、编译GEOS。参考http://www.cnblogs.com/denny402/p/4966558.html

2、编译GDAL. 参考http://www.cnblogs.com/sansan/p/3394636.html，注意要先编译GEOS，再编译GDAL，不然会提示找不到geos_c_i.lib这个库。

编译并设置好后，就可以开始代码测试了。

本例的数据及工具版本：

GDAL/OGR: 2.0.0

GEOS: 3.5.0

shp文件：中国国界和省界SHP文件 ，可点此下载

至于SHP文件的读取，可参考 http://www.cnblogs.com/denny402/p/4959867.html

测试代码：

#include "stdafx.h"
#include <iostream>
#include "ogrsf_frmts.h"
using namespace std;


int _tmain(int argc, _TCHAR* argv[])
{
    GDALAllRegister();
    GDALDataset   *poDS;
    CPLSetConfigOption("SHAPE_ENCODING","");  //解决中文乱码问题
    //读取shp文件
    poDS = (GDALDataset*) GDALOpenEx("d:/shp/province.shp", GDAL_OF_VECTOR, NULL, NULL, NULL );

    if( poDS == NULL )
    {
        printf( "Open failed.\n%s" );
        exit( 1 );
    }

    OGRLayer  *poLayer;
    poLayer = poDS->GetLayer(0); //读取层
    poLayer->ResetReading();
    OGRFeature *poFeature1,*poFeature2,*poFeature3;
    poFeature1=poLayer->GetFeature(205); //四川省
    poFeature2=poLayer->GetFeature(0); //黑龙江省
    poFeature3=poLayer->GetFeature(66); //青海省
    OGRGeometry *p1=poFeature1->GetGeometryRef();
    OGRGeometry *p3=poFeature2->GetGeometryRef();
    OGRGeometry *p2=poFeature3->GetGeometryRef();
    cout<<p1->IsEmpty()<<endl   //图形是否为空
        <<p1->IsSimple()<<endl  //是否是单个几何图形
        <<p1->getGeometryType()<<endl   //几何图形的类型，polygon返回3
        <<p1->getGeometryName()<<endl   //几何图形的名称
        <<p1->getDimension()<<endl     //图形的维度
        <<p1->getCoordinateDimension()<<endl   //坐标的维度
        <<p1->getSpatialReference()<<endl;    //空间参考
    if(p2->Disjoint(p1))
        cout<<"不相交"<<endl;
    else
    {
        if(p2->Touches(p1))
            cout<<"接触"<<endl;
        else if(p2->Overlaps(p1))
            cout<<"部分重叠"<<endl;
        else if(p2->Contains(p1))
            cout<<"包含"<<endl;
        else
            cout<<"unknown"<<endl;
    }

    system("pause");
    return 1;
}

测试结果：四川省和青海省接触