vtkreader.impl.hh 20.2 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1
2
#include <fstream>
#include <iterator>
3
#include <sstream>
Praetorius, Simon's avatar
Praetorius, Simon committed
4
5
#include <string>

6
#if HAVE_VTK_ZLIB
Praetorius, Simon's avatar
Praetorius, Simon committed
7
8
9
#include <zlib.h>
#endif

10
#include <dune/common/classname.hh>
11
#include <dune/common/version.hh>
12

13
#include "utility/errors.hh"
Praetorius, Simon's avatar
Praetorius, Simon committed
14
15
16
#include "utility/filesystem.hh"
#include "utility/string.hh"

17
namespace Dune {
Praetorius, Simon's avatar
Praetorius, Simon committed
18

19
template <class Grid, class Creator>
20
void VtkReader<Grid,Creator>::read (std::string const& filename, bool fillCreator)
Praetorius, Simon's avatar
Praetorius, Simon committed
21
22
{
  // check whether file exists!
Stenger, Florian's avatar
Stenger, Florian committed
23
  if (!Vtk::exists(filename))
Praetorius, Simon's avatar
Praetorius, Simon committed
24
25
26
    DUNE_THROW(IOError, "File " << filename << " does not exist!");

  std::ifstream input(filename, std::ios_base::in | std::ios_base::binary);
27
  assert(input.is_open());
Praetorius, Simon's avatar
Praetorius, Simon committed
28

29
  std::string ext = Vtk::Path(filename).extension().string();
30
  if (ext == ".vtu") {
31
    readSerialFileFromStream(input, fillCreator);
32
    pieces_.push_back(filename);
33
  } else if (ext == ".pvtu") {
34
    readParallelFileFromStream(input, comm().rank(), comm().size(), fillCreator);
35
  } else {
36
    DUNE_THROW(Dune::VtkError, "File has unknown file-extension '" << ext << "'. Allowed are only '.vtu' and '.pvtu'.");
37
38
39
40
41
  }
}


template <class Grid, class Creator>
42
void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bool fillCreator)
43
{
Praetorius, Simon's avatar
Praetorius, Simon committed
44
  clear();
45
46
  std::string compressor = "";
  std::string data_name = "", data_format = "";
Praetorius, Simon's avatar
Praetorius, Simon committed
47
48
  Vtk::DataTypes data_type = Vtk::UNKNOWN;
  std::size_t data_components = 0;
49
  std::uint64_t data_offset = 0;
Praetorius, Simon's avatar
Praetorius, Simon committed
50
51
52

  Sections section = NO_SECTION;
  for (std::string line; std::getline(input, line); ) {
Stenger, Florian's avatar
Stenger, Florian committed
53
    Vtk::ltrim(line);
Praetorius, Simon's avatar
Praetorius, Simon committed
54
55

    if (isSection(line, "VTKFile", section)) {
56
57
      bool closed = false;
      auto attr = parseXml(line, closed);
58
59

      if (!attr["type"].empty())
60
        VTK_ASSERT(attr["type"] == "UnstructuredGrid", "VtkReader supports UnstructuredGrid types");
Praetorius, Simon's avatar
Praetorius, Simon committed
61
      if (!attr["version"].empty())
62
        VTK_ASSERT(std::stod(attr["version"]) == 1.0, "File format must be 1.0");
63
      if (!attr["byte_order"].empty())
64
        VTK_ASSERT(attr["byte_order"] == "LittleEndian", "LittleEndian byte order supported");
65
      if (!attr["header_type"].empty())
66
        VTK_ASSERT(attr["header_type"] == "UInt64", "Header integer type must be UInt64");
67
      if (!attr["compressor"].empty()) {
Praetorius, Simon's avatar
Praetorius, Simon committed
68
        compressor = attr["compressor"];
69
        VTK_ASSERT(compressor == "vtkZLibDataCompressor", "Only ZLib compression supported");
70
      }
Praetorius, Simon's avatar
Praetorius, Simon committed
71
72
73
74
75
76
77
78
79
      section = VTK_FILE;
    }
    else if (isSection(line, "/VTKFile", section, VTK_FILE))
      section = NO_SECTION;
    else if (isSection(line, "UnstructuredGrid", section, VTK_FILE))
      section = UNSTRUCTURED_GRID;
    else if (isSection(line, "/UnstructuredGrid", section, UNSTRUCTURED_GRID))
      section = VTK_FILE;
    else if (isSection(line, "Piece", section, UNSTRUCTURED_GRID)) {
80
81
      bool closed = false;
      auto attr = parseXml(line, closed);
82

83
      VTK_ASSERT(attr.count("NumberOfPoints") > 0 && attr.count("NumberOfCells") > 0, "Number of points or cells in file must be > 0");
84
85
      numberOfPoints_ = std::stoul(attr["NumberOfPoints"]);
      numberOfCells_ = std::stoul(attr["NumberOfCells"]);
Praetorius, Simon's avatar
Praetorius, Simon committed
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
      section = PIECE;
    }
    else if (isSection(line, "/Piece", section, PIECE))
      section = UNSTRUCTURED_GRID;
    else if (isSection(line, "PointData", section, PIECE))
      section = POINT_DATA;
    else if (isSection(line, "/PointData", section, POINT_DATA))
      section = PIECE;
    else if (isSection(line, "CellData", section, PIECE))
      section = CELL_DATA;
    else if (isSection(line, "/CellData", section, CELL_DATA))
      section = PIECE;
    else if (isSection(line, "Points", section, PIECE))
      section = POINTS;
    else if (isSection(line, "/Points", section, POINTS))
      section = PIECE;
    else if (isSection(line, "Cells", section, PIECE))
      section = CELLS;
    else if (isSection(line, "/Cells", section, CELLS))
      section = PIECE;
    else if (line.substr(1,9) == "DataArray") {
107
108
      bool closed = false;
      auto attr = parseXml(line, closed);
Praetorius, Simon's avatar
Praetorius, Simon committed
109

110
      data_type = Vtk::Map::to_datatype[attr["type"]];
111
112

      if (!attr["Name"].empty())
Stenger, Florian's avatar
Stenger, Florian committed
113
        data_name = Vtk::to_lower(attr["Name"]);
114
115
116
117
      else if (section == POINTS)
        data_name = "points";

      data_components = 1;
Praetorius, Simon's avatar
Praetorius, Simon committed
118
119
120
121
      if (!attr["NumberOfComponents"].empty())
        data_components = std::stoul(attr["NumberOfComponents"]);

      // determine FormatType
Stenger, Florian's avatar
Stenger, Florian committed
122
      data_format = Vtk::to_lower(attr["format"]);
Praetorius, Simon's avatar
Praetorius, Simon committed
123
      if (data_format == "appended") {
124
        format_ = !compressor.empty() ? Vtk::COMPRESSED : Vtk::BINARY;
Praetorius, Simon's avatar
Praetorius, Simon committed
125
126
127
128
      } else {
        format_ = Vtk::ASCII;
      }

129
130
      // Offset makes sense in appended mode only
      data_offset = 0;
Praetorius, Simon's avatar
Praetorius, Simon committed
131
132
      if (!attr["offset"].empty()) {
        data_offset = std::stoul(attr["offset"]);
133
        VTK_ASSERT(data_format == "appended", "Attribute 'offset' only supported by appended mode");
Praetorius, Simon's avatar
Praetorius, Simon committed
134
135
      }

136
137
138
139
140
141
142
      // Store attributes of DataArray
      dataArray_[data_name] = {data_type, data_components, data_offset};

      // Skip section in appended mode
      if (data_format == "appended") {
        if (!closed) {
          while (std::getline(input, line)) {
Stenger, Florian's avatar
Stenger, Florian committed
143
            Vtk::ltrim(line);
144
            if (line.substr(1,10) == "/DataArray")
145
146
147
148
149
150
151
              break;
          }
        }
        continue;
      }

      if (section == POINT_DATA)
Praetorius, Simon's avatar
Praetorius, Simon committed
152
153
154
155
156
157
158
159
        section = PD_DATA_ARRAY;
      else if (section == POINTS)
        section = POINTS_DATA_ARRAY;
      else if (section == CELL_DATA)
        section = CD_DATA_ARRAY;
      else if (section == CELLS)
        section = CELLS_DATA_ARRAY;
      else
160
        DUNE_THROW(Dune::VtkError, "Wrong section for <DataArray>");
Praetorius, Simon's avatar
Praetorius, Simon committed
161
162
163
164
165
166
167
168
169
170
171
    }
    else if (line.substr(1,10) == "/DataArray") {
      if (section == PD_DATA_ARRAY)
        section = POINT_DATA;
      else if (section == POINTS_DATA_ARRAY)
        section = POINTS;
      else if (section == CD_DATA_ARRAY)
        section = CELL_DATA;
      else if (section == CELLS_DATA_ARRAY)
        section = CELLS;
      else
172
        DUNE_THROW(Dune::VtkError, "Wrong section for </DataArray>");
Praetorius, Simon's avatar
Praetorius, Simon committed
173
174
    }
    else if (isSection(line, "AppendedData", section, VTK_FILE)) {
175
176
      bool closed = false;
      auto attr = parseXml(line, closed);
177
      if (!attr["encoding"].empty())
178
        VTK_ASSERT(attr["encoding"] == "raw", "base64 encoding not supported");
179

180
181
      offset0_ = findAppendedDataPosition(input);
      if (dataArray_["points"].type == Vtk::FLOAT32)
182
183
184
185
186
187
        readPointsAppended<float>(input);
      else
        readPointsAppended<double>(input);

      readCellsAppended(input);
      section = NO_SECTION; // finish reading after appended section
Praetorius, Simon's avatar
Praetorius, Simon committed
188
189
190
191
192
193
    }
    else if (isSection(line, "/AppendedData", section, APPENDED_DATA))
      section = VTK_FILE;

    switch (section) {
      case PD_DATA_ARRAY:
194
195
        if (data_type == Vtk::FLOAT32) {
          std::vector<float> values;
196
197
          section = readPointData(input, values, data_name);
        } else if (data_type == Vtk::FLOAT64) {
198
          std::vector<double> values;
199
          section = readPointData(input, values, data_name);
200
        }
Praetorius, Simon's avatar
Praetorius, Simon committed
201
202
        break;
      case POINTS_DATA_ARRAY:
203
        section = readPoints(input, data_name);
Praetorius, Simon's avatar
Praetorius, Simon committed
204
205
        break;
      case CD_DATA_ARRAY:
206
207
        if (data_type == Vtk::FLOAT32) {
          std::vector<float> values;
208
209
          section = readCellData(input, values, data_name);
        } else if (data_type == Vtk::FLOAT64) {
210
          std::vector<double> values;
211
          section = readCellData(input, values, data_name);
212
        }
Praetorius, Simon's avatar
Praetorius, Simon committed
213
214
        break;
      case CELLS_DATA_ARRAY:
215
        section = readCells(input, data_name);
Praetorius, Simon's avatar
Praetorius, Simon committed
216
217
218
219
220
221
222
223
224
225
226
        break;
      default:
        // do nothing
        break;
    }

    if (section == NO_SECTION)
      break;
  }

  if (section != NO_SECTION)
227
    DUNE_THROW(Dune::VtkError, "VTK-File is incomplete. It must end with </VTKFile>!");
Praetorius, Simon's avatar
Praetorius, Simon committed
228

229
230
  if (fillCreator)
    fillGridCreator();
231
232
233
234
}


