X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=itk%2FclitkInvertVFFilter.txx;h=6c91dbb9ebf4fa14d72a411b5bf0fad28ef42986;hb=69dacfe0651ee24bef0f9e6b41171b9eec96fd2e;hp=ae6be89f8e85c424dbbc9c7a2c39bbc5f0947da3;hpb=2252420857d1554424d1f2c69ea68547bc759678;p=clitk.git

diff --git a/itk/clitkInvertVFFilter.txx b/itk/clitkInvertVFFilter.txx
index ae6be89..6c91dbb 100644
--- a/itk/clitkInvertVFFilter.txx
+++ b/itk/clitkInvertVFFilter.txx
@@ -47,7 +47,9 @@ public:
   //Typedefs
   typedef typename OutputImageType::PixelType        PixelType;
   typedef itk::Image<double, ImageDimension > WeightsImageType;
-  typedef itk::Image<itk::SimpleFastMutexLock, ImageDimension > MutexImageType;
+#if ITK_VERSION_MAJOR <= 4
+  typedef itk::Image<itk::SimpleFastMutexLock, ImageDimension> MutexImageType;
+#endif
 
   //===================================================================================
   //Set methods
@@ -55,11 +57,18 @@ public:
     m_Weights = input;
     this->Modified();
   }
+#if ITK_VERSION_MAJOR <= 4
   void SetMutexImage(const typename MutexImageType::Pointer input) {
     m_MutexImage=input;
     this->Modified();
     m_ThreadSafe=true;
   }
+#else
+  void SetMutexImage() {
+    this->Modified();
+    m_ThreadSafe=true;
+  }
+#endif
 
   //Get methods
   typename  WeightsImageType::Pointer GetWeights() {
@@ -74,14 +83,16 @@ protected:
   ~HelperClass1() {};
 
   //the actual processing
-  void GenerateInputRequestedRegion();
-  void GenerateOutputInformation();  
-  void BeforeThreadedGenerateData();
-  void ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, int threadId );
+  void BeforeThreadedGenerateData() ITK_OVERRIDE;
+  void ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId ) ITK_OVERRIDE;
 
   //member data
   typename  WeightsImageType::Pointer m_Weights;
+#if ITK_VERSION_MAJOR <= 4
   typename MutexImageType::Pointer m_MutexImage;
+#else
+  std::mutex m_Mutex;
+#endif
   bool m_ThreadSafe;
 
 };
@@ -101,33 +112,6 @@ HelperClass1<InputImageType, OutputImageType>::HelperClass1()
   m_ThreadSafe=false;
 }
 
-//=========================================================================================================================
-//Set input image params
-template <class InputImageType, class OutputImageType>
-void HelperClass1<InputImageType, OutputImageType>::GenerateInputRequestedRegion()
-{
-  //std::cout << "HelperClass1::GenerateInputRequestedRegion - IN" << std::endl;
-
-  typename InputImageType::Pointer inputPtr = const_cast< InputImageType * >(this->GetInput());
-  inputPtr->SetRequestedRegion(inputPtr->GetLargestPossibleRegion());
-  
-  //std::cout << "HelperClass1::GenerateInputRequestedRegion - OUT" << std::endl;
-}
-
-//=========================================================================================================================
-//Set output image params
-template<class InputImageType, class OutputImageType >
-void HelperClass1<InputImageType, OutputImageType>::GenerateOutputInformation()
-{
-  //std::cout << "HelperClass1::GenerateOutputInformation - IN" << std::endl;
-  Superclass::GenerateOutputInformation();
-
-  // it's the weight image that determines the size of the output
-  typename OutputImageType::Pointer outputPtr = this->GetOutput();
-  outputPtr->SetLargestPossibleRegion( m_Weights->GetLargestPossibleRegion() );
-  outputPtr->SetSpacing(m_Weights->GetSpacing());
-}
-
 //=========================================================================================================================
 //Before threaded data
 template<class InputImageType, class OutputImageType >
@@ -142,26 +126,21 @@ void HelperClass1<InputImageType, OutputImageType>::BeforeThreadedGenerateData()
 //=========================================================================================================================
 //update the output for the outputRegionForThread
 template<class InputImageType, class OutputImageType>
