diff --git a/CHANGELOG.md b/CHANGELOG.md index 53276a63..7bb204aa 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,4 +1,11 @@ # Changelog + +## [0.1.1] - 2022-09-02 +### Fixed +- Bad triangles at surface edges +- Polygons close to the boundary are defined as 'out of bounds' + + ## [0.1.0] - 2022-08-14 First release diff --git a/CITATION.bib b/CITATION.bib index f94e3b1e..7ec1655d 100644 --- a/CITATION.bib +++ b/CITATION.bib @@ -1,5 +1,5 @@ @article{Paden2022, - author={Paden, Ivan, Garc{\'i}a-S{\'a}nchez Clara and Ledoux, Hugo}, + author={Paden, Ivan and Garc{\'i}a-S{\'a}nchez, Clara and Ledoux, Hugo}, title={Towards Automatic Reconstruction of 3D City Models Tailored for Urban Flow Simulations}, journal={Frontiers in Built Environment}, volume={8}, diff --git a/README.md b/README.md index 1c963aff..36fcd45e 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ -[![docs](https://img.shields.io/badge/docs-Wiki-brightgreen)](https://github.com/tudelft3d/City4CFD/wiki) -[![GitHub license](https://img.shields.io/github/license/tudelft3d/City4CFD)](https://github.com/tudelft3d/City4CFD/blob/master/LICENSE) -[![DOI:10.3389/fbuil.2022.899332](http://img.shields.io/badge/DOI-10.3389/fbuil.2022.899332-B62030.svg)](https://www.frontiersin.org/articles/10.3389/fbuil.2022.899332) +[![docs](https://img.shields.io/badge/docs-Wiki-brightgreen?style=flat-square)](https://github.com/tudelft3d/City4CFD/wiki) +[![GitHub license](https://img.shields.io/github/license/tudelft3d/City4CFD?style=flat-square)](https://github.com/tudelft3d/City4CFD/blob/master/LICENSE) +[![DOI:10.3389/fbuil.2022.899332](http://img.shields.io/badge/DOI-10.3389/fbuil.2022.899332-B62030.svg?style=flat-square)](https://www.frontiersin.org/articles/10.3389/fbuil.2022.899332) # City4CFD diff --git a/src/Boundary.cpp b/src/Boundary.cpp index 0f31c230..b50dd428 100644 --- a/src/Boundary.cpp +++ b/src/Boundary.cpp @@ -59,7 +59,7 @@ void Boundary::set_bnd_poly(Polygon_2& bndPoly, Polygon_2& pcBndPoly, Polygon_2& } for (auto& pt : startBufferPoly) - pcBndPoly.push_back(pt - (pt - center) * 0.05 * std::abs(bufferLen)); + pcBndPoly.push_back(pt - (pt - center) * 0.001 * std::abs(bufferLen)); } else { startBufferPoly = bndPoly; pcBndPoly = startBufferPoly; @@ -67,7 +67,7 @@ void Boundary::set_bnd_poly(Polygon_2& bndPoly, Polygon_2& pcBndPoly, Polygon_2& } //-- Deactivate point cloud points that are out of bounds -void Boundary::set_bounds_to_pc(Point_set_3& pointCloud, const Polygon_2& pcBndPoly) { +void Boundary::set_bounds_to_buildings_pc(Point_set_3& pointCloud, const Polygon_2& pcBndPoly) { //-- Remove points out of the boundary region auto it = pointCloud.points().begin(); int count = 0; @@ -82,10 +82,10 @@ void Boundary::set_bounds_to_pc(Point_set_3& pointCloud, const Polygon_2& pcBndP pointCloud.collect_garbage(); // Free removed points from the memory } -void Boundary::set_bounds_to_terrain(Point_set_3& pointCloud, const Polygon_2& bndPoly, - const Polygon_2& pcBndPoly, const Polygon_2& startBufferPoly) { +void Boundary::set_bounds_to_terrain_pc(Point_set_3& pointCloud, const Polygon_2& bndPoly, + const Polygon_2& pcBndPoly, const Polygon_2& startBufferPoly) { //-- Remove points out of the boundary region - Boundary::set_bounds_to_pc(pointCloud, pcBndPoly); + Boundary::set_bounds_to_buildings_pc(pointCloud, pcBndPoly); //-- Add outer points to match the domain size to prescribed one SearchTree searchTree(pointCloud.points().begin(),pointCloud.points().end()); diff --git a/src/Boundary.h b/src/Boundary.h index ae2b50a4..b3f4ef2e 100644 --- a/src/Boundary.h +++ b/src/Boundary.h @@ -37,9 +37,9 @@ class Boundary : public TopoFeature { virtual ~Boundary(); static void set_bnd_poly(Polygon_2& bndPoly, Polygon_2& pcBndPoly, Polygon_2& startBufferPoly); - static void set_bounds_to_pc(Point_set_3& pointCloud, const Polygon_2& pcBndPoly); - static void set_bounds_to_terrain(Point_set_3& pointCloud, const Polygon_2& bndPoly, - const Polygon_2& pcBndPoly, const Polygon_2& startBufferPoly); + static void set_bounds_to_buildings_pc(Point_set_3& pointCloud, const Polygon_2& pcBndPoly); + static void set_bounds_to_terrain_pc(Point_set_3& pointCloud, const Polygon_2& bndPoly, + const Polygon_2& pcBndPoly, const Polygon_2& startBufferPoly); static std::vector get_domain_bbox(); virtual void reconstruct() = 0; diff --git a/src/Map3d.cpp b/src/Map3d.cpp index cff1ef4b..779c1656 100644 --- a/src/Map3d.cpp +++ b/src/Map3d.cpp @@ -248,13 +248,13 @@ void Map3d::set_bnd() { Boundary::set_bnd_poly(bndPoly, pcBndPoly, startBufferPoly); //-- Deactivate point cloud points that are out of bounds - Boundary::set_bounds_to_terrain(_pointCloud.get_terrain(), - bndPoly, pcBndPoly, startBufferPoly); - Boundary::set_bounds_to_pc(_pointCloud.get_buildings(), startBufferPoly); + Boundary::set_bounds_to_terrain_pc(_pointCloud.get_terrain(), + bndPoly, pcBndPoly, startBufferPoly); + Boundary::set_bounds_to_buildings_pc(_pointCloud.get_buildings(), startBufferPoly); //-- Check feature scope for surface layers now that the full domain is known for (auto& f: _surfaceLayers) { - f->check_feature_scope(pcBndPoly); + f->check_feature_scope(bndPoly); } this->clear_inactives(); } diff --git a/src/io.cpp b/src/io.cpp index 725c6029..d4786b3d 100644 --- a/src/io.cpp +++ b/src/io.cpp @@ -291,6 +291,7 @@ void IO::get_obj_pts(const Mesh& mesh, std::unordered_map& dPts) { for (auto& face : mesh.faces()) { + if (IO::is_degen(mesh, face)) continue; std::vector faceIdx; faceIdx.reserve(3); std::string fsTemp; std::string bsTemp; @@ -309,13 +310,7 @@ void IO::get_obj_pts(const Mesh& mesh, faceIdx.push_back(it->second); } } - - if (IO::not_small(faceIdx)) { - bs += "\nf"; - bs += bsTemp; -// } else { -// std::cerr << "Found duplicates!" << std::endl; - } + bs += "\nf" + bsTemp; } } @@ -324,15 +319,11 @@ void IO::get_stl_pts(Mesh& mesh, std::string& fs) { auto fnormals = mesh.add_property_map("f:normals", CGAL::NULL_VECTOR).first; PMP::compute_normals(mesh, vnormals, fnormals); for (auto& face : mesh.faces()) { + if (IO::is_degen(mesh, face)) continue; std::vector outputPts; for (auto index: CGAL::vertices_around_face(mesh.halfedge(face), mesh)) { outputPts.push_back(gen_key_bucket(mesh.point(index))); } - //-- Check for round off of small triangles - if (outputPts[0] == outputPts[1] - || outputPts[0] == outputPts[2] - || outputPts[1] == outputPts[2]) continue; - fs += "\nfacet normal " + gen_key_bucket(fnormals[face]); fs += "\n outer loop"; for (auto pt: outputPts) { @@ -349,6 +340,7 @@ void IO::get_cityjson_geom(const Mesh& mesh, nlohmann::json& g, std::unordered_m g["lod"] = Config::get().lod; g["boundaries"]; for (auto& face: mesh.faces()) { + if (IO::is_degen(mesh, face)) continue; std::vector faceIdx; faceIdx.reserve(3); std::vector tempPoly; @@ -367,23 +359,34 @@ void IO::get_cityjson_geom(const Mesh& mesh, nlohmann::json& g, std::unordered_m tempPoly.push_back(it->second); } } - - if (IO::not_small(faceIdx)) { - g["boundaries"].push_back({tempPoly}); -// } else { -// std::cerr << "Found duplicates!" << std::endl; - } + g["boundaries"].push_back({tempPoly}); } } -//-- Check for round off of small triangles -bool IO::not_small(std::vector idxLst) { +bool IO::not_same(std::vector idxLst) { std::sort(idxLst.begin(), idxLst.end()); auto it = std::unique(idxLst.begin(), idxLst.end()); return (it == idxLst.end()); } +bool IO::is_degen(const Mesh& mesh, Mesh::Face_index face) { + std::vector pts; pts.reserve(3); + for (auto index: CGAL::vertices_around_face(mesh.halfedge(face), mesh)) { + pts.push_back(mesh.point(index)); + } + //-- Precondition - check that the points are not the same + if (CGAL::squared_distance(pts[0], pts[1]) < 1e-3 || + CGAL::squared_distance(pts[0], pts[2]) < 1e-3 || + CGAL::squared_distance(pts[1], pts[2]) < 1e-3) { + return true; + } + if (sin(CGAL::approximate_angle(pts[0], pts[1], pts[2]) * M_PI / 180) < 0.00174) { + return true; + } + return false; +} + void IO::output_log() { if (!Config::get().outputLog) return; fs::current_path(Config::get().outputDir); @@ -429,7 +432,7 @@ bool IO::has_substr(const std::string& strMain, const std::string& subStr) { std::string IO::gen_key_bucket(const Point_2 p) { std::stringstream ss; - ss << std::fixed << std::setprecision(3) << p.x() << " " << p.y(); + ss << std::fixed << std::setprecision(6) << p.x() << " " << p.y(); return ss.str(); } @@ -437,7 +440,7 @@ std::string IO::gen_key_bucket(const Point_2 p) { template std::string IO::gen_key_bucket(const T& p) { std::stringstream ss; - ss << std::fixed << std::setprecision(3) << p.x() << " " << p.y() << " " << p.z(); + ss << std::fixed << std::setprecision(6) << p.x() << " " << p.y() << " " << p.z(); return ss.str(); } //- Explicit template instantiation diff --git a/src/io.h b/src/io.h index 13b24015..9b040653 100644 --- a/src/io.h +++ b/src/io.h @@ -50,7 +50,8 @@ namespace IO { void get_stl_pts(Mesh& mesh, std::string& fs); void get_cityjson_geom(const Mesh& mesh, nlohmann::json& g, std::unordered_map& dPts, std::string primitive); - bool not_small(std::vector idxLst); + bool not_same(std::vector idxLst); + bool is_degen(const Mesh& mesh, Mesh::Face_index face); void output_log(); bool has_substr(const std::string& strMain, const std::string& subStr); diff --git a/src/main.cpp b/src/main.cpp index 70e67373..339ccc72 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -29,7 +29,7 @@ #include "io.h" #include "Map3d.h" -std::string CITY4CFD_VERSION = "0.1.0"; +std::string CITY4CFD_VERSION = "0.1.1"; void printWelcome() { auto logo{