cmtkGroupwiseRegistrationFunctionalBase.h

Go to the documentation of this file.

/*
//
//  Copyright 1997-2009 Torsten Rohlfing
//
//  Copyright 2004-2011 SRI International
//
//  This file is part of the Computational Morphometry Toolkit.
//
//  http://www.nitrc.org/projects/cmtk/
//
//  The Computational Morphometry Toolkit is free software: you can
//  redistribute it and/or modify it under the terms of the GNU General Public
//  License as published by the Free Software Foundation, either version 3 of
//  the License, or (at your option) any later version.
//
//  The Computational Morphometry Toolkit is distributed in the hope that it
//  will be useful, but WITHOUT ANY WARRANTY; without even the implied
//  warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//  GNU General Public License for more details.
//
//  You should have received a copy of the GNU General Public License along
//  with the Computational Morphometry Toolkit.  If not, see
//  <http://www.gnu.org/licenses/>.
//
//  $Revision: 2752 $
//
//  $LastChangedDate: 2011-01-17 11:33:31 -0800 (Mon, 17 Jan 2011) $
//
//  $LastChangedBy: torstenrohlfing $
//
*/

#ifndef __cmtkGroupwiseRegistrationFunctionalBase_h_included_
#define __cmtkGroupwiseRegistrationFunctionalBase_h_included_

#include <cmtkconfig.h>

#include <Base/cmtkFunctional.h>
#include <Base/cmtkUniformVolume.h>
#include <Base/cmtkXform.h>

#include <System/cmtkSmartPtr.h>
#include <System/cmtkThreads.h>

#include <vector>

#ifdef CMTK_BUILD_MPI
#  include <mpi.h>
#endif

