From 0881e4b10e2b9d8803ebf16fc02a3414214eacfc Mon Sep 17 00:00:00 2001
From: Karsten Rink <karsten.rink@ufz.de>
Date: Thu, 13 Sep 2012 15:40:40 +0200
Subject: [PATCH] getSurfaceMesh is working now; getSurfaceNodes is working
now; meshes are now tested for unused nodes (such nodes are removed upon
loading)
---
FileIO/Legacy/MeshIO.cpp | 2 +-
FileIO/XmlIO/VTKInterface.cpp | 35 ++++-
FileIO/XmlIO/VTKInterface.h | 2 +
Gui/DataView/DirectConditionGenerator.cpp | 3 +-
Gui/DataView/MshEditDialog.cpp | 2 +-
Gui/DataView/MshLayerMapper.cpp | 3 +-
Gui/mainwindow.cpp | 11 +-
MeshLib/Elements/Element.h | 4 +-
MeshLib/Mesh.cpp | 79 ++++++----
MeshLib/Mesh.h | 3 +
MeshLib/MshEditor.cpp | 172 +++++++++++++---------
MeshLib/MshEditor.h | 13 +-
MeshLib/Node.h | 5 +-
13 files changed, 219 insertions(+), 115 deletions(-)
diff --git a/FileIO/Legacy/MeshIO.cpp b/FileIO/Legacy/MeshIO.cpp
index 85f8acc206c..082182c6771 100644
--- a/FileIO/Legacy/MeshIO.cpp
+++ b/FileIO/Legacy/MeshIO.cpp
@@ -36,7 +36,7 @@ MeshIO::MeshIO()
MeshLib::Mesh* MeshIO::loadMeshFromFile(const std::string& file_name)
{
- std::cout << "Reading OGS legacy mesh ... ";
+ std::cout << "Reading OGS legacy mesh ... " << std::endl;
std::ifstream in (file_name.c_str(),std::ios::in);
if (!in.is_open())
diff --git a/FileIO/XmlIO/VTKInterface.cpp b/FileIO/XmlIO/VTKInterface.cpp
index d706be51d0f..5a02297e5b7 100644
--- a/FileIO/XmlIO/VTKInterface.cpp
+++ b/FileIO/XmlIO/VTKInterface.cpp
@@ -67,8 +67,22 @@ MeshLib::Mesh* VTKInterface::readVTUFile(const std::string &file_name)
rapidxml::xml_node<>* root_node (doc.first_node());
if (isVTKUnstructuredGrid(root_node))
{
+ bool is_compressed(false);
//check if content is compressed
const rapidxml::xml_attribute<>* compressor (root_node->first_attribute("compressor"));
+ if (compressor )
+ {
+ if (std::string(compressor->value()).compare("vtkZLibDataCompressor") == 0)
+ {
+ is_compressed = true;
+ uncompressData(root_node);
+ }
+ else
+ {
+ std::cout << "VTKInterface::readVTUFile() - Unknown compression method." << std::endl;
+ return NULL;
+ }
+ }
//skip to <Piece>-tag and start parsing content
const rapidxml::xml_node<>* piece_node (doc.first_node()->first_node()->first_node("Piece"));
@@ -132,7 +146,14 @@ MeshLib::Mesh* VTKInterface::readVTUFile(const std::string &file_name)
if (connectivity_node && celltype_node)
{
- if (std::string(celltype_node->first_attribute("format")->value()).compare("ascii") == 0)
+ const std::string format = std::string(celltype_node->first_attribute("format")->value());
+ if (format.compare("ascii") == 0)
+ {
+ std::stringstream iss (celltype_node->value());
+ for(unsigned i=0; i<nElems; i++)
+ iss >> cell_types[i];
+ }
+ else if (format.compare("appended") == 0)
{
std::stringstream iss (celltype_node->value());
for(unsigned i=0; i<nElems; i++)
@@ -242,6 +263,18 @@ bool VTKInterface::isVTKUnstructuredGrid(const rapidxml::xml_node<>* vtk_root)
return false;
}
+unsigned char* VTKInterface::uncompressData(const rapidxml::xml_node<>* node)
+{
+ rapidxml::xml_node<>* data_node = node->first_node("AppendedData");
+ char* compressed_data (NULL);
+ if (data_node)
+ compressed_data = data_node->value();
+
+ return NULL;
+}
+
+
+
int VTKInterface::write(std::ostream& stream)
{
//if (this->_export_name.empty())
diff --git a/FileIO/XmlIO/VTKInterface.h b/FileIO/XmlIO/VTKInterface.h
index 1228d5791f9..6127937d165 100644
--- a/FileIO/XmlIO/VTKInterface.h
+++ b/FileIO/XmlIO/VTKInterface.h
@@ -76,6 +76,8 @@ private:
/// Construct an Element-object from the data given to the method and the data at the current stream position.
static MeshLib::Element* readElement(std::stringstream &iss, const std::vector<MeshLib::Node*> &nodes, unsigned material, unsigned type);
+ static unsigned char* uncompressData(const rapidxml::xml_node<>* node);
+
bool _use_compressor;
};
}
diff --git a/Gui/DataView/DirectConditionGenerator.cpp b/Gui/DataView/DirectConditionGenerator.cpp
index 59116f8fcba..5f3c217b2b0 100644
--- a/Gui/DataView/DirectConditionGenerator.cpp
+++ b/Gui/DataView/DirectConditionGenerator.cpp
@@ -33,7 +33,8 @@ const std::vector< std::pair<size_t,double> >& DirectConditionGenerator::directT
return _direct_values;
}
- const std::vector<GeoLib::PointWithID*> surface_nodes ( MeshLib::MshEditor::getSurfaceNodes(mesh) );
+ const double dir[3] = {0,0,1};
+ const std::vector<GeoLib::PointWithID*> surface_nodes ( MeshLib::MshEditor::getSurfaceNodes(mesh, dir) );
//std::vector<MeshLib::CNode*> nodes = mesh.nod_vector;
const size_t nNodes(surface_nodes.size());
_direct_values.reserve(nNodes);
diff --git a/Gui/DataView/MshEditDialog.cpp b/Gui/DataView/MshEditDialog.cpp
index 6cfcaae958b..bbbfe882122 100644
--- a/Gui/DataView/MshEditDialog.cpp
+++ b/Gui/DataView/MshEditDialog.cpp
@@ -115,7 +115,7 @@ void MshEditDialog::accept()
all_paths_set = false;
else
{
- for (unsigned i=1; i<_labels.size(); i++)
+ for (int i=1; i<_labels.size(); i++)
if (_edits[i]->text().length()==0)
all_paths_set = false;
}
diff --git a/Gui/DataView/MshLayerMapper.cpp b/Gui/DataView/MshLayerMapper.cpp
index 6716b972e47..ba67a7a6466 100644
--- a/Gui/DataView/MshLayerMapper.cpp
+++ b/Gui/DataView/MshLayerMapper.cpp
@@ -187,8 +187,7 @@ int MshLayerMapper::LayerMapping(MeshLib::Mesh* new_mesh, const std::string &ras
{
std::cout << "Warning: Removing " << noData_nodes.size()
<< " mesh nodes at NoData values ... " << std::endl;
- MeshLib::MshEditor msh_editor;
- MeshLib::Mesh* red_mesh = msh_editor.removeMeshNodes(new_mesh, noData_nodes);
+ MeshLib::Mesh* red_mesh = MeshLib::MshEditor::removeMeshNodes(new_mesh, noData_nodes);
if (new_mesh->getNElements() == 0)
{
delete new_mesh;
diff --git a/Gui/mainwindow.cpp b/Gui/mainwindow.cpp
index d2aae515435..78a96fa40b3 100644
--- a/Gui/mainwindow.cpp
+++ b/Gui/mainwindow.cpp
@@ -1111,7 +1111,16 @@ void MainWindow::showVisalizationPrefsDialog()
void MainWindow::FEMTestStart()
{
- //std::string name ("Test");
+ const double dir[3] = {0, 0, 1};
+ const MeshLib::Mesh* mesh = this->_project.getMesh("Ammer-Homogen100m-Final");
+ //_meshModels->addMesh( MeshLib::MshEditor::getMeshSurface(*mesh, dir) );
+ std::vector<GeoLib::PointWithID*> pnts (MeshLib::MshEditor::getSurfaceNodes(*mesh, dir));
+ std::vector<GeoLib::Point*> *sfcpnts = new std::vector<GeoLib::Point*>(pnts.size());
+ for (unsigned i=0; i<pnts.size(); ++i)
+ (*sfcpnts)[i] = pnts[i];
+ std::string name("SurfacePoints");
+ this->_geoModels->addPointVec(sfcpnts, name);
+ //std::string name ("test");
//_meshModels->addMesh(MshEditor::getMeshSurface(*_project.getMesh("Ammer-Homogen100m-Final")), name);
/*
diff --git a/MeshLib/Elements/Element.h b/MeshLib/Elements/Element.h
index 2a6b0fbce2a..7ac7218b746 100644
--- a/MeshLib/Elements/Element.h
+++ b/MeshLib/Elements/Element.h
@@ -26,10 +26,8 @@ class Node;
*/
class Element
{
- /* friend functions: */
+ /* friend classes */
friend class Mesh;//void Mesh::setElementInformationForNodes();
- //friend void Mesh::addElement(Element*);
-
public:
/// Compute the minimum and maximum squared edge length for this element
diff --git a/MeshLib/Mesh.cpp b/MeshLib/Mesh.cpp
index 94533bbc4ed..10b6c3532a3 100644
--- a/MeshLib/Mesh.cpp
+++ b/MeshLib/Mesh.cpp
@@ -39,6 +39,7 @@ Mesh::Mesh(const std::string &name, const std::vector<Node*> &nodes, const std::
//this->setNodesConnectedByEdges();
//this->setNodesConnectedByElements();
this->setElementsConnectedToElements();
+ this->removeUnusedMeshNodes();
}
Mesh::Mesh(const Mesh &mesh)
@@ -46,16 +47,16 @@ Mesh::Mesh(const Mesh &mesh)
{
const std::vector<Node*> nodes (mesh.getNodes());
const size_t nNodes (nodes.size());
- for (unsigned i=0; i<nNodes; i++)
+ for (unsigned i=0; i<nNodes; ++i)
_nodes[i] = new Node(*nodes[i]);
const std::vector<Element*> elements (mesh.getElements());
const size_t nElements (elements.size());
- for (unsigned i=0; i<nElements; i++)
+ for (unsigned i=0; i<nElements; ++i)
{
const size_t nElemNodes = elements[i]->getNNodes();
_elements[i] = elements[i]->clone();
- for (unsigned j=0; j<nElemNodes; j++)
+ for (unsigned j=0; j<nElemNodes; ++j)
_elements[i]->_nodes[j] = _nodes[elements[i]->getNode(j)->getID()];
}
@@ -69,11 +70,11 @@ Mesh::Mesh(const Mesh &mesh)
Mesh::~Mesh()
{
const size_t nElements (_elements.size());
- for (size_t i=0; i<nElements; i++)
+ for (size_t i=0; i<nElements; ++i)
delete _elements[i];
const size_t nNodes (_nodes.size());
- for (size_t i=0; i<nNodes; i++)
+ for (size_t i=0; i<nNodes; ++i)
delete _nodes[i];
}
@@ -84,10 +85,10 @@ void Mesh::makeNodesUnique()
//replace node pointers in elements
unsigned nElements (_elements.size());
- for (unsigned i=0; i<nElements; i++)
+ for (unsigned i=0; i<nElements; ++i)
{
unsigned nNodes (_elements[i]->getNNodes());
- for (unsigned j=0; j<nNodes; j++)
+ for (unsigned j=0; j<nNodes; ++j)
_elements[i]->getNodeIndex(j);
}
@@ -106,21 +107,21 @@ void Mesh::addElement(Element* elem)
// add element information to nodes
unsigned nNodes (elem->getNNodes());
- for (unsigned i=0; i<nNodes; i++)
+ for (unsigned i=0; i<nNodes; ++i)
elem->_nodes[i]->addElement(elem);
}
void Mesh::resetNodeIDs()
{
const size_t nNodes (this->_nodes.size());
- for (unsigned i=0; i<nNodes; i++)
+ for (unsigned i=0; i<nNodes; ++i)
_nodes[i]->setID(i);
}
void Mesh::setDimension()
{
const size_t nElements (_elements.size());
- for (unsigned i=0; i<nElements; i++)
+ for (unsigned i=0; i<nElements; ++i)
if (_elements[i]->getDimension() > _mesh_dimension)
_mesh_dimension = _elements[i]->getDimension();
}
@@ -128,19 +129,25 @@ void Mesh::setDimension()
void Mesh::setElementsConnectedToNodes()
{
const size_t nElements (_elements.size());
- for (unsigned i=0; i<nElements; i++)
+ for (unsigned i=0; i<nElements; ++i)
{
MeshLib::Element* element = _elements[i];
const unsigned nNodes (element->getNNodes());
- for (unsigned j=0; j<nNodes; j++)
+ for (unsigned j=0; j<nNodes; ++j)
element->_nodes[j]->addElement(element);
}
-#ifdef NDEBUG
+#ifndef NDEBUG
// search for nodes that are not part of any element
+ unsigned count(0);
const size_t nNodes (_nodes.size());
- for (unsigned i=0; i<nNodes; i++)
+ for (unsigned i=0; i<nNodes; ++i)
if (_nodes[i]->getNElements() == 0)
- WARN ("Warning: Node %d is not part of any element.", i);
+ {
+ WARN ("Node %d is not part of any element.", i);
+ ++count;
+ }
+ if (count)
+ WARN ("%d unused mesh nodes found.", count);
#endif
}
@@ -164,7 +171,7 @@ void Mesh::setElementsConnectedToElements()
#else
unsigned m(0);
#endif
- for (m=0; m<nElements; m++)
+ for (m=0; m<nElements; ++m)
{
// create vector with all elements connected to current element (includes lots of doubles!)
std::vector<Element*> neighbors;
@@ -172,7 +179,7 @@ void Mesh::setElementsConnectedToElements()
if (element->getType() != MshElemType::EDGE)
{
const size_t nNodes (element->getNNodes());
- for (unsigned n(0); n<nNodes; n++)
+ for (unsigned n(0); n<nNodes; ++n)
{
std::vector<Element*> const& conn_elems ((element->getNode(n)->getElements()));
neighbors.insert(neighbors.end(), conn_elems.begin(), conn_elems.end());
@@ -180,7 +187,7 @@ void Mesh::setElementsConnectedToElements()
const unsigned nNeighbors ( neighbors.size() );
- for (unsigned i(0); i<nNeighbors; i++)
+ for (unsigned i(0); i<nNeighbors; ++i)
{
if (element->addNeighbor(neighbors[i]) && neighbors[i]->getType() != MshElemType::EDGE)
{
@@ -194,21 +201,21 @@ void Mesh::setElementsConnectedToElements()
void Mesh::setNodesConnectedByEdges()
{
const size_t nNodes (this->_nodes.size());
- for (unsigned i=0; i<nNodes; i++)
+ for (unsigned i=0; i<nNodes; ++i)
{
MeshLib::Node* node (_nodes[i]);
std::vector<MeshLib::Node*> conn_set;
const std::vector<MeshLib::Element*> &conn_elems (node->getElements());
const size_t nConnElems (conn_elems.size());
- for (unsigned j=0; j<nConnElems; j++)
+ for (unsigned j=0; j<nConnElems; ++j)
{
const unsigned idx (conn_elems[j]->getNodeIDinElement(node));
const unsigned nElemNodes (conn_elems[j]->getNNodes());
- for (unsigned k(0); k<nElemNodes; k++)
+ for (unsigned k(0); k<nElemNodes; ++k)
{
bool is_in_vector (false);
const size_t nConnNodes (conn_set.size());
- for (unsigned l(0); l<nConnNodes; l++)
+ for (unsigned l(0); l<nConnNodes; ++l)
if (conn_elems[j]->getNode(k) == conn_set[l])
is_in_vector = true;
if (is_in_vector) continue;
@@ -223,22 +230,22 @@ void Mesh::setNodesConnectedByEdges()
void Mesh::setNodesConnectedByElements()
{
const size_t nNodes (this->_nodes.size());
- for (unsigned i=0; i<nNodes; i++)
+ for (unsigned i=0; i<nNodes; ++i)
{
MeshLib::Node* node (_nodes[i]);
std::vector<MeshLib::Node*> conn_vec;
const std::vector<MeshLib::Element*> &conn_elems (node->getElements());
const size_t nConnElems (conn_elems.size());
- for (unsigned j=0; j<nConnElems; j++)
+ for (unsigned j=0; j<nConnElems; ++j)
{
const unsigned nElemNodes (conn_elems[j]->getNNodes());
- for (unsigned k(0); k<nElemNodes; k++)
+ for (unsigned k(0); k<nElemNodes; ++k)
{
bool is_in_vector (false);
const MeshLib::Node* c_node (conn_elems[j]->getNode(k));
if (c_node == node) continue;
const size_t nConnNodes (conn_vec.size());
- for (unsigned l(0); l<nConnNodes; l++)
+ for (unsigned l(0); l<nConnNodes; ++l)
if (c_node == conn_vec[l])
is_in_vector = true;
if (!is_in_vector)
@@ -249,5 +256,25 @@ void Mesh::setNodesConnectedByElements()
}
}
+void Mesh::removeUnusedMeshNodes()
+{
+ unsigned count(0);
+ std::vector<MeshLib::Node*>::iterator it (this->_nodes.begin());
+ for (it; it != this->_nodes.end();)
+ {
+ if ((*it)->getNElements() == 0)
+ {
+ it = this->_nodes.erase(it);
+ ++count;
+ }
+ else ++it;
+ }
+ if (count)
+ {
+ std::cout << "Removed " << count << " unused mesh nodes." << std::endl;
+ this->resetNodeIDs();
+ }
+}
+
}
diff --git a/MeshLib/Mesh.h b/MeshLib/Mesh.h
index 179884d7b77..c499bf8aa3f 100644
--- a/MeshLib/Mesh.h
+++ b/MeshLib/Mesh.h
@@ -97,6 +97,9 @@ protected:
/// Checks the coordinates of all mesh nodes and removes identical nodes. Elements are adapted accordingly.
void makeNodesUnique();
+ /// Removes nodes that are not part of any element.
+ void removeUnusedMeshNodes();
+
/// Sets the dimension of the mesh.
void setDimension();
diff --git a/MeshLib/MshEditor.cpp b/MeshLib/MshEditor.cpp
index 6e56aaafd7b..9e91048d7d5 100644
--- a/MeshLib/MshEditor.cpp
+++ b/MeshLib/MshEditor.cpp
@@ -34,14 +34,14 @@ void MshEditor::getSurfaceAreaForNodes(const MeshLib::Mesh* mesh, std::vector<do
// for each node, a vector containing all the element idget every element
std::vector<MeshLib::Node*> nodes = mesh->getNodes();
const size_t nNodes ( mesh->getNNodes() );
- for (size_t n=0; n<nNodes; n++)
+ for (size_t n=0; n<nNodes; ++n)
{
double node_area (0);
std::vector<MeshLib::Element*> conn_elems = nodes[n]->getElements();
const size_t nConnElems (conn_elems.size());
- for (size_t i=0; i<nConnElems;i++)
+ for (size_t i=0; i<nConnElems; ++i)
{
const MeshLib::Element* elem (conn_elems[i]);
const unsigned nElemParts = (elem->getType() == MshElemType::TRIANGLE) ? 3 : 4;
@@ -67,12 +67,12 @@ MeshLib::Mesh* MshEditor::removeMeshNodes(MeshLib::Mesh* mesh,
// delete nodes and their connected elements and replace them with null pointers
const size_t delNodes = nodes.size();
std::vector<MeshLib::Node*> mesh_nodes = new_mesh->getNodes();
- for (size_t i = 0; i < delNodes; i++)
+ for (size_t i = 0; i < delNodes; ++i)
{
const MeshLib::Node* node = new_mesh->getNode(i);
std::vector<MeshLib::Element*> conn_elems = node->getElements();
- for (size_t j = 0; j < conn_elems.size(); j++)
+ for (size_t j = 0; j < conn_elems.size(); ++j)
{
delete conn_elems[j];
conn_elems[j] = NULL;
@@ -85,7 +85,7 @@ MeshLib::Mesh* MshEditor::removeMeshNodes(MeshLib::Mesh* mesh,
const size_t nNodes = new_mesh->getNNodes();
std::vector<int> id_map(nNodes, -1);
size_t count(0);
- for (size_t i = 0; i < nNodes; i++)
+ for (size_t i = 0; i < nNodes; ++i)
{
if (mesh_nodes[i])
{
@@ -115,81 +115,79 @@ MeshLib::Mesh* MshEditor::removeMeshNodes(MeshLib::Mesh* mesh,
return new_mesh;
}
-std::vector<GeoLib::PointWithID*> MshEditor::getSurfaceNodes(const MeshLib::Mesh &mesh)
+MeshLib::Mesh* MshEditor::getMeshSurface(const MeshLib::Mesh &mesh, const double* dir)
{
- /* TODO6
- INFO ("Extracting surface nodes...");
- // Sort points lexicographically
- size_t nNodes (mesh.nod_vector.size());
- std::vector<GeoLib::PointWithID*> nodes;
- std::vector<size_t> perm;
- for (size_t j(0); j<nNodes; j++)
- {
- nodes.push_back(new GeoLib::PointWithID(mesh.nod_vector[j]->getData(), j));
- perm.push_back(j);
- }
- Quicksort<GeoLib::PointWithID*> (nodes, 0, nodes.size(), perm);
-
- // Extract surface points
- double eps (std::numeric_limits<double>::epsilon());
- */
- std::vector<GeoLib::PointWithID*> surface_pnts;
- /* TODO6
- for (size_t k(1); k < nNodes; k++)
+ INFO ("Extracting mesh surface...");
+ const std::vector<MeshLib::Element*> all_elements (mesh.getElements());
+ const std::vector<MeshLib::Node*> all_nodes (mesh.getNodes());
+
+ std::vector<MeshLib::Element*> sfc_elements;
+ get2DSurfaceElements(all_elements, sfc_elements, dir, mesh.getDimension());
+
+ std::vector<MeshLib::Node*> sfc_nodes;
+ std::vector<unsigned> node_id_map(mesh.getNNodes());
+ get2DSurfaceNodes(all_nodes, sfc_nodes, sfc_elements, node_id_map);
+
+ // create new elements vector with newly created nodes
+ const size_t nNewElements (sfc_elements.size());
+ std::vector<MeshLib::Element*> new_elements(sfc_elements.size());
+ for (unsigned i=0; i<nNewElements; ++i)
{
- const GeoLib::PointWithID& p0 (*(nodes[k - 1]));
- const GeoLib::PointWithID& p1 (*(nodes[k]));
- if (fabs (p0[0] - p1[0]) > eps || fabs (p0[1] - p1[1]) > eps)
- surface_pnts.push_back (nodes[k - 1]);
+ MeshLib::Element* elem (sfc_elements[i]);
+ if (elem->getType() == MshElemType::TRIANGLE)
+ new_elements[i] = new MeshLib::Tri(sfc_nodes[node_id_map[elem->getNode(0)->getID()]],
+ sfc_nodes[node_id_map[elem->getNode(1)->getID()]],
+ sfc_nodes[node_id_map[elem->getNode(2)->getID()]]);
+ else
+ new_elements[i] = new MeshLib::Quad(sfc_nodes[node_id_map[elem->getNode(0)->getID()]],
+ sfc_nodes[node_id_map[elem->getNode(1)->getID()]],
+ sfc_nodes[node_id_map[elem->getNode(2)->getID()]],
+ sfc_nodes[node_id_map[elem->getNode(3)->getID()]]);
+ delete sfc_elements[i];
}
- // Add last point
- surface_pnts.push_back (nodes[nNodes - 1]);
- */
- return surface_pnts;
+ return new Mesh("SurfaceMesh", sfc_nodes, new_elements);
}
-MeshLib::Mesh* MshEditor::getMeshSurface(const MeshLib::Mesh &mesh, const double* dir)
+void MshEditor::get2DSurfaceElements(const std::vector<MeshLib::Element*> &all_elements, std::vector<MeshLib::Element*> &sfc_elements, const double* dir, unsigned mesh_dimension)
{
- INFO ("Extracting mesh surface...");
+ bool complete_surface (true);
+ if (dir)
+ complete_surface = ((dir[0]*dir[0] + dir[1]*dir[1] + dir[2]*dir[2]) == 0) ? true : false;
- const std::vector<MeshLib::Element*> elements = mesh.getElements();
- std::vector<MeshLib::Element*> new_elements;
- const size_t nElements (mesh.getNElements());
+ const size_t nElements (all_elements.size());
- bool complete_surface = ((dir[0]*dir[0] + dir[1]*dir[1] + dir[2]*dir[2]) == 0) ? true : false;
-
- // 2D meshes
- if (mesh.getDimension() == 2 )
+ if (mesh_dimension == 2 )
{
- if (complete_surface) return new MeshLib::Mesh(mesh); // simply copy the mesh
- else // check only surface normal directions of all elements
+ for (unsigned i=0; i<nElements; ++i)
{
- for (unsigned i=0; i<nElements; i++)
+ if (complete_surface)
+ sfc_elements.push_back(all_elements[i]);
+ else
{
- MeshLib::Face* face = dynamic_cast<MeshLib::Face*>(elements[i]);
+ MeshLib::Face* face = dynamic_cast<MeshLib::Face*>(all_elements[i]);
double normal[3];
face->getSurfaceNormal(normal);
if (MathLib::scpr(normal, dir, 3) > 0)
- new_elements.push_back(static_cast<MeshLib::Element*>(face));
+ sfc_elements.push_back(static_cast<MeshLib::Element*>(face));
}
}
}
- // 3D meshes
- else if (mesh.getDimension() == 3) //
+ else if (mesh_dimension == 3)
{
- for (unsigned i=0; i<nElements; i++)
+ for (unsigned i=0; i<nElements; ++i)
{
- if (const MeshLib::Cell* cell = dynamic_cast<MeshLib::Cell*>(elements[i]))
+ if (all_elements[i]->getDimension()==3)
{
+ const MeshLib::Cell* cell = static_cast<MeshLib::Cell*>(all_elements[i]);
if (cell->isOnSurface())
{
const unsigned nFaces (cell->getNFaces());
- for (unsigned j=0; j<nFaces; j++)
+ for (unsigned j=0; j<nFaces; ++j)
{
- if (cell->getNeighbor(i) == NULL)
+ if (cell->getNeighbor(j) == NULL)
{
- const MeshLib::Face* face = static_cast<const MeshLib::Face*>(cell->getFace(i));
+ const MeshLib::Face* face = static_cast<const MeshLib::Face*>(cell->getFace(j));
if (!complete_surface)
{
double normal[3];
@@ -199,35 +197,67 @@ MeshLib::Mesh* MshEditor::getMeshSurface(const MeshLib::Mesh &mesh, const double
}
if (face->getType() == MshElemType::TRIANGLE)
- new_elements.push_back(new MeshLib::Tri(*static_cast<const MeshLib::Tri*>(face)));
+ sfc_elements.push_back(new MeshLib::Tri(*static_cast<const MeshLib::Tri*>(face)));
else
- new_elements.push_back(new MeshLib::Quad(*static_cast<const MeshLib::Quad*>(face)));
+ sfc_elements.push_back(new MeshLib::Quad(*static_cast<const MeshLib::Quad*>(face)));
}
}
}
}
}
+ }
+}
- // now copy nodes
- const size_t nNewElements (new_elements.size());
- std::vector<const MeshLib::Node*> tmp_nodes(mesh.getNNodes(), NULL);
- const size_t nNodes (tmp_nodes.size());
- for (unsigned i=0; i<nNewElements; i++)
+void MshEditor::get2DSurfaceNodes(const std::vector<MeshLib::Node*> &all_nodes, std::vector<MeshLib::Node*> &sfc_nodes, const std::vector<MeshLib::Element*> &sfc_elements, std::vector<unsigned> &node_id_map)
+{
+ const size_t nNewElements (sfc_elements.size());
+ std::vector<const MeshLib::Node*> tmp_nodes(all_nodes.size(), NULL);
+ const size_t nNodes (tmp_nodes.size());
+ for (unsigned i=0; i<nNewElements; ++i)
+ {
+ const MeshLib::Element* elem (sfc_elements[i]);
+ for (unsigned j=0; j<elem->getNNodes(); ++j)
{
- const MeshLib::Element* elem (new_elements[i]);
- for (unsigned j=0; j<elem->getNNodes(); j++)
- {
- const MeshLib::Node* node (elem->getNode(i));
- tmp_nodes[node->getID()] = node;
- }
+ const MeshLib::Node* node (elem->getNode(j));
+ tmp_nodes[node->getID()] = node;
+ }
+ }
+ for (unsigned i=0; i<nNodes; ++i)
+ {
+ if (tmp_nodes[i])
+ {
+ node_id_map[i] = sfc_nodes.size();
+ sfc_nodes.push_back(new MeshLib::Node(tmp_nodes[i]->getCoords(), tmp_nodes[i]->getID()));
}
- std::vector<MeshLib::Node*> new_nodes;
- for (unsigned i=0; i<nNodes; i++)
- if (tmp_nodes[i])
- new_nodes.push_back(new MeshLib::Node(tmp_nodes[i]->getCoords()));
}
+}
+std::vector<GeoLib::PointWithID*> MshEditor::getSurfaceNodes(const MeshLib::Mesh &mesh, const double *dir)
+{
+ INFO ("Extracting surface nodes...");
+ const std::vector<MeshLib::Element*> all_elements (mesh.getElements());
+ const std::vector<MeshLib::Node*> all_nodes (mesh.getNodes());
+
+ std::vector<MeshLib::Element*> sfc_elements;
+ get2DSurfaceElements(all_elements, sfc_elements, dir, mesh.getDimension());
+
+ std::vector<MeshLib::Node*> sfc_nodes;
+ std::vector<unsigned> node_id_map(mesh.getNNodes());
+ get2DSurfaceNodes(all_nodes, sfc_nodes, sfc_elements, node_id_map);
+
+ const unsigned nElements (sfc_elements.size());
+ for (unsigned i=0; i<nElements; ++i)
+ delete sfc_elements[i];
+ const size_t nNodes (sfc_nodes.size());
+ std::vector<GeoLib::PointWithID*> surface_pnts(nNodes);
+ for (unsigned i=0; i<nNodes; ++i)
+ {
+ surface_pnts[i] = new GeoLib::PointWithID(sfc_nodes[i]->getCoords(), sfc_nodes[i]->getID());
+ delete sfc_nodes[i];
+ }
+ return surface_pnts;
}
+
} // end namespace MeshLib
diff --git a/MeshLib/MshEditor.h b/MeshLib/MshEditor.h
index 787593abd2f..c373feca762 100644
--- a/MeshLib/MshEditor.h
+++ b/MeshLib/MshEditor.h
@@ -25,6 +25,7 @@ namespace MeshLib {
// forward declarations
class Mesh;
class Element;
+class Node;
/**
* \brief A set of tools for manipulating existing meshes
@@ -32,9 +33,6 @@ class Element;
class MshEditor
{
public:
- MshEditor() {}
- ~MshEditor() {}
-
/// Returns the area assigned to each node on a surface mesh.
static void getSurfaceAreaForNodes(const MeshLib::Mesh* mesh, std::vector<double> &node_area_vec);
@@ -42,12 +40,17 @@ public:
static MeshLib::Mesh* removeMeshNodes(MeshLib::Mesh* mesh, const std::vector<std::size_t> &nodes);
/// Returns the surface nodes of a layered mesh.
- static std::vector<GeoLib::PointWithID*> getSurfaceNodes(const MeshLib::Mesh &mesh);
+ static std::vector<GeoLib::PointWithID*> getSurfaceNodes(const MeshLib::Mesh &mesh, const double* dir = NULL);
/// Returns the 2d-element mesh representing the surface of the given layered mesh.
- static MeshLib::Mesh* getMeshSurface(const MeshLib::Mesh &mesh, const double* dir);
+ static MeshLib::Mesh* getMeshSurface(const MeshLib::Mesh &mesh, const double* dir = NULL);
private:
+ /// Functionality needed for getSurfaceNodes() and getMeshSurface()
+ static void get2DSurfaceElements(const std::vector<MeshLib::Element*> &all_elements, std::vector<MeshLib::Element*> &sfc_elements, const double* dir, unsigned mesh_dimension);
+
+ /// Functionality needed for getSurfaceNodes() and getMeshSurface()
+ static void get2DSurfaceNodes(const std::vector<MeshLib::Node*> &all_nodes, std::vector<MeshLib::Node*> &sfc_nodes, const std::vector<MeshLib::Element*> &sfc_elements, std::vector<unsigned> &node_id_map);
};
} // end namespace MeshLib
diff --git a/MeshLib/Node.h b/MeshLib/Node.h
index e329d61fb28..fb9009b74ad 100644
--- a/MeshLib/Node.h
+++ b/MeshLib/Node.h
@@ -34,9 +34,8 @@ class Node : public GeoLib::PointWithID
{
/* friend functions: */
friend MeshLib::Mesh* MshEditor::removeMeshNodes(MeshLib::Mesh* mesh, const std::vector<std::size_t> &nodes);
- friend int MshLayerMapper::LayerMapping(MeshLib::Mesh* msh, const std::string &rasterfile,
- const std::size_t nLayers, const std::size_t layer_id,
- bool removeNoDataValues);
+ friend int MshLayerMapper::LayerMapping(MeshLib::Mesh* msh, const std::string &rasterfile, const std::size_t nLayers,
+ const std::size_t layer_id, bool removeNoDataValues);
/* friend classes: */
friend class Mesh;
friend class MeshCoarsener;
--
GitLab