Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

VTK surface format changes #2594

Merged
merged 5 commits into from
Mar 8, 2023
Merged

VTK surface format changes #2594

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
e4b3530
VTK surface format changes
Lestropie Mar 3, 2023
11d3bf2
VTK load: Change terminal message RE endianness
Lestropie Mar 5, 2023
2fa95de
VTK: Change binary output to single-precision
Lestropie Mar 6, 2023
95696c5
meshconvert: Updated tests for endianness support
Lestropie Mar 6, 2023
c49eb0d
Fixes to pass CI tests for VTK changes (#2594)
Lestropie Mar 6, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
166 changes: 104 additions & 62 deletions src/surface/mesh.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -92,6 +92,20 @@ namespace MR



namespace {
template<typename T>
void load_vtk_points_binary (std::ifstream& in, const size_t num_vertices, vector<Eigen::Matrix<T, 3, 1>>& out)
{
out.reserve (num_vertices);
Eigen::Matrix<T, 3, 1> v;
for (size_t i = 0; i != num_vertices; ++i) {
in.read (reinterpret_cast<char*>(v.data()), 3 * sizeof (T));
out.push_back (v);
}
}
}



void Mesh::load_vtk (const std::string& path)
{
Expand Down Expand Up @@ -136,6 +150,13 @@ namespace MR
in.seekg (offset);
}

// Won't know endianness of file when the vertex positions are read,
// only when the polygon information is encountered;
bool change_endianness = false;

// If both float and big-endian, need to store natively as floats and swap endianness later
vector<Eigen::Matrix<float, 3, 1>> vertices_float;

// From here, don't necessarily know which parts of the data will come first
while (!in.eof()) {

Expand Down Expand Up @@ -166,25 +187,22 @@ namespace MR
throw Exception ("Error in reading binary .vtk file: Unsupported datatype (\"" + line + "\"");

vertices.reserve (num_vertices);
Lestropie marked this conversation as resolved.
Show resolved Hide resolved
for (int i = 0; i != num_vertices; ++i) {
if (!is_double)
vertices_float.reserve (num_vertices);

Vertex v;
if (is_ascii) {
Vertex v;
Eigen::Matrix<float, 3, 1> v_float;

if (is_ascii) {
for (int i = 0; i != num_vertices; ++i) {
std::getline (in, line);
sscanf (line.c_str(), "%lf %lf %lf", &v[0], &v[1], &v[2]);
} else {
if (is_double) {
double data[3];
in.read (reinterpret_cast<char*>(&data[0]), 3 * sizeof (double));
v = { data[0], data[1], data[2] };
} else {
float data[3];
in.read (reinterpret_cast<char*>(&data[0]), 3 * sizeof (float));
v = { data[0], data[1], data[2] };
}
vertices.push_back (v);
}
vertices.push_back (v);

} else if (is_double) {
load_vtk_points_binary<double> (in, num_vertices, vertices);
} else {
load_vtk_points_binary<float> (in, num_vertices, vertices_float);
}

} else if (line.substr (0, 8) == "POLYGONS") {
Expand All @@ -206,8 +224,16 @@ namespace MR
in.read (reinterpret_cast<char*>(&vertex_count), sizeof (int));
}

if (!is_ascii) {
if (change_endianness) {
vertex_count = ByteOrder::swap (vertex_count);
} else if (vertex_count != 3 && vertex_count != 4) {
vertex_count = ByteOrder::swap (vertex_count);
change_endianness = true;
}
}
if (vertex_count != 3 && vertex_count != 4)
throw Exception ("Could not parse file \"" + path + "\"; only suppport 3- and 4-vertex polygons");
throw Exception ("Could not parse file \"" + path + "\": only support 3- and 4-vertex polygons");

vector<unsigned int> t (vertex_count, 0);

Expand Down Expand Up @@ -237,17 +263,50 @@ namespace MR
}
}

// TODO If reading a binary file, may want to test endianness of data
// There's no explicit flag for this, but just calculating the standard
// deviations of all vertex positions may be adequate
// (likely to be huge if the endianness is wrong)
// Alternatively, just test the polygon indices: if there's at least one that exceeds the
// number of vertices, it may be saved in big-endian format, so try flipping everything
// Actually, should pop up at the first polygon read: number of points in polygon won't be 3 or 4
#if MRTRIX_IS_BIG_ENDIAN
if (change_endianness) {
WARN("File \"" + path + "\" is little-endian, so is not format-compliant (may have been generated using an older MRtrix3 version); "
"imported contents will be converted to system big-endian");
} else {
INFO("File \"" + path + "\" is big-endian; no format conversion required as executing on big-endian system");
}
#else
if (change_endianness) {
INFO("Converting imported contents of file \"" + path + "\" to native little-endian");
} else {
WARN("File \"" + path + "\" already in native little-endian format, so no endianness conversion required; "
"but file is therefore not format-compliant (may have been generated using an older MRtrix3 version)");
}
#endif

if (change_endianness) {
for (auto& v : vertices) {
for (size_t i = 0; i != 3; ++i)
v[i] = ByteOrder::swap (v[i]);
}
for (auto& v : vertices_float) {
for (size_t i = 0; i != 3; ++i)
v[i] = ByteOrder::swap (v[i]);
}
for (auto& t : triangles) {
for (size_t i = 0; i != 3; ++i)
t[i] = ByteOrder::swap (t[i]);
}
for (auto& q : quads) {
for (size_t i = 0; i != 4; ++i)
q[i] = ByteOrder::swap (q[i]);
}
}

if (vertices_float.size()) {
assert (!vertices.size());
for (const auto& v : vertices_float)
vertices.emplace_back (Vertex (v.cast<double>()));
}

try {
verify_data();
} catch(Exception& e) {
} catch (Exception& e) {
throw Exception (e, "Error verifying surface data from VTK file \"" + path + "\"");
}
}
Expand Down Expand Up @@ -634,66 +693,49 @@ namespace MR
ProgressBar progress ("writing mesh to file", vertices.size() + triangles.size() + quads.size());
if (binary) {

// FIXME Binary VTK output _still_ not working (crashes ParaView)
// Can however export as binary then -reconvert to ASCII and al is well...?
// Changed to big-endian output, doesn't seem to have fixed...

out.close();
out.open (path, std::ios_base::out | std::ios_base::app | std::ios_base::binary);
const bool is_double = (sizeof(default_type) == 8);
const std::string str_datatype = is_double ? "double" : "float";
const std::string points_header ("POINTS " + str(vertices.size()) + " " + str_datatype + "\n");
const std::string points_header ("POINTS " + str(vertices.size()) + " float\n");
out.write (points_header.c_str(), points_header.size());
for (VertexList::const_iterator i = vertices.begin(); i != vertices.end(); ++i) {
//float temp[3];
//for (size_t id = 0; id != 3; ++id)
// MR::putBE ((*i)[id], &temp[id]);
if (is_double) {
const double temp[3] { double((*i)[0]), double((*i)[1]), double((*i)[2]) };
out.write (reinterpret_cast<const char*>(temp), 3 * sizeof(double));
} else {
const float temp[3] { float((*i)[0]), float((*i)[1]), float((*i)[2]) };
out.write (reinterpret_cast<const char*>(temp), 3 * sizeof(float));
}
std::array<float, 3> temp_vertex;
for (const auto& v : vertices) {
temp_vertex = { ByteOrder::BE (float(v[0])), ByteOrder::BE (float(v[1])), ByteOrder::BE (float(v[2])) };
out.write (reinterpret_cast<const char*>(&temp_vertex), 3 * sizeof(float));
++progress;
}
const std::string polygons_header ("POLYGONS " + str(triangles.size() + quads.size()) + " " + str(4*triangles.size() + 5*quads.size()) + "\n");
out.write (polygons_header.c_str(), polygons_header.size());
const uint32_t num_points_triangle = 3;
for (TriangleList::const_iterator i = triangles.begin(); i != triangles.end(); ++i) {
const uint32_t num_points_triangle = ByteOrder::BE (uint32_t(3));
std::array<uint32_t, 3> temp_triangle;
for (const auto& t : triangles) {
out.write (reinterpret_cast<const char*>(&num_points_triangle), sizeof(uint32_t));
//uint32_t temp[3];
//for (size_t id = 0; id != 3; ++id)
// MR::putBE ((*i)[id], &temp[id]);
const uint32_t temp[3] { (*i)[0], (*i)[1], (*i)[2] };
out.write (reinterpret_cast<const char*>(temp), 3 * sizeof(uint32_t));
temp_triangle = { ByteOrder::BE (t[0]), ByteOrder::BE (t[1]), ByteOrder::BE (t[2]) };
out.write (reinterpret_cast<const char*>(&temp_triangle), 3 * sizeof(uint32_t));
++progress;
}
const uint32_t num_points_quad = 4;
for (QuadList::const_iterator i = quads.begin(); i != quads.end(); ++i) {
const uint32_t num_points_quad = ByteOrder::BE (uint32_t(4));
std::array<uint32_t, 4> temp_quad;
for (const auto& q : quads) {
out.write (reinterpret_cast<const char*>(&num_points_quad), sizeof(uint32_t));
//uint32_t temp[4];
//for (size_t id = 0; id != 4; ++id)
// MR::putBE ((*i)[id], &temp[id]);
const uint32_t temp[4] { (*i)[0], (*i)[1], (*i)[2], (*i)[3] };
out.write (reinterpret_cast<const char*>(temp), 4 * sizeof(uint32_t));
temp_quad = { ByteOrder::BE (q[0]), ByteOrder::BE (q[1]), ByteOrder::BE (q[2]), ByteOrder::BE (q[3]) };
out.write (reinterpret_cast<const char*>(&temp_quad), 4 * sizeof(uint32_t));
++progress;
}

} else {

out << "POINTS " << str(vertices.size()) << " float\n";
for (VertexList::const_iterator i = vertices.begin(); i != vertices.end(); ++i) {
out << str((*i)[0]) << " " << str((*i)[1]) << " " << str((*i)[2]) << "\n";
for (const auto& v : vertices) {
out << str<float>(v[0]) << " " << str<float>(v[1]) << " " << str<float>(v[2]) << "\n";
++progress;
}
out << "POLYGONS " + str(triangles.size() + quads.size()) + " " + str(4*triangles.size() + 5*quads.size()) + "\n";
for (TriangleList::const_iterator i = triangles.begin(); i != triangles.end(); ++i) {
out << "3 " << str((*i)[0]) << " " << str((*i)[1]) << " " << str((*i)[2]) << "\n";
for (const auto& t : triangles) {
out << "3 " << str(t[0]) << " " << str(t[1]) << " " << str(t[2]) << "\n";
++progress;
}
for (QuadList::const_iterator i = quads.begin(); i != quads.end(); ++i) {
out << "4 " << str((*i)[0]) << " " << str((*i)[1]) << " " << str((*i)[2]) << " " << str((*i)[3]) << "\n";
for (const auto& q : quads) {
out << "4 " << str(q[0]) << " " << str(q[1]) << " " << str(q[2]) << " " << str(q[3]) << "\n";
++progress;
}

Expand Down
2 changes: 1 addition & 1 deletion testing/binaries/data
Open in desktop
Submodule data updated 3 files
+0 −0 meshconvert/in_ascii.vtk
+ − meshconvert/in_be.vtk
+ − meshconvert/in_le.vtk
2 changes: 1 addition & 1 deletion testing/binaries/tests/mesh2voxel
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1 +1 @@
mesh2voxel meshconvert/in.vtk meshconvert/image.mif.gz - | testing_diff_image - mesh2voxel/out.mif.gz -abs 1.5e-3
mesh2voxel meshconvert/in_ascii.vtk meshconvert/image.mif.gz - | testing_diff_image - mesh2voxel/out.mif.gz -abs 1.5e-3
22 changes: 12 additions & 10 deletions testing/binaries/tests/meshconvert
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,10 +1,12 @@
meshconvert meshconvert/in.vtk tmp.vtk -force && testing_diff_mesh tmp.vtk meshconvert/in.vtk 0.001
meshconvert meshconvert/in.vtk tmp.vtk -binary -force && testing_diff_mesh tmp.vtk meshconvert/in.vtk 0.001
meshconvert meshconvert/in.vtk tmp.obj -force && testing_diff_mesh tmp.obj meshconvert/in.vtk 0.001
meshconvert meshconvert/in.vtk tmp.obj -binary -force && testing_diff_mesh tmp.obj meshconvert/in.vtk 0.001
meshconvert meshconvert/in.vtk tmp.stl -force && testing_diff_mesh tmp.stl meshconvert/in.vtk 0.001
meshconvert meshconvert/in.vtk tmp.stl -binary -force && testing_diff_mesh tmp.stl meshconvert/in.vtk 0.001
meshconvert meshconvert/in.vtk tmp.vtk -transform real2first meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/first.vtk 0.001
meshconvert meshconvert/first.vtk tmp.vtk -transform first2real meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/in.vtk 0.001
meshconvert meshconvert/in.vtk tmp.vtk -transform real2voxel meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/voxel.vtk 0.001
meshconvert meshconvert/voxel.vtk tmp.vtk -transform voxel2real meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/in.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.vtk -force && testing_diff_mesh tmp.vtk meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.vtk -binary -force && testing_diff_mesh tmp.vtk meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.obj -force && testing_diff_mesh tmp.obj meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.obj -binary -force && testing_diff_mesh tmp.obj meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.stl -force && testing_diff_mesh tmp.stl meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.stl -binary -force && testing_diff_mesh tmp.stl meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.vtk -transform real2first meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/first.vtk 0.001
meshconvert meshconvert/first.vtk tmp.vtk -transform first2real meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_ascii.vtk tmp.vtk -transform real2voxel meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/voxel.vtk 0.001
meshconvert meshconvert/voxel.vtk tmp.vtk -transform voxel2real meshconvert/image.mif.gz -force && testing_diff_mesh tmp.vtk meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_le.vtk tmp.vtk -force && testing_diff_mesh tmp.vtk meshconvert/in_ascii.vtk 0.001
meshconvert meshconvert/in_be.vtk tmp.vtk -force && testing_diff_mesh tmp.vtk meshconvert/in_ascii.vtk 0.001