From 5f164a7720320f5356bc4d0efea85ae5702cd23c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ivan=20Pa=C4=91en?= <49401914+ipadjen@users.noreply.github.com> Date: Sat, 13 Jan 2024 10:55:27 +0100 Subject: [PATCH 1/7] Bump dev version --- CHANGELOG.md | 2 ++ src/main.cpp | 2 +- 2 files changed, 3 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 5955fc7e..4c00a215 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,5 +1,7 @@ # Changelog +## [Unreleased] + ## [0.4.4] - 2024-01-13 ### Changed - Improved fallbacks in case building import fails in late stages diff --git a/src/main.cpp b/src/main.cpp index 3b5ac231..c465b5e8 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -31,7 +31,7 @@ #include -std::string CITY4CFD_VERSION = "0.4.4"; +std::string CITY4CFD_VERSION = "0.4.4+dev"; void printWelcome() { auto logo{ From a505543bf3162fd3d7bb48680e73fcf594091a94 Mon Sep 17 00:00:00 2001 From: Ivan Paden Date: Wed, 6 Mar 2024 09:30:40 +0100 Subject: [PATCH 2/7] Improve checks for polygon import and failed reconstructions --- CMakeLists.txt | 2 +- src/Map3d.cpp | 35 ++++++++++++++++------------------- src/Map3d.h | 1 + src/PolyFeature.cpp | 26 +++++++++++++++----------- 4 files changed, 33 insertions(+), 31 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 2a5ac86a..f28898bb 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -4,7 +4,7 @@ project(city4cfd) #set( CMAKE_ALLOW_LOOSE_LOOP_CONSTRUCTS true ) -set(CMAKE_CXX_FLAGS "-O2") +set(CMAKE_CXX_FLAGS "-O3") set(CMAKE_BUILD_TYPE "Release") if (COMMAND cmake_policy) diff --git a/src/Map3d.cpp b/src/Map3d.cpp index 5f4280d8..41993b25 100644 --- a/src/Map3d.cpp +++ b/src/Map3d.cpp @@ -112,7 +112,7 @@ void Map3d::set_features() { _buildingsPtr.push_back(building); _allFeaturesPtr.push_back(building); } - if (Config::get().avoidBadPolys) this->clear_inactives(); // Remove buildings that potentially couldn't be imported + this->clear_inactives(); // Remove buildings that potentially couldn't be imported //- Imported buildings if (!_importedBuildingsJSON.empty()) { std::cout << "Importing CityJSON geometries" << std::endl; @@ -170,6 +170,7 @@ void Map3d::set_features() { _allFeaturesPtr.push_back(surfacePoly); } } + this->clear_inactives(); std::cout << "Polygons read: " << _allFeaturesPtr.size() << std::endl; //-- Set flat terrain or random thin terrain points @@ -349,15 +350,14 @@ void Map3d::reconstruct_buildings() { auto& b = _buildingsPtr[i]; if (b->is_active()) this->reconstruct_one_building(b); } - this->clear_inactives(); // in case of imported-reconstructed fallback + this->clear_inactives(); // renumber failed and in case of imported-reconstructed fallback // Gather failed reconstructions - int failed = 0; - for (auto& b : _buildingsPtr) if (b->has_failed_to_reconstruct()) ++failed; - std::cout << " Number of successfully reconstructed buildings: " << _buildingsPtr.size() - failed << std::endl; + std::cout << " Number of successfully reconstructed buildings: " << _buildingsPtr.size() + -_failedBuildingsPtr.size() << std::endl; Config::get().logSummary << "Building reconstruction summary: successfully reconstructed buildings: " - << _buildingsPtr.size() - failed << std::endl; + << _buildingsPtr.size() - _failedBuildingsPtr.size() << std::endl; Config::get().logSummary << " num of failed reconstructions: " - << failed << std::endl; + << _failedBuildingsPtr.size() << std::endl; } void Map3d::reconstruct_one_building(std::shared_ptr& building) { @@ -580,7 +580,7 @@ void Map3d::prep_cityjson_output() { // Temp impl, might change }; void Map3d::clear_inactives() { - for (unsigned long i = 0; i < _reconstructedBuildingsPtr.size();) { + for (int i = 0; i < _reconstructedBuildingsPtr.size();) { if (_reconstructedBuildingsPtr[i]->is_active()) ++i; else { _reconstructedBuildingsPtr.erase(_reconstructedBuildingsPtr.begin() + i); @@ -596,7 +596,7 @@ void Map3d::clear_inactives() { inactiveBuildingIdxs.push_back(importedBuilding->get_id()); } } - for (unsigned long i = 0; i < _importedBuildingsPtr.size();) { + for (int i = 0; i < _importedBuildingsPtr.size();) { auto it = std::find(inactiveBuildingIdxs.begin(), inactiveBuildingIdxs.end(), _importedBuildingsPtr[i]->get_id()); if (it == inactiveBuildingIdxs.end()) ++i; @@ -606,26 +606,27 @@ void Map3d::clear_inactives() { } } } else { - for (unsigned long i = 0; i < _importedBuildingsPtr.size();) { + for (int i = 0; i < _importedBuildingsPtr.size();) { if (_importedBuildingsPtr[i]->is_active()) ++i; else { _importedBuildingsPtr.erase(_importedBuildingsPtr.begin() + i); } } } - for (unsigned long i = 0; i < _buildingsPtr.size();) { - if (_buildingsPtr[i]->is_active() || _buildingsPtr[i]->has_failed_to_reconstruct()) ++i; + for (int i = 0; i < _buildingsPtr.size();) { + if (_buildingsPtr[i]->is_active()) ++i; else { + if (_buildingsPtr[i]->has_failed_to_reconstruct()) _failedBuildingsPtr.push_back(_buildingsPtr[i]); _buildingsPtr.erase(_buildingsPtr.begin() + i); } } - for (unsigned long i = 0; i < _surfaceLayersPtr.size();) { + for (int i = 0; i < _surfaceLayersPtr.size();) { if (_surfaceLayersPtr[i]->is_active()) ++i; else { _surfaceLayersPtr.erase(_surfaceLayersPtr.begin() + i); } } - for (unsigned long i = 0; i < _allFeaturesPtr.size();) { + for (int i = 0; i < _allFeaturesPtr.size();) { if (_allFeaturesPtr[i]->is_active()) ++i; else { _allFeaturesPtr.erase(_allFeaturesPtr.begin() + i); @@ -634,11 +635,7 @@ void Map3d::clear_inactives() { } BuildingsPtr Map3d::get_failed_buildings() const { - BuildingsPtr failedBuildings; - for (auto& b : _buildingsPtr) { - if (b->has_failed_to_reconstruct()) failedBuildings.push_back(b); - } - return failedBuildings; + return _failedBuildingsPtr; } //-- Templated functions diff --git a/src/Map3d.h b/src/Map3d.h index 607a7e14..a9cc253f 100644 --- a/src/Map3d.h +++ b/src/Map3d.h @@ -57,6 +57,7 @@ class Map3d { TerrainPtr _terrainPtr; BuildingsPtr _buildingsPtr; + BuildingsPtr _failedBuildingsPtr; ReconstructedBuildingsPtr _reconstructedBuildingsPtr; ImportedBuildingsPtr _importedBuildingsPtr; SurfaceLayersPtr _surfaceLayersPtr; diff --git a/src/PolyFeature.cpp b/src/PolyFeature.cpp index c5f850d7..1193eddd 100644 --- a/src/PolyFeature.cpp +++ b/src/PolyFeature.cpp @@ -89,12 +89,6 @@ PolyFeature::PolyFeature(const Polygon_with_attr& poly, const bool checkSimplici this->deactivate(); return; } - if (isOuterRing) { - if (tempPoly.is_clockwise_oriented()) tempPoly.reverse_orientation(); - isOuterRing = false; - } else { - if (tempPoly.is_counterclockwise_oriented()) tempPoly.reverse_orientation(); - } if (checkSimplicity) { if (!tempPoly.is_simple()) { if (Config::get().avoidBadPolys) { @@ -109,6 +103,12 @@ PolyFeature::PolyFeature(const Polygon_with_attr& poly, const bool checkSimplici } } } + if (isOuterRing) { + if (tempPoly.is_clockwise_oriented()) tempPoly.reverse_orientation(); + isOuterRing = false; + } else { + if (tempPoly.is_counterclockwise_oriented()) tempPoly.reverse_orientation(); + } _poly._rings.push_back(tempPoly); } } @@ -425,11 +425,10 @@ void PolyFeature::parse_json_poly(const nlohmann::json& poly, const bool checkSi } } geomutils::pop_back_if_equal_to_front(tempPoly); - - if (_poly._rings.empty()) { - if (tempPoly.is_clockwise_oriented()) tempPoly.reverse_orientation(); - } else { - if (tempPoly.is_counterclockwise_oriented()) tempPoly.reverse_orientation(); + if (tempPoly.size() < 3) { // Sanity check if it is even a polygon + std::cout << "WARNING: Skipping import of a zero-area polygon" << std::endl; + this->deactivate(); + return; } if (checkSimplicity) { if (!tempPoly.is_simple()) { @@ -445,6 +444,11 @@ void PolyFeature::parse_json_poly(const nlohmann::json& poly, const bool checkSi } } } + if (_poly._rings.empty()) { + if (tempPoly.is_clockwise_oriented()) tempPoly.reverse_orientation(); + } else { + if (tempPoly.is_counterclockwise_oriented()) tempPoly.reverse_orientation(); + } _poly._rings.push_back(tempPoly); } } \ No newline at end of file From d6ac8405b64a5d05e19f5025101ae46b648c701d Mon Sep 17 00:00:00 2001 From: Ivan Paden Date: Wed, 6 Mar 2024 09:32:35 +0100 Subject: [PATCH 3/7] Calcualte blockage ratio only when needed --- src/BoundingRegion.cpp | 28 +++++++++++++++------------- 1 file changed, 15 insertions(+), 13 deletions(-) diff --git a/src/BoundingRegion.cpp b/src/BoundingRegion.cpp index 3ecff697..43ac6fba 100644 --- a/src/BoundingRegion.cpp +++ b/src/BoundingRegion.cpp @@ -116,25 +116,27 @@ void BoundingRegion::calc_bnd_bpg(const Polygon_2& influRegionPoly, Config::get().topHeight = elevationMax + hMax * (Config::get().bpgDomainSize.back() - 1); //-- Blockage ratio handling - std::cout << "\nCalculating blockage ratio for flow direction (" << Config::get().flowDirection - << ")" << std::endl; + if (Config::get().bpgBlockageRatioFlag) { + std::cout << "\nCalculating blockage ratio for flow direction (" << Config::get().flowDirection + << ")" << std::endl; // double blockRatio = this->calc_blockage_ratio_from_chull(buildings, angle, localPoly); // double blockRatio = this->calc_blockage_ratio_from_ashape(buildings, angle, localPoly); - double blockRatio = this->calc_blockage_ratio_comb(buildings, angle, localPoly); + double blockRatio = this->calc_blockage_ratio_comb(buildings, angle, localPoly); // double blockRatio = this->calc_blockage_ratio_from_ashape_alt(buildings, angle, localPoly); // double blockRatio = this->calc_blockage_ratio_from_edges(buildings, angle, localPoly); - std::cout << " Blockage ratio is: " << blockRatio << std::endl; - if (Config::get().bpgBlockageRatioFlag && blockRatio > Config::get().bpgBlockageRatio) { - std::cout << "INFO: Blockage ratio is more than " << Config::get().bpgBlockageRatio * 100 - << "%. Expanding domain cross section to meet the guideline" - << std::endl; - double expRatio = std::sqrt(blockRatio / 0.03); - //-- Recalculate the bnd poly and height with new values - localPoly.clear(); - localPoly = this->calc_bnd_poly(candidatePts, hMax, angle, expRatio); - Config::get().topHeight = hMax * Config::get().bpgDomainSize.back() * expRatio; + std::cout << " Blockage ratio is: " << blockRatio << std::endl; + if (blockRatio > Config::get().bpgBlockageRatio) { + std::cout << "INFO: Blockage ratio is more than " << Config::get().bpgBlockageRatio * 100 + << "%. Expanding domain cross section to meet the guideline" + << std::endl; + double expRatio = std::sqrt(blockRatio / 0.03); + //-- Recalculate the bnd poly and height with new values + localPoly.clear(); + localPoly = this->calc_bnd_poly(candidatePts, hMax, angle, expRatio); + Config::get().topHeight = hMax * Config::get().bpgDomainSize.back() * expRatio; + } } //-- Return the points back to global coordinates for (auto& pt : localPoly) _boundingRegion.push_back(geomutils::rotate_pt(pt, angle)); From 32352a75dfdf317f087ffc6280bf97accfb4884b Mon Sep 17 00:00:00 2001 From: Ivan Paden Date: Wed, 6 Mar 2024 15:55:41 +0100 Subject: [PATCH 4/7] Fix polygon flattening with holes --- src/PolyFeature.cpp | 34 +++++++++++++++------------------- src/geomutils.cpp | 33 ++++++++++++++++++++++++++++++++- src/geomutils.h | 2 ++ 3 files changed, 49 insertions(+), 20 deletions(-) diff --git a/src/PolyFeature.cpp b/src/PolyFeature.cpp index 1193eddd..2c7860ab 100644 --- a/src/PolyFeature.cpp +++ b/src/PolyFeature.cpp @@ -237,9 +237,7 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, std::map> overlappingBuildings; //id-building map for (auto& pt3 : subsetPts) { Point_2 pt(pt3.x(), pt3.y()); - if (CGAL::bounded_side_2(_poly._rings.front().begin(), - _poly._rings.front().end(), - pt) != CGAL::ON_UNBOUNDED_SIDE) { + if (geomutils::point_in_poly_and_boundary(pt, _poly)) { auto itIdx = pointCloudConnectivity.find(pt3); auto pointSetIt = pointCloud.begin(); std::advance(pointSetIt, itIdx->second); @@ -277,7 +275,6 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, flattenCandidatePolys.emplace_back(geomutils::exact_poly_to_poly(flattenBndPoly)); } } - std::vector flattenBndPolys; if (!flattenCandidatePolys.empty()) { bool isFirst = true; for (auto& poly : flattenCandidatePolys) { @@ -285,10 +282,10 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, PolygonPtrVectorWH offset_poly = CGAL::create_interior_skeleton_and_offset_polygons_with_holes_2(offsetVal, poly.get_cgal_type()); - if (offset_poly.size() == 1) { // make a check whether the offset is successfully created or not - flattenBndPolys.push_back(offset_poly.front()->outer_boundary()); + if (offset_poly.size() > 0) { // make a check whether the offset is successfully created or not + poly = *(offset_poly.front()); } else { - std::cout << "Skeleton construction failed!" << std::endl; + std::cout << "WARNING: Skeleton construction failed! Some surfaces might not be flattened." << std::endl; return false; } if (isFirst) { @@ -300,15 +297,13 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, } } } else { - flattenBndPolys.push_back(_poly.outer_boundary()); + flattenCandidatePolys.push_back(_poly); } //-- Collect points that have not been already flattened for (auto& pt3 : subsetPts) { - for (auto& flattenBndPoly : flattenBndPolys) { + for (auto& flattenBndPoly : flattenCandidatePolys) { Point_2 pt(pt3.x(), pt3.y()); - if (CGAL::bounded_side_2(flattenBndPoly.begin(), - flattenBndPoly.end(), - pt) != CGAL::ON_UNBOUNDED_SIDE) { + if (geomutils::point_in_poly_and_boundary(pt, flattenBndPoly)) { auto itIdx = pointCloudConnectivity.find(pt3); auto it = flattenedPts.find(itIdx->second); @@ -331,16 +326,17 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, flattenedPts[i] = Point_3(pointCloud.point(i).x(), pointCloud.point(i).y(), avgHeight); } //-- Add additional new segments to constrain - if (!flattenBndPolys.empty()) { - for (auto& flattenBndPoly : flattenBndPolys) { - for (const auto& newEdge: flattenBndPoly.edges()) { - ePoint_3 pt1(newEdge.point(0).x(), newEdge.point(0).y(), avgHeight); - ePoint_3 pt2(newEdge.point(1).x(), newEdge.point(1).y(), avgHeight); - constrainedEdges.emplace_back(pt1, pt2); + if (!flattenCandidatePolys.empty()) { + for (auto& flattenBndPoly : flattenCandidatePolys) { + for (auto& ring: flattenBndPoly.rings()) { + for (const auto& newEdge: ring.edges()) { + ePoint_3 pt1(newEdge.point(0).x(), newEdge.point(0).y(), avgHeight); + ePoint_3 pt2(newEdge.point(1).x(), newEdge.point(1).y(), avgHeight); + constrainedEdges.emplace_back(pt1, pt2); + } } } } - return true; } diff --git a/src/geomutils.cpp b/src/geomutils.cpp index 03e75feb..09c3b8d3 100644 --- a/src/geomutils.cpp +++ b/src/geomutils.cpp @@ -272,7 +272,6 @@ template bool geomutils::point_in_poly(const Point_3& pt2, const Polygo template bool geomutils::point_in_poly(const T& pt2, const Polygon_with_holes_2& polygon) { Point_2 pt(pt2.x(), pt2.y()); - //-- Check if the point falls within the outer surface if (point_in_poly(pt, polygon.outer_boundary())) { // Check if the point falls within one of the holes @@ -289,6 +288,38 @@ bool geomutils::point_in_poly(const T& pt2, const Polygon_with_holes_2& polygon) template bool geomutils::point_in_poly(const Point_2& pt2, const Polygon_with_holes_2& polygon); template bool geomutils::point_in_poly(const Point_3& pt2, const Polygon_with_holes_2& polygon); +//-- Check if the point is inside a polygon on a 2D projection including polygon boundaries +template +bool geomutils::point_in_poly_and_boundary(const T& pt2, const Polygon_2& polygon) { + Point_2 pt(pt2.x(), pt2.y()); + if (CGAL::bounded_side_2(polygon.begin(), polygon.end(), pt) != CGAL::ON_UNBOUNDED_SIDE) { + return true; + } + return false; +} +//- Explicit template instantiation +template bool geomutils::point_in_poly_and_boundary(const Point_2& pt2, const Polygon_2& polygon); +template bool geomutils::point_in_poly_and_boundary(const Point_3& pt2, const Polygon_2& polygon); + +template +bool geomutils::point_in_poly_and_boundary(const T& pt2, const Polygon_with_holes_2& polygon) { + Point_2 pt(pt2.x(), pt2.y()); + //-- Check if the point falls within the outer surface + if (point_in_poly_and_boundary(pt, polygon.outer_boundary())) { + // Check if the point falls within one of the holes + for (auto it_hole = polygon.holes_begin(); it_hole != polygon.holes_end(); ++it_hole) { + if (point_in_poly_and_boundary(pt, *it_hole)) { + return false; + } + } + return true; + } + return false; +} +//- Explicit template instantiation +template bool geomutils::point_in_poly_and_boundary(const Point_2& pt2, const Polygon_with_holes_2& polygon); +template bool geomutils::point_in_poly_and_boundary(const Point_3& pt2, const Polygon_with_holes_2& polygon); + template void geomutils::make_round_poly(const Point_2& centre, double radius, T& poly) { const int nPts = 360; // Hardcoded diff --git a/src/geomutils.h b/src/geomutils.h index 29f6bc36..88ccbd42 100644 --- a/src/geomutils.h +++ b/src/geomutils.h @@ -50,6 +50,8 @@ namespace geomutils { //-- Templated functions template bool point_in_poly(const T& pt2, const Polygon_with_holes_2& polygon); template bool point_in_poly(const T& pt2, const Polygon_2& polygon); + template bool point_in_poly_and_boundary(const T& pt2, const Polygon_with_holes_2& polygon); + template bool point_in_poly_and_boundary(const T& pt2, const Polygon_2& polygon); template void make_round_poly(const Point_2& centre, double radius, T& poly); template void make_round_poly(const Point_2& centre, double radius1, double radius2, int nPts, double angInt, double ang, T& poly); From 86e9506fc5f6560c41166f1afccb4f33b2e5cc22 Mon Sep 17 00:00:00 2001 From: Ivan Paden Date: Wed, 6 Mar 2024 17:14:52 +0100 Subject: [PATCH 5/7] Fix log --- src/Map3d.cpp | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/src/Map3d.cpp b/src/Map3d.cpp index 41993b25..bfaf8fc7 100644 --- a/src/Map3d.cpp +++ b/src/Map3d.cpp @@ -352,10 +352,9 @@ void Map3d::reconstruct_buildings() { } this->clear_inactives(); // renumber failed and in case of imported-reconstructed fallback // Gather failed reconstructions - std::cout << " Number of successfully reconstructed buildings: " << _buildingsPtr.size() - -_failedBuildingsPtr.size() << std::endl; + std::cout << " Number of successfully reconstructed buildings: " << _buildingsPtr.size() << std::endl; Config::get().logSummary << "Building reconstruction summary: successfully reconstructed buildings: " - << _buildingsPtr.size() - _failedBuildingsPtr.size() << std::endl; + << _buildingsPtr.size() << std::endl; Config::get().logSummary << " num of failed reconstructions: " << _failedBuildingsPtr.size() << std::endl; } From ddfe53372db11c98d2a5bb6f4f5a8691731ab0ad Mon Sep 17 00:00:00 2001 From: Ivan Paden Date: Sat, 9 Mar 2024 10:13:24 +0100 Subject: [PATCH 6/7] Fix polyset handling in 2d boolean --- src/PolyFeature.cpp | 36 ++++++++++++++++++++++++++++++------ src/Terrain.cpp | 2 +- 2 files changed, 31 insertions(+), 7 deletions(-) diff --git a/src/PolyFeature.cpp b/src/PolyFeature.cpp index 2c7860ab..7cc2ba8a 100644 --- a/src/PolyFeature.cpp +++ b/src/PolyFeature.cpp @@ -222,6 +222,15 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, typedef boost::shared_ptr> PolygonPtrWH; typedef std::vector PolygonPtrVectorWH; +#ifdef CITY4CFD_POLYFEATURE_VERBOSE + std::cout << "\nINFO: Flattening a polygon..." << std::endl; + int nonSimpleRings = 0; + for (auto& ring : _poly.rings()) { + if (!ring.is_simple()) ++nonSimpleRings; + } + std::cout << "Non-simple rings: " << nonSimpleRings << std::endl; +#endif + std::vector indices; std::vector originalHeights; auto building_pt = pointCloud.property_map>("building_point").first; @@ -256,25 +265,32 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, //-- If next intersecting a building, clip poly with building std::vector flattenCandidatePolys; if (!overlappingBuildings.empty()) { +#ifdef CITY4CFD_POLYFEATURE_VERBOSE + std::cout << "Overlapping buildings: " << overlappingBuildings.size() << std::endl; +#endif // Add clipping polys to set CGAL::Polygon_set_2 polySet; Converter to_exact; for (auto& intersectBuilding : overlappingBuildings) { - polySet.insert(intersectBuilding.second->get_poly().get_exact_outer_boundary()); + polySet.join(intersectBuilding.second->get_poly().get_exact_outer_boundary()); } +#ifdef CITY4CFD_POLYFEATURE_VERBOSE + std::cout << "Polyset size? " << polySet.number_of_polygons_with_holes() << std::endl; +#endif + CGAL_assertion(polySet.is_valid()); // Clip with this polygon polySet.complement(); - polySet.intersection(this->get_poly().get_exact()); + polySet.intersection(_poly.get_exact()); // Store this as the new polygon std::vector> resPolys; polySet.polygons_with_holes(std::back_inserter(resPolys)); -#ifdef CITY4CFD_POLYFEATURE_VERBOSE - std::cout << "After flatten polygon clipping, total new polygons: " << resPolys.size() << std::endl; -#endif for (auto& flattenBndPoly : resPolys) { flattenCandidatePolys.emplace_back(geomutils::exact_poly_to_poly(flattenBndPoly)); } } +#ifdef CITY4CFD_POLYFEATURE_VERBOSE + std::cout << "After flatten polygon clipping, total new polygons: " << flattenCandidatePolys.size() << std::endl; +#endif if (!flattenCandidatePolys.empty()) { bool isFirst = true; for (auto& poly : flattenCandidatePolys) { @@ -282,7 +298,10 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, PolygonPtrVectorWH offset_poly = CGAL::create_interior_skeleton_and_offset_polygons_with_holes_2(offsetVal, poly.get_cgal_type()); - if (offset_poly.size() > 0) { // make a check whether the offset is successfully created or not +#ifdef CITY4CFD_POLYFEATURE_VERBOSE + std::cout << "Offset polygons created: " << offset_poly.size() << std::endl; +#endif + if (!offset_poly.empty()) { // make a check whether the offset is successfully created poly = *(offset_poly.front()); } else { std::cout << "WARNING: Skeleton construction failed! Some surfaces might not be flattened." << std::endl; @@ -297,6 +316,11 @@ bool PolyFeature::flatten_polygon_inner_points(const Point_set_3& pointCloud, } } } else { + if (!overlappingBuildings.empty()) { + std::cout << "WARNING: Polygon ID" << _id << " flattening failed due to unresolved" + " intersections with buildings!" << std::endl; + return false; + } flattenCandidatePolys.push_back(_poly); } //-- Collect points that have not been already flattened diff --git a/src/Terrain.cpp b/src/Terrain.cpp index acbf22ea..94eda108 100644 --- a/src/Terrain.cpp +++ b/src/Terrain.cpp @@ -102,7 +102,7 @@ void Terrain::constrain_features() { } IO::print_progress_bar(100); std::clog << std::endl; // extra edges to constrain when whole polygons couldn't be added - if (!_extraConstrainedEdges.empty()) std::cout << "\n Adding extra constrained edges" << std::endl; + if (!_extraConstrainedEdges.empty()) std::cout << "\n Inserting additional constrained edges" << std::endl; for (auto& extraEdge : _extraConstrainedEdges) { _cdt.insert_constraint(extraEdge.source(), extraEdge.target()); ++count; From 4246d242df9b588a0b1cc7d16eca20e5abf53ae0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ivan=20Pa=C4=91en?= <49401914+ipadjen@users.noreply.github.com> Date: Mon, 11 Mar 2024 11:27:45 +0100 Subject: [PATCH 7/7] Bump version and add changelog --- CHANGELOG.md | 6 +++++- src/main.cpp | 4 +--- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 9a66f884..816f83ff 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,6 +1,10 @@ # Changelog -## [Unreleased] +## [0.4.6] - 2024-03-11 +### Changed +- Improved robustness in handling bad polygons +### Fixed +- Code crashing when flattening specific cases ## [0.4.5] - 2024-02-12 ### Fixed diff --git a/src/main.cpp b/src/main.cpp index 5a1851e9..cfd6bade 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -1,6 +1,4 @@ /* - City4CFD - Copyright (c) 2021-2024, 3D Geoinformation Research Group, TU Delft This file is part of City4CFD. @@ -31,7 +29,7 @@ #include -std::string CITY4CFD_VERSION = "0.4.5+dev"; +std::string CITY4CFD_VERSION = "0.4.6"; void printWelcome() { auto logo{