+// -------------------------------------------------------------------------
+// @author Leonardo Florez-Valencia (florez-l@javeriana.edu.co)
+// -------------------------------------------------------------------------
+
+#ifndef __CPEXTENSIONS__ALGORITHMS__CPRFILTER__HXX__
+#define __CPEXTENSIONS__ALGORITHMS__CPRFILTER__HXX__
+
+#include <itkPoint.h>
+#include <itkMinimumMaximumImageCalculator.h>
+#include <itkNearestNeighborInterpolateImageFunction.h>
+#include <cpExtensions/Algorithms/BezierCurveFunction.h>
+
+// -------------------------------------------------------------------------
+template< class I, class S >
+const typename cpExtensions::Algorithms::CPRFilter< I, S >::
+TAxis* cpExtensions::Algorithms::CPRFilter< I, S >::
+GetAxis( ) const
+{
+}
+
+// -------------------------------------------------------------------------
+template< class I, class S >
+void cpExtensions::Algorithms::CPRFilter< I, S >::
+SetAxis( const TAxis* axis )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class I, class S >
+cpExtensions::Algorithms::CPRFilter< I, S >::
+CPRFilter( )
+ : Superclass( )
+{
+ this->m_Interpolator =
+ itk::NearestNeighborInterpolateImageFunction< I, S >::New( );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class S >
+cpExtensions::Algorithms::CPRFilter< I, S >::
+~CPRFilter( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class I, class S >
+void cpExtensions::Algorithms::CPRFilter< I, S >::
+GenerateData( )
+{
+ typedef itk::Point< S, VDimension > _P;
+ typedef typename _P::VectorType _V;
+ typedef cpExtensions::Algorithms::BezierCurveFunction< _V > _Bezier;
+ typedef typename _Bezier::TFrame _M;
+ typedef itk::MinimumMaximumImageCalculator< I > _MinMax;
+
+ auto image = this->GetInput( );
+ auto axis = this->GetAxis( );
+
+ // 0. Get image properties
+ typename _MinMax::Pointer minmax = _MinMax::New( );
+ minmax->SetImage( image );
+ minmax->Compute( );
+ auto background = minmax->GetMinimum( );
+
+ // 1. Construct Bezier curve
+ typename _Bezier::Pointer bezier = _Bezier::New( );
+ auto vertices = axis->GetVertexList( );
+ for( unsigned int i = 0; i < vertices->Size( ); ++i )
+ {
+ auto index = vertices->GetElement( i );
+ _P point;
+ image->TransformContinuousIndexToPhysicalPoint( index, point );
+ bezier->AddPoint( point.GetVectorFromOrigin( ) );
+
+ } // rof
+
+ // 2. Slice image and prepare join filter
+ this->m_Join = TJoin::New( );
+ unsigned int nSlices = this->m_NumberOfSlices;
+ if( nSlices == 0 )
+ nSlices = bezier->GetNumberOfPoints( );
+ _V vpre;
+ _M mpre;
+ S len = S( 0 );
+ for( unsigned int i = 0; i <= nSlices; ++i )
+ {
+ // Get geometrical data and correct Frenet frame
+ S u = S( i ) / S( nSlices );
+ _V vnex = bezier->Evaluate( u );
+ _M mnex = bezier->EvaluateFrenetFrame( u );
+ if( i > 0 )
+ {
+ len += ( vpre - vnex ).GetNorm( );
+
+ // TODO: Update Frenet frame orientation
+
+ } // fi
+
+ // Slice image
+ typename TSlicer::Pointer slicer = TSlicer::New( );
+ slicer->SetInput( image );
+ slicer->SetRotation( mnex );
+ slicer->SetTranslation( vnex );
+ slicer->SetSize( S( this->m_SliceRadius ) );
+ slicer->SpacingFromMinimumOn( );
+ slicer->SetDefaultValue( background );
+ slicer->SetInterpolator( this->m_Interpolator );
+
+ // Add output to join filter and keep track of the current slicer
+ this->m_Join->SetInput( i, slicer->GetOutput( ) );
+ this->m_Slices.push_back( slicer );
+
+ // Update values
+ vpre = vnex;
+ mpre = mnex;
+
+ } // rof
+
+ // 3. Finish join filter configuration
+ this->m_Join->SetSpacing( double( len / S( nSlices + 1 ) ) );
+ this->m_Join->SetOrigin( double( 0 ) );
+
+ // 4. Graft outputs
+ this->m_Join->GraftOutput( this->GetOutput( ) );
+ this->m_Join->UpdateLargestPossibleRegion( );
+ this->GraftOutput( this->m_Join->GetOutput( ) );
+}
+
+#endif // __CPEXTENSIONS__ALGORITHMS__CPRFILTER__HXX__
+
+// eof - $RCSfile$