freSingleValuedVnlSOMetricAdaptor.cxx

Go to the documentation of this file.

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

  Program:   F.R.E.E. - flexible registration evaluation engine
  Version:   v.1.0.0
  Date:      $Date: 2006/09/01 12:00:00 $
  Module:    $RCSfile: freSingleValuedVnlSOMetricAdaptor.cxx,v $
  Language:  C++



  Copyright (c) 2007 Ralf o Floca (Department of Medical Informatics,
  Institute for Medical Biometry and Informatics, University of Heidelberg,
  Germany). All rights reserved.
  See FREECopyright.txt or http://www.mi.med.uni-hd.de/free/copyright.htm
  for details.

     This software is distributed WITHOUT ANY WARRANTY; without even 
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 
     PURPOSE.  See the above copyright notices for more information.

=========================================================================*/
#include "freSingleValuedVnlSOMetricAdaptor.h"

namespace FREE
{
  
SingleValuedVnlSOMetricAdaptor 
::SingleValuedVnlSOMetricAdaptor(unsigned int spaceDimension):
  vnl_cost_function(spaceDimension) 
{ 
  m_ScalesInitialized = false;
  m_NegateCostFunction = false;
        m_Reporter = itk::Object::New();
}    

    
void
SingleValuedVnlSOMetricAdaptor
::SetScales(const ScalesType & scales)
{
        if (scales.size()!=this->get_number_of_unknowns()) throwExceptionMacro("Error: number of parameter scales ("<<scales.size()<<") doesn't equal the defined number of parameters ("<<this->get_number_of_unknowns()<<").")
  m_Scales = scales;
  m_ScalesInitialized = true;
}   

SingleValuedVnlSOMetricAdaptor::ParametersType
SingleValuedVnlSOMetricAdaptor::
ConvertInternalToParameter(const InternalParametersType & inparameters)
{
        // Use scales if they are provided
        ParametersType parameters(inparameters.size());
        if(m_ScalesInitialized)
        {
                for(unsigned int i=0;i<parameters.size();i++)
                {
                        parameters[i] = inparameters[i]/m_Scales[i];
                }
        }
        else
        {
                parameters.SetData(const_cast<double*>(inparameters.data_block()));
        }

        return parameters;
};

SingleValuedVnlSOMetricAdaptor::InternalParametersType
SingleValuedVnlSOMetricAdaptor::
ConvertParameterToInternal(const ParametersType & parameters)
{
        // Use scales if they are provided
        InternalParametersType inparameters(parameters.Size());
        
        for(unsigned int i=0;i<parameters.size();i++)
        {
                if(m_ScalesInitialized)
                {
                        inparameters[i] = parameters[i]*m_Scales[i];
                }
                else
                {
                        inparameters[i] = parameters[i];
                }
        }

        return parameters;
};

SingleValuedVnlSOMetricAdaptor::InternalMeasureType 
SingleValuedVnlSOMetricAdaptor
::f( const InternalParametersType & inparameters )
{
        if( !m_CostFunction )
        {
                itkGenericExceptionMacro(<<"Attempt to use a SingleValuedVnlSOMetricAdaptor without any CostFunction plugged in");
        }

        ParametersType parameters = ConvertInternalToParameter(inparameters);

        m_CachedValue = m_CostFunction->GetValue( parameters );
        m_CachedCurrentParameters = parameters;
        m_CachedDecomposedValue = this->m_CostFunction->GetCurrentDecomposedValue();

        if( m_NegateCostFunction )
        {
                m_CachedValue = -1.0*m_CachedValue;

                for (unsigned int index=0; index < m_CachedDecomposedValue.GetSize(); index++)
                {
                        m_CachedDecomposedValue[index] = -1.0*m_CachedDecomposedValue[index];
                }
        }

        InternalMeasureType value = static_cast<InternalMeasureType>(m_CachedValue);

        // Notify observers. This is used for overcoming the limitaion of VNL
        // optimizers of not providing callbacks per iteration.
        this->ReportIteration( itk::FunctionEvaluationIterationEvent() ); 

        return value;
}
  


void 
SingleValuedVnlSOMetricAdaptor
::gradf(  const InternalParametersType   & inparameters,
          InternalDerivativeType   & gradient       )
{
  if( !m_CostFunction )
    {
    itkGenericExceptionMacro("Attempt to use a SingleValuedVnlSOMetricAdaptor without any CostFunction plugged in");
    }

        ParametersType parameters = ConvertInternalToParameter(inparameters);
    
  m_CostFunction->GetDerivative( parameters, m_CachedDerivative );
  this->ConvertExternalToInternalGradient( m_CachedDerivative, gradient);
 
  // Notify observers. This is used for overcoming the limitaion of VNL
  // optimizers of not providing callbacks per iteration.
  // Note that m_CachedDerivative is already loaded in the GetDerivative() above.
  m_CachedCurrentParameters = parameters;
        this->ReportIteration( itk::GradientEvaluationIterationEvent() ); 
 
}
  


void 
SingleValuedVnlSOMetricAdaptor
::compute( const InternalParametersType   & x,
           InternalMeasureType      * f, 
           InternalDerivativeType   * g   ) 
{
  // delegate the computation to the CostFunction
  MeasureType    measure;

        ParametersType parameters = ConvertInternalToParameter(x);
  
  m_CostFunction->GetValueAndDerivative( parameters, measure, m_CachedDerivative );

  if( g ) // sometimes Vnl doesn't pass a valid pointer
    {
    this->ConvertExternalToInternalGradient( m_CachedDerivative, *g );
    }

        if( m_NegateCostFunction )
        {
                measure = -1.0*measure;
        }

        m_CachedValue = measure;
  m_CachedCurrentParameters = parameters;
        m_CachedDecomposedValue = this->m_CostFunction->GetCurrentDecomposedValue();

        if( m_NegateCostFunction )
        {
                m_CachedValue = -1.0*m_CachedValue;

                for (unsigned int index=0; index < m_CachedDecomposedValue.GetSize(); index++)
                {
                        m_CachedDecomposedValue[index] = -1.0*m_CachedDecomposedValue[index];
                }
        }

        if( f )
  {
    *f = static_cast<InternalMeasureType>( m_CachedValue );  
  }

        // Notify observers. This is used for overcoming the limitaion of VNL
  // optimizers of not providing callbacks per iteration.
  // Note that m_CachedDerivative is already loaded in the GetDerivative() above.
        this->ReportIteration( itk::FunctionAndGradientEvaluationIterationEvent() ); 
 
}
  

void 
SingleValuedVnlSOMetricAdaptor
::ConvertExternalToInternalGradient( const DerivativeType   & input,
                                     InternalDerivativeType & output ) const
{
  const unsigned int size = input.size();
  output = InternalDerivativeType(size);
  for( unsigned int i=0; i<size; i++ ) 
    {
    if( !m_NegateCostFunction )
      {
      output[i] = input[i];
      }
    else
      {
      output[i] = -input[i];
      }
    }
}


void 
SingleValuedVnlSOMetricAdaptor
::SetNegateCostFunction( bool flag )
{
  m_NegateCostFunction = flag;
}
                                   
bool 
SingleValuedVnlSOMetricAdaptor
::GetNegateCostFunction() const
{
  return m_NegateCostFunction;
}

void 
SingleValuedVnlSOMetricAdaptor
::ReportIteration( const itk::EventObject & event ) const
{
  this->m_Reporter->InvokeEvent( event );
}
 
unsigned long 
SingleValuedVnlSOMetricAdaptor
::AddObserver(const itk::EventObject & event, itk::Command * command) const
{
  return m_Reporter->AddObserver( event, command );
}


const SingleValuedVnlSOMetricAdaptor::MeasureType &
SingleValuedVnlSOMetricAdaptor
::GetCachedValue() const
{
  return m_CachedValue;
}

const SingleValuedVnlSOMetricAdaptor::DecomposedMeasureType &
SingleValuedVnlSOMetricAdaptor
::GetCachedDecomposedValue() const
{
  return m_CachedDecomposedValue;
}

const SingleValuedVnlSOMetricAdaptor::DerivativeType &
SingleValuedVnlSOMetricAdaptor
::GetCachedDerivative() const
{
  return m_CachedDerivative;
}

const SingleValuedVnlSOMetricAdaptor::ParametersType &
SingleValuedVnlSOMetricAdaptor
::GetCachedCurrentParameters() const
{
  return m_CachedCurrentParameters;
}

} // end namespace itk

---