template <class Grid, class Creator>
235
void VtkReader<Grid,Creator>::readParallelFileFromStream (std::ifstream& input, int commRank, int commSize, bool fillCreator)
236
{
Praetorius, Simon's avatar
Praetorius, Simon committed
237
  clear();
238
239
240

  Sections section = NO_SECTION;
  for (std::string line; std::getline(input, line); ) {
Stenger, Florian's avatar
Stenger, Florian committed
241
    Vtk::ltrim(line);
242
243
244
245
246
247

    if (isSection(line, "VTKFile", section)) {
      bool closed = false;
      auto attr = parseXml(line, closed);

      if (!attr["type"].empty())
248
        VTK_ASSERT(attr["type"] == "PUnstructuredGrid", "VtkReader supports PUnstructuredGrid types");
249
      if (!attr["version"].empty())
250
        VTK_ASSERT(std::stod(attr["version"]) == 1.0, "File format must be 1.0");
251
      if (!attr["byte_order"].empty())
252
        VTK_ASSERT(attr["byte_order"] == "LittleEndian", "LittleEndian byte order supported");
253
      if (!attr["header_type"].empty())
254
        VTK_ASSERT(attr["header_type"] == "UInt64", "Header integer type must be UInt64");
255
      if (!attr["compressor"].empty())
256
        VTK_ASSERT(attr["compressor"] == "vtkZLibDataCompressor", "Only ZLib compression supported");
257
258
259
260
261
262
263
264
265
266
267
268
      section = VTK_FILE;
    }
    else if (isSection(line, "/VTKFile", section, VTK_FILE))
      section = NO_SECTION;
    else if (isSection(line, "PUnstructuredGrid", section, VTK_FILE))
      section = UNSTRUCTURED_GRID;
    else if (isSection(line, "/PUnstructuredGrid", section, UNSTRUCTURED_GRID))
      section = VTK_FILE;
    else if (isSection(line, "Piece", section, UNSTRUCTURED_GRID)) {
      bool closed = false;
      auto attr = parseXml(line, closed);

269
      VTK_ASSERT(attr.count("Source") > , "No source files for partitions provided");
270
271
272
273
274
275
276
277
      pieces_.push_back(attr["Source"]);
    }

    if (section == NO_SECTION)
      break;
  }

  if (section != NO_SECTION)
278
    DUNE_THROW(Dune::VtkError, "VTK-File is incomplete. It must end with </VTKFile>!");
279

280
281
  if (fillCreator)
    fillGridCreator();
Praetorius, Simon's avatar
Praetorius, Simon committed
282
283
284
}


