Remove throw to avoid warning with c++11

[clitk.git] / itk / clitkBackProjectImageFilter.txx

diff --git a/itk/clitkBackProjectImageFilter.txx b/itk/clitkBackProjectImageFilter.txx
index 65ac8017d4cb98b2009b573f3de03afb07e694db..e6432266300f4e6f9017f557a94bda5accfd7b70 100644 (file)
--- a/itk/clitkBackProjectImageFilter.txx
+++ b/itk/clitkBackProjectImageFilter.txx
@@ -20,6 +20,7 @@
 #include "clitkBackProjectImageFilter.h"
 #include "itkContinuousIndex.h"
 #include "vnl/vnl_math.h"
+#include "itkLinearInterpolateImageFunction.h"
 
 namespace clitk
 {
@@ -39,6 +40,8 @@ namespace clitk
     this->m_SourceToAxis = 1000.0;
     this->m_EdgePaddingValue = itk::NumericTraits<OutputPixelType>::Zero;//density images
     this->m_RigidTransformMatrix.SetIdentity();
+    this->m_PanelShift[0] = 0.;
+    this->m_PanelShift[1] = 0.;
 
     //Parameters for output
     this->m_OutputSpacing.Fill(1);
@@ -187,7 +190,7 @@ namespace clitk
   template <class InputImageType, class OutputImageType>
   void 
   BackProjectImageFilter<InputImageType, OutputImageType>
-  ::Initialize( void ) throw (itk::ExceptionObject)
+  ::Initialize( void )
   {
     //Change the origin of the 2D input
     typename  InputImageType::ConstPointer inputPtr=this->GetInput();
@@ -220,7 +223,7 @@ namespace clitk
   BackProjectImageFilter<InputImageType, OutputImageType>
   ::CalculateProjectionMatrix( void )
   {
-    InputSpacingType inputSpacing=this->GetInput()->GetSpacing();
+    //InputSpacingType inputSpacing=this->GetInput()->GetSpacing();
     
     // Projection on YZ plane+pixelTrans
     itk::Matrix<double,3,4> temp;
@@ -299,13 +302,10 @@ namespace clitk
   //-----------------------------------------------------------------------
   template <class InputImageType, class OutputImageType>
   void 
-  BackProjectImageFilter<InputImageType, OutputImageType>
-  ::ThreadedGenerateData( const OutputImageRegionType & outputRegionForThread,  int threadId )
+  BackProjectImageFilter<InputImageType, OutputImageType>::ThreadedGenerateData( const OutputImageRegionType & outputRegionForThread,  itk::ThreadIdType threadId )
   {
     //Projection pointer
     InputImageConstPointer inputPtr=this->GetInput();
-    InputPixelType * beginPtr=const_cast<InputPixelType *>(this->GetInput()->GetBufferPointer());
-    InputPixelType * pp;
     
     //Volume pointer
     OutputImagePointer outputPTr= this->GetOutput();
@@ -324,9 +324,6 @@ namespace clitk
     iPoint.Fill(itk::NumericTraits<double>::Zero);
     OutputIndexType oIndex;
     ContinuousInputIndexType iIndex;
-    InputSizeType inputSize=inputPtr->GetLargestPossibleRegion().GetSize();
-    double dx,dy,dxm,dym;
-    int lx, ly;
 
     //Get the first output coordinate
     oIndex=iterator.GetIndex();//costly but only once a thread
@@ -337,9 +334,13 @@ namespace clitk
        
     //Compute the first input coordinate (invert Y/X)
     homInputPoint= (m_ProjectionMatrix * homOutputPoint);
-    iPoint[0]=-homInputPoint[0]/homInputPoint[2];
-    iPoint[1]=homInputPoint[1]/homInputPoint[2];
-    
+    iPoint[0]=-homInputPoint[0]/homInputPoint[2] + m_PanelShift[0];
+    iPoint[1]=homInputPoint[1]/homInputPoint[2] + m_PanelShift[1];
+
+    typedef itk::LinearInterpolateImageFunction< InputImageType, double > InterpolatorType;
+    typename InterpolatorType::Pointer interpolator = InterpolatorType::New();
+    interpolator->SetInputImage(this->GetInput());
+
     //Run over all output voxels
     for (unsigned int i=0; i<outputSizeForThread[2]; i++)
       {
@@ -347,25 +348,15 @@ namespace clitk
          {
            for (unsigned int k=0; k<outputSizeForThread[0]; k++)
              {
-               iPoint[0]=homInputPoint[0]/homInputPoint[2];
-               iPoint[1]=homInputPoint[1]/homInputPoint[2];
+               iPoint[0]=-homInputPoint[0]/homInputPoint[2] + m_PanelShift[0];
+               iPoint[1]=homInputPoint[1]/homInputPoint[2] + m_PanelShift[1];
 
                //Check wether inside, convert to index (use modified with correct origin)
-               if( m_ModifiedInput->TransformPhysicalPointToContinuousIndex(iPoint, iIndex) )
-                 {
-                   //Own (fast bilinear) interpolation   
-                   lx = (int)floor(iIndex[0]); dx = iIndex[0]-lx; dxm = 1.-dx;
-                   ly = (int)floor(iIndex[1]); dy = iIndex[1]-ly; dym = 1.-dy;
-                   pp = beginPtr + ly*inputSize[0]+lx;
-                   value =static_cast<OutputPixelType>( ( dxm  * dym*(double)(*pp) 
-                                                          + dx * dym*(double)(*(pp+1)) 
-                                                          + dxm* dy *(double)(*(pp + inputSize[0]))
-                                                          + dx * dy *(double)(*(pp + inputSize[0]+1))) );
-           
-                 }
+               if (m_ModifiedInput->TransformPhysicalPointToContinuousIndex(iPoint, iIndex) && interpolator->IsInsideBuffer(iIndex))
+                    value = interpolator->EvaluateAtContinuousIndex(iIndex);
                //Outside: padding value
-               else value=m_EdgePaddingValue;
-               
+               else
+                  value=m_EdgePaddingValue;
                //Set it
                iterator.Set(value);