read(const std::string &in, Format format = Default);
/**
- * Read grid data from an vtk file. Numerical data is ignored.
+ * Read grid data from a unstructured vtk file. The vtk file may contain
+ * the following VTK cell types: VTK_HEXAHEDRON, VTK_QUAD, and VTK_LINE.
+ *
+ * Depending on the template dimension, only some of the above are accepted.
+ *
+ * In particular, in three dimensions, this function expects the file to
+ * contain
+ *
+ * - VTK_HEXAHEDRON cell types
+ * - VTK_QUAD cell types, to specify optional boundary or interior quad faces
+ * - VTK_LINE cell types, to specify optional boundary or interior edges
+ *
+ * In two dimensions:
+ *
+ * - VTK_QUAD cell types
+ * - VTK_LINE cell types, to specify optional boundary or interior edges
+ *
+ * In one dimension
+ *
+ * - VTK_LINE cell types
+ *
+ * The input file may specify boundary ids, material ids, and manifold ids
+ * using the CELL_DATA section of the
+ * [VTK file format](http://www.vtk.org/VTK/img/file-formats.pdf).
+ *
+ * This function interprets two types of CELL_DATA contained in the input
+ * file: `SCALARS MaterialID`, used to specify the material_id of the cells,
+ * or the boundary_id of the faces and edges, and `SCALARS ManifoldID`, that
+ * can be used to specify the manifold id of any Triangulation object (cell,
+ * face, or edge).
+ *
+ * The companion GridOut::write_vtk function can be used to write VTK files
+ * compatible with this method.
*
* @author Mayank Sabharwal, Andreas Putz, 2013
+ * @author Luca Heltai, 2018
*/
void
read_vtk(std::istream &in);
void
GridIn<dim, spacedim>::read_vtk(std::istream &in)
{
- Assert((dim == 2) || (dim == 3), ExcNotImplemented());
std::string line;
// verify that the first, third and fourth lines match
unsigned int n_vertices;
in >> n_vertices;
- in.ignore(256,
- '\n'); // ignoring the number beside the total no. of points.
+ in >> keyword; // float, double, int, char, etc.
for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
{
ExcMessage(
"While reading VTK file, failed to find POINTS section"));
-
- //////////////////ignoring space between points and cells
- /// sections////////////////////
- std::string checkline;
- int no;
- in.ignore(256,
- '\n'); // this move pointer to the next line ignoring unwanted no.
- no = in.tellg();
- getline(in, checkline);
- if (checkline.compare("") != 0)
- {
- in.seekg(no);
- }
-
in >> keyword;
- ///////////////////Processing the CELLS section that contains cells(cells) and
- /// bound_quads(subcelldata)///////////////////////
+ unsigned int n_geometric_objects = 0;
+ unsigned int n_ints;
if (keyword == "CELLS")
{
- unsigned int n_geometric_objects;
in >> n_geometric_objects;
- in.ignore(256, '\n');
+ in >> n_ints; // Ignore this, since we don't need it.
if (dim == 3)
{
else if (type == 4)
{
+ // we assume that the file contains first all cells,
+ // then all faces, and finally all lines
+ AssertThrow(subcelldata.boundary_lines.size() == 0,
+ ExcNotImplemented());
+
subcelldata.boundary_quads.emplace_back();
for (unsigned int j = 0; j < type;
subcelldata.boundary_quads.back().material_id = 0;
}
+ else if (type == 2)
+ {
+ subcelldata.boundary_lines.emplace_back();
+
+ for (unsigned int j = 0; j < type;
+ j++) // loop to feed the data to the boundary
+ in >> subcelldata.boundary_lines.back().vertices[j];
+
+ subcelldata.boundary_lines.back().material_id = 0;
+ }
else
AssertThrow(
AssertThrow(
false,
ExcMessage(
- "While reading VTK file, unknown file type encountered"));
+ "While reading VTK file, unknown cell type encountered"));
+ }
+ }
+ else if (dim == 1)
+ {
+ for (unsigned int count = 0; count < n_geometric_objects; count++)
+ {
+ unsigned int type;
+ in >> type;
+
+ if (type == 2)
+ {
+ cells.emplace_back();
+
+ for (unsigned int j = 0; j < type; j++) // loop to feed data
+ in >> cells.back().vertices[j];
+
+ cells.back().material_id = 0;
+ }
+ else
+ AssertThrow(
+ false,
+ ExcMessage(
+ "While reading VTK file, unknown cell type encountered"));
}
}
else
if (keyword ==
"CELL_TYPES") // Entering the cell_types section and ignoring data.
{
- in.ignore(256, '\n');
-
- while (!in.fail() && !in.eof())
- {
- in >> keyword;
- if (keyword != "12" && keyword != "9")
- {
- break;
- }
- }
+ in >> n_ints;
+ int tmp_int;
+ for (unsigned int i = 0; i < n_ints; ++i)
+ in >> tmp_int;
}
+ else
+ AssertThrow(
+ false,
+ ExcMessage(std::string(
+ "While reading VTK file, missing CELL_TYPES section. Found <" +
+ keyword + "> instead.")));
+
+ // Ignore everything up to CELL_DATA
+ while (in >> keyword)
+ if (keyword == "CELL_DATA")
+ {
+ unsigned int n_ids;
+ in >> n_ids;
+
+ AssertThrow(n_ids == n_geometric_objects,
+ ExcMessage("The VTK reader found a CELL_DATA statement "
+ "that lists a total of " +
+ Utilities::int_to_string(n_ids) +
+ " cell data objects, but this needs to "
+ "equal the number of cells (which is " +
+ Utilities::int_to_string(cells.size()) +
+ ") plus the number of quads (" +
+ Utilities::int_to_string(
+ subcelldata.boundary_quads.size()) +
+ " in 3d or the number of lines (" +
+ Utilities::int_to_string(
+ subcelldata.boundary_lines.size()) +
+ ") in 2d."));
+
+ const std::vector<std::string> data_sets{"MaterialID",
+ "ManifoldID"};
+
+ for (unsigned int i = 0; i < data_sets.size(); ++i)
+ {
+ // Ignore everything until we get to a SCALARS data set
+ while (in >> keyword)
+ if (keyword == "SCALARS")
+ {
+ // Now see if we know about this type of data set,
+ // if not, just ignore everything till the next SCALARS
+ // keyword
+ std::string set = "";
+ in >> keyword;
+ for (auto set_cmp : data_sets)
+ if (keyword == set_cmp)
+ set = keyword;
+ if (set == "")
+ // keep ignoring everything until the next SCALARS
+ // keyword
+ continue;
+
+ // Now we got somewhere. Proceed from here.
+ // Ignore everything till the end of the line.
+ // SCALARS MaterialID 1 // this "1" is optional
+ in.ignore(256, '\n');
+
+ in >> keyword;
+ AssertThrow(
+ keyword == "LOOKUP_TABLE",
+ ExcMessage(
+ "While reading VTK file, missing keyword LOOKUP_TABLE"));
- ////////////////////////Processing the CELL_DATA
- /// section/////////////////////////////
-
- if (keyword == "CELL_DATA")
- {
- unsigned int n_ids;
- in >> n_ids;
-
- AssertThrow(
- n_ids == cells.size() +
- (dim == 3 ?
- subcelldata.boundary_quads.size() :
- (dim == 2 ? subcelldata.boundary_lines.size() : 0)),
- ExcMessage(
- "The VTK reader found a CELL_DATA statement "
- "that lists a total of " +
- Utilities::int_to_string(n_ids) +
- " cell data objects, but this needs to "
- "equal the number of cells (which is " +
- Utilities::int_to_string(cells.size()) +
- ") plus the number of quads (" +
- Utilities::int_to_string(subcelldata.boundary_quads.size()) +
- " in 3d or the number of lines (" +
- Utilities::int_to_string(subcelldata.boundary_lines.size()) +
- ") in 2d."));
-
-
- std::string linenew;
- std::string textnew[2];
- textnew[0] = "SCALARS MaterialID double";
- textnew[1] = "LOOKUP_TABLE default";
-
- in.ignore(256, '\n');
-
- for (unsigned int i = 0; i < 2; i++)
- {
- getline(in, linenew);
- if (i == 0)
- if (linenew.size() > textnew[0].size())
- linenew.resize(textnew[0].size());
-
- AssertThrow(linenew.compare(textnew[i]) == 0,
- ExcMessage(
- std::string(
- "While reading VTK file, failed to find <") +
- textnew[i] + "> section"));
- }
-
- // read material ids first for all cells, then for all
- // faces. the assumption that cells come before all faces
- // has been verified above via an assertion, so the order
- // used in the following blocks makes sense
- for (unsigned int i = 0; i < cells.size(); i++)
- {
- double id;
- in >> id;
- cells[i].material_id = id;
- }
-
- if (dim == 3)
- {
- for (unsigned int i = 0; i < subcelldata.boundary_quads.size();
- i++)
- {
- double id;
- in >> id;
- subcelldata.boundary_quads[i].material_id = id;
- }
- }
- else if (dim == 2)
- {
- for (unsigned int i = 0; i < subcelldata.boundary_lines.size();
- i++)
- {
- double id;
- in >> id;
- subcelldata.boundary_lines[i].material_id = id;
- }
- }
- }
+ in >> keyword;
+ AssertThrow(
+ keyword == "default",
+ ExcMessage(
+ "While reading VTK file, missing keyword default"));
+
+ // read material or manifold ids first for all cells,
+ // then for all faces, and finally for all lines. the
+ // assumption that cells come before all faces and
+ // lines has been verified above via an assertion, so
+ // the order used in the following blocks makes sense
+ for (unsigned int i = 0; i < cells.size(); i++)
+ {
+ double id;
+ in >> id;
+ if (set == "MaterialID")
+ cells[i].material_id = id;
+ else if (set == "ManifoldID")
+ cells[i].manifold_id = id;
+ else
+ Assert(false, ExcInternalError());
+ }
+
+ if (dim == 3)
+ {
+ for (unsigned int i = 0;
+ i < subcelldata.boundary_quads.size();
+ i++)
+ {
+ double id;
+ in >> id;
+ if (set == "MaterialID")
+ subcelldata.boundary_quads[i].material_id = id;
+ else if (set == "ManifoldID")
+ subcelldata.boundary_quads[i].manifold_id = id;
+ else
+ Assert(false, ExcInternalError());
+ }
+ for (unsigned int i = 0;
+ i < subcelldata.boundary_lines.size();
+ i++)
+ {
+ double id;
+ in >> id;
+ if (set == "MaterialID")
+ subcelldata.boundary_lines[i].material_id = id;
+ else if (set == "ManifoldID")
+ subcelldata.boundary_lines[i].manifold_id = id;
+ else
+ Assert(false, ExcInternalError());
+ }
+ }
+ else if (dim == 2)
+ {
+ for (unsigned int i = 0;
+ i < subcelldata.boundary_lines.size();
+ i++)
+ {
+ double id;
+ in >> id;
+ if (set == "MaterialID")
+ subcelldata.boundary_lines[i].material_id = id;
+ else if (set == "ManifoldID")
+ subcelldata.boundary_lines[i].manifold_id = id;
+ else
+ Assert(false, ExcInternalError());
+ }
+ }
+ }
+ }
+ }
Assert(subcelldata.check_consistency(dim), ExcInternalError());
}
-
template <int dim, int spacedim>
void
GridIn<dim, spacedim>::read_unv(std::istream &in)