segmentation/clitkExtractImageFilter.txx

   1 /*=========================================================================
   2   Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
   3
   4   Authors belong to:
   5   - University of LYON              http://www.universite-lyon.fr/
   6   - Léon Bérard cancer center       http://www.centreleonberard.fr
   7   - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
   8
   9   This software is distributed WITHOUT ANY WARRANTY; without even
  10   the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  11   PURPOSE.  See the copyright notices for more information.
  12
  13   It is distributed under dual licence
  14
  15   - BSD        See included LICENSE.txt file
  16   - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  17 ===========================================================================*/
  18
  19 #ifndef _clitkExtractImageFilter_txx
  20 #define _clitkExtractImageFilter_txx
  21
  22 #include "clitkExtractImageFilter.h"
  23 #include "itkImageRegionIterator.h"
  24 #include "itkImageRegionConstIterator.h"
  25 #include "itkObjectFactory.h"
  26 #include "itkExtractImageFilterRegionCopier.h"
  27 #include "itkProgressReporter.h"
  28
  29
  30 namespace clitk
  31 {
  32
  33 /**
  34  *
  35  */
  36 template <class TInputImage, class TOutputImage>
  37 ExtractImageFilter<TInputImage,TOutputImage>
  38 ::ExtractImageFilter()
  39 {
  40 }
  41
  42
  43 /**
  44  *
  45  */
  46 template <class TInputImage, class TOutputImage>
  47 void
  48 ExtractImageFilter<TInputImage,TOutputImage>
  49 ::PrintSelf(std::ostream& os, itk::Indent indent) const
  50 {
  51   Superclass::PrintSelf(os,indent);
  52
  53   os << indent << "ExtractionRegion: " << m_ExtractionRegion << std::endl;
  54   os << indent << "OutputImageRegion: " << m_OutputImageRegion << std::endl;
  55 }
  56
  57
  58 template<class TInputImage, class TOutputImage>
  59 void
  60 ExtractImageFilter<TInputImage,TOutputImage>
  61 ::CallCopyOutputRegionToInputRegion(InputImageRegionType &destRegion,
  62                                     const OutputImageRegionType &srcRegion)
  63 {
  64   ExtractImageFilterRegionCopierType extractImageRegionCopier;
  65   extractImageRegionCopier(destRegion, srcRegion, m_ExtractionRegion);
  66 }
  67
  68
  69 template <class TInputImage, class TOutputImage>
  70 void
  71 ExtractImageFilter<TInputImage,TOutputImage>
  72 ::SetExtractionRegion(InputImageRegionType extractRegion)
  73 {
  74   m_ExtractionRegion = extractRegion;
  75
  76   unsigned int nonzeroSizeCount = 0;
  77   InputImageSizeType inputSize = extractRegion.GetSize();
  78   OutputImageSizeType outputSize;
  79   OutputImageIndexType outputIndex;
  80
  81   /**
  82    * check to see if the number of non-zero entries in the extraction region
  83    * matches the number of dimensions in the output image.
  84    **/
  85   for (unsigned int i = 0; i < InputImageDimension; ++i)
  86     {
  87     if (inputSize[i])
  88       {
  89       outputSize[nonzeroSizeCount] = inputSize[i];
  90       outputIndex[nonzeroSizeCount] = extractRegion.GetIndex()[i];
  91       nonzeroSizeCount++;
  92       }
  93     }
  94
  95   if (nonzeroSizeCount != OutputImageDimension)
  96     {
  97     itkExceptionMacro("Extraction Region not consistent with output image");
  98     }
  99
 100   m_OutputImageRegion.SetSize(outputSize);
 101   m_OutputImageRegion.SetIndex(outputIndex);
 102   this->Modified();
 103 }
 104
 105
 106
 107 /**
 108  * ExtractImageFilter can produce an image which is a different resolution
 109  * than its input image.  As such, ExtractImageFilter needs to provide an
 110  * implementation for GenerateOutputInformation() in order to inform
 111  * the pipeline execution model.  The original documentation of this
 112  * method is below.
 113  *
 114  * \sa ProcessObject::GenerateOutputInformaton()
 115  */
 116 template <class TInputImage, class TOutputImage>
 117 void
 118 ExtractImageFilter<TInputImage,TOutputImage>
 119 ::GenerateOutputInformation()
 120 {
 121   // do not call the superclass' implementation of this method since
 122   // this filter allows the input and the output to be of different dimensions
 123
 124   // get pointers to the input and output
 125   typename Superclass::OutputImagePointer      outputPtr = this->GetOutput();
 126   typename Superclass::InputImageConstPointer  inputPtr  = this->GetInput();
 127
 128   if ( !outputPtr || !inputPtr)
 129     {
 130     return;
 131     }
 132
 133   // Set the output image size to the same value as the extraction region.
 134   outputPtr->SetLargestPossibleRegion( m_OutputImageRegion );
 135
 136   // Set the output spacing and origin
 137   const itk::ImageBase<InputImageDimension> *phyData;
 138
 139   phyData
 140     = dynamic_cast<const itk::ImageBase<InputImageDimension>*>(this->GetInput());
 141
 142   if (phyData)
 143     {
 144     // Copy what we can from the image from spacing and origin of the input
 145     // This logic needs to be augmented with logic that select which
 146     // dimensions to copy
 147
 148     unsigned int i;
 149     const typename InputImageType::SpacingType&
 150       inputSpacing = inputPtr->GetSpacing();
 151     const typename InputImageType::DirectionType&
 152       inputDirection = inputPtr->GetDirection();
 153     const typename InputImageType::PointType&
 154       inputOrigin = inputPtr->GetOrigin();
 155
 156     typename OutputImageType::SpacingType outputSpacing;
 157     typename OutputImageType::DirectionType outputDirection;
 158     typename OutputImageType::PointType outputOrigin;
 159
 160     if ( static_cast<unsigned int>(OutputImageDimension) >
 161          static_cast<unsigned int>(InputImageDimension )    )
 162       {
 163       // copy the input to the output and fill the rest of the
 164       // output with zeros.
 165       for (i=0; i < InputImageDimension; ++i)
 166         {
 167         outputSpacing[i] = inputSpacing[i];
 168         outputOrigin[i] = inputOrigin[i];
 169         for (unsigned int dim = 0; dim < InputImageDimension; ++dim)
 170           {
 171           outputDirection[i][dim] = inputDirection[i][dim];
 172           }
 173         }
 174       for (; i < OutputImageDimension; ++i)
 175         {
 176         outputSpacing[i] = 1.0;
 177         outputOrigin[i] = 0.0;
 178         for (unsigned int dim = 0; dim < InputImageDimension; ++dim)
 179           {
 180           outputDirection[i][dim] = 0.0;
 181           }
 182         outputDirection[i][i] = 1.0;
 183         }
 184       }
 185     else
 186       {
 187       // copy the non-collapsed part of the input spacing and origing to the output
 188       int nonZeroCount = 0;
 189       for (i=0; i < InputImageDimension; ++i)
 190         {
 191         if (m_ExtractionRegion.GetSize()[i])
 192           {
 193           outputSpacing[nonZeroCount] = inputSpacing[i];
 194           outputOrigin[nonZeroCount] = inputOrigin[i];
 195           int nonZeroCount2 = 0;
 196           for (unsigned int dim = 0; dim < InputImageDimension; ++dim)
 197             {
 198             if (m_ExtractionRegion.GetSize()[dim])
 199               {
 200               outputDirection[nonZeroCount][nonZeroCount2] =
 201                 inputDirection[i][dim];
 202               ++nonZeroCount2;
 203               }
 204             }
 205           nonZeroCount++;
 206           }
 207         }
 208       }
 209
 210     // set the spacing and origin
 211     outputPtr->SetSpacing( outputSpacing );
 212     outputPtr->SetDirection( outputDirection );
 213     outputPtr->SetOrigin( outputOrigin );
 214     outputPtr->SetNumberOfComponentsPerPixel(inputPtr->GetNumberOfComponentsPerPixel() );
 215     }
 216   else
 217     {
 218     // pointer could not be cast back down
 219     itkExceptionMacro(<< "itk::ExtractImageFilter::GenerateOutputInformation "
 220                       << "cannot cast input to "
 221                       << typeid(itk::ImageBase<InputImageDimension>*).name() );
 222     }
 223 }
 224
 225 /**
 226    * ExtractImageFilter can be implemented as a multithreaded filter.
 227    * Therefore, this implementation provides a ThreadedGenerateData()
 228    * routine which is called for each processing thread. The output
 229    * image data is allocated automatically by the superclass prior to
 230    * calling ThreadedGenerateData().  ThreadedGenerateData can only
 231    * write to the portion of the output image specified by the
 232    * parameter "outputRegionForThread"
 233    *
 234    * \sa ImageToImageFilter::ThreadedGenerateData(),
 235    *     ImageToImageFilter::GenerateData()
 236    */
 237 template <class TInputImage, class TOutputImage>
 238 void
 239 ExtractImageFilter<TInputImage,TOutputImage>
 240 ::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread,
 241                        int threadId)
 242 {
 243   itkDebugMacro(<<"Actually executing");
 244
 245   // Get the input and output pointers
 246   typename Superclass::InputImageConstPointer  inputPtr = this->GetInput();
 247   typename Superclass::OutputImagePointer outputPtr = this->GetOutput();
 248
 249   // support progress methods/callbacks
 250  itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
 251
 252   // Define the portion of the input to walk for this thread
 253   InputImageRegionType inputRegionForThread;
 254   this->CallCopyOutputRegionToInputRegion(inputRegionForThread, outputRegionForThread);
 255
 256   // Define the iterators.
 257   typedef itk::ImageRegionIterator<TOutputImage> OutputIterator;
 258   typedef itk::ImageRegionConstIterator<TInputImage> InputIterator;
 259
 260   OutputIterator outIt(outputPtr, outputRegionForThread);
 261   InputIterator inIt(inputPtr, inputRegionForThread);
 262
 263   // walk the output region, and sample the input image
 264   while( !outIt.IsAtEnd() )
 265     {
 266     // copy the input pixel to the output
 267     outIt.Set( static_cast<OutputImagePixelType>(inIt.Get()));
 268     ++outIt;
 269     ++inIt;
 270     progress.CompletedPixel();
 271     }
 272 }
 273
 274 } // end namespace clitk
 275
 276 #endif