Skip to content

Commit fbf70b4

Browse files
ktboltktbolt
authored
Read mesh domains from vtk file 363 (SimVascular#513)
* Add copying vtk cell data into an int vector. * Clean up data names, make VtkData interface consistent.
1 parent 8df8135 commit fbf70b4

6 files changed

Lines changed: 170 additions & 22 deletions

File tree

‎Code/Source/solver/VtkData.cpp‎

Lines changed: 74 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -629,6 +629,30 @@ Array<int> VtkVtpData::get_connectivity() const
629629
return conn;
630630
}
631631

632+
/// @brief Copy an array of cell data from an polydata mesh into the given Vector.
633+
//
634+
void VtkVtpData::copy_cell_data(const std::string& data_name, Vector<int>& mesh_data)
635+
{
636+
auto vtk_data = vtkIntArray::SafeDownCast(impl->vtk_polydata->GetCellData()->GetArray(data_name.c_str()));
637+
if (vtk_data == nullptr) {
638+
return;
639+
}
640+
641+
int num_data = vtk_data->GetNumberOfTuples();
642+
if (num_data == 0) {
643+
return;
644+
}
645+
646+
int num_comp = vtk_data->GetNumberOfComponents();
647+
648+
for (int i = 0; i < num_data; i++) {
649+
auto tuple = vtk_data->GetTuple(i);
650+
for (int j = 0; j < num_comp; j++) {
651+
mesh_data(i) = tuple[j];
652+
}
653+
}
654+
}
655+
632656
/// @brief Copy an array of point data from an polydata mesh into the given Array.
633657
//
634658
void VtkVtpData::copy_point_data(const std::string& data_name, Array<double>& mesh_data)
@@ -774,6 +798,19 @@ Array<double> VtkVtpData::get_points() const
774798
return points_array;
775799
}
776800

801+
bool VtkVtpData::has_cell_data(const std::string& data_name)
802+
{
803+
int num_arrays = impl->vtk_polydata->GetPointData()->GetNumberOfArrays();
804+
805+
for (int i = 0; i < num_arrays; i++) {
806+
if (!strcmp(impl->vtk_polydata->GetCellData()->GetArrayName(i), data_name.c_str())) {
807+
return true;
808+
}
809+
}
810+
811+
return false;
812+
}
813+
777814
bool VtkVtpData::has_point_data(const std::string& data_name)
778815
{
779816
int num_arrays = impl->vtk_polydata->GetPointData()->GetNumberOfArrays();
@@ -911,6 +948,30 @@ std::vector<std::string> VtkVtuData::get_point_data_names()
911948
return data_names;
912949
}
913950

951+
/// @brief Copy an array of cell data from an unstructured mesh into the given Vector.
952+
//
953+
void VtkVtuData::copy_cell_data(const std::string& data_name, Vector<int>& mesh_data)
954+
{
955+
auto vtk_data = vtkIntArray::SafeDownCast(impl->vtk_ugrid->GetCellData()->GetArray(data_name.c_str()));
956+
if (vtk_data == nullptr) {
957+
return;
958+
}
959+
960+
int num_data = vtk_data->GetNumberOfTuples();
961+
if (num_data == 0) {
962+
return;
963+
}
964+
965+
int num_comp = vtk_data->GetNumberOfComponents();
966+
967+
for (int i = 0; i < num_data; i++) {
968+
auto tuple = vtk_data->GetTuple(i);
969+
for (int j = 0; j < num_comp; j++) {
970+
mesh_data(i) = tuple[j];
971+
}
972+
}
973+
}
974+
914975
/// @brief Copy an array of point data from an unstructured grid into the given Array.
915976
//
916977
void VtkVtuData::copy_point_data(const std::string& data_name, Array<double>& mesh_data)
@@ -996,6 +1057,19 @@ void VtkVtuData::copy_points(Array<double>& points)
9961057
return;
9971058
}
9981059

