OGS 6.1.0-1699-ge946d4c5f
ComputeSurfaceNodeIDsInPolygonalRegion.cpp
Go to the documentation of this file.
1 
13 #include <algorithm>
14 #include <memory>
15 #include <string>
16 #include <vector>
17 
18 #include <tclap/CmdLine.h>
19 
22 #include "BaseLib/BuildInfo.h"
23 #include "BaseLib/Error.h"
24 #include "BaseLib/FileTools.h"
25 #include "GeoLib/GEOObjects.h"
26 #include "GeoLib/Polygon.h"
27 #include "MathLib/Vector3.h"
29 #include "MeshLib/Mesh.h"
31 #include "MeshLib/Node.h"
32 
33 void writeToFile(std::string const& id_area_fname, std::string const& csv_fname,
34  std::vector<std::pair<std::size_t, double>> const& ids_and_areas,
35  std::vector<MeshLib::Node*> const& mesh_nodes)
36 {
37  std::ofstream ids_and_area_out(id_area_fname);
38  if (!ids_and_area_out) {
39  OGS_FATAL("Unable to open the file '%s' - aborting.",
40  id_area_fname.c_str());
41  }
42  std::ofstream csv_out(csv_fname);
43  if (!csv_out) {
44  OGS_FATAL("Unable to open the file '%s' - aborting.",
45  csv_fname.c_str());
46  }
47 
48  ids_and_area_out << std::setprecision(20);
49  csv_out << std::setprecision(20);
50 
51  ids_and_area_out << ids_and_areas[0].first << " " << ids_and_areas[0].second;
52  csv_out << "ID x y z area node_id\n"; // CSV header
53  csv_out << 0 << " " << *mesh_nodes[ids_and_areas[0].first]
54  << ids_and_areas[0].second << " " << ids_and_areas[0].first;
55  for (std::size_t k(1); k<ids_and_areas.size(); k++) {
56  ids_and_area_out << "\n"
57  << ids_and_areas[k].first << " " << ids_and_areas[k].second;
58  csv_out << "\n" << k << " " << *mesh_nodes[ids_and_areas[k].first]
59  << ids_and_areas[k].second << " " << ids_and_areas[k].first;
60  }
61  ids_and_area_out << "\n";
62  csv_out << "\n";
63 }
64 
65 int main (int argc, char* argv[])
66 {
67  ApplicationsLib::LogogSetup logog_setup;
68 
69  TCLAP::CmdLine cmd(
70  "Computes ids of mesh nodes that are in polygonal regions and resides "
71  "on the top surface. The polygonal regions have to be given in a gml- "
72  "or gli-file. The found mesh nodes and the associated area are written "
73  "as txt and csv data. The documentation is available at "
74  "https://docs.opengeosys.org/docs/tools/model-preparation/"
75  "computesurfacenodeidsinpolygonalregion.\n\n"
76  "OpenGeoSys-6 software, version " +
78  ".\n"
79  "Copyright (c) 2012-2019, OpenGeoSys Community "
80  "(http://www.opengeosys.org)",
82  TCLAP::ValueArg<std::string> mesh_in("m", "mesh-input-file",
83  "the name of the file containing the input mesh", true,
84  "", "file name of input mesh");
85  cmd.add(mesh_in);
86  TCLAP::ValueArg<std::string> geo_in("g", "geo-file",
87  "the name of the gml file containing the polygons", true,
88  "", "file name of input geometry");
89  cmd.add(geo_in);
90 
91  cmd.parse(argc, argv);
92 
93  std::unique_ptr<MeshLib::Mesh const> mesh(MeshLib::IO::readMeshFromFile(mesh_in.getValue()));
94  INFO("Mesh read: %u nodes, %u elements.", mesh->getNumberOfNodes(), mesh->getNumberOfElements());
95 
96  GeoLib::GEOObjects geo_objs;
97  FileIO::readGeometryFromFile(geo_in.getValue(), geo_objs);
98  std::vector<std::string> geo_names;
99  geo_objs.getGeometryNames(geo_names);
100  INFO("Geometry '%s' read: %u points, %u polylines.",
101  geo_names[0].c_str(),
102  geo_objs.getPointVec(geo_names[0])->size(),
103  geo_objs.getPolylineVec(geo_names[0])->size());
104 
105  MathLib::Vector3 const dir(0.0, 0.0, -1.0);
106  double angle(90);
107 
108  auto computeElementTopSurfaceAreas = [](MeshLib::Mesh const& mesh,
109  MathLib::Vector3 const& d, double angle)
110  {
111  std::unique_ptr<MeshLib::Mesh> surface_mesh(
114  *surface_mesh.get());
115  };
116 
117  std::vector<double> areas(computeElementTopSurfaceAreas(*mesh, dir, angle));
118  std::vector<MeshLib::Node*> all_sfc_nodes(
120  );
121 
122  std::for_each(all_sfc_nodes.begin(), all_sfc_nodes.end(),
123  [](MeshLib::Node* p) { (*p)[2] = 0.0; });
124 
125  std::vector<MeshLib::Node*> const& mesh_nodes(mesh->getNodes());
126  GeoLib::PolylineVec const* ply_vec(
127  geo_objs.getPolylineVecObj(geo_names[0])
128  );
129  std::vector<GeoLib::Polyline*> const& plys(*(ply_vec->getVector()));
130 
131  for (std::size_t j(0); j<plys.size(); j++) {
132  if (! plys[j]->isClosed()) {
133  continue;
134  }
135  std::string polygon_name;
136  ply_vec->getNameOfElement(plys[j], polygon_name);
137  if (polygon_name.empty())
138  polygon_name = "Polygon-" + std::to_string(j);
139  // create Polygon from Polyline
140  GeoLib::Polygon const& polygon(*(plys[j]));
141  // ids of mesh nodes on surface that are within the given polygon
142  std::vector<std::pair<std::size_t, double>> ids_and_areas;
143  for (std::size_t k(0); k<all_sfc_nodes.size(); k++) {
144  MeshLib::Node const& pnt(*(all_sfc_nodes[k]));
145  if (polygon.isPntInPolygon(pnt[0], pnt[1], pnt[2])) {
146  ids_and_areas.push_back(std::make_pair(pnt.getID(), areas[k]));
147  }
148  }
149  if (ids_and_areas.empty()) {
150  ERR("Polygonal part of surface '%s' doesn't contains nodes. No "
151  "output will be generated.",
152  polygon_name.c_str());
153  continue;
154  }
155 
156  std::string const out_path(BaseLib::extractPath(geo_in.getValue()));
157  std::string id_and_area_fname(out_path + polygon_name);
158  std::string csv_fname(out_path + polygon_name);
159  id_and_area_fname += std::to_string(j) + ".txt";
160  csv_fname += std::to_string(j) + ".csv";
161  INFO(
162  "Polygonal part of surface '%s' contains %ul nodes. Writting to"
163  " files '%s' and '%s'.",
164  polygon_name.c_str(),
165  ids_and_areas.size(),
166  id_and_area_fname.c_str(),
167  csv_fname.c_str());
168  writeToFile(id_and_area_fname, csv_fname, ids_and_areas, mesh_nodes);
169  }
170 
171  std::for_each(all_sfc_nodes.begin(), all_sfc_nodes.end(),
172  std::default_delete<MeshLib::Node>());
173 
174  return 0;
175 }
void readGeometryFromFile(std::string const &fname, GeoLib::GEOObjects &geo_objs)
Container class for geometric objects.
Definition: GEOObjects.h:62
Definition of the Node class.
void writeToFile(std::string const &id_area_fname, std::string const &csv_fname, std::vector< std::pair< std::size_t, double >> const &ids_and_areas, std::vector< MeshLib::Node *> const &mesh_nodes)
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:105
bool isPntInPolygon(const GeoLib::Point &pnt) const
Definition: Polygon.cpp:67
Definition of the Mesh class.
MeshLib::Mesh * readMeshFromFile(const std::string &file_name)
std::size_t getID() const
Definition: Point3dWithID.h:63
Initialization and shutting down of the logog library.
Definition: LogogSetup.h:24
static std::vector< double > getSurfaceAreaForNodes(const MeshLib::Mesh &mesh)
Returns a vector of the areas assigned to each node on a surface mesh.
Definition of the MeshSurfaceExtraction class.
Definition of readMeshFromFile function.
BASELIB_EXPORT const std::string git_describe
static const double p
Definition of the Polygon class.
std::string extractPath(std::string const &pathname)
Definition: FileTools.cpp:158
static MeshLib::Mesh * getMeshSurface(const MeshLib::Mesh &subsfc_mesh, const MathLib::Vector3 &dir, double angle, std::string const &subsfc_node_id_prop_name="", std::string const &subsfc_element_id_prop_name="", std::string const &face_id_prop_name="")
The class TemplateVec takes a unique name and manages a std::vector of pointers to data elements of t...
Definition: TemplateVec.h:40
#define OGS_FATAL(fmt,...)
Definition: Error.h:71
Filename manipulation routines.
Definition of the Vector3 class. From: http://www.strout.net/info/coding/classlib/intro.html with modifications to derive from TemplatePoint.
int main(int argc, char *argv[])
static std::vector< MeshLib::Node * > getSurfaceNodes(const MeshLib::Mesh &mesh, const MathLib::Vector3 &dir, double angle)
Returns the surface nodes of a mesh.