00001 /*========================================================================= 00002 00003 Program: Insight Segmentation & Registration Toolkit 00004 Module: $RCSfile: itkQuadEdgeMeshDiscreteMeanCurvatureEstimator.h,v $ 00005 Language: C++ 00006 Date: $Date: 2010-01-14 18:19:08 $ 00007 Version: $Revision: 1.3 $ 00008 00009 Copyright (c) Insight Software Consortium. All rights reserved. 00010 See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details. 00011 00012 This software is distributed WITHOUT ANY WARRANTY; without even 00013 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 00014 PURPOSE. See the above copyright notices for more information. 00015 00016 =========================================================================*/ 00017 00018 #ifndef __itkQuadEdgeMeshDiscreteMeanCurvatureEstimator_h 00019 #define __itkQuadEdgeMeshDiscreteMeanCurvatureEstimator_h 00020 00021 #include "itkQuadEdgeMeshDiscreteCurvatureEstimator.h" 00022 #include "itkQuadEdgeMeshParamMatrixCoefficients.h" 00023 00024 namespace itk 00025 { 00034 template< class TInputMesh, class TOutputMesh > 00035 class QuadEdgeMeshDiscreteMeanCurvatureEstimator : 00036 public QuadEdgeMeshDiscreteCurvatureEstimator< TInputMesh, TOutputMesh > 00037 { 00038 public: 00039 typedef QuadEdgeMeshDiscreteMeanCurvatureEstimator Self; 00040 typedef SmartPointer< Self > Pointer; 00041 typedef SmartPointer< const Self > ConstPointer; 00042 typedef QuadEdgeMeshDiscreteCurvatureEstimator< TInputMesh, TOutputMesh > Superclass; 00043 00044 typedef typename Superclass::InputMeshType InputMeshType; 00045 typedef typename Superclass::InputMeshPointer InputMeshPointer; 00046 00047 typedef typename Superclass::OutputMeshType OutputMeshType; 00048 typedef typename Superclass::OutputMeshPointer OutputMeshPointer; 00049 typedef typename Superclass::OutputPointsContainerPointer OutputPointsContainerPointer; 00050 typedef typename Superclass::OutputPointsContainerIterator OutputPointsContainerIterator; 00051 typedef typename Superclass::OutputPointType OutputPointType; 00052 typedef typename Superclass::OutputVectorType OutputVectorType; 00053 typedef typename Superclass::OutputCoordType OutputCoordType; 00054 typedef typename Superclass::OutputPointIdentifier OutputPointIdentifier; 00055 typedef typename Superclass::OutputCellIdentifier OutputCellIdentifier; 00056 typedef typename Superclass::OutputQEType OutputQEType; 00057 typedef typename Superclass::OutputMeshTraits OutputMeshTraits; 00058 typedef typename Superclass::OutputCurvatureType OutputCurvatureType; 00059 00060 typedef typename Superclass::TriangleType TriangleType; 00061 00063 itkTypeMacro( QuadEdgeMeshDiscreteMeanCurvatureEstimator, QuadEdgeMeshDiscreteCurvatureEstimator ); 00064 00066 itkNewMacro( Self ); 00067 00068 typedef ConformalMatrixCoefficients< OutputMeshType > CoefficientType; 00069 00070 protected: 00071 QuadEdgeMeshDiscreteMeanCurvatureEstimator() {} 00072 ~QuadEdgeMeshDiscreteMeanCurvatureEstimator() {} 00073 00074 virtual OutputCurvatureType EstimateCurvature( const OutputPointType& iP ) 00075 { 00076 OutputMeshPointer output = this->GetOutput(); 00077 00078 OutputQEType* qe = iP.GetEdge( ); 00079 00080 OutputCurvatureType oH( 0. ); 00081 00082 OutputVectorType Laplace; 00083 Laplace.Fill( 0. ); 00084 00085 OutputCurvatureType area( 0. ); 00086 OutputVectorType normal; 00087 normal.Fill( 0. ); 00088 00089 if( qe != 0 ) 00090 { 00091 if( qe != qe->GetOnext() ) 00092 { 00093 CoefficientType coefficent; 00094 00095 OutputQEType* qe_it = qe; 00096 OutputQEType* qe_it2; 00097 00098 OutputCurvatureType temp_area; 00099 OutputCoordType temp_coeff; 00100 00101 OutputPointType q0, q1; 00102 OutputVectorType face_normal; 00103 00104 do 00105 { 00106 qe_it2 = qe_it->GetOnext(); 00107 q0 = output->GetPoint( qe_it->GetDestination() ); 00108 q1 = output->GetPoint( qe_it2->GetDestination() ); 00109 00110 temp_coeff = coefficent( output, qe_it ); 00111 Laplace += temp_coeff * ( iP - q0 ); 00112 00113 temp_area = ComputeMixedArea( qe_it, qe_it2 ); 00114 area += temp_area; 00115 00116 face_normal = TriangleType::ComputeNormal( q0, iP, q1 ); 00117 normal += face_normal; 00118 00119 qe_it = qe_it2; 00120 } while( qe_it != qe ); 00121 00122 if( area < 1e-6 ) 00123 { 00124 oH = 0.; 00125 } 00126 else 00127 { 00128 if( normal.GetSquaredNorm() > 0. ) 00129 { 00130 normal.Normalize(); 00131 Laplace *= 0.25 / area; 00132 oH = Laplace * normal; 00133 } 00134 else 00135 { 00136 oH = 0.; 00137 } 00138 } 00139 } 00140 } 00141 return oH; 00142 } 00143 00144 private: 00145 QuadEdgeMeshDiscreteMeanCurvatureEstimator( const Self& ); // purposely not implemented 00146 void operator = ( const Self& ); // purposely not implemented 00147 00148 }; 00149 00150 } 00151 #endif 00152