11#include <vtkAbstractArray.h>
12#include <vtkCellData.h>
13#include <vtkCellLocator.h>
14#include <vtkDataArray.h>
15#include <vtkDoubleArray.h>
16#include <vtkIntArray.h>
17#include <vtkLongLongArray.h>
32 std::array<double, 3>
const& ranges, std::array<double, 3>
const& cellsize)
34 if (cellsize[0] <= 0 || cellsize[1] <= 0 || cellsize[2] <= 0)
36 OGS_FATAL(
"A cellsize ({},{},{}) is not allowed to be <= 0",
37 cellsize[0], cellsize[1], cellsize[2]);
39 std::array<double, 3> numberOfVoxel = {ranges[0] / cellsize[0],
40 ranges[1] / cellsize[1],
41 ranges[2] / cellsize[2]};
43 if (ranges[0] < 0 || ranges[1] < 0 || ranges[2] < 0)
46 "The difference of max-min ({},{},{}) is not allowed to be < 0",
47 ranges[0], ranges[1], ranges[2]);
49 std::replace(numberOfVoxel.begin(), numberOfVoxel.end(), 0, 1);
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]))};
56std::vector<int>
assignCellIds(vtkSmartPointer<vtkUnstructuredGrid>
const& mesh,
58 std::array<std::size_t, 3>
const& dims,
59 std::array<double, 3>
const& cellsize)
61 vtkSmartPointer<vtkCellLocator> locator =
62 vtkSmartPointer<vtkCellLocator>::New();
63 locator->SetDataSet(mesh);
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]};
72 double const start[3] = {min[0] + cellsize[0] / 2.,
73 min[1] + cellsize[1] / 2.,
74 min[2] + cellsize[2] / 2.};
76 for (pnt[2] = start[2]; pnt[2] < grid_max[2]; pnt[2] += cellsize[2])
78 for (pnt[1] = start[1]; pnt[1] < grid_max[1]; pnt[1] += cellsize[1])
80 for (pnt[0] = start[0]; pnt[0] < grid_max[0]; pnt[0] += cellsize[0])
82 cell_ids.push_back(
static_cast<int>(locator->FindCell(pnt)));
90 std::unique_ptr<MeshLib::Mesh>& grid)
94 cell_id_name, 0,
static_cast<int>(mesh->GetNumberOfCells()),
true);
96 if (n_elems_marked == grid->getNumberOfElements())
98 ERR(
"No valid elements found. Aborting...");
111template <
typename T,
typename VTK_TYPE>
116 if (cell_ids ==
nullptr)
121 std::size_t
const n_elems = cell_ids->size();
124 vtk_arr->GetNumberOfComponents());
125 for (std::size_t j = 0; j < n_elems; ++j)
127 arr[j] = vtk_arr->GetValue((*cell_ids)[j]);
135 vtkSmartPointer<vtkCellData>
const cell_data,
136 char const*
const name)
138 using DataArrayType = vtkAOSDataArrayTemplate<T>;
139 vtkSmartPointer<DataArrayType>
const arr =
140 dynamic_cast<DataArrayType*
>(cell_data->GetArray(name));
150 std::unique_ptr<MeshLib::Mesh>
const& grid)
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)
157 auto const name = cell_data->GetArrayName(i);
164 WARN(
"Ignoring array '{:s}', array type {:s} not implemented...",
166 cell_data->GetArray(name)->GetDataTypeAsString());
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition of the Point3d class.
const std::vector< std::size_t > & getSearchedElementIDs() const
return marked elements
std::size_t searchByPropertyValueRange(std::string const &property_name, PROPERTY_TYPE const min_property_value, PROPERTY_TYPE const max_property_value, bool outside_of)
Properties & getProperties()
PropertyVector< T > * createNewPropertyVector(std::string_view name, MeshItemType mesh_item_type, std::size_t n_components=1)
PropertyVector< T > const * getPropertyVector(std::string_view name) const