ITK4: add otb::WarpImageFilter

Code/UtilitiesAdapters/ITKPendingPatches/otbWarpImageFilter.h

+/*=========================================================================
+
+  Program:   Insight Segmentation & Registration Toolkit
+  Module:    $RCSfile: itkWarpImageFilter.h,v $
+  Language:  C++
+  Date:      $Date: 2009-10-29 11:19:00 $
+  Version:   $Revision: 1.31 $
+
+  Copyright (c) Insight Software Consortium. All rights reserved.
+  See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
+
+     This software is distributed WITHOUT ANY WARRANTY; without even
+     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+     PURPOSE.  See the above copyright notices for more information.
+
+=========================================================================*/
+#ifndef __otbWarpImageFilter_h
+#define __otbWarpImageFilter_h
+#include "itkImageBase.h"
+#include "itkImageToImageFilter.h"
+#include "itkInterpolateImageFunction.h"
+#include "itkLinearInterpolateImageFunction.h"
+#include "itkPoint.h"
+#include "itkFixedArray.h"
+
+namespace otb
+{
+
+/** \class WarpImageFilter
+ * \brief Warps an image using an input deformation field.
+ *
+ * WarpImageFilter warps an existing image with respect to
+ * a given deformation field.
+ *
+ * A deformation field is represented as a image whose pixel type is some
+ * vector type with at least N elements, where N is the dimension of
+ * the input image. The vector type must support element access via operator
+ * [].
+ *
+ * The output image is produced by inverse mapping: the output pixels
+ * are mapped back onto the input image. This scheme avoids the creation of
+ * any holes and overlaps in the output image.
+ *
+ * Each vector in the deformation field represent the distance between
+ * a geometric point in the input space and a point in the output space such
+ * that:
+ *
+ * \f[ p_{in} = p_{out} + d \f]
+ *
+ * Typically the mapped position does not correspond to an integer pixel
+ * position in the input image. Interpolation via an image function
+ * is used to compute values at non-integer positions. The default
+ * interpolation typed used is the LinearInterpolateImageFunction.
+ * The user can specify a particular interpolation function via
+ * SetInterpolator(). Note that the input interpolator must derive
+ * from base class InterpolateImageFunction.
+ *
+ * Position mapped to outside of the input image buffer are assigned
+ * a edge padding value.
+ *
+ * The LargetPossibleRegion for the output is inherited
+ * from the input deformation field. The output image
+ * spacing, origin and orientation may be set via
+ * SetOutputSpacing, SetOutputOrigin and
+ * SetOutputDirection. The default are respectively a
+ * vector of 1's, a vector of 0's and an identity matrix.
+ *
+ * This class is templated over the type of the input image, the
+ * type of the output image and the type of the deformation field.
+ *
+ * The input image is set via SetInput. The input deformation field
+ * is set via SetDeformationField.
+ *
+ * This filter is implemented as a multithreaded filter.
+ *
+ * \warning This filter assumes that the input type, output type
+ * and deformation field type all have the same number of dimensions.
+ *
+ * \ingroup GeometricTransforms MultiThreaded Streamed
+ */
+template <
+  class TInputImage,
+  class TOutputImage,
+  class TDeformationField
+  >
+class ITK_EXPORT WarpImageFilter :
+    public itk::ImageToImageFilter<TInputImage, TOutputImage>
+{
+public:
+  /** Standard class typedefs. */
+  typedef WarpImageFilter                              Self;
+  typedef itk::ImageToImageFilter<TInputImage,TOutputImage> Superclass;
+  typedef itk::SmartPointer<Self>                           Pointer;
+  typedef itk::SmartPointer<const Self>                     ConstPointer;
+
+  /** Method for creation through the object factory. */
+  itkNewMacro(Self);
+
+  /** Run-time type information (and related methods) */
+  itkTypeMacro( WarpImageFilter, ImageToImageFilter );
+
+  /** Typedef to describe the output image region type. */
+  typedef typename TOutputImage::RegionType OutputImageRegionType;
+
+  /** Inherit some types from the superclass. */
+  typedef typename Superclass::InputImageType         InputImageType;
+  typedef typename Superclass::InputImagePointer      InputImagePointer;
+  typedef typename Superclass::OutputImageType        OutputImageType;
+  typedef typename Superclass::OutputImagePointer     OutputImagePointer;
+  typedef typename Superclass::InputImageConstPointer InputImageConstPointer;
+  typedef typename OutputImageType::IndexType         IndexType;
+  typedef typename OutputImageType::IndexValueType    IndexValueType;
+  typedef typename OutputImageType::SizeType          SizeType;
+  typedef typename OutputImageType::PixelType         PixelType;
+  typedef typename OutputImageType::SpacingType       SpacingType;
+
+  /** Determine the image dimension. */
+  itkStaticConstMacro(ImageDimension, unsigned int,
+                      TOutputImage::ImageDimension );
+  itkStaticConstMacro(InputImageDimension, unsigned int,
+                      TInputImage::ImageDimension );
+  itkStaticConstMacro(DeformationFieldDimension, unsigned int,
+                      TDeformationField::ImageDimension );
+  /** typedef for base image type at the current ImageDimension */
+  typedef itk::ImageBase<itkGetStaticConstMacro(ImageDimension)> ImageBaseType;
+
+  /** Deformation field typedef support. */
+  typedef TDeformationField                        DeformationFieldType;
+  typedef typename DeformationFieldType::Pointer   DeformationFieldPointer;
+  typedef typename DeformationFieldType::PixelType DisplacementType;
+
+  /** Interpolator typedef support. */
+  typedef double                                                CoordRepType;
+  typedef itk::InterpolateImageFunction<InputImageType,CoordRepType> InterpolatorType;
+  typedef typename InterpolatorType::Pointer                    InterpolatorPointer;
+  typedef itk::LinearInterpolateImageFunction<InputImageType,CoordRepType>
+                                                                DefaultInterpolatorType;
+
+  /** Point type */
+  typedef itk::Point<CoordRepType,itkGetStaticConstMacro(ImageDimension)> PointType;
+
+  /** Type for representing the direction of the output image */
+  typedef typename TOutputImage::DirectionType     DirectionType;
+
+  /** Set the deformation field. */
+  void SetDeformationField( const DeformationFieldType * field );
+
+  /** Get a pointer the deformation field. */
+  DeformationFieldType * GetDeformationField(void);
+
+  /** Set the interpolator function. */
+  itkSetObjectMacro( Interpolator, InterpolatorType );
+
+  /** Get a pointer to the interpolator function. */
+  itkGetObjectMacro( Interpolator, InterpolatorType );
+
+  /** Set the output image spacing. */
+  itkSetMacro(OutputSpacing, SpacingType);
+  virtual void SetOutputSpacing( const double* values);
+
+  /** Get the output image spacing. */
+  itkGetConstReferenceMacro(OutputSpacing, SpacingType);
+
+  /** Set the output image origin. */
+  itkSetMacro(OutputOrigin, PointType);
+  virtual void SetOutputOrigin( const double* values);
+
+  /** Get the output image origin. */
+  itkGetConstReferenceMacro(OutputOrigin, PointType);
+
+  /** Set/Get the direction (orientation) of the output image */
+  itkSetMacro(OutputDirection, DirectionType );
+  itkGetConstReferenceMacro(OutputDirection, DirectionType );
+
+  /** Helper method to set the output parameters based on this image */
+  void SetOutputParametersFromImage ( const ImageBaseType *image );
+
+  /** Set the start index of the output largest possible region.
+   * The default is an index of all zeros. */
+  itkSetMacro( OutputStartIndex, IndexType );
+
+  /** Get the start index of the output largest possible region. */
+  itkGetConstReferenceMacro( OutputStartIndex, IndexType );
+
+  /** Set the size of the output image. */
+  itkSetMacro( OutputSize, SizeType );
+
+  /** Get the size of the output image. */
+  itkGetConstReferenceMacro( OutputSize, SizeType );
+
+  /** Set the edge padding value */
+  itkSetMacro( EdgePaddingValue, PixelType );
+
+  /** Get the edge padding value */
+  itkGetConstMacro( EdgePaddingValue, PixelType );
+
+  /** WarpImageFilter produces an image which is a different
+   * size than its input image. As such, it needs to provide an
+   * implemenation for GenerateOutputInformation() which set
+   * the output information according the OutputSpacing, OutputOrigin
+   * and the deformation field's LargestPossibleRegion. */
+  virtual void GenerateOutputInformation();
+
+  /** It is difficult to compute in advance the input image region
+   * required to compute the requested output region. Thus the safest
+   * thing to do is to request for the whole input image.
+   *
+   * For the deformation field, the input requested region
+   * set to be the same as that of the output requested region. */
+  virtual void GenerateInputRequestedRegion();
+
+  /** This method is used to set the state of the filter before
+   * multi-threading. */
+  virtual void BeforeThreadedGenerateData();
+
+  /** This method is used to set the state of the filter after
+   * multi-threading. */
+  virtual void AfterThreadedGenerateData();
+
+#ifdef ITK_USE_CONCEPT_CHECKING
+  /** Begin concept checking */
+  itkConceptMacro(SameDimensionCheck1,
+    (itk::Concept::SameDimension<ImageDimension, InputImageDimension>));
+  itkConceptMacro(SameDimensionCheck2,
+    (itk::Concept::SameDimension<ImageDimension, DeformationFieldDimension>));
+  /** itkConceptMacro(InputHasNumericTraitsCheck,
+    (Concept::HasNumericTraits<typename TInputImage::PixelType>));*/
+  itkConceptMacro(DeformationFieldHasNumericTraitsCheck,
+                  (itk::Concept::HasNumericTraits<typename TDeformationField::PixelType::ValueType>));
+  /** End concept checking */
+#endif
+
+protected:
+  WarpImageFilter();
+  virtual ~WarpImageFilter() {};
+  void PrintSelf(std::ostream& os, itk::Indent indent) const;
+
+  /** WarpImageFilter is implemented as a multi-threaded filter.
+   * As such, it needs to provide and implementation for
+   * ThreadedGenerateData(). */
+  void ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread,
+                            int threadId );
+
+private:
+  WarpImageFilter(const Self&); //purposely not implemented
+  void operator=(const Self&); //purposely not implemented
+
+  /** This function should be in an interpolator but none of the ITK
+   * interpolators at this point handle edge conditions properly
+   */
+  DisplacementType EvaluateDeformationAtPhysicalPoint(const PointType &p);
+
+  PixelType                  m_EdgePaddingValue;
+  SpacingType                m_OutputSpacing;
+  PointType                  m_OutputOrigin;
+  DirectionType              m_OutputDirection;
+
+  InterpolatorPointer        m_Interpolator;
+  SizeType                   m_OutputSize;        // Size of the output image
+  IndexType                  m_OutputStartIndex;  // output image start index
+  bool                       m_DefFieldSizeSame;
+  // variables for deffield interpoator
+  IndexType m_StartIndex,m_EndIndex;
+
+
+};
+
+} // end namespace otb
+
+#ifndef OTB_MANUAL_INSTANTIATION
+#include "otbWarpImageFilter.txx"
+#endif
+
+#endif

