OGS
MeshToolsLib::MeshGenerator::VoxelGridFromMesh Namespace Reference

Functions

std::array< std::size_t, 3 > getNumberOfVoxelPerDimension (std::array< double, 3 > const &ranges, std::array< double, 3 > const &cellsize)
std::vector< int > assignCellIds (vtkSmartPointer< vtkUnstructuredGrid > const &mesh, MathLib::Point3d const &min, std::array< std::size_t, 3 > const &dims, std::array< double, 3 > const &cellsize)
bool removeUnusedGridCells (vtkSmartPointer< vtkUnstructuredGrid > const &mesh, std::unique_ptr< MeshLib::Mesh > &grid)
template<typename T, typename VTK_TYPE>
void mapArray (MeshLib::Mesh &grid, VTK_TYPE vtk_arr, std::string const &name)
template<typename T>
bool checkDyncast (MeshLib::Mesh &mesh, vtkSmartPointer< vtkCellData > const cell_data, char const *const name)
void mapMeshArraysOntoGrid (vtkSmartPointer< vtkUnstructuredGrid > const &mesh, std::unique_ptr< MeshLib::Mesh > const &grid)

Variables

static std::string const cell_id_name = "CellIds"

Function Documentation

◆ assignCellIds()

std::vector< int > MeshToolsLib::MeshGenerator::VoxelGridFromMesh::assignCellIds ( vtkSmartPointer< vtkUnstructuredGrid > const & mesh,
MathLib::Point3d const & min,
std::array< std::size_t, 3 > const & dims,
std::array< double, 3 > const & cellsize )

Definition at line 56 of file VoxelGridFromMesh.cpp.

60{
61 vtkSmartPointer<vtkCellLocator> locator =
62 vtkSmartPointer<vtkCellLocator>::New();
63 locator->SetDataSet(mesh);
64 locator->Update();
65
66 std::vector<int> cell_ids;
67 cell_ids.reserve(dims[0] * dims[1] * dims[2]);
68 std::array<double, 3> const grid_max = {min[0] + dims[0] * cellsize[0],
69 min[1] + dims[1] * cellsize[1],
70 min[2] + dims[2] * cellsize[2]};
71
72 double const start[3] = {min[0] + cellsize[0] / 2.,
73 min[1] + cellsize[1] / 2.,
74 min[2] + cellsize[2] / 2.};
75 double pnt[3];
76 for (pnt[2] = start[2]; pnt[2] < grid_max[2]; pnt[2] += cellsize[2])
77 {
78 for (pnt[1] = start[1]; pnt[1] < grid_max[1]; pnt[1] += cellsize[1])
79 {
80 for (pnt[0] = start[0]; pnt[0] < grid_max[0]; pnt[0] += cellsize[0])
81 {
82 cell_ids.push_back(static_cast<int>(locator->FindCell(pnt)));
83 }
84 }
85 }
86 return cell_ids;
87}

Referenced by Vtu2GridDialog::accept(), and main().

◆ checkDyncast()

template<typename T>
bool MeshToolsLib::MeshGenerator::VoxelGridFromMesh::checkDyncast ( MeshLib::Mesh & mesh,
vtkSmartPointer< vtkCellData > const cell_data,
char const *const name )

Definition at line 134 of file VoxelGridFromMesh.cpp.

137{
138 using DataArrayType = vtkAOSDataArrayTemplate<T>;
139 vtkSmartPointer<DataArrayType> const arr =
140 dynamic_cast<DataArrayType*>(cell_data->GetArray(name));
141 if (!arr)
142 {
143 return false;
144 }
146 return true;
147}
void mapArray(MeshLib::Mesh &grid, VTK_TYPE vtk_arr, std::string const &name)

References mapArray().

Referenced by mapMeshArraysOntoGrid().

◆ getNumberOfVoxelPerDimension()

std::array< std::size_t, 3 > MeshToolsLib::MeshGenerator::VoxelGridFromMesh::getNumberOfVoxelPerDimension ( std::array< double, 3 > const & ranges,
std::array< double, 3 > const & cellsize )

getNumberOfVoxelPerDimension is used to calculate how many voxel fit into a bounding box. For this calculation the difference of min and max point of the bounding box is divided by the cell size, for every dimension. The calculation is restricted to work only with positive values for the cell size. If the difference between min and max is zero, we assign one voxel for the respective dimension.

Definition at line 31 of file VoxelGridFromMesh.cpp.

33{
34 if (cellsize[0] <= 0 || cellsize[1] <= 0 || cellsize[2] <= 0)
35 {
36 OGS_FATAL("A cellsize ({},{},{}) is not allowed to be <= 0",
37 cellsize[0], cellsize[1], cellsize[2]);
38 }
39 std::array<double, 3> numberOfVoxel = {ranges[0] / cellsize[0],
40 ranges[1] / cellsize[1],
41 ranges[2] / cellsize[2]};
42
43 if (ranges[0] < 0 || ranges[1] < 0 || ranges[2] < 0)
44 {
46 "The difference of max-min ({},{},{}) is not allowed to be < 0",
47 ranges[0], ranges[1], ranges[2]);
48 }
49 std::replace(numberOfVoxel.begin(), numberOfVoxel.end(), 0, 1);
50
51 return {static_cast<std::size_t>(std::lround(numberOfVoxel[0])),
52 static_cast<std::size_t>(std::lround(numberOfVoxel[1])),
53 static_cast<std::size_t>(std::lround(numberOfVoxel[2]))};
54}
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

