diff --git a/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp b/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp index ced2bb7fa738c4b01172e2b3b71bdf61d91411d5..32606e29ef678a4c10ceb70c32b29f962405aab5 100644 --- a/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp +++ b/Applications/Utils/SimpleMeshCreation/createMeshElemPropertiesFromASCRaster.cpp @@ -30,7 +30,8 @@ #include "MathLib/MathTools.h" #include "MeshLib/MeshGenerators/VtkMeshConverter.h" -#include "MeshLib/MeshGenerators/ConvertRasterToMesh.h" +//#include "MeshLib/MeshGenerators/ConvertRasterToMesh.h" +#include "MeshLib/MeshGenerators/VtkMeshConverter.h" #include "MeshLib/Elements/Element.h" #include "MeshLib/Mesh.h" #include "MeshLib/MeshEditing/Mesh2MeshPropertyInterpolation.h" @@ -129,26 +130,26 @@ int main (int argc, char* argv[]) // read raster and if required manipulate it auto raster = std::unique_ptr<GeoLib::Raster>( FileIO::AsciiRasterInterface::getRasterFromASCFile(raster_arg.getValue())); + GeoLib::RasterHeader header (raster->getHeader()); if (refinement_arg.getValue() > 1) { raster->refineRaster(refinement_arg.getValue()); if (refinement_raster_output_arg.getValue()) { // write new asc file std::string new_raster_fname (BaseLib::dropFileExtension( raster_arg.getValue())); - new_raster_fname += "-" + std::to_string(raster->getNRows()) + "x" + - std::to_string(raster->getNCols()) + ".asc"; + new_raster_fname += "-" + std::to_string(header.n_rows) + "x" + + std::to_string(header.n_cols) + ".asc"; FileIO::AsciiRasterInterface::writeRasterAsASC(*raster, new_raster_fname); } } // put raster data in a std::vector GeoLib::Raster::const_iterator raster_it(raster->begin()); - std::size_t n_cols(raster->getNCols()), n_rows(raster->getNRows()); - std::size_t size(n_cols * n_rows); + std::size_t size(header.n_cols * header.n_rows); std::vector<double> src_properties(size); - for (unsigned row(0); row<n_rows; row++) { - for (unsigned col(0); col<n_cols; col++) { - src_properties[row * n_cols + col] = *raster_it; + for (unsigned row(0); row<header.n_rows; row++) { + for (unsigned col(0); col<header.n_cols; col++) { + src_properties[row * header.n_cols + col] = *raster_it; ++raster_it; } } @@ -160,10 +161,12 @@ int main (int argc, char* argv[]) INFO("Variance of source: %f.", var); } + std::unique_ptr<MeshLib::Mesh> src_mesh (MeshLib::VtkMeshConverter::convertRasterToMesh(*raster, MeshLib::MeshElemType::QUAD, MeshLib::UseIntensityAs::DATAVECTOR)); +/* std::unique_ptr<MeshLib::Mesh> src_mesh(MeshLib::ConvertRasterToMesh( *raster, MeshLib::MeshElemType::QUAD, MeshLib::UseIntensityAs::DATAVECTOR).execute()); - +*/ std::vector<std::size_t> src_perm(size); std::iota(src_perm.begin(), src_perm.end(), 0); BaseLib::quicksort<double>(src_properties, 0, size, src_perm); diff --git a/MeshLib/MeshGenerators/ConvertRasterToMesh.cpp b/MeshLib/MeshGenerators/ConvertRasterToMesh.cpp deleted file mode 100644 index f7c6c076aa9b4d146ce63a668d31f13fcea39d81..0000000000000000000000000000000000000000 --- a/MeshLib/MeshGenerators/ConvertRasterToMesh.cpp +++ /dev/null @@ -1,165 +0,0 @@ -/** - * @file ConvertRasterToMesh.cpp - * @author Thomas Fischer - * @date Nov 14, 2012 - * @brief Implementation of the ConvertRasterToMesh class. - * - * @copyright - * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org) - * Distributed under a Modified BSD License. - * See accompanying file LICENSE.txt or - * http://www.opengeosys.org/project/license - */ - -#include "ConvertRasterToMesh.h" - -#include "MeshLib/Node.h" -#include "MeshLib/Elements/Tri.h" -#include "MeshLib/Elements/Quad.h" - -namespace MeshLib { - -ConvertRasterToMesh::ConvertRasterToMesh(GeoLib::Raster const& raster, MeshElemType elem_type, - UseIntensityAs intensity_type) : - _raster(raster), _elem_type(elem_type), _intensity_type(intensity_type) -{} - -ConvertRasterToMesh::~ConvertRasterToMesh() -{} - -MeshLib::Mesh* ConvertRasterToMesh::execute() const -{ - GeoLib::RasterHeader const& header (_raster.getHeader()); - const std::size_t height(header.n_rows+1); - const std::size_t width(header.n_cols+1); - double* pix_vals(new double[height*width]); - bool* vis_nodes(new bool[height*width]); - - // determine a valid value for substitution of no data values - double substitution(getExistingValue(_raster.begin(), _raster.end())); - - // fill first row with non visual nodes - for (std::size_t j = 0; j < width; j++) { - pix_vals[j] = 0; - vis_nodes[j] = false; - } - - GeoLib::Raster::const_iterator raster_it(_raster.begin()); - for (std::size_t i = 0; i < height; ++i) { - for (std::size_t j = 0; j < width; ++j) { - const std::size_t index = (i+1) * width + j; - if (*raster_it == header.no_data) { - pix_vals[index] = substitution; - vis_nodes[index] = false; - } else { - pix_vals[index] = *raster_it; - vis_nodes[index] = true; - } - ++raster_it; - } - // fill last column with non-visual nodes - pix_vals[(i + 2) * width - 1] = 0; - vis_nodes[(i + 2) * width - 1] = false; - } - - MeshLib::Mesh* mesh = constructMesh(pix_vals, vis_nodes); - - delete [] pix_vals; - delete [] vis_nodes; - - return mesh; -} - -MeshLib::Mesh* ConvertRasterToMesh::constructMesh(const double* pix_vals, const bool* vis_nodes) const -{ - GeoLib::RasterHeader const& header (_raster.getHeader()); - const std::size_t height (header.n_rows+1); - const std::size_t width (header.n_cols+1); - std::size_t node_idx_count(0); - const double distance(header.cell_size); - - const std::size_t size(height*width); - int* node_idx_map(new int[size]); - for (std::size_t k(0); k<size; ++k) node_idx_map[k] = -1; - - std::vector<MeshLib::Node*> nodes; - std::vector<MeshLib::Element*> elements; - - for (std::size_t i = 0; i < height; i++) { - for (std::size_t j = 0; j < width; j++) { - const std::size_t index = i * width + j; - double zValue = (_intensity_type == UseIntensityAs::ELEVATION) ? - pix_vals[index] : header.origin[2]; - MeshLib::Node* node(new MeshLib::Node(header.origin[0] + (distance * j), - header.origin[1] + (distance * i), zValue)); - nodes.push_back(node); - node_idx_map[index] = node_idx_count; - node_idx_count++; - } - } - - std::vector<int> mat_ids; - // set mesh elements - for (std::size_t i = 0; i < height; i++) { - for (std::size_t j = 0; j < width; j++) { - const int index = i * width + j; - if ((node_idx_map[index] != -1) && (node_idx_map[index + 1] != -1) - && (node_idx_map[index + width] != -1) - && (node_idx_map[index + width + 1] != -1) - && (vis_nodes[index + width])) { - const int mat = (_intensity_type != UseIntensityAs::DATAVECTOR) ? 0 - : static_cast<int> (pix_vals[index + width]); - if (_elem_type == MeshElemType::TRIANGLE) { - MeshLib::Node** tri1_nodes = new MeshLib::Node*[3]; - tri1_nodes[0] = nodes[node_idx_map[index]]; - tri1_nodes[1] = nodes[node_idx_map[index + 1]]; - tri1_nodes[2] = nodes[node_idx_map[index + width]]; - - MeshLib::Node** tri2_nodes = new MeshLib::Node*[3]; - tri2_nodes[0] = nodes[node_idx_map[index + 1]]; - tri2_nodes[1] = nodes[node_idx_map[index + width + 1]]; - tri2_nodes[2] = nodes[node_idx_map[index + width]]; - - // upper left triangle - elements.push_back(new MeshLib::Tri(tri1_nodes)); - mat_ids.push_back(mat); - // lower right triangle - elements.push_back(new MeshLib::Tri(tri2_nodes)); - mat_ids.push_back(mat); - } - if (_elem_type == MeshElemType::QUAD) { - MeshLib::Node** quad_nodes = new MeshLib::Node*[4]; - quad_nodes[0] = nodes[node_idx_map[index]]; - quad_nodes[1] = nodes[node_idx_map[index + 1]]; - quad_nodes[2] = nodes[node_idx_map[index + width + 1]]; - quad_nodes[3] = nodes[node_idx_map[index + width]]; - elements.push_back(new MeshLib::Quad(quad_nodes)); - mat_ids.push_back(mat); - } - } - } - } - delete [] node_idx_map; - - // the name is only a temp-name, the name given in the dialog is set later - MeshLib::Properties properties; - boost::optional<MeshLib::PropertyVector<int>&> materials = - properties.createNewPropertyVector<int>("MaterialIDs", - MeshLib::MeshItemType::Cell); - assert(materials != boost::none); - materials->reserve(mat_ids.size()); - std::copy(mat_ids.cbegin(), mat_ids.cend(), std::back_inserter(*materials)); - return new MeshLib::Mesh("RasterDataMesh", nodes, elements, properties); -} - - -double ConvertRasterToMesh::getExistingValue(GeoLib::Raster::const_iterator first, GeoLib::Raster::const_iterator last) const -{ - double const no_data (_raster.getHeader().no_data); - for (GeoLib::Raster::const_iterator it(first); it != last; ++it) { - if (*it != no_data) - return *it; - } - return no_data; -} -} // end namespace MeshLib diff --git a/MeshLib/MeshGenerators/ConvertRasterToMesh.h b/MeshLib/MeshGenerators/ConvertRasterToMesh.h deleted file mode 100644 index 2655339e880f28fb46342c60c3d29d6e767fd757..0000000000000000000000000000000000000000 --- a/MeshLib/MeshGenerators/ConvertRasterToMesh.h +++ /dev/null @@ -1,46 +0,0 @@ -/** - * @file ConvertRasterToMesh.h - * @author Thomas Fischer - * @date Nov 14, 2012 - * @brief Definition of the ConvertRasterToMesh class. - * - * @copyright - * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org) - * Distributed under a Modified BSD License. - * See accompanying file LICENSE.txt or - * http://www.opengeosys.org/project/license - * - */ - -#ifndef CONVERTRASTERTOMESH_H_ -#define CONVERTRASTERTOMESH_H_ - -#include "GeoLib/Raster.h" -#include "MeshLib/MeshEnums.h" - -namespace MeshLib { - -// forward declaration -class Mesh; - - -/** - * Class to convert raster data into meshes. Based on Karsten Rinks algorithm. - */ -class ConvertRasterToMesh { -public: - ConvertRasterToMesh(GeoLib::Raster const& raster, MeshElemType elem_type, - UseIntensityAs intensity_type); - ~ConvertRasterToMesh(); - MeshLib::Mesh* execute() const; -private: - double getExistingValue(GeoLib::Raster::const_iterator beg, GeoLib::Raster::const_iterator last) const; - MeshLib::Mesh* constructMesh(const double* pix_vals, const bool* vis_nodes) const; - GeoLib::Raster const& _raster; - MeshElemType _elem_type; - UseIntensityAs _intensity_type; -}; - -} // end namespace MeshLib - -#endif /* CONVERTRASTERTOMESH_H_ */ diff --git a/MeshLib/MeshGenerators/VtkMeshConverter.cpp b/MeshLib/MeshGenerators/VtkMeshConverter.cpp index ae9f83c648536de2a81b4bf027f8d4e2afb6e8c3..29112e11924913554d9fc35c3e71b913d81e3e45 100644 --- a/MeshLib/MeshGenerators/VtkMeshConverter.cpp +++ b/MeshLib/MeshGenerators/VtkMeshConverter.cpp @@ -50,6 +50,13 @@ MeshLib::Element* createElementWithSameNodeOrder(const std::vector<MeshLib::Node } } +MeshLib::Mesh* VtkMeshConverter::convertRasterToMesh(GeoLib::Raster const& raster, + MeshElemType elem_type, + UseIntensityAs intensity_type) +{ + return convertImgToMesh(raster.begin(), raster.getHeader(), elem_type, intensity_type); +} + MeshLib::Mesh* VtkMeshConverter::convertImgToMesh(vtkImageData* img, const double origin[3], const double scalingFactor, diff --git a/MeshLib/MeshGenerators/VtkMeshConverter.h b/MeshLib/MeshGenerators/VtkMeshConverter.h index e98edd17f06a0a4f739be161c9455f9a3e89c8bd..422002c20fd093841f6838a5ff9729005095a7a0 100644 --- a/MeshLib/MeshGenerators/VtkMeshConverter.h +++ b/MeshLib/MeshGenerators/VtkMeshConverter.h @@ -45,6 +45,15 @@ class Properties; class VtkMeshConverter { public: + /** + * Converts greyscale image to a mesh + * \param elem_type defines if elements of the new mesh should be triangles or quads (or hexes for 3D) + * \param intensity_type defines how image intensities are interpreted + */ + static MeshLib::Mesh* convertRasterToMesh(GeoLib::Raster const& raster, + MeshElemType elem_type, + UseIntensityAs intensity_type); + /** * Converts greyscale image to a mesh * \param elem_type defines if elements of the new mesh should be triangles or quads (or hexes for 3D) diff --git a/Tests/MeshLib/TestMeshValidation.cpp b/Tests/MeshLib/TestMeshValidation.cpp index f89308a62049f85bd6da86d7f996eae5e2b08678..10ae99e4f94701faa843e7b862ca56810002df71 100644 --- a/Tests/MeshLib/TestMeshValidation.cpp +++ b/Tests/MeshLib/TestMeshValidation.cpp @@ -19,7 +19,8 @@ #include "MeshLib/Node.h" #include "MeshLib/MeshQuality/MeshValidation.h" #include "MeshLib/MeshGenerators/MeshGenerator.h" -#include "MeshLib/MeshGenerators/ConvertRasterToMesh.h" +#include "MeshLib/MeshGenerators/VtkMeshConverter.h" +//#include "MeshLib/MeshGenerators/ConvertRasterToMesh.h" #include "MeshLib/MeshEditing/DuplicateMeshComponents.h" #include "MeshLib/Elements/Element.h" @@ -47,8 +48,9 @@ TEST(MeshValidation, DetectHolesTri) GeoLib::RasterHeader const header = {4,3,MathLib::Point3d(std::array<double,3>{{0,0,0}}),1,-9999}; GeoLib::Raster const raster(header ,pix.begin(), pix.end()); - MeshLib::ConvertRasterToMesh conv(raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION); - auto mesh = std::unique_ptr<MeshLib::Mesh>{conv.execute()}; + //MeshLib::ConvertRasterToMesh conv(raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION); + //auto mesh = std::unique_ptr<MeshLib::Mesh>{conv.execute()}; + auto mesh = std::unique_ptr<MeshLib::Mesh>{MeshLib::VtkMeshConverter::convertRasterToMesh(raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION)}; ASSERT_EQ(0, MeshLib::MeshValidation::detectHoles(*mesh)); detectHoles(*mesh, {12}, 1); diff --git a/Tests/MeshLib/TestNodeSearch.cpp b/Tests/MeshLib/TestNodeSearch.cpp index 0510c90483136b17ca1c3dabc7b83853c5d5f4ef..a59e35f491dcfe6baaee918dc8bb0307da3a1887 100644 --- a/Tests/MeshLib/TestNodeSearch.cpp +++ b/Tests/MeshLib/TestNodeSearch.cpp @@ -14,7 +14,8 @@ #include "MeshLib/Node.h" #include "MeshLib/Elements/Element.h" #include "MeshLib/MeshSearch/NodeSearch.h" -#include "MeshLib/MeshGenerators/ConvertRasterToMesh.h" +#include "MeshLib/MeshGenerators/VtkMeshConverter.h" +//#include "MeshLib/MeshGenerators/ConvertRasterToMesh.h" #include "MeshLib/MeshEditing/DuplicateMeshComponents.h" #include "GeoLib/Raster.h" @@ -25,8 +26,9 @@ TEST(NodeSearch, UnusedNodes) GeoLib::RasterHeader const header = {4,3,MathLib::Point3d(std::array<double,3>{{0,0,0}}),1,-9999}; GeoLib::Raster const raster(header,pix.begin(), pix.end()); - MeshLib::ConvertRasterToMesh conv(raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION); - MeshLib::Mesh* mesh = conv.execute(); + MeshLib::Mesh* mesh (MeshLib::VtkMeshConverter::convertRasterToMesh(raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION)); + //MeshLib::ConvertRasterToMesh conv(raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION); + //MeshLib::Mesh* mesh = conv.execute(); MeshLib::NodeSearch ns(*mesh); ns.searchUnused(); std::vector<std::size_t> u_nodes = ns.getSearchedNodeIDs();