freVectorFieldSetVarianceCalculator.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: freVectorFieldSetVarianceCalculator.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 __freVectorFieldSetVarianceCalculator_txx
#define __freVectorFieldSetVarianceCalculator_txx

#include "freExceptions.h"
#include "freVectorFieldSetVarianceCalculator.h"
#include "freMeanVectorFieldGenerator.h"
#include "freVarianceFieldGenerator.h"
#include "freVectorSampleCharacteristicsCalculator.h"

#include "itkSubtractImageFilter.h"
#include "itkImageFileWriter.h"
#include "itkImageFileReader.h"
#include <stdio.h>

namespace FREE
{

template <unsigned int VImageDimension>
void
VectorFieldSetVarianceCalculator<VImageDimension>::
Compute()
{
        // Calculation of MVF
        MeanVectorFieldGenerator<VectorFieldType,VectorFieldType> meanGenerator;
        if (fnOnProgress.IsNotNull()) fnOnProgress->Execute(0, "Calculating mean VF", this);

        for (unsigned int i = 0; i<m_iSetSize; i++)
        {
                if (fnOnNextVectorField.IsNotNull()) fnOnNextVectorField->Execute(i, "", this);
                meanGenerator.AddField(GetVectorField(i));
        };

  m_smpMeanVectorField = meanGenerator.GetOutput();

        // Calculation of variance field and DVFs
        VarianceFieldGenerator<VectorFieldType,VarianceFieldType> varianceGenerator;
        if (fnOnProgress.IsNotNull()) fnOnProgress->Execute(3, "Calculating variance field", this);

        for (unsigned int iIndex = 0; iIndex<m_iSetSize; iIndex++)
        {
                if (fnOnNextVectorField.IsNotNull()) fnOnNextVectorField->Execute(iIndex, "", this);

                VectorFieldPointer smpDVF = varianceGenerator.AddField(GetVectorField(iIndex),m_smpMeanVectorField);

                if (fnOnNewDifferenceVectorField.IsNotNull()) fnOnNewDifferenceVectorField->Execute(iIndex, 0, smpDVF.GetPointer(), this);

                SetVectorField(smpDVF, iIndex);
        };

  m_smpVarianceDifferenceField = varianceGenerator.GetOutput();
};

template <unsigned int VImageDimension>
void
VectorFieldSetVarianceCalculator<VImageDimension>::
SetVectorField(const VectorFieldPointer& smpField, unsigned int iIndex)
{
        if (m_bFieldsOnDisk)
        {
                std::string sFileName = GetVectorFieldName(iIndex);
    try
    {
      typedef itk::ImageFileWriter< VectorFieldType >  FieldWriterType;
                        typename FieldWriterType::Pointer  fieldWriter = FieldWriterType::New();

                  fieldWriter->SetInput( smpField );
                        fieldWriter->SetFileName(sFileName.c_str());
                
                        fieldWriter->Update();
    }
                catchAllNPassMacro("Error; Saving vector field temporarily. File: " << sFileName);
        }
        else
        {
                if (m_VectorFields.size()>iIndex)
                {
                        m_VectorFields[iIndex] = smpField;
                }
                else
                {
                        m_VectorFields.push_back(smpField);
                }
        }
};

template <unsigned int VImageDimension>
typename VectorFieldSetVarianceCalculator<VImageDimension>::VectorFieldPointer
VectorFieldSetVarianceCalculator<VImageDimension>::
GetVectorField(unsigned int iIndex)
{
        VectorFieldPointer smpField;

        if (m_bFieldsOnDisk)
        {
                std::string sFileName = GetVectorFieldName(iIndex);
    try
    {
      typedef itk::ImageFileReader< VectorFieldType >  FieldReaderType;
                        typename FieldReaderType::Pointer  fieldReader = FieldReaderType::New();
                        
                        fieldReader->SetFileName(sFileName.c_str());
                        smpField = fieldReader->GetOutput();
                        
                        fieldReader->Update();
    }
                catchAllNPassMacro("Error; Loading vector field. File: " << sFileName);
        }
        else
        {
                if (m_VectorFields.size()>iIndex)
                {
                        smpField = m_VectorFields[iIndex];
                }
        }
        return smpField;
};

template <unsigned int VImageDimension>
void
VectorFieldSetVarianceCalculator<VImageDimension>::
ClearVectorFields()
{
        if (m_bFieldsOnDisk)
        {
                for (unsigned int i = 0; i<m_iSetSize; i++)
                {
                std::string sFileName = GetVectorFieldName(i);
                        if( remove( sFileName.c_str() ) == -1 ) throwExceptionMacro("Error; Deleting temporary vector field. File: " << sFileName);
                };
        }
        else
        {
                m_VectorFields.clear();
        }
};

template <unsigned int VImageDimension>
std::string
VectorFieldSetVarianceCalculator<VImageDimension>::
GetVectorFieldName(unsigned int iIndex)
{
        std::string result = m_sFieldsPath + "VFSDifference_Smpl_"+Convert::ToStr(iIndex)+".mhd";
        return result;
};

template <unsigned int VImageDimension>
void
VectorFieldSetVarianceCalculator<VImageDimension>::
AddVectorField(VectorFieldType* pVectorField)
{
        SetVectorField(pVectorField,m_iSetSize);
        m_iSetSize++;
};

template <unsigned int VImageDimension>
void
VectorFieldSetVarianceCalculator<VImageDimension>::
Reset()
{
        fnOnProgress = NULL;
        fnOnNextVectorField = NULL;
        fnOnNewDifferenceVectorField = NULL;

        m_sFieldsPath = GetGeneralFREEPath();

        m_smpMeanVectorField = NULL;

        m_smpVarianceDifferenceField = NULL;

        ClearVectorFields();
        m_iSetSize = 0;
        m_bFieldsOnDisk = false;
};

template <unsigned int VImageDimension>
VectorFieldSetVarianceCalculator<VImageDimension>::
VectorFieldSetVarianceCalculator()
{
        m_bFieldsOnDisk = false; //important for a first correct reseting
        Reset();
};

template <unsigned int VImageDimension>
VectorFieldSetVarianceCalculator<VImageDimension>::
~VectorFieldSetVarianceCalculator()
{
        Reset();
};

} //End of namespace FREE

#endif

---