-void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, int threadId )
+void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId )
 {
-  //std::cout << "HelperClass1::ThreadedGenerateData - IN" << std::endl;
+//   std::cout << "HelperClass1::ThreadedGenerateData - IN " << threadId << std::endl;
   //Get pointer to the input
   typename InputImageType::ConstPointer inputPtr = this->GetInput();
 
   //Get pointer to the output
   typename OutputImageType::Pointer outputPtr = this->GetOutput();
-/*  outputPtr->SetLargestPossibleRegion( m_Weights->GetLargestPossibleRegion() );
-  outputPtr->SetSpacing(m_Weights->GetSpacing());*/
-  typename OutputImageType::SizeType size=outputPtr->GetLargestPossibleRegion().GetSize();
-  //std::cout << outputPtr->GetSpacing() << size << std::endl;
+  //typename OutputImageType::SizeType size=outputPtr->GetLargestPossibleRegion().GetSize();
 
   //Iterator over input
-  //typedef itk::ImageRegionConstIteratorWithIndex<InputImageType> InputImageIteratorType;
-  typedef itk::ImageRegionConstIteratorWithIndex<OutputImageType> InputImageIteratorType;
+  typedef itk::ImageRegionConstIteratorWithIndex<InputImageType> InputImageIteratorType;
 
   //define them over the outputRegionForThread
-  //InputImageIteratorType outputIt(inputPtr, outputRegionForThread);
-  InputImageIteratorType outputIt(outputPtr, outputPtr->GetLargestPossibleRegion());
+  InputImageIteratorType inputIt(inputPtr, outputRegionForThread);
 
   //Initialize
   typename InputImageType::IndexType index;
@@ -170,14 +149,16 @@ void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const O
   typename OutputImageType::PointType opoint;
   typedef typename OutputImageType::PixelType DisplacementType;
   DisplacementType displacement;
-  outputIt.GoToBegin();
+  inputIt.GoToBegin();
+
+  typename OutputImageType::SizeType size = outputPtr->GetLargestPossibleRegion().GetSize();
 
   //define some temp variables
   signed long baseIndex[ImageDimension];
   double distance[ImageDimension];
   unsigned int dim, counter, upper;
   double totalOverlap,overlap;
-  typename OutputImageType::IndexType neighIndex, outIndex;
+  typename OutputImageType::IndexType neighIndex;
 
   //Find the number of neighbors
   unsigned int neighbors =  1 << ImageDimension;
@@ -185,17 +166,13 @@ void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const O
   //==================================================================================================
   //Loop over the output region and add the intensities from the input to the output and the weight to the weights
   //==================================================================================================
-  while( !outputIt.IsAtEnd() ) {
+  while( !inputIt.IsAtEnd() ) {
     // get the input image index
-    outIndex = outputIt.GetIndex();
-    outputPtr->TransformIndexToPhysicalPoint( outIndex,opoint );
-    for(unsigned int j = 0; j < ImageDimension; j++ ) ipoint[j] = opoint[j];
-    inputPtr->TransformPhysicalPointToIndex(ipoint, index);
+    index = inputIt.GetIndex();
     inputPtr->TransformIndexToPhysicalPoint( index,ipoint );
 
     // get the required displacement
-    //displacement = outputIt.Get();
-    displacement = inputPtr->GetPixel(index);
+    displacement = inputIt.Get();
 
     // compute the required output image point
     for(unsigned int j = 0; j < ImageDimension; j++ ) opoint[j] = ipoint[j] + (double)displacement[j];
@@ -205,7 +182,7 @@ void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const O
       for(dim = 0; dim < ImageDimension; dim++) {
         // The following  block is equivalent to the following line without
         // having to call floor. (Only for positive inputs, we already now that is in the image)
-        // baseIndex[dim] = (long) vcl_floor(contIndex[dim] );
+        // baseIndex[dim] = (long) std::floor(contIndex[dim] );
 
         baseIndex[dim] = (long) contIndex[dim];
         distance[dim] = contIndex[dim] - double( baseIndex[dim] );
@@ -227,8 +204,10 @@ void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const O
             overlap *= 1.0 - distance[dim];
           }
           upper >>= 1;
-        }
 
+          if (neighIndex[dim] >= size[dim])
+            neighIndex[dim] = size[dim] - 1;
+        }
 
 
         //Set neighbor value only if overlap is not zero
@@ -249,14 +228,22 @@ void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const O
 
           else {
             //Entering critilal section: shared memory
+#if ITK_VERSION_MAJOR <= 4
             m_MutexImage->GetPixel(neighIndex).Lock();
+#else
+            m_Mutex.lock();
+#endif
 
             //Set the pixel and weight at neighIndex
             outputPtr->SetPixel(neighIndex, outputPtr->GetPixel(neighIndex) - (displacement*overlap));
             m_Weights->SetPixel(neighIndex, m_Weights->GetPixel(neighIndex) + overlap);
 
             //Unlock
+#if ITK_VERSION_MAJOR <= 4
             m_MutexImage->GetPixel(neighIndex).Unlock();
+#else
+            m_Mutex.unlock();
+#endif
 
           }
           //Add to total overlap
@@ -270,10 +257,10 @@ void HelperClass1<InputImageType, OutputImageType>::ThreadedGenerateData(const O
       }
     }
 