285
286
287
288
289
290
291
// @{ implementation detail
/**
 * Read ASCII data from `input` stream into vector `values`
 * \param max_size  Upper bound for the number of values
 * \param section   Current XML section you are reading in
 * \param parent_section   XML Section to return when current `section` is finished.
 **/
Praetorius, Simon's avatar
Praetorius, Simon committed
292
template <class IStream, class T, class Sections>
293
294
Sections readDataArray (IStream& input, std::vector<T>& values, std::size_t max_size,
                        Sections section, Sections parent_section)
Praetorius, Simon's avatar
Praetorius, Simon committed
295
{
Praetorius, Simon's avatar
Praetorius, Simon committed
296
  values.reserve(max_size < std::size_t(-1) ? max_size : 0);
297
  using S = std::conditional_t<(sizeof(T) <= 1), std::uint16_t, T>; // problem when reading chars as ints
Praetorius, Simon's avatar
Praetorius, Simon committed
298
299

  std::size_t idx = 0;
300
  for (std::string line; std::getline(input, line);) {
Stenger, Florian's avatar
Stenger, Florian committed
301
    Vtk::trim(line);
302
    if (line.substr(1,10) == "/DataArray")
Praetorius, Simon's avatar
Praetorius, Simon committed
303
      return parent_section;
Praetorius, Simon's avatar
Praetorius, Simon committed
304
305
    if (line[0] == '<')
      break;
Praetorius, Simon's avatar
Praetorius, Simon committed
306
307
308
309
310

    std::istringstream stream(line);
    S value;
    for (; stream >> value; idx++)
      values.push_back(T(value));
Praetorius, Simon's avatar
Praetorius, Simon committed
311
312
313
314
315
316
317
318
319
320
321
    if (idx >= max_size)
      break;
  }

  return section;
}

