[geos-commits] r3676 - in branches/3.3: . include/geos/algorithm src/algorithm tests/unit tests/unit/capi
svn_geos at osgeo.org
svn_geos at osgeo.org
Fri Jun 22 08:02:06 PDT 2012
Author: strk
Date: 2012-06-22 08:02:06 -0700 (Fri, 22 Jun 2012)
New Revision: 3676
Added:
branches/3.3/tests/unit/capi/GEOSGetCentroidTest.cpp
Modified:
branches/3.3/NEWS
branches/3.3/include/geos/algorithm/CentroidArea.h
branches/3.3/src/algorithm/CentroidArea.cpp
branches/3.3/tests/unit/Makefile.am
Log:
Port robustness fix to CentroidArea (#559)
Modified: branches/3.3/NEWS
===================================================================
--- branches/3.3/NEWS 2012-06-22 11:47:08 UTC (rev 3675)
+++ branches/3.3/NEWS 2012-06-22 15:02:06 UTC (rev 3676)
@@ -3,6 +3,7 @@
- Bug fixes / improvements
- Correctly increment CAPI lib version from 3.3.3 (#558)
+ - Port robustness fix to CentroidArea (#559)
Changes in 3.3.4
2012-05-31
Modified: branches/3.3/include/geos/algorithm/CentroidArea.h
===================================================================
--- branches/3.3/include/geos/algorithm/CentroidArea.h 2012-06-22 11:47:08 UTC (rev 3675)
+++ branches/3.3/include/geos/algorithm/CentroidArea.h 2012-06-22 15:02:06 UTC (rev 3676)
@@ -12,6 +12,10 @@
* by the Free Software Foundation.
* See the COPYING file for more information.
*
+ **********************************************************************
+ *
+ * Last port: algorithm/CentroidArea.java r612
+ *
**********************************************************************/
#ifndef GEOS_ALGORITHM_CENTROIDAREA_H
@@ -54,7 +58,8 @@
CentroidArea()
:
basePt(0.0, 0.0),
- areasum2(0)
+ areasum2(0.0),
+ totalLength(0.0)
{}
~CentroidArea() {}
@@ -75,9 +80,10 @@
*/
void add(const geom::CoordinateSequence *ring);
+ // TODO: deprecate
geom::Coordinate* getCentroid() const;
- /// Return false if a centroid couldn't be computed
+ /// Return false if a centroid couldn't be computed ( empty polygon )
bool getCentroid(geom::Coordinate& ret) const;
private:
@@ -94,6 +100,10 @@
/// partial centroid sum
geom::Coordinate cg3;
+ // data for linear centroid computation, if needed
+ geom::Coordinate centSum;
+ double totalLength;
+
void setBasePoint(const geom::Coordinate &newbasePt);
void add(const geom::Polygon *poly);
@@ -111,6 +121,17 @@
static double area2(const geom::Coordinate &p1, const geom::Coordinate &p2,
const geom::Coordinate &p3);
+ /**
+ * Adds the linear segments defined by an array of coordinates
+ * to the linear centroid accumulators.
+ *
+ * This is done in case the polygon(s) have zero-area,
+ * in which case the linear centroid is computed instead.
+ *
+ * @param pts an array of {@link Coordinate}s
+ */
+ void addLinearSegments(const geom::CoordinateSequence& pts);
+
};
} // namespace geos::algorithm
@@ -118,11 +139,3 @@
#endif // GEOS_ALGORITHM_CENTROIDAREA_H
-
-/**********************************************************************
- * $Log$
- * Revision 1.1 2006/03/09 16:46:48 strk
- * geos::geom namespace definition, first pass at headers split
- *
- **********************************************************************/
-
Modified: branches/3.3/src/algorithm/CentroidArea.cpp
===================================================================
--- branches/3.3/src/algorithm/CentroidArea.cpp 2012-06-22 11:47:08 UTC (rev 3675)
+++ branches/3.3/src/algorithm/CentroidArea.cpp 2012-06-22 15:02:06 UTC (rev 3676)
@@ -14,6 +14,8 @@
*
**********************************************************************
*
+ * Last port: algorithm/CentroidArea.java r612
+ *
**********************************************************************/
#include <geos/algorithm/CentroidArea.h>
@@ -57,20 +59,25 @@
addShell(ring);
}
+/* TODO: deprecate this */
Coordinate*
CentroidArea::getCentroid() const
{
Coordinate *cent = new Coordinate();
- cent->x = cg3.x/3.0/areasum2;
- cent->y = cg3.y/3.0/areasum2;
+ getCentroid(*cent); // or return NULL on failure !
return cent;
}
bool
CentroidArea::getCentroid(Coordinate& ret) const
{
- if ( areasum2 == 0.0 ) return false;
- ret=Coordinate(cg3.x/3.0/areasum2, cg3.y/3.0/areasum2);
+ if ( areasum2 ) {
+ ret = Coordinate(cg3.x/3.0/areasum2, cg3.y/3.0/areasum2);
+ } else if ( totalLength ) {
+ ret = Coordinate(centSum.x/totalLength, centSum.y/totalLength);
+ } else {
+ return false;
+ }
return true;
}
@@ -94,11 +101,12 @@
CentroidArea::addShell(const CoordinateSequence *pts)
{
bool isPositiveArea=!CGAlgorithms::isCCW(pts);
- std::size_t const n=pts->getSize()-1;
+ std::size_t const n=pts->getSize()-1;
for(std::size_t i=0; i<n; ++i)
{
addTriangle(basePt, pts->getAt(i), pts->getAt(i+1), isPositiveArea);
}
+ addLinearSegments(*pts);
}
void
@@ -110,6 +118,7 @@
{
addTriangle(basePt, pts->getAt(i), pts->getAt(i+1), isPositiveArea);
}
+ addLinearSegments(*pts);
}
void
@@ -147,24 +156,20 @@
return (p2.x-p1.x)*(p3.y-p1.y)-(p3.x-p1.x)*(p2.y-p1.y);
}
+void
+CentroidArea::addLinearSegments(const geom::CoordinateSequence& pts)
+{
+ std::size_t const n = pts.size()-1;
+ for (std::size_t i = 0; i < n; ++i) {
+ double segmentLen = pts[i].distance(pts[i + 1]);
+ totalLength += segmentLen;
+
+ double midx = (pts[i].x + pts[i + 1].x) / 2;
+ centSum.x += segmentLen * midx;
+ double midy = (pts[i].y + pts[i + 1].y) / 2;
+ centSum.y += segmentLen * midy;
+ }
+}
+
} // namespace geos.algorithm
} //namespace geos
-
-/**********************************************************************
- * $Log$
- * Revision 1.24 2006/06/08 17:58:57 strk
- * Polygon::getNumInteriorRing() return size_t, Polygon::interiorRingN() takes size_t.
- *
- * Revision 1.23 2006/04/07 09:54:28 strk
- * Geometry::getNumGeometries() changed to return 'unsigned int'
- * rather then 'int'
- *
- * Revision 1.22 2006/03/21 11:12:23 strk
- * Cleanups: headers inclusion and Log section
- *
- * Revision 1.21 2006/03/09 16:46:45 strk
- * geos::geom namespace definition, first pass at headers split
- *
- **********************************************************************/
-
-
Modified: branches/3.3/tests/unit/Makefile.am
===================================================================
--- branches/3.3/tests/unit/Makefile.am 2012-06-22 11:47:08 UTC (rev 3675)
+++ branches/3.3/tests/unit/Makefile.am 2012-06-22 15:02:06 UTC (rev 3676)
@@ -102,6 +102,7 @@
capi/GEOSCoordSeqTest.cpp \
capi/GEOSGeomFromWKBTest.cpp \
capi/GEOSGeomToWKTTest.cpp \
+ capi/GEOSGetCentroidTest.cpp \
capi/GEOSContainsTest.cpp \
capi/GEOSIntersectsTest.cpp \
capi/GEOSWithinTest.cpp \
Added: branches/3.3/tests/unit/capi/GEOSGetCentroidTest.cpp
===================================================================
--- branches/3.3/tests/unit/capi/GEOSGetCentroidTest.cpp (rev 0)
+++ branches/3.3/tests/unit/capi/GEOSGetCentroidTest.cpp 2012-06-22 15:02:06 UTC (rev 3676)
@@ -0,0 +1,155 @@
+// $Id$
+//
+// Test Suite for C-API GEOSGetCentroid
+
+#include <tut.hpp>
+// geos
+#include <geos_c.h>
+// std
+#include <cstdarg>
+#include <cstdio>
+#include <cstdlib>
+#include <cstring>
+
+namespace tut
+{
+ //
+ // Test Group
+ //
+
+ // Common data used in test cases.
+ struct test_capicentroid_data
+ {
+ GEOSGeometry* geom1_;
+ GEOSGeometry* geom2_;
+ GEOSWKTWriter* wktw_;
+ char* wkt_;
+ double area_;
+
+ 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");
+ }
+
+ test_capicentroid_data()
+ : geom1_(0), geom2_(0), wkt_(0)
+ {
+ initGEOS(notice, notice);
+ wktw_ = GEOSWKTWriter_create();
+ GEOSWKTWriter_setTrim(wktw_, 1);
+ GEOSWKTWriter_setRoundingPrecision(wktw_, 8);
+ }
+
+ ~test_capicentroid_data()
+ {
+ GEOSGeom_destroy(geom1_);
+ GEOSGeom_destroy(geom2_);
+ GEOSWKTWriter_destroy(wktw_);
+ GEOSFree(wkt_);
+ geom1_ = 0;
+ geom2_ = 0;
+ wkt_ = 0;
+ finishGEOS();
+ }
+
+ };
+
+ typedef test_group<test_capicentroid_data> group;
+ typedef group::object object;
+
+ group test_capicentroid_group("capi::GEOSGetCentroid");
+
+ //
+ // Test Cases
+ //
+
+ // Single point
+ template<>
+ template<>
+ void object::test<1>()
+ {
+ geom1_ = GEOSGeomFromWKT("POINT(10 0)");
+
+ ensure( 0 != geom1_ );
+
+ geom2_ = GEOSGetCentroid(geom1_);
+
+ ensure( 0 != geom2_ );
+
+ wkt_ = GEOSWKTWriter_write(wktw_, geom2_);
+
+ ensure_equals(std::string(wkt_), std::string( "POINT (10 0)"));
+
+ }
+
+ // line
+ template<>
+ template<>
+ void object::test<2>()
+ {
+ geom1_ = GEOSGeomFromWKT("LINESTRING(0 0, 10 0)");
+
+ ensure( 0 != geom1_ );
+
+ geom2_ = GEOSGetCentroid(geom1_);
+
+ ensure( 0 != geom2_ );
+
+ wkt_ = GEOSWKTWriter_write(wktw_, geom2_);
+
+ ensure_equals(std::string(wkt_), std::string( "POINT (5 0)"));
+
+ }
+
+ // polygon
+ template<>
+ template<>
+ void object::test<3>()
+ {
+ geom1_ = GEOSGeomFromWKT("POLYGON((0 0, 10 0, 10 10, 0 10, 0 0))");
+
+ ensure( 0 != geom1_ );
+
+ geom2_ = GEOSGetCentroid(geom1_);
+
+ ensure( 0 != geom2_ );
+
+ wkt_ = GEOSWKTWriter_write(wktw_, geom2_);
+
+ ensure_equals(std::string(wkt_), std::string( "POINT (5 5)"));
+
+ }
+
+ // Tiny triangle, see http://trac.osgeo.org/geos/ticket/559
+ template<>
+ template<>
+ void object::test<4>()
+ {
+ geom1_ = GEOSGeomFromWKT(
+"POLYGON(( \
+56.528666666700 25.2101666667, \
+56.529000000000 25.2105000000, \
+56.528833333300 25.2103333333, \
+56.528666666700 25.2101666667))");
+
+ ensure( 0 != geom1_ );
+
+ geom2_ = GEOSGetCentroid(geom1_);
+
+ ensure( 0 != geom2_ );
+
+ wkt_ = GEOSWKTWriter_write(wktw_, geom2_);
+
+ ensure_equals(std::string(wkt_), std::string( "POINT (56.528833 25.210333)" ) );
+
+ }
+
+} // namespace tut
+
More information about the geos-commits
mailing list