-    ++outputIt;
+    ++inputIt;
   }
 
-  //std::cout << "HelperClass1::ThreadedGenerateData - OUT" << std::endl;
+//   std::cout << "HelperClass1::ThreadedGenerateData - OUT " << threadId << std::endl;
 }
 
 
@@ -323,8 +310,7 @@ protected:
 
 
   //the actual processing
-  void ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, int threadId );
-
+  void ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId ) ITK_OVERRIDE;
 
   //member data
   typename     WeightsImageType::Pointer m_Weights;
@@ -343,16 +329,19 @@ protected:
 //Empty constructor
 template<class InputImageType, class OutputImageType > HelperClass2<InputImageType, OutputImageType>::HelperClass2()
 {
-  m_EdgePaddingValue=itk::NumericTraits<PixelType>::Zero;
+  PixelType zero;
+  for(unsigned int i=0;i <PixelType::Dimension; i++) zero[i] = 0.0;
+  m_EdgePaddingValue=zero;
+  //m_EdgePaddingValue=itk::NumericTraits<PixelType>::Zero;
 }
 
 
 //=========================================================================================================================
 //update the output for the outputRegionForThread
-template<class InputImageType, class OutputImageType > void HelperClass2<InputImageType, OutputImageType>::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, int threadId )
+template<class InputImageType, class OutputImageType > void HelperClass2<InputImageType, OutputImageType>::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, itk::ThreadIdType threadId )
 {
-  //std::cout << "HelperClass2::ThreadedGenerateData - IN" << std::endl;
-  
+//   std::cout << "HelperClass2::ThreadedGenerateData - IN " << threadId << std::endl;
+
   //Get pointer to the input
   typename InputImageType::ConstPointer inputPtr = this->GetInput();
 
@@ -402,9 +391,9 @@ template<class InputImageType, class OutputImageType > void HelperClass2<InputIm
     ++inputIt;
 
   }//end while
-  
-  //std::cout << "HelperClass2::ThreadedGenerateData - OUT" << std::endl;
-  
+
+//   std::cout << "HelperClass2::ThreadedGenerateData - OUT " << threadId << std::endl;
+
 }//end member
 
 
@@ -424,65 +413,32 @@ namespace clitk
 template <class InputImageType, class OutputImageType>
 InvertVFFilter<InputImageType, OutputImageType>::InvertVFFilter()
 {
-  m_EdgePaddingValue=itk::NumericTraits<PixelType>::Zero; //no other reasonable value?
-  m_ThreadSafe=false;
-  m_Verbose=false;
-}
-
-//=========================================================================================================================
-//Set input image params
-template <class InputImageType, class OutputImageType>
-void InvertVFFilter<InputImageType, OutputImageType>::GenerateInputRequestedRegion()
-{
-  //std::cout << "InvertVFFilter::GenerateInputRequestedRegion - IN" << std::endl;
-  // call the superclass' implementation of this method
-  // Superclass::GenerateInputRequestedRegion();
-
-  typename InputImageType::Pointer inputPtr = const_cast< InputImageType * >(this->GetInput());
-  inputPtr->SetRequestedRegion(inputPtr->GetLargestPossibleRegion());
-  
-  //std::cout << "InvertVFFilter::GenerateInputRequestedRegion - OUT" << std::endl;
-}
 
