X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=registration%2FclitkBSplineDeformableTransform.txx;h=60a2136d7160c9bbfce3431e376bd2d9232dd08c;hb=afd2671211668e6106886aa59c7bb13a82c48a94;hp=db47430b5a82ab339d4103558a29d923bce7d6c8;hpb=573d80d0f7a17607d2ee883c21c940c0ba020282;p=clitk.git

diff --git a/registration/clitkBSplineDeformableTransform.txx b/registration/clitkBSplineDeformableTransform.txx
index db47430..60a2136 100644
--- a/registration/clitkBSplineDeformableTransform.txx
+++ b/registration/clitkBSplineDeformableTransform.txx
@@ -16,7 +16,7 @@
   - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
 ===========================================================================**/
 #ifndef __clitkBSplineDeformableTransform_txx
-#define __clitkBSplineDeformableTransform_tx
+#define __clitkBSplineDeformableTransform_txx
 #include "clitkBSplineDeformableTransform.h"
 
 //itk
@@ -31,8 +31,7 @@ namespace clitk
 
   // Constructor with default arguments
   template<class TCoordRep, unsigned int NInputDimensions, unsigned int NOutputDimensions>
-  BSplineDeformableTransform<TCoordRep, NInputDimensions, NOutputDimensions>
-  ::BSplineDeformableTransform():Superclass(OutputDimension,0)
+  BSplineDeformableTransform<TCoordRep, NInputDimensions, NOutputDimensions>::BSplineDeformableTransform():Superclass(0)
   {
     unsigned int i;
     
@@ -251,7 +250,7 @@ namespace clitk
 
   // Get the number of parameters
   template<class TCoordRep, unsigned int NInputDimensions, unsigned int NOutputDimensions>
-  unsigned int
+  typename BSplineDeformableTransform<TCoordRep, NInputDimensions, NOutputDimensions>::NumberOfParametersType
   BSplineDeformableTransform<TCoordRep, NInputDimensions, NOutputDimensions>
   ::GetNumberOfParameters(void) const
   {
@@ -576,10 +575,11 @@ namespace clitk
     m_CoefficientImage = m_WrappedImage;
     
     //JV Wrap jacobian into OutputDimension X Vectorial images
-    this->m_Jacobian.set_size( OutputDimension, this->GetNumberOfParameters() );
+    this->m_SharedDataBSplineJacobian.set_size( OutputDimension, this->GetNumberOfParameters() );
 
     // Use memset to set the memory
-    JacobianPixelType * jacobianDataPointer = reinterpret_cast<JacobianPixelType *>(this->m_Jacobian.data_block());
+    JacobianPixelType * jacobianDataPointer = reinterpret_cast<JacobianPixelType *>(this->m_SharedDataBSplineJacobian.data_block());
+
     memset(jacobianDataPointer, 0,  OutputDimension*numberOfPixels*sizeof(JacobianPixelType));
     m_LastJacobianIndex = m_ValidRegion.GetIndex();
     m_NeedResetJacobian = false;
@@ -905,19 +905,10 @@ namespace clitk
   // JV weights are identical as for transformpoint, could be done simultaneously in metric!!!!
   // Compute the Jacobian in one position 
   template<class TCoordRep, unsigned int NInputDimensions, unsigned int NOutputDimensions>
-  const 
-  typename BSplineDeformableTransform<TCoordRep, NInputDimensions, NOutputDimensions>
-  ::JacobianType & 
+  void
   BSplineDeformableTransform<TCoordRep, NInputDimensions, NOutputDimensions>
-  ::GetJacobian( const InputPointType & point ) const
+  ::ComputeJacobianWithRespectToParameters(const InputPointType & point, JacobianType & jacobian) const
   {
-    // Can only compute Jacobian if parameters are set via
-    // SetParameters or SetParametersByValue
-    //     if( m_InputParametersPointer == NULL )
-    //       {
-    // 	itkExceptionMacro( <<"Cannot compute Jacobian: parameters not set" );
-    //       }
-
     if (m_NeedResetJacobian)
       ResetJacobian();
 
@@ -926,30 +917,32 @@ namespace clitk
     //========================================================
 
     // Check if inside mask
-    if(m_Mask &&  !(m_Mask->IsInside(point) ) )
+    if (m_Mask &&  !(m_Mask->IsInside(point) ) )
       {
 	// Outside: no (deformable) displacement
-	return this->m_Jacobian;
-      }	
+        jacobian = m_SharedDataBSplineJacobian;
+	return;
+      }
 
-    //Get index   
+    //Get index
     this->TransformPointToContinuousIndex( point, m_Index );
 
     // NOTE: if the support region does not lie totally within the grid
     // we assume zero displacement and return the input point
     if ( !this->InsideValidRegion( m_Index ) )
       {
-	return this->m_Jacobian;
+        jacobian = m_SharedDataBSplineJacobian;
+	return;
       }
 
     //Compute interpolation weights
-    const WeightsDataType *weights=NULL;
+    const WeightsDataType *weights = NULL;
     m_VectorInterpolator->EvaluateWeightsAtContinuousIndex( m_Index, &weights, m_LastJacobianIndex);
     m_SupportRegion.SetIndex( m_LastJacobianIndex );
 
     //Reset the iterators
     unsigned int j = 0;
-    for ( j = 0; j < OutputDimension; j++ ) 
+    for ( j = 0; j < OutputDimension; j++ )
       m_Iterator[j] = IteratorType( m_JacobianImage[j], m_SupportRegion);
 
     // For each dimension, copy the weight to the support region
@@ -969,11 +962,12 @@ namespace clitk
     m_NeedResetJacobian = true;
 
     // Return the result
-    return this->m_Jacobian;
+    jacobian = m_SharedDataBSplineJacobian;
 
   }
 
 
+
   template<class TCoordRep, unsigned int NInputDimensions, unsigned int NOutputDimensions>
   void
   BSplineDeformableTransform<TCoordRep, NInputDimensions,NOutputDimensions>