namespace
cmtk
{


class GroupwiseRegistrationFunctionalBase : 
  public Functional
{
public:
  typedef Functional Superclass;

  typedef GroupwiseRegistrationFunctionalBase Self;

  typedef SmartPointer<Self> SmartPtr;

  class BadXform {};

  GroupwiseRegistrationFunctionalBase();

  virtual ~GroupwiseRegistrationFunctionalBase();

  virtual void SetFreeAndRereadImages( const bool flag = true )
  {
    this->m_FreeAndRereadImages = flag;
  }

  virtual void SetRepeatIntensityHistogramMatching( const bool flag = true )
  {
    this->m_RepeatIntensityHistogramMatching = flag;
    if ( flag )
      {
      this->SetFreeAndRereadImages( false );
      }
  }

  virtual void CreateTemplateGridFromTargets( const std::vector<UniformVolume::SmartPtr>& targets, const int downsample = 0 );

  virtual void CreateTemplateGrid( const DataGrid::IndexType& dims , const UniformVolume::CoordinateVectorType& deltas  );

  virtual void SetTemplateGrid( UniformVolume::SmartPtr& templateGrid , 
                                const int downsample = 1 , const bool useTemplateData = false  );

  virtual UniformVolume::SmartPtr& GetTemplateGrid() 
  { 
    return this->m_TemplateGrid; 
  }

  virtual const UniformVolume* GetTemplateGrid() const
  {
    return this->m_TemplateGrid; 
  }

  virtual void SetTargetImages( std::vector<UniformVolume::SmartPtr>& tImages );

  virtual void GetOriginalTargetImages( std::vector<UniformVolume::SmartPtr>& tImages )
  {
    tImages = this->m_OriginalImageVector;
  }  

  virtual std::vector<UniformVolume::SmartPtr>& GetOriginalTargetImages()
  {
    return this->m_OriginalImageVector;
  }  

  virtual size_t GetNumberOfTargetImages() const
  {
    return this->m_OriginalImageVector.size();
  }  

  virtual UniformVolume::SmartPtr GetOriginalTargetImage( const size_t imageIdx )
  {
    return this->m_OriginalImageVector[imageIdx];
  }  

  virtual const UniformVolume* GetOriginalTargetImage( const size_t imageIdx ) const
  {
    return this->m_OriginalImageVector[imageIdx];
  }  

  virtual void SetGaussianSmoothImagesSigma( const Types::Coordinate gaussianSmoothImagesSigma )
  {
    this->m_GaussianSmoothImagesSigma = gaussianSmoothImagesSigma;
  }

  virtual void SetUserBackgroundValue( const byte value = 0 )
  {
    this->m_UserBackgroundValue = value;
    this->m_UserBackgroundFlag = true;
  }

  virtual void UnsetUserBackgroundValue()
  {
    this->m_UserBackgroundFlag = false;
  }

  virtual void SetForceZeroSum( const bool forceZeroSum = true )
  {
    this->m_ForceZeroSum = forceZeroSum;
  }

  virtual void SetForceZeroSumFirstN( const size_t forceZeroSumFirstN )
  {
    this->m_ForceZeroSumFirstN = forceZeroSumFirstN;
    this->m_ForceZeroSum = this->m_ForceZeroSum || (this->m_ForceZeroSumFirstN>0);
  }

  virtual void SetActiveImagesFromTo( const size_t from, const size_t to )
  {
    this->m_ActiveImagesFrom = from;
    this->m_ActiveImagesTo = to;
  }

  virtual void SetActiveXformsFromTo( const size_t from, const size_t to )
  {
    this->m_ActiveXformsFrom = from;
    this->m_ActiveXformsTo = to;
  }

  virtual void SetProbabilisticSampleDensity( const float density )
  {
    this->m_ProbabilisticSampleDensity = density;
  }

  virtual void SetProbabilisticSampleUpdatesAfter( const int iterations )
  {
    this->m_ProbabilisticSampleUpdatesAfter = iterations;
    this->m_ProbabilisticSampleUpdatesSince = 0;
  }

  void SetXforms( const std::vector<Xform::SmartPtr>& xformVector )
  {
    this->m_XformVector.resize( xformVector.size() );
    for ( size_t i = 0; i < this->m_XformVector.size(); ++i )
      {
      this->m_XformVector[i] = xformVector[i];
      }
  }

  virtual const Xform* GetGenericXformByIndex( const size_t idx ) const
  {
    return this->m_XformVector[idx];
  }

  virtual Xform::SmartPtr GetGenericXformByIndex( const size_t idx )
  {
    return this->m_XformVector[idx];
  }

  virtual const Xform* GetGenericActiveXformByIndex( const size_t idx ) const
  {
    return this->m_XformVector[idx + this->m_ActiveXformsFrom];
  }

  virtual Xform::SmartPtr GetGenericActiveXformByIndex( const size_t idx )
  {
    return this->m_XformVector[idx + this->m_ActiveXformsFrom];
  }

  virtual Types::Coordinate GetParamStep( const size_t idx, const Types::Coordinate mmStep = 1 ) const 
  {
    const size_t xformIdx = idx / this->m_ParametersPerXform;
    if ( (xformIdx >= this->m_ActiveXformsFrom) && (xformIdx < this->m_ActiveXformsTo) )
      {
      return this->m_XformVector[xformIdx]->GetParamStep( idx % this->m_ParametersPerXform, this->m_ImageVector[xformIdx]->Size, mmStep );
      }
    else
      {
      return 0.0;
      }
  }
  
  virtual size_t ParamVectorDim() const 
  { 
    return this->m_ParametersPerXform * this->m_XformVector.size(); 
  }
  
  virtual size_t VariableParamVectorDim() const
  { 
    return this->ParamVectorDim();
  }
  
  virtual void GetParamVector( CoordinateVector& v );

  virtual void SetParamVector( CoordinateVector& v );

  virtual void SetParamVector( CoordinateVector& v, const size_t xformIdx );

  virtual void SetParameter( const size_t param, const Types::Coordinate value );

  virtual void SetParameter( const size_t xform, const size_t param, const Types::Coordinate value );

  virtual Self::ReturnType EvaluateAt ( CoordinateVector& v );

  virtual Self::ReturnType EvaluateWithGradient( CoordinateVector& v, CoordinateVector& g, const Types::Coordinate step = 1 );

  virtual void AllocateStorage();

  void DebugWriteImages();

protected:
  size_t m_NumberOfThreads;

  size_t m_NumberOfTasks;

  bool m_FreeAndRereadImages;

  bool m_ForceZeroSum;
  
  size_t m_ForceZeroSumFirstN;

  size_t m_ActiveImagesFrom;

  size_t m_ActiveImagesTo;

  size_t m_ActiveXformsFrom;

  size_t m_ActiveXformsTo;

  virtual void ForceZeroSumGradient( CoordinateVector& g ) const;

  size_t m_TemplateNumberOfPixels;

  size_t m_TemplateNumberOfSamples;

  UniformVolume::SmartPtr m_TemplateGrid;

  bool m_UseTemplateData;

  std::vector<byte> m_TemplateData;

  std::vector<UniformVolume::SmartPtr> m_ImageVector;

  std::vector<UniformVolume::SmartPtr> m_OriginalImageVector;

  std::vector<Xform::SmartPtr> m_XformVector;

  float m_ProbabilisticSampleDensity;

  std::vector<size_t> m_ProbabilisticSamples;

  int m_ProbabilisticSampleUpdatesAfter;

  int m_ProbabilisticSampleUpdatesSince;

  virtual void UpdateProbabilisticSamples();

  virtual void InterpolateAllImages();

  virtual void InterpolateImage( const size_t idx, byte* const destination ) { UNUSED(idx); UNUSED(destination); } // cannot make this pure virtual because we need to instantiate for affine initialization

  std::vector<byte*> m_Data;

  std::vector<byte> m_TempData;

  Types::Coordinate m_GaussianSmoothImagesSigma;

  static const byte m_PaddingValue = 255;

  byte m_UserBackgroundValue;
  
  bool m_UserBackgroundFlag;

  size_t m_ParametersPerXform;

  bool m_RepeatIntensityHistogramMatching;

  virtual bool Wiggle();

#ifdef CMTK_BUILD_MPI

  int m_RankMPI;

  int m_SizeMPI;
#endif

  virtual UniformVolume::SmartPtr PrepareSingleImage( UniformVolume::SmartPtr& image );

  virtual void PrepareTargetImages();

  virtual const UniformVolume::SmartPtr GetReformattedImage( const UniformVolume::SmartPtr& targetGrid, const size_t idx ) const;

private:
  void CopyTemplateData();

  friend class GroupwiseRegistrationFunctionalAffineInitializer;
};


} // namespace cmtk

#endif // #ifndef __cmtkGroupwiseRegistrationFunctionalBase_h_included_