foamVtkTools.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017-2020 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Namespace
27 Foam::vtk::Tools
28
29Description
30 A collection of static methods to assist converting OpenFOAM data
31 structures into VTK internal data structures.
32
33 Remapping of the symmTensor order is required in input or output
34 directions. OpenFOAM uses (XX, XY, XZ, YY, YZ, ZZ) order,
35 VTK uses (XX, YY, ZZ, XY, YZ, XZ) order.
36
37Note
38 The class is implemented as headers-only.
39
40SourceFiles
41 foamVtkToolsI.H
42 foamVtkToolsTemplates.C
43
44\*---------------------------------------------------------------------------*/
45
46#ifndef Foam_vtk_vtkTools_H
47#define Foam_vtk_vtkTools_H
48
49#include "faceList.H"
50#include "pointField.H"
51#include "symmTensor.H"
52#include "MinMax.H"
53
54// VTK includes
55#include "vtkCellArray.h"
56#include "vtkFloatArray.h"
57#include "vtkDoubleArray.h"
58#include "vtkIdTypeArray.h"
59#include "vtkSmartPointer.h"
60#include "vtkUnsignedCharArray.h"
61#include "vtkPoints.h"
62#include "vtkPolyData.h"
63
64#include <utility>
65
66// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
67
68// Forward Declarations
69class vtkDataSet;
70class vtkCellData;
71class vtkPointData;
72
73namespace Foam
74{
75namespace vtk
76{
77
78/*---------------------------------------------------------------------------*\
79 Class vtk::Caching Declaration
80\*---------------------------------------------------------------------------*/
81
82//- Bookkeeping for internal caching.
83// Retain an original copy of the geometry as well as a shallow copy
84// with the output fields.
85// The original copy is reused for different timesteps etc.
86template<class DataType>
87struct Caching
88{
89 typedef DataType dataType;
90
91 //- The geometry, without any cell/point data
92 vtkSmartPointer<dataType> vtkgeom;
93
94 //- The shallow-copy of geometry, plus additional data
95 vtkSmartPointer<dataType> dataset;
96
97 //- Number of points associated with the geometry
98 inline uint64_t nPoints() const
99 {
100 return vtkgeom ? vtkgeom->GetNumberOfPoints() : 0;
101 }
102
103 //- Clear geometry and dataset
105 {
106 vtkgeom = nullptr;
107 dataset = nullptr;
108 }
109
110 //- Return a shallow copy of vtkgeom for manipulation
111 vtkSmartPointer<dataType> getCopy() const
112 {
113 auto copy = vtkSmartPointer<dataType>::New();
114 if (vtkgeom)
115 {
116 copy->ShallowCopy(vtkgeom);
117 }
118 return copy;
119 }
120
121 //- Make a shallow copy of vtkgeom into dataset
122 void reuse()
123 {
124 dataset = vtkSmartPointer<dataType>::New();
125 if (vtkgeom)
126 {
127 dataset->ShallowCopy(vtkgeom);
128 }
129 }
130
131 //- Set the geometry and make a shallow copy to dataset
132 void set(vtkSmartPointer<dataType> geom)
133 {
134 vtkgeom = geom;
135 reuse();
136 }
137
138 //- Report basic information to output
139 void PrintSelf(std::ostream& os) const
140 {
141 os << "geom" << nl;
142 if (vtkgeom)
143 {
144 vtkgeom->PrintSelf(std::cout, vtkIndent(2));
145 }
146 else
147 {
148 os << "nullptr";
149 }
150 os << nl;
151
152 os << "copy" << nl;
153 if (dataset)
154 {
155 dataset->PrintSelf(std::cout, vtkIndent(2));
156 }
157 else
158 {
159 os << "nullptr";
160 }
161 os << nl;
162 }
163};
164
165
166/*---------------------------------------------------------------------------*\
167 Namespace vtk::Tools
168\*---------------------------------------------------------------------------*/
169
170namespace Tools
171{
172 //- Wrap vtkUnsignedCharArray as a UList
174 (
175 vtkUnsignedCharArray* array,
176 const label size
177 );
178
179 //- Wrap vtkIdTypeArray as a UList
181 (
182 vtkIdTypeArray* array,
183 const label size
184 );
185
186 //- Return a list of points as vtkPoints
187 inline vtkSmartPointer<vtkPoints> Points
188 (
189 const UList<point>& pts
190 );
191
192 //- Return an indirect list of points as vtkPoints
193 inline vtkSmartPointer<vtkPoints> Points
194 (
195 const UList<point>& pts,
196 const labelUList& addr
197 );
198
199 //- Convert a list of faces (or triFaces) to vtk polygon cells
200 template<class Face>
201 vtkSmartPointer<vtkCellArray> Faces(const UList<Face>& faces);
202
203 //- Return vtkPolyData of vertices for each point
204 inline vtkSmartPointer<vtkPolyData> Vertices
205 (
206 const UList<point>& pts
207 );
208
209 //- Return vtkPolyData of vertices for each point
210 inline vtkSmartPointer<vtkPolyData> Vertices
211 (
212 const UList<point>& pts,
213 const labelUList& addr
214 );
215
216 //- Min/Max of scalar, or mag() of non-scalars. Includes nullptr check.
217 inline scalarMinMax rangeOf(vtkDataArray* data);
218
219
220 //- Convert OpenFOAM patch to vtkPolyData
221 struct Patch
222 {
223 //- Return local patch points as vtkPoints
224 template<class PatchType>
225 static vtkSmartPointer<vtkPoints> points(const PatchType& p);
226
227 //- Convert patch faces to vtk polygon cells
228 template<class PatchType>
229 static vtkSmartPointer<vtkCellArray> faces(const PatchType& p);
230
231 //- Convert patch points/faces to vtkPolyData
232 template<class PatchType>
233 static vtkSmartPointer<vtkPolyData> mesh(const PatchType& p);
234
235 //- Convert patch face normals to vtkFloatArray
236 template<class PatchType>
237 static vtkSmartPointer<vtkFloatArray> faceNormals(const PatchType& p);
238
239 //- Return patch face centres as vtkPoints
240 template<class PatchType>
241 static vtkSmartPointer<vtkPoints> faceCentres(const PatchType& p);
242
243 //- Convert points/faces component to vtkPolyData
244 template<class Face>
245 static vtkSmartPointer<vtkPolyData> mesh
246 (
247 const UList<point>& pts,
248 const UList<Face>& fcs
249 );
250 };
251
252
253 //- Remapping for some OpenFOAM data types (eg, symmTensor)
254 // \param data[in,out] The data to be remapped in-place
255 template<class Type>
256 inline void remapTuple(float data[]) {}
257
258 //- Template specialization for symmTensor ordering
259 template<>
260 inline void remapTuple<symmTensor>(float data[]);
261
262 //- Remapping for some OpenFOAM data types (eg, symmTensor)
263 // \param data[in,out] The data to be remapped in-place
264 template<class Type>
265 inline void remapTuple(double data[]) {}
266
267 //- Template specialization for symmTensor ordering
268 template<>
269 inline void remapTuple<symmTensor>(double data[]);
270
271 //- Copy/transcribe OpenFOAM data types to VTK format
272 // This allows a change of data type (float vs double) as well as
273 // addressing any swapping issues (eg, symmTensor)
274 //
275 // \param output[out] The output scratch space. Must be long enough
276 // to hold the result.
277 // \param val[in] The input data to copy/transcribe
278 template<class Type>
279 inline void foamToVtkTuple(float output[], const Type& val);
280
281 //- Copy/transcribe OpenFOAM data types to VTK format
282 // This allows a change of data type (float vs double) as well as
283 // addressing any swapping issues (eg, symmTensor)
284 //
285 // \param output[out] The output scratch space. Must be long enough
286 // to hold the result.
287 // \param val[in] The input data to copy/transcribe
288 template<class Type>
289 inline void foamToVtkTuple(double output[], const Type& val);
290
291
292 // Field Conversion Functions
293
294 //- Copy list to pre-allocated vtk array.
295 // \return number of input items copied
296 template<class Type>
298 (
299 vtkFloatArray* array,
300 const UList<Type>& input,
301 vtkIdType start = 0
302 );
303
304 //- Create named field initialized to zero
305 template<class Type>
306 vtkSmartPointer<vtkFloatArray> zeroField
307 (
308 const word& name,
309 const label size
310 );
311
312 //- Convert field data to a vtkFloatArray
313 template<class Type>
314 vtkSmartPointer<vtkFloatArray> convertFieldToVTK
315 (
316 const word& name,
317 const UList<Type>& fld
318 );
319
320 //- An identity list of VTK_VERTEX
321 inline vtkSmartPointer<vtkCellArray> identityVertices
322 (
323 const label size
324 );
325
326} // End namespace Tools
327
328// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
329
330} // End namespace vtk
331} // End namespace Foam
332
333
334// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
335
336// Specializations
337
338//- Template specialization for symmTensor ordering
339template<>
340void Foam::vtk::Tools::remapTuple<Foam::symmTensor>(float data[])
341{
342 std::swap(data[1], data[3]); // swap XY <-> YY
343 std::swap(data[2], data[5]); // swap XZ <-> ZZ
344}
345
346
347//- Template specialization for symmTensor ordering
348template<>
349void Foam::vtk::Tools::remapTuple<Foam::symmTensor>(double data[])
350{
351 std::swap(data[1], data[3]); // swap XY <-> YY
352 std::swap(data[2], data[5]); // swap XZ <-> ZZ
353}
354
355
356// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
357
358#include "foamVtkToolsI.H"
359
360// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
361
362#ifdef NoRepository
363 #include "foamVtkToolsTemplates.C"
364#endif
365
366// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
367
368#endif
369
370// ************************************************************************* //
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: UList.H:94
Database for solution data, solver performance and other reduced data.
Definition: data.H:58
A class for handling words, derived from Foam::string.
Definition: word.H:68
A class representing the concept of a field of 0 used to avoid unnecessary manipulations for objects ...
Definition: zeroField.H:56
volScalarField & p
dynamicFvMesh & mesh
OBJstream os(runTime.globalPath()/outputName)
const pointField & points
void foamToVtkTuple(float output[], const Type &val)
Copy/transcribe OpenFOAM data types to VTK format.
UList< uint8_t > asUList(vtkUnsignedCharArray *array, const label size)
Wrap vtkUnsignedCharArray as a UList.
Definition: foamVtkToolsI.H:33
scalarMinMax rangeOf(vtkDataArray *data)
Min/Max of scalar, or mag() of non-scalars. Includes nullptr check.
vtkSmartPointer< vtkPoints > Points(const UList< point > &pts)
Return a list of points as vtkPoints.
Definition: foamVtkToolsI.H:59
label transcribeFloatData(vtkFloatArray *array, const UList< Type > &input, vtkIdType start=0)
Copy list to pre-allocated vtk array.
vtkSmartPointer< vtkFloatArray > convertFieldToVTK(const word &name, const UList< Type > &fld)
Convert field data to a vtkFloatArray.
void remapTuple< symmTensor >(float data[])
Template specialization for symmTensor ordering.
void remapTuple(float data[])
Remapping for some OpenFOAM data types (eg, symmTensor)
Definition: foamVtkTools.H:256
vtkSmartPointer< vtkCellArray > identityVertices(const label size)
An identity list of VTK_VERTEX.
vtkSmartPointer< vtkCellArray > Faces(const UList< Face > &faces)
Convert a list of faces (or triFaces) to vtk polygon cells.
vtkSmartPointer< vtkPolyData > Vertices(const UList< point > &pts)
Return vtkPolyData of vertices for each point.
Definition: foamVtkToolsI.H:93
Namespace for OpenFOAM.
static Istream & input(Istream &is, IntRange< T > &range)
Definition: IntRanges.C:55
static Ostream & output(Ostream &os, const IntRange< T > &range)
Definition: IntRanges.C:66
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Bookkeeping for internal caching.
Definition: foamVtkTools.H:88
void set(vtkSmartPointer< dataType > geom)
Set the geometry and make a shallow copy to dataset.
Definition: foamVtkTools.H:132
void PrintSelf(std::ostream &os) const
Report basic information to output.
Definition: foamVtkTools.H:139
void clearGeom()
Clear geometry and dataset.
Definition: foamVtkTools.H:104
vtkSmartPointer< dataType > vtkgeom
The geometry, without any cell/point data.
Definition: foamVtkTools.H:92
void reuse()
Make a shallow copy of vtkgeom into dataset.
Definition: foamVtkTools.H:122
vtkSmartPointer< dataType > dataset
The shallow-copy of geometry, plus additional data.
Definition: foamVtkTools.H:95
uint64_t nPoints() const
Number of points associated with the geometry.
Definition: foamVtkTools.H:98
vtkSmartPointer< dataType > getCopy() const
Return a shallow copy of vtkgeom for manipulation.
Definition: foamVtkTools.H:111
Convert OpenFOAM patch to vtkPolyData.
Definition: foamVtkTools.H:222
static vtkSmartPointer< vtkCellArray > faces(const PatchType &p)
Convert patch faces to vtk polygon cells.
static vtkSmartPointer< vtkFloatArray > faceNormals(const PatchType &p)
Convert patch face normals to vtkFloatArray.
static vtkSmartPointer< vtkPoints > faceCentres(const PatchType &p)
Return patch face centres as vtkPoints.