freVarianceFieldGenerator.txx

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: freVarianceFieldGenerator.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 __freVarianceFieldGenerator_txx
#define __freVarianceFieldGenerator_txx

#include "freVarianceFieldGenerator.h"

#include "itkCastImageFilter.h"
#include "itkNumericTraits.h"
#include "itkSubtractImageFilter.h"

namespace FREE
{

template <class TInputVectorField,class TOutputField>
void
VarianceFieldGenerator<TInputVectorField, TOutputField>::
ResizeSumField(const InputImageType* elementField)
{
        if (!m_iFieldCount)
        {
                if (!m_bUserRegion)
                {
                        m_OutputRegion = elementField->GetLargestPossibleRegion();
                        m_OutputOrigin = elementField->GetOrigin();
                }

                m_SumField = OutputImageType::New();
                m_SumField->SetOrigin(m_OutputOrigin);
                m_SumField->SetSpacing(elementField->GetSpacing());
                m_SumField->SetRegions(m_OutputRegion);
                m_SumField->Allocate();

                //Reset the variance field to no variance (all values are 0)
    m_SumField->FillBuffer(itk::NumericTraits<OutputValueType>::Zero);
        }
};

template <class TInputVectorField,class TOutputField>
typename VarianceFieldGenerator<TInputVectorField, TOutputField>::InputImagePointer
VarianceFieldGenerator<TInputVectorField, TOutputField>::
AddField(const InputImageType* pField, const InputImageType* pMeanField)
{
        ResizeSumField(pField);

        typedef itk::SubtractImageFilter< InputImageType, InputImageType, InputImageType > SubtractFilterType;
        typename SubtractFilterType::Pointer smpSubFilter = SubtractFilterType::New();
                
        smpSubFilter->SetInput1(pField);
        smpSubFilter->SetInput2(pMeanField);
                
        InputImagePointer smpDVF = smpSubFilter->GetOutput();
        smpSubFilter->Update();

        OutputRegionType newRegion = smpDVF->GetLargestPossibleRegion();

        if (newRegion.Crop(m_OutputRegion))
        { //new field is inside the output region, so add values to the sum, otherwise safe time
                // and skip;
                typedef itk::ImageRegionIterator<OutputImageType> OutputIteratorType;
                typedef itk::ImageRegionConstIterator<InputImageType> InputConstIteratorType;

    OutputIteratorType sumIterator(m_SumField,newRegion);
    InputConstIteratorType iterator(smpDVF,newRegion);

    for (sumIterator.GoToBegin(), iterator.GoToBegin(); !(sumIterator.IsAtEnd()); ++sumIterator, ++iterator)
    {
                        InputValueType dif = iterator.Get();
                        OutputValueType result = dif*dif;
      sumIterator.Set(sumIterator.Get()+result);  
                }
        }

        m_iFieldCount++;

        return smpDVF;
};
  
        /* Function creates an output vector field containing
         * the mean of every vector.*/
template <class TInputVectorField,class TOutputField>
typename VarianceFieldGenerator<TInputVectorField, TOutputField>::OutputImagePointer
VarianceFieldGenerator<TInputVectorField, TOutputField>::
GetOutput (void) const
{
        OutputImagePointer outputField = OutputImageType::New();
        outputField->SetOrigin(m_SumField->GetOrigin());
        outputField->SetSpacing(m_SumField->GetSpacing());
        outputField->SetRegions(m_SumField->GetLargestPossibleRegion());
        outputField->Allocate();

        typedef itk::ImageRegionIterator<OutputImageType> OutputIteratorType;
        typedef itk::ImageRegionConstIterator<OutputImageType> OutputConstIteratorType;

  OutputConstIteratorType sumIterator(m_SumField,m_SumField->GetLargestPossibleRegion());
  OutputIteratorType outIterator(outputField,m_SumField->GetLargestPossibleRegion());

        OutputValueType factor = 0;
        if (m_iFieldCount>1) factor = 1.0/(m_iFieldCount - 1.0);

        for (sumIterator.GoToBegin(),outIterator.GoToBegin(); !(sumIterator.IsAtEnd()); ++sumIterator, ++outIterator)
        {
                outIterator.Set( (sumIterator.Get())*factor);
        }

        return outputField;
};

template <class TInputVectorField,class TOutputField>
void
VarianceFieldGenerator<TInputVectorField, TOutputField>::
Reset()
{
        m_iFieldCount = 0;
        m_bUserRegion = false;
};

template <class TInputVectorField,class TOutputField>
VarianceFieldGenerator<TInputVectorField, TOutputField>::
VarianceFieldGenerator()
{
        Reset();
};

template <class TInputVectorField,class TOutputField>
VarianceFieldGenerator<TInputVectorField, TOutputField>::
~VarianceFieldGenerator()
{
};

template <class TInputVectorField,class TOutputField>
typename VarianceFieldGenerator<TInputVectorField, TOutputField>::OutputRegionType
VarianceFieldGenerator<TInputVectorField, TOutputField>::
GetOutputRegion() const
{
        return m_OutputRegion;
};

template <class TInputVectorField,class TOutputField>
void
VarianceFieldGenerator<TInputVectorField, TOutputField>::
SetOutputGeometry(const OutputRegionType& outputRegion, const OutputPointType& origin)
{
        m_OutputRegion = outputRegion;
        m_OutputOrigin = origin;
        m_bUserRegion = true;
};

}//End of Namespace free

#endif

---