Compute Normals of a Mesh#

Synopsis#

Compute normals of a mesh

Results#

Input mesh

Input mesh#

Interactive input mesh

Output:

Index * Point * Normal
0 * [0.461602, -0.018446, -0.115989] * [0.964656, 0.0898966, -0.247704]
1 * [0.202784, 0.359481, -0.217858] * [0.398458, 0.53089, -0.747922]
2 * [-0.481417, -0.046787, 0.167623] * [-0.987132, 0.0342443, 0.1562]
3 * [-0.028388, -0.193602, -0.329748] * [-0.73488, -0.637902, -0.230287]
4 * [-0.087584, 0.116513, -0.399884] * [-0.396325, 0.718368, -0.571729]
5 * [0.310978, -0.239208, -0.176146] * [-0.256754, 0.55162, -0.793595]
6 * [-0.303314, 0.043035, -0.515252] * [-0.26465, 0.684736, -0.679041]
7 * [-0.256717, -0.39897, 0.177429] * [0.0907732, -0.987892, 0.125819]
8 * [0.393919, -0.388691, -0.207638] * [0.290631, -0.955747, -0.0456223]
9 * [-0.225987, -0.368634, -0.367823] * [0.154093, 0.365541, -0.917952]
10 * [-0.142556, -0.187471, -0.613449] * [0.902758, -0.0825592, -0.422153]
11 * [0.405923, -0.12839, -0.433595] * [0.432278, -0.901416, 0.0241831]
12 * [0.107094, 0.05873, -0.53405] * [0.452665, 0.664453, -0.594639]
13 * [0.179737, -0.358912, 0.264171] * [0.622873, -0.779647, 0.0646554]
14 * [0.270291, 0.167529, -0.267412] * [0.282452, 0.630941, -0.722588]
15 * [0.190508, -0.139941, -0.30917] * [0.642572, -0.758423, 0.109063]
16 * [-0.405901, 0.088723, -0.347001] * [-0.736406, 0.557068, -0.383904]
17 * [0.150462, -0.314033, -0.695425] * [-0.171163, -0.38232, -0.908039]
18 * [-0.134383, -0.384493, 0.060186] * [-0.0256952, -0.994633, 0.100226]
19 * [-0.333505, -0.116041, -0.567548] * [-0.847824, -0.186998, -0.496212]
20 * [-0.133643, -0.267196, -0.238155] * [0.123997, 0.550013, -0.8259]
21 * [-0.274406, 0.185852, -0.366349] * [-0.0630695, 0.682696, -0.727976]
22 * [-0.34352, -0.372556, 0.17102] * [-0.565292, -0.779407, 0.270129]
23 * [-0.468522, 0.005617, -0.307273] * [-0.945988, 0.0866878, -0.312397]
24 * [0.012837, 0.243987, -0.417788] * [-0.288839, 0.714998, -0.636671]
25 * [0.004712, 0.118755, -0.486604] * [-0.451072, 0.542366, -0.708783]

[...]

1300 * [0.429046, 0.119465, 0.002968] * [0.838127, 0.413114, -0.356201]
1301 * [-0.346539, -0.163042, 0.316913] * [-0.323662, 0.716134, 0.618381]
1302 * [-0.111388, 0.551166, -0.050984] * [-0.558618, 0.686741, -0.465115]
1303 * [0.376285, 0.211351, 0.347429] * [0.754427, -0.0746173, 0.652129]
1304 * [-0.026005, 0.583073, 0.081228] * [0.60426, 0.298334, 0.738828]
1305 * [-0.196168, -0.079437, 0.372088] * [-0.859283, 0.463078, 0.217237]
1306 * [-0.000569, 0.431162, 0.415514] * [0.683671, 0.70426, 0.191341]

Code#

C++#

#include "itkVector.h"
#include "itkQuadEdgeMesh.h"
#include "itkMeshFileReader.h"
#include "itkQuadEdgeMeshExtendedTraits.h"
#include "itkNormalQuadEdgeMeshFilter.h"

