Implement vtkPolygonalMeshPointNormals

vtkTriangleMeshPointNormals is a filter that computes point normals for
a triangle mesh to enable high-performance rendering. It is a fast-path
version of the vtkPolyDataNormals filter in order to be able to compute
normals for triangle meshes deforming rapidly.

The computed normals (a vtkFloatArray) are set to be the active normals
(using SetNormals()) of the PointData. The array name is "Normals", so
they can be retrieved either with `output->GetPointData()->GetNormals()`
or with `output->GetPointData()->GetArray("Normals")`.

The algorithm works by determining normals for each triangle and adding
these vectors to the triangle points. The resulting vectors at each
point are then normalized.

@warning
Normals are computed only for triangular polygons: the filter can not
handle meshes with other types of cells (Verts, Lines, Strips) or Polys
with the wrong number of components (not equal to 3).

@warning
Unlike the vtkPolyDataNormals filter, this filter does not apply any
splitting nor checks for cell orientation consistency in order to speed
up the computation. Moreover, normals are not calculated the exact same
way as the vtkPolyDataNormals filter since the triangle normals are not
normalized before being added to the point normals: those cell normals
are therefore weighted by the triangle area. This is not more nor less
correct than normalizing them before adding them, but it is much faster.

@sa
If you do not need to do high-performance rendering, you should use
vtkPolyDataNormals if your mesh is not only triangular, if you need to
split vertices at sharp edges, or if you need to check that the cell
orientations are consistent to flip inverted normals.

@sa
If you still need high-performance rendering but your input polydata is
not a triangular mesh and/or does not have consistent cell orientations
(causing inverted normals), you can still use this filter by using
vtkTriangleFilter and/or vtkCleanPolyData respectively beforehand. If
your mesh is deforming rapidly, you should be deforming the output mesh
of those two filters instead in order to only run them once.


This filter was implemented for a real-time simulation application
where a mesh vertices are deforming at every frame. With 110k points,
the framerates were the following:
- 100 fps with no normals computation
- 5 fps with vtkPolyDataNormals with default options
- 15 fps with vtkPolyDataNormals with SplittingOff & ConsistencyOff
- 80 fps with vtkPolygonalMeshPointNormals
This is an improvement of 5 to 16 times faster for that use case.

Note: for high-performance rendering of deformable meshes, two
other improvements can be made to boost the performance:
- skip frustrum culling
```
renderer->RemoveCuller(renderer->GetCullers()->GetLastItem());
```
- skip auto shift & scale of VBO:
```
openglmapper->SetVBOShiftScaleMethod(vtkOpenGLVertexBufferObject::DISABLE_SHIFT_SCALE);
```
Those two changes allow to skip calling `GetBounds()`, which requires
to recompute the bounds at each frame. This brings the framerate from
80 fps to 200 fps.

Filters/Core/CMakeLists.txt

