All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
vtkWritePolyDataTemplates.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2021-2023 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "vtkWritePolyData.H"
27 #include "OFstream.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 template<class Type, class DataType>
33 (
34  std::ostream& os,
35  const bool binary,
36  const wordList& fieldNames,
37  const boolList& fieldIsPointValues,
38  const UPtrList<const Field<Type>>& fieldTypeValues,
39  const bool writePointValues
40 )
41 {
42  forAll(fieldNames, fieldi)
43  {
44  if
45  (
46  fieldIsPointValues[fieldi] == writePointValues
47  && fieldTypeValues.set(fieldi)
48  )
49  {
50  const label nCmpt = pTraits<Type>::nComponents;
51 
52  os << fieldNames[fieldi] << ' ' << pTraits<Type>::nComponents
53  << ' ' << fieldTypeValues[fieldi].size() << ' '
54  << (std::is_integral<DataType>::value ? "int" : "float") << nl;
55 
56  List<DataType> data(nCmpt*fieldTypeValues[fieldi].size());
57  label i = 0;
58  forAll(fieldTypeValues[fieldi], fieldValuei)
59  {
60  for (direction cmpt = 0; cmpt < nCmpt; ++ cmpt)
61  {
62  data[i ++] =
64  (
65  component
66  (
67  fieldTypeValues[fieldi][fieldValuei],
68  cmpt
69  )
70  );
71  }
72  }
73 
74  vtkWriteOps::write(os, binary, data);
75  }
76  }
77 }
78 
79 
80 template<class PointField, class VertexList, class LineList, class FaceList>
82 (
83  const fileName& file,
84  const word& title,
85  const bool binary,
86  const PointField& points,
87  const VertexList& vertices,
88  const LineList& lines,
89  const FaceList& faces,
90  const wordList& fieldNames,
91  const boolList& fieldIsPointValues,
92  const UPtrList<const Field<label>>& fieldLabelValues
93  #define FieldTypeValuesConstArg(Type, nullArg) \
94  , const UPtrList<const Field<Type>>& field##Type##Values
97 )
98 {
99  // Open the file
100  std::ofstream os(file, std::ios::binary);
101 
102  // Write the header
103  vtkWriteOps::writeHeader(os, binary, title);
104  os << "DATASET POLYDATA" << nl;
105 
106  // Write the points
107  {
108  os << "POINTS " << points.size() << " float" << nl;
109  List<floatScalar> coordinates(points.size()*3);
110  forAll(points, pointi)
111  {
112  const point& p = points[pointi];
113  forAll(p, i)
114  {
115  coordinates[3*pointi + i] = vtkWriteOps::cast(p[i]);
116  }
117  }
118  vtkWriteOps::write(os, binary, coordinates);
119  }
120 
121  // Write the vertices
122  if (vertices.size())
123  {
124  os << "VERTICES " << vertices.size() << ' '
125  << 2*vertices.size() << nl;
126  labelList data(2*vertices.size());
127  forAll(vertices, vertexi)
128  {
129  data[2*vertexi] = 1;
130  data[2*vertexi + 1] = vertices[vertexi];
131  }
132  vtkWriteOps::write(os, binary, data);
133  }
134 
135  // Write the lines
136  if (lines.size())
137  {
138  label nLineNodes = 0;
139  forAll(lines, facei)
140  {
141  nLineNodes += lines[facei].size();
142  }
143  os << "LINES " << lines.size() << ' '
144  << lines.size() + nLineNodes << nl;
145  labelList data(lines.size() + nLineNodes);
146  label i = 0;
147  forAll(lines, linei)
148  {
149  data[i ++] = lines[linei].size();
150  forAll(lines[linei], linePointi)
151  {
152  data[i ++] = lines[linei][linePointi];
153  }
154  }
155  vtkWriteOps::write(os, binary, data);
156  }
157 
158  // Write the faces
159  if (faces.size())
160  {
161  label nFaceNodes = 0;
162  forAll(faces, facei)
163  {
164  nFaceNodes += faces[facei].size();
165  }
166  os << "POLYGONS " << faces.size() << ' '
167  << faces.size() + nFaceNodes << nl;
168  labelList data(faces.size() + nFaceNodes);
169  label i = 0;
170  forAll(faces, facei)
171  {
172  data[i ++] = faces[facei].size();
173  forAll(faces[facei], facePointi)
174  {
175  data[i ++] = faces[facei][facePointi];
176  }
177  }
178  vtkWriteOps::write(os, binary, data);
179  }
180 
181  // Write the fields
182  const label nPointFields = count(fieldIsPointValues, true);
183  const label nFaceFields = count(fieldIsPointValues, false);
184  if (nPointFields > 0)
185  {
186  os << "POINT_DATA " << points.size() << nl
187  << "FIELD attributes " << nPointFields << nl;
188  writeFieldTypeValues<label, label>
189  (
190  os,
191  binary,
192  fieldNames,
193  fieldIsPointValues,
194  fieldLabelValues,
195  true
196  );
197  #define WriteFieldTypeValues(Type, nullArg) \
198  writeFieldTypeValues<Type, floatScalar> \
199  ( \
200  os, \
201  binary, \
202  fieldNames, \
203  fieldIsPointValues, \
204  field##Type##Values, \
205  true \
206  );
208  #undef WriteFieldTypeValues
209  }
210  if (nFaceFields > 0)
211  {
212  os << "CELL_DATA "
213  << vertices.size() + lines.size() + faces.size() << nl
214  << "FIELD attributes " << nFaceFields << nl;
215  writeFieldTypeValues<label, label>
216  (
217  os,
218  binary,
219  fieldNames,
220  fieldIsPointValues,
221  fieldLabelValues,
222  false
223  );
224  #define WriteFieldTypeValues(Type, nullArg) \
225  writeFieldTypeValues<Type, floatScalar> \
226  ( \
227  os, \
228  binary, \
229  fieldNames, \
230  fieldIsPointValues, \
231  field##Type##Values, \
232  false \
233  );
235  #undef WriteFieldTypeValues
236  }
237 }
238 
239 
240 template
241 <
242  class PointField,
243  class VertexList,
244  class LineList,
245  class FaceList,
246  class ... Args
247 >
249 (
250  const fileName& file,
251  const word& title,
252  const bool binary,
253  const PointField& points,
254  const VertexList& vertices,
255  const LineList& lines,
256  const FaceList& faces,
257  const Args& ... args
258 )
259 {
260  const label nFields = sizeof...(Args)/3;
261 
262  wordList fieldNames(nFields);
263  boolList fieldIsPointValues(nFields);
264  UPtrList<const Field<label>> fieldLabelValues(nFields);
265  #define DeclareFieldTypeValues(Type, nullArg) \
266  UPtrList<const Field<Type>> field##Type##Values(nFields);
268  #undef DeclareFieldTypeValues
269 
271  (
272  fieldNames,
273  fieldIsPointValues,
274  fieldLabelValues
275  #define FieldTypeValuesParameter(Type, nullArg) , field##Type##Values
278  args ...
279  );
280 
281  write
282  (
283  file,
284  title,
285  binary,
286  points,
287  vertices,
288  lines,
289  faces,
290  fieldNames,
291  fieldIsPointValues,
292  fieldLabelValues
293  #define FieldTypeValuesParameter(Type, nullArg) , field##Type##Values
296  );
297 }
298 
299 
300 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static List< word > fieldNames
Definition: globalFoam.H:46
const pointField & points
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
void writeHeader(std::ostream &, const bool isBinary, const std::string &title)
Write header.
label cast(const label &x)
Cast an integer to the corresponding VTK write type. Does nothing.
void write(const fileName &file, const word &title, const bool binary, const PointField &points, const VertexList &vertices, const LineList &lines, const FaceList &faces, const wordList &fieldNames, const boolList &fieldIsPointValues, const UPtrList< const Field< label >> &fieldLabelValues #define FieldTypeValuesConstArg(Type, nullArg))
Write VTK polygonal data to a file. Takes a PtrList of fields of labels and.
void writeFieldTypeValues(std::ostream &os, const bool binary, const wordList &fieldNames, const boolList &fieldIsPointValues, const UPtrList< const Field< Type >> &fieldTypeValues, const bool writePointValues)
Write the field values out for a type.
void unpackFieldTypeValues(wordList &fieldNames, boolList &fieldIsPointValues, UPtrList< const Field< label >> &fieldLabelValues #define FieldTypeValuesNonConstArg(Type, nullArg))
Helper for templated write.
List< word > wordList
A List of words.
Definition: fileName.H:54
pointField vertices(const blockVertexList &bvl)
List< label > labelList
A List of labels.
Definition: labelList.H:56
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
vector point
Point is a vector.
Definition: point.H:41
GeometricField< Type, pointPatchField, pointMesh > PointField
FOR_ALL_FIELD_TYPES(makeFieldSourceTypedef)
label count(const ListType &l, typename ListType::const_reference x)
Count the number of occurrences of a value in a list.
static const char nl
Definition: Ostream.H:266
uint8_t direction
Definition: direction.H:45
Foam::argList args(argc, argv)
volScalarField & p
#define FieldTypeValuesConstArg(Type, nullArg)
#define FieldTypeValuesParameter(Type, nullArg)
#define WriteFieldTypeValues(Type, nullArg)
#define DeclareFieldTypeValues(Type, nullArg)