int
main(int argc, char * argv[])
{
  if (argc < 2)
  {
    std::cerr << "Usage:" << std::endl;
    std::cerr << argv[0] << "<InputFileName> <WeightType>" << std::endl;
    std::cerr << "Weight type: " << std::endl;
    std::cerr << "  * 0:  GOURAUD" << std::endl;
    std::cerr << "  * 1:  THURMER" << std::endl;
    std::cerr << "  * 2:  AREA" << std::endl;
    return EXIT_FAILURE;
  }

  constexpr unsigned int Dimension = 3;
  using CoordType = double;

  using InputMeshType = itk::QuadEdgeMesh<CoordType, Dimension>;

  using VectorType = itk::Vector<CoordType, Dimension>;

  using Traits =
    itk::QuadEdgeMeshExtendedTraits<VectorType, Dimension, 2, CoordType, CoordType, VectorType, bool, bool>;

  using ReaderType = itk::MeshFileReader<InputMeshType>;
  auto reader = ReaderType::New();
  reader->SetFileName(argv[1]);
  reader->Update();

  using OutputMeshType = itk::QuadEdgeMesh<VectorType, Dimension, Traits>;

  using NormalFilterType = itk::NormalQuadEdgeMeshFilter<InputMeshType, OutputMeshType>;
  NormalFilterType::WeightEnum weight_type;

  int param = std::stoi(argv[2]);

  if ((param < 0) || (param > 2))
  {
    std::cerr << "Weight type must be either: " << std::endl;
    std::cerr << "   * 0:  GOURAUD" << std::endl;
    std::cerr << "   * 1:  THURMER" << std::endl;
    std::cerr << "   * 2:  AREA" << std::endl;
    return EXIT_FAILURE;
  }
  else
  {
    switch (param)
    {
      default:
      case 0:
        weight_type = itk::NormalQuadEdgeMeshFilterEnums::Weight::GOURAUD;
        break;
      case 1:
        weight_type = itk::NormalQuadEdgeMeshFilterEnums::Weight::THURMER;
        break;
      case 2:
        weight_type = itk::NormalQuadEdgeMeshFilterEnums::Weight::AREA;
        break;
    }
  }

  auto normals = NormalFilterType::New();
  normals->SetInput(reader->GetOutput());
  normals->SetWeight(weight_type);
  normals->Update();

  OutputMeshType::Pointer output = normals->GetOutput();

  OutputMeshType::PointsContainerPointer  points = output->GetPoints();
  OutputMeshType::PointsContainerIterator p_it = points->Begin();

  OutputMeshType::PointDataContainerPointer  container = output->GetPointData();
  OutputMeshType::PointDataContainerIterator d_it = container->Begin();

  std::cout << "Index * Point * Normal" << std::endl;

  while (p_it != points->End())
  {
    std::cout << p_it.Index() << " * ";
    std::cout << p_it.Value() << " * ";
    std::cout << d_it.Value() << std::endl;

    ++p_it;
    ++d_it;
  }

  return EXIT_SUCCESS;
}

Classes demonstrated#

template<typename TInputMesh, typename TOutputMesh>
class NormalQuadEdgeMeshFilter : public itk::QuadEdgeMeshToQuadEdgeMeshFilter<TInputMesh, TOutputMesh>

Filter which computes normals to faces and vertices and store it in the output mesh. Normals to face are first computed, then normals to vertices are computed as linear combination of neighbor face normals, i.e.

n_v = \frac{\sum_{i=0}^{N_f} \omega_i \cdot n_i}{\| \sum_{k=0}^{N_f} \omega_x \cdot n_k\|}

The difference between each method relies in the definition of the weight \omega_i that you can specify by the method SetWeight.

  • GOURAUD \omega_i = 1 [1]

  • THURMER \omega_i = Angle of the considered triangle at the given vertex} [2]

  • AREA \omega_i = Area(t_i) [3]

These weights are defined in the literature [45], [56] and [127].

Todo:

Fix run-time issues regarding the difference between the Traits of TInputMesh and the one of TOutputMesh. Right now, it only works if TInputMesh::MeshTraits == TOutputMesh::MeshTraits (and of course it requires that the output have some itk::Vector for point data and cell data.

ITK Sphinx Examples:

Note

By default the weight is set to the THURMER weight.

See itk::NormalQuadEdgeMeshFilter for additional documentation.