foamVtkToolsI.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
26\*---------------------------------------------------------------------------*/
27
28#include <numeric>
29
30// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31
33(
34 vtkUnsignedCharArray* array,
35 const label size
36)
37{
38 array->SetNumberOfComponents(1);
39 array->SetNumberOfTuples(size);
40
41 return UList<uint8_t>(array->WritePointer(0, size), size);
42}
43
44
46(
47 vtkIdTypeArray* array,
48 const label size
49)
50{
51 array->SetNumberOfComponents(1);
52 array->SetNumberOfTuples(size);
53
54 return UList<vtkIdType>(array->WritePointer(0, size), size);
55}
56
57
58inline vtkSmartPointer<vtkPoints>
60{
61 auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
62
63 vtkpoints->SetNumberOfPoints(pts.size());
64
65 vtkIdType pointId = 0;
66 for (const point& pt : pts)
67 {
68 vtkpoints->SetPoint(pointId++, pt.v_);
69 }
70
71 return vtkpoints;
72}
73
74
75inline vtkSmartPointer<vtkPoints>
77{
78 auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
79
80 vtkpoints->SetNumberOfPoints(addr.size());
81
82 vtkIdType pointId = 0;
83 for (const label pointi : addr)
84 {
85 vtkpoints->SetPoint(pointId++, pts[pointi].v_);
86 }
87
88 return vtkpoints;
89}
90
91
92inline vtkSmartPointer<vtkPolyData>
94{
95 auto vtkmesh = vtkSmartPointer<vtkPolyData>::New();
96
97 vtkmesh->SetPoints(Tools::Points(pts));
98 vtkmesh->SetVerts(Tools::identityVertices(pts.size()));
99
100 return vtkmesh;
101}
102
103
104inline vtkSmartPointer<vtkPolyData>
106{
107 auto vtkmesh = vtkSmartPointer<vtkPolyData>::New();
108
109 vtkmesh->SetPoints(Tools::Points(pts, addr));
110 vtkmesh->SetVerts(Tools::identityVertices(addr.size()));
111
112 return vtkmesh;
113}
114
115
117{
118 double range[2]{GREAT, -GREAT};
119
120 if (data)
121 {
122 if (data->GetNumberOfComponents() == 1)
123 {
124 data->GetRange(range);
125 }
126 else
127 {
128 // Mag
129 data->GetRange(range, -1);
130 }
131 }
132
133 return scalarMinMax(range[0], range[1]);
134}
135
136
137inline vtkSmartPointer<vtkCellArray> Foam::vtk::Tools::identityVertices
138(
139 const label size
140)
141{
142 auto cells = vtkSmartPointer<vtkCellArray>::New();
143
144 #ifdef VTK_CELL_ARRAY_V2
145
146 // Offsets
147 // [0, n1, n1+n2, n1+n2+n3... ]
148
149 auto offsets = vtkSmartPointer<vtkIdTypeArray>::New();
150
151 {
152 const vtkIdType nOffsets(size+1);
153
154 offsets->SetNumberOfTuples(nOffsets);
155
156 vtkIdType* iter = offsets->WritePointer(0, nOffsets);
157
158 std::iota(iter, (iter + nOffsets), vtkIdType(0));
159 }
160
161
162 auto connect = vtkSmartPointer<vtkIdTypeArray>::New();
163
164 // Connectivity
165 {
166 const vtkIdType nConnect(size);
167
168 connect->SetNumberOfTuples(nConnect);
169
170 vtkIdType* iter = connect->WritePointer(0, nConnect);
171
172 std::iota(iter, (iter + nConnect), vtkIdType(0));
173 }
174
175 // Move into a vtkCellArray
176
177 cells->SetData(offsets, connect);
178
179 #else
180
181 // In VTK-8.2.0 and older,
182 // sizes are interwoven (prefixed) in the connectivity
183
184 // Connectivity size, with prefixed size information
185
186 // Cell connectivity for vertex
187 // [size, ids.., size, ids...] -> therefore [1, id, 1, id, ...]
188
189 const vtkIdType nElem(size);
190 const vtkIdType nConnect(2*size);
191
192 {
193 cells->GetData()->SetNumberOfTuples(nConnect);
194
195 vtkIdType* iter = cells->WritePointer(nElem, nConnect);
196
197 // Fill in the connectivity array, with prefixed size information
198
199 for (vtkIdType id = 0; id < nElem; ++id)
200 {
201 *(iter++) = 1;
202 *(iter++) = id;
203 }
204 }
205
206 #endif
207
208 return cells;
209};
210
211
212template<class Type>
214(
215 float output[],
216 const Type& val
217)
218{
219 for (direction cmpt=0; cmpt < pTraits<Type>::nComponents; ++cmpt)
220 {
221 output[cmpt] = component(val, cmpt);
222 }
223 remapTuple<Type>(output);
224}
225
226
227template<class Type>
229(
230 double output[],
231 const Type& val
232)
233{
234 for (direction cmpt=0; cmpt < pTraits<Type>::nComponents; ++cmpt)
235 {
236 output[cmpt] = component(val, cmpt);
237 }
238 remapTuple<Type>(output);
239}
240
241
242// ************************************************************************* //
scalar range
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
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Database for solution data, solver performance and other reduced data.
Definition: data.H:58
const cellShapeList & cells
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
vtkSmartPointer< vtkCellArray > identityVertices(const label size)
An identity list of VTK_VERTEX.
vtkSmartPointer< vtkPolyData > Vertices(const UList< point > &pts)
Return vtkPolyData of vertices for each point.
Definition: foamVtkToolsI.H:93
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
uint8_t direction
Definition: direction.H:56