Loading [MathJax]/extensions/MathMenu.js
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 57 of file VoxelGridFromMesh.cpp.

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

◆ checkDyncast()

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

Definition at line 135 of file VoxelGridFromMesh.cpp.

138{
139 using DataArrayType = vtkAOSDataArrayTemplate<T>;
140 vtkSmartPointer<DataArrayType> const arr =
141 dynamic_cast<DataArrayType*>(cell_data->GetArray(name));
142 if (!arr)
143 {
144 return false;
145 }
146 mapArray<T, vtkSmartPointer<DataArrayType>>(mesh, arr, name);
147 return true;
148}

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 32 of file VoxelGridFromMesh.cpp.

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

References OGS_FATAL.

◆ 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 113 of file VoxelGridFromMesh.cpp.

114{
115 auto const* cell_ids = grid.getProperties().getPropertyVector<int>(
116 cell_id_name, MeshLib::MeshItemType::Cell, 1);
117 if (cell_ids == nullptr)
118 {
119 // Error message
120 return;
121 }
122 std::size_t const n_elems = cell_ids->size();
123 auto& arr = *grid.getProperties().createNewPropertyVector<T>(
124 name, MeshLib::MeshItemType::Cell, n_elems,
125 vtk_arr->GetNumberOfComponents());
126 for (std::size_t j = 0; j < n_elems; ++j)
127 {
128 arr[j] = vtk_arr->GetValue((*cell_ids)[j]);
129 }
130}
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

References MeshLib::Cell, cell_id_name, MeshLib::Properties::createNewPropertyVector(), MeshLib::Mesh::getProperties(), and MeshLib::Properties::getPropertyVector().

Referenced by checkDyncast().

◆ mapMeshArraysOntoGrid()

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

Definition at line 150 of file VoxelGridFromMesh.cpp.

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

References checkDyncast(), and WARN().

◆ removeUnusedGridCells()

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

Definition at line 90 of file VoxelGridFromMesh.cpp.

92{
93 MeshLib::ElementSearch search(*grid);
94 std::size_t const n_elems_marked = search.searchByPropertyValueRange<int>(
95 cell_id_name, 0, static_cast<int>(mesh->GetNumberOfCells()), true);
96
97 if (n_elems_marked == grid->getNumberOfElements())
98 {
99 ERR("No valid elements found. Aborting...");
100 return false;
101 }
102
103 if (n_elems_marked)
104 {
106 *grid, search.getSearchedElementIDs(), "trimmed_grid"));
107 }
108 return true;
109}
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:45
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().

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 mapArray(), and removeUnusedGridCells().