00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049 #ifndef GEOS_GEOM_BINARYOP_H
00050 #define GEOS_GEOM_BINARYOP_H
00051
00052 #include <geos/algorithm/BoundaryNodeRule.h>
00053 #include <geos/geom/Geometry.h>
00054 #include <geos/geom/GeometryCollection.h>
00055 #include <geos/geom/Polygon.h>
00056 #include <geos/geom/Lineal.h>
00057 #include <geos/geom/PrecisionModel.h>
00058 #include <geos/geom/GeometryFactory.h>
00059 #include <geos/precision/CommonBitsRemover.h>
00060 #include <geos/precision/SimpleGeometryPrecisionReducer.h>
00061 #include <geos/precision/GeometryPrecisionReducer.h>
00062
00063 #include <geos/operation/overlay/snap/GeometrySnapper.h>
00064
00065 #include <geos/simplify/TopologyPreservingSimplifier.h>
00066 #include <geos/operation/IsSimpleOp.h>
00067 #include <geos/operation/valid/IsValidOp.h>
00068 #include <geos/operation/valid/TopologyValidationError.h>
00069 #include <geos/util/TopologyException.h>
00070 #include <geos/util.h>
00071
00072 #include <memory>
00073
00074
00075
00076 #ifdef GEOS_DEBUG_BINARYOP
00077 # include <iostream>
00078 # include <iomanip>
00079 # include <sstream>
00080 #endif
00081
00082
00083
00084
00085
00086 #ifndef USE_ORIGINAL_INPUT
00087 # define USE_ORIGINAL_INPUT 1
00088 #endif
00089
00090
00091
00092
00093
00094
00095 #ifndef USE_PRECISION_REDUCTION_POLICY
00096 # define USE_PRECISION_REDUCTION_POLICY 1
00097 #endif
00098
00099
00100
00101
00102
00103
00104 #ifndef USE_TP_SIMPLIFY_POLICY
00105
00106 #endif
00107
00108
00109
00110
00111
00112
00113 #ifndef USE_COMMONBITS_POLICY
00114 # define USE_COMMONBITS_POLICY 1
00115 #endif
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128 #define GEOS_CHECK_COMMONBITS_VALIDITY 1
00129
00130
00131
00132
00133 #ifndef USE_SNAPPING_POLICY
00134 # define USE_SNAPPING_POLICY 1
00135 #endif
00136
00137 namespace geos {
00138 namespace geom {
00139
00140 inline bool
00141 check_valid(const Geometry& g, const std::string& label, bool doThrow=false, bool validOnly=false)
00142 {
00143 if ( dynamic_cast<const Lineal*>(&g) ) {
00144 if ( ! validOnly ) {
00145 operation::IsSimpleOp sop(g, algorithm::BoundaryNodeRule::getBoundaryEndPoint());
00146 if ( ! sop.isSimple() )
00147 {
00148 if ( doThrow ) {
00149 throw geos::util::TopologyException(
00150 label + " is not simple");
00151 }
00152 return false;
00153 }
00154 }
00155 } else {
00156 operation::valid::IsValidOp ivo(&g);
00157 if ( ! ivo.isValid() )
00158 {
00159 using operation::valid::TopologyValidationError;
00160 TopologyValidationError* err = ivo.getValidationError();
00161 #ifdef GEOS_DEBUG_BINARYOP
00162 std::cerr << label << " is INVALID: "
00163 << err->toString()
00164 << " (" << std::setprecision(20)
00165 << err->getCoordinate() << ")"
00166 << std::endl;
00167 #endif
00168 if ( doThrow ) {
00169 throw geos::util::TopologyException(
00170 label + " is invalid: " + err->toString(),
00171 err->getCoordinate());
00172 }
00173 return false;
00174 }
00175 }
00176 return true;
00177 }
00178
00179
00180
00181
00182
00183
00184
00185 inline std::auto_ptr<Geometry>
00186 fix_self_intersections(std::auto_ptr<Geometry> g, const std::string& label)
00187 {
00188 #ifdef GEOS_DEBUG_BINARYOP
00189 std::cerr << label << " fix_self_intersection (UnaryUnion)" << std::endl;
00190 #endif
00191
00192
00193 if ( ! dynamic_cast<const GeometryCollection*>(g.get()) ) return g;
00194
00195 using operation::valid::IsValidOp;
00196
00197 IsValidOp ivo(g.get());
00198
00199
00200 if ( ivo.isValid() ) return g;
00201
00202
00203
00204 using operation::valid::TopologyValidationError;
00205 TopologyValidationError* err = ivo.getValidationError();
00206 switch ( err->getErrorType() ) {
00207 case TopologyValidationError::eRingSelfIntersection:
00208 case TopologyValidationError::eTooFewPoints:
00209 #ifdef GEOS_DEBUG_BINARYOP
00210 std::cerr << label << " ATTEMPT_TO_FIX: " << err->getErrorType() << ": " << *g << std::endl;
00211 #endif
00212 g = g->Union();
00213 #ifdef GEOS_DEBUG_BINARYOP
00214 std::cerr << label << " ATTEMPT_TO_FIX succeeded.. " << std::endl;
00215 #endif
00216 return g;
00217 case TopologyValidationError::eSelfIntersection:
00218
00219 default:
00220 #ifdef GEOS_DEBUG_BINARYOP
00221 std::cerr << label << " invalidity is: " << err->getErrorType() << std::endl;
00222 #endif
00223 return g;
00224 }
00225 }
00226
00227
00233 template <class BinOp>
00234 std::auto_ptr<Geometry>
00235 SnapOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00236 {
00237 typedef std::auto_ptr<Geometry> GeomPtr;
00238
00239 #define CBR_BEFORE_SNAPPING 1
00240
00241
00242 using geos::operation::overlay::snap::GeometrySnapper;
00243
00244
00245
00246 double snapTolerance = GeometrySnapper::computeOverlaySnapTolerance(*g0, *g1);
00247 #if GEOS_DEBUG_BINARYOP
00248 std::cerr<< std::setprecision(20) << "Computed snap tolerance: "<<snapTolerance<<std::endl;
00249 #endif
00250
00251
00252 #if CBR_BEFORE_SNAPPING
00253
00254 geos::precision::CommonBitsRemover cbr;
00255 cbr.add(g0); cbr.add(g1);
00256 #if GEOS_DEBUG_BINARYOP
00257 std::cerr<<"Computed common bits: "<<cbr.getCommonCoordinate()<<std::endl;
00258 #endif
00259
00260
00261 GeomPtr rG0( cbr.removeCommonBits(g0->clone()) );
00262 GeomPtr rG1( cbr.removeCommonBits(g1->clone()) );
00263
00264 #if GEOS_DEBUG_BINARYOP
00265 check_valid(*rG0, "CBR: removed-bits geom 0");
00266 check_valid(*rG1, "CBR: removed-bits geom 1");
00267 #endif
00268
00269 const Geometry& operand0 = *rG0;
00270 const Geometry& operand1 = *rG1;
00271 #else // don't CBR before snapping
00272 const Geometry& operand0 = *g0;
00273 const Geometry& operand1 = *g1;
00274 #endif
00275
00276
00277 GeometrySnapper snapper0( operand0 );
00278 GeomPtr snapG0( snapper0.snapTo(operand1, snapTolerance) );
00279
00280
00281
00282 GeometrySnapper snapper1( operand1 );
00283 GeomPtr snapG1( snapper1.snapTo(*snapG0, snapTolerance) );
00284
00285
00286
00287 GeomPtr result( _Op(snapG0.get(), snapG1.get()) );
00288
00289 #if GEOS_DEBUG_BINARYOP
00290 check_valid(*result, "SNAP: result (before common-bits addition");
00291 #endif
00292
00293 #if CBR_BEFORE_SNAPPING
00294
00295 cbr.addCommonBits( result.get() );
00296
00297
00298 check_valid(*result, "CBR: result (after common-bits addition)", true);
00299
00300 #endif
00301
00302 return result;
00303 }
00304
00305 template <class BinOp>
00306 std::auto_ptr<Geometry>
00307 BinaryOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00308 {
00309 typedef std::auto_ptr<Geometry> GeomPtr;
00310
00311 GeomPtr ret;
00312 geos::util::TopologyException origException;
00313
00314 #ifdef USE_ORIGINAL_INPUT
00315
00316 try
00317 {
00318 #if GEOS_DEBUG_BINARYOP
00319 std::cerr << "Trying with original input." << std::endl;
00320 #endif
00321 ret.reset(_Op(g0, g1));
00322 return ret;
00323 }
00324 catch (const geos::util::TopologyException& ex)
00325 {
00326 origException=ex;
00327 #if GEOS_DEBUG_BINARYOP
00328 std::cerr << "Original exception: " << ex.what() << std::endl;
00329 #endif
00330 }
00331
00332 check_valid(*g0, "Input geom 0", true, true);
00333 check_valid(*g1, "Input geom 1", true, true);
00334
00335 #if GEOS_DEBUG_BINARYOP
00336
00337 check_valid(*g0, "Input geom 0");
00338 check_valid(*g1, "Input geom 1");
00339 #endif
00340
00341 #endif // USE_ORIGINAL_INPUT
00342
00343 #ifdef USE_COMMONBITS_POLICY
00344
00345
00346
00347
00348
00349
00350
00351 try
00352 {
00353 GeomPtr rG0;
00354 GeomPtr rG1;
00355 precision::CommonBitsRemover cbr;
00356
00357 #if GEOS_DEBUG_BINARYOP
00358 std::cerr << "Trying with Common Bits Remover (CBR)" << std::endl;
00359 #endif
00360
00361 cbr.add(g0);
00362 cbr.add(g1);
00363
00364 rG0.reset( cbr.removeCommonBits(g0->clone()) );
00365 rG1.reset( cbr.removeCommonBits(g1->clone()) );
00366
00367 #if GEOS_DEBUG_BINARYOP
00368 check_valid(*rG0, "CBR: geom 0 (after common-bits removal)");
00369 check_valid(*rG1, "CBR: geom 1 (after common-bits removal)");
00370 #endif
00371
00372 ret.reset( _Op(rG0.get(), rG1.get()) );
00373
00374 #if GEOS_DEBUG_BINARYOP
00375 check_valid(*ret, "CBR: result (before common-bits addition)");
00376 #endif
00377
00378 cbr.addCommonBits( ret.get() );
00379
00380 check_valid(*ret, "CBR: result (after common-bits addition)", true);
00381
00382 #if GEOS_CHECK_COMMONBITS_VALIDITY
00383
00384
00385 using operation::valid::IsValidOp;
00386 using operation::valid::TopologyValidationError;
00387 IsValidOp ivo(ret.get());
00388 if ( ! ivo.isValid() )
00389 {
00390 TopologyValidationError* e = ivo.getValidationError();
00391 throw geos::util::TopologyException(
00392 "Result of overlay became invalid "
00393 "after re-addin common bits of operand "
00394 "coordinates: " + e->toString(),
00395 e->getCoordinate());
00396 }
00397 #endif // GEOS_CHECK_COMMONBITS_VALIDITY
00398
00399 return ret;
00400 }
00401 catch (const geos::util::TopologyException& ex)
00402 {
00403 ::geos::ignore_unused_variable_warning(ex);
00404 #if GEOS_DEBUG_BINARYOP
00405 std::cerr << "CBR: " << ex.what() << std::endl;
00406 #endif
00407 }
00408 #endif
00409
00410
00411
00412
00413
00414
00415
00416
00417 #if USE_SNAPPING_POLICY
00418
00419 #if GEOS_DEBUG_BINARYOP
00420 std::cerr << "Trying with snapping " << std::endl;
00421 #endif
00422
00423 try {
00424 ret = SnapOp(g0, g1, _Op);
00425 #if GEOS_DEBUG_BINARYOP
00426 std::cerr << "SnapOp succeeded" << std::endl;
00427 #endif
00428 return ret;
00429
00430 }
00431 catch (const geos::util::TopologyException& ex)
00432 {
00433 ::geos::ignore_unused_variable_warning(ex);
00434 #if GEOS_DEBUG_BINARYOP
00435 std::cerr << "SNAP: " << ex.what() << std::endl;
00436 #endif
00437 }
00438
00439 #endif // USE_SNAPPING_POLICY }
00440
00441
00442 #if USE_PRECISION_REDUCTION_POLICY
00443
00444
00445
00446 try
00447 {
00448 long unsigned int g0scale = g0->getFactory()->getPrecisionModel()->getScale();
00449 long unsigned int g1scale = g1->getFactory()->getPrecisionModel()->getScale();
00450
00451 #if GEOS_DEBUG_BINARYOP
00452 std::cerr << "Original input scales are: "
00453 << g0scale
00454 << " and "
00455 << g1scale
00456 << std::endl;
00457 #endif
00458
00459 double maxScale = 1e16;
00460
00461
00462 if ( g0scale && g0scale < maxScale ) maxScale = g0scale;
00463 if ( g1scale && g1scale < maxScale ) maxScale = g1scale;
00464
00465
00466 for (double scale=maxScale; scale >= 1; scale /= 10)
00467 {
00468 PrecisionModel pm(scale);
00469 GeometryFactory gf(&pm);
00470 #if GEOS_DEBUG_BINARYOP
00471 std::cerr << "Trying with scale " << scale << std::endl;
00472 #endif
00473
00474 precision::GeometryPrecisionReducer reducer( gf );
00475 GeomPtr rG0( reducer.reduce(*g0) );
00476 GeomPtr rG1( reducer.reduce(*g1) );
00477
00478 try
00479 {
00480 ret.reset( _Op(rG0.get(), rG1.get()) );
00481
00482 if ( g0->getFactory()->getPrecisionModel()->compareTo( g1->getFactory()->getPrecisionModel() ) < 0 ) {
00483 ret.reset( g0->getFactory()->createGeometry(ret.get()) );
00484 }
00485 else {
00486 ret.reset( g1->getFactory()->createGeometry(ret.get()) );
00487 }
00488 return ret;
00489 }
00490 catch (const geos::util::TopologyException& ex)
00491 {
00492 if ( scale == 1 ) throw ex;
00493 #if GEOS_DEBUG_BINARYOP
00494 std::cerr << "Reduced with scale (" << scale << "): "
00495 << ex.what() << std::endl;
00496 #endif
00497 }
00498
00499 }
00500
00501 }
00502 catch (const geos::util::TopologyException& ex)
00503 {
00504 #if GEOS_DEBUG_BINARYOP
00505 std::cerr << "Reduced: " << ex.what() << std::endl;
00506 #endif
00507 }
00508
00509 #endif
00510
00511
00512
00513
00514
00515
00516
00517 #if USE_TP_SIMPLIFY_POLICY
00518
00519
00520 try
00521 {
00522
00523 double maxTolerance = 0.04;
00524 double minTolerance = 0.01;
00525 double tolStep = 0.01;
00526
00527 for (double tol = minTolerance; tol <= maxTolerance; tol += tolStep)
00528 {
00529 #if GEOS_DEBUG_BINARYOP
00530 std::cerr << "Trying simplifying with tolerance " << tol << std::endl;
00531 #endif
00532
00533 GeomPtr rG0( simplify::TopologyPreservingSimplifier::simplify(g0, tol) );
00534 GeomPtr rG1( simplify::TopologyPreservingSimplifier::simplify(g1, tol) );
00535
00536 try
00537 {
00538 ret.reset( _Op(rG0.get(), rG1.get()) );
00539 return ret;
00540 }
00541 catch (const geos::util::TopologyException& ex)
00542 {
00543 if ( tol >= maxTolerance ) throw ex;
00544 #if GEOS_DEBUG_BINARYOP
00545 std::cerr << "Simplified with tolerance (" << tol << "): "
00546 << ex.what() << std::endl;
00547 #endif
00548 }
00549
00550 }
00551
00552 return ret;
00553
00554 }
00555 catch (const geos::util::TopologyException& ex)
00556 {
00557 #if GEOS_DEBUG_BINARYOP
00558 std::cerr << "Simplified: " << ex.what() << std::endl;
00559 #endif
00560 }
00561
00562 #endif
00563
00564
00565 throw origException;
00566 }
00567
00568
00569 }
00570 }
00571
00572 #endif // GEOS_GEOM_BINARYOP_H