24#ifndef DUMUX_POINTCLOUD_VTK_WRITER_HH
25#define DUMUX_POINTCLOUD_VTK_WRITER_HH
30#include <dune/common/fvector.hh>
31#include <dune/common/exceptions.hh>
32#include <dune/common/path.hh>
33#include <dune/grid/io/file/vtk/common.hh>
46template<
class Scalar,
class GlobalPosition>
51 using DimWorldVector = Dune::FieldVector<Scalar, GlobalPosition::size()>;
53 static constexpr unsigned int precision = 6;
54 static constexpr unsigned int numBeforeLineBreak = 15;
59 template<
class ContainerType>
70 VTKFunction(
const ContainerType& data,
const std::string& name,
const int numComponents) : data_(data), name_(name), numComponents_(numComponents)
76 const std::string& name()
const
84 int numComponents()
const
86 return numComponents_;
94 auto& operator() (
const int idx)
const {
return data_[idx]; }
96 decltype(
auto) begin()
const
101 decltype(
auto) end()
const
115 const ContainerType& data_;
116 const std::string name_;
117 const int numComponents_;
135 void write(
const std::string& name, Dune::VTK::OutputType type = Dune::VTK::ascii)
142 writeCoordinates_(file, coordinates_);
143 writeDataInfo_(file);
145 for (
auto&& data : scalarPointData_)
146 writeData_(file, data);
148 for (
auto&& data :vectorPointData_)
149 writeData_(file, data);
151 if (!scalarPointData_.empty() || !vectorPointData_.empty())
152 file <<
"</PointData>\n";
154 file <<
"</Piece>\n";
155 file <<
"</PolyData>\n";
156 file <<
"</VTKFile>";
177 void pwrite(
const std::string & name,
const std::string & path,
const std::string & extendpath,
178 Dune::VTK::OutputType type = Dune::VTK::ascii)
180 DUNE_THROW(Dune::NotImplemented,
"parallel point cloud vtk output not supported yet");
189 void addPointData(
const std::vector<Scalar>& v,
const std::string &name)
191 assert(v.size() == coordinates_.size());
201 void addPointData(
const std::vector<DimWorldVector>& v,
const std::string &name)
203 assert(v.size() == coordinates_.size());
212 scalarPointData_.clear();
213 vectorPointData_.clear();
228 const std::string& path)
const
230 static const std::string extension =
".vtp";
232 return Dune::concatPaths(path, name+extension);
248 const std::string& path,
251 std::ostringstream s;
252 if(path.size() > 0) {
254 if(path[path.size()-1] !=
'/')
257 s <<
's' << std::setw(4) << std::setfill(
'0') << commSize <<
'-';
268 void writeHeader_(std::ostream& file)
270 std::string header =
"<?xml version=\"1.0\"?>\n";
271 header +=
"<VTKFile type=\"PolyData\" version=\"0.1\" byte_order=\"LittleEndian\">\n";
272 header +=
"<PolyData>\n";
273 header +=
"<Piece NumberOfLines=\"0\" NumberOfPoints=\"" + std::to_string(coordinates_.size()) +
"\">\n";
280 void writeDataInfo_(std::ostream& file)
282 std::string scalarName;
283 std::string vectorName;
284 bool foundScalar =
false;
285 bool foundVector =
false;
287 for(
auto&& data : scalarPointData_)
289 if(data.numComponents() == 1 && !foundScalar)
291 scalarName = data.name();
296 if(data.numComponents() > 1 && !foundVector)
298 vectorName = data.name();
303 for(
auto&& data : vectorPointData_)
305 if(data.numComponents() > 1 && !foundVector)
307 vectorName = data.name();
314 file <<
"<PointData Scalars=\"" << scalarName <<
"\" Vectors=\"" << vectorName <<
"\">\n";
316 file <<
"<PointData Scalars=\"" << scalarName <<
"\">\n";
318 file <<
"<PointData Vectors=\"" << vectorName <<
"\">\n";
329 void writeCoordinates_(std::ostream& file,
const std::vector<GlobalPosition>& positions)
332 file <<
"<Points>\n";
333 file <<
"<DataArray type=\"Float32\" Name=\"Coordinates\" NumberOfComponents=\"3\" format=\"ascii\">\n";
335 for(
auto&& x : positions)
347 if((++counter) > numBeforeLineBreak)
353 file <<
"\n</DataArray>\n";
354 file <<
"</Points>\n";
364 void writeData_(std::ostream& file,
const T& data)
366 file <<
"<DataArray type=\"Float32\" Name=\"" << data.name() <<
"\" NumberOfComponents=\"" << data.numComponents() <<
"\" format=\"ascii\">\n";
368 for(
auto&& value : data)
371 writeToFile_(file, value);
374 if((++counter) > numBeforeLineBreak)
380 file <<
"\n</DataArray>\n";
389 void writeToFile_(std::ostream& file,
const Scalar& s)
400 void writeToFile_(std::ostream& file,
const DimWorldVector& g)
402 assert(g.size() > 1 && g.size() < 4);
409 const std::vector<GlobalPosition>& coordinates_;
410 std::list<ScalarFunction> scalarPointData_;
411 std::list<VectorFunction> vectorPointData_;
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
A VTK output module to simplify writing dumux simulation data to VTK format.
Definition: pointcloudvtkwriter.hh:48
std::string getSerialPieceName(const std::string &name, const std::string &path) const
Return name of a serial header file.
Definition: pointcloudvtkwriter.hh:227
VTKFunction< std::vector< Scalar > > ScalarFunction
Definition: pointcloudvtkwriter.hh:122
std::string getParallelHeaderName(const std::string &name, const std::string &path, int commSize) const
Return name of a parallel header file.
Definition: pointcloudvtkwriter.hh:247
void addPointData(const std::vector< DimWorldVector > &v, const std::string &name)
Add a vector of vector data that live on arbitrary points to the visualization.
Definition: pointcloudvtkwriter.hh:201
void addPointData(const std::vector< Scalar > &v, const std::string &name)
Add a vector of scalar data that live on arbitrary points to the visualization.
Definition: pointcloudvtkwriter.hh:189
VTKFunction< std::vector< DimWorldVector > > VectorFunction
Definition: pointcloudvtkwriter.hh:123
PointCloudVtkWriter(const std::vector< GlobalPosition > &coordinates)
Definition: pointcloudvtkwriter.hh:126
void write(const std::string &name, Dune::VTK::OutputType type=Dune::VTK::ascii)
Create an output file.
Definition: pointcloudvtkwriter.hh:135
void clear()
Clears the data lists.
Definition: pointcloudvtkwriter.hh:210
void pwrite(const std::string &name, const std::string &path, const std::string &extendpath, Dune::VTK::OutputType type=Dune::VTK::ascii)
Create an output file in parallel.
Definition: pointcloudvtkwriter.hh:177