Referenced by Vtu2GridDialog::accept(), and main().

◆ mapArray()

template<typename T, typename VTK_TYPE>
void MeshToolsLib::MeshGenerator::VoxelGridFromMesh::mapArray ( MeshLib::Mesh & grid,
VTK_TYPE vtk_arr,
std::string const & name )

Definition at line 112 of file VoxelGridFromMesh.cpp.

113{
114 auto const* cell_ids = grid.getProperties().getPropertyVector<int>(
116 if (cell_ids == nullptr)
117 {
118 // Error message
119 return;
120 }
121 std::size_t const n_elems = cell_ids->size();
122 auto& arr = *grid.getProperties().createNewPropertyVector<T>(
123 name, MeshLib::MeshItemType::Cell, n_elems,
124 vtk_arr->GetNumberOfComponents());
125 for (std::size_t j = 0; j < n_elems; ++j)
126 {
127 arr[j] = vtk_arr->GetValue((*cell_ids)[j]);
128 }
129}
Properties & getProperties()
Definition Mesh.h:136
PropertyVector< T > * createNewPropertyVector(std::string_view name, MeshItemType mesh_item_type, std::size_t n_components=1)
Definition Properties.h:17
PropertyVector< T > const * getPropertyVector(std::string_view name) const
constexpr std::size_t size() const

References MeshLib::Cell, cell_id_name, MeshLib::Properties::createNewPropertyVector(), MeshLib::Mesh::getProperties(), MeshLib::Properties::getPropertyVector(), and MeshLib::PropertyVector< PROP_VAL_TYPE >::size().

Referenced by checkDyncast().

◆ mapMeshArraysOntoGrid()

void MeshToolsLib::MeshGenerator::VoxelGridFromMesh::mapMeshArraysOntoGrid ( vtkSmartPointer< vtkUnstructuredGrid > const & mesh,
std::unique_ptr< MeshLib::Mesh > const & grid )

Definition at line 149 of file VoxelGridFromMesh.cpp.

151{
152 assert(mesh != nullptr);
153 assert(grid != nullptr);
154 vtkSmartPointer<vtkCellData> const cell_data = mesh->GetCellData();
155 for (int i = 0; i < cell_data->GetNumberOfArrays(); ++i)
156 {
157 auto const name = cell_data->GetArrayName(i);
158
159 if (!(checkDyncast<double>(*grid, cell_data, name) ||
160 checkDyncast<long>(*grid, cell_data, name) ||
161 checkDyncast<long long>(*grid, cell_data, name) ||
162 checkDyncast<int>(*grid, cell_data, name)))
163 {
164 WARN("Ignoring array '{:s}', array type {:s} not implemented...",
165 name,
166 cell_data->GetArray(name)->GetDataTypeAsString());
167 }
168 }
169}
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:42
bool checkDyncast(MeshLib::Mesh &mesh, vtkSmartPointer< vtkCellData > const cell_data, char const *const name)

References checkDyncast(), and WARN().

Referenced by Vtu2GridDialog::accept(), and main().

◆ removeUnusedGridCells()

bool MeshToolsLib::MeshGenerator::VoxelGridFromMesh::removeUnusedGridCells ( vtkSmartPointer< vtkUnstructuredGrid > const & mesh,
std::unique_ptr< MeshLib::Mesh > & grid )

Definition at line 89 of file VoxelGridFromMesh.cpp.

91{
92 MeshLib::ElementSearch search(*grid);
93 std::size_t const n_elems_marked = search.searchByPropertyValueRange<int>(
94 cell_id_name, 0, static_cast<int>(mesh->GetNumberOfCells()), true);
95
96 if (n_elems_marked == grid->getNumberOfElements())
97 {
98 ERR("No valid elements found. Aborting...");
99 return false;
100 }
101
102 if (n_elems_marked)
103 {
105 *grid, search.getSearchedElementIDs(), "trimmed_grid"));
106 }
107 return true;
108}
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:48
Element search class.
MeshLib::Mesh * removeElements(const MeshLib::Mesh &mesh, const std::vector< std::size_t > &removed_element_ids, const std::string &new_mesh_name)

References cell_id_name, ERR(), MeshLib::ElementSearch::getSearchedElementIDs(), MeshToolsLib::removeElements(), and MeshLib::ElementSearch::searchByPropertyValueRange().

Referenced by Vtu2GridDialog::accept(), and main().

Variable Documentation

◆ cell_id_name

std::string const MeshToolsLib::MeshGenerator::VoxelGridFromMesh::cell_id_name = "CellIds"
static

Definition at line 30 of file VoxelGridFromMesh.h.

Referenced by Vtu2GridDialog::accept(), main(), mapArray(), and removeUnusedGridCells().