ITK
4.1.0
Insight Segmentation and Registration Toolkit
|
00001 /*========================================================================= 00002 * 00003 * Copyright Insight Software Consortium 00004 * 00005 * Licensed under the Apache License, Version 2.0 (the "License"); 00006 * you may not use this file except in compliance with the License. 00007 * You may obtain a copy of the License at 00008 * 00009 * http://www.apache.org/licenses/LICENSE-2.0.txt 00010 * 00011 * Unless required by applicable law or agreed to in writing, software 00012 * distributed under the License is distributed on an "AS IS" BASIS, 00013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 00014 * See the License for the specific language governing permissions and 00015 * limitations under the License. 00016 * 00017 *=========================================================================*/ 00018 #ifndef __itkShapeLabelObject_h 00019 #define __itkShapeLabelObject_h 00020 00021 #include "itkLabelObject.h" 00022 #include "itkLabelMap.h" 00023 #include "itkAffineTransform.h" 00024 00025 namespace itk 00026 { 00041 template< class TLabel, unsigned int VImageDimension > 00042 class ITK_EXPORT ShapeLabelObject:public LabelObject< TLabel, VImageDimension > 00043 { 00044 public: 00046 typedef ShapeLabelObject Self; 00047 typedef LabelObject< TLabel, VImageDimension > Superclass; 00048 typedef typename Superclass::LabelObjectType LabelObjectType; 00049 typedef SmartPointer< Self > Pointer; 00050 typedef SmartPointer< const Self > ConstPointer; 00051 typedef WeakPointer< const Self > ConstWeakPointer; 00052 00054 itkNewMacro(Self); 00055 00057 itkTypeMacro(ShapeLabelObject, LabelObject); 00058 00059 typedef LabelMap< Self > LabelMapType; 00060 00061 itkStaticConstMacro(ImageDimension, unsigned int, VImageDimension); 00062 00063 typedef typename Superclass::IndexType IndexType; 00064 00065 typedef TLabel LabelType; 00066 00067 typedef typename Superclass::LineType LineType; 00068 00069 typedef typename Superclass::LengthType LengthType; 00070 00071 typedef typename Superclass::AttributeType AttributeType; 00072 00074 itkStaticConstMacro(NUMBER_OF_PIXELS, AttributeType, 100); 00075 00079 itkStaticConstMacro(PHYSICAL_SIZE, AttributeType, 101); 00080 00084 itkStaticConstMacro(CENTROID, AttributeType, 104); 00085 00086 itkStaticConstMacro(BOUNDING_BOX, AttributeType, 105); 00087 00094 itkStaticConstMacro(NUMBER_OF_PIXELS_ON_BORDER, AttributeType, 106); 00095 00103 itkStaticConstMacro(PERIMETER_ON_BORDER, AttributeType, 107); 00104 00108 itkStaticConstMacro(FERET_DIAMETER, AttributeType, 108); 00109 00111 itkStaticConstMacro(PRINCIPAL_MOMENTS, AttributeType, 109); 00112 00114 itkStaticConstMacro(PRINCIPAL_AXES, AttributeType, 110); 00115 00119 itkStaticConstMacro(ELONGATION, AttributeType, 111); 00120 00122 itkStaticConstMacro(PERIMETER, AttributeType, 112); 00123 00124 itkStaticConstMacro(ROUNDNESS, AttributeType, 113); 00125 00129 itkStaticConstMacro(EQUIVALENT_SPHERICAL_RADIUS, AttributeType, 114); 00130 00134 itkStaticConstMacro(EQUIVALENT_SPHERICAL_PERIMETER, AttributeType, 115); 00135 00139 itkStaticConstMacro(EQUIVALENT_ELLIPSOID_DIAMETER, AttributeType, 116); 00140 00141 itkStaticConstMacro(FLATNESS, AttributeType, 117); 00142 00143 itkStaticConstMacro(PERIMETER_ON_BORDER_RATIO, AttributeType, 118); 00144 00145 static AttributeType GetAttributeFromName(const std::string & s) 00146 { 00147 if ( s == "NumberOfPixels" ) 00148 { 00149 return NUMBER_OF_PIXELS; 00150 } 00151 else if ( s == "PhysicalSize" ) 00152 { 00153 return PHYSICAL_SIZE; 00154 } 00155 else if ( s == "Centroid" ) 00156 { 00157 return CENTROID; 00158 } 00159 else if ( s == "BoundingBox" ) 00160 { 00161 return BOUNDING_BOX; 00162 } 00163 else if ( s == "NumberOfPixelsOnBorder" ) 00164 { 00165 return NUMBER_OF_PIXELS_ON_BORDER; 00166 } 00167 else if ( s == "PerimeterOnBorder" ) 00168 { 00169 return PERIMETER_ON_BORDER; 00170 } 00171 else if ( s == "FeretDiameter" ) 00172 { 00173 return FERET_DIAMETER; 00174 } 00175 else if ( s == "PrincipalMoments" ) 00176 { 00177 return PRINCIPAL_MOMENTS; 00178 } 00179 else if ( s == "PrincipalAxes" ) 00180 { 00181 return PRINCIPAL_AXES; 00182 } 00183 else if ( s == "Elongation" ) 00184 { 00185 return ELONGATION; 00186 } 00187 else if ( s == "Perimeter" ) 00188 { 00189 return PERIMETER; 00190 } 00191 else if ( s == "Roundness" ) 00192 { 00193 return ROUNDNESS; 00194 } 00195 else if ( s == "EquivalentSphericalRadius" ) 00196 { 00197 return EQUIVALENT_SPHERICAL_RADIUS; 00198 } 00199 else if ( s == "EquivalentSphericalPerimeter" ) 00200 { 00201 return EQUIVALENT_SPHERICAL_PERIMETER; 00202 } 00203 else if ( s == "EquivalentEllipsoidDiameter" ) 00204 { 00205 return EQUIVALENT_ELLIPSOID_DIAMETER; 00206 } 00207 else if ( s == "Flatness" ) 00208 { 00209 return FLATNESS; 00210 } 00211 else if ( s == "PerimeterOnBorderRatio" ) 00212 { 00213 return PERIMETER_ON_BORDER_RATIO; 00214 } 00215 // can't recognize the name 00216 return Superclass::GetAttributeFromName(s); 00217 } 00218 00219 static std::string GetNameFromAttribute(const AttributeType & a) 00220 { 00221 std::string name; 00222 switch ( a ) 00223 { 00224 case NUMBER_OF_PIXELS: 00225 name = "NumberOfPixels"; 00226 break; 00227 case PHYSICAL_SIZE: 00228 name = "PhysicalSize"; 00229 break; 00230 case CENTROID: 00231 name = "Centroid"; 00232 break; 00233 case BOUNDING_BOX: 00234 name = "BoundingBox"; 00235 break; 00236 case NUMBER_OF_PIXELS_ON_BORDER: 00237 name = "NumberOfPixelsOnBorder"; 00238 break; 00239 case PERIMETER_ON_BORDER: 00240 name = "PerimeterOnBorder"; 00241 break; 00242 case FERET_DIAMETER: 00243 name = "FeretDiameter"; 00244 break; 00245 case PRINCIPAL_MOMENTS: 00246 name = "PrincipalMoments"; 00247 break; 00248 case PRINCIPAL_AXES: 00249 name = "PrincipalAxes"; 00250 break; 00251 case ELONGATION: 00252 name = "Elongation"; 00253 break; 00254 case PERIMETER: 00255 name = "Perimeter"; 00256 break; 00257 case ROUNDNESS: 00258 name = "Roundness"; 00259 break; 00260 case EQUIVALENT_SPHERICAL_RADIUS: 00261 name = "EquivalentSphericalRadius"; 00262 break; 00263 case EQUIVALENT_SPHERICAL_PERIMETER: 00264 name = "EquivalentSphericalPerimeter"; 00265 break; 00266 case EQUIVALENT_ELLIPSOID_DIAMETER: 00267 name = "EquivalentEllipsoidDiameter"; 00268 break; 00269 case FLATNESS: 00270 name = "Flatness"; 00271 break; 00272 case PERIMETER_ON_BORDER_RATIO: 00273 name = "PerimeterOnBorderRatio"; 00274 break; 00275 default: 00276 // can't recognize the name 00277 name = Superclass::GetNameFromAttribute(a); 00278 break; 00279 } 00280 return name; 00281 } 00282 00283 typedef ImageRegion< VImageDimension > RegionType; 00284 00285 typedef Point< double, VImageDimension > CentroidType; 00286 00287 typedef Matrix< double, VImageDimension, VImageDimension > MatrixType; 00288 00289 typedef Vector< double, VImageDimension > VectorType; 00290 00291 const RegionType & GetBoundingBox() const 00292 { 00293 return m_BoundingBox; 00294 } 00295 00296 void SetBoundingBox(const RegionType & v) 00297 { 00298 m_BoundingBox = v; 00299 } 00300 00301 const double & GetPhysicalSize() const 00302 { 00303 return m_PhysicalSize; 00304 } 00305 00306 void SetPhysicalSize(const double & v) 00307 { 00308 m_PhysicalSize = v; 00309 } 00310 00311 const SizeValueType & GetNumberOfPixels() const 00312 { 00313 return m_NumberOfPixels; 00314 } 00315 00316 void SetNumberOfPixels(const SizeValueType & v) 00317 { 00318 m_NumberOfPixels = v; 00319 } 00320 00321 const CentroidType & GetCentroid() const 00322 { 00323 return m_Centroid; 00324 } 00325 00326 void SetCentroid(const CentroidType & centroid) 00327 { 00328 m_Centroid = centroid; 00329 } 00330 00331 const SizeValueType & GetNumberOfPixelsOnBorder() const 00332 { 00333 return m_NumberOfPixelsOnBorder; 00334 } 00335 00336 void SetNumberOfPixelsOnBorder(const SizeValueType & v) 00337 { 00338 m_NumberOfPixelsOnBorder = v; 00339 } 00340 00341 const double & GetPerimeterOnBorder() const 00342 { 00343 return m_PerimeterOnBorder; 00344 } 00345 00346 void SetPerimeterOnBorder(const double & v) 00347 { 00348 m_PerimeterOnBorder = v; 00349 } 00350 00351 const double & GetFeretDiameter() const 00352 { 00353 return m_FeretDiameter; 00354 } 00355 00356 void SetFeretDiameter(const double & v) 00357 { 00358 m_FeretDiameter = v; 00359 } 00360 00361 const VectorType & GetPrincipalMoments() const 00362 { 00363 return m_PrincipalMoments; 00364 } 00365 00366 void SetPrincipalMoments(const VectorType & v) 00367 { 00368 m_PrincipalMoments = v; 00369 } 00370 00371 const MatrixType & GetPrincipalAxes() const 00372 { 00373 return m_PrincipalAxes; 00374 } 00375 00376 void SetPrincipalAxes(const MatrixType & v) 00377 { 00378 m_PrincipalAxes = v; 00379 } 00380 00381 const double & GetElongation() const 00382 { 00383 return m_Elongation; 00384 } 00385 00386 void SetElongation(const double & v) 00387 { 00388 m_Elongation = v; 00389 } 00390 00391 const double & GetPerimeter() const 00392 { 00393 return m_Perimeter; 00394 } 00395 00396 void SetPerimeter(const double & v) 00397 { 00398 m_Perimeter = v; 00399 } 00400 00401 const double & GetRoundness() const 00402 { 00403 return m_Roundness; 00404 } 00405 00406 void SetRoundness(const double & v) 00407 { 00408 m_Roundness = v; 00409 } 00410 00411 const double & GetEquivalentSphericalRadius() const 00412 { 00413 return m_EquivalentSphericalRadius; 00414 } 00415 00416 void SetEquivalentSphericalRadius(const double & v) 00417 { 00418 m_EquivalentSphericalRadius = v; 00419 } 00420 00421 const double & GetEquivalentSphericalPerimeter() const 00422 { 00423 return m_EquivalentSphericalPerimeter; 00424 } 00425 00426 void SetEquivalentSphericalPerimeter(const double & v) 00427 { 00428 m_EquivalentSphericalPerimeter = v; 00429 } 00430 00431 const VectorType & GetEquivalentEllipsoidDiameter() const 00432 { 00433 return m_EquivalentEllipsoidDiameter; 00434 } 00435 00436 void SetEquivalentEllipsoidDiameter(const VectorType & v) 00437 { 00438 m_EquivalentEllipsoidDiameter = v; 00439 } 00440 00441 const double & GetFlatness() const 00442 { 00443 return m_Flatness; 00444 } 00445 00446 void SetFlatness(const double & v) 00447 { 00448 m_Flatness = v; 00449 } 00450 00451 const double & GetPerimeterOnBorderRatio() const 00452 { 00453 return m_PerimeterOnBorderRatio; 00454 } 00455 00456 void SetPerimeterOnBorderRatio(const double & v) 00457 { 00458 m_PerimeterOnBorderRatio = v; 00459 } 00460 00461 // some helper methods - not really required, but really useful! 00462 00464 typedef AffineTransform< double, VImageDimension > AffineTransformType; 00465 typedef typename AffineTransformType::Pointer AffineTransformPointer; 00466 00470 AffineTransformPointer GetPrincipalAxesToPhysicalAxesTransform() const 00471 { 00472 typename AffineTransformType::MatrixType matrix; 00473 typename AffineTransformType::OffsetType offset; 00474 for ( unsigned int i = 0; i < VImageDimension; i++ ) 00475 { 00476 offset[i] = m_Centroid[i]; 00477 for ( unsigned int j = 0; j < VImageDimension; j++ ) 00478 { 00479 matrix[j][i] = m_PrincipalAxes[i][j]; // Note the transposition 00480 } 00481 } 00483 00484 AffineTransformPointer result = AffineTransformType::New(); 00485 00486 result->SetMatrix(matrix); 00487 result->SetOffset(offset); 00488 00489 return result; 00490 } 00491 00496 AffineTransformPointer GetPhysicalAxesToPrincipalAxesTransform(void) const 00497 { 00498 typename AffineTransformType::MatrixType matrix; 00499 typename AffineTransformType::OffsetType offset; 00500 for ( unsigned int i = 0; i < VImageDimension; i++ ) 00501 { 00502 offset[i] = m_Centroid[i]; 00503 for ( unsigned int j = 0; j < VImageDimension; j++ ) 00504 { 00505 matrix[j][i] = m_PrincipalAxes[i][j]; // Note the transposition 00506 } 00507 } 00509 00510 AffineTransformPointer result = AffineTransformType::New(); 00511 result->SetMatrix(matrix); 00512 result->SetOffset(offset); 00513 00514 AffineTransformPointer inverse = AffineTransformType::New(); 00515 result->GetInverse(inverse); 00516 00517 return inverse; 00518 } 00519 00520 virtual void CopyAttributesFrom(const LabelObjectType *lo) 00521 { 00522 Superclass::CopyAttributesFrom(lo); 00523 00524 // copy the data of the current type if possible 00525 const Self *src = dynamic_cast< const Self * >( lo ); 00526 if ( src == NULL ) 00527 { 00528 return; 00529 } 00530 m_BoundingBox = src->m_BoundingBox; 00531 m_NumberOfPixels = src->m_NumberOfPixels; 00532 m_PhysicalSize = src->m_PhysicalSize; 00533 m_Centroid = src->m_Centroid; 00534 m_NumberOfPixelsOnBorder = src->m_NumberOfPixelsOnBorder; 00535 m_PerimeterOnBorder = src->m_PerimeterOnBorder; 00536 m_FeretDiameter = src->m_FeretDiameter; 00537 m_PrincipalMoments = src->m_PrincipalMoments; 00538 m_PrincipalAxes = src->m_PrincipalAxes; 00539 m_Elongation = src->m_Elongation; 00540 m_Perimeter = src->m_Perimeter; 00541 m_Roundness = src->m_Roundness; 00542 m_EquivalentSphericalRadius = src->m_EquivalentSphericalRadius; 00543 m_EquivalentSphericalPerimeter = src->m_EquivalentSphericalPerimeter; 00544 m_EquivalentEllipsoidDiameter = src->m_EquivalentEllipsoidDiameter; 00545 m_Flatness = src->m_Flatness; 00546 m_PerimeterOnBorderRatio = src->m_PerimeterOnBorderRatio; 00547 } 00548 00549 protected: 00550 ShapeLabelObject() 00551 { 00552 m_NumberOfPixels = 0; 00553 m_PhysicalSize = 0; 00554 m_Centroid.Fill(0); 00555 m_NumberOfPixelsOnBorder = 0; 00556 m_PerimeterOnBorder = 0; 00557 m_FeretDiameter = 0; 00558 m_PrincipalMoments.Fill(0); 00559 m_PrincipalAxes.Fill(0); 00560 m_Elongation = 0; 00561 m_Perimeter = 0; 00562 m_Roundness = 0; 00563 m_EquivalentSphericalRadius = 0; 00564 m_EquivalentSphericalPerimeter = 0; 00565 m_EquivalentEllipsoidDiameter.Fill(0); 00566 m_Flatness = 0; 00567 m_PerimeterOnBorderRatio = 0; 00568 } 00569 00570 void PrintSelf(std::ostream & os, Indent indent) const 00571 { 00572 Superclass::PrintSelf(os, indent); 00573 00574 os << indent << "NumberOfPixels: " << m_NumberOfPixels << std::endl; 00575 os << indent << "PhysicalSize: " << m_PhysicalSize << std::endl; 00576 os << indent << "Perimeter: " << m_Perimeter << std::endl; 00577 os << indent << "NumberOfPixelsOnBorder: " << m_NumberOfPixelsOnBorder << std::endl; 00578 os << indent << "PerimeterOnBorder: " << m_PerimeterOnBorder << std::endl; 00579 os << indent << "PerimeterOnBorderRatio: " << m_PerimeterOnBorderRatio << std::endl; 00580 os << indent << "Elongation: " << m_Elongation << std::endl; 00581 os << indent << "Flatness: " << m_Flatness << std::endl; 00582 os << indent << "Roundness: " << m_Roundness << std::endl; 00583 os << indent << "Centroid: " << m_Centroid << std::endl; 00584 os << indent << "BoundingBox: "; 00585 m_BoundingBox.Print(os, indent); 00586 os << indent << "EquivalentSphericalRadius: " << m_EquivalentSphericalRadius << std::endl; 00587 os << indent << "EquivalentSphericalPerimeter: " << m_EquivalentSphericalPerimeter << std::endl; 00588 os << indent << "EquivalentEllipsoidDiameter: " << m_EquivalentEllipsoidDiameter << std::endl; 00589 os << indent << "PrincipalMoments: " << m_PrincipalMoments << std::endl; 00590 os << indent << "PrincipalAxes: " << std::endl << m_PrincipalAxes; 00591 os << indent << "FeretDiameter: " << m_FeretDiameter << std::endl; 00592 } 00593 00594 private: 00595 ShapeLabelObject(const Self &); //purposely not implemented 00596 void operator=(const Self &); //purposely not implemented 00597 00598 RegionType m_BoundingBox; 00599 SizeValueType m_NumberOfPixels; 00600 double m_PhysicalSize; 00601 CentroidType m_Centroid; 00602 SizeValueType m_NumberOfPixelsOnBorder; 00603 double m_PerimeterOnBorder; 00604 double m_FeretDiameter; 00605 VectorType m_PrincipalMoments; 00606 MatrixType m_PrincipalAxes; 00607 double m_Elongation; 00608 double m_Perimeter; 00609 double m_Roundness; 00610 double m_EquivalentSphericalRadius; 00611 double m_EquivalentSphericalPerimeter; 00612 VectorType m_EquivalentEllipsoidDiameter; 00613 double m_Flatness; 00614 double m_PerimeterOnBorderRatio; 00615 }; 00616 } // end namespace itk 00617 00618 #endif 00619