template <class IStream, class Sections>
Sections skipRestOfDataArray (IStream& input, Sections section, Sections parent_section)
{
  for (std::string line; std::getline(input, line);) {
Stenger, Florian's avatar
Stenger, Florian committed
322
    Vtk::ltrim(line);
Praetorius, Simon's avatar
Praetorius, Simon committed
323
324
    if (line.substr(1,10) == "/DataArray")
      return parent_section;
Praetorius, Simon's avatar
Praetorius, Simon committed
325
326
327
328
  }

  return section;
}
329
// @}
Praetorius, Simon's avatar
Praetorius, Simon committed
330
331


332
333
template <class Grid, class Creator>
typename VtkReader<Grid,Creator>::Sections
334
VtkReader<Grid,Creator>::readPoints (std::ifstream& input, std::string name)
Praetorius, Simon's avatar
Praetorius, Simon committed
335
336
{
  using T = typename GlobalCoordinate::value_type;
337
338
339
  assert(numberOfPoints_ > 0);
  assert(dataArray_["points"].components == 3u);

340
341
  Sections sec;

Praetorius, Simon's avatar
Praetorius, Simon committed
342
  std::vector<T> point_values;
343
  sec = readDataArray(input, point_values, 3*numberOfPoints_, POINTS_DATA_ARRAY, POINTS);
Praetorius, Simon's avatar
Praetorius, Simon committed
344
345
  if (sec != POINTS)
    sec = skipRestOfDataArray(input, POINTS_DATA_ARRAY, POINTS);
Praetorius, Simon's avatar
Praetorius, Simon committed
346
  assert(sec == POINTS);
347
  assert(point_values.size() == 3*numberOfPoints_);
Praetorius, Simon's avatar
Praetorius, Simon committed
348
349
350

  // extract points from continuous values
  GlobalCoordinate p;
351
  vec_points.reserve(numberOfPoints_);
Praetorius, Simon's avatar
Praetorius, Simon committed
352
  std::size_t idx = 0;
353
  for (std::size_t i = 0; i < numberOfPoints_; ++i) {
Praetorius, Simon's avatar
Praetorius, Simon committed
354
355
356
    for (std::size_t j = 0; j < p.size(); ++j)
      p[j] = point_values[idx++];
    idx += (3u - p.size());
357
    vec_points.push_back(p);
Praetorius, Simon's avatar
Praetorius, Simon committed
358
359
360
361
362
363
  }

  return sec;
}


