diff --git a/src/BuildTiles/Clipper/Makefile.am b/src/BuildTiles/Clipper/Makefile.am index 29352cfa..6557c986 100644 --- a/src/BuildTiles/Clipper/Makefile.am +++ b/src/BuildTiles/Clipper/Makefile.am @@ -1,7 +1,7 @@ noinst_LIBRARIES = libClipper.a libClipper_a_SOURCES = clipper.cxx clipper.hxx priorities.cxx priorities.hxx - +dist_pkgdata_DATA = default_priorities.txt noinst_PROGRAMS = testclipper testclipper_SOURCES = testclipper.cxx @@ -19,3 +19,4 @@ testclipper_LDADD = \ -lgenpolyclip -lz INCLUDES = -I$(top_srcdir)/src/Lib -I$(top_srcdir)/src/BuildTiles + diff --git a/src/BuildTiles/Clipper/clipper.cxx b/src/BuildTiles/Clipper/clipper.cxx index 96615061..73a677eb 100644 --- a/src/BuildTiles/Clipper/clipper.cxx +++ b/src/BuildTiles/Clipper/clipper.cxx @@ -73,7 +73,7 @@ bool TGClipper::load_polys(const string& path) { bool poly3d = false; string first_line; string poly_name; - AreaType poly_type = DefaultArea; + AreaType poly_type; int contours, count, i, j; int hole_flag; double startx, starty, startz, x, y, z, lastx, lasty, lastz; @@ -172,17 +172,17 @@ bool TGClipper::load_polys(const string& path) { } } } - + + int area = (int)poly_type; + + add_poly(area, poly); + + // FILE *ofp= fopen("outfile", "w"); + // gpc_write_polygon(ofp, &polys.landuse); + in >> skipcomment; } - int area = (int)poly_type; - - add_poly(area, poly); - - // FILE *ofp= fopen("outfile", "w"); - // gpc_write_polygon(ofp, &polys.landuse); - return true; } @@ -192,7 +192,7 @@ bool TGClipper::load_polys(const string& path) { bool TGClipper::load_osgb36_polys(const string& path) { // cout << "Loading osgb36 poly\n"; string poly_name; - AreaType poly_type = DefaultArea; + AreaType poly_type; int contours, count, i, j; int hole_flag; double startx, starty, x, y, lastx, lasty; @@ -275,26 +275,16 @@ bool TGClipper::load_osgb36_polys(const string& path) { // gpc_add_contour( poly, &v_list, hole_flag ); } + int area = (int)poly_type; + + add_poly( area, poly); + + // FILE *ofp= fopen("outfile", "w"); + // gpc_write_polygon(ofp, &polys.landuse); + in >> skipcomment; } - int area = (int)poly_type; - - // if ( area == OceanArea ) { - // TEST - Ignore - // } else - - if ( area < TG_MAX_AREA_TYPES ) { - polys_in.polys[area].push_back(poly); - } else { - SG_LOG( SG_CLIPPER, SG_ALERT, "Polygon type out of range = " - << (int)poly_type); - exit(-1); - } - - // FILE *ofp= fopen("outfile", "w"); - // gpc_write_polygon(ofp, &polys.landuse); - return true; } @@ -391,7 +381,7 @@ void TGClipper::merge_slivers( TGPolyList& clipped, TGPolygon& slivers ) { for ( area = 0; area < TG_MAX_AREA_TYPES && !done; ++area ) { - if ( area == HoleArea ) { + if ( is_hole_area( area ) ) { // don't merge a non-hole sliver in with a hole continue; } @@ -438,31 +428,6 @@ void TGClipper::merge_slivers( TGPolyList& clipped, TGPolygon& slivers ) { } -static bool -is_water_area (AreaType type) -{ - switch (type) { - case PondArea: - case LakeArea: - case DryLakeArea: - case IntLakeArea: - case ReservoirArea: - case IntReservoirArea: - case StreamArea: - case IntStreamArea: - case CanalArea: - case OceanArea: - case BogArea: - case MarshArea: - case LittoralArea: - case WaterBodyCover: - return true; - default: - return false; - } -} - - // Clip all the polygons against each other in a priority scheme based // on order of the polygon type in the polygon type enum. bool TGClipper::clip_all(const point2d& min, const point2d& max) { @@ -490,33 +455,30 @@ bool TGClipper::clip_all(const point2d& min, const point2d& max) { // best representation of land vs. ocean. If we have other less // accurate data that spills out into the ocean, we want to just // clip it. - TGPolygon land_mask; + // also set up a mask for all water and islands + TGPolygon land_mask, water_mask, island_mask; land_mask.erase(); - for ( i = 0; i < (int)polys_in.polys[DefaultArea].size(); ++i ) { - land_mask = - tgPolygonUnion( land_mask, polys_in.polys[DefaultArea][i] ); - } - - // set up a mask for all water. - TGPolygon water_mask; water_mask.erase(); + island_mask.erase(); for ( i = 0; i < TG_MAX_AREA_TYPES; i++ ) { - if (is_water_area(AreaType(i))) { + if ( is_landmass_area( i ) ) { + for ( unsigned j = 0; j < (int)polys_in.polys[i].size(); ++j ) { + land_mask = + tgPolygonUnion( land_mask, polys_in.polys[i][j] ); + } + } else if ( is_water_area( i ) ) { for (unsigned int j = 0; j < polys_in.polys[i].size(); j++) { water_mask = tgPolygonUnion( water_mask, polys_in.polys[i][j] ); } + } else if ( is_water_area( i ) ) { + for (unsigned int j = 0; j < polys_in.polys[i].size(); j++) { + island_mask = + tgPolygonUnion( island_mask, polys_in.polys[i][j] ); + } } } - // set up island mask, for cutting holes in lakes - TGPolygon island_mask; - island_mask.erase(); - for ( i = 0; i < (int)polys_in.polys[IslandArea].size(); ++i ) { - island_mask = - tgPolygonUnion( island_mask, polys_in.polys[IslandArea][i] ); - } - // process polygons in priority order for ( i = 0; i < TG_MAX_AREA_TYPES; ++i ) { cout << "num polys of type (" << i << ") = " @@ -529,7 +491,7 @@ bool TGClipper::clip_all(const point2d& min, const point2d& max) { tmp = current; // if not a hole, clip the area to the land_mask - if ( i != HoleArea ) { + if ( ! is_hole_area( i ) ) { tmp = tgPolygonInt( tmp, land_mask ); } @@ -549,7 +511,7 @@ bool TGClipper::clip_all(const point2d& min, const point2d& max) { // } // if a water area, cut out potential islands - if ( is_water_area(AreaType(i)) ) { + if ( is_water_area( i ) ) { // clip against island mask tmp = tgPolygonDiff( tmp, island_mask ); } @@ -615,7 +577,7 @@ bool TGClipper::clip_all(const point2d& min, const point2d& max) { } if ( remains.contours() > 0 ) { - polys_clipped.polys[(int)OceanArea].push_back(remains); + polys_clipped.polys[(int)get_sliver_target_area_type()].push_back(remains); } } diff --git a/src/BuildTiles/Clipper/default_priorities.txt b/src/BuildTiles/Clipper/default_priorities.txt new file mode 100644 index 00000000..4bb644a2 --- /dev/null +++ b/src/BuildTiles/Clipper/default_priorities.txt @@ -0,0 +1,68 @@ +Default # Default area which can be overridden by + # raster landcover data (e.g. USGS) +Ocean # collect slivers as ocean + +# Area types in order of descending priority +SomeSort other +Hole hole # Leave area completely empty +Airport other +Freeway road +Road road +Railroad road +Pond lake +Lake lake +DryLake lake +IntLake lake +Reservoir lake +IntReservoir lake +Stream stream +IntStream stream +Canal stream +Glacier other # Solid ice/snow +PackIce other # Water with ice packs +PolarIce other +Ocean ocean +Urban other # Densely-populated city or large town +Town other # Small town or village +FloodLand other # Land subject to flooding +Bog other # Bog +Marsh other # Marshland or swamp +Sand other # Sand-covered area +Littoral other # Tidal, Sand-covered area +Lava other # Lava-covered area + +# USGS Land Covers +# These are low-priority, since known polygons should always win. + +BuiltUpCover other # Urban and Built-Up Land +DryCropPastureCover other # Dryland Cropland and Pasture +IrrCropPastureCover other # Irrigated Cropland and Pasture +MixedCropPastureCover other # Mixed Dryland/Irrigated Cropland and Pasture +CropGrassCover other # Cropland/Grassland Mosaic +CropWoodCover other # Cropland/Woodland Mosaic +GrassCover other # Grassland +ShrubCover other # Shrubland +ShrubGrassCover other # Mixed Shrubland/Grassland +SavannaCover other # Savanna +DeciduousBroadCover other # Deciduous Broadleaf Forest +DeciduousNeedleCover other # Deciduous Needleleaf Forest +EvergreenBroadCover other # Evergreen Broadleaf Forest +EvergreenNeedleCover other # Evergreen Needleleaf Forest +MixedForestCover other # Mixed Forest +WaterBodyCover lake # Water Bodies +HerbWetlandCover other # Herbaceous Wetland +WoodedWetlandCover other # Wooded Wetland +BarrenCover other # Barren or Sparsely Vegetated +HerbTundraCover other # Herbaceous Tundra +WoodedTundraCover other # Wooded Tundra +MixedTundraCover other # Mixed Tundra +BareTundraCover other # Bare Ground Tundra +SnowCover other # Snow or Ice + +Island island # any island area not covered otherwise +Default landmass # any land area not covered otherwise + +Void other +Null other +Unknown other + diff --git a/src/BuildTiles/Clipper/priorities.cxx b/src/BuildTiles/Clipper/priorities.cxx index 445c3943..2face93c 100644 --- a/src/BuildTiles/Clipper/priorities.cxx +++ b/src/BuildTiles/Clipper/priorities.cxx @@ -17,133 +17,165 @@ // You should have received a copy of the GNU General Public License // along with this program; if not, write to the Free Software // Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. - + #include #include +#include +#include +#include #include #include +#include + #include "priorities.hxx" +using std::fstream; using std::string; using std::map; +using std::vector; -typedef map area_type_map; +typedef enum { + Hole, + Landmass, + Island, + Ocean, + Lake, + Stream, + Road, + Other +} AreaKind; + +typedef struct { + string name; + AreaKind kind; +} area_type_descriptor; + +typedef vector area_type_list; typedef map area_name_map; -static area_type_map area_types; +static area_type_list area_types; static area_name_map area_names; +static AreaType default_area_type; +static AreaType sliver_target_area_type; - -inline static void set_area (const string &name, AreaType type) -{ - area_types[type] = name; - area_names[name] = type; +int load_area_types( const std::string& filename ) { + fstream in ( filename.c_str() ); + + if ( ! in ) { + SG_LOG(SG_GENERAL, SG_ALERT, "Unable to open file " << filename); + return 0; + } + + in >> skipcomment; + string sliver_area_name, default_area_name; + + in >> default_area_name; + in >> skipcomment; + in >> sliver_area_name; + in >> skipcomment; + while ( !in.eof() ) { + area_type_descriptor descriptor; + string type; + descriptor.kind = Other; + in >> descriptor.name; + in >> type; + if ( type=="hole" ) { + descriptor.kind = Hole; + } else if ( type=="landmass" ) { + descriptor.kind = Landmass; + } else if ( type=="island" ) { + descriptor.kind = Island; + } else if ( type=="ocean" ) { + descriptor.kind = Ocean; + } else if ( type=="lake" ) { + descriptor.kind = Lake; + } else if ( type=="stream" ) { + descriptor.kind = Stream; + } else if ( type=="road" ) { + descriptor.kind = Road; + } else if ( type=="other" ) { + descriptor.kind = Other; + } + AreaType index = (AreaType)area_types.size(); + area_types.push_back(descriptor); + area_names[descriptor.name]=index; + SG_LOG(SG_GENERAL, SG_INFO, " " << descriptor.name << " " << descriptor.kind); + in >> skipcomment; + } + + in.close(); + + sliver_target_area_type = get_area_type( sliver_area_name ); + default_area_type = get_area_type( default_area_name ); + + return 1; } - - -static bool _initialized = false; - - -inline static void init () -{ - if (_initialized) - return; - - set_area("SomeSort", SomeSortOfArea); - set_area("Hole", HoleArea); - set_area("Airport", AirportArea); - set_area("Island", IslandArea); - set_area("Pond", PondArea); - set_area("Swamp or Marsh", MarshArea); - set_area("Marsh", MarshArea); - set_area("Littoral", LittoralArea); - set_area("Bog", BogArea); - set_area("Sand", SandArea); - set_area("Lava", LavaArea); - set_area("FloodLand", FloodLandArea); - set_area("Lake", LakeArea); - set_area("Lake Dry", DryLakeArea); - set_area("DryLake", DryLakeArea); - set_area("Lake Intermittent", IntLakeArea); - set_area("IntermittentLake", IntLakeArea); - set_area("Reservoir", ReservoirArea); - set_area("Reservoir Intermittent", IntReservoirArea); - set_area("IntermittentReservoir", IntReservoirArea); - set_area("Freeway", FreewayArea); - set_area("Road", RoadArea); - set_area("Railroad", RailroadArea); - set_area("Stream", StreamArea); - set_area("IntermittentStream", IntStreamArea); - set_area("Canal", CanalArea); - set_area("Glacier", GlacierArea); - set_area("PackIce", PackIceArea); - set_area("PolarIce", PolarIceArea); - set_area("Urban", UrbanArea); - set_area("Town", TownArea); - set_area("BuiltUpCover", BuiltUpCover); - set_area("DryCropPastureCover", DryCropPastureCover); - set_area("IrrCropPastureCover", IrrCropPastureCover); - set_area("MixedCropPastureCover", MixedCropPastureCover); - set_area("CropGrassCover", CropGrassCover); - set_area("CropWoodCover", CropWoodCover); - set_area("GrassCover", GrassCover); - set_area("ShrubCover", ShrubCover); - set_area("ShrubGrassCover", ShrubGrassCover); - set_area("SavannaCover", SavannaCover); - set_area("DeciduousBroadCover", DeciduousBroadCover); - set_area("DeciduousNeedleCover", DeciduousNeedleCover); - set_area("EvergreenBroadCover", EvergreenBroadCover); - set_area("EvergreenNeedleCover", EvergreenNeedleCover); - set_area("MixedForestCover", MixedForestCover); - set_area("WaterBodyCover", WaterBodyCover); - set_area("HerbWetlandCover", HerbWetlandCover); - set_area("WoodedWetlandCover", WoodedWetlandCover); - set_area("BarrenCover", BarrenCover); - set_area("HerbTundraCover", HerbTundraCover); - set_area("WoodedTundraCover", WoodedTundraCover); - set_area("MixedTundraCover", MixedTundraCover); - set_area("BareTundraCover", BareTundraCover); - set_area("SnowCover", SnowCover); - set_area("Default", DefaultArea); - set_area("Bay Estuary or Ocean", OceanArea); - set_area("Ocean", OceanArea); - set_area("Void Area", VoidArea); - set_area("Null", NullArea); - - _initialized = true; -} - - // return area type from text name AreaType get_area_type (const string &area) { - init(); area_name_map::const_iterator it = area_names.find(area); if (it != area_names.end()) { return it->second; } else { - SG_LOG(SG_GENERAL, SG_WARN, "unknown area = '" << area << "'"); - // SG_LOG(SG_GENERAL, SG_DEBUG, "area = " << area); - // for ( int i = 0; i < area.length(); i++ ) { - // SG_LOG(SG_GENERAL, SG_DEBUG, i << ") " << (int)area[i]); - // } - return UnknownArea; + SG_LOG(SG_GENERAL, SG_ALERT, "unknown area = '" << area << "'"); + exit(-1); } } +static area_type_descriptor& get_area_descriptor( AreaType area ) { + if ( 0<=area || area < area_types.size() ) { + return area_types[area]; + } else { + SG_LOG(SG_GENERAL, SG_ALERT, "unknown area code = " << (int)area); + exit(-1); + } +} + // return text from of area name string get_area_name( AreaType area ) { - init(); - area_type_map::const_iterator it = area_types.find(area); - if (it != area_types.end()) { - return it->second; - } else { - SG_LOG(SG_GENERAL, SG_WARN, "unknown area code = " << (int)area); - return "Unknown"; - } + return get_area_descriptor( area ).name; } +bool is_hole_area( AreaType area ) { + return get_area_descriptor( area ).kind==Hole; +} +bool is_water_area( AreaType area ) { + const AreaKind kind = get_area_descriptor( area ).kind; + return (kind==Ocean || kind==Lake || kind==Stream); +} + +bool is_landmass_area( AreaType area ) { + return get_area_descriptor( area ).kind==Landmass; +} + +bool is_lake_area( AreaType area ) { + const AreaKind kind = get_area_descriptor( area ).kind; + return (kind==Lake); +} + +bool is_stream_area( AreaType area ) { + const AreaKind kind = get_area_descriptor( area ).kind; + return (kind==Stream); +} + +bool is_road_area( AreaType area ) { + const AreaKind kind = get_area_descriptor( area ).kind; + return (kind==Road); +} + +bool is_ocean_area( AreaType area ) { + const AreaKind kind = get_area_descriptor( area ).kind; + return (kind==Ocean); +} + +AreaType get_sliver_target_area_type() { + return sliver_target_area_type; +} + +AreaType get_default_area_type() { + return default_area_type; +} diff --git a/src/BuildTiles/Clipper/priorities.hxx b/src/BuildTiles/Clipper/priorities.hxx index d57649b2..7c8c570f 100644 --- a/src/BuildTiles/Clipper/priorities.hxx +++ b/src/BuildTiles/Clipper/priorities.hxx @@ -29,74 +29,17 @@ #include -// Posible shape file types. Note the order of these is important and -// defines the priority of these shapes if they should intersect. The -// smaller the number, the higher the priority. -enum AreaType { - SomeSortOfArea = 0, - HoleArea, // Leave area completely empty - AirportArea, - FreewayArea, - RoadArea, - RailroadArea, - PondArea, - LakeArea, - DryLakeArea, - IntLakeArea, - ReservoirArea, - IntReservoirArea, - StreamArea, - IntStreamArea, - CanalArea, - GlacierArea, // Solid ice/snow - PackIceArea, // Water with ice packs - PolarIceArea, - OceanArea, - UrbanArea, // Densely-populated city or large town - TownArea, // Small town or village - FloodLandArea, // Land subject to flooding - BogArea, // Bog - MarshArea, // Marshland or swamp - SandArea, // Sand-covered area - LittoralArea, // Tidal, Sand-covered area - LavaArea, // Lava-covered area - - // USGS Land Covers - // These are low-priority, since known polygons should always win. - - BuiltUpCover, // Urban and Built-Up Land - DryCropPastureCover, // Dryland Cropland and Pasture - IrrCropPastureCover, // Irrigated Cropland and Pasture - MixedCropPastureCover, // Mixed Dryland/Irrigated Cropland and Pasture - CropGrassCover, // Cropland/Grassland Mosaic - CropWoodCover, // Cropland/Woodland Mosaic - GrassCover, // Grassland - ShrubCover, // Shrubland - ShrubGrassCover, // Mixed Shrubland/Grassland - SavannaCover, // Savanna - DeciduousBroadCover, // Deciduous Broadleaf Forest - DeciduousNeedleCover, // Deciduous Needleleaf Forest - EvergreenBroadCover, // Evergreen Broadleaf Forest - EvergreenNeedleCover, // Evergreen Needleleaf Forest - MixedForestCover, // Mixed Forest - WaterBodyCover, // Water Bodies - HerbWetlandCover, // Herbaceous Wetland - WoodedWetlandCover, // Wooded Wetland - BarrenCover, // Barren or Sparsely Vegetated - HerbTundraCover, // Herbaceous Tundra - WoodedTundraCover, // Wooded Tundra - MixedTundraCover, // Mixed Tundra - BareTundraCover, // Bare Ground Tundra - SnowCover, // Snow or Ice - - IslandArea, // any island area not covered otherwise - DefaultArea, // any land area not covered otherwise - - VoidArea, - NullArea, - UnknownArea -}; - +typedef unsigned int AreaType; +int load_area_types( const std::string& filename ); +bool is_hole_area(AreaType areaType); +bool is_landmass_area(AreaType areaType); +bool is_water_area(AreaType areaType); +bool is_lake_area(AreaType areaType); +bool is_stream_area(AreaType areaType); +bool is_road_area(AreaType areaType); +bool is_ocean_area(AreaType areaType); +AreaType get_sliver_target_area_type(); +AreaType get_default_area_type(); // return area type from text name AreaType get_area_type( const std::string &area ); diff --git a/src/BuildTiles/Main/Makefile.am b/src/BuildTiles/Main/Makefile.am index 24c7ba1a..04a3ab3e 100644 --- a/src/BuildTiles/Main/Makefile.am +++ b/src/BuildTiles/Main/Makefile.am @@ -26,9 +26,11 @@ fgfs_construct_LDADD = \ fgfs_master_SOURCES = master.cxx fgfs_master_LDADD = -lsgbucket -lsgmisc -lsgdebug -lsgxml -AM_CPPFLAGS = -DDEFAULT_USGS_MAPFILE="\"$(pkgdatadir)/usgsmap.txt\"" +AM_CPPFLAGS = -DDEFAULT_USGS_MAPFILE="\"$(pkgdatadir)/usgsmap.txt\"" \ + -DDEFAULT_PRIORITIES_FILE="\"$(pkgdatadir)/default_priorities.txt\"" INCLUDES = \ -I$(top_srcdir)/src \ -I$(top_srcdir)/src/Lib \ -I$(top_srcdir)/src/BuildTiles + diff --git a/src/BuildTiles/Main/main.cxx b/src/BuildTiles/Main/main.cxx index adb772d9..c3ac5ade 100644 --- a/src/BuildTiles/Main/main.cxx +++ b/src/BuildTiles/Main/main.cxx @@ -154,7 +154,7 @@ static AreaType get_area_type (const LandCover &cover, int cover_value = cover.getValue(xpos, ypos); AreaType area = translateUSGSCover(cover_value); - if (area != DefaultArea) { + if ( area != get_default_area_type() ) { // Non-default area is fine. return area; } else { @@ -165,7 +165,7 @@ static AreaType get_area_type (const LandCover &cover, if (x != xpos || y != ypos) { cover_value = cover.getValue(x, y); area = translateUSGSCover(cover_value); - if (area != DefaultArea) { + if (area != get_default_area_type() ) { return area; } } @@ -174,7 +174,7 @@ static AreaType get_area_type (const LandCover &cover, } // OK, give up and return default - return DefaultArea; + return get_default_area_type(); } @@ -232,7 +232,7 @@ static void make_area( const LandCover &cover, const TGArray &array, AreaType area = get_area_type( cover, x1 + half_dx, y1 + half_dy, x2 - x1, y2 - y1 ); - if (area != DefaultArea) { + if ( area != get_default_area_type() ) { // Create a square polygon and merge it into the list. TGPolygon poly; poly.erase(); @@ -512,7 +512,7 @@ static void fix_point_heights( TGConstruct& c, const TGArray& array ) // It might be better to eventually iterate, and allow // some flexibility in elevations to handle rivers and // things like that. - if ( (a == LakeArea) || (a == ReservoirArea) ) { + if ( is_lake_area( a ) ) { e1 = raw_nodes[n1].z(); e2 = raw_nodes[n2].z(); e3 = raw_nodes[n3].z(); @@ -525,8 +525,7 @@ static void fix_point_heights( TGConstruct& c, const TGArray& array ) raw_nodes[n1].setz( min ); raw_nodes[n2].setz( min ); raw_nodes[n3].setz( min ); - } else if ( (a == StreamArea) || (a == IntStreamArea) - || (a == CanalArea) ) { + } else if ( is_stream_area( a ) ) { e1 = raw_nodes[n1].z(); e2 = raw_nodes[n2].z(); e3 = raw_nodes[n3].z(); @@ -546,8 +545,7 @@ static void fix_point_heights( TGConstruct& c, const TGArray& array ) if ( max1 < e1 ) { raw_nodes[n1].setz( max1 ); } if ( max2 < e2 ) { raw_nodes[n2].setz( max2 ); } if ( max3 < e3 ) { raw_nodes[n3].setz( max3 ); } - } else if ( (a == RoadArea) || (a == FreewayArea) - || (a == RailroadArea) ) { + } else if ( is_road_area( a ) ) { e1 = raw_nodes[n1].z(); e2 = raw_nodes[n2].z(); e3 = raw_nodes[n3].z(); @@ -581,7 +579,7 @@ static void fix_point_heights( TGConstruct& c, const TGArray& array ) n3 = t.get_n3(); a = (AreaType)((int)(t.get_attribute())); - if ( a == OceanArea ) { + if ( is_ocean_area( a ) ) { raw_nodes[n1].setz( 0.0 ); raw_nodes[n2].setz( 0.0 ); raw_nodes[n3].setz( 0.0 ); @@ -611,7 +609,7 @@ static void fix_land_cover_assignments( TGConstruct& c ) { cout << " Total triangles = " << tri_elements.size() << endl; for ( unsigned int i = 0; i < tri_elements.size(); ++i ) { TGTriEle t = tri_elements[i]; - if ( t.get_attribute() == DefaultArea ) { + if ( t.get_attribute() == get_default_area_type() ) { Point3D p1 = geod_nodes[t.get_n1()]; Point3D p2 = geod_nodes[t.get_n2()]; Point3D p3 = geod_nodes[t.get_n3()]; @@ -1043,6 +1041,7 @@ static void usage( const string name ) { cout << " --xdist=" << endl; cout << " --ydist=" << endl; cout << " --nudge=" << endl; + cout << " --priorities=" << endl; cout << " --usgs-map=" << endl; cout << " --useUKgrid" << endl; cout << " --no-write-shared-edges" << endl; @@ -1056,6 +1055,7 @@ int main(int argc, char **argv) { string output_dir = "."; string work_dir = "."; string cover = ""; + string priorities_file = DEFAULT_PRIORITIES_FILE; string usgs_map_file = DEFAULT_USGS_MAPFILE; double lon = -110.664244; // P13 double lat = 33.352890; @@ -1102,6 +1102,8 @@ int main(int argc, char **argv) { nudge = atof(arg.substr(8).c_str())*SG_EPSILON; } else if (arg.find("--cover=") == 0) { cover = arg.substr(8); + } else if (arg.find("--priorities=") == 0) { + priorities_file = arg.substr(13); } else if (arg.find("--usgs-map=") == 0) { usgs_map_file = arg.substr(11); } else if (arg.find("--useUKgrid") == 0) { @@ -1132,6 +1134,11 @@ int main(int argc, char **argv) { load_dirs.push_back(argv[i]); cout << "Load directory: " << argv[i] << endl; } + cout << "Priorities file is " << priorities_file << endl; + if ( ! load_area_types( priorities_file ) ) { + SG_LOG(SG_GENERAL, SG_ALERT, "Failed to load priorities file " << priorities_file); + exit(-1); + } cout << "USGS Map file is " << usgs_map_file << endl; if ( ! load_usgs_map( usgs_map_file ) ) { SG_LOG(SG_GENERAL, SG_ALERT, "Failed to load USGS map file " << usgs_map_file); diff --git a/src/BuildTiles/Main/usgs.cxx b/src/BuildTiles/Main/usgs.cxx index 92665e27..fd464af1 100644 --- a/src/BuildTiles/Main/usgs.cxx +++ b/src/BuildTiles/Main/usgs.cxx @@ -59,9 +59,9 @@ int load_usgs_map( const std::string& filename ) { AreaType translateUSGSCover (int usgs_value) { if ( 0contour[j].vertex[0].x, - // gpc_poly->contour[j].vertex[0].y); - // fclose(junkfp); - } - - /* if ( i == OceanArea ) { - cout << "temporary exit point" << endl; - exit(-1); - } */ - - // for each contour, calculate a point inside (but not - // also inside any interior contours - - // new way - - // try to make sure our polygons aren't goofy -#if 0 + sprintf(junkn, "c%d", j); + gpc_poly.write_contour( j, junkn ); + */ + + for ( int k = 0; k < gpc_poly.contour_size( j ); k++ ) { + Point3D p = gpc_poly.get_pt( j, k ); + index = in_nodes.unique_add( p ); + // junkp = in_nodes.get_node( index ); + // fprintf(junkfp, "%.4f %.4f\n", junkp.x(), junkp.y()); + // cout << " - " << index << endl; + } + // fprintf(junkfp, "%.4f %.4f\n", + // gpc_poly->contour[j].vertex[0].x, + // gpc_poly->contour[j].vertex[0].y); + // fclose(junkfp); + } + + /* if ( i == OceanArea ) { + cout << "temporary exit point" << endl; + exit(-1); + } */ + + // for each contour, calculate a point inside (but not + // also inside any interior contours + + // new way + + // try to make sure our polygons aren't goofy + #if 0 // CLO 09/18/2001: if we snap polygons including holes // will this screw up the edge matching when objects are // inserted into their holes? - gpc_poly = snap(gpc_poly, 0.000001); -#endif - gpc_poly = remove_dups( gpc_poly ); - gpc_poly = reduce_degeneracy( gpc_poly ); - gpc_poly = reduce_degeneracy( gpc_poly ); // can happen multiple time - gpc_poly = remove_dups( gpc_poly ); - gpc_poly = remove_bad_contours( gpc_poly ); - gpc_poly = remove_cycles( gpc_poly ); - - cout << "after sanity checks, contours = " - << gpc_poly.contours() << endl; - - /* + gpc_poly = snap(gpc_poly, 0.000001); + #endif + gpc_poly = remove_dups( gpc_poly ); + gpc_poly = reduce_degeneracy( gpc_poly ); + gpc_poly = reduce_degeneracy( gpc_poly ); // can happen multiple time + gpc_poly = remove_dups( gpc_poly ); + gpc_poly = remove_bad_contours( gpc_poly ); + gpc_poly = remove_cycles( gpc_poly ); + + cout << "after sanity checks, contours = " + << gpc_poly.contours() << endl; + + /* for ( j = 0; j < gpc_poly.contours(); ++j ) { - cout << " contour " << j << " size = " - << gpc_poly.contour_size( j ) << endl; - char junkn[256]; - sprintf(junkn, "d%d", j); - gpc_poly.write_contour( j, junkn ); + cout << " contour " << j << " size = " + << gpc_poly.contour_size( j ) << endl; + char junkn[256]; + sprintf(junkn, "d%d", j); + gpc_poly.write_contour( j, junkn ); } - */ - - cout << "before calc_points_inside()" << endl; - calc_points_inside( gpc_poly ); - cout << "after calc_points_inside()" << endl; - -#if 0 - // old way - Point3D inside_pt; - for ( j = 0; j < gpc_poly.contours(); ++j ) { - inside_pt = calc_point_inside( gpc_poly, j, in_nodes ); - gpc_poly.set_point_inside( j, inside_pt ); - } -#endif - - polylist[i].push_back( gpc_poly ); - -#if 0 - // temporary ... write out hole/polygon info for debugging - for ( j = 0; j < (int)gpc_poly.contours(); ++j ) { - char pname[256]; - sprintf(pname, "poly%02d-%02d-%02d", i, debug_counter, j); - cout << "writing to " << pname << endl; - FILE *fp = fopen( pname, "w" ); - Point3D point; - for ( int k = 0; k < gpc_poly.contour_size( j ); ++k ) { - point = gpc_poly.get_pt( j, k ); - fprintf( fp, "%.6f %.6f\n", point.x(), point.y() ); - } - point = gpc_poly.get_pt( j, 0 ); - fprintf( fp, "%.6f %.6f\n", point.x(), point.y() ); - fclose(fp); - - char hname[256]; - sprintf(hname, "hole%02d-%02d-%02d", i, debug_counter, j); - FILE *fh = fopen( hname, "w" ); - point = gpc_poly.get_point_inside( j ); - fprintf( fh, "%.6f %.6f\n", point.x(), point.y() ); - fclose(fh); - } - - // cout << "type a letter + enter to continue: "; - // string input; - // cin >> input; -#endif - - ++debug_counter; - } + */ + + cout << "before calc_points_inside()" << endl; + calc_points_inside( gpc_poly ); + cout << "after calc_points_inside()" << endl; + + #if 0 + // old way + Point3D inside_pt; + for ( j = 0; j < gpc_poly.contours(); ++j ) { + inside_pt = calc_point_inside( gpc_poly, j, in_nodes ); + gpc_poly.set_point_inside( j, inside_pt ); + } + #endif + + polylist[i].push_back( gpc_poly ); + + #if 0 + // temporary ... write out hole/polygon info for debugging + for ( j = 0; j < (int)gpc_poly.contours(); ++j ) { + char pname[256]; + sprintf(pname, "poly%02d-%02d-%02d", i, debug_counter, j); + cout << "writing to " << pname << endl; + FILE *fp = fopen( pname, "w" ); + Point3D point; + for ( int k = 0; k < gpc_poly.contour_size( j ); ++k ) { + point = gpc_poly.get_pt( j, k ); + fprintf( fp, "%.6f %.6f\n", point.x(), point.y() ); + } + point = gpc_poly.get_pt( j, 0 ); + fprintf( fp, "%.6f %.6f\n", point.x(), point.y() ); + fclose(fp); + + char hname[256]; + sprintf(hname, "hole%02d-%02d-%02d", i, debug_counter, j); + FILE *fh = fopen( hname, "w" ); + point = gpc_poly.get_point_inside( j ); + fprintf( fh, "%.6f %.6f\n", point.x(), point.y() ); + fclose(fh); + } + + // cout << "type a letter + enter to continue: "; + // string input; + // cin >> input; + #endif + + ++debug_counter; + } } - + // last, do the rest of the height nodes for ( i = 0; i < (int)fit_list.size(); ++i ) { - index = in_nodes.course_add( fit_list[i] ); + index = in_nodes.course_add( fit_list[i] ); } - + for ( i = 0; i < TG_MAX_AREA_TYPES; ++i ) { - if ( polylist[i].size() ) { - cout << get_area_name((AreaType)i) << " = " - << polylist[i].size() << endl; - } + if ( polylist[i].size() ) { + cout << get_area_name((AreaType)i) << " = " + << polylist[i].size() << endl; + } } - + // traverse the polygon lists and build the segment (edge) list // that is used by the "Triangle" lib. - + cout << "building segment list" << endl; int i1, i2; Point3D p1, p2; point_list node_list = in_nodes.get_node_list(); TGPolygon poly; - + for ( i = 0; i < TG_MAX_AREA_TYPES; ++i ) { - cout << "area type = " << i << endl; - poly_list_iterator tp_current, tp_last; - tp_current = polylist[i].begin(); - tp_last = polylist[i].end(); - - // process each polygon in list - for ( ; tp_current != tp_last; ++tp_current ) { - poly = *tp_current; - cout << " processing a polygon with contours = " - << poly.contours() << endl; - for ( int j = 0; j < (int)poly.contours(); ++j) { - for ( int k = 0; k < (int)(poly.contour_size(j) - 1); ++k ) { - p1 = poly.get_pt( j, k ); - p2 = poly.get_pt( j, k + 1 ); - i1 = in_nodes.find( p1 ); - i2 = in_nodes.find( p2 ); - // calc_line_params(i1, i2, &m, &b); - if ( i == (int)HoleArea ) { - // mark as a boundary - in_segs.unique_divide_and_add( node_list, - TGTriSeg(i1, i2, 1) ); - } else { - // non boundary - in_segs.unique_divide_and_add( node_list, - TGTriSeg(i1, i2, 0) ); - } - } - p1 = poly.get_pt( j, 0 ); - p2 = poly.get_pt( j, poly.contour_size(j) - 1 ); - i1 = in_nodes.find( p1 ); - i2 = in_nodes.find( p2 ); - // calc_line_params(i1, i2, &m, &b); - if ( i == (int)HoleArea ) { - // mark as a boundary - in_segs.unique_divide_and_add( node_list, - TGTriSeg(i1, i2, 1) ); - } else { - // non boundary - in_segs.unique_divide_and_add( node_list, - TGTriSeg(i1, i2, 0) ); - } - } - } + cout << "area type = " << i << endl; + poly_list_iterator tp_current, tp_last; + tp_current = polylist[i].begin(); + tp_last = polylist[i].end(); + + // process each polygon in list + for ( ; tp_current != tp_last; ++tp_current ) { + poly = *tp_current; + cout << " processing a polygon with contours = " + << poly.contours() << endl; + for ( int j = 0; j < (int)poly.contours(); ++j) { + for ( int k = 0; k < (int)(poly.contour_size(j) - 1); ++k ) { + p1 = poly.get_pt( j, k ); + p2 = poly.get_pt( j, k + 1 ); + i1 = in_nodes.find( p1 ); + i2 = in_nodes.find( p2 ); + // calc_line_params(i1, i2, &m, &b); + if ( is_hole_area(i) ) { + // mark as a boundary + in_segs.unique_divide_and_add( node_list, + TGTriSeg(i1, i2, 1) ); + } else { + // non boundary + in_segs.unique_divide_and_add( node_list, + TGTriSeg(i1, i2, 0) ); + } + } + p1 = poly.get_pt( j, 0 ); + p2 = poly.get_pt( j, poly.contour_size(j) - 1 ); + i1 = in_nodes.find( p1 ); + i2 = in_nodes.find( p2 ); + // calc_line_params(i1, i2, &m, &b); + if ( is_hole_area(i) ) { + // mark as a boundary + in_segs.unique_divide_and_add( node_list, + TGTriSeg(i1, i2, 1) ); + } else { + // non boundary + in_segs.unique_divide_and_add( node_list, + TGTriSeg(i1, i2, 0) ); + } + } + } } - + return 0; } @@ -282,10 +282,10 @@ TGTriangle::build( const point_list& corner_list, int TGTriangle::rebuild( TGConstruct& c ) { in_nodes.clear(); in_segs.clear(); - + in_nodes = c.get_tri_nodes(); in_segs = c.get_tri_segs(); - + return 0; } @@ -303,118 +303,126 @@ int TGTriangle::run_triangulate( double angle, const int pass ) { struct triangulateio in, out, vorout; int counter; int i, j; - + // point list point_list node_list = in_nodes.get_node_list(); in.numberofpoints = node_list.size(); in.pointlist = (REAL *) malloc(in.numberofpoints * 2 * sizeof(REAL)); - + for ( i = 0; i < in.numberofpoints; ++i ) { - in.pointlist[2*i] = node_list[i].x(); - in.pointlist[2*i + 1] = node_list[i].y(); + in.pointlist[2*i] = node_list[i].x(); + in.pointlist[2*i + 1] = node_list[i].y(); } - + in.numberofpointattributes = 1; in.pointattributelist = (REAL *) malloc(in.numberofpoints * - in.numberofpointattributes * - sizeof(REAL)); + in.numberofpointattributes * + sizeof(REAL)); for ( i = 0; i < in.numberofpoints * in.numberofpointattributes; ++i) { - in.pointattributelist[i] = node_list[i].z(); + in.pointattributelist[i] = node_list[i].z(); } - + in.pointmarkerlist = (int *) malloc(in.numberofpoints * sizeof(int)); for ( i = 0; i < in.numberofpoints; ++i) { - in.pointmarkerlist[i] = 0; + in.pointmarkerlist[i] = 0; } - + // triangle list in.numberoftriangles = 0; - + // segment list triseg_list seg_list = in_segs.get_seg_list(); in.numberofsegments = seg_list.size(); in.segmentlist = (int *) malloc(in.numberofsegments * 2 * sizeof(int)); in.segmentmarkerlist = (int *) malloc(in.numberofsegments * sizeof(int)); - + triseg_list_iterator s_current, s_last; s_current = seg_list.begin(); s_last = seg_list.end(); counter = 0; for ( ; s_current != s_last; ++s_current ) { - in.segmentlist[counter++] = s_current->get_n1(); - in.segmentlist[counter++] = s_current->get_n2(); + in.segmentlist[counter++] = s_current->get_n1(); + in.segmentlist[counter++] = s_current->get_n2(); } s_current = seg_list.begin(); s_last = seg_list.end(); counter = 0; for ( ; s_current != s_last; ++s_current ) { - in.segmentmarkerlist[counter++] = s_current->get_boundary_marker(); + in.segmentmarkerlist[counter++] = s_current->get_boundary_marker(); } - + // hole list (make holes for airport ignore areas) poly_list_iterator h_current, h_last; in.numberofholes = 0; - h_current = polylist[(int)HoleArea].begin(); - h_last = polylist[(int)HoleArea].end(); - for ( ; h_current != h_last; ++h_current ) { - poly = *h_current; - for ( j = 0; j < poly.contours(); ++j ) { - in.numberofholes++; + for ( i = 0; i < TG_MAX_AREA_TYPES; i++) { + if ( is_hole_area( i ) ) { + h_current = polylist[i].begin(); + h_last = polylist[i].end(); + for ( ; h_current != h_last; ++h_current ) { + poly = *h_current; + for ( j = 0; j < poly.contours(); ++j ) { + in.numberofholes++; + } + } } } - + in.holelist = (REAL *) malloc(in.numberofholes * 2 * sizeof(REAL)); - - h_current = polylist[(int)HoleArea].begin(); - h_last = polylist[(int)HoleArea].end(); - counter = 0; - for ( ; h_current != h_last; ++h_current ) { - poly = *h_current; - for ( j = 0; j < poly.contours(); ++j ) { - p = poly.get_point_inside( j ); - in.holelist[counter++] = p.x(); - in.holelist[counter++] = p.y(); - } + + for ( i = 0; i < TG_MAX_AREA_TYPES; i++) { + if ( is_hole_area( i ) ) { + h_current = polylist[i].begin(); + h_last = polylist[i].end(); + counter = 0; + for ( ; h_current != h_last; ++h_current ) { + poly = *h_current; + for ( j = 0; j < poly.contours(); ++j ) { + p = poly.get_point_inside( j ); + in.holelist[counter++] = p.x(); + in.holelist[counter++] = p.y(); + } + } + } } - + // region list in.numberofregions = 0; for ( i = 0; i < TG_MAX_AREA_TYPES; ++i ) { - poly_list_iterator h_current, h_last; - h_current = polylist[i].begin(); - h_last = polylist[i].end(); - for ( ; h_current != h_last; ++h_current ) { - poly = *h_current; - for ( j = 0; j < poly.contours(); ++j ) { - if ( ! poly.get_hole_flag( j ) ) { - ++in.numberofregions; - } - } - } + poly_list_iterator h_current, h_last; + h_current = polylist[i].begin(); + h_last = polylist[i].end(); + for ( ; h_current != h_last; ++h_current ) { + poly = *h_current; + for ( j = 0; j < poly.contours(); ++j ) { + if ( ! poly.get_hole_flag( j ) ) { + ++in.numberofregions; + } + } + } } - + in.regionlist = (REAL *) malloc(in.numberofregions * 4 * sizeof(REAL)); counter = 0; for ( i = 0; i < TG_MAX_AREA_TYPES; ++i ) { - poly_list_iterator h_current, h_last; - h_current = polylist[(int)i].begin(); - h_last = polylist[(int)i].end(); - for ( ; h_current != h_last; ++h_current ) { - poly = *h_current; - for ( j = 0; j < poly.contours(); ++j ) { - if ( ! poly.get_hole_flag( j ) ) { - p = poly.get_point_inside( j ); - cout << "Region point = " << p << endl; - in.regionlist[counter++] = p.x(); // x coord - in.regionlist[counter++] = p.y(); // y coord - in.regionlist[counter++] = i; // region attribute - in.regionlist[counter++] = -1.0; // area constraint - // (unused) - } - } - } + poly_list_iterator h_current, h_last; + h_current = polylist[(int)i].begin(); + h_last = polylist[(int)i].end(); + for ( ; h_current != h_last; ++h_current ) { + poly = *h_current; + for ( j = 0; j < poly.contours(); ++j ) { + if ( ! poly.get_hole_flag( j ) ) { + p = poly.get_point_inside( j ); + cout << "Region point = " << p << endl; + in.regionlist[counter++] = p.x(); // x coord + in.regionlist[counter++] = p.y(); // y coord + in.regionlist[counter++] = i; // region attribute + in.regionlist[counter++] = -1.0; // area constraint + // (unused) + } + } + } } - + // prep the output structures out.pointlist = (REAL *) NULL; // Not needed if -N switch used. // Not needed if -N switch used or number of point attributes is zero: @@ -430,7 +438,7 @@ int TGTriangle::run_triangulate( double angle, const int pass ) { out.segmentmarkerlist = (int *) NULL; out.edgelist = (int *) NULL; // Needed only if -e switch used. out.edgemarkerlist = (int *) NULL; // Needed if -e used and -B not used. - + vorout.pointlist = (REAL *) NULL; // Needed only if -v switch used. // Needed only if -v switch used and number of attributes is not zero: vorout.pointattributelist = (REAL *) NULL; @@ -439,81 +447,81 @@ int TGTriangle::run_triangulate( double angle, const int pass ) { // TEMPORARY // write_tri_data(&in); exit(1); - + // Triangulate the points. Switches are chosen to read and write // a PSLG (p), preserve the convex hull (c), number everything // from zero (z), assign a regional attribute to each element (A), // and produce an edge list (e), and a triangle neighbor list (n). - + string tri_options; if ( pass == 1 ) { - // use a quality value of 10 (q10) meaning no interior - // triangle angles less than 10 degrees - // tri_options = "pczAen"; - if ( angle < 0.00001 ) { - tri_options = "pczAen"; - } else { + // use a quality value of 10 (q10) meaning no interior + // triangle angles less than 10 degrees + // tri_options = "pczAen"; + if ( angle < 0.00001 ) { + tri_options = "pczAen"; + } else { char angle_str[256]; sprintf( angle_str, "%.2f", angle ); - tri_options = "pczq"; + tri_options = "pczq"; tri_options += angle_str; tri_options += "Aen"; - } - // // string tri_options = "pzAen"; - // // string tri_options = "pczq15S400Aen"; + } + // // string tri_options = "pzAen"; + // // string tri_options = "pczq15S400Aen"; } else if ( pass == 2 ) { - // no new points on boundary (Y), no internal segment - // splitting (YY), no quality refinement (q) - tri_options = "pczYYAen"; + // no new points on boundary (Y), no internal segment + // splitting (YY), no quality refinement (q) + tri_options = "pczYYAen"; } else { - cout << "unknown pass number = " << pass - << " in TGTriangle::run_triangulate()" << endl; - exit(-1); + cout << "unknown pass number = " << pass + << " in TGTriangle::run_triangulate()" << endl; + exit(-1); } cout << "Triangulation with options = " << tri_options << endl; - + triangulate( (char *)tri_options.c_str(), &in, &out, &vorout ); - + // TEMPORARY // write_tri_data(&out); - + // now copy the results back into the corresponding TGTriangle // structures - + // nodes out_nodes.clear(); for ( i = 0; i < out.numberofpoints; ++i ) { - Point3D p( out.pointlist[2*i], out.pointlist[2*i + 1], - out.pointattributelist[i] ); - out_nodes.simple_add( p ); + Point3D p( out.pointlist[2*i], out.pointlist[2*i + 1], + out.pointattributelist[i] ); + out_nodes.simple_add( p ); } - + // segments out_segs.clear(); for ( i = 0; i < out.numberofsegments; ++i ) { - out_segs.unique_add( TGTriSeg( out.segmentlist[2*i], - out.segmentlist[2*i+1], - out.segmentmarkerlist[i] ) ); + out_segs.unique_add( TGTriSeg( out.segmentlist[2*i], + out.segmentlist[2*i+1], + out.segmentmarkerlist[i] ) ); } - + // triangles elelist.clear(); int n1, n2, n3; double attribute; for ( i = 0; i < out.numberoftriangles; ++i ) { - n1 = out.trianglelist[i * 3]; - n2 = out.trianglelist[i * 3 + 1]; - n3 = out.trianglelist[i * 3 + 2]; - if ( out.numberoftriangleattributes > 0 ) { - attribute = out.triangleattributelist[i]; - } else { - attribute = 0.0; - } - // cout << "triangle = " << n1 << " " << n2 << " " << n3 << endl; - - elelist.push_back( TGTriEle( n1, n2, n3, attribute ) ); + n1 = out.trianglelist[i * 3]; + n2 = out.trianglelist[i * 3 + 1]; + n3 = out.trianglelist[i * 3 + 2]; + if ( out.numberoftriangleattributes > 0 ) { + attribute = out.triangleattributelist[i]; + } else { + attribute = 0.0; + } + // cout << "triangle = " << n1 << " " << n2 << " " << n3 << endl; + + elelist.push_back( TGTriEle( n1, n2, n3, attribute ) ); } - + // free mem allocated to the "Triangle" structures free(in.pointlist); free(in.pointattributelist); @@ -534,7 +542,7 @@ int TGTriangle::run_triangulate( double angle, const int pass ) { free(vorout.pointattributelist); free(vorout.edgelist); free(vorout.normlist); - + return 0; }