1060+
bool VtkVtuData::has_cell_data(const std::string& data_name)
1061+
{
1062+
int num_arrays = impl->vtk_ugrid->GetPointData()->GetNumberOfArrays();
1063+
1064+
for (int i = 0; i < num_arrays; i++) {
1065+
if (!strcmp(impl->vtk_ugrid->GetCellData()->GetArrayName(i), data_name.c_str())) {
1066+
return true;
1067+
}
1068+
}
1069+
1070+
return false;
1071+
}
1072+
9991073
bool VtkVtuData::has_point_data(const std::string& data_name)
10001074
{
10011075
int num_arrays = impl->vtk_ugrid->GetPointData()->GetNumberOfArrays();

‎Code/Source/solver/VtkData.h‎

Lines changed: 37 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -32,11 +32,20 @@ class VtkData {
3232
virtual void set_points(const Array<double>& points) = 0;
3333
virtual void set_connectivity(const int nsd, const Array<int>& conn, const int pid = 0) = 0;
3434

35+
virtual bool has_cell_data(const std::string& data_name) = 0;
3536
virtual bool has_point_data(const std::string& data_name) = 0;
3637

3738
virtual void copy_points(Array<double>& points) = 0;
39+
3840
virtual void copy_point_data(const std::string& data_name, Array<double>& mesh_data) = 0;
3941
virtual void copy_point_data(const std::string& data_name, Vector<double>& mesh_data) = 0;
42+
virtual void copy_point_data(const std::string& data_name, Vector<int>& mesh_data) = 0;
43+
44+
virtual void copy_cell_data(const std::string& data_name, Vector<int>& mesh_data) = 0;
45+
46+
virtual Array<double> get_point_data(const std::string& data_name) = 0;
47+
virtual std::vector<std::string> get_point_data_names() = 0;
48+
4049
virtual void write() = 0;
4150

4251
static VtkData* create_reader(const std::string& file_name);
@@ -64,13 +73,20 @@ class VtkVtpData : public VtkData {
6473
virtual int num_points() const override;
6574
virtual void read_file(const std::string& file_name) override;
6675

67-
void copy_points(Array<double>& points) override;
68-
void copy_point_data(const std::string& data_name, Array<double>& mesh_data) override;
69-
void copy_point_data(const std::string& data_name, Vector<double>& mesh_data) override;
70-
void copy_point_data(const std::string& data_name, Vector<int>& mesh_data);
71-
Array<double> get_point_data(const std::string& data_name);
72-
std::vector<std::string> get_point_data_names();
73-
bool has_point_data(const std::string& data_name) override;
76+
virtual void copy_points(Array<double>& points) override;
77+
78+
virtual void copy_point_data(const std::string& data_name, Array<double>& mesh_data) override;
79+
virtual void copy_point_data(const std::string& data_name, Vector<double>& mesh_data) override;
80+
virtual void copy_point_data(const std::string& data_name, Vector<int>& mesh_data) override;
81+
82+
virtual void copy_cell_data(const std::string& data_name, Vector<int>& mesh_data) override;
83+
84+
virtual Array<double> get_point_data(const std::string& data_name) override;
85+
virtual std::vector<std::string> get_point_data_names() override;
86+
87+
virtual bool has_cell_data(const std::string& data_name) override;
88+
virtual bool has_point_data(const std::string& data_name) override;
89+
7490
virtual void set_connectivity(const int nsd, const Array<int>& conn, const int pid = 0) override;
7591

7692
virtual void set_element_data(const std::string& data_name, const Array<double>& data) override;
@@ -101,15 +117,22 @@ class VtkVtuData : public VtkData {
101117
virtual int num_points() const override;
102118
virtual void read_file(const std::string& file_name) override;
103119

104-
void copy_points(Array<double>& points) override;
105-
void copy_point_data(const std::string& data_name, Array<double>& mesh_data) override;
106-
void copy_point_data(const std::string& data_name, Vector<double>& mesh_data) override;
107-
void copy_point_data(const std::string& data_name, Vector<int>& mesh_data);
120+
virtual void copy_points(Array<double>& points) override;
121+
122+
virtual void copy_point_data(const std::string& data_name, Array<double>& mesh_data) override;
123+
virtual void copy_point_data(const std::string& data_name, Vector<double>& mesh_data) override;
124+
virtual void copy_point_data(const std::string& data_name, Vector<int>& mesh_data) override;
125+
126+
virtual void copy_cell_data(const std::string& data_name, Vector<int>& mesh_data) override;
127+
128+
virtual Array<double> get_point_data(const std::string& data_name) override;
129+
virtual std::vector<std::string> get_point_data_names() override;
108130

109-
Array<double> get_point_data(const std::string& data_name);
110-
std::vector<std::string> get_point_data_names();
111131
virtual Array<double> get_points() const override;
112-
bool has_point_data(const std::string& data_name) override;
132+
133+
virtual bool has_cell_data(const std::string& data_name) override;
134+
virtual bool has_point_data(const std::string& data_name) override;
135+
113136
virtual void set_connectivity(const int nsd, const Array<int>& conn, const int pid = 0) override;
114137

115138
virtual void set_element_data(const std::string& data_name, const Array<double>& data) override;

‎Code/Source/solver/read_msh.cpp‎

Lines changed: 23 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -2086,16 +2086,32 @@ void set_dmn_id_ff(Simulation* simulation, mshType& lM, const std::string& file_
20862086
}
20872087
}
20882088

2089-
/// @brief Read mesh domains from a vtu/vtp file.
2089+
/// @brief Read mesh domains from a vtk vtu/vtp file.
20902090
///
2091-
/// \todo [NOTE] Not implemented.
2092-
//
2093-
void set_dmn_id_vtk(Simulation* simulation, mshType& mesh, const std::string& file_name, const std::string& kwrd)
2091+
void set_dmn_id_vtk(Simulation* simulation, mshType& lM, const std::string& file_name, const std::string& data_name)
20942092
{
2095-
int btSiz = std::numeric_limits<int>::digits;
2093+
#define n_debug_set_dmn_id_vtk
2094+
#ifdef debug_set_dmn_id_vtk
2095+
DebugMsg dmsg(__func__, simulation->com_mod.cm.idcm());
2096+
dmsg.banner();
2097+
dmsg << "file_name: " << file_name;
2098+
dmsg << "lM.gnEl: " << lM.gnEl;
2099+
#endif
2100+
2101+
if (lM.eId.size() == 0) {
2102+
lM.eId.resize(lM.gnEl);
2103+
}
2104+
2105+
Vector<int> element_data(lM.gnEl);
2106+
vtk_xml::read_element_data(lM, file_name, data_name, element_data);
2107+
2108+
for (int a = 0; a < lM.gnEl; a++) {
2109+
lM.eId(a) = lM.eId(a) | (1UL << element_data(a));
2110+
}
2111+
20962112
}
20972113

2098-
/// @brief This routines associates two faces with each other and sets gN.
2114+
/// @brief Associate two faces with each other and set gN.
20992115
///
21002116
/// Data set
21012117
/// \code {.cpp}
@@ -2232,4 +2248,4 @@ void set_projector(Simulation* simulation, utils::stackType& avNds)
22322248
}
22332249
}
22342250

2235-
};
2251+
};