@@ -65,6 +65,7 @@ set(Module_SRCS
   vtkThresholdPoints.cxx
   vtkTransposeTable.cxx
   vtkTriangleFilter.cxx
+  vtkTriangleMeshPointNormals.cxx
   vtkTubeFilter.cxx
   vtkUnstructuredGridQuadricDecimation.cxx
   vtkVectorDot.cxx

Filters/Core/Testing/Cxx/CMakeLists.txt

@@ -41,6 +41,7 @@ vtk_add_test_cxx(${vtk-module}CxxTests tests
   TestThreshold.cxx,NO_VALID
   TestThresholdPoints.cxx,NO_VALID
   TestTransposeTable.cxx,NO_VALID
+  TestTriangleMeshPointNormals.cxx
   TestTubeFilter.cxx
   TestUnstructuredGridQuadricDecimation.cxx,NO_VALID
   UnitTestMaskPoints.cxx,NO_VALID

Filters/Core/Testing/Cxx/TestTriangleMeshPointNormals.cxx

+/*=========================================================================
+
+  Program:   Visualization Toolkit
+  Module:    TestTriangleMeshPointNormals.cxx
+
+  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
+  All rights reserved.
+  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
+
+     This software is distributed WITHOUT ANY WARRANTY; without even
+     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+     PURPOSE.  See the above copyright notice for more information.
+=========================================================================*/
+#include "vtkActor.h"
+#include "vtkArrowSource.h"
+#include "vtkCamera.h"
+#include "vtkCleanPolyData.h"
+#include "vtkGlyph3D.h"
+#include "vtkInteractorStyleTrackballCamera.h"
+#include "vtkNew.h"
+#include "vtkPolyDataMapper.h"
+#include "vtkRegressionTestImage.h"
+#include "vtkRenderer.h"
+#include "vtkRenderWindow.h"
+#include "vtkRenderWindowInteractor.h"
+#include "vtkTesting.h"
+#include "vtkTriangleFilter.h"
+#include "vtkTriangleMeshPointNormals.h"
+#include "vtkXMLPolyDataReader.h"
+
+int TestTriangleMeshPointNormals(int argc, char* argv[])
+{
+  vtkSmartPointer<vtkTesting> testHelper =
+    vtkSmartPointer<vtkTesting>::New();
+  testHelper->AddArguments(argc, argv);
+  if (!testHelper->IsFlagSpecified("-D"))
+  {
+    std::cerr << "Error: -D /path/to/data was not specified.";
+    return EXIT_FAILURE;
+  }
+
+  std::string dataRoot = testHelper->GetDataRoot();
+  std::string fileName = dataRoot + "/Data/cow.vtp";
+  std::cout << fileName << std::endl;
+
+  // reader
+  vtkNew<vtkXMLPolyDataReader> reader;
+  reader->SetFileName(fileName.c_str());
+
+  // triangle filter
+  vtkNew<vtkTriangleFilter> triFilter;
+  triFilter->SetInputConnection(reader->GetOutputPort());
+
+  // cleaning filter
+  vtkNew<vtkCleanPolyData> cleanFilter;
+  cleanFilter->SetInputConnection(triFilter->GetOutputPort());
+
+  // normals
+  vtkNew<vtkTriangleMeshPointNormals> normFilter;
+  normFilter->SetInputConnection(cleanFilter->GetOutputPort());
+
+  // mapper, actor
+  vtkNew<vtkPolyDataMapper> mapper;
+  mapper->SetInputConnection(normFilter->GetOutputPort());
+  vtkNew<vtkActor> actor;
+  actor->SetMapper(mapper.GetPointer());
+
+  // glyphs
+  vtkNew<vtkArrowSource> glyphsource;
+  vtkNew<vtkGlyph3D> glyph;
+  glyph->SetInputConnection(normFilter->GetOutputPort());
+  glyph->SetSourceConnection(glyphsource->GetOutputPort());
+  glyph->SetVectorModeToUseNormal();
+  glyph->SetColorModeToColorByVector();
+  glyph->SetScaleModeToScaleByVector();
+  glyph->SetScaleFactor(0.5);
+  vtkNew<vtkPolyDataMapper> glyphmapper;
+  glyphmapper->SetInputConnection(glyph->GetOutputPort());
+  vtkNew<vtkActor> glyphactor;
+  glyphactor->SetMapper(glyphmapper.GetPointer());
+
+  // renderer
+  vtkNew<vtkRenderer> renderer;
+  renderer->AddActor(actor.GetPointer());
+  renderer->AddActor(glyphactor.GetPointer());
+  renderer->SetBackground(0.0, 0.0, 0.0);
+  renderer->ResetCamera();
+
+  // renderwindow, interactor
+  vtkNew<vtkRenderWindow> renWin ;
+  renWin->AddRenderer(renderer.GetPointer());
+  renWin->SetSize(300,300);
+  renWin->SetMultiSamples(0);
+  vtkNew<vtkRenderWindowInteractor> iren;
+  iren->SetRenderWindow(renWin.GetPointer());
+
+  iren->Initialize();
+  renWin->Render();
+
+  int retVal = vtkRegressionTestImage(renWin.GetPointer());
+  if (retVal == vtkRegressionTester::DO_INTERACTOR)
+  {
+    vtkNew<vtkInteractorStyleTrackballCamera> iStyle;
+    iren->SetInteractorStyle(iStyle.GetPointer());
+    renWin->SetSize(1000,1000);
+    iren->Start();
+  }
+
+  return !retVal;
+}

Filters/Core/Testing/Data/Baseline/TestTriangleMeshPointNormals.png.md5

+341887c94cfd82ba024bd626b0093f97

Filters/Core/vtkPolyDataNormals.h

@@ -46,6 +46,12 @@
  * @warning
  * Triangle strips are broken up into triangle polygons. You may want to
  * restrip the triangles.
+ *
+ * @sa
+ * For high-performance rendering, you could use vtkPolygonalMeshPointNormals
+ * if you know that you have a polygonal mesh which does not require splitting
+ * nor consistency check on the cell orientations.
+ *
 */
 
 #ifndef vtkPolyDataNormals_h

Filters/Core/vtkTriangleMeshPointNormals.cxx

+/*=========================================================================
+
+  Program:   Visualization Toolkit
+  Module:    vtkTriangleMeshPointNormals.cxx
+
+  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
+  All rights reserved.
+  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
+
+     This software is distributed WITHOUT ANY WARRANTY; without even
+     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+     PURPOSE.  See the above copyright notice for more information.
+
+=========================================================================*/
+#include "vtkTriangleMeshPointNormals.h"
+
+#include "vtkCellArray.h"
+#include "vtkCellData.h"
+#include "vtkFloatArray.h"
+#include "vtkMath.h"
+#include "vtkInformation.h"
+#include "vtkInformationVector.h"
+#include "vtkObjectFactory.h"
+#include "vtkPointData.h"
+#include "vtkPolyData.h"
+
+vtkStandardNewMacro(vtkTriangleMeshPointNormals);
+
+namespace
+{
+template<typename ptDataType>
+const char* computeNormalsDirection(vtkPolyData* mesh, float *n)
+{
+  ptDataType *points = reinterpret_cast<ptDataType*>
+    (mesh->GetPoints()->GetData()->GetVoidPointer(0));
+  vtkIdType *cells = reinterpret_cast<vtkIdType*>
+    (mesh->GetPolys()->GetData()->GetVoidPointer(0));
+  vtkIdType v0Offset, v1Offset, v2Offset;
+  ptDataType *p0, *p1, *p2;
+  float a[3], b[3], tn[3];
+  for (vtkIdType i = 0; i < mesh->GetNumberOfPolys(); ++i)
+  {
+    // First value in cellArray indicates number of points in cell.
+    // We need 3 to compute normals.
+    if (*cells == 3)
+    {
+      // vertex offsets
+      v0Offset = 3 * cells[1];
+      v1Offset = 3 * cells[2];
+      v2Offset = 3 * cells[3];
+      // pointer to each vertex
+      p0 = points + v0Offset;
+      p1 = points + v1Offset;
+      p2 = points + v2Offset;
+      // two vectors
+      a[0] = p2[0] - p1[0];
+      a[1] = p2[1] - p1[1];
+      a[2] = p2[2] - p1[2];
+      b[0] = p0[0] - p1[0];
+      b[1] = p0[1] - p1[1];
+      b[2] = p0[2] - p1[2];
+      // cell normals by cross-product
+      // (no need to normalize those + it's faster not to)
+      tn[0] = (a[1] * b[2] - a[2] * b[1]);
+      tn[1] = (a[2] * b[0] - a[0] * b[2]);
+      tn[2] = (a[0] * b[1] - a[1] * b[0]);
+      // append triangle normals to point normals
+      *(n + v0Offset) += tn[0];
+      *(n + v0Offset + 1) += tn[1];
+      *(n + v0Offset + 2) += tn[2];
+      *(n + v1Offset) += tn[0];
+      *(n + v1Offset + 1) += tn[1];
+      *(n + v1Offset + 2) += tn[2];
+      *(n + v2Offset) += tn[0];
+      *(n + v2Offset + 1) += tn[1];
+      *(n + v2Offset + 2) += tn[2];
+      // move to next cells
+      cells += (*cells + 1); // current cell nbr of pts + 1
+    }
+    // If degenerate cell
+    else if (*cells < 3)
+    {
+      return "Some cells are degenerate (less than 3 points). "
+             "Use vtkCleanPolyData beforehand to correct this.";
+    }
+    // If cell not triangle
+    else
+    {
+      return "Some cells have too many points (more than 3 points). "
+             "Use vtkTriangulate to correct this.";
+    }
+  }
+  return NULL;
+}
+}
+
+// Generate normals for polygon meshes
+int vtkTriangleMeshPointNormals::RequestData(
+  vtkInformation *vtkNotUsed(request),
+  vtkInformationVector **inputVector,
+  vtkInformationVector *outputVector)
+{
+  // get the info objects
+  vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
+  vtkInformation *outInfo = outputVector->GetInformationObject(0);
+
+  // get the input and output
+  vtkPolyData *input = vtkPolyData::SafeDownCast(
+    inInfo->Get(vtkDataObject::DATA_OBJECT()));
+  vtkPolyData *output = vtkPolyData::SafeDownCast(
+    outInfo->Get(vtkDataObject::DATA_OBJECT()));
+
+  vtkDebugMacro(<<"Generating surface normals");
+
+  vtkIdType numPts = input->GetNumberOfPoints(); // nbr of points from input
+  if ( numPts < 1 )
+  {
+    vtkDebugMacro(<<"No data to generate normals for!");
+    return 1;
+  }
+
+  if (input->GetVerts()->GetNumberOfCells() != 0 ||
+      input->GetLines()->GetNumberOfCells() != 0 ||
+      input->GetStrips()->GetNumberOfCells() != 0)
+  {
+    vtkErrorMacro(<< "Can not compute normals for a mesh with Verts, Lines or Strips, as it will "
+                  << "corrupt the number of points used during the normals computation."
+                  << "Make sure your input PolyData only has triangles (Polys with 3 components)).");
+    return 0;
+  }
+
+  // Copy structure and cell data
+  output->CopyStructure(input);
+  output->GetCellData()->PassData(input->GetCellData());
+
+  // If there is nothing to do, pass the point data through
+  if (input->GetNumberOfPolys() < 1)
+  {
+    output->GetPointData()->PassData(input->GetPointData());
+    return 1;
+  }
+  // Else pass everything but normals
+  output->GetPointData()->CopyNormalsOff();
+  output->GetPointData()->PassData(input->GetPointData());
+
+  // Prepare array for normals
+  vtkFloatArray *normals = vtkFloatArray::New();
+  normals->SetNumberOfComponents(3);
+  normals->SetNumberOfTuples(numPts);
+  normals->SetName("Normals");
+  normals->FillValue(0.0);
+  output->GetPointData()->SetNormals(normals);
+  normals->Delete();
+
+  this->UpdateProgress(0.1);
+
+  // Compute normals direction
+  float *n = reinterpret_cast<float*>(normals->GetVoidPointer(0));
+  switch (output->GetPoints()->GetDataType())
+  {
+    vtkTemplateMacro(
+      const char *warning = computeNormalsDirection<VTK_TT>(output, n);
+      if(warning)
+      {
+        vtkWarningMacro( << warning);
+      }
+    );
+  }
+  this->UpdateProgress(0.5);
+
+  // Normalize point normals
+  float l;
+  unsigned int i3;
+  for (vtkIdType i = 0; i < numPts; ++i)
+  {
+    i3 = i * 3;
+    if ((l = sqrt(n[i3] * n[i3] +
+        n[i3 + 1] * n[i3 + 1] +
+        n[i3 + 2] * n[i3 + 2])) != 0.0)
+    {
+        n[i3] /= l;
+        n[i3 + 1] /= l;
+        n[i3 + 2] /= l;
+    }
+  }
+  this->UpdateProgress(0.9);
+
+  // Update modified time
+  normals->Modified();
+
+  return 1;
+}
+
+void vtkTriangleMeshPointNormals::PrintSelf(ostream& os, vtkIndent indent)
+{
+  this->Superclass::PrintSelf(os,indent);
+}

Filters/Core/vtkTriangleMeshPointNormals.h

+/*=========================================================================
+
+  Program:   Visualization Toolkit
+  Module:    vtkTriangleMeshPointNormals.h
+
+  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
+  All rights reserved.
+  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
+
+     This software is distributed WITHOUT ANY WARRANTY; without even
+     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+     PURPOSE.  See the above copyright notice for more information.
+
+=========================================================================*/
+/**
+ * @class   vtkTriangleMeshPointNormals
+ * @brief   compute point normals for triangle mesh
+ *
+ * vtkTriangleMeshPointNormals is a filter that computes point normals for
+ * a triangle mesh to enable high-performance rendering. It is a fast-path
+ * version of the vtkPolyDataNormals filter in order to be able to compute
+ * normals for triangle meshes deforming rapidly.
+ *
+ * The computed normals (a vtkFloatArray) are set to be the active normals
+ * (using SetNormals()) of the PointData. The array name is "Normals", so
+ * they can be retrieved either with `output->GetPointData()->GetNormals()`
+ * or with `output->GetPointData()->GetArray("Normals")`.
+ *
+ * The algorithm works by determining normals for each triangle and adding
+ * these vectors to the triangle points. The resulting vectors at each
+ * point are then normalized.
+ *
+ * @warning
+ * Normals are computed only for triangular polygons: the filter can not
+ * handle meshes with other types of cells (Verts, Lines, Strips) or Polys
+ * with the wrong number of components (not equal to 3).
+ *
+ * @warning
+ * Unlike the vtkPolyDataNormals filter, this filter does not apply any
+ * splitting nor checks for cell orientation consistency in order to speed
+ * up the computation. Moreover, normals are not calculated the exact same
+ * way as the vtkPolyDataNormals filter since the triangle normals are not
+ * normalized before being added to the point normals: those cell normals
+ * are therefore weighted by the triangle area. This is not more nor less
+ * correct than normalizing them before adding them, but it is much faster.
+ *
+ * @sa
+ * If you do not need to do high-performance rendering, you should use
+ * vtkPolyDataNormals if your mesh is not only triangular, if you need to
+ * split vertices at sharp edges, or if you need to check that the cell
+ * orientations are consistent to flip inverted normals.
+ *
+ * @sa
+ * If you still need high-performance rendering but your input polydata is
+ * not a triangular mesh and/or does not have consistent cell orientations
+ * (causing inverted normals), you can still use this filter by using
+ * vtkTriangleFilter and/or vtkCleanPolyData respectively beforehand. If
+ * your mesh is deforming rapidly, you should be deforming the output mesh
+ * of those two filters instead in order to only run them once.
+ *
+*/
+
+#ifndef vtkTriangleMeshPointNormals_h
+#define vtkTriangleMeshPointNormals_h
+
+#include "vtkFiltersCoreModule.h" // For export macro
+#include "vtkPolyDataAlgorithm.h"
+
+class vtkPolyData;
+
+
+class VTKFILTERSCORE_EXPORT vtkTriangleMeshPointNormals : public vtkPolyDataAlgorithm
+{
+public:
+  vtkTypeMacro(vtkTriangleMeshPointNormals,vtkPolyDataAlgorithm);
+  void PrintSelf(ostream& os, vtkIndent indent) VTK_OVERRIDE;
+
+  static vtkTriangleMeshPointNormals *New();
+
+protected:
+  vtkTriangleMeshPointNormals() {}
+  ~vtkTriangleMeshPointNormals() VTK_OVERRIDE {}
+
+  // Usual data generation method
+  int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *) VTK_OVERRIDE;
+
+private:
+  vtkTriangleMeshPointNormals(const vtkTriangleMeshPointNormals&) VTK_DELETE_FUNCTION;
+  void operator=(const vtkTriangleMeshPointNormals&) VTK_DELETE_FUNCTION;
+};
+
+#endif