Medical Imaging Interaction Toolkit  2016.11.0
Medical Imaging Interaction Toolkit
mitkAnisotropicRegistrationCommon.cpp
Go to the documentation of this file.
1 /*===================================================================
2 
3 The Medical Imaging Interaction Toolkit (MITK)
4 
5 Copyright (c) German Cancer Research Center,
6 Division of Medical and Biological Informatics.
7 All rights reserved.
8 
9 This software is distributed WITHOUT ANY WARRANTY; without
10 even the implied warranty of MERCHANTABILITY or FITNESS FOR
11 A PARTICULAR PURPOSE.
12 
13 See LICENSE.txt or http://www.mitk.org for details.
14 
15 ===================================================================*/
16 
18 #include <mitkPointSet.h>
19 #include <vtkPoints.h>
20 
22  const CovarianceMatrix &sigma_X, const CovarianceMatrix &sigma_Y)
23 {
24  WeightMatrix returnValue;
25 
26  WeightMatrix sum = sigma_X + sigma_Y;
27  vnl_svd<double> svd(sum.GetVnlMatrix());
28 
29  WeightMatrix diag;
30  diag.Fill(0.0);
31  diag[0][0] = 1.0 / sqrt(svd.W(0));
32  diag[1][1] = 1.0 / sqrt(svd.W(1));
33  diag[2][2] = 1.0 / sqrt(svd.W(2));
34 
35  WeightMatrix V; // convert vnl matrix to itk matrix...
36  for (unsigned int i = 0; i < 3; ++i)
37  for (unsigned int j = 0; j < 3; ++j)
38  V[i][j] = svd.V()[i][j];
39 
40  // add weighting matrix for point j1 (corresponding to identity transform)
41  returnValue = V * diag * V.GetTranspose();
42 
43  return returnValue;
44 }
45 
47  vtkPoints *dst,
48  const Rotation &rotation,
49  const Translation &translation)
50 {
51 #pragma omp parallel for
52  for (vtkIdType i = 0; i < src->GetNumberOfPoints(); ++i)
53  {
54  double p_in[3];
55  double p_out[3];
56  src->GetPoint(i, p_in);
57 
58  for (unsigned int j = 0; j < 3; ++j)
59  {
60  p_out[j] = p_in[0] * rotation[j][0] + p_in[1] * rotation[j][1] + p_in[2] * rotation[j][2] + translation[j];
61  }
62 
63  dst->SetPoint(i, p_out);
64  }
65 }
66 
68  MatrixList &dst,
69  const Rotation &rotation)
70 {
71  const vnl_matrix_fixed<double, 3, 3> rotationT = rotation.GetTranspose();
72 
73 #pragma omp parallel for
74  for (int i = 0; i < static_cast<int>(src.size()); ++i)
75  {
76  dst[i] = rotation * src[i] * rotationT;
77  }
78 }
79 
81  const mitk::PointSet *fixedTargets,
82  const Rotation &rotation,
83  const Translation &translation)
84 {
85  double tre = 0.0;
86 
87  for (int i = 0; i < movingTargets->GetSize(); ++i)
88  {
89  mitk::Point3D pm = movingTargets->GetPoint(i);
90  mitk::Point3D ps = fixedTargets->GetPoint(i);
91 
92  // transform point
93  pm = rotation * pm + translation;
94 
95  const double dist =
96  (ps[0] - pm[0]) * (ps[0] - pm[0]) + (ps[1] - pm[1]) * (ps[1] - pm[1]) + (ps[2] - pm[2]) * (ps[2] - pm[2]);
97 
98  tre += dist;
99  }
100 
101  tre /= movingTargets->GetSize();
102 
103  return sqrt(tre);
104 }
virtual int GetSize(unsigned int t=0) const
returns the current size of the point-list
PointType GetPoint(PointIdentifier id, int t=0) const
Get the point with ID id in world coordinates.
static Matrix3D rotation
static double ComputeTargetRegistrationError(const mitk::PointSet *movingTargets, const mitk::PointSet *fixedTargets, const Rotation &rotation, const Translation &translation)
Compute the target registration error between two point sets.
static void PropagateMatrices(const MatrixList &src, MatrixList &dst, const Rotation &rotation)
Propagate a list of matrices with a rotation matrix.
Data structure which stores a set of points. Superclass of mitk::Mesh.
Definition: mitkPointSet.h:79
static void TransformPoints(vtkPoints *src, vtkPoints *dst, const Rotation &rotation, const Translation &translation)
Transforms a point cloud with a Rotation and Translation.
static WeightMatrix CalculateWeightMatrix(const CovarianceMatrix &sigma_X, const CovarianceMatrix &sigma_Y)
Method that computes a WeightMatrix with two CovarianceMatrices.