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 ::geos::ignore_unused_variable_warning(label);
00189 #ifdef GEOS_DEBUG_BINARYOP
00190 std::cerr << label << " fix_self_intersection (UnaryUnion)" << std::endl;
00191 #endif
00192
00193
00194 if ( ! dynamic_cast<const GeometryCollection*>(g.get()) ) return g;
00195
00196 using operation::valid::IsValidOp;
00197
00198 IsValidOp ivo(g.get());
00199
00200
00201 if ( ivo.isValid() ) return g;
00202
00203
00204
00205 using operation::valid::TopologyValidationError;
00206 TopologyValidationError* err = ivo.getValidationError();
00207 switch ( err->getErrorType() ) {
00208 case TopologyValidationError::eRingSelfIntersection:
00209 case TopologyValidationError::eTooFewPoints:
00210 #ifdef GEOS_DEBUG_BINARYOP
00211 std::cerr << label << " ATTEMPT_TO_FIX: " << err->getErrorType() << ": " << *g << std::endl;
00212 #endif
00213 g = g->Union();
00214 #ifdef GEOS_DEBUG_BINARYOP
00215 std::cerr << label << " ATTEMPT_TO_FIX succeeded.. " << std::endl;
00216 #endif
00217 return g;
00218 case TopologyValidationError::eSelfIntersection:
00219
00220 default:
00221 #ifdef GEOS_DEBUG_BINARYOP
00222 std::cerr << label << " invalidity is: " << err->getErrorType() << std::endl;
00223 #endif
00224 return g;
00225 }
00226 }
00227
00228
00234 template <class BinOp>
00235 std::auto_ptr<Geometry>
00236 SnapOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00237 {
00238 typedef std::auto_ptr<Geometry> GeomPtr;
00239
00240 #define CBR_BEFORE_SNAPPING 1
00241
00242
00243 using geos::operation::overlay::snap::GeometrySnapper;
00244
00245
00246
00247 double snapTolerance = GeometrySnapper::computeOverlaySnapTolerance(*g0, *g1);
00248 #if GEOS_DEBUG_BINARYOP
00249 std::cerr<< std::setprecision(20) << "Computed snap tolerance: "<<snapTolerance<<std::endl;
00250 #endif
00251
00252
00253 #if CBR_BEFORE_SNAPPING
00254
00255 geos::precision::CommonBitsRemover cbr;
00256 cbr.add(g0); cbr.add(g1);
00257 #if GEOS_DEBUG_BINARYOP
00258 std::cerr<<"Computed common bits: "<<cbr.getCommonCoordinate()<<std::endl;
00259 #endif
00260
00261
00262 GeomPtr rG0( cbr.removeCommonBits(g0->clone()) );
00263 GeomPtr rG1( cbr.removeCommonBits(g1->clone()) );
00264
00265 #if GEOS_DEBUG_BINARYOP
00266 check_valid(*rG0, "CBR: removed-bits geom 0");
00267 check_valid(*rG1, "CBR: removed-bits geom 1");
00268 #endif
00269
00270 const Geometry& operand0 = *rG0;
00271 const Geometry& operand1 = *rG1;
00272 #else // don't CBR before snapping
00273 const Geometry& operand0 = *g0;
00274 const Geometry& operand1 = *g1;
00275 #endif
00276
00277
00278 GeometrySnapper snapper0( operand0 );
00279 GeomPtr snapG0( snapper0.snapTo(operand1, snapTolerance) );
00280
00281
00282
00283 GeometrySnapper snapper1( operand1 );
00284 GeomPtr snapG1( snapper1.snapTo(*snapG0, snapTolerance) );
00285
00286
00287
00288 GeomPtr result( _Op(snapG0.get(), snapG1.get()) );
00289
00290 #if GEOS_DEBUG_BINARYOP
00291 check_valid(*result, "SNAP: result (before common-bits addition");
00292 #endif
00293
00294 #if CBR_BEFORE_SNAPPING
00295
00296 cbr.addCommonBits( result.get() );
00297
00298
00299 check_valid(*result, "CBR: result (after common-bits addition)", true);
00300
00301 #endif
00302
00303 return result;
00304 }
00305
00306 template <class BinOp>
00307 std::auto_ptr<Geometry>
00308 BinaryOp(const Geometry* g0, const Geometry *g1, BinOp _Op)
00309 {
00310 typedef std::auto_ptr<Geometry> GeomPtr;
00311
00312 GeomPtr ret;
00313 geos::util::TopologyException origException;
00314
00315 #ifdef USE_ORIGINAL_INPUT
00316
00317 try
00318 {
00319 #if GEOS_DEBUG_BINARYOP
00320 std::cerr << "Trying with original input." << std::endl;
00321 #endif
00322 ret.reset(_Op(g0, g1));
00323 return ret;
00324 }
00325 catch (const geos::util::TopologyException& ex)
00326 {
00327 origException=ex;
00328 #if GEOS_DEBUG_BINARYOP
00329 std::cerr << "Original exception: " << ex.what() << std::endl;
00330 #endif
00331 }
00332
00333 check_valid(*g0, "Input geom 0", true, true);
00334 check_valid(*g1, "Input geom 1", true, true);
00335
00336 #if GEOS_DEBUG_BINARYOP
00337
00338 check_valid(*g0, "Input geom 0");
00339 check_valid(*g1, "Input geom 1");
00340 #endif
00341
00342 #endif // USE_ORIGINAL_INPUT
00343
00344 #ifdef USE_COMMONBITS_POLICY
00345
00346
00347
00348
00349
00350
00351
00352 try
00353 {
00354 GeomPtr rG0;
00355 GeomPtr rG1;
00356 precision::CommonBitsRemover cbr;
00357
00358 #if GEOS_DEBUG_BINARYOP
00359 std::cerr << "Trying with Common Bits Remover (CBR)" << std::endl;
00360 #endif
00361
00362 cbr.add(g0);
00363 cbr.add(g1);
00364
00365 rG0.reset( cbr.removeCommonBits(g0->clone()) );
00366 rG1.reset( cbr.removeCommonBits(g1->clone()) );
00367
00368 #if GEOS_DEBUG_BINARYOP
00369 check_valid(*rG0, "CBR: geom 0 (after common-bits removal)");
00370 check_valid(*rG1, "CBR: geom 1 (after common-bits removal)");
00371 #endif
00372
00373 ret.reset( _Op(rG0.get(), rG1.get()) );
00374
00375 #if GEOS_DEBUG_BINARYOP
00376 check_valid(*ret, "CBR: result (before common-bits addition)");
00377 #endif
00378
00379 cbr.addCommonBits( ret.get() );
00380
00381 check_valid(*ret, "CBR: result (after common-bits addition)", true);
00382
00383 #if GEOS_CHECK_COMMONBITS_VALIDITY
00384
00385
00386 using operation::valid::IsValidOp;
00387 using operation::valid::TopologyValidationError;
00388 IsValidOp ivo(ret.get());
00389 if ( ! ivo.isValid() )
00390 {
00391 TopologyValidationError* e = ivo.getValidationError();
00392 throw geos::util::TopologyException(
00393 "Result of overlay became invalid "
00394 "after re-addin common bits of operand "
00395 "coordinates: " + e->toString(),
00396 e->getCoordinate());
00397 }
00398 #endif // GEOS_CHECK_COMMONBITS_VALIDITY
00399
00400 return ret;
00401 }
00402 catch (const geos::util::TopologyException& ex)
00403 {
00404 ::geos::ignore_unused_variable_warning(ex);
00405 #if GEOS_DEBUG_BINARYOP
00406 std::cerr << "CBR: " << ex.what() << std::endl;
00407 #endif
00408 }
00409 #endif
00410
00411
00412
00413
00414
00415
00416
00417
00418 #if USE_SNAPPING_POLICY
00419
00420 #if GEOS_DEBUG_BINARYOP
00421 std::cerr << "Trying with snapping " << std::endl;
00422 #endif
00423
00424 try {
00425 ret = SnapOp(g0, g1, _Op);
00426 #if GEOS_DEBUG_BINARYOP
00427 std::cerr << "SnapOp succeeded" << std::endl;
00428 #endif
00429 return ret;
00430
00431 }
00432 catch (const geos::util::TopologyException& ex)
00433 {
00434 ::geos::ignore_unused_variable_warning(ex);
00435 #if GEOS_DEBUG_BINARYOP
00436 std::cerr << "SNAP: " << ex.what() << std::endl;
00437 #endif
00438 }
00439
00440 #endif // USE_SNAPPING_POLICY }
00441
00442
00443 #if USE_PRECISION_REDUCTION_POLICY
00444
00445
00446
00447 try
00448 {
00449 long unsigned int g0scale =
00450 static_cast<long unsigned int>(g0->getFactory()->getPrecisionModel()->getScale());
00451 long unsigned int g1scale =
00452 static_cast<long unsigned int>(g1->getFactory()->getPrecisionModel()->getScale());
00453
00454 #if GEOS_DEBUG_BINARYOP
00455 std::cerr << "Original input scales are: "
00456 << g0scale
00457 << " and "
00458 << g1scale
00459 << std::endl;
00460 #endif
00461
00462 double maxScale = 1e16;
00463
00464
00465 if ( g0scale && g0scale < maxScale ) maxScale = g0scale;
00466 if ( g1scale && g1scale < maxScale ) maxScale = g1scale;
00467
00468
00469 for (double scale=maxScale; scale >= 1; scale /= 10)
00470 {
00471 PrecisionModel pm(scale);
00472 GeometryFactory::unique_ptr gf = GeometryFactory::create(&pm);
00473 #if GEOS_DEBUG_BINARYOP
00474 std::cerr << "Trying with scale " << scale << std::endl;
00475 #endif
00476
00477 precision::GeometryPrecisionReducer reducer( *gf );
00478 GeomPtr rG0( reducer.reduce(*g0) );
00479 GeomPtr rG1( reducer.reduce(*g1) );
00480
00481 try
00482 {
00483 ret.reset( _Op(rG0.get(), rG1.get()) );
00484
00485 if ( g0->getFactory()->getPrecisionModel()->compareTo( g1->getFactory()->getPrecisionModel() ) < 0 ) {
00486 ret.reset( g0->getFactory()->createGeometry(ret.get()) );
00487 }
00488 else {
00489 ret.reset( g1->getFactory()->createGeometry(ret.get()) );
00490 }
00491 return ret;
00492 }
00493 catch (const geos::util::TopologyException& ex)
00494 {
00495 if ( scale == 1 ) throw ex;
00496 #if GEOS_DEBUG_BINARYOP
00497 std::cerr << "Reduced with scale (" << scale << "): "
00498 << ex.what() << std::endl;
00499 #endif
00500 }
00501
00502 }
00503
00504 }
00505 catch (const geos::util::TopologyException& ex)
00506 {
00507 #if GEOS_DEBUG_BINARYOP
00508 std::cerr << "Reduced: " << ex.what() << std::endl;
00509 #endif
00510 ::geos::ignore_unused_variable_warning(ex);
00511 }
00512
00513 #endif
00514
00515
00516
00517
00518
00519
00520
00521 #if USE_TP_SIMPLIFY_POLICY
00522
00523
00524 try
00525 {
00526
00527 double maxTolerance = 0.04;
00528 double minTolerance = 0.01;
00529 double tolStep = 0.01;
00530
00531 for (double tol = minTolerance; tol <= maxTolerance; tol += tolStep)
00532 {
00533 #if GEOS_DEBUG_BINARYOP
00534 std::cerr << "Trying simplifying with tolerance " << tol << std::endl;
00535 #endif
00536
00537 GeomPtr rG0( simplify::TopologyPreservingSimplifier::simplify(g0, tol) );
00538 GeomPtr rG1( simplify::TopologyPreservingSimplifier::simplify(g1, tol) );
00539
00540 try
00541 {
00542 ret.reset( _Op(rG0.get(), rG1.get()) );
00543 return ret;
00544 }
00545 catch (const geos::util::TopologyException& ex)
00546 {
00547 if ( tol >= maxTolerance ) throw ex;
00548 #if GEOS_DEBUG_BINARYOP
00549 std::cerr << "Simplified with tolerance (" << tol << "): "
00550 << ex.what() << std::endl;
00551 #endif
00552 }
00553
00554 }
00555
00556 return ret;
00557
00558 }
00559 catch (const geos::util::TopologyException& ex)
00560 {
00561 #if GEOS_DEBUG_BINARYOP
00562 std::cerr << "Simplified: " << ex.what() << std::endl;
00563 #endif
00564 }
00565
00566 #endif
00567
00568
00569 throw origException;
00570 }
00571
00572
00573 }
00574 }
00575
00576 #endif // GEOS_GEOM_BINARYOP_H