File: PointSetWithCovariantVectors.cxx

package info (click to toggle)
insighttoolkit4 4.13.3withdata-dfsg1-4
  • links: PTS, VCS
  • area: main
  • in suites: bullseye
  • size: 489,260 kB
  • sloc: cpp: 557,342; ansic: 146,850; fortran: 34,788; python: 16,572; sh: 2,187; lisp: 2,070; tcl: 993; java: 362; perl: 200; makefile: 129; csh: 81; pascal: 69; xml: 19; ruby: 10
file content (194 lines) | stat: -rw-r--r-- 7,129 bytes parent folder | download | duplicates (5)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
/*=========================================================================
 *
 *  Copyright Insight Software Consortium
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *         http://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/

//  Software Guide : BeginLatex
//
//  It is common to represent geometric objects by using points on their surfaces
//  and normals associated with those points.  This structure can be easily
//  instantiated with the \doxygen{PointSet} class.
//
//  The natural class for representing normals to surfaces and
//  gradients of functions is the \doxygen{CovariantVector}. A
//  covariant vector differs from a vector in the way it behaves
//  under affine transforms, in particular under anisotropic
//  scaling. If a covariant vector represents the gradient of a
//  function, the transformed covariant vector will still be the valid
//  gradient of the transformed function, a property which would not
//  hold with a regular vector.
//
//  \index{itk::PointSet!itk::CovariantVector}
//  \index{itk::CovariantVector!itk::PointSet}
//
//  The following example demonstrates how a \code{CovariantVector} can
//  be used as the \code{PixelType} for the \code{PointSet} class.  The
//  example illustrates how a deformable model could move under
//  the influence of the gradient of a potential function.
//
//  In order to use the CovariantVector class it is necessary to
//  include its header file along with the header of the point set.
//
//  \index{itk::CovariantVector!Header}
//
//  Software Guide : EndLatex

// Software Guide : BeginCodeSnippet
#include "itkCovariantVector.h"
#include "itkPointSet.h"
// Software Guide : EndCodeSnippet

int main(int, char *[])
{
  //  Software Guide : BeginLatex
  //
  //  The CovariantVector class is templated over the type used to
  //  represent the spatial coordinates and over the space dimension.  Since
  //  the PixelType is independent of the PointType, we are free to select any
  //  dimension for the covariant vectors to be used as pixel type. However, we
  //  want to illustrate here the spirit of a deformable model. It is then
  //  required for the vectors representing gradients to be of the same
  //  dimension as the points in space.
  //
  //  \index{itk::CovariantVector!Instantiation}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  const unsigned int Dimension = 3;
  typedef itk::CovariantVector< float, Dimension >    PixelType;
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  Then we use the PixelType (which are actually CovariantVectors) to
  //  instantiate the PointSet type and subsequently create a PointSet object.
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  typedef itk::PointSet< PixelType, Dimension > PointSetType;
  PointSetType::Pointer  pointSet = PointSetType::New();
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  The following code generates a circle and assigns gradient values to
  //  the points. The components of the CovariantVectors in this example are
  //  computed to represent the normals to the circle.
  //
  //  \index{itk::PointSet!SetPoint()}
  //  \index{itk::PointSet!SetPointData()}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  PointSetType::PixelType   gradient;
  PointSetType::PointType   point;

  unsigned int pointId =  0;
  const double radius = 300.0;

  for(unsigned int i=0; i<360; i++)
    {
    const double angle = i * std::atan(1.0) / 45.0;
    point[0] = radius * std::sin( angle );
    point[1] = radius * std::cos( angle );
    point[2] = 1.0;   // flat on the Z plane
    gradient[0] =  std::sin(angle);
    gradient[1] =  std::cos(angle);
    gradient[2] = 0.0;  // flat on the Z plane
    pointSet->SetPoint( pointId, point );
    pointSet->SetPointData( pointId, gradient );
    pointId++;
    }
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  We can now visit all the points and use the vector on the pixel values
  //  to apply a deformation on the points by following the gradient of the
  //  function. This is along the spirit of what a deformable model could do
  //  at each one of its iterations. To be more formal we should use the
  //  function gradients as forces and multiply them by local stress tensors
  //  in order to obtain local deformations.  The resulting deformations
  //  would finally be used to apply displacements on the points.  However,
  //  to shorten the example, we will ignore this complexity for the moment.
  //
  //  \index{itk::PointSet!PointDataIterator}
  //
  //  Software Guide : EndLatex


  // Software Guide : BeginCodeSnippet
  typedef  PointSetType::PointDataContainer::ConstIterator PointDataIterator;
  PointDataIterator pixelIterator = pointSet->GetPointData()->Begin();
  PointDataIterator pixelEnd      = pointSet->GetPointData()->End();

  typedef  PointSetType::PointsContainer::Iterator     PointIterator;
  PointIterator pointIterator = pointSet->GetPoints()->Begin();
  PointIterator pointEnd      = pointSet->GetPoints()->End();

  while( pixelIterator != pixelEnd  && pointIterator != pointEnd )
    {
    point    = pointIterator.Value();
    gradient = pixelIterator.Value();
    for(unsigned int i=0; i<Dimension; i++)
      {
      point[i] += gradient[i];
      }
    pointIterator.Value() = point;
    ++pixelIterator;
    ++pointIterator;
    }
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  The CovariantVector class does not overload the \code{+}
  //  operator with the \doxygen{Point}. In other words, CovariantVectors can
  //  not be added to points in order to get new points. Further, since we
  //  are ignoring physics in the example, we are also forced to do the
  //  illegal addition manually between the components of the gradient and
  //  the coordinates of the points.
  //
  //  Note that the absence of some basic operators on the ITK geometry classes
  //  is completely intentional with the aim of preventing the  incorrect use
  //  of the mathematical concepts they represent.
  //
  //  \index{itk::CovariantVector}
  //
  //  Software Guide : EndLatex
  //


  //  We can finally visit all the points and print out the new values.
  //
  pointIterator = pointSet->GetPoints()->Begin();
  pointEnd      = pointSet->GetPoints()->End();
  while( pointIterator != pointEnd )
    {
    std::cout << pointIterator.Value() << std::endl;
    ++pointIterator;
    }


  return EXIT_SUCCESS;
}