[geos-commits] r4021 - in trunk: . capi include/geos/operation include/geos/operation/intersection src/geom src/operation src/operation/intersection tests/unit tests/unit/capi tests/unit/operation tests/unit/operation/intersection
svn_geos at osgeo.org
svn_geos at osgeo.org
Thu Sep 25 03:46:21 PDT 2014
Author: strk
Date: 2014-09-25 03:46:20 -0700 (Thu, 25 Sep 2014)
New Revision: 4021
Added:
trunk/include/geos/operation/intersection/
trunk/include/geos/operation/intersection/Makefile.am
trunk/include/geos/operation/intersection/Rectangle.h
trunk/include/geos/operation/intersection/RectangleIntersection.h
trunk/include/geos/operation/intersection/RectangleIntersectionBuilder.h
trunk/src/operation/intersection/
trunk/src/operation/intersection/Makefile.am
trunk/src/operation/intersection/Rectangle.cpp
trunk/src/operation/intersection/RectangleIntersection.cpp
trunk/src/operation/intersection/RectangleIntersectionBuilder.cpp
trunk/tests/unit/capi/GEOSClipByRectTest.cpp
trunk/tests/unit/operation/intersection/
trunk/tests/unit/operation/intersection/RectangleIntersectionTest.cpp
Modified:
trunk/NEWS
trunk/capi/geos_c.cpp
trunk/capi/geos_c.h.in
trunk/capi/geos_ts_c.cpp
trunk/configure.ac
trunk/include/geos/operation/Makefile.am
trunk/src/geom/Geometry.cpp
trunk/src/operation/Makefile.am
trunk/tests/unit/Makefile.am
Log:
Add optimized RectangleIntersection functionality
Includes:
C++ API, with tests
C-API GEOSClipByRect, with tests
Initial C++ code provided by Mika Heiskanen.
Modified by me to work with arbitrarily ordered polygon ring vertices.
See #699 for background
Modified: trunk/NEWS
===================================================================
--- trunk/NEWS 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/NEWS 2014-09-25 10:46:20 UTC (rev 4021)
@@ -6,6 +6,7 @@
- PHP: Geometry->normalize method
- GEOS_USE_ONLY_R_API macro support (#695)
- PHP: WKBReader->read() & WKBWriter::write() methods (Benjamin Morel)
+ - RectangleIntersection, GEOSClipByRect (Mika Heiskanen, Sandro Santilli)
...
- Improvements:
- Speed-up intersection and difference between geometries
Modified: trunk/capi/geos_c.cpp
===================================================================
--- trunk/capi/geos_c.cpp 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/capi/geos_c.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -485,6 +485,11 @@
}
+Geometry *
+GEOSClipByRect(const Geometry *g, double xmin, double ymin, double xmax, double ymax)
+{
+ return GEOSClipByRect_r( handle, g, xmin, ymin, xmax, ymax );
+}
Modified: trunk/capi/geos_c.h.in
===================================================================
--- trunk/capi/geos_c.h.in 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/capi/geos_c.h.in 2014-09-25 10:46:20 UTC (rev 4021)
@@ -511,6 +511,12 @@
const GEOSGeometry* g);
extern GEOSGeometry GEOS_DLL *GEOSNode_r(GEOSContextHandle_t handle,
const GEOSGeometry* g);
+/* Fast, non-robust intersection between an arbitrary geometry and
+ * a rectangle. The returned geometry may be invalid. */
+extern GEOSGeometry GEOS_DLL *GEOSClipByRect_r(GEOSContextHandle_t handle,
+ const GEOSGeometry* g,
+ double xmin, double ymin,
+ double xmax, double ymax);
/*
* all arguments remain ownership of the caller
@@ -1366,6 +1372,7 @@
extern GEOSGeometry GEOS_DLL *GEOSPointOnSurface(const GEOSGeometry* g);
extern GEOSGeometry GEOS_DLL *GEOSGetCentroid(const GEOSGeometry* g);
extern GEOSGeometry GEOS_DLL *GEOSNode(const GEOSGeometry* g);
+extern GEOSGeometry GEOS_DLL *GEOSClipByRect(const GEOSGeometry* g, double xmin, double ymin, double xmax, double ymax);
/*
* all arguments remain ownership of the caller
Modified: trunk/capi/geos_ts_c.cpp
===================================================================
--- trunk/capi/geos_ts_c.cpp 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/capi/geos_ts_c.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -54,6 +54,8 @@
#include <geos/operation/linemerge/LineMerger.h>
#include <geos/operation/overlay/OverlayOp.h>
#include <geos/operation/overlay/snap/GeometrySnapper.h>
+#include <geos/operation/intersection/Rectangle.h>
+#include <geos/operation/intersection/RectangleIntersection.h>
#include <geos/operation/polygonize/Polygonizer.h>
#include <geos/operation/relate/RelateOp.h>
#include <geos/operation/sharedpaths/SharedPathsOp.h>
@@ -2158,6 +2160,48 @@
return NULL;
}
+Geometry *
+GEOSClipByRect_r(GEOSContextHandle_t extHandle, const Geometry *g, double xmin, double ymin, double xmax, double ymax)
+{
+ if ( 0 == extHandle )
+ {
+ return NULL;
+ }
+
+ GEOSContextHandleInternal_t *handle = 0;
+ handle = reinterpret_cast<GEOSContextHandleInternal_t*>(extHandle);
+ if ( 0 == handle->initialized )
+ {
+ return NULL;
+ }
+
+ try
+ {
+ using geos::operation::intersection::Rectangle;
+ using geos::operation::intersection::RectangleIntersection;
+ Rectangle rect(xmin, ymin, xmax, ymax);
+ std::auto_ptr<Geometry> g3 = RectangleIntersection::clip(*g, rect);
+ return g3.release();
+ }
+ catch (const std::exception &e)
+ {
+#if VERBOSE_EXCEPTIONS
+ std::ostringstream s;
+ s << "Exception on GEOSClipByRect with following inputs:" << std::endl;
+ s << "A: "<<g1->toString() << std::endl;
+ s << "B: "<<g2->toString() << std::endl;
+ handle->NOTICE_MESSAGE("%s", s.str().c_str());
+#endif // VERBOSE_EXCEPTIONS
+ handle->ERROR_MESSAGE("%s", e.what());
+ }
+ catch (...)
+ {
+ handle->ERROR_MESSAGE("Unknown exception thrown");
+ }
+
+ return NULL;
+}
+
//-------------------------------------------------------------------
// memory management functions
//------------------------------------------------------------------
Modified: trunk/configure.ac
===================================================================
--- trunk/configure.ac 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/configure.ac 2014-09-25 10:46:20 UTC (rev 4021)
@@ -453,6 +453,7 @@
include/geos/operation/Makefile
include/geos/operation/buffer/Makefile
include/geos/operation/distance/Makefile
+ include/geos/operation/intersection/Makefile
include/geos/operation/linemerge/Makefile
include/geos/operation/overlay/Makefile
include/geos/operation/overlay/snap/Makefile
@@ -484,6 +485,7 @@
src/operation/Makefile
src/operation/buffer/Makefile
src/operation/distance/Makefile
+ src/operation/intersection/Makefile
src/operation/linemerge/Makefile
src/operation/overlay/Makefile
src/operation/polygonize/Makefile
Modified: trunk/include/geos/operation/Makefile.am
===================================================================
--- trunk/include/geos/operation/Makefile.am 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/include/geos/operation/Makefile.am 2014-09-25 10:46:20 UTC (rev 4021)
@@ -4,6 +4,7 @@
SUBDIRS = \
buffer \
distance \
+ intersection \
linemerge \
overlay \
polygonize \
Added: trunk/include/geos/operation/intersection/Makefile.am
===================================================================
--- trunk/include/geos/operation/intersection/Makefile.am (rev 0)
+++ trunk/include/geos/operation/intersection/Makefile.am 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,13 @@
+#
+# This file is part of project GEOS (http://trac.osgeo.org/geos/)
+#
+#SUBDIRS =
+
+#EXTRA_DIST =
+
+geosdir = $(includedir)/geos/operation/intersection
+
+geos_HEADERS = \
+ Rectangle.h \
+ RectangleIntersection.h \
+ RectangleIntersectionBuilder.h
Added: trunk/include/geos/operation/intersection/Rectangle.h
===================================================================
--- trunk/include/geos/operation/intersection/Rectangle.h (rev 0)
+++ trunk/include/geos/operation/intersection/Rectangle.h 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,209 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2014 Mika Heiskanen <mika.heiskanen at fmi.fi>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#ifndef GEOS_OP_INTERSECTION_RECTANGLE_H
+#define GEOS_OP_INTERSECTION_RECTANGLE_H
+
+#include <geos/export.h>
+
+#ifdef _MSC_VER
+#pragma warning(push)
+#pragma warning(disable: 4251) // warning C4251: needs to have dll-interface to be used by clients of class
+#endif
+
+// Forward declarations
+namespace geos {
+ namespace geom {
+ class GeometryFactory;
+ class Geometry;
+ class Polygon;
+ class LinearRing;
+ }
+}
+
+namespace geos {
+namespace operation { // geos::operation
+namespace intersection { // geos::operation::intersection
+
+/**
+ * \brief Clipping rectangle
+ *
+ * A clipping rectangle defines the boundaries of the rectangle
+ * by defining the limiting x- and y-coordinates. The clipping
+ * rectangle must be non-empty. In addition, methods are provided
+ * for specifying the location of a given coordinate with respect
+ * to the clipping rectangle similarly to the Cohen-Sutherland
+ * clipping algorithm.
+ *
+ */
+
+class GEOS_DLL Rectangle
+{
+ public:
+
+ /**
+ * \brief Construct a clipping rectangle
+ *
+ * @param x1 x-coordinate of the left edge
+ * @param y1 y-coordinate of the bottom edge
+ * @param x2 x-coordinate of the right edge
+ * @param y2 y-coordinate of the right edge
+ * @throws IllegalArgumentException if the rectangle is empty
+ */
+
+ Rectangle(double x1, double y1, double x2, double y2);
+
+ /**
+ * \@return the minimum x-coordinate of the rectangle
+ */
+ double xmin() const { return xMin; }
+
+ /**
+ * \@return the minimum y-coordinate of the rectangle
+ */
+
+ double ymin() const { return yMin; }
+
+
+ /**
+ * \@return the maximum x-coordinate of the rectangle
+ */
+
+ double xmax() const { return xMax; }
+
+
+ /**
+ * \@return the maximum y-coordinate of the rectangle
+ */
+
+ double ymax() const { return yMax; }
+
+ /**
+ * \@return the rectangle as a polygon geometry
+ *
+ * Ownership transferred to caller
+ */
+ geom::Polygon* toPolygon(const geom::GeometryFactory &f) const;
+
+ geom::LinearRing* toLinearRing(const geom::GeometryFactory &f) const;
+
+ /**
+ * @brief Position with respect to a clipping rectangle
+ */
+
+ enum Position
+ {
+ Inside = 1,
+ Outside = 2,
+
+ Left = 4,
+ Top = 8,
+ Right = 16,
+ Bottom = 32,
+
+ TopLeft = Top|Left, // 12
+ TopRight = Top|Right, // 24
+ BottomLeft = Bottom|Left, // 36
+ BottomRight = Bottom|Right // 48
+ };
+
+ /**
+ * @brief Test if the given position is on a {@link Rectangle] edge
+ * @param pos {@link Rectangle} {@link Position}
+ * @return true, if the position is on an edge
+ */
+
+ static bool onEdge(Position pos)
+ {
+ return (pos > Outside);
+ }
+
+ /**
+ * @brief Test if the given positions are on the same {@link Rectangle} edge
+ * @param pos1 {@link Rectangle} {@link Position} of first coordinate
+ * @param pos2 {@link Rectangle} {@link Position} of second coordinate
+ * @return true, if the positions are on the same edge
+ */
+
+ static bool onSameEdge(Position pos1, Position pos2)
+ {
+ return onEdge(Position(pos1 & pos2));
+ }
+
+ /**
+ * @brief Establish position of coordinate with respect to a {@link Rectangle}
+ * @param x x-coordinate
+ * @param y y-coordinate
+ * @return {@link Position} of the coordinate
+ */
+
+ Position position(double x, double y) const
+ {
+ // We assume the point to be inside and test it first
+ if(x>xMin && x<xMax && y>yMin && y<yMax)
+ return Inside;
+ // Next we assume the point to be outside and test it next
+ if(x<xMin || x>xMax || y<yMin || y>yMax)
+ return Outside;
+ // Slower cases
+ unsigned int pos = 0;
+ if(x==xMin)
+ pos |= Left;
+ else if(x==xMax)
+ pos |= Right;
+ if(y==yMin)
+ pos |= Bottom;
+ else if(y==yMax)
+ pos |= Top;
+ return Position(pos);
+ }
+
+ /**
+ * @brief Next edge in clock-wise direction
+ * @param pos {@link Rectangle} {@link Position}
+ * @return next {@Rectangle} {@link Position} in clock-wise direction
+ */
+
+ static Position nextEdge(Position pos)
+ {
+ switch(pos)
+ {
+ case BottomLeft:
+ case Left: return Top;
+ case TopLeft:
+ case Top: return Right;
+ case TopRight:
+ case Right: return Bottom;
+ case BottomRight:
+ case Bottom: return Left;
+ /* silences compiler warnings, Inside & Outside are not handled explicitly */
+ default: return pos;
+ }
+ }
+
+ private:
+
+ Rectangle();
+ double xMin;
+ double yMin;
+ double xMax;
+ double yMax;
+
+}; // class RectangleIntersection
+
+} // namespace geos::operation::intersection
+} // namespace geos::operation
+} // namespace geos
+
+#endif // GEOS_OP_INTERSECTION_RECTANGLE_H
Added: trunk/include/geos/operation/intersection/RectangleIntersection.h
===================================================================
--- trunk/include/geos/operation/intersection/RectangleIntersection.h (rev 0)
+++ trunk/include/geos/operation/intersection/RectangleIntersection.h 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,179 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2014 Mika Heiskanen <mika.heiskanen at fmi.fi>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#ifndef GEOS_OP_RECTANGLE_INTERSECTION_H
+#define GEOS_OP_RECTANGLE_INTERSECTION_H
+
+#include <geos/export.h>
+
+#include <memory>
+
+#ifdef _MSC_VER
+#pragma warning(push)
+#pragma warning(disable: 4251) // warning C4251: needs to have dll-interface to be used by clients of class
+#endif
+
+// Forward declarations
+namespace geos {
+ namespace geom {
+ class Point;
+ class MultiPoint;
+ class Polygon;
+ class MultiPolygon;
+ class LineString;
+ class MultiLineString;
+ class Geometry;
+ class GeometryCollection;
+ class GeometryFactory;
+ class CoordinateSequenceFactory;
+ }
+ namespace operation {
+ namespace intersection {
+ class Rectangle;
+ class RectangleIntersectionBuilder;
+ }
+ }
+}
+
+namespace geos {
+namespace operation { // geos::operation
+namespace intersection { // geos::operation::intersection
+
+/**
+ * \brief Speed-optimized clipping of a {@link Geometry} with a rectangle.
+ *
+ * Two different methods are provided. The first performs normal
+ * clipping, the second clips the boundaries of polygons, not
+ * the polygons themselves. In the first case a polygon will remain
+ * a polygon or is completely cut out. In the latter case polygons
+ * will be converted to polylines if any vertex is outside the clipping
+ * rectangle, or will be cut out completely.
+ *
+ * The algorithm works best when the number of intersections is very low.
+ * For example, if the geometry is completely to the left of the
+ * clipping rectangle, only the x-coordinate of the geometry is ever
+ * tested and is only compared with the x-coordinate of the left edge
+ * of the rectangle. Hence clipping may be faster than calculating
+ * the envelope of the geometry for trivial overlap tests.
+ *
+ * The input geometry must be valid. In particular all {@link LinearRing}s must
+ * be properly closed, or the algorithm may not terminate.
+ *
+ */
+class GEOS_DLL RectangleIntersection
+{
+ public:
+
+ /**
+ * \brief Clip geometry with a rectangle
+ *
+ * @param geom a {@link Geometry}
+ * @param rect a {@link Rectangle}
+ * @return the clipped geometry
+ * @return NULL if the geometry is outside the {@link Rectangle}
+ */
+ static std::auto_ptr<geom::Geometry> clip(const geom::Geometry & geom,
+ const Rectangle & rect);
+
+ /**
+ * \brief Clip boundary of a geometry with a rectangle
+ *
+ *
+ * Any polygon which intersects the rectangle will be converted to
+ * a polyline or a multipolyline - including the holes.
+ *
+ * @param geom a {@link Geometry}
+ * @param rect a {@link Rectangle}
+ * @return the clipped geometry
+ * @return NULL if the geometry is outside the {@link Rectangle}
+ */
+ static std::auto_ptr<geom::Geometry> clipBoundary(const geom::Geometry & geom,
+ const Rectangle & rect);
+
+private:
+
+ RectangleIntersection(const geom::Geometry& geom, const Rectangle& rect);
+
+ std::auto_ptr<geom::Geometry> clipBoundary();
+
+ std::auto_ptr<geom::Geometry> clip();
+
+ const geom::Geometry &_geom;
+ const Rectangle &_rect;
+ const geom::GeometryFactory *_gf;
+ const geom::CoordinateSequenceFactory *_csf;
+
+ void clip_geom(const geom::Geometry * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons);
+
+ void clip_point(const geom::Point * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect);
+
+ void clip_multipoint(const geom::MultiPoint * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect);
+
+ void clip_linestring(const geom::LineString * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect);
+
+ void clip_multilinestring(const geom::MultiLineString * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect);
+
+ void clip_polygon(const geom::Polygon * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons);
+
+ void clip_multipolygon(const geom::MultiPolygon * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons);
+
+ void clip_geometrycollection(
+ const geom::GeometryCollection * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons);
+
+ void clip_polygon_to_linestrings(const geom::Polygon * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect);
+
+ void clip_polygon_to_polygons(const geom::Polygon * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect);
+
+
+ /**
+ * \brief Clip geometry.
+ *
+ * Returns true if the geometry was fully inside, and does not output
+ * anything to RectangleIntersectionBuilder.
+ */
+ bool clip_linestring_parts(const geom::LineString * gi,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect);
+
+}; // class RectangleIntersection
+
+} // namespace geos::operation::intersection
+} // namespace geos::operation
+} // namespace geos
+
+#endif // GEOS_OP_RECTANGLE_INTERSECTION_H
Added: trunk/include/geos/operation/intersection/RectangleIntersectionBuilder.h
===================================================================
--- trunk/include/geos/operation/intersection/RectangleIntersectionBuilder.h (rev 0)
+++ trunk/include/geos/operation/intersection/RectangleIntersectionBuilder.h 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,161 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2014 Mika Heiskanen <mika.heiskanen at fmi.fi>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#ifndef GEOS_OP_INTERSECTION_RECTANGLEINTERSECTIONBUILDER_H
+#define GEOS_OP_INTERSECTION_RECTANGLEINTERSECTIONBUILDER_H
+
+#include <geos/export.h>
+
+#ifdef _MSC_VER
+#pragma warning(push)
+#pragma warning(disable: 4251) // warning C4251: needs to have dll-interface to be used by clients of class
+#endif
+
+#include <memory>
+#include <list>
+#include <vector>
+
+
+// Forward declarations
+namespace geos {
+ namespace geom {
+ class Coordinate;
+ class Geometry;
+ class GeometryFactory;
+ class Polygon;
+ class LineString;
+ class Point;
+ }
+ namespace operation {
+ namespace intersection {
+ class Rectangle;
+ }
+ }
+}
+
+namespace geos {
+namespace operation { // geos::operation
+namespace intersection { // geos::operation::intersection
+
+/**
+ * \brief Rebuild geometries from subpaths left by clipping with a rectangle
+ *
+ * The RectangleIntersectionBuilder is used to maintain lists of polygons,
+ * linestrings and points left from clipping a {@link Geometry} with a {@link Rectangle}.
+ * Once all clipping has been done, the class builds a valid {@link Geometry} from
+ * the components.
+ *
+ * This is a utility class needed by {@link RectangleIntersection}, and is not
+ * intended for public use.
+ */
+
+class GEOS_DLL RectangleIntersectionBuilder
+{
+ // Regular users are not supposed to use this utility class.
+ friend class RectangleIntersection;
+
+public:
+
+ ~RectangleIntersectionBuilder();
+
+private:
+
+ /**
+ * \brief Build the result geometry from partial results and clean up
+ */
+ std::auto_ptr<geom::Geometry> build();
+
+ /**
+ * \brief Build polygons from parts left by clipping one
+ *
+ * 1. Build exterior ring(s) from lines
+ * 2. Attach polygons as holes to the exterior ring(s)
+ */
+ void reconnectPolygons(const Rectangle & rect);
+
+ /**
+ * Reconnect disjointed parts
+ *
+ * When we clip a LinearRing we may get multiple linestrings.
+ * Often the first and last ones can be reconnected to simplify
+ * output.
+ *
+ * Sample clip with a rectangle 0,0 --> 10,10 without reconnecting:
+ *
+ * Input: POLYGON ((5 10,0 0,10 0,5 10))
+ * Output: MULTILINESTRING ((5 10,0 0),(10 0,5 10))
+ * Desired: LINESTRING (10 0,5 10,0 0)
+ *
+ * TODO: If there is a very sharp spike from inside the rectangle
+ * outside, and then back in, it is possible that the
+ * intersection points at the edge are equal. In this
+ * case we could reconnect the linestrings. The task is
+ * the same we're already doing for the 1st/last linestrings,
+ * we'd just do it for any adjacent pair as well.
+ */
+ void reconnect();
+
+ void reverseLines();
+
+ /**
+ * Export parts to another container
+ */
+ void release(RectangleIntersectionBuilder & parts);
+
+ // Adding Geometry components
+ void add(geom::Polygon * g);
+ void add(geom::LineString * g);
+ void add(geom::Point * g);
+
+ // Trivial methods
+ bool empty() const;
+ void clear();
+
+ // Added components
+ std::list<geom::Polygon *> polygons;
+ std::list<geom::LineString *> lines;
+ std::list<geom::Point *> points;
+
+ /**
+ * \brief Close a ring clockwise along rectangle edges
+ *
+ * Only the 4 corners and x1,y1 need to be considered. The possible
+ * cases are:
+ *
+ * x1,y1
+ * corner1 x1,y1
+ * corner1 corner2 x1,y1
+ * corner1 corner2 corner3 x1,y1
+ * corner1 corner2 corner3 corner4 x1,y1
+ */
+ void close_boundary(
+ const Rectangle & rect,
+ std::vector<geom::Coordinate> * ring,
+ double x1, double y1,
+ double x2, double y2);
+
+ void close_ring(const Rectangle & rect, std::vector<geom::Coordinate> * ring);
+
+ RectangleIntersectionBuilder(const geom::GeometryFactory& f)
+ : _gf(f) {}
+
+ const geom::GeometryFactory &_gf;
+
+}; // class RectangleIntersectionBuilder
+
+} // namespace geos::operation::intersection
+} // namespace geos::operation
+} // namespace geos
+
+#endif // GEOS_OP_INTERSECTION_RECTANGLEINTERSECTIONBUILDER_H
Modified: trunk/src/geom/Geometry.cpp
===================================================================
--- trunk/src/geom/Geometry.cpp 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/src/geom/Geometry.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -39,6 +39,8 @@
#include <geos/algorithm/InteriorPointLine.h>
#include <geos/algorithm/InteriorPointArea.h>
#include <geos/algorithm/ConvexHull.h>
+#include <geos/operation/intersection/Rectangle.h>
+#include <geos/operation/intersection/RectangleIntersection.h>
#include <geos/operation/predicate/RectangleContains.h>
#include <geos/operation/predicate/RectangleIntersects.h>
#include <geos/operation/relate/RelateOp.h>
@@ -67,6 +69,7 @@
#endif
#define SHORTCIRCUIT_PREDICATES 1
+//#define USE_RECTANGLE_INTERSECTION 1
using namespace std;
using namespace geos::algorithm;
@@ -500,8 +503,8 @@
Geometry::intersection(const Geometry *other) const
{
/**
- * TODO: MD - add optimization for P-A case using Point-In-Polygon
- */
+ * TODO: MD - add optimization for P-A case using Point-In-Polygon
+ */
// special case: if one input is empty ==> empty
if (isEmpty() || other->isEmpty() )
@@ -509,6 +512,24 @@
return getFactory()->createGeometryCollection();
}
+#ifdef USE_RECTANGLE_INTERSECTION
+ // optimization for rectangle arguments
+ using operation::intersection::Rectangle;
+ using operation::intersection::RectangleIntersection;
+ if ( isRectangle() ) {
+ const Envelope* env = getEnvelopeInternal();
+ Rectangle rect(env->getMinX(), env->getMinY(),
+ env->getMaxX(), env->getMaxY());
+ return RectangleIntersection::clip(*other, rect).release();
+ }
+ if (other->isRectangle()) {
+ const Envelope* env = other->getEnvelopeInternal();
+ Rectangle rect(env->getMinX(), env->getMinY(),
+ env->getMaxX(), env->getMaxY());
+ return RectangleIntersection::clip(*this, rect).release();
+ }
+#endif
+
return BinaryOp(this, other, overlayOp(OverlayOp::opINTERSECTION)).release();
}
Modified: trunk/src/operation/Makefile.am
===================================================================
--- trunk/src/operation/Makefile.am 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/src/operation/Makefile.am 2014-09-25 10:46:20 UTC (rev 4021)
@@ -4,6 +4,7 @@
SUBDIRS = \
buffer \
distance \
+ intersection \
linemerge \
overlay \
polygonize \
@@ -24,6 +25,7 @@
liboperation_la_LIBADD = \
buffer/libopbuffer.la \
distance/libopdistance.la \
+ intersection/libopintersection.la \
linemerge/liboplinemerge.la \
overlay/libopoverlay.la \
polygonize/liboppolygonize.la \
Added: trunk/src/operation/intersection/Makefile.am
===================================================================
--- trunk/src/operation/intersection/Makefile.am (rev 0)
+++ trunk/src/operation/intersection/Makefile.am 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,15 @@
+#
+# This file is part of project GEOS (http://trac.osgeo.org/geos/)
+#
+SUBDIRS =
+
+noinst_LTLIBRARIES = libopintersection.la
+
+AM_CPPFLAGS = -I$(top_srcdir)/include
+
+libopintersection_la_SOURCES = \
+ Rectangle.cpp \
+ RectangleIntersection.cpp \
+ RectangleIntersectionBuilder.cpp
+
+libopintersection_la_LIBADD =
Added: trunk/src/operation/intersection/Rectangle.cpp
===================================================================
--- trunk/src/operation/intersection/Rectangle.cpp (rev 0)
+++ trunk/src/operation/intersection/Rectangle.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,65 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2014 Mika Heiskanen <mika.heiskanen at fmi.fi>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#include <geos/operation/intersection/Rectangle.h>
+#include <geos/util/IllegalArgumentException.h>
+#include <geos/geom/Polygon.h>
+#include <geos/geom/GeometryFactory.h>
+#include <geos/geom/CoordinateSequenceFactory.h>
+#include <geos/geom/CoordinateSequence.h>
+#include <geos/geom/Coordinate.h>
+
+namespace geos {
+namespace operation { // geos::operation
+namespace intersection { // geos::operation::intersection
+
+ /*
+ * Create a clipping rectangle
+ */
+
+ Rectangle::Rectangle(double x1, double y1, double x2, double y2)
+ : xMin(x1)
+ , yMin(y1)
+ , xMax(x2)
+ , yMax(y2)
+ {
+ if(xMin >= xMax || yMin >= yMax)
+ {
+ throw util::IllegalArgumentException("Clipping rectangle must be non-empty");
+ }
+ }
+
+ geom::Polygon*
+ Rectangle::toPolygon(const geom::GeometryFactory &f) const
+ {
+ geom::LinearRing* ls = toLinearRing(f);
+ return f.createPolygon(ls, 0);
+ }
+
+ geom::LinearRing*
+ Rectangle::toLinearRing(const geom::GeometryFactory &f) const
+ {
+ const geom::CoordinateSequenceFactory *csf = f.getCoordinateSequenceFactory();
+ geom::CoordinateSequence *seq = csf->create(5, 2);
+ seq->setAt(geom::Coordinate(xMin, yMin), 0);
+ seq->setAt(geom::Coordinate(xMin, yMax), 1);
+ seq->setAt(geom::Coordinate(xMax, yMax), 2);
+ seq->setAt(geom::Coordinate(xMax, yMin), 3);
+ seq->setAt(seq->getAt(0), 4); // close
+ return f.createLinearRing(seq);
+ }
+
+} // namespace geos::operation::intersection
+} // namespace geos::operation
+} // namespace geos
Added: trunk/src/operation/intersection/RectangleIntersection.cpp
===================================================================
--- trunk/src/operation/intersection/RectangleIntersection.cpp (rev 0)
+++ trunk/src/operation/intersection/RectangleIntersection.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,704 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2014 Mika Heiskanen <mika.heiskanen at fmi.fi>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#include <geos/algorithm/CGAlgorithms.h>
+#include <geos/operation/intersection/RectangleIntersection.h>
+#include <geos/operation/intersection/Rectangle.h>
+#include <geos/operation/intersection/RectangleIntersectionBuilder.h>
+#include <geos/operation/predicate/RectangleIntersects.h>
+#include <geos/geom/GeometryFactory.h>
+#include <geos/geom/CoordinateSequenceFactory.h>
+#include <geos/geom/CoordinateSequence.h>
+#include <geos/geom/Polygon.h>
+#include <geos/geom/MultiPolygon.h>
+#include <geos/geom/Point.h>
+#include <geos/geom/Geometry.h>
+#include <geos/geom/LineString.h>
+#include <geos/geom/LinearRing.h>
+#include <geos/geom/Location.h>
+#include <geos/util/UnsupportedOperationException.h>
+#include <list>
+#include <stdexcept>
+
+using geos::operation::intersection::Rectangle;
+using geos::operation::intersection::RectangleIntersectionBuilder;
+using namespace geos::geom;
+namespace geos {
+namespace operation { // geos::operation
+namespace intersection { // geos::operation::intersection
+
+/**
+ * \brief Test if two coordinates are different
+ */
+
+inline
+bool different(double x1, double y1, double x2, double y2)
+{
+ return !(x1==x2 && y1==y2);
+}
+
+/**
+ * \brief Calculate a line intersection point
+ *
+ * Note:
+ * - Calling this with x1,y1 and x2,y2 swapped cuts the other end of the line
+ * - Calling this with x and y swapped cuts in y-direction instead
+ * - Calling with 1<->2 and x<->y swapped works too
+ */
+
+inline
+void clip_one_edge(double & x1, double & y1, double x2, double y2, double limit)
+{
+ if(x1 != x2)
+ {
+ y1 += (y2-y1)*(limit-x1)/(x2-x1);
+ x1 = limit;
+ }
+}
+
+/**
+ * \brief Start point is outside, endpoint is definitely inside
+ *
+ * Note: Even though one might think using >= etc operators would produce
+ * the same result, that is not the case. We rely on the fact
+ * that nothing is clipped unless the point is truly outside the
+ * rectangle! Without this handling lines ending on the edges of
+ * the rectangle would be very difficult.
+ */
+
+void clip_to_edges(double & x1, double & y1,
+ double x2, double y2,
+ const Rectangle & rect)
+{
+ if(x1 < rect.xmin())
+ clip_one_edge(x1,y1,x2,y2,rect.xmin());
+ else if(x1 > rect.xmax())
+ clip_one_edge(x1,y1,x2,y2,rect.xmax());
+
+ if(y1<rect.ymin())
+ clip_one_edge(y1,x1,y2,x2,rect.ymin());
+ else if(y1>rect.ymax())
+ clip_one_edge(y1,x1,y2,x2,rect.ymax());
+}
+
+
+/**
+ * \brief Clip geometry
+ *
+ * Here outGeom may also be a MultiPoint
+ */
+
+void
+RectangleIntersection::clip_point(const geom::Point * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect)
+{
+ if(g == NULL)
+ return;
+
+ double x = g->getX();
+ double y = g->getY();
+
+ if(rect.position(x,y) == Rectangle::Inside)
+ parts.add(dynamic_cast<geom::Point*>(g->clone()));
+}
+
+bool
+RectangleIntersection::clip_linestring_parts(const geom::LineString * gi,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect)
+{
+ using namespace geos::geom;
+
+ int n = gi->getNumPoints();
+
+ if(gi == NULL || n<1)
+ return false;
+
+ // For shorthand code
+
+ std::vector<Coordinate> cs;
+ gi->getCoordinatesRO()->toVector(cs);
+ //const geom::CoordinateSequence &cs = *(gi->getCoordinatesRO());
+
+ // Keep a record of the point where a line segment entered
+ // the rectangle. If the boolean is set, we must insert
+ // the point to the beginning of the linestring which then
+ // continues inside the rectangle.
+
+ double x0 = 0;
+ double y0 = 0;
+ bool add_start = false;
+
+ // Start iterating
+
+ int i = 0;
+
+ while(i<n)
+ {
+ // Establish initial position
+
+ double x = cs[i].x;
+ double y = cs[i].y;
+ Rectangle::Position pos = rect.position(x,y);
+
+ if(pos == Rectangle::Outside)
+ {
+ // Skip points as fast as possible until something has to be checked
+ // in more detail.
+
+ ++i; // we already know it is outside
+
+ if(x < rect.xmin())
+ while(i < n && cs[i].x < rect.xmin())
+ ++i;
+
+ else if(x > rect.xmax())
+ while(i < n && cs[i].x > rect.xmax())
+ ++i;
+
+ else if(y < rect.ymin())
+ while(i < n && cs[i].y < rect.ymin())
+ ++i;
+
+ else if(y > rect.ymax())
+ while(i < n && cs[i].y > rect.ymax())
+ ++i;
+
+ if(i >= n)
+ return false;
+
+ // Establish new position
+ x = cs[i].x;
+ y = cs[i].y;
+ pos = rect.position(x,y);
+
+ // Handle all possible cases
+ x0 = cs[i-1].x;
+ y0 = cs[i-1].y;
+ clip_to_edges(x0,y0,x,y,rect);
+
+ if(pos == Rectangle::Inside)
+ {
+ add_start = true; // x0,y0 must have clipped the rectangle
+ // Main loop will enter the Inside/Edge section
+
+ }
+ else if(pos == Rectangle::Outside)
+ {
+ // From Outside to Outside. We need to check whether
+ // we created a line segment inside the box. In any
+ // case, we will continue the main loop after this
+ // which will then enter the Outside section.
+
+ // Clip the other end too
+ clip_to_edges(x,y,x0,y0,rect);
+
+ Rectangle::Position prev_pos = rect.position(x0,y0);
+ pos = rect.position(x,y);
+
+ if( different(x0,y0,x,y) && // discard corners etc
+ Rectangle::onEdge(prev_pos) && // discard if does not intersect rect
+ Rectangle::onEdge(pos) &&
+ !Rectangle::onSameEdge(prev_pos,pos) // discard if travels along edge
+ )
+ {
+ std::vector<Coordinate> *coords = new std::vector<Coordinate>(2);
+ (*coords)[0] = Coordinate(x0,y0);
+ (*coords)[1] = Coordinate(x,y);
+ CoordinateSequence *seq = _csf->create(coords);
+ geom::LineString * line = _gf->createLineString(seq);
+ parts.add(line);
+ }
+
+ // Continue main loop outside the rect
+
+ }
+ else
+ {
+ // From outside to edge. If the edge we hit first when
+ // following the line is not the edge we end at, then
+ // clearly we must go through the rectangle and hence
+ // a start point must be set.
+
+ Rectangle::Position newpos = rect.position(x0,y0);
+ if(!Rectangle::onSameEdge(pos,newpos))
+ {
+ add_start = true;
+ }
+ else
+ {
+ // we ignore the travel along the edge and continue
+ // the main loop at the last edge point
+ }
+ }
+ }
+
+ else
+ {
+ // The point is now stricly inside or on the edge.
+ // Keep iterating until the end or the point goes
+ // outside. We may have to output partial linestrings
+ // while iterating until we go strictly outside
+
+ int start_index = i; // 1st valid original point
+ bool go_outside = false;
+
+ while(!go_outside && ++i<n)
+ {
+ x = cs[i].x;
+ y = cs[i].y;
+
+ Rectangle::Position prev_pos = pos;
+ pos = rect.position(x,y);
+
+ if(pos == Rectangle::Inside)
+ {
+ // Just keep going
+ }
+ else if(pos == Rectangle::Outside)
+ {
+ go_outside = true;
+
+ // Clip the outside point to edges
+ clip_to_edges(x, y, cs[i-1].x, cs[i-1].y, rect);
+ pos = rect.position(x,y);
+
+ // Does the line exit through the inside of the box?
+
+ bool through_box = (different(x,y,cs[i].x,cs[i].y) &&
+ !Rectangle::onSameEdge(prev_pos,pos));
+
+ // Output a LineString if it at least one segment long
+
+ if(start_index < i-1 || add_start || through_box)
+ {
+ std::vector<Coordinate> *coords = new std::vector<Coordinate>();
+ if(add_start)
+ {
+ coords->push_back(Coordinate(x0, y0));
+ add_start = false;
+ }
+ //line->addSubLineString(&g, start_index, i-1);
+ coords->insert(coords->end(), cs.begin()+start_index, cs.begin()+i);
+
+ if(through_box) coords->push_back(Coordinate(x,y));
+
+ CoordinateSequence *seq = _csf->create(coords);
+ geom::LineString * line = _gf->createLineString(seq);
+ parts.add(line);
+ }
+ // And continue main loop on the outside
+ }
+ else
+ {
+ // on same edge?
+ if(Rectangle::onSameEdge(prev_pos,pos))
+ {
+ // Nothing to output if we haven't been elsewhere
+ if(start_index < i-1 || add_start)
+ {
+ std::vector<Coordinate> *coords = new std::vector<Coordinate>();
+ //geom::LineString * line = new geom::LineString();
+ if(add_start)
+ {
+ //line->addPoint(x0,y0);
+ coords->push_back(Coordinate(x0, y0));
+ add_start = false;
+ }
+ //line->addSubLineString(&g, start_index, i-1);
+ coords->insert(coords->end(), cs.begin()+start_index, cs.begin()+i);
+
+ CoordinateSequence *seq = _csf->create(coords);
+ geom::LineString * line = _gf->createLineString(seq);
+ parts.add(line);
+ }
+ start_index = i;
+ }
+ else
+ {
+ // On different edge. Must have gone through the box
+ // then - keep collecting points that generate inside
+ // line segments
+ }
+ }
+ }
+
+ // Was everything in? If so, generate no output but return true in this case only
+ if(start_index == 0 && i >= n)
+ return true;
+
+ // Flush the last line segment if data ended and there is something to flush
+
+ if(!go_outside && // meaning data ended
+ (start_index < i-1 || add_start)) // meaning something has to be generated
+ {
+ std::vector<Coordinate> *coords = new std::vector<Coordinate>();
+ //geom::LineString * line = new geom::LineString();
+ if(add_start)
+ {
+ //line->addPoint(x0,y0);
+ coords->push_back(Coordinate(x0, y0));
+ add_start = false;
+ }
+ //line->addSubLineString(&g, start_index, i-1);
+ coords->insert(coords->end(), cs.begin()+start_index, cs.begin()+i);
+
+ CoordinateSequence *seq = _csf->create(coords);
+ geom::LineString * line = _gf->createLineString(seq);
+ parts.add(line);
+ }
+
+ }
+ }
+
+ return false;
+
+}
+
+/**
+ * \brief Clip polygon, do not close clipped ones
+ */
+
+void
+RectangleIntersection::clip_polygon_to_linestrings(const geom::Polygon * g,
+ RectangleIntersectionBuilder & toParts,
+ const Rectangle & rect)
+{
+ if(g == NULL || g->isEmpty())
+ return;
+
+ // Clip the exterior first to see what's going on
+
+ RectangleIntersectionBuilder parts(*_gf);
+
+ // If everything was in, just clone the original
+
+ if(clip_linestring_parts(g->getExteriorRing(), parts, rect))
+ {
+ toParts.add(dynamic_cast<geom::Polygon *>(g->clone()));
+ return;
+ }
+
+ // Now, if parts is empty, our rectangle may be inside the polygon
+ // If not, holes are outside too
+
+ if(parts.empty())
+ {
+ // We could now check whether the rectangle is inside the outer
+ // ring to avoid checking the holes. However, if holes are much
+ // smaller than the exterior ring just checking the holes
+ // separately could be faster.
+
+ if(g->getNumInteriorRing() == 0)
+ return;
+
+ }
+ else
+ {
+ // The exterior must have been clipped into linestrings.
+ // Move them to the actual parts collector, clearing parts.
+
+ parts.reconnect();
+ parts.release(toParts);
+ }
+
+ // Handle the holes now:
+ // - Clipped ones become linestrings
+ // - Intact ones become new polygons without holes
+
+ for(int i=0, n=g->getNumInteriorRing(); i<n; ++i)
+ {
+ if(clip_linestring_parts(g->getInteriorRingN(i), parts, rect))
+ {
+ // clones
+ LinearRing *hole = dynamic_cast<LinearRing*>(g->getInteriorRingN(i)->clone());
+ // becomes exterior
+ Polygon *poly = _gf->createPolygon(hole, 0);
+ toParts.add(poly);
+ }
+ else if(!parts.empty())
+ {
+ parts.reconnect();
+ parts.release(toParts);
+ }
+ }
+}
+
+/**
+ * \brief Clip polygon, close clipped ones
+ */
+
+void
+RectangleIntersection::clip_polygon_to_polygons(const geom::Polygon * g,
+ RectangleIntersectionBuilder & toParts,
+ const Rectangle & rect)
+{
+ if(g == NULL || g->isEmpty())
+ return;
+
+ // Clip the exterior first to see what's going on
+
+ RectangleIntersectionBuilder parts(*_gf);
+
+ // If everything was in, just clone the original
+
+ const LineString *shell = g->getExteriorRing();
+ if(clip_linestring_parts(shell, parts, rect))
+ {
+ toParts.add(dynamic_cast<geom::Polygon *>(g->clone()));
+ return;
+ }
+
+ // If there were no intersections, the outer ring might be
+ // completely outside.
+
+ using geos::algorithm::CGAlgorithms;
+ if( parts.empty() )
+ {
+ const Coordinate rectCorner(rect.xmin(), rect.ymin());
+ if ( CGAlgorithms::locatePointInRing(rectCorner,
+ *g->getExteriorRing()->getCoordinatesRO())
+ == Location::EXTERIOR )
+ {
+ return;
+ }
+ }
+ else
+ {
+ // TODO: make CCW checking part of clip_linestring_parts ?
+ if ( CGAlgorithms::isCCW(shell->getCoordinatesRO()) ) {
+ parts.reverseLines();
+ }
+ }
+
+ // Must do this to make sure all end points are on the edges
+
+ parts.reconnect();
+
+ // Handle the holes now:
+ // - Clipped ones become part of the exterior
+ // - Intact ones become holes in new polygons formed by exterior parts
+
+
+ for(int i=0, n=g->getNumInteriorRing(); i<n; ++i)
+ {
+ RectangleIntersectionBuilder holeparts(*_gf);
+ const LineString *hole = g->getInteriorRingN(i);
+ if(clip_linestring_parts(hole, holeparts, rect))
+ {
+ // becomes exterior
+ LinearRing *cloned = dynamic_cast<LinearRing*>(hole->clone());
+ Polygon *poly = _gf->createPolygon(cloned, 0);
+ parts.add(poly);
+ }
+ else
+ {
+ if(!holeparts.empty())
+ {
+ // TODO: make CCW checking part of clip_linestring_parts ?
+ if ( ! CGAlgorithms::isCCW(hole->getCoordinatesRO()) ) {
+ holeparts.reverseLines();
+ }
+ holeparts.reconnect();
+ holeparts.release(parts);
+ }
+ else
+ {
+ using geos::algorithm::CGAlgorithms;
+ if( CGAlgorithms::isPointInRing(Coordinate(rect.xmin(), rect.ymin()),
+ g->getInteriorRingN(i)->getCoordinatesRO()) )
+ {
+ // Completely inside the hole
+ return;
+ }
+ }
+ }
+ }
+
+ parts.reconnectPolygons(rect);
+ parts.release(toParts);
+
+}
+
+/**
+ * \brief Clip geometry
+ */
+
+void
+RectangleIntersection::clip_polygon(const geom::Polygon * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons)
+{
+ if(keep_polygons)
+ clip_polygon_to_polygons(g, parts, rect);
+ else
+ clip_polygon_to_linestrings(g, parts, rect);
+}
+
+/**
+ * \brief Clip geometry
+ */
+
+void
+RectangleIntersection::clip_linestring(const geom::LineString * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect)
+{
+ if(g == NULL || g->isEmpty())
+ return;
+
+ // If everything was in, just clone the original
+
+ if(clip_linestring_parts(g, parts, rect))
+ parts.add(dynamic_cast<geom::LineString *>(g->clone()));
+
+}
+
+void
+RectangleIntersection::clip_multipoint(const geom::MultiPoint * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect)
+{
+ if(g == NULL || g->isEmpty())
+ return;
+ for(int i=0, n=g->getNumGeometries(); i<n; ++i)
+ {
+ clip_point(dynamic_cast<const geom::Point *>(g->getGeometryN(i)),
+ parts, rect);
+ }
+}
+
+void
+RectangleIntersection::clip_multilinestring(const geom::MultiLineString * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect)
+{
+ if(g == NULL || g->isEmpty())
+ return;
+
+ for(int i=0, n=g->getNumGeometries(); i<n; ++i)
+ {
+ clip_linestring(dynamic_cast<const geom::LineString *>(g->getGeometryN(i)),
+ parts, rect);
+ }
+}
+
+void
+RectangleIntersection::clip_multipolygon(const geom::MultiPolygon * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons)
+{
+ if(g == NULL || g->isEmpty())
+ return;
+
+ for(int i=0, n=g->getNumGeometries(); i<n; ++i)
+ {
+ clip_polygon(dynamic_cast<const geom::Polygon *>(g->getGeometryN(i)),
+ parts, rect, keep_polygons);
+ }
+}
+
+void
+RectangleIntersection::clip_geometrycollection(
+ const geom::GeometryCollection * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons)
+{
+ if(g == NULL || g->isEmpty())
+ return;
+
+ for(int i=0, n=g->getNumGeometries(); i<n; ++i)
+ {
+ clip_geom(g->getGeometryN(i),
+ parts, rect, keep_polygons);
+ }
+}
+
+void
+RectangleIntersection::clip_geom(const geom::Geometry * g,
+ RectangleIntersectionBuilder & parts,
+ const Rectangle & rect,
+ bool keep_polygons)
+{
+ if ( const Point* p = dynamic_cast<const geom::Point *>(g) )
+ return clip_point(p, parts, rect);
+ else if ( const MultiPoint* p = dynamic_cast<const geom::MultiPoint *>(g) )
+ return clip_multipoint(p, parts, rect);
+ else if ( const LineString* p = dynamic_cast<const geom::LineString *>(g) )
+ return clip_linestring(p, parts, rect);
+ else if ( const MultiLineString* p = dynamic_cast<const geom::MultiLineString *>(g) )
+ return clip_multilinestring(p, parts, rect);
+ else if ( const Polygon* p = dynamic_cast<const geom::Polygon *>(g) )
+ return clip_polygon(p, parts, rect, keep_polygons);
+ else if ( const MultiPolygon* p = dynamic_cast<const geom::MultiPolygon *>(g) )
+ return clip_multipolygon(p, parts, rect, keep_polygons);
+ else if ( const GeometryCollection* p = dynamic_cast<const geom::GeometryCollection *>(g) )
+ return clip_geometrycollection(p, parts, rect, keep_polygons);
+ else {
+ throw util::UnsupportedOperationException("Encountered an unknown geometry component when clipping polygons");
+ }
+}
+
+/* public static */
+std::auto_ptr<geom::Geometry>
+RectangleIntersection::clipBoundary(const geom::Geometry & g, const Rectangle & rect)
+{
+ RectangleIntersection ri(g,rect);
+ return ri.clipBoundary();
+}
+
+std::auto_ptr<geom::Geometry>
+RectangleIntersection::clipBoundary()
+{
+ RectangleIntersectionBuilder parts(*_gf);
+
+ bool keep_polygons = false;
+ clip_geom(&_geom, parts, _rect, keep_polygons);
+
+ return parts.build();
+}
+
+/* public static */
+std::auto_ptr<geom::Geometry>
+RectangleIntersection::clip(const geom::Geometry & g, const Rectangle & rect)
+{
+ RectangleIntersection ri(g,rect);
+ return ri.clip();
+}
+
+std::auto_ptr<geom::Geometry>
+RectangleIntersection::clip()
+{
+ RectangleIntersectionBuilder parts(*_gf);
+
+ bool keep_polygons = true;
+ clip_geom(&_geom, parts, _rect, keep_polygons);
+
+ return parts.build();
+}
+
+RectangleIntersection::RectangleIntersection(const geom::Geometry& geom, const Rectangle& rect)
+ : _geom(geom), _rect(rect),
+ _gf(geom.getFactory())
+{
+ _csf = _gf->getCoordinateSequenceFactory();
+}
+
+} // namespace geos::operation::intersection
+} // namespace geos::operation
+} // namespace geos
Added: trunk/src/operation/intersection/RectangleIntersectionBuilder.cpp
===================================================================
--- trunk/src/operation/intersection/RectangleIntersectionBuilder.cpp (rev 0)
+++ trunk/src/operation/intersection/RectangleIntersectionBuilder.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,518 @@
+/**********************************************************************
+ *
+ * GEOS - Geometry Engine Open Source
+ * http://geos.osgeo.org
+ *
+ * Copyright (C) 2014 Mika Heiskanen <mika.heiskanen at fmi.fi>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU Lesser General Public Licence as published
+ * by the Free Software Foundation.
+ * See the COPYING file for more information.
+ *
+ **********************************************************************/
+
+#include <geos/operation/intersection/Rectangle.h>
+#include <geos/operation/intersection/RectangleIntersectionBuilder.h>
+#include <geos/geom/CoordinateSequenceFactory.h>
+#include <geos/geom/GeometryFactory.h>
+#include <geos/geom/GeometryCollection.h>
+#include <geos/geom/Polygon.h>
+#include <geos/geom/Point.h>
+#include <geos/geom/LineString.h>
+#include <geos/geom/LinearRing.h>
+#include <geos/algorithm/CGAlgorithms.h>
+
+namespace geos {
+namespace operation { // geos::operation
+namespace intersection { // geos::operation::intersection
+
+using namespace geos::geom;
+
+RectangleIntersectionBuilder::~RectangleIntersectionBuilder()
+{
+ for (std::list<geom::Polygon *>::iterator i=polygons.begin(), e=polygons.end(); i!=e; ++i)
+ delete *i;
+ for (std::list<geom::LineString *>::iterator i=lines.begin(), e=lines.end(); i!=e; ++i)
+ delete *i;
+ for (std::list<geom::Point *>::iterator i=points.begin(), e=points.end(); i!=e; ++i)
+ delete *i;
+}
+
+void
+RectangleIntersectionBuilder::reconnect()
+{
+ // Nothing to reconnect if there aren't at least two lines
+ if(lines.size() < 2)
+ return;
+
+ geom::LineString * line1 = lines.front();
+ const geom::CoordinateSequence &cs1 = *line1->getCoordinatesRO();
+
+ geom::LineString * line2 = lines.back();
+ const geom::CoordinateSequence &cs2 = *line2->getCoordinatesRO();
+
+ const int n1 = cs1.size();
+ const int n2 = cs2.size();
+
+ // Safety check against bad input to prevent segfaults
+ if(n1==0 || n2==0)
+ return;
+
+ if (cs1[0] != cs2[n2-1]) return;
+
+ // Merge the two linestrings
+
+ CoordinateSequence *ncs = CoordinateSequence::removeRepeatedPoints(&cs2);
+ ncs->add(&cs1, false, true);
+
+ delete line1;
+ delete line2;
+
+ LineString * nline = _gf.createLineString(ncs);
+ lines.pop_front();
+ lines.pop_back();
+
+ lines.push_front(nline);
+}
+
+
+void RectangleIntersectionBuilder::release(RectangleIntersectionBuilder & theParts)
+{
+ for (std::list<geom::Polygon *>::iterator i=polygons.begin(), e=polygons.end(); i!=e; ++i)
+ theParts.add(*i);
+
+ for (std::list<geom::LineString *>::iterator i=lines.begin(), e=lines.end(); i!=e; ++i)
+ theParts.add(*i);
+
+ for (std::list<geom::Point *>::iterator i=points.begin(), e=points.end(); i!=e; ++i)
+ theParts.add(*i);
+
+ clear();
+}
+
+/**
+ * \brief Clear the parts having transferred ownership elsewhere
+ */
+
+void RectangleIntersectionBuilder::clear()
+{
+ polygons.clear();
+ lines.clear();
+ points.clear();
+}
+
+/**
+ * \brief Test if there are no parts at all
+ */
+
+bool RectangleIntersectionBuilder::empty() const
+{
+ return polygons.empty() && lines.empty() && points.empty();
+}
+
+/**
+ * \brief Add intermediate Polygon
+ */
+
+void RectangleIntersectionBuilder::add(geom::Polygon * thePolygon)
+{
+ polygons.push_back(thePolygon);
+}
+
+/**
+ * \brief Add intermediate LineString
+ */
+
+void RectangleIntersectionBuilder::add(geom::LineString * theLine)
+{
+ lines.push_back(theLine);
+}
+
+/**
+ * \brief Add intermediate Point
+ */
+
+void RectangleIntersectionBuilder::add(geom::Point * thePoint)
+{
+ points.push_back(thePoint);
+}
+
+std::auto_ptr<geom::Geometry>
+RectangleIntersectionBuilder::build()
+{
+ // Total number of objects
+
+ std::size_t n = polygons.size() + lines.size() + points.size();
+
+ if(n == 0)
+ return std::auto_ptr<Geometry>(_gf.createGeometryCollection());
+
+ std::vector<Geometry *> *geoms = new std::vector<Geometry *>;
+ geoms->reserve(n);
+
+ for (std::list<geom::Polygon *>::iterator i=polygons.begin(), e=polygons.end(); i!=e; ++i)
+ geoms->push_back(*i);
+ polygons.clear();
+
+ for (std::list<geom::LineString *>::iterator i=lines.begin(), e=lines.end(); i!=e; ++i)
+ geoms->push_back(*i);
+ lines.clear();
+
+ for (std::list<geom::Point *>::iterator i=points.begin(), e=points.end(); i!=e; ++i)
+ geoms->push_back(*i);
+ points.clear();
+
+ return std::auto_ptr<Geometry>(
+ (*geoms)[0]->getFactory()->buildGeometry(geoms)
+ );
+}
+
+/**
+ * Distance of 1st point of linestring to last point of linearring
+ * along rectangle edges
+ */
+
+double distance(const Rectangle & rect,
+ double x1, double y1,
+ double x2, double y2)
+{
+ double dist = 0;
+
+ Rectangle::Position pos = rect.position(x1,y1);
+ Rectangle::Position endpos = rect.position(x2,y2);
+
+ while(true)
+ {
+ // Close up when we have the same edge and the
+ // points are in the correct clockwise order
+ if((pos & endpos) != 0 &&
+ (
+ (x1==rect.xmin() && y2>=y1) ||
+ (y1==rect.ymax() && x2>=x1) ||
+ (x1==rect.xmax() && y2<=y1) ||
+ (y1==rect.ymin() && x2<=x1))
+ )
+ {
+ dist += fabs(x2-x1) + fabs(y2-y1);
+ break;
+ }
+
+ pos = Rectangle::nextEdge(pos);
+ if(pos & Rectangle::Left)
+ {
+ dist += x1-rect.xmin();
+ x1 = rect.xmin();
+ }
+ else if(pos & Rectangle::Top)
+ {
+ dist += rect.ymax()-y1;
+ y1 = rect.ymax();
+ }
+ else if(pos & Rectangle::Right)
+ {
+ dist += rect.xmax()-x1;
+ x1 = rect.xmax();
+ }
+ else
+ {
+ dist += y1-rect.ymin();
+ y1 = rect.ymin();
+ }
+ }
+ return dist;
+}
+
+double distance(const Rectangle & rect,
+ const std::vector<Coordinate> &ring,
+ const geom::LineString * line)
+{
+ double nr = ring.size();
+ const Coordinate &c1 = ring[nr-1];
+
+ const CoordinateSequence * linecs = line->getCoordinatesRO();
+ const Coordinate &c2 = linecs->getAt(0);
+
+ return distance(rect,c1.x,c1.y,c2.x,c2.y);
+}
+
+double distance(const Rectangle & rect,
+ const std::vector<Coordinate> &ring)
+{
+ double nr = ring.size();
+ const Coordinate& c1 = ring[nr-1]; // TODO: ring.back() ?
+ const Coordinate& c2 = ring[0]; // TODO: ring.front() ?
+ return distance(rect, c1.x, c1.y, c2.x, c2.y);
+}
+
+/**
+ * \brief Reverse given segment in a coordinate vector
+ */
+void reverse_points(std::vector<Coordinate> &v, int start, int end)
+{
+ geom::Coordinate p1;
+ geom::Coordinate p2;
+ while(start < end)
+ {
+ p1 = v[start];
+ p2 = v[end];
+ v[start] = p2;
+ v[end] = p1;
+ ++start;
+ --end;
+ }
+}
+
+/**
+ * \brief Normalize a ring into lexicographic order
+ */
+void
+normalize_ring(std::vector<Coordinate> &ring)
+{
+ if(ring.empty())
+ return;
+
+ // Find the "smallest" coordinate
+
+ int best_pos = 0;
+ int n = ring.size();
+ for(int pos = 0; pos<n; ++pos)
+ {
+ // TODO: use CoordinateLessThan ?
+ if(ring[pos].x < ring[best_pos].x)
+ best_pos = pos;
+ else if(ring[pos].x == ring[best_pos].x &&
+ ring[pos].y < ring[best_pos].y )
+ best_pos = pos;
+ }
+
+ // Quick exit if the ring is already normalized
+ if(best_pos == 0)
+ return;
+
+ // Flip hands -algorithm to the part without the
+ // duplicate last coordinate at n-1:
+
+ reverse_points(ring,0,best_pos-1);
+ reverse_points(ring,best_pos,n-2);
+ reverse_points(ring,0,n-2);
+
+ // And make sure the ring is valid by duplicating the first coordinate
+ // at the end:
+
+ geom::Coordinate c;
+ c = ring[0];
+ ring[n-1] = c;
+}
+
+void
+RectangleIntersectionBuilder::close_boundary(
+ const Rectangle & rect,
+ std::vector<Coordinate> * ring,
+ double x1, double y1,
+ double x2, double y2)
+{
+ Rectangle::Position endpos = rect.position(x2,y2);
+ Rectangle::Position pos = rect.position(x1,y1);
+
+ while(true)
+ {
+ // Close up when we have the same edge and the
+ // points are in the correct clockwise order
+ if((pos & endpos) != 0 &&
+ (
+ (x1==rect.xmin() && y2>=y1) ||
+ (y1==rect.ymax() && x2>=x1) ||
+ (x1==rect.xmax() && y2<=y1) ||
+ (y1==rect.ymin() && x2<=x1))
+ )
+ {
+ if(x1!=x2 || y1!=y2) // the polygon may have started at a corner
+ ring->push_back(Coordinate(x2,y2));
+ break;
+ }
+
+ pos = Rectangle::nextEdge(pos);
+ if(pos & Rectangle::Left)
+ x1 = rect.xmin();
+ else if(pos & Rectangle::Top)
+ y1 = rect.ymax();
+ else if(pos & Rectangle::Right)
+ x1 = rect.xmax();
+ else
+ y1 = rect.ymin();
+
+ ring->push_back(Coordinate(x1,y1));
+
+ }
+}
+
+void
+RectangleIntersectionBuilder::close_ring(const Rectangle & rect,
+ std::vector<Coordinate> * ring)
+{
+ double nr = ring->size();
+ Coordinate &c2 = (*ring)[0];
+ Coordinate &c1 = (*ring)[nr-1];
+ double x2 = c2.x;
+ double y2 = c2.y;
+ double x1 = c1.x;
+ double y1 = c1.y;
+
+ close_boundary(rect, ring, x1, y1, x2, y2);
+
+}
+
+
+void
+RectangleIntersectionBuilder::reconnectPolygons(const Rectangle & rect)
+{
+ // Build the exterior rings first
+
+ typedef std::vector< geom::Geometry *> LinearRingVect;
+ typedef std::pair< geom::LinearRing *, LinearRingVect * > ShellAndHoles;
+ typedef std::list< ShellAndHoles > ShellAndHolesList;
+
+ ShellAndHolesList exterior;
+
+ const CoordinateSequenceFactory &_csf = *_gf.getCoordinateSequenceFactory();
+
+ // If there are no lines, the rectangle must have been
+ // inside the exterior ring.
+
+ if(lines.empty())
+ {
+ geom::LinearRing * ring = rect.toLinearRing(_gf);
+ exterior.push_back(make_pair(ring, new LinearRingVect()));
+ }
+ else
+ {
+ // Reconnect all lines into one or more linearrings
+ // using box boundaries if necessary
+
+ std::vector<Coordinate> *ring = NULL;
+
+ while(!lines.empty() || ring != NULL)
+ {
+ if(ring == NULL)
+ {
+ ring = new std::vector<Coordinate>();
+ LineString *line = lines.front();
+ lines.pop_front();
+ line->getCoordinatesRO()->toVector(*ring);
+ delete line;
+ }
+
+ // Distance to own endpoint
+ double own_distance = distance(rect, *ring);
+
+ // Find line to connect to
+ // TODO: should we use LineMerge op ?
+ double best_distance = -1;
+ std::list<LineString*>::iterator best_pos = lines.begin();
+ for(std::list<LineString*>::iterator iter=lines.begin(); iter!=lines.end(); ++iter)
+ {
+ double d = distance(rect, *ring, *iter);
+ if(best_distance < 0 || d<best_distance)
+ {
+ best_distance = d;
+ best_pos = iter;
+ }
+ }
+
+ // If own end point is closest, close the ring and continue
+ if(best_distance < 0 || own_distance < best_distance)
+ {
+ close_ring(rect,ring);
+ normalize_ring(*ring);
+ geom::CoordinateSequence *shell_cs = _csf.create(ring);
+ geom::LinearRing *shell = _gf.createLinearRing(shell_cs);
+ exterior.push_back(make_pair(shell, new LinearRingVect()));
+ ring = NULL;
+ }
+ else
+ {
+ LineString * line = *best_pos;
+ int nr = ring->size();
+ const CoordinateSequence& cs = *line->getCoordinatesRO();
+ close_boundary(rect, ring,
+ (*ring)[nr-1].x,
+ (*ring)[nr-1].y,
+ cs[0].x,
+ cs[0].y);
+ // above function adds the 1st point
+ for (size_t i=1; i<cs.size(); ++i)
+ ring->push_back(cs[i]);
+ //ring->addSubLineString(line,1);
+ delete line;
+ lines.erase(best_pos);
+ }
+ }
+ }
+
+ // Attach holes to polygons
+
+ for (std::list<geom::Polygon *>::iterator i=polygons.begin(), e=polygons.end(); i!=e; ++i)
+ {
+ geom::Polygon *poly = *i;
+ const geom::LineString *hole = poly->getExteriorRing();
+
+ if(exterior.size() == 1)
+ {
+ exterior.front().second->push_back( hole->clone() );
+ }
+ else
+ {
+ using geos::algorithm::CGAlgorithms;
+ geom::Coordinate c;
+ hole->getCoordinatesRO()->getAt(0, c);
+ for (ShellAndHolesList::iterator i=exterior.begin(), e=exterior.end(); i!=e; ++i)
+ {
+ ShellAndHoles &p = *i;
+ const CoordinateSequence *shell_cs = p.first->getCoordinatesRO();
+ if( CGAlgorithms::isPointInRing(c, shell_cs) )
+ {
+ // add hole to shell
+ p.second->push_back(hole->clone());
+ break;
+ }
+ }
+ }
+
+ delete poly;
+ }
+
+ // Build the result polygons
+
+ std::list<geom::Polygon *> new_polygons;
+ for (ShellAndHolesList::iterator i=exterior.begin(), e=exterior.end(); i!=e; ++i)
+ {
+ ShellAndHoles &p = *i;
+ geom::Polygon * poly = _gf.createPolygon(p.first, p.second);
+ new_polygons.push_back(poly);
+ }
+
+ clear();
+ polygons = new_polygons;
+}
+
+void
+RectangleIntersectionBuilder::reverseLines()
+{
+ std::list<geom::LineString *> new_lines;
+ for (std::list<geom::LineString *>::reverse_iterator i=lines.rbegin(), e=lines.rend(); i!=e; ++i)
+ {
+ LineString *ol = *i;
+ new_lines.push_back(dynamic_cast<LineString*>(ol->reverse()));
+ delete ol;
+ }
+ lines = new_lines;
+#if GEOS_DEBUG
+ std::cout << "After lines reverse, parts are " << *this << std::endl;
+#endif
+}
+
+
+} // namespace geos::operation::intersection
+} // namespace geos::operation
+} // namespace geos
Modified: trunk/tests/unit/Makefile.am
===================================================================
--- trunk/tests/unit/Makefile.am 2014-09-24 07:54:58 UTC (rev 4020)
+++ trunk/tests/unit/Makefile.am 2014-09-25 10:46:20 UTC (rev 4021)
@@ -89,6 +89,7 @@
operation/buffer/BufferOpTest.cpp \
operation/buffer/BufferParametersTest.cpp \
operation/distance/DistanceOpTest.cpp \
+ operation/intersection/RectangleIntersectionTest.cpp \
operation/IsSimpleOpTest.cpp \
operation/linemerge/LineMergerTest.cpp \
operation/linemerge/LineSequencerTest.cpp \
@@ -114,6 +115,7 @@
triangulate/DelaunayTest.cpp \
triangulate/VoronoiTest.cpp \
util/UniqueCoordinateArrayFilterTest.cpp \
+ capi/GEOSClipByRectTest.cpp \
capi/GEOSCoordSeqTest.cpp \
capi/GEOSDelaunayTriangulationTest.cpp \
capi/GEOSVoronoiDiagramTest.cpp \
Added: trunk/tests/unit/capi/GEOSClipByRectTest.cpp
===================================================================
--- trunk/tests/unit/capi/GEOSClipByRectTest.cpp (rev 0)
+++ trunk/tests/unit/capi/GEOSClipByRectTest.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,191 @@
+//
+// Test Suite for C-API GEOSClipByRect
+
+#include <tut.hpp>
+// geos
+#include <geos_c.h>
+// std
+#include <cstdarg>
+#include <cstdio>
+#include <cstdlib>
+#include <memory>
+
+namespace tut
+{
+ //
+ // Test Group
+ //
+
+ // Common data used in test cases.
+ struct test_capigeosclipbyrect_data
+ {
+ GEOSGeometry* geom1_;
+ GEOSGeometry* geom2_;
+ GEOSGeometry* geom3_;
+ GEOSWKTWriter* w_;
+
+ static void notice(const char *fmt, ...)
+ {
+ std::fprintf( stdout, "NOTICE: ");
+
+ va_list ap;
+ va_start(ap, fmt);
+ std::vfprintf(stdout, fmt, ap);
+ va_end(ap);
+
+ std::fprintf(stdout, "\n");
+ }
+
+ void isEqual(GEOSGeom g, const char *exp_wkt)
+ {
+ geom3_ = GEOSGeomFromWKT(exp_wkt);
+ bool eq = GEOSEquals(geom3_, g);
+ if ( ! eq ) {
+ std::printf("EXP: %s\n", exp_wkt);
+ char *obt_wkt = GEOSWKTWriter_write(w_, g);
+ std::printf("OBT: %s\n", obt_wkt);
+ free(obt_wkt);
+ }
+ ensure(eq);
+ }
+
+ test_capigeosclipbyrect_data()
+ : geom1_(0), geom2_(0), geom3_(0), w_(0)
+ {
+ initGEOS(notice, notice);
+ w_ = GEOSWKTWriter_create();
+ GEOSWKTWriter_setTrim(w_, 1);
+ GEOSWKTWriter_setRoundingPrecision(w_, 8);
+ }
+
+ ~test_capigeosclipbyrect_data()
+ {
+ GEOSGeom_destroy(geom1_);
+ GEOSGeom_destroy(geom2_);
+ GEOSGeom_destroy(geom3_);
+ GEOSWKTWriter_destroy(w_);
+ geom1_ = 0;
+ geom2_ = 0;
+ geom3_ = 0;
+ finishGEOS();
+ }
+
+ };
+
+ typedef test_group<test_capigeosclipbyrect_data> group;
+ typedef group::object object;
+
+ group test_capigeosclipbyrect_group("capi::GEOSClipByRect");
+
+ //
+ // Test Cases
+ //
+
+ /// Point outside
+ template<> template<> void object::test<1>()
+ {
+ geom1_ = GEOSGeomFromWKT("POINT(0 0)");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "POINT EMPTY");
+ }
+
+ /// Point inside
+ template<> template<> void object::test<2>()
+ {
+ geom1_ = GEOSGeomFromWKT("POINT(15 15)");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "POINT(15 15)");
+ }
+
+ /// Point on boundary
+ template<> template<> void object::test<3>()
+ {
+ geom1_ = GEOSGeomFromWKT("POINT(15 10)");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "POINT EMPTY");
+ }
+
+ /// Line outside
+ template<> template<> void object::test<4>()
+ {
+ geom1_ = GEOSGeomFromWKT("LINESTRING(0 0, -5 5)");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "LINESTRING EMPTY");
+ }
+
+ /// Line inside
+ template<> template<> void object::test<5>()
+ {
+ geom1_ = GEOSGeomFromWKT("LINESTRING(15 15, 16 15)");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "LINESTRING(15 15, 16 15)");
+ }
+
+ /// Line on boundary
+ template<> template<> void object::test<6>()
+ {
+ geom1_ = GEOSGeomFromWKT("LINESTRING(10 15, 10 10, 15 10)");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "LINESTRING EMPTY");
+ }
+
+ /// Line splitting rectangle
+ template<> template<> void object::test<7>()
+ {
+ geom1_ = GEOSGeomFromWKT("LINESTRING(10 5, 25 20)");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "LINESTRING (15 10, 20 15)");
+ }
+
+ /// Polygon shell (CCW) fully on rectangle boundary
+ template<> template<> void object::test<8>()
+ {
+ geom1_ = GEOSGeomFromWKT("POLYGON((10 10, 20 10, 20 20, 10 20, 10 10))");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "POLYGON((10 10, 20 10, 20 20, 10 20, 10 10))");
+ }
+
+ /// Polygon shell (CW) fully on rectangle boundary
+ template<> template<> void object::test<9>()
+ {
+ geom1_ = GEOSGeomFromWKT("POLYGON((10 10, 10 20, 20 20, 20 10, 10 10))");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "POLYGON((10 10, 20 10, 20 20, 10 20, 10 10))");
+ }
+
+ /// Polygon hole (CCW) fully on rectangle boundary
+ template<> template<> void object::test<10>()
+ {
+ geom1_ = GEOSGeomFromWKT("POLYGON((0 0, 0 30, 30 30, 30 0, 0 0),(10 10, 20 10, 20 20, 10 20, 10 10))");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "POLYGON EMPTY");
+ }
+
+ /// Polygon hole (CW) fully on rectangle boundary
+ template<> template<> void object::test<11>()
+ {
+ geom1_ = GEOSGeomFromWKT("POLYGON((0 0, 0 30, 30 30, 30 0, 0 0),(10 10, 10 20, 20 20, 20 10, 10 10))");
+ geom2_ = GEOSClipByRect(geom1_, 10, 10, 20, 20);
+ isEqual(geom2_, "POLYGON EMPTY");
+ }
+
+ /// Polygon fully within rectangle
+ template<> template<> void object::test<12>()
+ {
+ const char *wkt = "POLYGON((1 1, 1 30, 30 30, 30 1, 1 1),(10 10, 20 10, 20 20, 10 20, 10 10))";
+ geom1_ = GEOSGeomFromWKT(wkt);
+ geom2_ = GEOSClipByRect(geom1_, 0, 0, 40, 40);
+ isEqual(geom2_, wkt);
+ }
+
+ /// Polygon overlapping rectangle
+ template<> template<> void object::test<13>()
+ {
+ const char *wkt = "POLYGON((0 0, 0 30, 30 30, 30 0, 0 0),(10 10, 20 10, 20 20, 10 20, 10 10))";
+ geom1_ = GEOSGeomFromWKT(wkt);
+ geom2_ = GEOSClipByRect(geom1_, 5, 5, 15, 15);
+ isEqual(geom2_, "POLYGON ((5 5, 5 15, 10 15, 10 10, 15 10, 15 5, 5 5))");
+ }
+
+} // namespace tut
+
Added: trunk/tests/unit/operation/intersection/RectangleIntersectionTest.cpp
===================================================================
--- trunk/tests/unit/operation/intersection/RectangleIntersectionTest.cpp (rev 0)
+++ trunk/tests/unit/operation/intersection/RectangleIntersectionTest.cpp 2014-09-25 10:46:20 UTC (rev 4021)
@@ -0,0 +1,1506 @@
+//
+// Test Suite for geos::operation::intersection::RectangleIntersection class.
+
+// tut
+#include <tut.hpp>
+// geos
+#include <geos/operation/intersection/Rectangle.h>
+#include <geos/operation/intersection/RectangleIntersection.h>
+#include <geos/geom/GeometryFactory.h>
+#include <geos/geom/Geometry.h>
+#include <geos/geom/Polygon.h>
+#include <geos/geom/Point.h>
+#include <geos/io/WKTReader.h>
+#include <geos/io/WKTWriter.h>
+// std
+#include <memory>
+#include <string>
+#include <vector>
+#include <iostream>
+
+namespace tut
+{
+ //
+ // Test Group
+ //
+
+ // Common data used by tests
+ struct test_rectangleintersectiontest_data
+ {
+ geos::geom::GeometryFactory gf;
+ geos::io::WKTReader wktreader;
+ geos::io::WKTWriter wktwriter;
+
+ typedef geos::geom::Geometry::AutoPtr GeomPtr;
+ typedef geos::geom::Geometry Geom;
+ typedef geos::operation::intersection::Rectangle Rectangle;
+ typedef geos::operation::intersection::RectangleIntersection RectangleIntersection;
+
+ test_rectangleintersectiontest_data()
+ : gf(),
+ wktreader(&gf)
+ {
+ wktwriter.setTrim(true);
+ }
+
+ GeomPtr readWKT(const std::string& inputWKT)
+ {
+ return GeomPtr(wktreader.read(inputWKT));
+ }
+
+ static GeomPtr normalize(const Geom& g)
+ {
+ GeomPtr g2 ( g.clone() );
+ g2->normalize();
+ return g2;
+ }
+
+ bool isEqual(const Geom& a, const Geom& b, double tolerance=0)
+ {
+ using std::cout;
+ using std::endl;
+ bool eq;
+ // TODO: use HausdorffDistance ?
+ GeomPtr a2 = normalize(a);
+ GeomPtr b2 = normalize(b);
+ eq = a2->equalsExact(b2.get(), tolerance);
+ if ( ! eq ) {
+ cout << "OBTAINED: " << wktwriter.write(b2.get()) << endl;
+ }
+ return eq;
+ }
+
+ void doLineClipTest(const char* inputWKT, const std::string& expectedWKT, const Rectangle& rect, double tolerance=0)
+ {
+ GeomPtr g = readWKT(inputWKT);
+ ensure(g.get());
+ GeomPtr obt = RectangleIntersection::clipBoundary(*g,rect);
+ ensure(obt.get());
+ bool ok = isEqual(*readWKT(expectedWKT), *obt, tolerance);
+ ensure(ok);
+ }
+
+ void doClipTest(const char* inputWKT, const std::string& expectedWKT, const Rectangle& rect, double tolerance=0)
+ {
+ GeomPtr g = readWKT(inputWKT);
+ ensure(g.get());
+ GeomPtr obt = RectangleIntersection::clip(*g,rect);
+ ensure(obt.get());
+ bool ok = isEqual(*readWKT(expectedWKT), *obt, tolerance);
+ ensure(ok);
+// Compare with GEOSIntersection output
+#if 0
+ GeomPtr g2 ( rect.toPolygon(*g->getFactory()) );
+ obt.reset(g->intersection(g2.get()));
+ ensure(obt.get());
+ ok = isEqual(*readWKT(expectedWKT), *obt, tolerance);
+ ensure(ok);
+#endif
+ }
+
+ };
+
+ typedef test_group<test_rectangleintersectiontest_data, 256> group;
+ typedef group::object object;
+
+ group test_rectangleintersectiontest_group("geos::operation::intersection::RectangleIntersection");
+
+ // inside
+ template<> template<> void object::test<1>()
+ {
+ doLineClipTest(
+ "LINESTRING (1 1,1 9,9 9,9 1)",
+ "LINESTRING (1 1,1 9,9 9,9 1)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // outside
+ template<> template<> void object::test<2>()
+ {
+ doLineClipTest(
+ "LINESTRING (-1 -9,-1 11,9 11)",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go in from left
+ template<> template<> void object::test<3>()
+ {
+ doLineClipTest(
+ "LINESTRING (-1 5,5 5,9 9)",
+ "LINESTRING (0 5,5 5,9 9)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go out from right
+ template<> template<> void object::test<4>()
+ {
+ doLineClipTest(
+ "LINESTRING (5 5,8 5,12 5)",
+ "LINESTRING (5 5,8 5,10 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go in and out
+ template<> template<> void object::test<5>()
+ {
+ doLineClipTest(
+ "LINESTRING (5 -1,5 5,1 2,-3 2,1 6)",
+ "MULTILINESTRING ((5 0,5 5,1 2,0 2),(0 5,1 6))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go along left edge
+ template<> template<> void object::test<6>()
+ {
+ doLineClipTest(
+ "LINESTRING (0 3,0 5,0 7)",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go out from left edge
+ template<> template<> void object::test<7>()
+ {
+ doLineClipTest(
+ "LINESTRING (0 3,0 5,-1 7)",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go in from left edge
+ template<> template<> void object::test<8>()
+ {
+ doLineClipTest(
+ "LINESTRING (0 3,0 5,2 7)",
+ "LINESTRING (0 5,2 7)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // triangle corner at bottom left corner
+ template<> template<> void object::test<9>()
+ {
+ doLineClipTest(
+ "LINESTRING (2 1,0 0,1 2)",
+ "LINESTRING (2 1,0 0,1 2)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go from in to edge and back in
+ template<> template<> void object::test<10>()
+ {
+ doLineClipTest(
+ "LINESTRING (3 3,0 3,0 5,2 7)",
+ "MULTILINESTRING ((3 3,0 3),(0 5,2 7))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go from in to edge and then straight out
+ template<> template<> void object::test<11>()
+ {
+ doLineClipTest(
+ "LINESTRING (5 5,10 5,20 5)",
+ "LINESTRING (5 5,10 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // triangle corner at left edge
+ template<> template<> void object::test<12>()
+ {
+ doLineClipTest(
+ "LINESTRING (3 3,0 6,3 9)",
+ "LINESTRING (3 3,0 6,3 9)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon completely inside
+ template<> template<> void object::test<13>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 5,5 6,6 6,6 5,5 5))",
+ "POLYGON ((5 5,5 6,6 6,6 5,5 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon completely outside
+ template<> template<> void object::test<14>()
+ {
+ doLineClipTest(
+ "POLYGON ((15 15,15 16,16 16,16 15,15 15))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon surrounds the rectangle
+ template<> template<> void object::test<15>()
+ {
+ doLineClipTest(
+ "POLYGON ((-1 -1,-1 11,11 11,11 -1,-1 -1))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon cuts the rectangle
+ template<> template<> void object::test<16>()
+ {
+ doLineClipTest(
+ "POLYGON ((-1 -1,-1 5,5 5,5 -1,-1 -1))",
+ "LINESTRING (0 5,5 5,5 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon with hole cuts the rectangle
+ template<> template<> void object::test<17>()
+ {
+ doLineClipTest(
+ "POLYGON ((-2 -2,-2 5,5 5,5 -2,-2 -2), (3 3,4 4,4 2,3 3))",
+ "GEOMETRYCOLLECTION (POLYGON ((3 3,4 4,4 2,3 3)),LINESTRING (0 5,5 5,5 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // rectangle cuts both the polygon and the hole
+ template<> template<> void object::test<18>()
+ {
+ doLineClipTest(
+ "POLYGON ((-2 -2,-2 5,5 5,5 -2,-2 -2), (-1 -1,3 1,3 3,-1 -1))",
+ "MULTILINESTRING ((0 5,5 5,5 0),(1 0,3 1,3 3,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle at two corners and one edge
+ template<> template<> void object::test<19>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 0,10 0,5 10,0 0))",
+ "LINESTRING (10 0,5 10,0 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Same triangle with another starting point
+ template<> template<> void object::test<20>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 10,0 0,10 0,5 10))",
+ "LINESTRING (10 0,5 10,0 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle intersection at corner and edge
+ template<> template<> void object::test<21>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 -5,5 5,5 -5,-5 -5))",
+ "LINESTRING (0 0,5 5,5 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle intersection at adjacent edges
+ template<> template<> void object::test<22>()
+ {
+ doLineClipTest(
+ "POLYGON ((-1 5,5 11,5 5,-1 5))",
+ "MULTILINESTRING ((0 6,4 10),(5 10,5 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // One triangle intersection and one inside edge
+ template<> template<> void object::test<23>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 5,5 10,5 5,-5 5))",
+ "LINESTRING (0.0 7.5,5 10,5 5,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle intersection at center and end of the same edge
+ template<> template<> void object::test<24>()
+ {
+ doLineClipTest(
+ "POLYGON ((-10 5,10 10,10 5,-10 5))",
+ "MULTILINESTRING ((0.0 7.5,10 10),(10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Two different edges clips
+ template<> template<> void object::test<25>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 5,15 15,15 5,-5 5))",
+ "MULTILINESTRING ((0.0 7.5,5 10),(10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Inside triangle with all corners at edges
+ template<> template<> void object::test<26>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 5,5 10,10 5,0 5))",
+ "POLYGON ((0 5,5 10,10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Inside triangle whose base is one of the edges
+ template<> template<> void object::test<27>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 0,5 5,10 0,0 0))",
+ "LINESTRING (0 0,5 5,10 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle touching two corners on the outside
+ template<> template<> void object::test<28>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 5,5 15,15 5,-5 5))",
+ "LINESTRING (10 5,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle with a diagonal and sharing two edges
+ template<> template<> void object::test<29>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 0,10 10,10 0,0 0))",
+ "LINESTRING (0 0,10 10)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle exits edge at a corner
+ template<> template<> void object::test<30>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 0,5 10,5 0,-5 0))",
+ "LINESTRING (0 5,5 10,5 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle enters edge at a corner
+ template<> template<> void object::test<31>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 10,5 10,0 0,-5 10))",
+ "LINESTRING (5 10,0 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle enters and exits the same edge
+ template<> template<> void object::test<32>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 0,5 10,15 0,-5 0))",
+ "LINESTRING (0 5,5 10,10 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle enters at a corner and exits at another
+ template<> template<> void object::test<33>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 -5,15 15,15 -5,-5 -5))",
+ "LINESTRING (0 0,10 10)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // From outside to nearest edge etc
+ template<> template<> void object::test<34>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 -5,0 5,5 0,-5 -5))",
+ "LINESTRING (0 5,5 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // From outside to opposite edge etc
+ template<> template<> void object::test<35>()
+ {
+ doLineClipTest(
+ "POLYGON ((-10 5,10 5,0 -5,-10 5))",
+ "LINESTRING (0 5,10 5,5 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Drew all combinations I could imagine on paper, and added the following.
+ // All triangles fully inside
+ template<> template<> void object::test<36>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 0,0 10,10 10,0 0))",
+ "LINESTRING (10 10,0 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<37>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 5,0 10,10 10,0 5))",
+ "LINESTRING (10 10,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<38>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,10 10,5 0,0 10))",
+ "LINESTRING (10 10,5 0,0 10)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<39>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,10 10,5 5,0 10))",
+ "LINESTRING (10 10,5 5,0 10)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<40>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,5 10,0 5,0 10))",
+ "LINESTRING (5 10,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<41>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,10 5,0 5,0 10))",
+ "LINESTRING (0 10,10 5,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<42>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,10 0,0 5,0 10))",
+ "LINESTRING (0 10,10 0,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<43>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,5 0,0 5,0 10))",
+ "LINESTRING (0 10,5 0,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<44>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,5 5,0 5,0 10))",
+ "LINESTRING (0 10,5 5,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<45>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,7 7,3 3,0 10))",
+ "POLYGON ((0 10,7 7,3 3,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<46>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,5 5,5 0,0 10))",
+ "POLYGON ((0 10,5 5,5 0,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<47>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 10,10 5,5 0,0 10))",
+ "POLYGON ((0 10,10 5,5 0,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<48>()
+ {
+ doLineClipTest(
+ "POLYGON ((2 5,5 7,7 5,2 5))",
+ "POLYGON ((2 5,5 7,7 5,2 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<49>()
+ {
+ doLineClipTest(
+ "POLYGON ((2 5,5 10,7 5,2 5))",
+ "POLYGON ((2 5,5 10,7 5,2 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<50>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 5,5 10,5 5,0 5))",
+ "POLYGON ((0 5,5 10,5 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<51>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 5,5 10,10 5,0 5))",
+ "POLYGON ((0 5,5 10,10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<52>()
+ {
+ doLineClipTest(
+ "POLYGON ((0 5,5 7,10 5,0 5))",
+ "POLYGON ((0 5,5 7,10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // No points inside, one intersection
+ template<> template<> void object::test<53>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 10,0 15,0 10,-5 10))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<54>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 10,0 5,-5 0,-5 10))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // No points inside, two intersections
+ template<> template<> void object::test<55>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 5,0 10,0 0,-5 5))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<56>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 5,0 10,0 5,-5 5))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<57>()
+ {
+ doLineClipTest(
+ "POLYGON ((-5 5,0 7,0 3,-5 5))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // One point inside
+ template<> template<> void object::test<58>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 5,-5 0,-5 10,5 5))",
+ "LINESTRING (0.0 7.5,5 5,0.0 2.5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<59>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 0,-5 0,-5 10,5 0))",
+ "LINESTRING (0 5,5 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<60>()
+ {
+ doLineClipTest(
+ "POLYGON((10 0,-10 0,-10 10,10 0))",
+ "LINESTRING (0 5,10 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<61>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 0,-5 5,-5 10,5 0))",
+ "LINESTRING (0 5,5 0,0.0 2.5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<62>()
+ {
+ doLineClipTest(
+ "POLYGON ((10 5,-10 0,-10 10,10 5))",
+ "LINESTRING (0.0 7.5,10 5,0.0 2.5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<63>()
+ {
+ doLineClipTest(
+ "POLYGON ((10 10,-10 0,-10 5,10 10))",
+ "LINESTRING (0.0 7.5,10 10,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<64>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 5,-5 -5,-5 15,5 5))",
+ "LINESTRING (0 10,5 5,0 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<65>()
+ {
+ doLineClipTest(
+ "POLYGON ((10 5,-10 -5,-10 15,10 5))",
+ "LINESTRING (0 10,10 5,0 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<66>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 0,-5 0,-5 20,5 0))",
+ "LINESTRING (0 10,5 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<67>()
+ {
+ doLineClipTest(
+ "POLYGON ((10 0,-10 0,-10 20,10 0))",
+ "LINESTRING (0 10,10 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<68>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 5,-10 5,0 15,5 5))",
+ "LINESTRING (2.5 10.0,5 5,0 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<69>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 5,-5 -5,0 15,5 5))",
+ "LINESTRING (2.5 10.0,5 5,0 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<70>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 5,-15 -20,-15 30,5 5))",
+ "LINESTRING (1 10,5 5,1 0)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Two points inside
+ template<> template<> void object::test<71>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 7,5 3,-5 5,5 7))",
+ "LINESTRING (0 6,5 7,5 3,0 4)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<72>()
+ {
+ doLineClipTest(
+ "POLYGON ((5 7,5 3,-5 13,5 7))",
+ "LINESTRING (0 10,5 7,5 3,0 8)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<73>()
+ {
+ doLineClipTest(
+ "POLYGON ((6 6,4 4,-4 14,6 6))",
+ "LINESTRING (1.0 10.0,6 6,4 4,0 9)",
+ Rectangle(0,0,10,10), 1e-12
+ );
+ }
+
+ // Polygon with hole which surrounds the rectangle
+ template<> template<> void object::test<74>()
+ {
+ doLineClipTest(
+ "POLYGON ((-2 -2,-2 12,12 12,12 -2,-2 -2),(-1 -1,11 -1,11 11,-1 11,-1 -1))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Polygon surrounding the rect, but with a hole inside the rect
+ template<> template<> void object::test<75>()
+ {
+ doLineClipTest(
+ "POLYGON ((-2 -2,-2 12,12 12,12 -2,-2 -2),(1 1,9 1,9 9,1 9,1 1))",
+ "POLYGON ((1 1,9 1,9 9,1 9,1 1))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+
+ // inside
+ template<> template<> void object::test<76>()
+ {
+ doClipTest(
+ "LINESTRING (1 1,1 9,9 9,9 1)",
+ "LINESTRING (1 1,1 9,9 9,9 1)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // outside
+ template<> template<> void object::test<77>()
+ {
+ doClipTest(
+ "LINESTRING (-1 -9,-1 11,9 11)",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go in from left
+ template<> template<> void object::test<78>()
+ {
+ doClipTest(
+ "LINESTRING (-1 5,5 5,9 9)",
+ "LINESTRING (0 5,5 5,9 9)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go out from right
+ template<> template<> void object::test<79>()
+ {
+ doClipTest(
+ "LINESTRING (5 5,8 5,12 5)",
+ "LINESTRING (5 5,8 5,10 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go in and out
+ template<> template<> void object::test<80>()
+ {
+ doClipTest(
+ "LINESTRING (5 -1,5 5,1 2,-3 2,1 6)",
+ "MULTILINESTRING ((5 0,5 5,1 2,0 2),(0 5,1 6))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go along left edge
+ template<> template<> void object::test<81>()
+ {
+ doClipTest(
+ "LINESTRING (0 3,0 5,0 7)",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go out from left edge
+ template<> template<> void object::test<82>()
+ {
+ doClipTest(
+ "LINESTRING (0 3,0 5,-1 7)",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go in from left edge
+ template<> template<> void object::test<83>()
+ {
+ doClipTest(
+ "LINESTRING (0 3,0 5,2 7)",
+ "LINESTRING (0 5,2 7)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // triangle corner at bottom left corner
+ template<> template<> void object::test<84>()
+ {
+ doClipTest(
+ "LINESTRING (2 1,0 0,1 2)",
+ "LINESTRING (2 1,0 0,1 2)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go from in to edge and back in
+ template<> template<> void object::test<85>()
+ {
+ doClipTest(
+ "LINESTRING (3 3,0 3,0 5,2 7)",
+ "MULTILINESTRING ((3 3,0 3),(0 5,2 7))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // go from in to edge and then straight out
+ template<> template<> void object::test<86>()
+ {
+ doClipTest(
+ "LINESTRING (5 5,10 5,20 5)",
+ "LINESTRING (5 5,10 5)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // triangle corner at left edge
+ template<> template<> void object::test<87>()
+ {
+ doClipTest(
+ "LINESTRING (3 3,0 6,3 9)",
+ "LINESTRING (3 3,0 6,3 9)",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon completely inside
+ template<> template<> void object::test<88>()
+ {
+ doClipTest(
+ "POLYGON ((5 5,5 6,6 6,6 5,5 5))",
+ "POLYGON ((5 5,5 6,6 6,6 5,5 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon completely outside
+ template<> template<> void object::test<89>()
+ {
+ doClipTest(
+ "POLYGON ((15 15,15 16,16 16,16 15,15 15))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon surrounds the rectangle
+ template<> template<> void object::test<90>()
+ {
+ doClipTest(
+ "POLYGON ((-1 -1,-1 11,11 11,11 -1,-1 -1))",
+ "POLYGON ((0 0,0 10,10 10,10 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon cuts the rectangle
+ template<> template<> void object::test<91>()
+ {
+ doClipTest(
+ "POLYGON ((-1 -1,-1 5,5 5,5 -1,-1 -1))",
+ "POLYGON ((0 0,0 5,5 5,5 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // polygon with hole cuts the rectangle
+ template<> template<> void object::test<92>()
+ {
+ doClipTest(
+ "POLYGON ((-2 -2,-2 5,5 5,5 -2,-2 -2), (3 3,4 4,4 2,3 3))",
+ "POLYGON ((0 0,0 5,5 5,5 0,0 0),(3 3,4 4,4 2,3 3))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // rectangle cuts both the polygon and the hole
+ template<> template<> void object::test<93>()
+ {
+ doClipTest(
+ "POLYGON ((-2 -2,-2 5,5 5,5 -2,-2 -2), (-1 -1,3 1,3 3,-1 -1))",
+ "POLYGON ((0 0,0 5,5 5,5 0,1 0,3 1,3 3,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle at two corners and one edge
+ template<> template<> void object::test<94>()
+ {
+ doClipTest(
+ "POLYGON ((0 0,10 0,5 10,0 0))",
+ "POLYGON ((0 0,10 0,5 10,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Same triangle with another starting point
+ template<> template<> void object::test<95>()
+ {
+ doClipTest(
+ "POLYGON ((5 10,0 0,10 0,5 10))",
+ "POLYGON ((0 0,10 0,5 10,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Triangle intersection at corner and edge
+ template<> template<> void object::test<96>()
+ {
+ doClipTest(
+ "POLYGON ((-5 -5,5 5,5 -5,-5 -5))",
+ "POLYGON ((0 0,5 5,5 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // All triangles fully inside
+ template<> template<> void object::test<97>()
+ {
+ doClipTest(
+ "POLYGON ((0 0,0 10,10 10,0 0))",
+ "POLYGON ((0 0,0 10,10 10,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<98>()
+ {
+ doClipTest(
+ "POLYGON ((0 5,0 10,10 10,0 5))",
+ "POLYGON ((0 5,0 10,10 10,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<99>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,10 10,5 0,0 10))",
+ "POLYGON ((0 10,10 10,5 0,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<100>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,10 10,5 5,0 10))",
+ "POLYGON ((0 10,10 10,5 5,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<101>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,5 10,0 5,0 10))",
+ "POLYGON ((0 5,0 10,5 10,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<102>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,10 5,0 5,0 10))",
+ "POLYGON ((0 5,0 10,10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<103>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,10 0,0 5,0 10))",
+ "POLYGON ((0 5,0 10,10 0,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<104>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,5 0,0 5,0 10))",
+ "POLYGON ((0 5,0 10,5 0,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<105>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,5 5,0 5,0 10))",
+ "POLYGON ((0 5,0 10,5 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<106>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,7 7,3 3,0 10))",
+ "POLYGON ((0 10,7 7,3 3,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<107>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,5 5,5 0,0 10))",
+ "POLYGON ((0 10,5 5,5 0,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<108>()
+ {
+ doClipTest(
+ "POLYGON ((0 10,10 5,5 0,0 10))",
+ "POLYGON ((0 10,10 5,5 0,0 10))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<109>()
+ {
+ doClipTest(
+ "POLYGON ((2 5,5 7,7 5,2 5))",
+ "POLYGON ((2 5,5 7,7 5,2 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<110>()
+ {
+ doClipTest(
+ "POLYGON ((2 5,5 10,7 5,2 5))",
+ "POLYGON ((2 5,5 10,7 5,2 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<111>()
+ {
+ doClipTest(
+ "POLYGON ((0 5,5 10,5 5,0 5))",
+ "POLYGON ((0 5,5 10,5 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<112>()
+ {
+ doClipTest(
+ "POLYGON ((0 5,5 10,10 5,0 5))",
+ "POLYGON ((0 5,5 10,10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<113>()
+ {
+ doClipTest(
+ "POLYGON ((0 5,5 7,10 5,0 5))",
+ "POLYGON ((0 5,5 7,10 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // No points inside, one intersection
+ template<> template<> void object::test<114>()
+ {
+ doClipTest(
+ "POLYGON ((-5 10,0 15,0 10,-5 10))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<115>()
+ {
+ doClipTest(
+ "POLYGON ((-5 10,0 5,-5 0,-5 10))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // No points inside, two intersections
+ template<> template<> void object::test<116>()
+ {
+ doClipTest(
+ "POLYGON ((-5 5,0 10,0 0,-5 5))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<117>()
+ {
+ doClipTest(
+ "POLYGON ((-5 5,0 10,0 5,-5 5))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<118>()
+ {
+ doClipTest(
+ "POLYGON ((-5 5,0 7,0 3,-5 5))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // One point inside
+ template<> template<> void object::test<119>()
+ {
+ doClipTest(
+ "POLYGON ((5 5,-5 0,-5 10,5 5))",
+ "POLYGON ((0.0 2.5,0.0 7.5,5 5,0.0 2.5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<120>()
+ {
+ doClipTest(
+ "POLYGON ((5 0,-5 0,-5 10,5 0))",
+ "POLYGON ((0 0,0 5,5 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<121>()
+ {
+ doClipTest(
+ "POLYGON ((10 0,-10 0,-10 10,10 0))",
+ "POLYGON ((0 0,0 5,10 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<122>()
+ {
+ doClipTest(
+ "POLYGON ((5 0,-5 5,-5 10,5 0))",
+ "POLYGON ((0.0 2.5,0 5,5 0,0.0 2.5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<123>()
+ {
+ doClipTest(
+ "POLYGON ((10 5,-10 0,-10 10,10 5))",
+ "POLYGON ((0.0 2.5,0.0 7.5,10 5,0.0 2.5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<124>()
+ {
+ doClipTest(
+ "POLYGON ((10 10,-10 0,-10 5,10 10))",
+ "POLYGON ((0 5,0.0 7.5,10 10,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<125>()
+ {
+ doClipTest(
+ "POLYGON ((5 5,-5 -5,-5 15,5 5))",
+ "POLYGON ((0 0,0 10,5 5,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<126>()
+ {
+ doClipTest(
+ "POLYGON ((10 5,-10 -5,-10 15,10 5))",
+ "POLYGON ((0 0,0 10,10 5,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<127>()
+ {
+ doClipTest(
+ "POLYGON ((5 0,-5 0,-5 20,5 0))",
+ "POLYGON ((0 0,0 10,5 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<128>()
+ {
+ doClipTest(
+ "POLYGON ((10 0,-10 0,-10 20,10 0))",
+ "POLYGON ((0 0,0 10,10 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<129>()
+ {
+ doClipTest(
+ "POLYGON ((5 5,-10 5,0 15,5 5))",
+ "POLYGON ((0 5,0 10,2.5 10.0,5 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<130>()
+ {
+ doClipTest(
+ "POLYGON ((5 5,-5 -5,0 15,5 5))",
+ "POLYGON ((0 0,0 10,2.5 10.0,5 5,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<131>()
+ {
+ doClipTest(
+ "POLYGON ((5 5,-15 -20,-15 30,5 5))",
+ "POLYGON ((0 0,0 10,1 10,5 5,1 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Two points inside
+ template<> template<> void object::test<132>()
+ {
+ doClipTest(
+ "POLYGON ((5 7,5 3,-5 5,5 7))",
+ "POLYGON ((0 4,0 6,5 7,5 3,0 4))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<133>()
+ {
+ doClipTest(
+ "POLYGON ((5 7,5 3,-5 13,5 7))",
+ "POLYGON ((0 8,0 10,5 7,5 3,0 8))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ template<> template<> void object::test<134>()
+ {
+ doClipTest(
+ "POLYGON ((6 6,4 4,-4 14,6 6))",
+ "POLYGON ((0 9,0 10,1.0 10.0,6 6,4 4,0 9))",
+ Rectangle(0,0,10,10), 1e-12
+ );
+ }
+
+ // Polygon with hole which surrounds the rectangle
+ template<> template<> void object::test<135>()
+ {
+ doClipTest(
+ "POLYGON ((-2 -2,-2 12,12 12,12 -2,-2 -2),(-1 -1,11 -1,11 11,-1 11,-1 -1))",
+ "GEOMETRYCOLLECTION EMPTY",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Polygon surrounding the rect, but with a hole inside the rect
+ template<> template<> void object::test<136>()
+ {
+ doClipTest(
+ "POLYGON ((-2 -2,-2 12,12 12,12 -2,-2 -2),(1 1,9 1,9 9,1 9,1 1))",
+ "POLYGON ((0 0,0 10,10 10,10 0,0 0),(1 1,9 1,9 9,1 9,1 1))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Polygon with hole cut at the right corner
+ template<> template<> void object::test<137>()
+ {
+ doClipTest(
+ "POLYGON ((5 5,15 5,15 -5,5 -5,5 5),(8 1,8 -1,9 -1,9 1,8 1))",
+ "POLYGON ((5 0,5 5,10 5,10 0,9 0,9 1,8 1,8 0,5 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Polygon going around a corner
+ template<> template<> void object::test<138>()
+ {
+ doClipTest(
+ "POLYGON ((-6 5,5 5,5 -6,-6 5))",
+ "POLYGON ((0 0,0 5,5 5,5 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Hole in a corner
+ template<> template<> void object::test<139>()
+ {
+ doClipTest(
+ "POLYGON ((-15 -15,-15 15,15 15,15 -15,-15 -15),(-5 5,-5 -5,5 -5,5 5,-5 5))",
+ "POLYGON ((0 5,0 10,10 10,10 0,5 0,5 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Hole going around a corner
+ template<> template<> void object::test<140>()
+ {
+ doClipTest(
+ "POLYGON ((-15 -15,-15 15,15 15,15 -15,-15 -15),(-6 5,5 -6,5 5,-6 5))",
+ "POLYGON ((0 5,0 10,10 10,10 0,5 0,5 5,0 5))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Surround the rectangle, hole outside rectangle
+ template<> template<> void object::test<141>()
+ {
+ doClipTest(
+ "POLYGON ((-15 -15,-15 15,15 15,15 -15,-15 -15),(-5 5,-6 5,-6 6,-5 6,-5 5))",
+ "POLYGON ((0 0,0 10,10 10,10 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+ // Surround the rectangle, hole outside rectangle but shares edge
+ template<> template<> void object::test<142>()
+ {
+ doClipTest(
+ "POLYGON ((-15 -15,-15 15,15 15,15 -15,-15 -15),(0 5,-1 5,-1 6,0 6,0 5))",
+ "POLYGON ((0 0,0 10,10 10,10 0,0 0))",
+ Rectangle(0,0,10,10)
+ );
+ }
+
+//-- NEW TESTS HERE
+
+ // Simple overlap, counter-clockwise shell
+ // Found in TestFunctionAA.xml case1 test1
+ template<> template<> void object::test<204>()
+ {
+ Rectangle r(10,10,100,100);
+ const char *inp =
+ "POLYGON ((50 50,200 50,200 200,50 200,50 50))"; // CCW
+ //"POLYGON((50 50,50 200,200 200,200 50,50 50))"; // CW
+ const char *exp =
+ "POLYGON ((50 50, 50 100, 100 100, 100 50, 50 50))";
+ doClipTest(inp, exp, r);
+ }
+
+ // Clockwise shell, clockwise hole (both clipped)
+ template<> template<> void object::test<205>()
+ {
+ Rectangle r(0,0,10,10);
+ const char *inp =
+ "POLYGON ("
+ "(-10 2,-10 8,8 8,8 2,-10 2)," // CW
+ "(-5 6,5 6,5 4,-5 4,-5 6)" // CW
+ ")";
+ const char *exp =
+ "POLYGON((0 8,8 8, 8 2, 0 2, 0 4, 5 4, 5 6, 0 6, 0 8))";
+ doClipTest(inp, exp, r);
+ }
+
+ // Counterclockwise shell, clockwise hole (both clipped)
+ template<> template<> void object::test<206>()
+ {
+ Rectangle r(0,0,10,10);
+ const char *inp =
+ "POLYGON ("
+ "(-10 2,8 2,8 8,-10 8,-10 2)," // CCW
+ "(-5 6,5 6,5 4,-5 4,-5 6)" // CW
+ ")";
+ const char *exp =
+ "POLYGON((0 8,8 8, 8 2, 0 2, 0 4, 5 4, 5 6, 0 6, 0 8))";
+ doClipTest(inp, exp, r);
+ }
+}
More information about the geos-commits
mailing list