foamVtkTools.H 8.88 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2017-2018 OpenCFD Ltd.
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

24
Namespace
25 26 27 28 29 30
    Foam::vtk::Tools

Description
    A collection of static methods to assist converting OpenFOAM data
    structures into VTK internal data structures.

31 32 33 34
    Remapping of the symmTensor order is required in input or output
    directions. OpenFOAM uses (XX, XY, XZ, YY, YZ, ZZ) order,
    VTK uses (XX, YY, ZZ, XY, YZ, XZ) order.

35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51
Note
    The class is implemented as headers-only.

SourceFiles
    foamVtkToolsI.H
    foamVtkToolsTemplates.C

\*---------------------------------------------------------------------------*/

#ifndef foamVtkTools_H
#define foamVtkTools_H

#include "labelList.H"
#include "faceList.H"
#include "pointField.H"
#include "symmTensor.H"

52
// VTK includes
53 54 55 56 57 58 59 60
#include "vtkCellArray.h"
#include "vtkFloatArray.h"
#include "vtkDoubleArray.h"
#include "vtkIdTypeArray.h"
#include "vtkSmartPointer.h"
#include "vtkUnsignedCharArray.h"
#include "vtkPoints.h"
#include "vtkPolyData.h"
61 62 63 64 65

#include <utility>

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

66 67 68 69 70 71 72
// Forward declarations
class vtkDataSet;
class vtkCellData;
class vtkPointData;

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166
namespace Foam
{
namespace vtk
{

/*---------------------------------------------------------------------------*\
                        Class vtk::Caching Declaration
\*---------------------------------------------------------------------------*/

//- Bookkeeping for internal caching.
//  Retain an original copy of the geometry as well as a shallow copy
//  with the output fields.
//  The original copy is reused for different timesteps etc.
template<class DataType>
struct Caching
{
    typedef DataType dataType;

    //- The geometry, without any cell/point data
    vtkSmartPointer<dataType> vtkgeom;

    //- The shallow-copy of geometry, plus additional data
    vtkSmartPointer<dataType> dataset;

    //- Number of points associated with the geometry
    inline uint64_t nPoints() const
    {
        return vtkgeom ? vtkgeom->GetNumberOfPoints() : 0;
    }

    //- Clear geometry and dataset
    void clearGeom()
    {
        vtkgeom = nullptr;
        dataset = nullptr;
    }

    //- Return a shallow copy of vtkgeom for manipulation
    vtkSmartPointer<dataType> getCopy() const
    {
        auto copy = vtkSmartPointer<dataType>::New();
        if (vtkgeom)
        {
            copy->ShallowCopy(vtkgeom);
        }
        return copy;
    }

    //- Make a shallow copy of vtkgeom into dataset
    void reuse()
    {
        dataset = vtkSmartPointer<dataType>::New();
        if (vtkgeom)
        {
            dataset->ShallowCopy(vtkgeom);
        }
    }

    //- Set the geometry and make a shallow copy to dataset
    void set(vtkSmartPointer<dataType> geom)
    {
        vtkgeom = geom;
        reuse();
    }