Code/UtilitiesAdapters/ITKPendingPatches/otbWarpImageFilter.txx

+/*=========================================================================
+
+  Program:   Insight Segmentation & Registration Toolkit
+  Module:    $RCSfile: itkWarpImageFilter.txx,v $
+  Language:  C++
+  Date:      $Date: 2009-10-29 11:19:10 $
+  Version:   $Revision: 1.34 $
+
+  Copyright (c) Insight Software Consortium. All rights reserved.
+  See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
+
+     This software is distributed WITHOUT ANY WARRANTY; without even
+     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+     PURPOSE.  See the above copyright notices for more information.
+
+=========================================================================*/
+#ifndef __otbWarpImageFilter_txx
+#define __otbWarpImageFilter_txx
+#include "otbWarpImageFilter.h"
+
+#include "itkImageRegionIterator.h"
+#include "itkImageRegionIteratorWithIndex.h"
+#include "itkNumericTraits.h"
+#include "itkProgressReporter.h"
+#include "itkContinuousIndex.h"
+#include "vnl/vnl_math.h"
+
+#include "itkVariableLengthVector.h"
+#include "otbPixelBuilder.h"
+
+namespace otb
+{
+
+/**
+ * Default constructor.
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::WarpImageFilter()
+{
+  // Setup the number of required inputs
+  this->SetNumberOfRequiredInputs( 2 );
+
+  // Setup default values
+  m_OutputSpacing.Fill( 1.0 );
+  m_OutputOrigin.Fill( 0.0 );
+  m_OutputDirection.SetIdentity();
+  m_OutputSize.Fill(0);
+  PixelBuilder<PixelType>::Zero(m_EdgePaddingValue,1);
+  m_OutputStartIndex.Fill(0);
+  // Setup default interpolator
+  typename DefaultInterpolatorType::Pointer interp =
+    DefaultInterpolatorType::New();
+
+  m_Interpolator =
+    static_cast<InterpolatorType*>( interp.GetPointer() );
+}
+
+/**
+ * Standard PrintSelf method.
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::PrintSelf(std::ostream& os, itk::Indent indent) const
+{
+
+  Superclass::PrintSelf(os, indent);
+
+  os << indent << "OutputSpacing: " << m_OutputSpacing << std::endl;
+  os << indent << "OutputOrigin: " << m_OutputOrigin << std::endl;
+  os << indent << "OutputDirection: " << m_OutputDirection << std::endl;
+  os << indent << "OutputSize: " << m_OutputSize << std::endl;
+  os << indent << "OutputStartIndex: " << m_OutputStartIndex << std::endl;
+  os << indent << "EdgePaddingValue: "
+     << static_cast<typename itk::NumericTraits<PixelType>::PrintType>(m_EdgePaddingValue)
+     << std::endl;
+  os << indent << "Interpolator: " << m_Interpolator.GetPointer() << std::endl;
+
+}
+
+
+/**
+ * Set the output image spacing.
+ *
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::SetOutputSpacing(
+  const double* spacing)
+{
+  SpacingType s(spacing);
+  this->SetOutputSpacing( s );
+}
+
+
+/**
+ * Set the output image origin.
+ *
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::SetOutputOrigin(
+  const double* origin)
+{
+  PointType p(origin);
+  this->SetOutputOrigin(p);
+}
+
+/** Helper method to set the output parameters based on this image */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::SetOutputParametersFromImage ( const ImageBaseType * image )
+{
+  this->SetOutputOrigin ( image->GetOrigin() );
+  this->SetOutputSpacing ( image->GetSpacing() );
+  this->SetOutputDirection ( image->GetDirection() );
+  this->SetOutputStartIndex ( image->GetLargestPossibleRegion().GetIndex() );
+  this->SetOutputSize ( image->GetLargestPossibleRegion().GetSize() );
+}
+
+/**
+ * Set deformation field as Inputs[1] for this ProcessObject.
+ *
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::SetDeformationField(
+  const DeformationFieldType * field )
+{
+  // const cast is needed because the pipeline is not const-correct.
+  DeformationFieldType * input =
+       const_cast< DeformationFieldType * >( field );
+  this->itk::ProcessObject::SetNthInput( 1, input );
+}
+
+
+/**
+ * Return a pointer to the deformation field.
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+typename WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::DeformationFieldType *
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::GetDeformationField(void)
+{
+  return static_cast<DeformationFieldType *>
+    ( this->itk::ProcessObject::GetInput( 1 ));
+}
+
+
+/**
+ * Setup state of filter before multi-threading.
+ * InterpolatorType::SetInputImage is not thread-safe and hence
+ * has to be setup before ThreadedGenerateData
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::BeforeThreadedGenerateData()
+{
+
+  if( !m_Interpolator )
+    {
+    itkExceptionMacro(<< "Interpolator not set");
+    }
+  DeformationFieldPointer fieldPtr = this->GetDeformationField();
+
+  // Connect input image to interpolator
+  m_Interpolator->SetInputImage( this->GetInput() );
+  typename DeformationFieldType::RegionType defRegion =
+    fieldPtr->GetLargestPossibleRegion();
+  typename OutputImageType::RegionType outRegion =
+    this->GetOutput()->GetLargestPossibleRegion();
+  m_DefFieldSizeSame = outRegion == defRegion;
+  if(!m_DefFieldSizeSame)
+    {
+    m_StartIndex = fieldPtr->GetBufferedRegion().GetIndex();
+    for(unsigned i = 0; i < ImageDimension; i++)
+      {
+      m_EndIndex[i] = m_StartIndex[i] +
+        fieldPtr->GetBufferedRegion().GetSize()[i] - 1;
+      }
+    }
+}
+
+/**
+ * Setup state of filter after multi-threading.
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::AfterThreadedGenerateData()
+{
+  // Disconnect input image from interpolator
+  m_Interpolator->SetInputImage( NULL );
+}
+
+
+template <class TInputImage,class TOutputImage,class TDeformationField>
+typename WarpImageFilter<TInputImage,
+                         TOutputImage,
+                         TDeformationField>::DisplacementType
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::EvaluateDeformationAtPhysicalPoint(const PointType &point)
+{
+  DeformationFieldPointer fieldPtr = this->GetDeformationField();
+  itk::ContinuousIndex<double,ImageDimension> index;
+  fieldPtr->TransformPhysicalPointToContinuousIndex(point,index);
+  unsigned int dim;  // index over dimension
+  /**
+   * Compute base index = closest index below point
+   * Compute distance from point to base index
+   */
+  IndexType baseIndex;
+  IndexType neighIndex;
+  double distance[ImageDimension];
+
+  for( dim = 0; dim < ImageDimension; dim++ )
+    {
+    baseIndex[dim] = itk::Math::Floor< IndexValueType >( index[dim] );
+
+    if( baseIndex[dim] >=  m_StartIndex[dim] )
+      {
+      if( baseIndex[dim] <  m_EndIndex[dim] )
+        {
+        distance[dim] = index[dim] - static_cast<double>( baseIndex[dim] );
+        }
+      else
+        {
+        baseIndex[dim] = m_EndIndex[dim];
+        distance[dim] = 0.0;
+        }
+      }
+    else
+      {
+      baseIndex[dim] = m_StartIndex[dim];
+      distance[dim] = 0.0;
+      }
+    }
+
+  /**
+   * Interpolated value is the weight some of each of the surrounding
+   * neighbors. The weight for each neighbour is the fraction overlap
+   * of the neighbor pixel with respect to a pixel centered on point.
+   */
+  DisplacementType output;
+  output.Fill(0);
+
+  double totalOverlap = 0.0;
+  unsigned int numNeighbors(1 << TInputImage::ImageDimension);
+
+  for( unsigned int counter = 0; counter < numNeighbors; counter++ )
+    {
+    double overlap = 1.0;          // fraction overlap
+    unsigned int upper = counter;  // each bit indicates upper/lower neighbour
+
+    // get neighbor index and overlap fraction
+    for( dim = 0; dim < ImageDimension; dim++ )
+      {
+
+      if( upper & 1 )
+        {
+        neighIndex[dim] = baseIndex[dim] + 1;
+        overlap *= distance[dim];
+        }
+      else
+        {
+        neighIndex[dim] = baseIndex[dim];
+        overlap *= 1.0 - distance[dim];
+        }
+
+      upper >>= 1;
+
+      }
+
+    // get neighbor value only if overlap is not zero
+    if( overlap )
+      {
+      const DisplacementType input =
+        fieldPtr->GetPixel( neighIndex );
+      for(unsigned int k = 0; k < PixelSizeFinder(input); k++ )
+        {
+        output[k] += overlap * static_cast<double>( input[k] );
+        }
+      totalOverlap += overlap;
+      }
+
+    if( totalOverlap == 1.0 )
+      {
+      // finished
+      break;
+      }
+
+    }
+  return ( output );
+}
+
+template <class PixelType> unsigned int PixelSizeFinder(itk::VariableLengthVector<PixelType> pix)
+{
+  return pix.Size();
+}
+template <class PixelType> unsigned int PixelSizeFinder(PixelType pix)
+{
+  return PixelType::Dimension;
+}
+
+/**
+ * Compute the output for the region specified by outputRegionForThread.
+ */
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::ThreadedGenerateData(
+  const OutputImageRegionType& outputRegionForThread,
+  int threadId )
+{
+
+  InputImageConstPointer inputPtr = this->GetInput();
+  OutputImagePointer outputPtr = this->GetOutput();
+  DeformationFieldPointer fieldPtr = this->GetDeformationField();
+
+  // support progress methods/callbacks
+  itk::ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
+
+  // iterator for the output image
+  itk::ImageRegionIteratorWithIndex<OutputImageType> outputIt(
+    outputPtr, outputRegionForThread );
+  IndexType index;
+  PointType point;
+  DisplacementType displacement;
+  if(this->m_DefFieldSizeSame)
+    {
+    // iterator for the deformation field
+    itk::ImageRegionIterator<DeformationFieldType>
+      fieldIt(fieldPtr, outputRegionForThread );
+
+    while( !outputIt.IsAtEnd() )
+      {
+      // get the output image index
+      index = outputIt.GetIndex();
+      outputPtr->TransformIndexToPhysicalPoint( index, point );
+
+      // get the required displacement
+      displacement = fieldIt.Get();
+
+      // compute the required input image point
+      for(unsigned int j = 0; j < ImageDimension; j++ )
+        {
+        point[j] += displacement[j];
+        }
+
+      // get the interpolated value
+      if( m_Interpolator->IsInsideBuffer( point ) )
+        {
+        PixelType value =
+          static_cast<PixelType>(m_Interpolator->Evaluate( point ) );
+        outputIt.Set( value );
+        }
+      else
+        {
+        outputIt.Set( m_EdgePaddingValue );
+        }
+      ++outputIt;
+      ++fieldIt;
+      progress.CompletedPixel();
+      }
+    }
+  else
+    {
+    while( !outputIt.IsAtEnd() )
+      {
+      // get the output image index
+      index = outputIt.GetIndex();
+      outputPtr->TransformIndexToPhysicalPoint( index, point );
+
+      displacement = this->EvaluateDeformationAtPhysicalPoint(point);
+      // compute the required input image point
+      for(unsigned int j = 0; j < ImageDimension; j++ )
+        {
+        point[j] += displacement[j];
+        }
+
+      // get the interpolated value
+      if( m_Interpolator->IsInsideBuffer( point ) )
+        {
+        PixelType value =
+          static_cast<PixelType>(m_Interpolator->Evaluate( point ) );
+        outputIt.Set( value );
+        }
+      else
+        {
+        outputIt.Set( m_EdgePaddingValue );
+        }
+      ++outputIt;
+      progress.CompletedPixel();
+      }
+    }
+}
+
+
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::GenerateInputRequestedRegion()
+{
+
+  // call the superclass's implementation
+  Superclass::GenerateInputRequestedRegion();
+
+  // request the largest possible region for the input image
+  InputImagePointer inputPtr =
+    const_cast< InputImageType * >( this->GetInput() );
+
+  if( inputPtr )
+    {
+    inputPtr->SetRequestedRegionToLargestPossibleRegion();
+    }
+
+  // just propagate up the output requested region for the
+  // deformation field.
+  DeformationFieldPointer fieldPtr = this->GetDeformationField();
+  OutputImagePointer outputPtr = this->GetOutput();
+  if(fieldPtr.IsNotNull() )
+    {
+    fieldPtr->SetRequestedRegion( outputPtr->GetRequestedRegion() );
+    if(!fieldPtr->VerifyRequestedRegion())
+      {
+      fieldPtr->SetRequestedRegion(fieldPtr->GetLargestPossibleRegion());
+      }
+    }
+}
+
+
+template <class TInputImage,class TOutputImage,class TDeformationField>
+void
+WarpImageFilter<TInputImage,TOutputImage,TDeformationField>
+::GenerateOutputInformation()
+{
+  // call the superclass's implementation of this method
+  Superclass::GenerateOutputInformation();
+
+  OutputImagePointer outputPtr = this->GetOutput();
+
+  outputPtr->SetSpacing( m_OutputSpacing );
+  outputPtr->SetOrigin( m_OutputOrigin );
+  outputPtr->SetDirection( m_OutputDirection );
+
+  DeformationFieldPointer fieldPtr = this->GetDeformationField();
+  if( this->m_OutputSize[0] == 0 &&
+      fieldPtr.IsNotNull())
+    {
+    outputPtr->SetLargestPossibleRegion( fieldPtr->
+                                         GetLargestPossibleRegion() );
+    }
+  else
+    {
+    OutputImageRegionType region;
+    region.SetSize(this->m_OutputSize);
+    region.SetIndex(this->m_OutputStartIndex);
+    outputPtr->SetLargestPossibleRegion(region);
+    }
+}
+
+
+} // end namespace otb
+
+#endif