Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #ifndef __itkSimplexMeshAdaptTopologyFilter_h
00018 #define __itkSimplexMeshAdaptTopologyFilter_h
00019
00020 #include "itkMesh.h"
00021 #include "itkPolygonCell.h"
00022 #include "itkMapContainer.h"
00023 #include "itkCellInterfaceVisitor.h"
00024
00025 #include "itkSimplexMesh.h"
00026 #include "itkSimplexMeshGeometry.h"
00027 #include "itkMeshToMeshFilter.h"
00028 #include "itkVectorContainer.h"
00029 #include <itkCovariantVector.h>
00030
00031 #include <vxl_version.h>
00032 #if VXL_VERSION_DATE_FULL > 20040406
00033 # include <vnl/vnl_cross.h>
00034 # define itk_cross_3d vnl_cross_3d
00035 #else
00036 # define itk_cross_3d cross_3d
00037 #endif
00038
00039 namespace itk
00040 {
00051 template <class TInputMesh, class TOutputMesh>
00052 class ITK_EXPORT SimplexMeshAdaptTopologyFilter : public MeshToMeshFilter<TInputMesh, TOutputMesh>
00053 {
00054
00055 public:
00057 typedef SimplexMeshAdaptTopologyFilter Self;
00058
00060 typedef MeshToMeshFilter<TInputMesh, TOutputMesh> Superclass;
00061
00063 typedef SmartPointer<Self> Pointer;
00064
00066 typedef SmartPointer<const Self> ConstPointer;
00067
00069 itkNewMacro(Self);
00070
00072 itkTypeMacro(SimplexMeshAdaptTopologyFilter,MeshToMeshFilter);
00073
00074 typedef TInputMesh InputMeshType;
00075 typedef typename InputMeshType::Pointer InputMeshPointer;
00076 typedef typename InputMeshType::PointType InputPointType;
00077 typedef typename InputMeshType::VectorType InputVectorType;
00078 typedef typename InputMeshType::PixelType InputPixelType;
00079 typedef typename InputMeshType::MeshTraits::CellTraits InputCellTraitsType;
00080 typedef typename InputMeshType::CellType InputCellType;
00081 typedef typename InputCellType::PointIdIterator InputCellPointIdIterator;
00082 typedef typename InputCellType::CellAutoPointer InputCellAutoPointer;
00083 typedef typename InputMeshType::CellAutoPointer CellAutoPointer;
00084 typedef itk::PolygonCell<InputCellType> InputPolygonType;
00085 typedef typename InputPolygonType::PointIdIterator InputPolygonPointIdIterator;
00086 typedef CovariantVector<
00087 typename InputVectorType::ValueType, 3 > CovariantVectorType;
00088 typedef TOutputMesh OutputMeshType;
00089 typedef typename OutputMeshType::Pointer OutputMeshPointer;
00090 typedef typename OutputMeshType::CellType OutputCellType;
00091 typedef itk::PolygonCell<OutputCellType> OutputPolygonType;
00092
00093 typedef itk::MapContainer<unsigned long, double> DoubleValueMapType;
00094 typedef typename DoubleValueMapType::Iterator DoubleContainerIterator;
00095
00096
00103 class SimplexCellVisitor
00104 {
00105
00106 public:
00107 InputMeshPointer mesh;
00108 double totalArea;
00109 double totalCurvature;
00110 double minCellSize;
00111 double maxCellSize;
00112 DoubleValueMapType::Pointer areaMap;
00113 DoubleValueMapType::Pointer curvatureMap;
00114
00115 double minCurvature;
00116 double maxCurvature;
00117
00118 SimplexCellVisitor()
00119 {
00120 areaMap = DoubleValueMapType::New();
00121 curvatureMap = DoubleValueMapType::New();
00122 totalArea = 0;
00123 totalCurvature = 0;
00124 minCellSize = NumericTraits<double>::max();
00125 maxCellSize = 0;
00126 minCurvature = NumericTraits<double>::max();
00127 maxCurvature = 0;
00128 }
00129
00133 void Visit(unsigned long cellId, InputPolygonType * poly)
00134 {
00135 typename InputPolygonType::PointIdIterator it = poly->PointIdsBegin();
00136
00137 double meanCurvature = 0;
00138 unsigned long refPoint = *it;
00139 double val = mesh->GetMeanCurvature(*it++);
00140 meanCurvature += vcl_abs(val);
00141
00142 unsigned long id1 = *it;
00143 val = mesh->GetMeanCurvature(*it++);
00144 meanCurvature += vcl_abs(val);
00145
00146 unsigned long id2;
00147
00148 double area = 0;
00149
00150 int cnt = 0;
00151
00152 while ( it != poly->PointIdsEnd() )
00153 {
00154 id2 = *it;
00155 area += ComputeArea(refPoint,id1,id2);
00156 id1 = id2;
00157 val = mesh->GetMeanCurvature(*it);
00158 meanCurvature += vcl_abs(val);
00159 cnt++;
00160 it++;
00161 }
00162
00163 meanCurvature /= (double)cnt;
00164 totalArea += area;
00165 totalCurvature += meanCurvature;
00166
00167 areaMap->InsertElement(cellId, area);
00168 curvatureMap->InsertElement(cellId, meanCurvature);
00169
00170 if (area > maxCellSize ) maxCellSize = area;
00171 if (area < minCellSize ) minCellSize = area;
00172 if (meanCurvature > maxCurvature ) maxCurvature = meanCurvature;
00173 if (meanCurvature < minCurvature ) minCurvature = meanCurvature;
00174 }
00175
00176 double ComputeArea(unsigned long p1,unsigned long p2, unsigned long p3)
00177 {
00178 InputPointType v1,v2,v3;
00179 v1.Fill(0);
00180 v2.Fill(0);
00181 v3.Fill(0);
00182
00183 mesh->GetPoint(p1, &v1);
00184 mesh->GetPoint(p2, &v2);
00185 mesh->GetPoint(p3, &v3);
00186 return vcl_abs (itk_cross_3d((v2-v1).GetVnlVector(), (v3-v1).GetVnlVector()).two_norm() /2.0);
00187 }
00188
00189 DoubleValueMapType::Pointer GetAreaMap()
00190 {
00191 return areaMap;
00192 }
00193
00194 DoubleValueMapType::Pointer GetCurvatureMap()
00195 {
00196 return curvatureMap;
00197 }
00198
00199 double GetTotalMeshArea()
00200 {
00201 return totalArea;
00202 }
00203
00204 double GetTotalMeanCurvature()
00205 {
00206 return totalCurvature/(curvatureMap->Size());
00207 }
00208
00209 double GetMaximumCellSize()
00210 {
00211 return maxCellSize;
00212 }
00213
00214 double GetMinimumCellSize()
00215 {
00216 return minCellSize;
00217 }
00218
00219 double GetMaximumCurvature()
00220 {
00221 return maxCurvature;
00222 }
00223
00224 double GetMinimumCurvature()
00225 {
00226 return minCurvature;
00227 }
00228 };
00229
00230
00231 typedef itk::CellInterfaceVisitorImplementation<InputPixelType,
00232 InputCellTraitsType,
00233 InputPolygonType,
00234 SimplexCellVisitor>
00235 SimplexVisitorInterfaceType;
00236
00237 typedef typename SimplexVisitorInterfaceType::Pointer SimplexVisitorInterfacePointer;
00238 typedef typename InputCellType::MultiVisitor CellMultiVisitorType;
00239 typedef typename CellMultiVisitorType::Pointer CellMultiVisitorPointer;
00240
00241
00242 itkSetMacro(Threshold, double);
00243 itkGetConstMacro(Threshold, double);
00244
00245 itkSetMacro(SelectionMethod, int);
00246 itkGetConstMacro(SelectionMethod, int);
00247
00248 itkGetConstMacro(ModifiedCount, int);
00249
00250
00251 protected:
00252
00253 SimplexMeshAdaptTopologyFilter();
00254 ~SimplexMeshAdaptTopologyFilter();
00255 SimplexMeshAdaptTopologyFilter(const Self&) {}
00256
00257 void operator=(const Self&) {}
00258
00259 void PrintSelf(std::ostream& os, Indent indent) const;
00260
00261 virtual void GenerateData();
00262
00263
00267 void Initialize();
00268
00274 void ComputeCellParameters();
00275
00279 void InsertNewCells();
00280
00286 void ModifyNeighborCells(unsigned long id1, unsigned long id2, unsigned long insertPointId);
00287
00291 InputPointType ComputeCellCenter(InputCellAutoPointer &simplexCell);
00292
00296 unsigned long m_IdOffset;
00297
00302 double m_Threshold;
00303
00307 int m_SelectionMethod;
00308
00313 int m_ModifiedCount;
00314
00319 OutputMeshPointer m_Output;
00320
00321 InputCellAutoPointer m_NewSimplexCellPointer;
00322
00323 };
00324
00325 }
00326
00327 #ifndef ITK_MANUAL_INSTANTIATION
00328 #include "itkSimplexMeshAdaptTopologyFilter.txx"
00329 #endif
00330
00331 #endif // __itkSimplexMeshAdaptTopologyFilter_h
00332