364
template <class Grid, class Creator>
Praetorius, Simon's avatar
Praetorius, Simon committed
365
  template <class T>
366
void VtkReader<Grid,Creator>::readPointsAppended (std::ifstream& input)
Praetorius, Simon's avatar
Praetorius, Simon committed
367
{
368
369
  assert(numberOfPoints_ > 0);
  assert(dataArray_["points"].components == 3u);
Praetorius, Simon's avatar
Praetorius, Simon committed
370
  std::vector<T> point_values;
371
372
  readAppended(input, point_values, dataArray_["points"].offset);
  assert(point_values.size() == 3*numberOfPoints_);
Praetorius, Simon's avatar
Praetorius, Simon committed
373
374
375

  // extract points from continuous values
  GlobalCoordinate p;
376
  vec_points.reserve(numberOfPoints_);
Praetorius, Simon's avatar
Praetorius, Simon committed
377
  std::size_t idx = 0;
378
  for (std::size_t i = 0; i < numberOfPoints_; ++i) {
Praetorius, Simon's avatar
Praetorius, Simon committed
379
380
381
    for (std::size_t j = 0; j < p.size(); ++j)
      p[j] = T(point_values[idx++]);
    idx += (3u - p.size());
382
    vec_points.push_back(p);
Praetorius, Simon's avatar
Praetorius, Simon committed
383
384
385
386
  }
}


