zxhMetricFFDBase.h

/*=========================================================================

  Program:   ZXH Registration Software
  Author:    Xiahai Zhuang
  Module:    $RCSfle: zxhMetricFFDBase.h    $
  Language:  C++
  Date:      $Date: From  2004-01 $
  Version:   $Revision: 1.0, 2.0, 2.1 $

=========================================================================*/

#ifndef zxhMetricFFDBase_h
#define zxhMetricFFDBase_h

#include "zxhTransformFFDBase.h"
#include "zxhMetricBase.h"
#include "zxhImageModelingLinear.h"


class zxhMetricFFDBase 
{
public:
    zxhMetricFFDBase(void);
    ~zxhMetricFFDBase(void);
    virtual zxhMetricFFDBase*Clone(zxhMetricFFDBase*&p);
    virtual std::string GetMetricType() {return "_FFD_"; } ;

    virtual zxhlfloat  GetBendingEnergy(zxhTransformFFDBase*pTrans) ;
    virtual bool ComputeCtrPntOnTargetForceForBendingEnergy(zxhTransformFFDBase*pTrans, zxhTransformFFDBase * gradient) ;

  
    virtual int GetValueByFFDGrid(zxhlfloat&fValue,int aiCtrPnt[ImageDimensionMax])=0;
    
    virtual int ComputeConstanceForFFDGradient()=0  ;

    virtual int GetDiffByFFDGrid( float * afPD, // output partial differentiation 
                                  zxhTransformFFDBase * pFFDGrid,   
                                  int * aiCtrPnt,
                                  zxhTransformFFDBase * pDifferential ) ;
    
    virtual int GetDiffByFFDGridOnFFDPhysicalVector( float * afPD, // output partial differentiation 
                                          zxhTransformFFDBase * pFFDGrid,
                                          int * aiCtrPnt,
                                          zxhTransformFFDBase * pDifferential,
                                          float * vector); 

    virtual int GetBendingEnergyDerivative( float * afPD, int * aiCtrPnt)
    {
        if( m_pCtrPntDerivativeBendingEnergyConst==0 ||
            m_pCtrPntDerivativeBendingEnergyConst->GetCtrPntsWithCheck(aiCtrPnt[0], aiCtrPnt[1], aiCtrPnt[2], aiCtrPnt[3])==0)
            return -1 ;
        for( int id=0; id<m_pCtrPntDerivativeBendingEnergyConst->GetDimension(); ++id )
        {
            afPD[id] = m_pCtrPntDerivativeBendingEnergyConst->GetCtrPntsWithoutCheck(aiCtrPnt[0], aiCtrPnt[1], aiCtrPnt[2], aiCtrPnt[3])[id];
        }
        return 1;
    };

    virtual void SetPreComputeFFDsCurrPosition( float * afCurrCtrPntWorld, float *afFiniteDiffPosition) ;
    virtual void SetWeightFFDBendingEnergy(float f)     {m_fWeightFFDBendingEnergy=f;};
    virtual float GetWeightFFDBendingEnergy(void)       {return m_fWeightFFDBendingEnergy;};

    virtual void    SetUsingFFDsSubVolume(bool b)           {m_bFFDsSubVolume=b;};
    virtual bool    GetUsingFFDsSubVolume(void)             {return m_bFFDsSubVolume; };

    void GetCtrPntDerivativeConst(zxhTransformFFDBase&ret)  const { zxhTransformBase*p=&ret ; m_pCtrPntDerivativeConst->Clone(p) ; } ;
    virtual zxhTransformFFDBase * GetCtrPntDerivativeConst()        {return m_pCtrPntDerivativeConst ; } ;
protected:

    virtual void    GetFFDControlPointSupportVolumeTest(float *afCtrVolumeTest, int aiCtrPnt[] ) ;