-//=========================================================================================================================
-//Set output image params
-template<class InputImageType, class OutputImageType >
-void InvertVFFilter<InputImageType, OutputImageType>::GenerateOutputInformation()
-{
-  //std::cout << "InvertVFFilter::GenerateOutputInformation - IN" << std::endl;
-  Superclass::GenerateOutputInformation();
+  //m_EdgePaddingValue=itk::NumericTraits<PixelType>::Zero; //no other reasonable value?
+  PixelType zero;
+  for(unsigned int i=0;i <PixelType::Dimension; i++) zero[i] = 0.0;
+  m_EdgePaddingValue=zero; //no other reasonable value?
 
-  typename InputImageType::ConstPointer inputPtr=this->GetInput();
-  typename WeightsImageType::RegionType region = inputPtr->GetLargestPossibleRegion();
-  region.SetSize(m_OutputSize);
-
-  typename OutputImageType::Pointer outputPtr = this->GetOutput();
-  outputPtr->SetRegions( region );
-  outputPtr->SetSpacing(m_OutputSpacing);
-  outputPtr->SetOrigin(inputPtr->GetOrigin());
-
-  //std::cout << "InvertVFFilter::GenerateOutputInformation - OUT" << std::endl;
+  m_ThreadSafe=false;
+  m_Verbose=false;
 }
 
 //=========================================================================================================================
 //Update
 template <class InputImageType, class OutputImageType> void InvertVFFilter<InputImageType, OutputImageType>::GenerateData()
 {
-  //std::cout << "InvertVFFilter::GenerateData - IN" << std::endl;
+  // std::cout << "InvertVFFilter::GenerateData - IN" << std::endl;
 
   //Get the properties of the input
   typename InputImageType::ConstPointer inputPtr=this->GetInput();
   typename WeightsImageType::RegionType region = inputPtr->GetLargestPossibleRegion();
-  //region.SetSize(size);
-  region.SetSize(m_OutputSize);
 
   //Allocate the weights
   typename WeightsImageType::Pointer weights=WeightsImageType::New();
   weights->SetOrigin(inputPtr->GetOrigin());
   weights->SetRegions(region);
-  weights->SetSpacing(m_OutputSpacing);
   weights->Allocate();
-  //weights->SetSpacing(inputPtr->GetSpacing());
+  weights->SetSpacing(inputPtr->GetSpacing());
 
   //===========================================================================
   //Inversion is divided in in two loops, for each we will call a threaded helper class
@@ -503,19 +459,28 @@ template <class InputImageType, class OutputImageType> void InvertVFFilter<Input
   typename HelperClass1Type::Pointer helper1=HelperClass1Type::New();
 
   //Set input
-  if(m_NumberOfThreadsIsGiven)helper1->SetNumberOfThreads(m_NumberOfThreads);
+  if(m_NumberOfThreadsIsGiven) {
+#if ITK_VERSION_MAJOR <= 4
+    helper1->SetNumberOfThreads(m_NumberOfThreads);
+#else
+    helper1->SetNumberOfWorkUnits(m_NumberOfWorkUnits);
+#endif
+  }
   helper1->SetInput(inputPtr);
   helper1->SetWeights(weights);
 
   //Threadsafe?
   if(m_ThreadSafe) {
     //Allocate the mutex image
+#if ITK_VERSION_MAJOR <= 4
     typename MutexImageType::Pointer mutex=InvertVFFilter::MutexImageType::New();
     mutex->SetRegions(region);
     mutex->Allocate();
-    //mutex->SetSpacing(inputPtr->GetSpacing());
-    mutex->SetSpacing(m_OutputSpacing);
+    mutex->SetSpacing(inputPtr->GetSpacing());
     helper1->SetMutexImage(mutex);
+#else
+    helper1->SetMutexImage();
+#endif
     if (m_Verbose) std::cout <<"Inverting using a thread-safe algorithm" <<std::endl;
   } else  if(m_Verbose)std::cout <<"Inverting using a thread-unsafe algorithm" <<std::endl;
 
@@ -534,8 +499,13 @@ template <class InputImageType, class OutputImageType> void InvertVFFilter<Input
   typename HelperClass2Type::Pointer helper2=HelperClass2Type::New();
 
   //Set temporary output as input
-  //std::cout << temp->GetSpacing() << temp->GetLargestPossibleRegion().GetSize() << std::endl;
-  //std::cout << weights->GetSpacing() << weights->GetLargestPossibleRegion().GetSize() << std::endl;
+  if(m_NumberOfThreadsIsGiven) {
+#if ITK_VERSION_MAJOR <= 4
+    helper2->SetNumberOfThreads(m_NumberOfThreads);
+#else
+    helper2->SetNumberOfWorkUnits(m_NumberOfWorkUnits);
+#endif
+  }
   helper2->SetInput(temp);
   helper2->SetWeights(weights);
   helper2->SetEdgePaddingValue(m_EdgePaddingValue);
@@ -546,6 +516,8 @@ template <class InputImageType, class OutputImageType> void InvertVFFilter<Input
 
   //Set the output
   this->SetNthOutput(0, helper2->GetOutput());
+
+  //std::cout << "InvertVFFilter::GenerateData - OUT" << std::endl;
 }