#ifndef __cpExtensions__Algorithms__CPRImageFilter__hxx__
#define __cpExtensions__Algorithms__CPRImageFilter__hxx__

#include <itkMinimumMaximumImageCalculator.h>

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
_TCurve* cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
GetCurve( )
{
  return(
    dynamic_cast< _TCurve* >( this->itk::ProcessObject::GetInput( 1 ) )
    );
}

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
const _TCurve* cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
GetCurve( ) const
{
  return(
    dynamic_cast< const _TCurve* >( this->itk::ProcessObject::GetInput( 1 ) )
    );
}

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
void cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
SetCurve( _TCurve* curve )
{
  this->itk::ProcessObject::SetNthInput( 1, curve );
}

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
CPRImageFilter( )
  : Superclass( ),
    m_SliceRadius( 0 )
{
  this->Superclass::SetNumberOfRequiredInputs( 2 );
}

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
~CPRImageFilter( )
{
}

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
void cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
GenerateOutputInformation( )
{
}

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
void cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
GenerateInputRequestedRegion( )
{
  TImage* input = const_cast< TImage* >( this->GetInput( ) );
  if( input != NULL )
    input->SetRequestedRegionToLargestPossibleRegion( );
}

// -------------------------------------------------------------------------
template< class _TImage, class _TCurve >
void cpExtensions::Algorithms::CPRImageFilter< _TImage, _TCurve >::
GenerateData( )
{
  typedef itk::MinimumMaximumImageCalculator< _TImage > _TMinMax;

  auto input = this->GetInput( );
  auto curve = this->GetCurve( );
  auto output = this->GetOutput( );

  // Compute image intensity range
  typename _TMinMax::Pointer minmax = _TMinMax::New( );
  minmax->SetImage( input );
  minmax->Compute( );

  // Main loop
  this->m_Join = TJoinFilter::New( );
  this->m_Slicers.clear( );
  unsigned long N = curve->GetNumberOfPoints( );
  double len = 0;
  typename _TCurve::TPoint p, q;
  for( unsigned long n = 0; n < N; ++n )
  {
    p = curve->GetPoint( n );
    if( n > 0 )
      len += p.EuclideanDistanceTo( q );
    q = p;

    // Prepare slicer
    typename TSlicer::Pointer slicer = TSlicer::New( );
    slicer->SetInput( input );
    slicer->SetDefaultValue( minmax->GetMinimum( ) );
    slicer->SpacingFromMinimumOn( );
    if( this->m_Interpolator.IsNotNull( ) )
      slicer->SetInterpolator( this->m_Interpolator );
    slicer->SizeFromMaximumOff( );
    if( this->m_SliceRadius > double( 0 ) )
    {
      slicer->SizeFromMinimumOff( );
      slicer->SetSize( this->m_SliceRadius );
    }
    else
      slicer->SizeFromMinimumOn( );
    slicer->SetTranslation( p );
    slicer->SetRotation( curve->GetFrame( n ) );
    slicer->Update( );

    this->m_Join->SetInput( n, slicer->GetOutput( ) );
    this->m_Slicers.push_back( slicer );

  } // rof
  this->m_Join->SetSpacing( len / double( N ) );
  this->m_Join->Update( );
  this->GetOutput( )->Graft( this->m_Join->GetOutput( ) );
}

#endif // __cpExtensions__Algorithms__CPRImageFilter__hxx__

// eof - $RCSfile$