叠加分析算法对比研究（vector map overlay algorithm）

叠加分析算法对比研究（vector map overlay algorithm）

 ---   by  wangsh

 

1.      介绍：

叠加分析是 GIS中的一项非常重要的空间分析功能。是指在统一空间参考系统下，通过对两个数据进行的一系列集合运算，产生新数据的过程。这里提到的数据可以是图层对应的数据集，也可以是地物对象。叠加分析的目的就是通过对空间数据的加工或分析，提取用户需要的新的空间几何信息。在叠加分析中至少涉及到三个数据，其中一个数据的类型可以是点、线、面等，被称作输入数据（Demo中称作被操作图层/地物）；另一个数据是面数据被称作叠加数据（Demo中称作操作图层/地物）；还有一个数据就是叠加结果数据，包含叠加后数据的几何信息和属性信息。

叠加分析类型包括：

（参见http://i.supermap.com.cn:7080/demo/help/OverlayAnalystHelp_GL.htm介绍）

• 裁剪：裁剪运算是用一个叠加数据从一个源数据中抽取部分特征（点、线、面）集合的运算。
• 擦除：擦除运算是用来擦除掉源数据中与叠加数据中多边形（集合）相重叠部分的操作。
• 同一：同一运算是对两个数据进行相交计算。如果源数据为点数据或者线数据，则新生成的数据中保留源数据的所有对象；如果源数据为面数据，则保留源数据中所有的面对象，并且将与叠加数据相交的对象在相交的地方分割成多个对象。
• 求交：求交运算是求两个数据的交集的操作，两个数据中共同的部分将被输出到结果数据中，其余部分将被排除。
• 对称差：对称差运算是在两个面数据之间进行的交集取反操作。运算的结果会保留两个数据中所有的对象，同时对那些存在于不同的数据中并且相交的面对象进行操作。操作的结果是，对于每一个面对象，去掉其与另一个数据中的几何对象相交的部分，而保留下的剩下的部分。
• 合并：合并运算是求两个数据集的并集的操作，只限于两个面数据之间进行

2.      说明：本文介绍了叠加分析的主要算法进展，没有做过对比试验，仅仅列出开源项目的代码位置和相关的项目信息。

3.      开源项目介绍

 

参考资料1给出了拥有overlay开源项目总结，但是我不同意Blog主的想法，从代码风格来看CGAL固然好，但是GEOS和GRASS的代码也是各有千秋，值得大家学习的。

 

JTS（见参考资料2）是一个成熟的开源GIS项目，Java语言实现的拓扑分析的开发套件。（源码位置参见com.vividsolutions.jts.operation.overlay）

 

GEOS（见参考资料3）是一个成熟的开源GIS项目，是JTS的C++版本，现在已经到了3.1版本，为QGIS等开源GIS平台使用（代码位置geos31/source/geom/Geometry.cpp），参见：

代码参见：

 

Geometry*Geometry::intersection(constGeometry *other)const

{

    returnOverlayOp::overlayOp(this,other,OverlayOp::opINTERSECTION);

}

 

Geometry*Geometry::Union(constGeometry *other)const

    //throw(TopologyException *, IllegalArgumentException *)

{

    Geometry *out =NULL;

 

#ifdefSHORTCIRCUIT_PREDICATES

    // if envelopes are disjoint return a MULTI geom or

    // a geometrycollection

    if ( !getEnvelopeInternal()->intersects(other->getEnvelopeInternal()) )

    {

//cerr<<"SHORTCIRCUITED-UNION engaged"<<endl;

        constGeometryCollection *coll;

        size_tngeoms, i;

        vector<Geometry *> *v =newvector<Geometry *>();

 

        if (NULL != (coll =dynamic_cast<constGeometryCollection *>(this)) )

        {

             ngeoms =coll->getNumGeometries();

             for (i=0;i<ngeoms;i++)

                  v->push_back(coll->getGeometryN(i)->clone());

        }else {

             v->push_back(this->clone());

        }

 

        if (NULL != (coll =dynamic_cast<constGeometryCollection *>(other)) )

        {

             ngeoms =coll->getNumGeometries();

             for (i=0;i<ngeoms;i++)

                  v->push_back(coll->getGeometryN(i)->clone());

        }else {

             v->push_back(other->clone());

        }

 

        out =factory->buildGeometry(v);

        returnout;

    }

#endif

 

    returnOverlayOp::overlayOp(this,other,OverlayOp::opUNION);

}

 