    //- Report basic information to output
    void PrintSelf(std::ostream& os) const
    {
        os << "geom" << nl;
        if (vtkgeom)
        {
            vtkgeom->PrintSelf(std::cout, vtkIndent(2));
        }
        else
        {
            os << "nullptr";
        }
        os << nl;

        os << "copy" << nl;
        if (dataset)
        {
            dataset->PrintSelf(std::cout, vtkIndent(2));
        }
        else
        {
            os << "nullptr";
        }
        os << nl;
    }
};


/*---------------------------------------------------------------------------*\
167
                             Namespace vtk::Tools
168 169
\*---------------------------------------------------------------------------*/

170
namespace Tools
171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208
{
    //- Wrap vtkUnsignedCharArray as a UList
    inline static UList<uint8_t> asUList
    (
        vtkUnsignedCharArray* array,
        const label size
    );

    //- Wrap vtkIdTypeArray as a UList
    inline static UList<vtkIdType> asUList
    (
        vtkIdTypeArray* array,
        const label size
    );

    //- Wrap vtkCellArray as a UList
    inline static UList<vtkIdType> asUList
    (
        vtkCellArray* cells,
        const label nCells,
        const label size
    );


    //- Convert OpenFOAM patch to vtkPolyData
    struct Patch
    {
        //- Convert local patch points to vtkPoints
        template<class PatchType>
        static vtkSmartPointer<vtkPoints> points(const PatchType& p);

        //- Convert patch faces to vtk polygon cells
        template<class PatchType>
        static vtkSmartPointer<vtkCellArray> faces(const PatchType& p);

        //- Convert patch points/faces to vtkPolyData
        template<class PatchType>
        static vtkSmartPointer<vtkPolyData> mesh(const PatchType& p);
209 210 211 212

        //- Convert patch face normals to vtkFloatArray
        template<class PatchType>
        static vtkSmartPointer<vtkFloatArray> faceNormals(const PatchType& p);
213 214 215
    };


216
    //- Remapping for some OpenFOAM data types (eg, symmTensor)
217
    //  \param data[in,out] The data to be remapped in-place
218
    template<class Type>
219
    inline static void remapTuple(float data[]) {}
220

221
    //- Remapping for some OpenFOAM data types (eg, symmTensor)
222
    //  \param data[in,out] The data to be remapped in-place
223
    template<class Type>
224
    inline static void remapTuple(double data[]) {}
225

226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246
    //- Copy/transcribe OpenFOAM data types to VTK format
    //  This allows a change of data type (float vs double) as well as
    //  addressing any swapping issues (eg, symmTensor)
    //
    //  \param output[out] The output scratch space. Must be long enough
    //     to hold the result.
    //  \param val[in] The input data to copy/transcribe
    template<class Type>
    inline static void foamToVtkTuple(float output[], const Type& val);

    //- Copy/transcribe OpenFOAM data types to VTK format
    //  This allows a change of data type (float vs double) as well as
    //  addressing any swapping issues (eg, symmTensor)
    //
    //  \param output[out] The output scratch space. Must be long enough
    //     to hold the result.
    //  \param val[in] The input data to copy/transcribe
    template<class Type>
    inline static void foamToVtkTuple(double output[], const Type& val);


247 248 249 250 251 252 253 254 255
    // Field Conversion Functions

    //- Copy list to pre-allocated vtk array.
    //  \return number of input items copied
    template<class Type>
    static label transcribeFloatData
    (
        vtkFloatArray* array,
        const UList<Type>& input,
256
        vtkIdType start = 0         //!< The write offset into output array
257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282
    );

    //- Create named field initialized to zero
    template<class Type>
    static vtkSmartPointer<vtkFloatArray> zeroField
    (
        const word& name,
        const label size
    );

    //- Convert field data to a vtkFloatArray
    template<class Type>
    static vtkSmartPointer<vtkFloatArray> convertFieldToVTK
    (
        const word& name,
        const UList<Type>& fld
    );

    //- An identity list of VTK_VERTEX
    static inline vtkSmartPointer<vtkCellArray> identityVertices
    (
        const label size
    );
};


283
//- Template specialization for symmTensor ordering
284
template<>
285
inline void Foam::vtk::Tools::remapTuple<Foam::symmTensor>(float data[])
286
{
287 288
    std::swap(data[1], data[3]);    // swap XY <-> YY
    std::swap(data[2], data[5]);    // swap XZ <-> ZZ
289 290 291
}


292
//- Template specialization for symmTensor ordering
293
template<>
294
inline void Foam::vtk::Tools::remapTuple<Foam::symmTensor>(double data[])
295
{
296 297
    std::swap(data[1], data[3]);    // swap XY <-> YY
    std::swap(data[2], data[5]);    // swap XZ <-> ZZ
298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320
}


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace vtk
} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#include "foamVtkToolsI.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#ifdef NoRepository
    #include "foamVtkToolsTemplates.C"
#endif

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //