...

diff --git a/lib/fpa/Common/SliceBySliceRandomWalker.hxx b/lib/fpa/Common/SliceBySliceRandomWalker.hxx
index 54b1fbcfcebdc2735270a05d7754dc3b97fd0f98..d1fbd9cfdf495e39a5036db6fdbc1e2ba96c9c12 100644 (file)
--- a/lib/fpa/Common/SliceBySliceRandomWalker.hxx
+++ b/lib/fpa/Common/SliceBySliceRandomWalker.hxx
@@ -11,6 +11,8 @@
 #include <itkImageRegionIterator.h>
 #include <itkJoinSeriesImageFilter.h>
 #include <itkMinimumMaximumImageCalculator.h>
+#include <fpa/Common/RandomWalker.h>
+#include <fpa/Functors/Dijkstra/Image/Gaussian.h>
 
 // -------------------------------------------------------------------------
 template< class _TImage, class _TLabels, class _TScalarImage >
@@ -51,6 +53,8 @@ GenerateData( )
   typename TImage::SizeType sliceSize = region.GetSize( );
   int numSlices = sliceSize[ 2 ];
   sliceSize[ 2 ] = 0;
+  this->m_VesselnessValue =
+    this->m_VesselnessThreshold * this->m_VesselnessMax / TScalar( 100 );
 
   // Composite image
   typename TScalarImage::Pointer composite = TScalarImage::New( );
@@ -172,155 +176,76 @@ _GoDown(
   typedef typename _TSliceBinaryPtr::ObjectType _TSliceBinary;
   typedef typename _TSliceData3DPtr::ObjectType _TSliceData3D;
   typedef typename _TSliceFusionPtr::ObjectType _TSliceFusion;
+  typedef typename _TSliceFusion::PixelType     _TSliceScalar;
 
   int z = -1;
   while( z < numSlices - 2 )
   {
     z++;
     _TSliceBinaryPtr tmp = binaryTree[ z ];
-    _TSliceBinaryPtr next = binaryTree[ z  + 1];
+    _TSliceBinaryPtr next = binaryTree[ z  + 1 ];
     _TSliceData3DPtr data = data3D[ z + 1 ];
     _TSliceFusionPtr vess = fusion[ z + 1 ];
+    typename _TSliceBinary::RegionType region = tmp->GetRequestedRegion( );
 
-    // Create seeds
+    // Create seeds and background
     _TSliceBinaryPtr seeds = _TSliceBinary::New( );
     seeds->CopyInformation( tmp );
-    seeds->SetRequestedRegion( tmp->GetRequestedRegion( ) );
+    seeds->SetRequestedRegion( region );
     seeds->Allocate( );
     seeds->FillBuffer( 0 );
 
+    itk::ImageRegionConstIterator< _TSliceBinary > tIt( tmp, region );
+    itk::ImageRegionConstIterator< _TSliceFusion > vIt( vess, region );
+    itk::ImageRegionIterator< _TSliceBinary > sIt( seeds, region );
+    tIt.GoToBegin( );
+    vIt.GoToBegin( );
+    sIt.GoToBegin( );
+    int len = 0;
+    int numO = 0;
+    int numB = 0;
+    for( ; !tIt.IsAtEnd( ); ++tIt, ++vIt, ++sIt )
+    {
+      if( tIt.Get( ) != 0 )
+      {
+        len++;
+        if( vIt.Get( ) > this->m_VesselnessValue )
+        {
+          sIt.Set( 1 );
+          numO++;
+
+        } // fi
+
+      } // fi
+
+      if( vIt.Get( ) < this->m_VesselnessValue )
+      {
+        sIt.Set( 2 );
+        numB++;
+
+      } // fi
+
+    } // rof
+
+    if( len == 0 )
+      continue;
+
+    // Random walker function
+    typedef fpa::Functors::Dijkstra::Image::Gaussian< _TSliceData3D, _TSliceScalar > _TFunction;
+    typename TFunction::Pointer edge = TFunction::New( );
+    edge->SetBeta( this->m_Beta );
+    edge->SetEpsilon( this->m_Eps );
+    edge->TreatAsWeightOff( );
+
+    // Real random walker
+    typedef fpa::Common::RandomWalker< _TSliceData3D, _TSliceBinary, _TSliceScalar > _TRandomWalker;
+    typename _TRandomWalker::Pointer rw = _TRandomWalker::New( );
+    rw->SetInput( data );
+    rw->SetInputLabels( seeds );
+    rw->SetEdgeFunction( edge );
+
+    // Keep maximum
     /* TODO
-
-
-       SliceType::RegionType      sliceRegion = tmp->GetLargestPossibleRegion();
-       SliceType::SizeType                size = sliceRegion.GetSize();
-       SliceType::SpacingType     spacing = tmp->GetSpacing();
-       SliceType::PointType       origin = tmp->GetOrigin();
-       SliceType::IndexType       start = sliceRegion.GetIndex();
-
-       SliceType::Pointer seedsObject = SliceType::New();
-       seedsObject->SetRegions(sliceRegion);
-       seedsObject->SetOrigin(origin);
-       seedsObject->SetSpacing(spacing);
-       seedsObject->Allocate();
-       seedsObject->FillBuffer(0);
-
-       SliceType::Pointer seedsBackground = SliceType::New();
-       seedsBackground->SetRegions(sliceRegion);
-       seedsBackground->SetOrigin(origin);
-       seedsBackground->SetSpacing(spacing);
-       seedsBackground->Allocate();
-       seedsBackground->FillBuffer(0);
-
-       typedef itk::ImageRegionIteratorWithIndex<SliceType> SliceIteratorType;
-       SliceIteratorType itSliceMori(tmp, tmp->GetLargestPossibleRegion());
-
-       int len = 0;
-       int numO = 0;
-       int numB = 0;
-
-       itSliceMori.GoToBegin();
-
-       while (!itSliceMori.IsAtEnd())
-       {
-       if (itSliceMori.Get() != 0)
-       len++;
-
-       typename TImage::PixelType vessVal = vess->GetPixel(itSliceMori.GetIndex());
-       typename TImage::PixelType imgVal = data->GetPixel(itSliceMori.GetIndex());
-       typename TImage::PixelType bw = tmp->GetPixel(itSliceMori.GetIndex());
-
-       if (bw != 0 && vessVal > this->m_VessTh * this->m_vessMax / 100)
-       {
-       seedsObject->SetPixel(itSliceMori.GetIndex(), 1);
-       numO++;
-       }
-
-       itSliceMori++;
-       }
-
-       if (len == 0)
-       {
-       std::cout << "Slice " << z << ": - brak pikseli zakwalifikowanych do drzewa" << std::endl;
-       continue;
-       }
-
-       SliceIteratorType itSliceVess(vess, vess->GetLargestPossibleRegion());
-       itSliceVess.GoToBegin();
-
-       while (!itSliceVess.IsAtEnd())
-       {
-       if (itSliceVess.Get() < this->m_VessTh * this->m_vessMax / 100)
-       {
-       seedsBackground->SetPixel(itSliceVess.GetIndex(), 1);
-       numB++;
-       }
-       itSliceVess++;
-
-       }
-
-       std::cout << "slice: " << z << "seeds obj: " << numO << " seeds bkg " << numB << std::endl;
-       // przygotowania do Random Walkera
-
-       typename TImage::PixelType  *bufferImg = 0;
-       typename TImage::PixelType  *bufferSeedsObj = 0;
-       typename TImage::PixelType  *bufferSeedsBkg = 0;
-
-       int dims[2], dims2[2], dims3[2];
-       double spac[2], orig[2], spac2[2], orig2[2], spac3[2], orig3[2];
-
-       image2buffer(bufferImg, dims, spac, orig, 2, data);
-       image2buffer(bufferSeedsObj, dims2, spac2, orig2, 2, seedsObject);
-       image2buffer(bufferSeedsBkg, dims3, spac3, orig3, 2, seedsBackground);
-
-       // TO DO HERE
-       int width = dims[0];
-       int height = dims[1];
-
-       PixelType *mask = new PixelType[width*height];
-       double *probs = new double[2 * width*height];
-
-       int i, j, no_seeds = 0;
-
-       for (i = 0; i < width*height; i++)
-
-       if (bufferSeedsObj[i] == 1 || bufferSeedsBkg[i] == 1)
-       no_seeds++;
-
-       int *seed_indexes = new int[no_seeds];
-       int *seed_labels = new int[no_seeds];
-
-       for (j = 0, i = 0; i < width*height; i++)
-       {
-       if (bufferSeedsObj[i] == 1)
-       {
-       seed_indexes[j] = i + 1;//=i+1 (i)
-       seed_labels[j] = this->m_replaceValue;
-       j++;
-       }
-       else if (bufferSeedsBkg[i] == 1)
-       {
-       seed_indexes[j] = i + 1;//=i+1 (i)
-       seed_labels[j] = 0;
-       j++;
-       }
-       }
-
-       call_walker(mask, probs, bufferImg, height, width, seed_indexes, no_seeds, seed_labels, 2, this->m_beta);
-
-       delete[] bufferSeedsObj;
-       delete[] bufferSeedsBkg;
-       delete[] seed_indexes;
-       delete[] seed_labels;
-
-       SliceType::Pointer img2 = buffer2image(mask, dims, spac, orig, 2);
-
-       SliceType::PixelType *probs_img = new PixelType[width*height];
-       for (i = 0; i < width*height; i++)
-       probs_img[i] = floor(255 * probs[i]);
-
-       SliceType::Pointer img3 = buffer2image(probs_img, dims, spac, orig, 2);
-
        typedef itk::MaximumImageFilter<SliceType> MaxFilterType;
        MaxFilterType::Pointer maxFilter = MaxFilterType::New();
        maxFilter->SetInput(0, img2);
@@ -329,7 +254,6 @@ _GoDown(
 
        SliceType::Pointer sth = maxFilter->GetOutput();
        ImageAlgorithm::Copy(sth.GetPointer(), binaryTree[z + 1].GetPointer(), sth->GetRequestedRegion(), binaryTree[z + 1]->GetRequestedRegion());
-
        delete[] probs;
     */
   }