‎Code/Source/solver/read_msh.h‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -52,7 +52,7 @@ namespace read_msh_ns {
5252
void read_msh(Simulation* simulation);
5353

5454
void set_dmn_id_ff(Simulation* simulation, mshType& mesh, const std::string& file_name);
55-
void set_dmn_id_vtk(Simulation* simulation, mshType& mesh, const std::string& file_name, const std::string& kwrd);
55+
void set_dmn_id_vtk(Simulation* simulation, mshType& lM, const std::string& file_name, const std::string& kwrd);
5656
void set_projector(Simulation* simulation, utils::stackType& avNds);
5757
void set_ris_projector(Simulation* simulation);
5858
void set_uris_meshes(Simulation* simulation);

‎Code/Source/solver/vtk_xml.cpp‎

Lines changed: 33 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -454,6 +454,39 @@ void read_vtp(const std::string& file_name, faceType& face)
454454
}
455455
}
456456

457+
/// @brief Read in element data from an integer vector array data
458+
/// named `data_name` stored in VTK VTU file into 'element_data'.
459+
///
460+
void read_element_data(const mshType& mesh, const std::string& file_name, const std::string& data_name, Vector<int>& element_data)
461+
{
462+
if (FILE *file = fopen(file_name.c_str(), "r")) {
463+
fclose(file);
464+
} else {
465+
throw std::runtime_error("The VTK element data file '" + file_name + "' can't be read.");
466+
}
467+
468+
auto vtk_data = VtkData::create_reader(file_name);
469+
if (vtk_data == nullptr) {
470+
throw std::runtime_error("Failed to create VTK reader for file: " + file_name);
471+
}
472+
473+
int num_elems = vtk_data->num_elems();
474+
int num_points = vtk_data->num_points();
475+
476+
if (num_elems != mesh.gnEl) {
477+
throw std::runtime_error("The number of elements (" + std::to_string(num_elems) +
478+
") in the file '" + file_name + "' is not equal to the number of elements ("
479+
+ std::to_string(mesh.gnEl) + ") for the mesh named '" + mesh.name + "'.");
480+
}
481+
482+
if (!vtk_data->has_cell_data(data_name)) {
483+
throw std::runtime_error("No CellData DataArray named '" + data_name +
484+
"' found in the VTK file '" + file_name + "' for the mesh named '" + mesh.name + "'.");
485+
}
486+
487+
vtk_data->copy_cell_data(data_name, element_data);
488+
}
489+
457490
//----------------
458491
// read_vtp_pdata
459492
//----------------

‎Code/Source/solver/vtk_xml.h‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -22,6 +22,8 @@ void read_vtu(const std::string& file_name, mshType& mesh);
2222

2323
void read_precomputed_solution_vtu(const std::string& file_name, const std::string& field_name, mshType& mesh);
2424

25+
void read_element_data(const mshType& mesh, const std::string& fName, const std::string& kwrd, Vector<int>& tmpR);
26+
2527
void read_vtu_pdata(const std::string& fName, const std::string& kwrd, const int nsd, const int m, const int idx, mshType& mesh);
2628

2729
void read_vtus(Simulation* simulation, Array<double>& lA, Array<double>& lY, Array<double>& lD, const std::string& fName);

0 commit comments

Comments
 (0)