1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102
|
/*=========================================================================
Program: Insight Segmentation & Registration Toolkit
Module: $RCSfile: itkQuaternionRigidTransformGradientDescentOptimizer.cxx,v $
Language: C++
Date: $Date: 2007-03-29 19:37:01 $
Version: $Revision: 1.16 $
Copyright (c) Insight Software Consortium. All rights reserved.
See ITKCopyright.txt or http://www.itk.org/HTML/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.
=========================================================================*/
#ifndef _itkQuaternionRigidTransformGradientDescentOptimizer_txx
#define _itkQuaternionRigidTransformGradientDescentOptimizer_txx
#include "itkQuaternionRigidTransformGradientDescentOptimizer.h"
#include "vnl/vnl_quaternion.h"
#include "itkEventObject.h"
namespace itk
{
/**
* Advance one Step following the gradient direction
*/
void
QuaternionRigidTransformGradientDescentOptimizer
::AdvanceOneStep( void )
{
double direction;
if( m_Maximize )
{
direction = 1.0;
}
else
{
direction = -1.0;
}
ScalesType scales = this->GetScales();
const unsigned int spaceDimension = m_CostFunction->GetNumberOfParameters();
// Make sure the scales have been set
if (scales.size() != spaceDimension)
{
itkExceptionMacro(<< "The size of Scales is "
<< scales.size()
<< ", but the NumberOfParameters is "
<< spaceDimension
<< ".");
}
DerivativeType transformedGradient( spaceDimension);
for ( unsigned int i=0; i< spaceDimension; i++)
{
transformedGradient[i] = m_Gradient[i] / scales[i];
}
ParametersType currentPosition = this->GetCurrentPosition();
// compute new quaternion value
vnl_quaternion<double> newQuaternion;
for ( unsigned int j=0; j < 4; j++ )
{
newQuaternion[j] = currentPosition[j] + direction * m_LearningRate *
transformedGradient[j];
}
newQuaternion.normalize();
ParametersType newPosition( spaceDimension );
// update quaternion component of currentPosition
for ( unsigned int j=0; j < 4; j++ )
{
newPosition[j] = newQuaternion[j];
}
// update the translation component
for (unsigned int j=4; j< spaceDimension; j++)
{
newPosition[j] = currentPosition[j] +
direction * m_LearningRate * transformedGradient[j];
}
// First invoke the event, so the current position
// still corresponds to the metric values.
this->InvokeEvent( IterationEvent() );
this->SetCurrentPosition( newPosition );
}
} // end namespace itk
#endif
|