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

BinaryOp.h

00001 /**********************************************************************
00002  * $Id: BinaryOp.h 3606 2012-04-10 15:02:23Z 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   // TODO: check for presence of self-intersections first ?
00160   return g->Union();
00161 }
00162 
00163 
00169 template <class BinOp>
00170 std::auto_ptr<Geometry>
00171 SnapOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00172 {
00173         typedef std::auto_ptr<Geometry> GeomPtr;
00174 
00175 #define CBR_BEFORE_SNAPPING 1
00176 
00177         //using geos::precision::GeometrySnapper;
00178         using geos::operation::overlay::snap::GeometrySnapper;
00179 
00180         // Snap tolerance must be computed on the original
00181         // (not commonbits-removed) geoms
00182         double snapTolerance = GeometrySnapper::computeOverlaySnapTolerance(*g0, *g1);
00183 #if GEOS_DEBUG_BINARYOP
00184         std::cerr<< std::setprecision(20) << "Computed snap tolerance: "<<snapTolerance<<std::endl;
00185 #endif
00186 
00187 
00188 #if CBR_BEFORE_SNAPPING
00189         // Compute common bits
00190         geos::precision::CommonBitsRemover cbr;
00191         cbr.add(g0); cbr.add(g1);
00192 #if GEOS_DEBUG_BINARYOP
00193         std::cerr<<"Computed common bits: "<<cbr.getCommonCoordinate()<<std::endl;
00194 #endif
00195 
00196         // Now remove common bits
00197         GeomPtr rG0( cbr.removeCommonBits(g0->clone()) );
00198         GeomPtr rG1( cbr.removeCommonBits(g1->clone()) );
00199 
00200 #if GEOS_DEBUG_BINARYOP
00201         check_valid(*rG0, "CBR: removed-bits geom 0");
00202         check_valid(*rG1, "CBR: removed-bits geom 1");
00203 #endif
00204 
00205         const Geometry& operand0 = *rG0;
00206         const Geometry& operand1 = *rG1;
00207 #else // don't CBR before snapping
00208         const Geometry& operand0 = *g0;
00209         const Geometry& operand1 = *g1;
00210 #endif
00211 
00212 
00213         GeometrySnapper snapper0( operand0 );
00214         GeomPtr snapG0( snapper0.snapTo(operand1, snapTolerance) );
00215         snapG0 = fix_self_intersections(snapG0, "SNAP: snapped geom 0");
00216 
00217         // NOTE: second geom is snapped on the snapped first one
00218         GeometrySnapper snapper1( operand1 );
00219         GeomPtr snapG1( snapper1.snapTo(*snapG0, snapTolerance) );
00220         snapG1 = fix_self_intersections(snapG1, "SNAP: snapped geom 1");
00221 
00222         // Run the binary op
00223         GeomPtr result( _Op(snapG0.get(), snapG1.get()) );
00224 
00225 #if GEOS_DEBUG_BINARYOP
00226         check_valid(*result, "SNAP: result (before common-bits addition");
00227 #endif
00228 
00229 #if CBR_BEFORE_SNAPPING
00230         // Add common bits back in
00231         cbr.addCommonBits( result.get() );
00232         result = fix_self_intersections(result, "SNAP: result (after common-bits addition)");
00233 #endif
00234 
00235 #if GEOS_DEBUG_BINARYOP
00236         check_valid(*result, "SNAP: result (after common-bits addition");
00237 #endif
00238 
00239         return result;
00240 }
00241 
00242 template <class BinOp>
00243 std::auto_ptr<Geometry>
00244 BinaryOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00245 {
00246         typedef std::auto_ptr<Geometry> GeomPtr;
00247 
00248         GeomPtr ret;
00249         geos::util::TopologyException origException;
00250 
00251 #ifdef USE_ORIGINAL_INPUT
00252         // Try with original input
00253         try
00254         {
00255 #if GEOS_DEBUG_BINARYOP
00256                 std::cerr << "Trying with original input." << std::endl;
00257 #endif
00258                 ret.reset(_Op(g0, g1));
00259                 return ret;
00260         }
00261         catch (const geos::util::TopologyException& ex)
00262         {
00263                 origException=ex;
00264 #if GEOS_DEBUG_BINARYOP
00265                 std::cerr << "Original exception: " << ex.what() << std::endl;
00266 #endif
00267         }
00268 
00269 //#if GEOS_DEBUG_BINARYOP
00270   // Should we just give up here ?
00271   check_valid(*g0, "Input geom 0");
00272   check_valid(*g1, "Input geom 1");
00273 //#endif
00274 
00275 #endif // USE_ORIGINAL_INPUT
00276 
00277 
00278 #ifdef USE_COMMONBITS_POLICY
00279         // Try removing common bits (possibly obsoleted by snapping below)
00280         //
00281         // NOTE: this policy was _later_ implemented 
00282         //       in JTS as EnhancedPrecisionOp
00283         // TODO: consider using the now-ported EnhancedPrecisionOp
00284         //       here too
00285         // 
00286         try
00287         {
00288                 GeomPtr rG0;
00289                 GeomPtr rG1;
00290                 precision::CommonBitsRemover cbr;
00291 
00292 #if GEOS_DEBUG_BINARYOP
00293                 std::cerr << "Trying with Common Bits Remover (CBR)" << std::endl;
00294 #endif
00295 
00296                 cbr.add(g0);
00297                 cbr.add(g1);
00298 
00299                 rG0.reset( cbr.removeCommonBits(g0->clone()) );
00300                 rG1.reset( cbr.removeCommonBits(g1->clone()) );
00301 
00302 #if GEOS_DEBUG_BINARYOP
00303                 check_valid(*rG0, "CBR: geom 0 (after common-bits removal)");
00304                 check_valid(*rG1, "CBR: geom 1 (after common-bits removal)");
00305 #endif
00306 
00307                 ret.reset( _Op(rG0.get(), rG1.get()) );
00308 
00309 #if GEOS_DEBUG_BINARYOP
00310                 check_valid(*ret, "CBR: result (before common-bits addition)");
00311 #endif
00312 
00313                 cbr.addCommonBits( ret.get() ); 
00314 
00315 #if GEOS_DEBUG_BINARYOP
00316                 check_valid(*ret, "CBR: result (after common-bits addition)");
00317 #endif
00318 
00319                 // Common-bits removal policy could introduce self-intersections
00320                 ret = fix_self_intersections(ret, "CBR: result (after common-bits addition)");
00321 
00322 #if GEOS_DEBUG_BINARYOP
00323                 check_valid(*ret, "CBR: result (after common-bits addition and fix_self_intersections)");
00324 #endif
00325 
00326 #if GEOS_CHECK_COMMONBITS_VALIDITY
00327                 // check that result is a valid geometry after the
00328                 // reshift to orginal precision (see EnhancedPrecisionOp)
00329                 using operation::valid::IsValidOp;
00330                 using operation::valid::TopologyValidationError;
00331                 IsValidOp ivo(ret.get());
00332                 if ( ! ivo.isValid() )
00333                 {
00334                         TopologyValidationError* e = ivo.getValidationError();
00335                         throw geos::util::TopologyException(
00336                                 "Result of overlay became invalid "
00337                                 "after re-addin common bits of operand "
00338                                 "coordinates: " + e->toString(),
00339                                 e->getCoordinate());
00340                 }
00341 #endif // GEOS_CHECK_COMMONBITS_VALIDITY
00342 
00343                 return ret;
00344         }
00345         catch (const geos::util::TopologyException& ex)
00346         {
00347         ::geos::ignore_unused_variable_warning(ex);
00348 #if GEOS_DEBUG_BINARYOP
00349                 std::cerr << "CBR: " << ex.what() << std::endl;
00350 #endif
00351         }
00352 #endif
00353 
00354         // Try with snapping
00355         //
00356         // TODO: possible optimization would be reusing the
00357         //       already common-bit-removed inputs and just
00358         //       apply geometry snapping, whereas the current
00359         //       SnapOp function does both.
00360 // {
00361 #if USE_SNAPPING_POLICY
00362 
00363 #if GEOS_DEBUG_BINARYOP
00364         std::cerr << "Trying with snapping " << std::endl;
00365 #endif
00366 
00367         try {
00368                 ret = SnapOp(g0, g1, _Op);
00369 #if GEOS_DEBUG_BINARYOP
00370         std::cerr << "SnapOp succeeded" << std::endl;
00371 #endif
00372                 return ret;
00373                 
00374         }
00375         catch (const geos::util::TopologyException& ex)
00376         {
00377         ::geos::ignore_unused_variable_warning(ex);
00378 #if GEOS_DEBUG_BINARYOP
00379                 std::cerr << "SNAP: " << ex.what() << std::endl;
00380 #endif
00381         }
00382 
00383 #endif // USE_SNAPPING_POLICY }
00384 
00385 
00386 
00387 // {
00388 #if USE_PRECISION_REDUCTION_POLICY
00389 
00390 
00391         // Try reducing precision
00392         try
00393         {
00394                 int maxPrecision=25;
00395 
00396                 for (int precision=maxPrecision; precision; --precision)
00397                 {
00398                         std::auto_ptr<PrecisionModel> pm(new PrecisionModel(precision));
00399 #if GEOS_DEBUG_BINARYOP
00400                         std::cerr << "Trying with precision " << precision << std::endl;
00401 #endif
00402 
00403                         precision::SimpleGeometryPrecisionReducer reducer( pm.get() );
00404                         GeomPtr rG0( reducer.reduce(g0) );
00405                         GeomPtr rG1( reducer.reduce(g1) );
00406 
00407                         try
00408                         {
00409                                 ret.reset( _Op(rG0.get(), rG1.get()) );
00410                                 return ret;
00411                         }
00412                         catch (const geos::util::TopologyException& ex)
00413                         {
00414                                 if ( precision == 1 ) throw ex;
00415 #if GEOS_DEBUG_BINARYOP
00416                                 std::cerr << "Reduced with precision (" << precision << "): "
00417                                           << ex.what() << std::endl;
00418 #endif
00419                         }
00420 
00421                 }
00422 
00423         }
00424         catch (const geos::util::TopologyException& ex)
00425         {
00426 #if GEOS_DEBUG_BINARYOP
00427                 std::cerr << "Reduced: " << ex.what() << std::endl;
00428 #endif
00429         }
00430 
00431 #endif
00432 // USE_PRECISION_REDUCTION_POLICY }
00433 
00434 // {
00435 #if USE_TP_SIMPLIFY_POLICY 
00436 
00437         // Try simplifying
00438         try
00439         {
00440 
00441                 double maxTolerance = 0.04;
00442                 double minTolerance = 0.01;
00443                 double tolStep = 0.01;
00444 
00445                 for (double tol = minTolerance; tol <= maxTolerance; tol += tolStep)
00446                 {
00447 #if GEOS_DEBUG_BINARYOP
00448                         std::cerr << "Trying simplifying with tolerance " << tol << std::endl;
00449 #endif
00450 
00451                         GeomPtr rG0( simplify::TopologyPreservingSimplifier::simplify(g0, tol) );
00452                         GeomPtr rG1( simplify::TopologyPreservingSimplifier::simplify(g1, tol) );
00453 
00454                         try
00455                         {
00456                                 ret.reset( _Op(rG0.get(), rG1.get()) );
00457                                 return ret;
00458                         }
00459                         catch (const geos::util::TopologyException& ex)
00460                         {
00461                                 if ( tol >= maxTolerance ) throw ex;
00462 #if GEOS_DEBUG_BINARYOP
00463                                 std::cerr << "Simplified with tolerance (" << tol << "): "
00464                                           << ex.what() << std::endl;
00465 #endif
00466                         }
00467 
00468                 }
00469 
00470                 return ret;
00471 
00472         }
00473         catch (const geos::util::TopologyException& ex)
00474         {
00475 #if GEOS_DEBUG_BINARYOP
00476                 std::cerr << "Simplified: " << ex.what() << std::endl;
00477 #endif
00478         }
00479 
00480 #endif
00481 // USE_TP_SIMPLIFY_POLICY }
00482 
00483         throw origException;
00484 }
00485 
00486 
00487 } // namespace geos::geom
00488 } // namespace geos
00489 
00490 #endif // GEOS_GEOM_BINARYOP_H

Generated on Tue Jun 5 14:58:07 2012 for GEOS by  doxygen 1.3.9.1