X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=itk%2FclitkSegmentationUtils.txx;h=4990b247c2b85dac1e342b77dc64901c80c553f6;hb=f4836557301d6fafe12dd9f767f911796992f855;hp=42ddc4759867b1231016ca847ade9c65d403068b;hpb=88b91aa3b3d43bb047dc2fe9637e4361405b3d38;p=clitk.git

diff --git a/itk/clitkSegmentationUtils.txx b/itk/clitkSegmentationUtils.txx
index 42ddc47..4990b24 100644
--- a/itk/clitkSegmentationUtils.txx
+++ b/itk/clitkSegmentationUtils.txx
@@ -3,7 +3,7 @@
 
   Authors belong to: 
   - University of LYON              http://www.universite-lyon.fr/
-  - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
+  - Léon Bérard cancer center       http://www.centreleonberard.fr
   - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
 
   This software is distributed WITHOUT ANY WARRANTY; without even
@@ -32,6 +32,8 @@
 #include <itkBinaryDilateImageFilter.h>
 #include <itkConstantPadImageFilter.h>
 #include <itkImageSliceIteratorWithIndex.h>
+#include <itkBinaryMorphologicalOpeningImageFilter.h>
+#include <itkImageDuplicator.h>
 
 namespace clitk {
 
@@ -297,7 +299,8 @@ namespace clitk {
   {
     typename ImageType::Pointer working_image;
     working_image = Labelize<ImageType>(input, BG, isFullyConnected, minimalComponentSize);
-    working_image = RemoveLabels<ImageType>(working_image, BG, param->GetLabelsToRemove());
+    if (param->GetLabelsToRemove().size() != 0)
+      working_image = RemoveLabels<ImageType>(working_image, BG, param->GetLabelsToRemove());
     working_image = KeepLabels<ImageType>(working_image, 
                                           BG, FG, 
                                           param->GetFirstKeep(), 
@@ -336,7 +339,8 @@ namespace clitk {
                                std::string orientation, 
                                bool uniqueConnectedComponent, 
                                double spacing, 
-                               bool inverseflag) 
+                               bool autocropFlag, 
+                               bool singleObjectCCL) 
   {
     typedef SliceBySliceRelativePositionFilter<MaskImageType> SliceRelPosFilterType;
     typename SliceRelPosFilterType::Pointer sliceRelPosFilter = SliceRelPosFilterType::New();
@@ -350,8 +354,10 @@ namespace clitk {
     sliceRelPosFilter->SetIntermediateSpacingFlag((spacing != -1));
     sliceRelPosFilter->SetIntermediateSpacing(spacing);
     sliceRelPosFilter->SetUniqueConnectedComponentBySlice(uniqueConnectedComponent);
-    sliceRelPosFilter->SetInverseOrientationFlag(inverseflag);
-    //  sliceRelPosFilter->SetAutoCropFlag(true); ??
+    sliceRelPosFilter->SetUseASingleObjectConnectedComponentBySliceFlag(singleObjectCCL);
+    //    sliceRelPosFilter->SetInverseOrientationFlag(inverseflag); 
+    sliceRelPosFilter->SetAutoCropFlag(autocropFlag); 
+    sliceRelPosFilter->IgnoreEmptySliceObjectFlagOn();
     sliceRelPosFilter->Update();
     return sliceRelPosFilter->GetOutput();
   }
@@ -384,7 +390,7 @@ namespace clitk {
     /*
       loop over input pixels, store the index in the fg that is max
       according to the given direction. 
-    */
+    */    
     typedef itk::ImageRegionConstIteratorWithIndex<ImageType> IteratorType;
     IteratorType iter(input, input->GetLargestPossibleRegion());
     iter.GoToBegin();
@@ -416,7 +422,7 @@ namespace clitk {
   //--------------------------------------------------------------------
   template<class ImageType>
   typename ImageType::Pointer
-  CropImageAbove(const ImageType * image, 
+  CropImageRemoveGreaterThan(const ImageType * image, 
                  int dim, double min, bool autoCrop,
                  typename ImageType::PixelType BG) 
   {
@@ -431,7 +437,7 @@ namespace clitk {
   //--------------------------------------------------------------------
   template<class ImageType>
   typename ImageType::Pointer
-  CropImageBelow(const ImageType * image, 
+  CropImageRemoveLowerThan(const ImageType * image, 
                  int dim, double max, bool autoCrop,
                  typename ImageType::PixelType BG) 
   {
@@ -502,12 +508,79 @@ namespace clitk {
     statFilter->Update();
     typename LabelMapType::Pointer labelMap = statFilter->GetOutput();
 
+    centroids.clear();
+    typename ImageType::PointType dummy;
+    centroids.push_back(dummy); // label 0 -> no centroid, use dummy point for BG 
+    //DS FIXME (not useful ! to change ..)
+    for(uint i=0; i<labelMap->GetNumberOfLabelObjects(); i++) {
+      int label = labelMap->GetLabels()[i];
+      centroids.push_back(labelMap->GetLabelObject(label)->GetCentroid());
+    } 
+  }
+  //--------------------------------------------------------------------
+
+
+  //--------------------------------------------------------------------
+  template<class ImageType, class LabelType>
+  typename itk::LabelMap< itk::ShapeLabelObject<LabelType, ImageType::ImageDimension> >::Pointer
+  ComputeLabelMap(const ImageType * image, 
+                  typename ImageType::PixelType BG, 
+                  bool computePerimeterFlag) 
+  {
+    static const unsigned int Dim = ImageType::ImageDimension;
+    typedef itk::ShapeLabelObject< LabelType, Dim > LabelObjectType;
+    typedef itk::LabelMap< LabelObjectType > LabelMapType;
+    typedef itk::LabelImageToLabelMapFilter<ImageType, LabelMapType> ImageToMapFilterType;
+    typename ImageToMapFilterType::Pointer imageToLabelFilter = ImageToMapFilterType::New(); 
+    typedef itk::ShapeLabelMapFilter<LabelMapType, ImageType> ShapeFilterType; 
+    typename ShapeFilterType::Pointer statFilter = ShapeFilterType::New();
+    imageToLabelFilter->SetBackgroundValue(BG);
+    imageToLabelFilter->SetInput(image);
+    statFilter->SetInput(imageToLabelFilter->GetOutput());
+    statFilter->SetComputePerimeter(computePerimeterFlag);
+    statFilter->Update();
+    return statFilter->GetOutput();
+  }
+  //--------------------------------------------------------------------
+
+
+  //--------------------------------------------------------------------
+  template<class ImageType>
+  void
+  ComputeCentroids2(const ImageType * image, 
+                   typename ImageType::PixelType BG, 
+                   std::vector<typename ImageType::PointType> & centroids) 
+  {
+    typedef long LabelType;
+    static const unsigned int Dim = ImageType::ImageDimension;
+    typedef itk::ShapeLabelObject< LabelType, Dim > LabelObjectType;
+    typedef itk::LabelMap< LabelObjectType > LabelMapType;
+    typedef itk::LabelImageToLabelMapFilter<ImageType, LabelMapType> ImageToMapFilterType;
+    typename ImageToMapFilterType::Pointer imageToLabelFilter = ImageToMapFilterType::New(); 
+    typedef itk::ShapeLabelMapFilter<LabelMapType, ImageType> ShapeFilterType; 
+    typename ShapeFilterType::Pointer statFilter = ShapeFilterType::New();
+    imageToLabelFilter->SetBackgroundValue(BG);
+    imageToLabelFilter->SetInput(image);
+    statFilter->SetInput(imageToLabelFilter->GetOutput());
+    statFilter->Update();
+    typename LabelMapType::Pointer labelMap = statFilter->GetOutput();
+
     centroids.clear();
     typename ImageType::PointType dummy;
     centroids.push_back(dummy); // label 0 -> no centroid, use dummy point
     for(uint i=1; i<labelMap->GetNumberOfLabelObjects()+1; i++) {
       centroids.push_back(labelMap->GetLabelObject(i)->GetCentroid());
     } 
+    
+    for(uint i=1; i<labelMap->GetNumberOfLabelObjects()+1; i++) {
+      DD(labelMap->GetLabelObject(i)->GetBinaryPrincipalAxes());
+      DD(labelMap->GetLabelObject(i)->GetBinaryFlatness());
+      DD(labelMap->GetLabelObject(i)->GetRoundness ());      
+
+      // search for the point on the boundary alog PA
+
+    }
+
   }
   //--------------------------------------------------------------------
 
@@ -531,27 +604,6 @@ namespace clitk {
   //--------------------------------------------------------------------
 
 
-  //--------------------------------------------------------------------
-  template<class ImageType>
-  typename ImageType::Pointer
-  JoinSlices(std::vector<typename itk::Image<typename ImageType::PixelType, 
-                                             ImageType::ImageDimension-1>::Pointer > & slices, 
-             const ImageType * input, 
-             int direction) {
-    typedef typename itk::Image<typename ImageType::PixelType, ImageType::ImageDimension-1> SliceType;
-    typedef itk::JoinSeriesImageFilter<SliceType, ImageType> JoinSeriesFilterType;
-    typename JoinSeriesFilterType::Pointer joinFilter = JoinSeriesFilterType::New();
-    joinFilter->SetOrigin(input->GetOrigin()[direction]);
-    joinFilter->SetSpacing(input->GetSpacing()[direction]);
-    for(unsigned int i=0; i<slices.size(); i++) {
-      joinFilter->PushBackInput(slices[i]);
-    }
-    joinFilter->Update();
-    return joinFilter->GetOutput();
-  }
-  //--------------------------------------------------------------------
-
-
   //--------------------------------------------------------------------
   template<class ImageType>
   void
@@ -575,7 +627,7 @@ namespace clitk {
   //--------------------------------------------------------------------
   template<class ImageType>
   void 
-  PointsUtils<ImageType>::Convert2DTo3DList(const MapPoint2DType & map, 
+  PointsUtils<ImageType>::Convert2DMapTo3DList(const MapPoint2DType & map, 
                                             const ImageType * image, 
                                             VectorPoint3DType & list)
   {
@@ -589,6 +641,24 @@ namespace clitk {
   }
   //--------------------------------------------------------------------
 
+
+  //--------------------------------------------------------------------
+  template<class ImageType>
+  void 
+  PointsUtils<ImageType>::Convert2DListTo3DList(const VectorPoint2DType & p2D, 
+                                                int slice,
+                                                const ImageType * image, 
+                                                VectorPoint3DType & list) 
+  {
+    for(uint i=0; i<p2D.size(); i++) {
+      PointType3D p;
+      Convert2DTo3D(p2D[i], image, slice, p);
+      list.push_back(p);
+    }
+  }
+  //--------------------------------------------------------------------
+
+
   //--------------------------------------------------------------------
   template<class ImageType>
   void 
@@ -685,6 +755,34 @@ namespace clitk {
   //--------------------------------------------------------------------
 
 
+  //--------------------------------------------------------------------
+  template<class ImageType>
+  typename ImageType::Pointer 
+  Opening(const ImageType * image, typename ImageType::SizeType radius,
+         typename ImageType::PixelType BG,
+         typename ImageType::PixelType FG)
+  {
+    // Kernel 
+    typedef itk::BinaryBallStructuringElement<typename ImageType::PixelType, 
+                                              ImageType::ImageDimension> KernelType;    
+    KernelType structuringElement;
+    structuringElement.SetRadius(radius);
+    structuringElement.CreateStructuringElement();
+    
+    // Filter
+    typedef itk::BinaryMorphologicalOpeningImageFilter<ImageType, ImageType , KernelType> OpeningFilterType;
+    typename OpeningFilterType::Pointer open = OpeningFilterType::New();
+    open->SetInput(image);
+    open->SetBackgroundValue(BG);
+    open->SetForegroundValue(FG);
+    open->SetKernel(structuringElement);
+    open->Update();
+    return open->GetOutput();
+  }
+  //--------------------------------------------------------------------
+
+
+
   //--------------------------------------------------------------------
   template<class ValueType, class VectorType>
   void ConvertOption(std::string optionName, uint given, 
@@ -837,6 +935,7 @@ namespace clitk {
     typedef itk::BinaryThresholdImageFilter<ImageType, ImageType> BinaryThresholdFilterType;
     typename BinaryThresholdFilterType::Pointer binarizeFilter = BinaryThresholdFilterType::New();
     binarizeFilter->SetInput(input);
+    binarizeFilter->InPlaceOff();
     binarizeFilter->SetLowerThreshold(lower);
     binarizeFilter->SetUpperThreshold(upper);
     binarizeFilter->SetInsideValue(FG);
@@ -919,5 +1018,101 @@ namespace clitk {
   //--------------------------------------------------------------------
 
 
+  //--------------------------------------------------------------------
+  /* Consider an input object, for each slice, find the extrema
+     position according to a given direction and build a line segment
+     passing throught this point in a given direction.  Output is a
+     vector of line (from point A to B), for each slice;
+   */
+  template<class ImageType>
+  void 
+  SliceBySliceBuildLineSegmentAccordingToExtremaPosition(const ImageType * input, 
+                                                         typename ImageType::PixelType BG, 
+                                                         int sliceDimension, 
+                                                         int extremaDirection, 
+                                                         bool extremaOppositeFlag, 
+                                                         int lineDirection,
+                                                         double margin,
+                                                         std::vector<typename ImageType::PointType> & A, 
+                                                         std::vector<typename ImageType::PointType> & B)
+  {
+    // Type of a slice
+    typedef typename itk::Image<typename ImageType::PixelType, ImageType::ImageDimension-1> SliceType;
+    
+    // Build the list of slices
+    std::vector<typename SliceType::Pointer> slices;
+    clitk::ExtractSlices<ImageType>(input, sliceDimension, slices);
+
+    // Build the list of 2D points
+    std::map<int, typename SliceType::PointType> position2D;
+    for(uint i=0; i<slices.size(); i++) {
+      typename SliceType::PointType p;
+      bool found = 
+        clitk::FindExtremaPointInAGivenDirection<SliceType>(slices[i], BG, 
+                                                            extremaDirection, extremaOppositeFlag, p);
+      if (found) {
+        position2D[i] = p;
+      }    
+    }
+    
+    // Convert 2D points in slice into 3D points
+    clitk::PointsUtils<ImageType>::Convert2DMapTo3DList(position2D, input, A);
+    
+    // Create additional point just right to the previous ones, on the
+    // given lineDirection, in order to create a horizontal/vertical line.
+    for(uint i=0; i<A.size(); i++) {
+      typename ImageType::PointType p = A[i];
+      p[lineDirection] += 10;
+      B.push_back(p);
+      // Margins ?
+      A[i][1] += margin;
+      B[i][1] += margin;
+    }
+
+  }
+  //--------------------------------------------------------------------
+
+
+  //--------------------------------------------------------------------
+  template<class ImageType>
+  typename ImageType::Pointer
+  SliceBySliceKeepMainCCL(const ImageType * input, 
+                          typename ImageType::PixelType BG,
+                          typename ImageType::PixelType FG)  {
+    
+    // Extract slices
+    const int d = ImageType::ImageDimension-1;
+    typedef typename itk::Image<typename ImageType::PixelType, d> SliceType;
+    std::vector<typename SliceType::Pointer> slices;
+    clitk::ExtractSlices<ImageType>(input, d, slices);
+    
+    // Labelize and keep the main one
+    std::vector<typename SliceType::Pointer> o;
+    for(uint i=0; i<slices.size(); i++) {
+      o.push_back(clitk::Labelize<SliceType>(slices[i], BG, false, 1));
+      o[i] = clitk::KeepLabels<SliceType>(o[i], BG, FG, 1, 1, true);
+    }
+    
+    // Join slices
+    typename ImageType::Pointer output;
+    output = clitk::JoinSlices<ImageType>(o, input, d);
+    return output;
+  }
+  //--------------------------------------------------------------------
+
+
+  //--------------------------------------------------------------------
+  template<class ImageType>
+  typename ImageType::Pointer
+  Clone(const ImageType * input) {
+    typedef itk::ImageDuplicator<ImageType> DuplicatorType;
+    typename DuplicatorType::Pointer duplicator = DuplicatorType::New();
+    duplicator->SetInputImage(input);
+    duplicator->Update();
+    return duplicator->GetOutput();
+  }
+  //--------------------------------------------------------------------
+
+
 } // end of namespace