387
388
389
template <class Grid, class Creator>
typename VtkReader<Grid,Creator>::Sections
VtkReader<Grid,Creator>::readCells (std::ifstream& input, std::string name)
Praetorius, Simon's avatar
Praetorius, Simon committed
390
391
392
{
  Sections sec = CELLS_DATA_ARRAY;

393
  assert(numberOfCells_ > 0);
Praetorius, Simon's avatar
Praetorius, Simon committed
394
  if (name == "types") {
395
396
    sec = readDataArray(input, vec_types, numberOfCells_, CELLS_DATA_ARRAY, CELLS);
    assert(vec_types.size() == numberOfCells_);
Praetorius, Simon's avatar
Praetorius, Simon committed
397
  } else if (name == "offsets") {
398
399
    sec = readDataArray(input, vec_offsets, numberOfCells_, CELLS_DATA_ARRAY, CELLS);
    assert(vec_offsets.size() == numberOfCells_);
Praetorius, Simon's avatar
Praetorius, Simon committed
400
  } else if (name == "connectivity") {
Praetorius, Simon's avatar
Praetorius, Simon committed
401
    sec = readDataArray(input, vec_connectivity, std::size_t(-1), CELLS_DATA_ARRAY, CELLS);
402
403
404
  } else if (name == "global_point_ids") {
    sec = readDataArray(input, vec_point_ids, numberOfPoints_, CELLS_DATA_ARRAY, CELLS);
    assert(vec_point_ids.size() == numberOfPoints_);
Praetorius, Simon's avatar
Praetorius, Simon committed
405
406
407
408
409
410
  }

  return sec;
}


411
412
template <class Grid, class Creator>
void VtkReader<Grid,Creator>::readCellsAppended (std::ifstream& input)
Praetorius, Simon's avatar
Praetorius, Simon committed
413
{
414
415
  assert(numberOfCells_ > 0);
  auto types_data = dataArray_["types"];
416
  auto offsets_data = dataArray_["offsets"];
417
  auto connectivity_data = dataArray_["connectivity"];
Praetorius, Simon's avatar
Praetorius, Simon committed
418

419
420
421
  assert(types_data.type == Vtk::UINT8);
  readAppended(input, vec_types, types_data.offset);
  assert(vec_types.size() == numberOfCells_);
Praetorius, Simon's avatar
Praetorius, Simon committed
422

423
424
425
426
427
428
429
430
431
  if (offsets_data.type == Vtk::INT64)
    readAppended(input, vec_offsets, offsets_data.offset);
  else if (offsets_data.type == Vtk::INT32) {
    std::vector<std::int32_t> offsets;
    readAppended(input, offsets, offsets_data.offset);
    vec_offsets.resize(offsets.size());
    std::copy(offsets.begin(), offsets.end(), vec_offsets.begin());
  }
  else { DUNE_THROW(Dune::NotImplemented, "Unsupported DataType in Cell offsets."); }
432
  assert(vec_offsets.size() == numberOfCells_);
Praetorius, Simon's avatar
Praetorius, Simon committed
433

434
435
436
437
438
439
440
441
442
  if (connectivity_data.type == Vtk::INT64)
    readAppended(input, vec_connectivity, connectivity_data.offset);
  else if (connectivity_data.type == Vtk::INT32) {
    std::vector<std::int32_t> connectivity;
    readAppended(input, connectivity, connectivity_data.offset);
    vec_connectivity.resize(connectivity.size());
    std::copy(connectivity.begin(), connectivity.end(), vec_connectivity.begin());
  }
  else { DUNE_THROW(Dune::NotImplemented, "Unsupported DataType in Cell connectivity."); }
443
  assert(vec_connectivity.size() == std::size_t(vec_offsets.back()));
444
445
446
447
448
449
450

  if (dataArray_.count("global_point_ids") > 0) {
    auto point_id_data = dataArray_["global_point_ids"];
    assert(point_id_data.type == Vtk::UINT64);
    readAppended(input, vec_point_ids, point_id_data.offset);
    assert(vec_point_ids.size() == numberOfPoints_);
  }
Praetorius, Simon's avatar
Praetorius, Simon committed
451
452
453
}


