diff --git a/src/Lib/terragear/tg_chopper.cxx b/src/Lib/terragear/tg_chopper.cxx index 092678bf..3fbf8721 100644 --- a/src/Lib/terragear/tg_chopper.cxx +++ b/src/Lib/terragear/tg_chopper.cxx @@ -10,7 +10,9 @@ #include "tg_chopper.hxx" #include "tg_shapefile.hxx" +#if HAS_CLIP_ROW unsigned int strip_clip = 0; +#endif tgPolygon tgChopper::Clip( const tgPolygon& subject, const std::string& type, @@ -66,6 +68,7 @@ tgPolygon tgChopper::Clip( const tgPolygon& subject, return result; } +#if HAS_CLIP_ROW // Pass in the center lat for clipping buckets from the row. // We can't rely on sgBucketOffset, as rounding error sometimes causes it to look like there are 2 rows // (the first being a sliver) @@ -85,6 +88,9 @@ void tgChopper::ClipRow( const tgPolygon& subject, const double& center_lat, con Clip( subject, type, b_cur ); } } +#endif + +#if HAS_CLIP_ROW void tgChopper::Add( const tgPolygon& subject, const std::string& type ) { @@ -157,6 +163,125 @@ void tgChopper::Add( const tgPolygon& subject, const std::string& type ) } } +#else + +void tgChopper::Add( const tgPolygon& subject, const std::string& type ) +{ + tgRectangle bb; + SGGeod p; + + // bail out immediately if polygon is empty + if ( subject.Contours() == 0 ) + return; + + bb = subject.GetBoundingBox(); + + SG_LOG( SG_GENERAL, SG_DEBUG, " min = " << bb.getMin() << " max = " << bb.getMax() ); + + // find buckets for min, and max points of convex hull. + // note to self: self, you should think about checking for + // polygons that span the date line + SGBucket b_min( bb.getMin() ); + SGBucket b_max( bb.getMax() ); + SG_LOG( SG_GENERAL, SG_DEBUG, " Bucket min = " << b_min ); + SG_LOG( SG_GENERAL, SG_DEBUG, " Bucket max = " << b_max ); + + if ( b_min == b_max ) { + // shape entirely contained in a single bucket, write and bail + Clip( subject, type, b_min ); + } else { + SGBucket b_cur; + int dx, dy; + + sgBucketDiff(b_min, b_max, &dx, &dy); + SG_LOG( SG_GENERAL, SG_DEBUG, " polygon spans tile boundaries" ); + SG_LOG( SG_GENERAL, SG_DEBUG, " dx = " << dx << " dy = " << dy ); + + if ( (dx > 2880) || (dy > 1440) ) + throw sg_exception("something is really wrong in split_polygon()!!!!"); + + if ( dy <= 1 ) { + // we are down to at most two rows, write each column and then bail + double min_center_lat = b_min.get_center_lat(); + double min_center_lon = b_min.get_center_lon(); + for ( int j = 0; j <= dy; ++j ) { + for ( int i = 0; i <= dx; ++i ) { + b_cur = sgBucketOffset(min_center_lon, min_center_lat, i, j); + Clip( subject, type, b_cur ); + } + } + return; + } + + // we have two or more rows left, split in half (along a + // horizontal dividing line) and recurse with each half + + // find mid point (integer math) + int mid = (dy + 1) / 2 - 1; + + // determine horizontal clip line + SGBucket b_clip = sgBucketOffset( bb.getMin().getLongitudeDeg(), bb.getMin().getLatitudeDeg(), 0, mid); + double clip_line = b_clip.get_center_lat(); + if ( (clip_line >= -90.0 + SG_HALF_BUCKET_SPAN) + && (clip_line < 90.0 - SG_HALF_BUCKET_SPAN) ) + clip_line += SG_HALF_BUCKET_SPAN; + else if ( clip_line < -89.0 ) + clip_line = -89.0; + else if ( clip_line >= 89.0 ) + clip_line = 90.0; + else { + SG_LOG( SG_GENERAL, SG_ALERT, "Out of range latitude in clip_and_write_poly() = " << clip_line ); + } + + { + // + // Crop bottom area (hopefully by putting this in it's own + // scope we can shorten the life of some really large data + // structures to reduce memory use) + // + + SG_LOG( SG_GENERAL, SG_DEBUG, "Generating bottom half (" << bb.getMin().getLatitudeDeg() << "-" << clip_line << ")" ); + + tgPolygon bottom, bottom_clip; + + bottom.AddNode( 0, SGGeod::fromDeg(-180.0, bb.getMin().getLatitudeDeg()) ); + bottom.AddNode( 0, SGGeod::fromDeg( 180.0, bb.getMin().getLatitudeDeg()) ); + bottom.AddNode( 0, SGGeod::fromDeg( 180.0, clip_line) ); + bottom.AddNode( 0, SGGeod::fromDeg(-180.0, clip_line) ); + + bottom_clip = tgPolygon::Intersect( subject, bottom ); + + if ( (bottom_clip.TotalNodes() > 0) && (bottom_clip.TotalNodes() != subject.TotalNodes()) ) { + Add( bottom_clip, type ); + } + } + + { + // + // Crop top area (hopefully by putting this in it's own scope + // we can shorten the life of some really large data + // structures to reduce memory use) + // + + SG_LOG( SG_GENERAL, SG_DEBUG, "Generating top half (" << clip_line << "-" << bb.getMax().getLatitudeDeg() << ")" ); + + tgPolygon top, top_clip; + + top.AddNode( 0, SGGeod::fromDeg(-180.0, clip_line) ); + top.AddNode( 0, SGGeod::fromDeg( 180.0, clip_line) ); + top.AddNode( 0, SGGeod::fromDeg( 180.0, bb.getMax().getLatitudeDeg()) ); + top.AddNode( 0, SGGeod::fromDeg(-180.0, bb.getMax().getLatitudeDeg()) ); + + top_clip = tgPolygon::Intersect( subject, top ); + + if ( (top_clip.TotalNodes() > 0) && (top_clip.TotalNodes() != subject.TotalNodes()) ) { + Add( top_clip, type ); + } + } + } +} +#endif + long int tgChopper::GenerateIndex( std::string path ) { std::string index_file = path + "/chop.idx"; diff --git a/src/Lib/terragear/tg_chopper.hxx b/src/Lib/terragear/tg_chopper.hxx index 3c526ff5..c46b2f1d 100644 --- a/src/Lib/terragear/tg_chopper.hxx +++ b/src/Lib/terragear/tg_chopper.hxx @@ -2,6 +2,8 @@ #include "tg_polygon.hxx" +#define HAS_CLIP_ROW (0) + // for ogr-decode : generate a bunch of polygons, mapped by bucket id typedef std::map bucket_polys_map; typedef bucket_polys_map::iterator bucket_polys_map_interator; @@ -18,7 +20,11 @@ public: private: long int GenerateIndex( std::string path ); + +#if HAS_CLIP_ROW void ClipRow( const tgPolygon& subject, const double& center_lat, const std::string& type ); +#endif + tgPolygon Clip( const tgPolygon& subject, const std::string& type, SGBucket& b ); void Chop( const tgPolygon& subject, const std::string& type ); diff --git a/src/Lib/terragear/tg_contour.cxx b/src/Lib/terragear/tg_contour.cxx index 35ba70c6..f79a53b4 100644 --- a/src/Lib/terragear/tg_contour.cxx +++ b/src/Lib/terragear/tg_contour.cxx @@ -73,6 +73,7 @@ double tgContour::GetArea( void ) const return fabs(area * 0.5); } +#if HAS_IS_INSIDE bool tgContour::IsInside( const tgContour& inside, const tgContour& outside ) { // first contour is inside second if the intersection of first with second is == first @@ -89,6 +90,7 @@ bool tgContour::IsInside( const tgContour& inside, const tgContour& outside ) return isInside; } +#endif bool tgContour::RemoveCycles( const tgContour& subject, tgcontour_list& result ) { @@ -133,6 +135,7 @@ bool tgContour::RemoveCycles( const tgContour& subject, tgcontour_list& result ) second.AddNode( subject.GetNode(n) ); } +#if HAS_IS_INSIDE // determine hole vs boundary if ( IsInside( first, second ) ) { SG_LOG(SG_GENERAL, SG_DEBUG, "first contur is within second contour " ); @@ -153,6 +156,7 @@ bool tgContour::RemoveCycles( const tgContour& subject, tgcontour_list& result ) first.SetHole( subject.GetHole() ); second.SetHole( subject.GetHole() ); } +#endif SG_LOG(SG_GENERAL, SG_DEBUG, "remove first: size " << first.GetSize() ); first.SetHole( subject.GetHole() ); @@ -498,6 +502,7 @@ tgPolygon tgContour::Diff( const tgContour& subject, tgPolygon& clip ) return result; } +#if HAS_INTERSECT tgPolygon tgContour::Intersect( const tgContour& subject, const tgContour& clip ) { tgPolygon result; @@ -523,6 +528,7 @@ tgPolygon tgContour::Intersect( const tgContour& subject, const tgContour& clip return result; } +#endif static bool FindIntermediateNode( const SGGeod& start, const SGGeod& end, const std::vector& nodes, SGGeod& result, diff --git a/src/Lib/terragear/tg_contour.hxx b/src/Lib/terragear/tg_contour.hxx index 88a2bc90..8c6f7656 100644 --- a/src/Lib/terragear/tg_contour.hxx +++ b/src/Lib/terragear/tg_contour.hxx @@ -13,6 +13,9 @@ #include "tg_rectangle.hxx" #include "clipper.hpp" +#define HAS_IS_INSIDE (0) +#define HAS_INTERSECT (0) + /* forward declarations */ class tgPolygon; typedef std::vector tgpolygon_list; @@ -104,9 +107,15 @@ public: static tgPolygon Union( const tgContour& subject, tgPolygon& clip ); static tgPolygon Diff( const tgContour& subject, tgPolygon& clip ); - static tgPolygon Intersect( const tgContour& subject, const tgContour& clip ); +#if HAS_INTERSECT + static tgPolygon Intersect( const tgContour& subject, const tgContour& clip ); +#endif + +#if HAS_IS_INSIDE static bool IsInside( const tgContour& inside, const tgContour& outside ); +#endif + static tgContour AddColinearNodes( const tgContour& subject, UniqueSGGeodSet& nodes ); static tgContour AddColinearNodes( const tgContour& subject, std::vector& nodes ); static bool FindColinearLine( const tgContour& subject, const SGGeod& node, SGGeod& start, SGGeod& end ); diff --git a/src/Prep/OGRDecode/ogr-decode.cxx b/src/Prep/OGRDecode/ogr-decode.cxx index eaa8cef3..1c28eee0 100644 --- a/src/Prep/OGRDecode/ogr-decode.cxx +++ b/src/Prep/OGRDecode/ogr-decode.cxx @@ -44,6 +44,8 @@ /* stretch endpoints to reduce slivers in linear data ~.1 meters */ #define EP_STRETCH (0.1) +#define HAS_SIMPLIFY (0) + using std::string; // scope? @@ -201,7 +203,10 @@ void Decoder::processPolygon(OGRPolygon* poGeometry, const string& area_type ) // first add the outer ring tgPolygon shape = tgShapefile::ToPolygon( poGeometry ); + +#if HAS_SIMPLIFY shape = tgPolygon::Simplify( shape ); +#endif if ( max_segment_length > 0 ) { shape = tgPolygon::SplitLongEdges( shape, max_segment_length );