Commit b2b4517c authored by Mark Olesen's avatar Mark Olesen
Browse files

ENH: code cleanup for vtk transcription

- use cellCentres() instead of volField equivalent for vtk conversion

- make looping variables more consistent

- centralize the transcription of OpenFOAM -> vtk tuples
parent fbf476d6
......@@ -40,11 +40,11 @@ InClass
#include "areaFaMesh.H"
#include "areaFields.H"
// vtk includes
#include "vtkFloatArray.h"
#include "vtkCellData.h"
#include "vtkPointData.h"
#include "vtkSmartPointer.h"
// VTK includes
#include <vtkFloatArray.h>
#include <vtkCellData.h>
#include <vtkPointData.h>
#include <vtkSmartPointer.h>
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
//
......@@ -68,10 +68,8 @@ void Foam::vtkPVFoam::convertVolField
autoPtr<GeometricField<Type, pointPatchField, pointMesh>> ptfPtr;
if (interpField)
{
if (debug)
{
Info<< "convertVolField interpolating:" << fld.name() << nl;
}
DebugInfo
<< "convertVolField interpolating:" << fld.name() << nl;
ptfPtr.reset
(
......@@ -102,7 +100,7 @@ void Foam::vtkPVFoam::convertVolField
foamVtpData& vtpData = iter.object();
auto dataset = vtpData.dataset;
const labelList& patchIds = vtpData.additionalIds();
const labelUList& patchIds = vtpData.additionalIds();
if (patchIds.empty())
{
......@@ -291,13 +289,13 @@ void Foam::vtkPVFoam::convertVolFields
catch (Foam::IOerror& ioErr)
{
ioErr.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
catch (Foam::error& err)
{
// Bit of trickery to get the original message
err.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
// Restore previous exception throwing state
......@@ -373,13 +371,13 @@ void Foam::vtkPVFoam::convertDimFields
catch (Foam::IOerror& ioErr)
{
ioErr.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
catch (Foam::error& err)
{
// Bit of trickery to get the original message
err.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
// Restore previous exception throwing state
......@@ -501,13 +499,13 @@ void Foam::vtkPVFoam::convertAreaFields
catch (Foam::IOerror& ioErr)
{
ioErr.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
catch (Foam::error& err)
{
// Bit of trickery to get the original message
err.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
// Restore previous exception throwing state
......@@ -545,10 +543,8 @@ void Foam::vtkPVFoam::convertPointFields
continue;
}
if (debug)
{
Info<< "convertPointFields : " << fieldName << nl;
}
DebugInfo
<< "convertPointFields : " << fieldName << nl;
// Throw FatalError, FatalIOError as exceptions
const bool throwingError = FatalError.throwExceptions();
......@@ -581,7 +577,7 @@ void Foam::vtkPVFoam::convertPointFields
foamVtpData& vtpData = iter.object();
auto dataset = vtpData.dataset;
const labelList& patchIds = vtpData.additionalIds();
const labelUList& patchIds = vtpData.additionalIds();
if (patchIds.size() != 1)
{
continue;
......@@ -645,13 +641,13 @@ void Foam::vtkPVFoam::convertPointFields
catch (Foam::IOerror& ioErr)
{
ioErr.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
catch (Foam::error& err)
{
// Bit of trickery to get the original message
err.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
// Restore previous exception throwing state
......@@ -706,7 +702,6 @@ vtkSmartPointer<vtkFloatArray> Foam::vtkPVFoam::convertPointField
const foamVtuData& vtuData
)
{
const int nComp(pTraits<Type>::nComponents);
const labelUList& addPointCellLabels = vtuData.additionalIds();
const labelUList& pointMap = vtuData.pointMap();
......@@ -714,7 +709,7 @@ vtkSmartPointer<vtkFloatArray> Foam::vtkPVFoam::convertPointField
const label nPoints = (pointMap.size() ? pointMap.size() : pfld.size());
auto data = vtkSmartPointer<vtkFloatArray>::New();
data->SetNumberOfComponents(nComp);
data->SetNumberOfComponents(static_cast<int>(pTraits<Type>::nComponents));
data->SetNumberOfTuples(nPoints + addPointCellLabels.size());
// Note: using the name of the original volField
......@@ -729,75 +724,54 @@ vtkSmartPointer<vtkFloatArray> Foam::vtkPVFoam::convertPointField
data->SetName(pfld.name().c_str());
}
if (debug)
{
Info<< "convert Point field: "
<< pfld.name()
<< " size=" << (nPoints + addPointCellLabels.size())
<< " (" << nPoints << " + " << addPointCellLabels.size()
<< ") nComp=" << nComp << nl;
}
DebugInfo
<< "Convert point field: " << pfld.name()
<< " size=" << (nPoints + addPointCellLabels.size())
<< " (" << nPoints << " + " << addPointCellLabels.size()
<< ") nComp=" << static_cast<int>(pTraits<Type>::nComponents) << nl;
float vec[nComp];
float scratch[pTraits<Type>::nComponents];
label pointi = 0;
vtkIdType pointi = 0;
if (pointMap.size())
{
forAll(pointMap, i)
for (const label meshPointi : pointMap)
{
const Type& t = pfld[pointMap[i]];
for (direction d=0; d<nComp; ++d)
{
vec[d] = component(t, d);
}
remapTuple<Type>(vec);
data->SetTuple(pointi++, vec);
vtk::Tools::foamToVtkTuple(scratch, pfld[meshPointi]);
data->SetTuple(pointi++, scratch);
}
}
else
{
forAll(pfld, i)
for (const Type& val : pfld)
{
const Type& t = pfld[i];
for (direction d=0; d<nComp; ++d)
{
vec[d] = component(t, d);
}
remapTuple<Type>(vec);
data->SetTuple(pointi++, vec);
vtk::Tools::foamToVtkTuple(scratch, val);
data->SetTuple(pointi++, scratch);
}
}
// Continue with additional points
// - correspond to cell centres
if (notNull(vfld))
{
forAll(addPointCellLabels, apI)
for (const label meshCelli : addPointCellLabels)
{
const Type& t = vfld[addPointCellLabels[apI]];
for (direction d=0; d<nComp; ++d)
{
vec[d] = component(t, d);
}
remapTuple<Type>(vec);
data->SetTuple(pointi++, vec);
vtk::Tools::foamToVtkTuple(scratch, vfld[meshCelli]);
data->SetTuple(pointi++, scratch);
}
}
else
{
forAll(addPointCellLabels, apI)
for (const label meshCelli : addPointCellLabels)
{
Type t = interpolatePointToCell(pfld, addPointCellLabels[apI]);
for (direction d=0; d<nComp; ++d)
{
vec[d] = component(t, d);
}
remapTuple<Type>(vec);
data->SetTuple(pointi++, vec);
vtk::Tools::foamToVtkTuple
(
scratch,
interpolatePointToCell(pfld, meshCelli)
);
data->SetTuple(pointi++, scratch);
}
}
......@@ -846,13 +820,13 @@ void Foam::vtkPVFoam::convertLagrangianFields
catch (Foam::IOerror& ioErr)
{
ioErr.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
catch (Foam::error& err)
{
// Bit of trickery to get the original message
err.write(Warning, false);
Info << nl << endl;
Info<< nl << endl;
}
// Restore previous exception throwing state
......@@ -876,53 +850,44 @@ Foam::vtkPVFoam::convertFaceFieldToVTK
const labelUList& faceLabels
) const
{
if (debug)
{
Info<< "convert face field: "
<< fld.name()
<< " size=" << faceLabels.size()
<< " nComp=" << int(pTraits<Type>::nComponents) << nl;
}
DebugInfo
<< "Convert face field: " << fld.name()
<< " size=" << faceLabels.size()
<< " nComp=" << static_cast<int>(pTraits<Type>::nComponents) << nl;
const fvMesh& mesh = fld.mesh();
const int nComp(pTraits<Type>::nComponents);
const label nInternalFaces = mesh.nInternalFaces();
const labelList& faceOwner = mesh.faceOwner();
const labelList& faceNeigh = mesh.faceNeighbour();
const labelUList& faceOwner = mesh.faceOwner();
const labelUList& faceNeigh = mesh.faceNeighbour();
auto data = vtkSmartPointer<vtkFloatArray>::New();
data->SetName(fld.name().c_str());
data->SetNumberOfComponents(nComp);
data->SetNumberOfComponents(static_cast<int>(pTraits<Type>::nComponents));
data->SetNumberOfTuples(faceLabels.size());
float scratch[nComp];
// Interior faces: average owner/neighbour
// Boundary faces: the owner value
forAll(faceLabels, idx)
float scratch[pTraits<Type>::nComponents];
vtkIdType faceId = 0;
for (const label meshFacei : faceLabels)
{
const label faceNo = faceLabels[idx];
if (faceNo < nInternalFaces)
if (meshFacei < nInternalFaces)
{
Type t = 0.5*(fld[faceOwner[faceNo]] + fld[faceNeigh[faceNo]]);
Type val =
0.5*(fld[faceOwner[meshFacei]] + fld[faceNeigh[meshFacei]]);
for (direction d=0; d<nComp; ++d)
{
scratch[d] = component(t, d);
}
vtk::Tools::foamToVtkTuple(scratch, val);
}
else
{
const Type& t = fld[faceOwner[faceNo]];
for (direction d=0; d<nComp; ++d)
{
scratch[d] = component(t, d);
}
const Type& val = fld[faceOwner[meshFacei]];
vtk::Tools::foamToVtkTuple(scratch, val);
}
remapTuple<Type>(scratch);
data->SetTuple(idx, scratch);
data->SetTuple(faceId++, scratch);
}
return data;
......@@ -937,35 +902,26 @@ Foam::vtkPVFoam::convertVolFieldToVTK
const foamVtuData& vtuData
) const
{
const int nComp(pTraits<Type>::nComponents);
const labelUList& cellMap = vtuData.cellMap();
auto data = vtkSmartPointer<vtkFloatArray>::New();
data->SetName(fld.name().c_str());
data->SetNumberOfComponents(nComp);
data->SetNumberOfComponents(static_cast<int>(pTraits<Type>::nComponents));
data->SetNumberOfTuples(cellMap.size());
if (debug)
{
Info<< "convert volField: "
<< fld.name()
<< " size=" << cellMap.size()
<< " (" << fld.size() << " + "
<< (cellMap.size() - fld.size())
<< ") nComp=" << nComp << nl;
}
DebugInfo
<< "Convert volField: " << fld.name() << " size="
<< cellMap.size() << " (field "
<< fld.size() << ") nComp="
<< static_cast<int>(pTraits<Type>::nComponents) << nl;
float scratch[nComp];
forAll(cellMap, idx)
{
const Type& t = fld[cellMap[idx]];
for (direction d=0; d<nComp; ++d)
{
scratch[d] = component(t, d);
}
remapTuple<Type>(scratch);
float scratch[pTraits<Type>::nComponents];
data->SetTuple(idx, scratch);
vtkIdType celli = 0;
for (const label meshCelli : cellMap)
{
vtk::Tools::foamToVtkTuple(scratch, fld[meshCelli]);
data->SetTuple(celli++, scratch);
}
return data;
......
......@@ -172,7 +172,7 @@ void Foam::vtkPVFoam::convertMeshPatches()
{
// Point movement on single patch is OK
const labelList& patchIds = vtpData.additionalIds();
const labelUList& patchIds = vtpData.additionalIds();
if (patchIds.size() == 1)
{
vtkgeom = vtpData.getCopy();
......
......@@ -49,23 +49,21 @@ vtkSmartPointer<vtkPoints> Foam::vtkPVFoam::foamVtuData::points
const pointField& pts = mesh.points();
// Additional cell centres
const labelList& addPoints = this->additionalIds();
const labelUList& addPoints = this->additionalIds();
vtkpoints->SetNumberOfPoints(pts.size() + addPoints.size());
// Normal points
label pointId = 0;
vtkIdType pointId = 0;
for (const point& p : pts)
{
vtkpoints->SetPoint(pointId, p.v_);
++pointId;
vtkpoints->SetPoint(pointId++, p.v_);
}
// Cell centres
for (const label ptId : addPoints)
for (const label meshCelli : addPoints)
{
vtkpoints->SetPoint(pointId, mesh.C()[ptId].v_);
++pointId;
vtkpoints->SetPoint(pointId++, mesh.cellCentres()[meshCelli].v_);
}
return vtkpoints;
......@@ -85,23 +83,21 @@ vtkSmartPointer<vtkPoints> Foam::vtkPVFoam::foamVtuData::points
const pointField& pts = mesh.points();
// Additional cell centres
const labelList& addPoints = this->additionalIds();
const labelUList& addPoints = this->additionalIds();
vtkpoints->SetNumberOfPoints(pointMap.size() + addPoints.size());
// Normal points
label pointId = 0;
for (const label ptId : pointMap)
vtkIdType pointId = 0;
for (const label meshPointi : pointMap)
{
vtkpoints->SetPoint(pointId, pts[ptId].v_);
++pointId;
vtkpoints->SetPoint(pointId++, pts[meshPointi].v_);
}
// Cell centres
for (const label ptId : addPoints)
for (const label meshCelli : addPoints)
{
vtkpoints->SetPoint(pointId, mesh.C()[ptId].v_);
++pointId;
vtkpoints->SetPoint(pointId++, mesh.cellCentres()[meshCelli].v_);
}
return vtkpoints;
......
......@@ -49,6 +49,7 @@ SourceFiles
#include "pointField.H"
#include "symmTensor.H"
// VTK includes
#include <vtkCellArray.h>
#include <vtkFloatArray.h>
#include <vtkDoubleArray.h>
......@@ -209,13 +210,36 @@ public:
//- Remapping for some OpenFOAM data types (eg, symmTensor)
// \param data[in,out] The data to be remapped in-place
template<class Type>
inline static void remapTuple(float data[]) {}
//- Remapping for some OpenFOAM data types (eg, symmTensor)
// \param data[in,out] The data to be remapped in-place
template<class Type>
inline static void remapTuple(double data[]) {}
//- 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);
// Field Conversion Functions
//- Copy list to pre-allocated vtk array.
......@@ -225,7 +249,7 @@ public:
(
vtkFloatArray* array,
const UList<Type>& input,
const label start = 0
vtkIdType start = 0 //!< The write offset into output array
);
//- Create named field initialized to zero
......@@ -270,7 +294,6 @@ inline void Foam::vtk::Tools::remapTuple<Foam::symmTensor>(double data[])
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace vtk
......
......@@ -87,4 +87,34 @@ inline vtkSmartPointer<vtkCellArray> Foam::vtk::Tools::identityVertices
};
template<class Type>
inline void Foam::vtk::Tools::foamToVtkTuple
(
float output[],
const Type& val
)
{
for (direction cmpt=0; cmpt < pTraits<Type>::nComponents; ++cmpt)
{
output[cmpt] = component(val, cmpt);
}
remapTuple<Type>(output);
}
template<class Type>
inline void Foam::vtk::Tools::foamToVtkTuple
(
double output[],
const Type& val
)
{
for (direction cmpt=0; cmpt < pTraits<Type>::nComponents; ++cmpt)
{
output[cmpt] = component(val, cmpt);
}
remapTuple<Type>(output);
}
// ************************************************************************* //
......@@ -50,12 +50,11 @@ Foam::vtk::Tools::Patch::points(const PatchType& p)
auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
vtkpoints->SetNumberOfPoints(pts.size());
vtkIdType pointId = 0;
vtkIdType pointId = 0;
for (const point& p : pts)
{
vtkpoints->SetPoint(pointId, p.v_);
++pointId;
vtkpoints->SetPoint(pointId++, p.v_);
}
return vtkpoints;
......@@ -127,8 +126,7 @@ Foam::vtk::Tools::Patch::faceNormals(const PatchType& p)
vtkIdType faceId = 0;
for (const vector& n : norms)
{
array->SetTuple(faceId, n.v_);
++faceId;
array->SetTuple(faceId++, n.v_);
}
return array;
......@@ -145,7 +143,7 @@ Foam::label Foam::vtk::Tools::transcribeFloatData
(
vtkFloatArray* array,
const UList<Type>& input,
const label start
vtkIdType start
)
{
const int nComp(pTraits<Type>::nComponents);
......@@ -162,7 +160,7 @@ Foam::label Foam::vtk::Tools::transcribeFloatData
}
const vtkIdType maxSize = array->GetNumberOfTuples();
const vtkIdType endPos = vtkIdType(start) + vtkIdType(input.size());
const vtkIdType endPos = start + vtkIdType(input.size());
if (!maxSize)
{
......@@ -174,7 +172,7 @@ Foam::label Foam::vtk::Tools::transcribeFloatData