ITK  4.2.0
Insight Segmentation and Registration Toolkit
itkShapeLabelObject.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 #ifndef __itkShapeLabelObject_h
19 #define __itkShapeLabelObject_h
20 
21 #include "itkLabelObject.h"
22 #include "itkLabelMap.h"
23 #include "itkAffineTransform.h"
24 
25 namespace itk
26 {
41 template< class TLabel, unsigned int VImageDimension >
42 class ITK_EXPORT ShapeLabelObject:public LabelObject< TLabel, VImageDimension >
43 {
44 public:
48  typedef typename Superclass::LabelObjectType LabelObjectType;
52 
54  itkNewMacro(Self);
55 
57  itkTypeMacro(ShapeLabelObject, LabelObject);
58 
60 
61  itkStaticConstMacro(ImageDimension, unsigned int, VImageDimension);
62 
63  typedef typename Superclass::IndexType IndexType;
64 
65  typedef TLabel LabelType;
66 
67  typedef typename Superclass::LineType LineType;
68 
69  typedef typename Superclass::LengthType LengthType;
70 
72 
74  itkStaticConstMacro(NUMBER_OF_PIXELS, AttributeType, 100);
75 
79  itkStaticConstMacro(PHYSICAL_SIZE, AttributeType, 101);
80 
84  itkStaticConstMacro(CENTROID, AttributeType, 104);
85 
86  itkStaticConstMacro(BOUNDING_BOX, AttributeType, 105);
87 
94  itkStaticConstMacro(NUMBER_OF_PIXELS_ON_BORDER, AttributeType, 106);
95 
103  itkStaticConstMacro(PERIMETER_ON_BORDER, AttributeType, 107);
104 
108  itkStaticConstMacro(FERET_DIAMETER, AttributeType, 108);
109 
111  itkStaticConstMacro(PRINCIPAL_MOMENTS, AttributeType, 109);
112 
114  itkStaticConstMacro(PRINCIPAL_AXES, AttributeType, 110);
115 
119  itkStaticConstMacro(ELONGATION, AttributeType, 111);
120 
122  itkStaticConstMacro(PERIMETER, AttributeType, 112);
123 
124  itkStaticConstMacro(ROUNDNESS, AttributeType, 113);
125 
129  itkStaticConstMacro(EQUIVALENT_SPHERICAL_RADIUS, AttributeType, 114);
130 
134  itkStaticConstMacro(EQUIVALENT_SPHERICAL_PERIMETER, AttributeType, 115);
135 
139  itkStaticConstMacro(EQUIVALENT_ELLIPSOID_DIAMETER, AttributeType, 116);
140 
141  itkStaticConstMacro(FLATNESS, AttributeType, 117);
142 
143  itkStaticConstMacro(PERIMETER_ON_BORDER_RATIO, AttributeType, 118);
144 
145  static AttributeType GetAttributeFromName(const std::string & s)
146  {
147  if ( s == "NumberOfPixels" )
148  {
149  return NUMBER_OF_PIXELS;
150  }
151  else if ( s == "PhysicalSize" )
152  {
153  return PHYSICAL_SIZE;
154  }
155  else if ( s == "Centroid" )
156  {
157  return CENTROID;
158  }
159  else if ( s == "BoundingBox" )
160  {
161  return BOUNDING_BOX;
162  }
163  else if ( s == "NumberOfPixelsOnBorder" )
164  {
165  return NUMBER_OF_PIXELS_ON_BORDER;
166  }
167  else if ( s == "PerimeterOnBorder" )
168  {
169  return PERIMETER_ON_BORDER;
170  }
171  else if ( s == "FeretDiameter" )
172  {
173  return FERET_DIAMETER;
174  }
175  else if ( s == "PrincipalMoments" )
176  {
177  return PRINCIPAL_MOMENTS;
178  }
179  else if ( s == "PrincipalAxes" )
180  {
181  return PRINCIPAL_AXES;
182  }
183  else if ( s == "Elongation" )
184  {
185  return ELONGATION;
186  }
187  else if ( s == "Perimeter" )
188  {
189  return PERIMETER;
190  }
191  else if ( s == "Roundness" )
192  {
193  return ROUNDNESS;
194  }
195  else if ( s == "EquivalentSphericalRadius" )
196  {
197  return EQUIVALENT_SPHERICAL_RADIUS;
198  }
199  else if ( s == "EquivalentSphericalPerimeter" )
200  {
201  return EQUIVALENT_SPHERICAL_PERIMETER;
202  }
203  else if ( s == "EquivalentEllipsoidDiameter" )
204  {
205  return EQUIVALENT_ELLIPSOID_DIAMETER;
206  }
207  else if ( s == "Flatness" )
208  {
209  return FLATNESS;
210  }
211  else if ( s == "PerimeterOnBorderRatio" )
212  {
213  return PERIMETER_ON_BORDER_RATIO;
214  }
215  // can't recognize the name
216  return Superclass::GetAttributeFromName(s);
217  }
218 
219  static std::string GetNameFromAttribute(const AttributeType & a)
220  {
221  std::string name;
222  switch ( a )
223  {
224  case NUMBER_OF_PIXELS:
225  name = "NumberOfPixels";
226  break;
227  case PHYSICAL_SIZE:
228  name = "PhysicalSize";
229  break;
230  case CENTROID:
231  name = "Centroid";
232  break;
233  case BOUNDING_BOX:
234  name = "BoundingBox";
235  break;
236  case NUMBER_OF_PIXELS_ON_BORDER:
237  name = "NumberOfPixelsOnBorder";
238  break;
239  case PERIMETER_ON_BORDER:
240  name = "PerimeterOnBorder";
241  break;
242  case FERET_DIAMETER:
243  name = "FeretDiameter";
244  break;
245  case PRINCIPAL_MOMENTS:
246  name = "PrincipalMoments";
247  break;
248  case PRINCIPAL_AXES:
249  name = "PrincipalAxes";
250  break;
251  case ELONGATION:
252  name = "Elongation";
253  break;
254  case PERIMETER:
255  name = "Perimeter";
256  break;
257  case ROUNDNESS:
258  name = "Roundness";
259  break;
260  case EQUIVALENT_SPHERICAL_RADIUS:
261  name = "EquivalentSphericalRadius";
262  break;
263  case EQUIVALENT_SPHERICAL_PERIMETER:
264  name = "EquivalentSphericalPerimeter";
265  break;
266  case EQUIVALENT_ELLIPSOID_DIAMETER:
267  name = "EquivalentEllipsoidDiameter";
268  break;
269  case FLATNESS:
270  name = "Flatness";
271  break;
272  case PERIMETER_ON_BORDER_RATIO:
273  name = "PerimeterOnBorderRatio";
274  break;
275  default:
276  // can't recognize the name
277  name = Superclass::GetNameFromAttribute(a);
278  break;
279  }
280  return name;
281  }
282 
284 
286 
288 
290 
291  const RegionType & GetBoundingBox() const
292  {
293  return m_BoundingBox;
294  }
295 
296  void SetBoundingBox(const RegionType & v)
297  {
298  m_BoundingBox = v;
299  }
300 
301  const double & GetPhysicalSize() const
302  {
303  return m_PhysicalSize;
304  }
305 
306  void SetPhysicalSize(const double & v)
307  {
308  m_PhysicalSize = v;
309  }
310 
311  const SizeValueType & GetNumberOfPixels() const
312  {
313  return m_NumberOfPixels;
314  }
315 
316  void SetNumberOfPixels(const SizeValueType & v)
317  {
318  m_NumberOfPixels = v;
319  }
320 
321  const CentroidType & GetCentroid() const
322  {
323  return m_Centroid;
324  }
325 
326  void SetCentroid(const CentroidType & centroid)
327  {
328  m_Centroid = centroid;
329  }
330 
331  const SizeValueType & GetNumberOfPixelsOnBorder() const
332  {
333  return m_NumberOfPixelsOnBorder;
334  }
335 
336  void SetNumberOfPixelsOnBorder(const SizeValueType & v)
337  {
338  m_NumberOfPixelsOnBorder = v;
339  }
340 
341  const double & GetPerimeterOnBorder() const
342  {
343  return m_PerimeterOnBorder;
344  }
345 
346  void SetPerimeterOnBorder(const double & v)
347  {
348  m_PerimeterOnBorder = v;
349  }
350 
351  const double & GetFeretDiameter() const
352  {
353  return m_FeretDiameter;
354  }
355 
356  void SetFeretDiameter(const double & v)
357  {
358  m_FeretDiameter = v;
359  }
360 
361  const VectorType & GetPrincipalMoments() const
362  {
363  return m_PrincipalMoments;
364  }
365 
366  void SetPrincipalMoments(const VectorType & v)
367  {
368  m_PrincipalMoments = v;
369  }
370 
371  const MatrixType & GetPrincipalAxes() const
372  {
373  return m_PrincipalAxes;
374  }
375 
376  void SetPrincipalAxes(const MatrixType & v)
377  {
378  m_PrincipalAxes = v;
379  }
380 
381  const double & GetElongation() const
382  {
383  return m_Elongation;
384  }
385 
386  void SetElongation(const double & v)
387  {
388  m_Elongation = v;
389  }
390 
391  const double & GetPerimeter() const
392  {
393  return m_Perimeter;
394  }
395 
396  void SetPerimeter(const double & v)
397  {
398  m_Perimeter = v;
399  }
400 
401  const double & GetRoundness() const
402  {
403  return m_Roundness;
404  }
405 
406  void SetRoundness(const double & v)
407  {
408  m_Roundness = v;
409  }
410 
411  const double & GetEquivalentSphericalRadius() const
412  {
413  return m_EquivalentSphericalRadius;
414  }
415 
416  void SetEquivalentSphericalRadius(const double & v)
417  {
418  m_EquivalentSphericalRadius = v;
419  }
420 
421  const double & GetEquivalentSphericalPerimeter() const
422  {
423  return m_EquivalentSphericalPerimeter;
424  }
425 
426  void SetEquivalentSphericalPerimeter(const double & v)
427  {
428  m_EquivalentSphericalPerimeter = v;
429  }
430 
431  const VectorType & GetEquivalentEllipsoidDiameter() const
432  {
433  return m_EquivalentEllipsoidDiameter;
434  }
435 
436  void SetEquivalentEllipsoidDiameter(const VectorType & v)
437  {
438  m_EquivalentEllipsoidDiameter = v;
439  }
440 
441  const double & GetFlatness() const
442  {
443  return m_Flatness;
444  }
445 
446  void SetFlatness(const double & v)
447  {
448  m_Flatness = v;
449  }
450 
451  const double & GetPerimeterOnBorderRatio() const
452  {
453  return m_PerimeterOnBorderRatio;
454  }
455 
456  void SetPerimeterOnBorderRatio(const double & v)
457  {
458  m_PerimeterOnBorderRatio = v;
459  }
460 
461  // some helper methods - not really required, but really useful!
462 
466 
470  AffineTransformPointer GetPrincipalAxesToPhysicalAxesTransform() const
471  {
472  typename AffineTransformType::MatrixType matrix;
473  typename AffineTransformType::OffsetType offset;
474  for ( unsigned int i = 0; i < VImageDimension; i++ )
475  {
476  offset[i] = m_Centroid[i];
477  for ( unsigned int j = 0; j < VImageDimension; j++ )
478  {
479  matrix[j][i] = m_PrincipalAxes[i][j]; // Note the transposition
480  }
481  }
483 
484  AffineTransformPointer result = AffineTransformType::New();
485 
486  result->SetMatrix(matrix);
487  result->SetOffset(offset);
488 
489  return result;
490  }
491 
496  AffineTransformPointer GetPhysicalAxesToPrincipalAxesTransform(void) const
497  {
498  typename AffineTransformType::MatrixType matrix;
499  typename AffineTransformType::OffsetType offset;
500  for ( unsigned int i = 0; i < VImageDimension; i++ )
501  {
502  offset[i] = m_Centroid[i];
503  for ( unsigned int j = 0; j < VImageDimension; j++ )
504  {
505  matrix[j][i] = m_PrincipalAxes[i][j]; // Note the transposition
506  }
507  }
509 
510  AffineTransformPointer result = AffineTransformType::New();
511  result->SetMatrix(matrix);
512  result->SetOffset(offset);
513 
514  AffineTransformPointer inverse = AffineTransformType::New();
515  result->GetInverse(inverse);
516 
517  return inverse;
518  }
519 
520  virtual void CopyAttributesFrom(const LabelObjectType *lo)
521  {
522  Superclass::CopyAttributesFrom(lo);
523 
524  // copy the data of the current type if possible
525  const Self *src = dynamic_cast< const Self * >( lo );
526  if ( src == NULL )
527  {
528  return;
529  }
530  m_BoundingBox = src->m_BoundingBox;
531  m_NumberOfPixels = src->m_NumberOfPixels;
532  m_PhysicalSize = src->m_PhysicalSize;
533  m_Centroid = src->m_Centroid;
534  m_NumberOfPixelsOnBorder = src->m_NumberOfPixelsOnBorder;
535  m_PerimeterOnBorder = src->m_PerimeterOnBorder;
536  m_FeretDiameter = src->m_FeretDiameter;
537  m_PrincipalMoments = src->m_PrincipalMoments;
538  m_PrincipalAxes = src->m_PrincipalAxes;
539  m_Elongation = src->m_Elongation;
540  m_Perimeter = src->m_Perimeter;
541  m_Roundness = src->m_Roundness;
542  m_EquivalentSphericalRadius = src->m_EquivalentSphericalRadius;
543  m_EquivalentSphericalPerimeter = src->m_EquivalentSphericalPerimeter;
544  m_EquivalentEllipsoidDiameter = src->m_EquivalentEllipsoidDiameter;
545  m_Flatness = src->m_Flatness;
546  m_PerimeterOnBorderRatio = src->m_PerimeterOnBorderRatio;
547  }
548 
549 protected:
551  {
552  m_NumberOfPixels = 0;
553  m_PhysicalSize = 0;
554  m_Centroid.Fill(0);
555  m_NumberOfPixelsOnBorder = 0;
556  m_PerimeterOnBorder = 0;
557  m_FeretDiameter = 0;
558  m_PrincipalMoments.Fill(0);
559  m_PrincipalAxes.Fill(0);
560  m_Elongation = 0;
561  m_Perimeter = 0;
562  m_Roundness = 0;
563  m_EquivalentSphericalRadius = 0;
564  m_EquivalentSphericalPerimeter = 0;
565  m_EquivalentEllipsoidDiameter.Fill(0);
566  m_Flatness = 0;
567  m_PerimeterOnBorderRatio = 0;
568  }
569 
570  void PrintSelf(std::ostream & os, Indent indent) const
571  {
572  Superclass::PrintSelf(os, indent);
573 
574  os << indent << "NumberOfPixels: " << m_NumberOfPixels << std::endl;
575  os << indent << "PhysicalSize: " << m_PhysicalSize << std::endl;
576  os << indent << "Perimeter: " << m_Perimeter << std::endl;
577  os << indent << "NumberOfPixelsOnBorder: " << m_NumberOfPixelsOnBorder << std::endl;
578  os << indent << "PerimeterOnBorder: " << m_PerimeterOnBorder << std::endl;
579  os << indent << "PerimeterOnBorderRatio: " << m_PerimeterOnBorderRatio << std::endl;
580  os << indent << "Elongation: " << m_Elongation << std::endl;
581  os << indent << "Flatness: " << m_Flatness << std::endl;
582  os << indent << "Roundness: " << m_Roundness << std::endl;
583  os << indent << "Centroid: " << m_Centroid << std::endl;
584  os << indent << "BoundingBox: ";
585  m_BoundingBox.Print(os, indent);
586  os << indent << "EquivalentSphericalRadius: " << m_EquivalentSphericalRadius << std::endl;
587  os << indent << "EquivalentSphericalPerimeter: " << m_EquivalentSphericalPerimeter << std::endl;
588  os << indent << "EquivalentEllipsoidDiameter: " << m_EquivalentEllipsoidDiameter << std::endl;
589  os << indent << "PrincipalMoments: " << m_PrincipalMoments << std::endl;
590  os << indent << "PrincipalAxes: " << std::endl << m_PrincipalAxes;
591  os << indent << "FeretDiameter: " << m_FeretDiameter << std::endl;
592  }
593 
594 private:
595  ShapeLabelObject(const Self &); //purposely not implemented
596  void operator=(const Self &); //purposely not implemented
597 
607  double m_Elongation;
608  double m_Perimeter;
609  double m_Roundness;
613  double m_Flatness;
615 };
616 } // end namespace itk
617 
618 #endif
619