454
455
456
457
458
459
460
461
// @{ implementation detail
/**
 * Read compressed data into `buffer_in`, uncompress it and store the result in
 * the concrete-data-type `buffer`
 * \param bs     Size of the uncompressed data
 * \param cbs    Size of the compressed data
 * \param input  Stream to read from.
 **/
Praetorius, Simon's avatar
Praetorius, Simon committed
462
463
464
465
template <class T, class IStream>
void read_compressed (T* buffer, unsigned char* buffer_in,
                      std::uint64_t bs, std::uint64_t cbs, IStream& input)
{
466
#if HAVE_VTK_ZLIB
Praetorius, Simon's avatar
Praetorius, Simon committed
467
468
469
470
471
472
  uLongf uncompressed_space = uLongf(bs);
  uLongf compressed_space = uLongf(cbs);

  Bytef* compressed_buffer = reinterpret_cast<Bytef*>(buffer_in);
  Bytef* uncompressed_buffer = reinterpret_cast<Bytef*>(buffer);

473
  input.read((char*)(compressed_buffer), compressed_space);
474
  assert(uLongf(input.gcount()) == compressed_space);
Praetorius, Simon's avatar
Praetorius, Simon committed
475
476
477
478
479

  if (uncompress(uncompressed_buffer, &uncompressed_space, compressed_buffer, compressed_space) != Z_OK) {
    std::cerr << "Zlib error while uncompressing data.\n";
    std::abort();
  }
480
  assert(uLongf(bs) == uncompressed_space);
Praetorius, Simon's avatar
Praetorius, Simon committed
481
482
483
484
485
#else
  std::cerr << "Can not call read_compressed without compression enabled!\n";
  std::abort();
#endif
}
486
// @}
Praetorius, Simon's avatar
Praetorius, Simon committed
487
488


489
template <class Grid, class Creator>
Praetorius, Simon's avatar
Praetorius, Simon committed
490
  template <class T>
491
void VtkReader<Grid,Creator>::readAppended (std::ifstream& input, std::vector<T>& values, std::uint64_t offset)
Praetorius, Simon's avatar
Praetorius, Simon committed
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
{
  input.seekg(offset0_ + offset);

  std::uint64_t size = 0;

  std::uint64_t num_blocks = 0;
  std::uint64_t block_size = 0;
  std::uint64_t last_block_size = 0;
  std::vector<std::uint64_t> cbs; // compressed block sizes

  // read total size / block-size(s)
  if (format_ == Vtk::COMPRESSED) {
    input.read((char*)&num_blocks, sizeof(std::uint64_t));
    input.read((char*)&block_size, sizeof(std::uint64_t));
    input.read((char*)&last_block_size, sizeof(std::uint64_t));

    assert(block_size % sizeof(T) == 0);

    // total size of the uncompressed data
    size = block_size * (num_blocks-1) + last_block_size;

    // size of the compressed blocks
    cbs.resize(num_blocks);
    input.read((char*)cbs.data(), num_blocks*sizeof(std::uint64_t));
  } else {
    input.read((char*)&size, sizeof(std::uint64_t));
  }
519
  assert(size > 0 && (size % sizeof(T)) == 0);
Praetorius, Simon's avatar
Praetorius, Simon committed
520
521
522
523
524
525
526
527
528
529
530
531
532
533
  values.resize(size / sizeof(T));

  if (format_ == Vtk::COMPRESSED) {
    // upper bound for compressed block-size
    std::uint64_t compressed_block_size = block_size + (block_size + 999)/1000 + 12;
    // number of values in the full blocks
    std::size_t num_values = block_size / sizeof(T);

    std::vector<unsigned char> buffer_in(compressed_block_size);
    for (std::size_t i = 0; i < std::size_t(num_blocks); ++i) {
      std::uint64_t bs = i < std::size_t(num_blocks-1) ? block_size : last_block_size;
      read_compressed(values.data() + i*num_values, buffer_in.data(), bs, cbs[i], input);
    }
  } else {
534
    input.read((char*)(values.data()), size);
535
    assert(input.gcount() == std::streamsize(size));
Praetorius, Simon's avatar
Praetorius, Simon committed
536
537
538
539
  }
}


