freFieldPrecisionSOMetric.txx

Go to the documentation of this file.

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

  Program:   F.R.E.E. - flexible registration evaluation engine
  Version:   v.0.7.2
  Date:      $Date: 2006/09/01 12:00:00 $
  Module:    $RCSfile: frePrecisionSOMetric.txx,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.

=========================================================================*/
#ifndef __frePrecisionSOMetric_txx
#define __frePrecisionSOMetric_txx

#include "frePrecisionSOMetric.h"
#include "freGenericSetupToImageAdaptor.h"
#include "freVectorFieldSetVarianceCalculator.h"
#include "freVectorToNormImageFilter.h"
#include "freImageSampleCharacteristicsCalculator.h"
#include "freSessionProcessor.h"

#include "itkNumericTraits.h"
#include "itkAddImageFilter.h"

namespace FREE
{
 
  template <unsigned int VImageDimension>
  void
        PrecisionSOMetric<VImageDimension>::
  SetResultFieldPath(const IDPath& path)
        {
    itkDebugMacro("setting ResultFieldPath to " << path.ToStr() );
    this->m_ResultFieldPath = path;
    this->Modified(); \
        };

  template <unsigned int VImageDimension>
  void
        PrecisionSOMetric<VImageDimension>::
  SetInterimFieldPath(const std::string& path)
        {
    itkDebugMacro("setting InterimFieldPath to " << path );
    this->m_InterimFieldPath = path;
    this->Modified();
        };

  template <unsigned int VImageDimension>
  typename PrecisionSOMetric<VImageDimension>::DecomposedMeasureType
        PrecisionSOMetric<VImageDimension>::
	ComputeDecomposedValue( const ParametersType & parameters ) const
        {
                if (!this->m_pAdaptations) throwExceptionMacro("Error. Cannot calculate metric value, adaptiation list is null.");
                if (this->m_Transform.IsNull()) throwExceptionMacro("Error. Cannot calculate metric value, no setup transform set.");

                if (this->m_Transform->GetNumberOfParameters() != parameters.size()) throwExceptionMacro("Error. Cannot calculate metric value, number of passed parameters is not equal to the parameter size of the transform. Passed parameter count: "<< parameters.size());

                m_dMinError = itk::NumericTraits< double >::max();
                m_dMaxError = 0.0;
                m_lFailedRegistrations = 0;
                m_lActAdaptation = 0;
                m_dErrorMean = 0.0;
                m_dErrorVar = 0.0;

                m_Transform->SetParameters(parameters);
                m_smpSetup = m_Transform->GenerateTransformedSetup();

    typedef FREE::VectorFieldSetVarianceCalculator<VImageDimension> EvaluatorType;
    EvaluatorType evaluator;

    evaluator.SetSaveFieldsOnDisk(m_SaveInterimOnDisc);
          evaluator.SetFieldsPath(m_InterimFieldPath);

    GenericSetupToImageAdaptor adaptor;
                adaptor.SetTemplateSetup(m_smpSetup);
                adaptor.AddAdaptations(*m_pAdaptations);

                while (!adaptor.EndOfAdaptation())
                {
      typename Setup::Pointer pSetup = adaptor.AdaptNextSetup();

                        if (this->fnOnNextAdaptation.IsNotNull()) this->fnOnNextAdaptation->Execute(adaptor.GetCurrentAdaptationID(),"",(void*)this);
                        this->InvokeEvent( NextAdaptationObserverEvent() );

                        m_lActAdaptation = adaptor.GetCurrentAdaptationID();

      for (unsigned int iIndex = 0; iIndex<m_SampleSize; iIndex++)
      {
        SessionProcessor processor;

        processor.SetSetup(pSetup);

        processor.DefineOutput("resultField", this->m_ResultFieldPath);

        try
        {
          processor.InitializeSession();
      
          typename TransformationFieldType::Pointer smpResult = processor.GetCastedOutput<TransformationFieldType>("resultField");

                          evaluator.AddVectorField(smpResult);
        }
              catch( FREE::ExceptionBase & err ) 
              { 
                      if (this->fnOnEvaluationFailed.IsNotNull()) this->fnOnEvaluationFailed->Execute(m_lActAdaptation,std::string(err.what()),(void*)this);
                      m_lFailedRegistrations++;
              } 
              catch(...)
              {
                      if (this->fnOnEvaluationFailed.IsNotNull()) this->fnOnEvaluationFailed->Execute(m_lActAdaptation,"unkown error",(void*)this);
                      m_lFailedRegistrations++;
              }
      }

      if (fnOnEvaluationDone.IsNotNull()) fnOnEvaluationDone->Execute(adaptor.GetCurrentAdaptationID(),"",(void*)this);
      this->InvokeEvent( EvaluationDoneObserverEvent() );
                };

                evaluator.Compute();
                typename EvaluatorType::VarianceFieldType::Pointer smpVarField = evaluator.GetVarianceField();

    FREE::ImageSampleCharacteristicsCalculator<typename EvaluatorType::VarianceFieldType> varCalculator;
                varCalculator.SetImage(smpVarField);
                varCalculator.Compute();

                m_dErrorMean = varCalculator.GetMean();
                m_dErrorVar = varCalculator.GetVariance();

          //checking for minimum or maximum
          typedef itk::ImageRegionConstIterator<typename  EvaluatorType::VarianceFieldType > NormIterator;
          NormIterator it(smpVarField,smpVarField->GetBufferedRegion());
          for ( it.GoToBegin(); !it.IsAtEnd(); ++it )
          {
                  if (it.Get()<m_dMinError) m_dMinError = it.Get();
                  if (it.Get()>m_dMaxError) m_dMaxError = it.Get();
          }

          //saving the computed results as measure
    DecomposedMeasureType results(5);
                results.SetElement(0,m_dErrorMean);
                results.SetElement(1,std::sqrt(m_dErrorVar));
                results.SetElement(2,m_dMinError);
                results.SetElement(3,m_dMaxError);
                results.SetElement(4,m_lFailedRegistrations);

                return results;
        };

template <unsigned int VImageDimension>
PrecisionSOMetric<VImageDimension>::
PrecisionSOMetric()
{
  m_IterationEvent = FREE::ProgressEvent<Self>::New(this,&PrecisionSOMetric<VImageDimension>::OnEvaluationProgress);

  m_SaveInterimOnDisc = true;
  m_InterimFieldPath = "";
};

template <unsigned int VImageDimension>
void
PrecisionSOMetric<VImageDimension>::
OnEvaluationProgress(const long status, const std::string& sComment, void* pSender, long threadID)
{
        if (fnOnEvaluationProgress.IsNotNull()) fnOnEvaluationProgress->Execute(status, "", pSender, threadID);
        this->InvokeEvent( itk::IterationEvent() );
}

} // end namespace FREE

#endif

---