From 6cdb80b24cc95a952e00b3bc643796616fe48566 Mon Sep 17 00:00:00 2001 From: Mark Olesen <Mark.Olesen@esi-group.com> Date: Fri, 1 Feb 2019 10:37:15 +0100 Subject: [PATCH] ENH: add vtk::Tools::Patch faceCentres() method --- src/conversion/vtk/adaptor/foamVtkTools.H | 6 +++- .../vtk/adaptor/foamVtkToolsTemplates.C | 35 ++++++++++++++++++- 2 files changed, 39 insertions(+), 2 deletions(-) diff --git a/src/conversion/vtk/adaptor/foamVtkTools.H b/src/conversion/vtk/adaptor/foamVtkTools.H index 3d35938abf..91ee3f178b 100644 --- a/src/conversion/vtk/adaptor/foamVtkTools.H +++ b/src/conversion/vtk/adaptor/foamVtkTools.H @@ -195,7 +195,7 @@ namespace Tools //- Convert OpenFOAM patch to vtkPolyData struct Patch { - //- Convert local patch points to vtkPoints + //- Return local patch points as vtkPoints template<class PatchType> static vtkSmartPointer<vtkPoints> points(const PatchType& p); @@ -210,6 +210,10 @@ namespace Tools //- Convert patch face normals to vtkFloatArray template<class PatchType> static vtkSmartPointer<vtkFloatArray> faceNormals(const PatchType& p); + + //- Return patch face centres as vtkPoints + template<class PatchType> + static vtkSmartPointer<vtkPoints> faceCentres(const PatchType& p); }; diff --git a/src/conversion/vtk/adaptor/foamVtkToolsTemplates.C b/src/conversion/vtk/adaptor/foamVtkToolsTemplates.C index ce7abd026a..0032964635 100644 --- a/src/conversion/vtk/adaptor/foamVtkToolsTemplates.C +++ b/src/conversion/vtk/adaptor/foamVtkToolsTemplates.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2017-2018 OpenCFD Ltd. + \\ / A nd | Copyright (C) 2017-2019 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -139,6 +139,39 @@ Foam::vtk::Tools::Patch::faceNormals(const PatchType& p) } +template<class PatchType> +vtkSmartPointer<vtkPoints> +Foam::vtk::Tools::Patch::faceCentres(const PatchType& p) +{ + auto vtkpoints = vtkSmartPointer<vtkPoints>::New(); + + vtkpoints->SetNumberOfPoints(p.size()); + + // Use cached values if available or loop over faces + // (avoid triggering cache) + + vtkIdType pointId = 0; + + if (p.hasFaceCentres()) + { + for (const point& pt : p.faceCentres()) + { + vtkpoints->SetPoint(pointId++, pt.v_); + } + } + else + { + for (const auto& f : p) + { + const point pt(f.centre(p.points())); + vtkpoints->SetPoint(pointId++, pt.v_); + } + } + + return vtkpoints; +} + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // // Low-Level conversions -- GitLab