    virtual void TransformPointToWorldPlusDisplacementAdvanceByBsplineFFDs( zxhImageDataT<float>* aPreComputeTransformPointToWorld, int dimension,
                                                      int x,int y,int z,int t,  float fux,float fuy,float fuz, float *pWorldRef )
    {
        float v= zxh::BSpline(fux)*zxh::BSpline(fuy)*zxh::BSpline(fuz);

        for( int id=0; id<dimension; ++id )
            pWorldRef[id] = aPreComputeTransformPointToWorld[id].GetPixelGreyscaleClosest( x,y,z,t )
                         +  v*m_afFiniteDiffPosition[id]; //m_afFiniteDiffPosition has been projected from physical to world coordinate
    };
    // Mainly for non-MI metric where sample points can be off image-grid
    virtual void TransformPointToWorldPlusDisplacementAdvanceByBsplineFFDs( zxhImageDataT<float>* aPreComputeTransformPointToWorld, int dimension,
                                         float x,float y,float z,float t,   float fux,float fuy,float fuz, float *pWorldRef )
    {
        float v=zxh::BSpline(fux)*zxh::BSpline(fuy)*zxh::BSpline(fuz);
        zxhImageModelingLinearT<float,float> mod ;
        for( int id=0; id<dimension; ++id )
        {
            mod.SetImage( &aPreComputeTransformPointToWorld[id] ) ;
            pWorldRef[id] = mod.GetPixelFloatValueWithCheck( x,y,z,t ) + v*m_afFiniteDiffPosition[id];//m_afFiniteDiffPosition has been projected from physical to world coordinate
        }
    };

    virtual void DisplacementAdvanceByBsplineFFDsUsingPreComputedWorld( int dimension, float ux,float uy,float uz, float *pAdWorld )
    {
        float v=zxh::BSpline(ux)*zxh::BSpline(uy)*zxh::BSpline(uz);
        for( int id=0; id<dimension; ++id )
            pAdWorld[id] = v *m_afFiniteDiffPosition[id];//m_afFiniteDiffPosition has been projected from physical to world coordinate
    };
protected:
    virtual int AddFiniteDiffPenaltyByFFDGrid(float * afPD, // output partial differentiation
                                      zxhTransformFFDBase * pFFDGrid,
                                      int * aiCtrPnt,
                                      zxhTransformFFDBase * pDifferential );
    virtual bool AddAllCtrPntOnTargetForceForBendingEnergy() ;

    virtual void InitMetricFFDConstant()
    {
        // new the constant variable for storage
        if( m_pCtrPntDerivativeConst ) { delete m_pCtrPntDerivativeConst ;m_pCtrPntDerivativeConst=0;};
        zxhTransformBase*ptmptrans = m_pCtrPntDerivativeConst ;
        m_pCtrPntDerivativeConst = (zxhTransformFFDBase*)((zxhTransformFFDBase*)m_pMetricBase->GetTransform())->zxhTransformFFDBase::Clone( ptmptrans ) ;
        m_pCtrPntDerivativeConst->SetTransformPara( 0 ) ;
        if( m_pCtrPntDerivativeBendingEnergyConst ){delete m_pCtrPntDerivativeBendingEnergyConst ;m_pCtrPntDerivativeBendingEnergyConst=0;};
        ptmptrans = m_pCtrPntDerivativeBendingEnergyConst ;
        m_pCtrPntDerivativeBendingEnergyConst = (zxhTransformFFDBase*)((zxhTransformFFDBase*)m_pMetricBase->GetTransform())->zxhTransformFFDBase::Clone( ptmptrans ) ;
        m_pCtrPntDerivativeBendingEnergyConst->SetTransformPara( 0 ) ;
    }

    /*  variables */
protected:
    virtual void    SetFFDMetricBase()=0;
    zxhMetricBase*          m_pMetricBase;

    zxhTransformFFDBase*    m_pCtrPntDerivativeConst;
    zxhTransformFFDBase*    m_pCtrPntDerivativeBendingEnergyConst;
    float   m_fWeightFFDBendingEnergy ;

    float                   m_afCurrFFDsCtrPntWorld[ImageDimensionMax];
    float                   m_afFiniteDiffPosition[ImageDimensionMax] ;
    bool                    m_bFFDsSubVolume ;
};

#endif