Geometry*Geometry::difference(constGeometry *other)const

    //throw(IllegalArgumentException *)

{

    returnOverlayOp::overlayOp(this,other,OverlayOp::opDIFFERENCE);

}

 

Geometry*Geometry::symDifference(constGeometry *other)const

{

    returnOverlayOp::overlayOp(this,other,OverlayOp::opSYMDIFFERENCE);

}

 

CGAL（见参考资料4）是一个经典的计算几何开源库，重在计算精度，效率一般，模板设计比较好，代码风格一致，是一个计算几何库的经典之作。

 

代码

#include <CGAL/Arr_overlay.h>

template<typename Traits, typename Dcel1, typename Dcel2, typename ResDcel, typenameOverlayTraits>

void

overlay (

Arrangement_2<Traits,Dcel1> arr1,
Arrangement_2<Traits,Dcel2> arr2,
Arrangement_2<Traits,ResDcel>& res,
OverlayTraits& ovl_tr)

 

 

GRASS（见参考资料5）是一个开源的GIS平台项目，是用c语言实现的项目，C实现，具有图层间overlay的能力，代码很复杂，代码位置：grass-6.4.0RC6/lib/vector/Vlib/overlay.c

 

代码：

/*!

  /brief Overlay 2 vector maps and create new one

  /param AMap vector map A

  /param AList unused ?

  /param AAList unused ?

  /param BMap vector map B

  /param BList unused ?

  /param BAList unused ?

  /param operator operator code

  /param OMap output vector map

  /return 0 on success

 */

intVect_overlay(structMap_info *AMap,intatype, struct ilist *AList,struct ilist *AAList,/* map A */

         structMap_info *BMap,int btype, struct ilist *BList,structilist *BAList, /* map B */

         intoperator, struct Map_info *OMap)

{                 /* output map */

   switch (operator) {

   caseGV_O_AND:

    Vect_overlay_and(AMap,atype,AList, AAList, BMap, btype, BList,

              BAList,OMap);

    break;

   default:

    G_fatal_error("Vect_overlay(): %s",_("unknown operator"));

   }

   return 0;

}

      LEDA  是优秀的计算几何库（参考资料6），可惜现在已经不开源了，我这仅有3.3版本，这里不贴代码，只是强调这是一个高效的计算几何源码库。

      Wykobi是一个简单的计算几何开源版本，对于计算小数据量计算几何算子演示非常出色，参见参考资料7.

 

     

      众多与叠加分析算法相关计算几何的优秀源码参考资料8.

 

4.      参考资料

a)        拥有overlay开源项目http://www.scrumo.com/ricepig/?p=93

b)       JTShttp://www.vividsolutions.com/JTS/javadoc/overview-summary.html

c)       GEOShttp://download.osgeo.org/geos/doxygen/classgeos_1_1operation_1_1overlay_1_1OverlayOp.html

d)       CGALhttp://www.cgal.org/Manual/3.3/doc_html/cgal_manual/Arrangement_2_ref/Function_overlay.html

e)        GRASS http://grass.itc.it/grass57/tutorial/vectoroverlay.html

f)         LEDA  http://www.algorithmic-solutions.com/index.htm

g)       Wykobihttp://www.wykobi.com/

http://www.codeproject.com/KB/recipes/Wykobi.aspx

h)       计算几何源码http://compgeom.cs.uiuc.edu/~jeffe/compgeom/code.html

（转载注明来源：http://blog.csdn.net/wsh6759）