KWStyle - itkPoint.txx
 
Matrix View
Description
1 /*=========================================================================
2
3
4 Program:   Insight Segmentation & Registration Toolkit
5 HRD   Module:    $RCSfile: itkPoint.txx.html,v $
6 HRD   Language:  C++
7 HRD   Date:      $Date: 2006/01/17 19:15:43 $
8   Version:   $Revision: 1.4 $
9
10 HRD   Copyright (c) Insight Software Consortium. All rights reserved.
11   See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
12
13 HRD      This software is distributed WITHOUT ANY WARRANTY; without even 
14 HRD      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 
15      PURPOSE.  See the above copyright notices for more information.
16 HRD
17 =========================================================================*/
18 DEF #ifndef _itkPoint_txx
19 DEF #define _itkPoint_txx
20 #include "itkPoint.h" 
21 #include <vnl/vnl_math.h>
22 #include "itkObject.h"
23
24
25 namespace itk
26 {
27
28
29 /*
30 IND ** Assignment Operator
31 IND **/
32 template<class T, unsigned int TPointDimension>
33 Point<T, TPointDimension>&
34 Point<T, TPointDimension>
35 ::operator= (const Self& r)
36 {
37   BaseArray::operator=(r);
38   return *this;
39 }
40
41
42 /*
43 IND ** Assignemt from a plain array
44 IND **/
45 template<class T, unsigned int TPointDimension>
46 Point<T, TPointDimension>&
47 Point<T, TPointDimension>
48 ::operator= (const ValueType r[TPointDimension])
49 {
50   BaseArray::operator=(r);
51   return *this;
52 }
53
54
55 /*
56 IND ** In place increment by a vector
57 IND **/
58 template<class T, unsigned int TPointDimension>
59 const Point<T, TPointDimension> &
60 Point<T, TPointDimension>
61 ::operator+=( const VectorType & vec )
62 {
63   for( unsigned int i=0; i<TPointDimension; i++) 
64     {
65     (*this)[i] += vec[i];
66     }
67   return *this;
68 }
69
70  
71 /*
72  * In place subtract a vector
73  */
74 template<class T, unsigned int TPointDimension>
75 const Point<T, TPointDimension> &
76 Point<T, TPointDimension>
77 ::operator-=( const VectorType & vec )
78 {
79   for( unsigned int i=0; i<TPointDimension; i++) 
80     {
81     (*this)[i] -= vec[i];
82     }
83   return *this;
84 }
85
86  
87
88 /*
89  * Add operator with a vector
90  */
91 template<class T, unsigned int TPointDimension>
92 Point<T, TPointDimension> 
93 Point<T, TPointDimension>
94 ::operator+( const VectorType & vec ) const
95 {
96   Self result;
97   for( unsigned int i=0; i<TPointDimension; i++) 
98     {
99     result[i] = (*this)[i] + vec[i];
100     }
101   return result;
102 }
103
104
105 EML
106 /*
107  * Subtract a vector, return a point
108  */
109 template<class T, unsigned int TPointDimension>
110 Point<T, TPointDimension> 
111 Point<T, TPointDimension>
112 ::operator-( const VectorType & vec )  const
113 {
114   Self result;
115   for( unsigned int i=0; i<TPointDimension; i++) 
116     {
117     result[i] = (*this)[i] - vec[i];
118     }
119   return result;
120 }
121
122
123 EML
124 /*
125  * Difference between two points
126  */
127 template<class T, unsigned int TPointDimension>
128 Vector<T, TPointDimension> 
129 Point<T, TPointDimension>
130 ::operator-( const Self & pnt )  const
131 {
132   VectorType result;
133   for( unsigned int i=0; i<TPointDimension; i++) 
134     {
135     result[i] = (*this)[i] - pnt[i];
136     }
137   return result;
138 }
139
140 /*
141  * Return a vnl_vector_ref
142  */
143 template<class T, unsigned int TPointDimension >
144 vnl_vector_ref< T >
145 Point<T, TPointDimension>
146 ::GetVnlVector( void ) 
147 {
148   return vnl_vector_ref< T >( TPointDimension, this->GetDataPointer());
149 }
150
151 /**
152  * Return a vnl_vector const
153  */
154 template<class T, unsigned int TPointDimension>
155 vnl_vector< T >
156 Point<T, TPointDimension>
157 ::GetVnlVector( void ) const 
158 {
159   // Return a vector_ref<>.  This will be automatically converted to a
160   // vnl_vector<>.  We have to use a const_cast<> which would normally
161   // be prohibited in a const method, but it is safe to do here
162   // because the cast to vnl_vector<> will ultimately copy the data.
163   return vnl_vector_ref<T>( TPointDimension,
164                             const_cast<T*>(this->GetDataPointer()));
165 }
166
167 /*
168 IND ** Return a vnl_vector_ref
169 IND **/
170 template<class T, unsigned int TPointDimension >
171 vnl_vector_ref< T >
172 Point<T, TPointDimension>
173 ::Get_vnl_vector( void ) 
174 {
175   return vnl_vector_ref< T >( TPointDimension, this->GetDataPointer());
176 }
177
178 /**
179  * Return a vnl_vector const
180  */
181 template<class T, unsigned int TPointDimension>
182 vnl_vector< T >
183 Point<T, TPointDimension>
184 ::Get_vnl_vector( void ) const 
185 {
186   // Return a vector_ref<>.  This will be automatically converted to a
187   // vnl_vector<>.  We have to use a const_cast<> which would normally
188   // be prohibited in a const method, but it is safe to do here
189   // because the cast to vnl_vector<> will ultimately copy the data.
190   return vnl_vector_ref<T>( TPointDimension,
191                             const_cast<T*>(this->GetDataPointer()));
192 }
193
194 template<class T, unsigned int TPointDimension >
195 typename Point<T, TPointDimension>::VectorType
196 Point<T, TPointDimension>
197 ::GetVectorFromOrigin() const
198 {
199   // VectorType knows how to construct from ValueType*.
200   return &(*this)[0];
201 }
202
203
204 EML
205 /*
206 IND ** Set the point to the median point of the two arguments
207 IND **/
208 template<class T, unsigned int TPointDimension>
209 void
210 Point<T, TPointDimension>
211 ::SetToMidPoint( const Self & A, const Self & B )  
212 {
213   for( unsigned int i=0; i<TPointDimension; i++) 
214     {
215     (*this)[i] = ( A[i] + B[i] ) /2.0;
216     }
217 }
218
219
220 EML
221 EML
222 EML
223 /*
224  * Set the point to the barycentric combination of two points
225  */
226 template<class T, unsigned int TPointDimension>
227 void
228 Point<T, TPointDimension>
229 ::SetToBarycentricCombination( const Self & A,
230                                const Self & B,
231                                double alpha )  
232 {
233   const double wa = alpha;
234   const double wb = 1.0 - alpha;
235   for( unsigned int i=0; i<TPointDimension; i++) 
236     {
237     (*this)[i] =  wa * A[i] + wb * B[i];
238     }
239 }
240
241
242 EML
243 /*
244  * Set the point to the barycentric combination of three points
245  */
246 template<class T, unsigned int TPointDimension>
247 void
248 Point<T, TPointDimension>
249 ::SetToBarycentricCombination( const Self & A,
250                                const Self & B,
251                                const Self & C,
252                                double weightForA, 
253                                double weightForB  )
254 {
255   const double weightForC = 1.0 - weightForA - weightForB;
256   for( unsigned int i=0; i<TPointDimension; i++) 
257     {
258     (*this)[i] =  weightForA * A[i] + weightForB * B[i] + weightForC * C[i];
259     }
260 }
261
262
263 EML
264 /*
265  * Set the point to the barycentric combination of N points
266  */
267 template<class T, unsigned int TPointDimension>
268 void
269 Point<T, TPointDimension>
270 ::SetToBarycentricCombination( const Self * P,
271                                const double * weights, unsigned int N )
272 {
273   Fill( NumericTraits<T>::Zero ); // put this point to null
274   double weightSum = 0.0;
275   for( unsigned int j=0; j<N-1; j++) 
276     {
277     const double weight = weights[j];
278     weightSum += weight;
279     for( unsigned int i=0; i<TPointDimension; i++) 
280       {
281       (*this)[i] +=  weight * P[j][i];
282       }
283     }
284   const double weight = ( 1.0 - weightSum );
285   for( unsigned int i=0; i<TPointDimension; i++) 
286     {
287     (*this)[i] +=  weight * P[N-1][i];
288     }
289
290 }
291
292
293 EML
294 EML
295 /*
296  * Set the point to the barycentric combination of N points in a Container
297  */
298 template< class TPointContainer, class TWeightContainer >
299 typename BarycentricCombination<TPointContainer,TWeightContainer>
300 ::PointType 
301 BarycentricCombination<TPointContainer,TWeightContainer>
302 ::Evaluate( 
303   const PointContainerPointer & points, 
304   const WeightContainerType & weights ) 
305 {
306   
307   typedef typename TPointContainer::Element   PointType;
308   typedef typename PointType::ValueType       ValueType;
309   PointType barycentre;
310   barycentre.Fill( NumericTraits< ValueType >::Zero ); // set to null
311
312   typename TPointContainer::Iterator point = points->Begin();
313   typename TPointContainer::Iterator final = points->End();
314   typename TPointContainer::Iterator last  = final;
315   last--; // move to the (N)th point
316   
317   double weightSum = 0.0;
318
319   const unsigned int PointDimension = PointType::PointDimension;
320   unsigned int j = 0;
321
322   while( point != last )
323     {
324     const double weight = weights[j++];
325     weightSum += weight;
326     for( unsigned int i=0; i<PointDimension; i++) 
327       {
328       barycentre[i] +=  weight * (point->Value())[i];
329       }
330     point++;
331     }
332
333   // Compute directly the last one
334   // to make sure that the weights sum 1
335   const double weight = ( 1.0 - weightSum );
336   for( unsigned int i=0; i<PointDimension; i++) 
337     {
338     barycentre[i] +=  weight * (last->Value())[i];
339     }
340
341   return barycentre;
342
343 }
344
345
346 EML
347 EML
348 EML
349 /*
350 IND ** Print content to an ostream
351 IND **/
352 template<class T, unsigned int TPointDimension>
353 std::ostream &
354 operator<<(std::ostream& os,const Point<T,TPointDimension> & vct ) 
355 {
356   os << "[";
357   if ( TPointDimension == 1)
358     {
359     os << vct[0];
360     }
361   else
362     {
363     for( unsigned int i=0; i+1<TPointDimension; i++)
364       {
365       os <<  vct[i] << ", ";
366       }
367     os <<  vct[TPointDimension-1];
368     }
369   os << "]";
370   return os;
371 }
372
373
374 EML
375 /*
376  * Read content from an istream
377  */
378 template<class T, unsigned int TPointDimension>
379 std::istream &
380 operator>>(std::istream& is, Point<T,TPointDimension> & vct ) 
381 {
382   for( unsigned int i=0; i<TPointDimension; i++)
383     {
384     is >>  vct[i];
385     }
386   return is;
387 }
388
389
390 EML
391 } // end namespace itk
392
393
394 #endif
395
Generated by KWStyle 1.0b on Tuesday January,17 at 02:14:39PM
© Kitware Inc.