vtkreader.hh 11.9 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1
2
3
4
#pragma once

#include <iosfwd>
#include <map>
5
6
7
#include <memory>
#include <vector>

Praetorius, Simon's avatar
Praetorius, Simon committed
8
#include <dune/common/shared_ptr.hh>
9
#include <dune/common/typelist.hh>
10
#include <dune/common/typeutilities.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
11

Praetorius, Simon's avatar
Praetorius, Simon committed
12
#include <dune/vtk/filereader.hh>
Stenger, Florian's avatar
Stenger, Florian committed
13
#include <dune/vtk/types.hh>
14
#include <dune/vtk/utility/errors.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
15
16
17

// default GridCreator
#include <dune/vtk/gridcreators/continuousgridcreator.hh>
18
#include <dune/vtk/gridfunctions/common.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
19

20
namespace Dune
Praetorius, Simon's avatar
Praetorius, Simon committed
21
{
22
  /// File-Reader for Vtk unstructured .vtu files
23
24
25
26
  /**
   * Reads .vtu files and constructs a grid from the cells stored in the file
   * Additionally, stored data can be read.
   *
27
28
   * NOTE: Assumption on the file structure: Each XML tag must be on a separate line.
   *
Praetorius, Simon's avatar
Praetorius, Simon committed
29
30
31
32
   * \tparam Grid       The type of the grid to construct.
   * \tparam GC         GridCreator policy type to control what to pass to a grid factory with
   *                    data given from the file. [ContinuousGridCreator]
   * \tparam FieldType  Type of the components of the data to extract from the file [default: double]
33
   **/
Praetorius, Simon's avatar
Praetorius, Simon committed
34
  template <class Grid, class GC = Vtk::ContinuousGridCreator<Grid>, class FieldType = double>
Praetorius, Simon's avatar
Praetorius, Simon committed
35
  class VtkReader
Praetorius, Simon's avatar
Praetorius, Simon committed
36
      : public Vtk::FileReader<Grid, VtkReader<Grid, GC>>
Praetorius, Simon's avatar
Praetorius, Simon committed
37
  {
38
    // Sections visited during the xml parsing
Praetorius, Simon's avatar
Praetorius, Simon committed
39
40
41
42
43
    enum Sections {
      NO_SECTION = 0, VTK_FILE, UNSTRUCTURED_GRID, PIECE, POINT_DATA, PD_DATA_ARRAY, CELL_DATA, CD_DATA_ARRAY,
      POINTS, POINTS_DATA_ARRAY, CELLS, CELLS_DATA_ARRAY, APPENDED_DATA, XML_NAME, XML_NAME_ASSIGN, XML_VALUE
    };

44
    // Type storing information about read data
45
46
    struct DataArrayAttributes
    {
47
      std::string name;
48
      Vtk::DataTypes type;
49
      unsigned int components = 1;
50
      std::uint64_t offset = 0;
51
      Sections section = NO_SECTION;
52
53
    };

54
    // Type of global world coordinates
Praetorius, Simon's avatar
Praetorius, Simon committed
55
    using GlobalCoordinate = typename GC::GlobalCoordinate;
Praetorius, Simon's avatar
Praetorius, Simon committed
56

57
58
    // Template representing a grid-function that is created in getPointData() and getCellData()
    // with Context either Vtk::PointContext or Vek::CellContext, respectively.
Praetorius, Simon's avatar
Praetorius, Simon committed
59
    // To each GridCreator a GridFunction is associated, see, e.g. Vtk::ContinuousGridFunction
60
61
    // or Vtk::LagrangeGridFunction.
    template <class Context>
Praetorius, Simon's avatar
Praetorius, Simon committed
62
63
64
65
66
67
68
69
70
71
    using GridFunction = typename Vtk::AssociatedGridFunction<GC, FieldType, Context>::type;

  public:
    using GridCreator = GC;

    /// GridFunction representing the data stored on the points in the file
    using PointGridFunction = GridFunction<Vtk::PointContext>;

    /// GridFunction representing the data stored on the cells in the file
    using CellGridFunction = GridFunction<Vtk::CellContext>;
72

Praetorius, Simon's avatar
Praetorius, Simon committed
73
  public:
Praetorius, Simon's avatar
Praetorius, Simon committed
74
    /// Constructor. Creates a new GridCreator with the passed factory
Praetorius, Simon's avatar
Praetorius, Simon committed
75
76
77
78
79
80
    /**
     * \param args... Either pass a GridFactory by reference or shared_ptr, or a list of arguments
     *                passed to the constructor of a Dune::GridFactory (typically and empty parameter
     *                list). See the constructor of \ref GridCreatorInterface and the GridCreator
     *                passed to this reader.
     **/
81
    template <class... Args,
Praetorius, Simon's avatar
Praetorius, Simon committed
82
      std::enable_if_t<std::is_constructible<GridCreator, Args...>::value,int> = 0>
Praetorius, Simon's avatar
Praetorius, Simon committed
83
    explicit VtkReader (Args&&... args)
84
      : VtkReader(std::make_shared<GridCreator>(std::forward<Args>(args)...))
Praetorius, Simon's avatar
Praetorius, Simon committed
85
86
    {}

87
88
89
    /// Constructor. Stores the references in a non-destroying shared_ptr
    explicit VtkReader (GridCreator& creator)
      : VtkReader(stackobject_to_shared_ptr(creator))
Praetorius, Simon's avatar
Praetorius, Simon committed
90
    {}
Praetorius, Simon's avatar
Praetorius, Simon committed
91

92
93
94
95
    /// Constructor. Stores the shared_ptr
    explicit VtkReader (std::shared_ptr<GridCreator> creator)
      : creator_(std::move(creator))
    {}
96

97
    /// Read the grid from file with `filename` into the GridCreator
98
    /**
99
100
101
102
103
104
     * This function fills internal data containers representing the information from the
     * passed file.
     *
     * \param filename     The name of the input file
     * \param fillCreator  If `false`, only fill internal data structures, if `true`, pass
     *                     the internal data to the GridCreator. [true]
105
     **/
106
    void read (std::string const& filename, bool fillCreator = true);
107

108
    /// Obtains the creator of the reader
Praetorius, Simon's avatar
Praetorius, Simon committed
109
110
111
112
113
    GridCreator& gridCreator ()
    {
      return *creator_;
    }

114
    /// Construct the actual grid using the GridCreator
Praetorius, Simon's avatar
Praetorius, Simon committed
115
    /// [[expects: read_ == true]]
116
117
118
119
120
    std::unique_ptr<Grid> createGrid () const
    {
      return creator_->createGrid();
    }

121
    /// Construct a grid-function representing the point-data with the given name
Praetorius, Simon's avatar
Praetorius, Simon committed
122
    /// [[expects: read_ == true]]
123
124
    GridFunction<Vtk::PointContext> getPointData (std::string const& name) const
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
125
126
127
128
      auto data_it = dataArray_.find("PointData." + name);
      auto point_it = pointData_.find("PointData." + name);
      VTK_ASSERT_MSG(data_it != dataArray_.end() && point_it != pointData_.end(),
        "The data to extract is not found in point-data. Try `getCellData()` instead!");
129
      VTK_ASSERT(data_it->second.section == POINT_DATA);
130

131
132
      return {*creator_, point_it->second, data_it->second.name, data_it->second.components,
              data_it->second.type, vec_types, vec_offsets, vec_connectivity};
133
134
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
135
136
137
138
139
140
141
    /// Return a vector of DataArrayAttributes for all POINT_DATA blocks
    /// [[expects: read_ == true]]
    std::vector<DataArrayAttributes> getPointDataAttributes () const
    {
      std::vector<DataArrayAttributes> attributes;
      attributes.reserve(pointData_.size());
      for (auto const& da : dataArray_) {
142
        if (da.second.section == POINT_DATA)
Praetorius, Simon's avatar
Praetorius, Simon committed
143
144
145
146
147
          attributes.push_back(da.second);
      }
      return attributes;
    }

148
    /// Construct a grid-function representing the cell-data with the given name
Praetorius, Simon's avatar
Praetorius, Simon committed
149
    /// [[expects: read_ == true]]
150
151
    GridFunction<Vtk::CellContext> getCellData (std::string const& name) const
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
152
153
154
155
156
      auto data_it = dataArray_.find("CellData." + name);
      auto cell_it = cellData_.find("CellData." + name);
      VTK_ASSERT_MSG(data_it != dataArray_.end() && cell_it != cellData_.end(),
        "The data to extract is not found in cell-data. Try `getPointData()` instead!");
      VTK_ASSERT(data_it->second.section == CELL_DATA);
157

158
159
      return {*creator_, cell_it->second, data_it->second.name, data_it->second.components,
              data_it->second.type, vec_types, vec_offsets, vec_connectivity};
160
161
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
162
163
164
165
166
167
168
    /// Return a vector of DataArrayAttributes for all CELL_DATA blocks
    /// [[expects: read_ == true]]
    std::vector<DataArrayAttributes> getCellDataAttributes () const
    {
      std::vector<DataArrayAttributes> attributes;
      attributes.reserve(cellData_.size());
      for (auto const& da : dataArray_) {
169
        if (da.second.section == CELL_DATA)
Praetorius, Simon's avatar
Praetorius, Simon committed
170
171
172
173
174
          attributes.push_back(da.second);
      }
      return attributes;
    }

175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
    /// Advanced read methods
    /// @{

    /// Read the grid from an input stream, referring to a .vtu file, into the GridFactory \ref factory_
    /**
     * \param input   A STL input stream to read the VTK file from.
     * \param create  If `false`, only fill internal data structures, if `true`, also create the grid. [true]
     **/
    void readSerialFileFromStream (std::ifstream& input, bool create = true);

    /// Read the grid from and input stream, referring to a .pvtu file, into the GridFactory \ref factory_
    /**
     * \param input   A STL input stream to read the VTK file from.
     * \param create  If `false`, only fill internal data structures, if `true`, also create the grid. [true]
     **/
    void readParallelFileFromStream (std::ifstream& input, int rank, int size, bool create = true);

192
    /// Insert all internal data to the GridCreator
193
    /// NOTE: requires an aforegoing call to \ref read()
194
    void fillGridCreator (bool insertPieces = true);
Praetorius, Simon's avatar
Praetorius, Simon committed
195

196
197
    /// @}

198
199
200
201
202
203
    /// Return the filenames of parallel pieces
    std::vector<std::string> const& pieces () const
    {
      return pieces_;
    }

204
205
206
207
208
209
210
211
212
#ifndef DOXYGEN
    // Implementation of the FileReader interface
    static void fillFactoryImpl (GridFactory<Grid>& factory, std::string const& filename)
    {
      VtkReader reader{factory};
      reader.read(filename);
    }
#endif

Praetorius, Simon's avatar
Praetorius, Simon committed
213
  private:
214
215
    // Read values stored on the cells with ID `id`
    Sections readCellData (std::ifstream& input, std::string id);
216

217
218
    template <class F, class H>
    void readCellDataAppended (MetaType<F>, MetaType<H>, std::ifstream& input, std::string id);
219

220
221
    // Read values stored on the points with ID `id`
    Sections readPointData (std::ifstream& input, std::string id);
Praetorius, Simon's avatar
Praetorius, Simon committed
222

223
224
    template <class F, class H>
    void readPointDataAppended (MetaType<F>, MetaType<H>, std::ifstream& input, std::string id);
225

Praetorius, Simon's avatar
Praetorius, Simon committed
226

227
    // Read vertex coordinates from `input` stream and store in into `factory`
228
    Sections readPoints (std::ifstream& input, std::string id);
Praetorius, Simon's avatar
Praetorius, Simon committed
229

230
231
    template <class F, class H>
    void readPointsAppended (MetaType<F>, MetaType<H>, std::ifstream& input);
Praetorius, Simon's avatar
Praetorius, Simon committed
232

233

234
    // Read cell type, cell offsets and connectivity from `input` stream
235
    Sections readCells (std::ifstream& input, std::string id);
Praetorius, Simon's avatar
Praetorius, Simon committed
236

237
238
    template <class H>
    void readCellsAppended (MetaType<H>, std::ifstream& input);
239

240
    // Read data from appended section in vtk file, starting from `offset`
241
242
    template <class FloatType, class HeaderType>
    void readAppended (std::ifstream& input, std::vector<FloatType>& values, HeaderType offset);
Praetorius, Simon's avatar
Praetorius, Simon committed
243

244
245
246
247
    // Test whether line belongs to section
    bool isSection (std::string line,
                    std::string key,
                    Sections current,
248
                    Sections parent = NO_SECTION) const
Praetorius, Simon's avatar
Praetorius, Simon committed
249
250
251
252
253
254
255
    {
      bool result = line.substr(1, key.length()) == key;
      if (result && current != parent)
        DUNE_THROW(Exception , "<" << key << "> in wrong section." );
      return result;
    }

256
257
258
    // Convert a section into a string
    std::string toString (Sections s) const;

259
260
261
    // Find beginning of appended binary data
    std::uint64_t findAppendedDataPosition (std::ifstream& input) const;

262
    // Read attributes from current xml tag
Praetorius, Simon's avatar
Praetorius, Simon committed
263
    std::map<std::string, std::string> parseXml (std::string const& line, bool& closed);
Praetorius, Simon's avatar
Praetorius, Simon committed
264

Praetorius, Simon's avatar
Praetorius, Simon committed
265
266
    // clear all vectors
    void clear ();
Praetorius, Simon's avatar
Praetorius, Simon committed
267

268
269
270
271
272
    auto comm () const
    {
      return MPIHelper::getCollectiveCommunication();
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
273
  private:
274
    std::shared_ptr<GridCreator> creator_;
275
276

    /// Data format, i.e. ASCII, BINARY or COMPRESSED. Read from xml attributes.
Praetorius, Simon's avatar
Praetorius, Simon committed
277
278
    Vtk::FormatTypes format_;

279
280
281
    /// Type of compression algorithm used for binary data
    Vtk::CompressorTypes compressor_;

282
    // Temporary data to construct the grid elements
283
    std::vector<GlobalCoordinate> vec_points;
Praetorius, Simon's avatar
Praetorius, Simon committed
284
285
286
    std::vector<std::uint64_t> vec_point_ids;   //< Global unique vertex ID
    std::vector<std::uint8_t> vec_types;        //< VTK cell type ID
    std::vector<std::int64_t> vec_offsets;      //< offset of vertices of cell
287
    std::vector<std::int64_t> vec_connectivity; //< vertex indices of cell
Praetorius, Simon's avatar
Praetorius, Simon committed
288

Praetorius, Simon's avatar
Praetorius, Simon committed
289
290
    std::size_t numberOfCells_ = 0;   //< Number of cells in the grid
    std::size_t numberOfPoints_ = 0;  //< Number of vertices in the grid
Praetorius, Simon's avatar
Praetorius, Simon committed
291

292
    // offset information for appended data
293
    // map: id -> {Name,DataType,NumberOfComponents,Offset}
294
    std::map<std::string, DataArrayAttributes> dataArray_;
Praetorius, Simon's avatar
Praetorius, Simon committed
295

296
    // storage for internal point and cell data
297
    // map: id -> vector of data entries per point/cell
298
299
300
    std::map<std::string, std::vector<FieldType>> pointData_;
    std::map<std::string, std::vector<FieldType>> cellData_;

301
302
303
    // vector of filenames of parallel pieces
    std::vector<std::string> pieces_;

304
    /// Offset of beginning of appended data
Praetorius, Simon's avatar
Praetorius, Simon committed
305
    std::uint64_t offset0_ = 0;
Praetorius, Simon's avatar
Praetorius, Simon committed
306
307

    bool read_ = false;
Praetorius, Simon's avatar
Praetorius, Simon committed
308
309
  };

310
311
312
  // deduction guides
  template <class Grid>
  VtkReader (GridFactory<Grid>&)
Stenger, Florian's avatar
Stenger, Florian committed
313
    -> VtkReader<Grid, Vtk::ContinuousGridCreator<Grid>>;
314
315
316
317
318
319
320
321
322
323
324

  template <class GridCreator,
    class = std::void_t<typename GridCreator::Grid>>
  VtkReader (GridCreator&)
    -> VtkReader<typename GridCreator::Grid, GridCreator>;

  template <class GridCreator,
    class = std::void_t<typename GridCreator::Grid>>
  VtkReader (std::shared_ptr<GridCreator>)
    -> VtkReader<typename GridCreator::Grid, GridCreator>;

325
} // end namespace Dune
Praetorius, Simon's avatar
Praetorius, Simon committed
326

327
#include "vtkreader.impl.hh"