Main Page | Namespace List | Class Hierarchy | Class List | File List | Namespace Members | Class Members | Related Pages

BinaryOp.h

00001 /**********************************************************************
00002  * $Id: BinaryOp.h 3713 2012-09-10 08:05:28Z strk $
00003  *
00004  * GEOS - Geometry Engine Open Source
00005  * http://geos.refractions.net
00006  *
00007  * Copyright (C) 2006 Refractions Research Inc.
00008  *
00009  * This is free software; you can redistribute and/or modify it under
00010  * the terms of the GNU Lesser General Public Licence as published
00011  * by the Free Software Foundation. 
00012  * See the COPYING file for more information.
00013  *
00014  **********************************************************************
00015  *
00016  * Last port: ORIGINAL WORK
00017  *
00018  **********************************************************************
00019  *
00020  * This file provides a single templated function, taking two
00021  * const Geometry pointers, applying a binary operator to them
00022  * and returning a result Geometry in an auto_ptr<>.
00023  *
00024  * The binary operator is expected to take two const Geometry pointers
00025  * and return a newly allocated Geometry pointer, possibly throwing
00026  * a TopologyException to signal it couldn't succeed due to robustness
00027  * issues.
00028  *
00029  * This function will catch TopologyExceptions and try again with
00030  * slightly modified versions of the input. The following heuristic
00031  * is used:
00032  *
00033  *      - Try with original input.
00034  *      - Try removing common bits from input coordinate values
00035  *      - Try snaping input geometries to each other
00036  *      - Try snaping input coordinates to a increasing grid (size from 1/25 to 1)
00037  *      - Try simplifiying input with increasing tolerance (from 0.01 to 0.04)
00038  *
00039  * If none of the step succeeds the original exception is thrown.
00040  *
00041  * Note that you can skip Grid snapping, Geometry snapping and Simplify policies
00042  * by a compile-time define when building geos.
00043  * See USE_TP_SIMPLIFY_POLICY, USE_PRECISION_REDUCTION_POLICY and
00044  * USE_SNAPPING_POLICY macros below.
00045  *
00046  *
00047  **********************************************************************/
00048 
00049 #ifndef GEOS_GEOM_BINARYOP_H
00050 #define GEOS_GEOM_BINARYOP_H
00051 
00052 #include <geos/geom/Geometry.h>
00053 #include <geos/geom/PrecisionModel.h>
00054 #include <geos/precision/CommonBitsRemover.h>
00055 #include <geos/precision/SimpleGeometryPrecisionReducer.h>
00056 
00057 #include <geos/operation/overlay/snap/GeometrySnapper.h>
00058 
00059 #include <geos/simplify/TopologyPreservingSimplifier.h>
00060 #include <geos/operation/valid/IsValidOp.h>
00061 #include <geos/operation/valid/TopologyValidationError.h>
00062 #include <geos/util/TopologyException.h>
00063 #include <geos/util.h>
00064 
00065 #include <memory> // for auto_ptr
00066 
00067 //#define GEOS_DEBUG_BINARYOP 1
00068 
00069 #ifdef GEOS_DEBUG_BINARYOP
00070 # include <iostream>
00071 # include <iomanip>
00072 # include <sstream>
00073 #endif
00074 
00075 
00076 /*
00077  * Always try original input first
00078  */
00079 #ifndef USE_ORIGINAL_INPUT
00080 # define USE_ORIGINAL_INPUT 1
00081 #endif
00082 
00083 /*
00084  * Define this to use PrecisionReduction policy
00085  * in an attempt at by-passing binary operation
00086  * robustness problems (handles TopologyExceptions)
00087  */
00088 #ifndef USE_PRECISION_REDUCTION_POLICY
00089 //# define USE_PRECISION_REDUCTION_POLICY 1
00090 #endif
00091 
00092 /*
00093  * Define this to use TopologyPreserving simplification policy
00094  * in an attempt at by-passing binary operation
00095  * robustness problems (handles TopologyExceptions)
00096  */
00097 #ifndef USE_TP_SIMPLIFY_POLICY 
00098 //# define USE_TP_SIMPLIFY_POLICY 1
00099 #endif
00100 
00101 /*
00102  * Use common bits removal policy.
00103  * If enabled, this would be tried /before/
00104  * Geometry snapping.
00105  */
00106 #ifndef USE_COMMONBITS_POLICY
00107 # define USE_COMMONBITS_POLICY 1
00108 #endif
00109 
00110 /*
00111  * Check validity of operation performed
00112  * by common bits removal policy.
00113  *
00114  * This matches what EnhancedPrecisionOp does in JTS
00115  * and fixes 5 tests of invalid outputs in our testsuite
00116  * (stmlf-cases-20061020-invalid-output.xml)
00117  * and breaks 1 test (robustness-invalid-output.xml) so much
00118  * to prevent a result.
00119  *
00120  */
00121 #define GEOS_CHECK_COMMONBITS_VALIDITY 1
00122 
00123 /*
00124  * Use snapping policy
00125  */
00126 #ifndef USE_SNAPPING_POLICY
00127 # define USE_SNAPPING_POLICY 1
00128 #endif
00129 
00130 namespace geos {
00131 namespace geom { // geos::geom
00132 
00133 inline bool
00134 check_valid(const Geometry& g, const std::string& label)
00135 {
00136         operation::valid::IsValidOp ivo(&g);
00137         if ( ! ivo.isValid() )
00138         {
00139 #ifdef GEOS_DEBUG_BINARYOP
00140                 using operation::valid::TopologyValidationError;
00141                 TopologyValidationError* err = ivo.getValidationError();
00142                 std::cerr << label << " is INVALID: "
00143                         << err->toString()
00144                         << " (" << std::setprecision(20)
00145                         << err->getCoordinate() << ")"
00146                         << std::endl;
00147 #else
00148         (void)label;
00149 #endif
00150                 return false;
00151         } 
00152         return true;
00153 }
00154 
00155 /* A single component may become a multi component */
00156 inline std::auto_ptr<Geometry>
00157 fix_self_intersections(std::auto_ptr<Geometry> g, const std::string& label)
00158 {
00159 #ifdef GEOS_DEBUG_BINARYOP
00160         std::cerr << label << " fix_self_intersection (UnaryUnion)" << std::endl;
00161 #endif
00162   // TODO: check for the sole self-intersections case ?
00163   if ( ! check_valid(*g, label) ) return g->Union();
00164   else return std::auto_ptr<Geometry>(g->clone());
00165 }
00166 
00167 
00173 template <class BinOp>
00174 std::auto_ptr<Geometry>
00175 SnapOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00176 {
00177         typedef std::auto_ptr<Geometry> GeomPtr;
00178 
00179 #define CBR_BEFORE_SNAPPING 1
00180 
00181         //using geos::precision::GeometrySnapper;
00182         using geos::operation::overlay::snap::GeometrySnapper;
00183 
00184         // Snap tolerance must be computed on the original
00185         // (not commonbits-removed) geoms
00186         double snapTolerance = GeometrySnapper::computeOverlaySnapTolerance(*g0, *g1);
00187 #if GEOS_DEBUG_BINARYOP
00188         std::cerr<< std::setprecision(20) << "Computed snap tolerance: "<<snapTolerance<<std::endl;
00189 #endif
00190 
00191 
00192 #if CBR_BEFORE_SNAPPING
00193         // Compute common bits
00194         geos::precision::CommonBitsRemover cbr;
00195         cbr.add(g0); cbr.add(g1);
00196 #if GEOS_DEBUG_BINARYOP
00197         std::cerr<<"Computed common bits: "<<cbr.getCommonCoordinate()<<std::endl;
00198 #endif
00199 
00200         // Now remove common bits
00201         GeomPtr rG0( cbr.removeCommonBits(g0->clone()) );
00202         GeomPtr rG1( cbr.removeCommonBits(g1->clone()) );
00203 
00204 #if GEOS_DEBUG_BINARYOP
00205         check_valid(*rG0, "CBR: removed-bits geom 0");
00206         check_valid(*rG1, "CBR: removed-bits geom 1");
00207 #endif
00208 
00209         const Geometry& operand0 = *rG0;
00210         const Geometry& operand1 = *rG1;
00211 #else // don't CBR before snapping
00212         const Geometry& operand0 = *g0;
00213         const Geometry& operand1 = *g1;
00214 #endif
00215 
00216 
00217         GeometrySnapper snapper0( operand0 );
00218         GeomPtr snapG0( snapper0.snapTo(operand1, snapTolerance) );
00219         snapG0 = fix_self_intersections(snapG0, "SNAP: snapped geom 0");
00220 
00221         // NOTE: second geom is snapped on the snapped first one
00222         GeometrySnapper snapper1( operand1 );
00223         GeomPtr snapG1( snapper1.snapTo(*snapG0, snapTolerance) );
00224         snapG1 = fix_self_intersections(snapG1, "SNAP: snapped geom 1");
00225 
00226         // Run the binary op
00227         GeomPtr result( _Op(snapG0.get(), snapG1.get()) );
00228 
00229 #if GEOS_DEBUG_BINARYOP
00230         check_valid(*result, "SNAP: result (before common-bits addition");
00231 #endif
00232 
00233 #if CBR_BEFORE_SNAPPING
00234         // Add common bits back in
00235         cbr.addCommonBits( result.get() );
00236         result = fix_self_intersections(result, "SNAP: result (after common-bits addition)");
00237 #endif
00238 
00239 #if GEOS_DEBUG_BINARYOP
00240         check_valid(*result, "SNAP: result (after common-bits addition");
00241 #endif
00242 
00243         return result;
00244 }
00245 
00246 template <class BinOp>
00247 std::auto_ptr<Geometry>
00248 BinaryOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00249 {
00250         typedef std::auto_ptr<Geometry> GeomPtr;
00251 
00252         GeomPtr ret;
00253         geos::util::TopologyException origException;
00254 
00255 #ifdef USE_ORIGINAL_INPUT
00256         // Try with original input
00257         try
00258         {
00259 #if GEOS_DEBUG_BINARYOP
00260                 std::cerr << "Trying with original input." << std::endl;
00261 #endif
00262                 ret.reset(_Op(g0, g1));
00263                 return ret;
00264         }
00265         catch (const geos::util::TopologyException& ex)
00266         {
00267                 origException=ex;
00268 #if GEOS_DEBUG_BINARYOP
00269                 std::cerr << "Original exception: " << ex.what() << std::endl;
00270 #endif
00271         }
00272 
00273 //#if GEOS_DEBUG_BINARYOP
00274   // Should we just give up here ?
00275   check_valid(*g0, "Input geom 0");
00276   check_valid(*g1, "Input geom 1");
00277 //#endif
00278 
00279 #endif // USE_ORIGINAL_INPUT
00280 
00281 
00282 #ifdef USE_COMMONBITS_POLICY
00283         // Try removing common bits (possibly obsoleted by snapping below)
00284         //
00285         // NOTE: this policy was _later_ implemented 
00286         //       in JTS as EnhancedPrecisionOp
00287         // TODO: consider using the now-ported EnhancedPrecisionOp
00288         //       here too
00289         // 
00290         try
00291         {
00292                 GeomPtr rG0;
00293                 GeomPtr rG1;
00294                 precision::CommonBitsRemover cbr;
00295 
00296 #if GEOS_DEBUG_BINARYOP
00297                 std::cerr << "Trying with Common Bits Remover (CBR)" << std::endl;
00298 #endif
00299 
00300                 cbr.add(g0);
00301                 cbr.add(g1);
00302 
00303                 rG0.reset( cbr.removeCommonBits(g0->clone()) );
00304                 rG1.reset( cbr.removeCommonBits(g1->clone()) );
00305 
00306 #if GEOS_DEBUG_BINARYOP
00307                 check_valid(*rG0, "CBR: geom 0 (after common-bits removal)");
00308                 check_valid(*rG1, "CBR: geom 1 (after common-bits removal)");
00309 #endif
00310 
00311                 ret.reset( _Op(rG0.get(), rG1.get()) );
00312 
00313 #if GEOS_DEBUG_BINARYOP
00314                 check_valid(*ret, "CBR: result (before common-bits addition)");
00315 #endif
00316 
00317                 cbr.addCommonBits( ret.get() ); 
00318 
00319 #if GEOS_DEBUG_BINARYOP
00320                 check_valid(*ret, "CBR: result (after common-bits addition)");
00321 #endif
00322 
00323                 // Common-bits removal policy could introduce self-intersections
00324                 ret = fix_self_intersections(ret, "CBR: result (after common-bits addition)");
00325 
00326 #if GEOS_DEBUG_BINARYOP
00327                 check_valid(*ret, "CBR: result (after common-bits addition and fix_self_intersections)");
00328 #endif
00329 
00330 #if GEOS_CHECK_COMMONBITS_VALIDITY
00331                 // check that result is a valid geometry after the
00332                 // reshift to orginal precision (see EnhancedPrecisionOp)
00333                 using operation::valid::IsValidOp;
00334                 using operation::valid::TopologyValidationError;
00335                 IsValidOp ivo(ret.get());
00336                 if ( ! ivo.isValid() )
00337                 {
00338                         TopologyValidationError* e = ivo.getValidationError();
00339                         throw geos::util::TopologyException(
00340                                 "Result of overlay became invalid "
00341                                 "after re-addin common bits of operand "
00342                                 "coordinates: " + e->toString(),
00343                                 e->getCoordinate());
00344                 }
00345 #endif // GEOS_CHECK_COMMONBITS_VALIDITY
00346 
00347                 return ret;
00348         }
00349         catch (const geos::util::TopologyException& ex)
00350         {
00351         ::geos::ignore_unused_variable_warning(ex);
00352 #if GEOS_DEBUG_BINARYOP
00353                 std::cerr << "CBR: " << ex.what() << std::endl;
00354 #endif
00355         }
00356 #endif
00357 
00358         // Try with snapping
00359         //
00360         // TODO: possible optimization would be reusing the
00361         //       already common-bit-removed inputs and just
00362         //       apply geometry snapping, whereas the current
00363         //       SnapOp function does both.
00364 // {
00365 #if USE_SNAPPING_POLICY
00366 
00367 #if GEOS_DEBUG_BINARYOP
00368         std::cerr << "Trying with snapping " << std::endl;
00369 #endif
00370 
00371         try {
00372                 ret = SnapOp(g0, g1, _Op);
00373 #if GEOS_DEBUG_BINARYOP
00374         std::cerr << "SnapOp succeeded" << std::endl;
00375 #endif
00376                 return ret;
00377                 
00378         }
00379         catch (const geos::util::TopologyException& ex)
00380         {
00381         ::geos::ignore_unused_variable_warning(ex);
00382 #if GEOS_DEBUG_BINARYOP
00383                 std::cerr << "SNAP: " << ex.what() << std::endl;
00384 #endif
00385         }
00386 
00387 #endif // USE_SNAPPING_POLICY }
00388 
00389 
00390 
00391 // {
00392 #if USE_PRECISION_REDUCTION_POLICY
00393 
00394 
00395         // Try reducing precision
00396         try
00397         {
00398                 int maxPrecision=25;
00399 
00400                 for (int precision=maxPrecision; precision; --precision)
00401                 {
00402                         std::auto_ptr<PrecisionModel> pm(new PrecisionModel(precision));
00403 #if GEOS_DEBUG_BINARYOP
00404                         std::cerr << "Trying with precision " << precision << std::endl;
00405 #endif
00406 
00407                         precision::SimpleGeometryPrecisionReducer reducer( pm.get() );
00408                         GeomPtr rG0( reducer.reduce(g0) );
00409                         GeomPtr rG1( reducer.reduce(g1) );
00410 
00411                         try
00412                         {
00413                                 ret.reset( _Op(rG0.get(), rG1.get()) );
00414                                 return ret;
00415                         }
00416                         catch (const geos::util::TopologyException& ex)
00417                         {
00418                                 if ( precision == 1 ) throw ex;
00419 #if GEOS_DEBUG_BINARYOP
00420                                 std::cerr << "Reduced with precision (" << precision << "): "
00421                                           << ex.what() << std::endl;
00422 #endif
00423                         }
00424 
00425                 }
00426 
00427         }
00428         catch (const geos::util::TopologyException& ex)
00429         {
00430 #if GEOS_DEBUG_BINARYOP
00431                 std::cerr << "Reduced: " << ex.what() << std::endl;
00432 #endif
00433         }
00434 
00435 #endif
00436 // USE_PRECISION_REDUCTION_POLICY }
00437 
00438 // {
00439 #if USE_TP_SIMPLIFY_POLICY 
00440 
00441         // Try simplifying
00442         try
00443         {
00444 
00445                 double maxTolerance = 0.04;
00446                 double minTolerance = 0.01;
00447                 double tolStep = 0.01;
00448 
00449                 for (double tol = minTolerance; tol <= maxTolerance; tol += tolStep)
00450                 {
00451 #if GEOS_DEBUG_BINARYOP
00452                         std::cerr << "Trying simplifying with tolerance " << tol << std::endl;
00453 #endif
00454 
00455                         GeomPtr rG0( simplify::TopologyPreservingSimplifier::simplify(g0, tol) );
00456                         GeomPtr rG1( simplify::TopologyPreservingSimplifier::simplify(g1, tol) );
00457 
00458                         try
00459                         {
00460                                 ret.reset( _Op(rG0.get(), rG1.get()) );
00461                                 return ret;
00462                         }
00463                         catch (const geos::util::TopologyException& ex)
00464                         {
00465                                 if ( tol >= maxTolerance ) throw ex;
00466 #if GEOS_DEBUG_BINARYOP
00467                                 std::cerr << "Simplified with tolerance (" << tol << "): "
00468                                           << ex.what() << std::endl;
00469 #endif
00470                         }
00471 
00472                 }
00473 
00474                 return ret;
00475 
00476         }
00477         catch (const geos::util::TopologyException& ex)
00478         {
00479 #if GEOS_DEBUG_BINARYOP
00480                 std::cerr << "Simplified: " << ex.what() << std::endl;
00481 #endif
00482         }
00483 
00484 #endif
00485 // USE_TP_SIMPLIFY_POLICY }
00486 
00487         throw origException;
00488 }
00489 
00490 
00491 } // namespace geos::geom
00492 } // namespace geos
00493 
00494 #endif // GEOS_GEOM_BINARYOP_H

Generated on Fri Nov 16 16:52:52 2012 for GEOS by  doxygen 1.3.9.1