540
template <class Grid, class Creator>
541
void VtkReader<Grid,Creator>::fillGridCreator (bool insertPieces)
Praetorius, Simon's avatar
Praetorius, Simon committed
542
{
543
544
545
  assert(vec_points.size() == numberOfPoints_);
  assert(vec_types.size() == numberOfCells_);
  assert(vec_offsets.size() == numberOfCells_);
Praetorius, Simon's avatar
Praetorius, Simon committed
546

Praetorius, Simon's avatar
Praetorius, Simon committed
547
  if (!vec_points.empty())
548
    creator_->insertVertices(vec_points, vec_point_ids);
Praetorius, Simon's avatar
Praetorius, Simon committed
549
  if (!vec_types.empty())
550
    creator_->insertElements(vec_types, vec_offsets, vec_connectivity);
551
  if (insertPieces)
552
    creator_->insertPieces(pieces_);
553
}
Praetorius, Simon's avatar
Praetorius, Simon committed
554

555
556
557
// Assume input already read the line <AppendedData ...>
template <class Grid, class Creator>
std::uint64_t VtkReader<Grid,Creator>::findAppendedDataPosition (std::ifstream& input) const
558
559
560
561
562
563
564
565
566
567
568
569
{
  char c;
  while (input.get(c) && std::isblank(c)) { /*do nothing*/ }

  std::uint64_t offset = input.tellg();
  if (c != '_')
    --offset; // if char is not '_', assume it is part of the data.

  return offset;
}


570
571
template <class Grid, class Creator>
std::map<std::string, std::string> VtkReader<Grid,Creator>::parseXml (std::string const& line, bool& closed)
Praetorius, Simon's avatar
Praetorius, Simon committed
572
{
573
  closed = false;
Praetorius, Simon's avatar
Praetorius, Simon committed
574
575
576
577
578
579
580
581
582
583
584
585
586
587
  std::map<std::string, std::string> attr;

  Sections sec = NO_SECTION;
  bool escape = false;

  std::string name = "";
  std::string value = "";
  for (auto c : line) {
    switch (sec) {
    case NO_SECTION:
      if (std::isalpha(c) || c == '_') {
        name.clear();
        sec = XML_NAME;
        name.push_back(c);
588
589
      } else if (c == '/') {
        closed = true;
Praetorius, Simon's avatar
Praetorius, Simon committed
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
      }
      break;
    case XML_NAME:
      if (std::isalpha(c) || c == '_')
        name.push_back(c);
      else
        sec = (c == '=' ? XML_NAME_ASSIGN : NO_SECTION);
      break;
    case XML_NAME_ASSIGN:
      value.clear();
      escape = false;
      assert( c == '"' && "Format error!" );
      sec = XML_VALUE;
      break;
    case XML_VALUE:
      if (c == '"' && !escape) {
        attr[name] = value;
        sec = NO_SECTION;
      } else if (c == '\\' && !escape) {
        escape = true;
      }  else {
        value.push_back(c);
        escape = false;
      }
      break;
    default:
      assert(false && "Format error!");
    }
  }

  return attr;
}

Praetorius, Simon's avatar
Praetorius, Simon committed
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639

template <class Grid, class Creator>
void VtkReader<Grid,Creator>::clear ()
{
  vec_points.clear();
  vec_point_ids.clear();
  vec_types.clear();
  vec_offsets.clear();
  vec_connectivity.clear();
  dataArray_.clear();
  pieces_.clear();

  numberOfCells_ = 0;
  numberOfPoints_ = 0;
  offset0_ = 0;
}

640
} // end namespace Dune