Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes

itk::TensorDerivedMeasurementsFilter< TPixel > Class Template Reference

#include <itkTensorDerivedMeasurementsFilter.h>

List of all members.

Public Types

enum  Measure {
  FA, RA, AD, RD,
  CA, L2, L3
}
typedef itk::DiffusionTensor3D
< TPixel > 
TensorType
typedef itk::Image< TensorType, 3 > TensorImageType
typedef itk::Image< TPixel, 3 > OutputImageType
typedef
TensorDerivedMeasurementsFilter 
Self
typedef SmartPointer< SelfPointer
typedef SmartPointer< const SelfConstPointer
typedef ImageToImageFilter
< TensorType, OutputImageType
SuperClass

Public Member Functions

virtual const char * GetClassName () const
virtual void SetMeasure (Measure _arg)

Static Public Member Functions

static Pointer New ()

Protected Member Functions

 TensorDerivedMeasurementsFilter ()
 ~TensorDerivedMeasurementsFilter ()
void GenerateData ()

Protected Attributes

Measure m_Measure

Detailed Description

template<class TPixel>
class itk::TensorDerivedMeasurementsFilter< TPixel >

Definition at line 30 of file itkTensorDerivedMeasurementsFilter.h.


Member Typedef Documentation

template<class TPixel >
typedef SmartPointer<const Self> itk::TensorDerivedMeasurementsFilter< TPixel >::ConstPointer

Definition at line 54 of file itkTensorDerivedMeasurementsFilter.h.

template<class TPixel >
typedef itk::Image< TPixel, 3 > itk::TensorDerivedMeasurementsFilter< TPixel >::OutputImageType

Definition at line 50 of file itkTensorDerivedMeasurementsFilter.h.

template<class TPixel >
typedef SmartPointer<Self> itk::TensorDerivedMeasurementsFilter< TPixel >::Pointer

Definition at line 53 of file itkTensorDerivedMeasurementsFilter.h.

Definition at line 52 of file itkTensorDerivedMeasurementsFilter.h.

template<class TPixel >
typedef ImageToImageFilter< TensorType, OutputImageType > itk::TensorDerivedMeasurementsFilter< TPixel >::SuperClass

Definition at line 57 of file itkTensorDerivedMeasurementsFilter.h.

template<class TPixel >
typedef itk::Image< TensorType, 3 > itk::TensorDerivedMeasurementsFilter< TPixel >::TensorImageType

Definition at line 49 of file itkTensorDerivedMeasurementsFilter.h.

template<class TPixel >
typedef itk::DiffusionTensor3D<TPixel> itk::TensorDerivedMeasurementsFilter< TPixel >::TensorType

Definition at line 48 of file itkTensorDerivedMeasurementsFilter.h.


Member Enumeration Documentation

template<class TPixel >
enum itk::TensorDerivedMeasurementsFilter::Measure
Enumerator:
FA 
RA 
AD 
RD 
CA 
L2 
L3 

Definition at line 37 of file itkTensorDerivedMeasurementsFilter.h.

  {
    FA,// Fractional anisotropy
    RA,// Relative anisotropy
    AD,// Axial diffusivity
    RD,// Radial diffusivity
    CA, // Clustering anisotropy 1-(λ2+λ3)/(2*λ1)
    L2,
    L3
  };

Constructor & Destructor Documentation

template<class TPixel >
itk::TensorDerivedMeasurementsFilter< TPixel >::TensorDerivedMeasurementsFilter (  ) [protected]

Definition at line 12 of file itkTensorDerivedMeasurementsFilter.txx.

                                                                               : m_Measure(AD)
  {
  }
template<class TPixel >
itk::TensorDerivedMeasurementsFilter< TPixel >::~TensorDerivedMeasurementsFilter (  ) [inline, protected]

Definition at line 70 of file itkTensorDerivedMeasurementsFilter.h.

{};

Member Function Documentation

template<class TPixel >
void itk::TensorDerivedMeasurementsFilter< TPixel >::GenerateData (  ) [protected]

Definition at line 17 of file itkTensorDerivedMeasurementsFilter.txx.

  {
    typename TensorImageType::Pointer tensorImage = static_cast< TensorImageType * >( this->ProcessObject::GetInput(0) );
    typedef ImageRegionConstIterator< TensorImageType > TensorImageIteratorType;
    typedef ImageRegionIterator< OutputImageType > OutputImageIteratorType;

    
    typename OutputImageType::Pointer outputImage = 
        static_cast< OutputImageType * >(this->ProcessObject::GetOutput(0));

    typename TensorImageType::RegionType region = tensorImage->GetLargestPossibleRegion();


    outputImage->SetSpacing( tensorImage->GetSpacing() );   // Set the image spacing
    outputImage->SetOrigin( tensorImage->GetOrigin() );     // Set the image origin
    outputImage->SetDirection( tensorImage->GetDirection() );  // Set the image direction
    outputImage->SetRegions( tensorImage->GetLargestPossibleRegion());
    outputImage->Allocate();

    TensorImageIteratorType tensorIt(tensorImage, tensorImage->GetLargestPossibleRegion());
    OutputImageIteratorType outputIt(outputImage, outputImage->GetLargestPossibleRegion());

    tensorIt.GoToBegin();
    outputIt.GoToBegin();

    while(!tensorIt.IsAtEnd() && !outputIt.IsAtEnd()){

      TensorType tensor = tensorIt.Get();

      switch(m_Measure)
      {
      case FA:
        {
          TPixel diffusionIndex = tensor.GetFractionalAnisotropy();
          outputIt.Set(diffusionIndex);
          break;
        }
      case RA:
        {
          TPixel diffusionIndex = tensor.GetRelativeAnisotropy();
          outputIt.Set(diffusionIndex);
          break;
        }
      case AD:
        {
          // eigenvalues are sorted in ascending order by default because the
          // itk::SymmetricEigenAnalysis defaults are not touched in the tensor implementation
          typename TensorType::EigenValuesArrayType evs;
          tensor.ComputeEigenValues(evs);
          outputIt.Set(evs[2]);
          break;
        }
      case RD:
        {
          // eigenvalues are sorted in ascending order by default because the
          // itk::SymmetricEigenAnalysis defaults are not touched in the tensor implementation
          typename TensorType::EigenValuesArrayType evs;
          tensor.ComputeEigenValues(evs);
          outputIt.Set((evs[0]+evs[1])/2);
          
          break;
        }
      case CA:
        {
          // eigenvalues are sorted in ascending order by default because the
          // itk::SymmetricEigenAnalysis defaults are not touched in the tensor implementation
          typename TensorType::EigenValuesArrayType evs;
          tensor.ComputeEigenValues(evs);
          if (evs[2] == 0)
          {
            outputIt.Set(0);
            break;
          }
          outputIt.Set(1-(evs[0]+evs[1])/(2*evs[2]));
          break;
        }
      case L2:
        {
          // eigenvalues are sorted in ascending order by default because the
          // itk::SymmetricEigenAnalysis defaults are not touched in the tensor implementation
          typename TensorType::EigenValuesArrayType evs;
          tensor.ComputeEigenValues(evs);
          outputIt.Set(evs[1]);
          break;
        }
      case L3:
        {
          // eigenvalues are sorted in ascending order by default because the
          // itk::SymmetricEigenAnalysis defaults are not touched in the tensor implementation
          typename TensorType::EigenValuesArrayType evs;
          tensor.ComputeEigenValues(evs);
          outputIt.Set(evs[0]);
          break;
        }
      }

      ++tensorIt;
      ++outputIt;
    }

    
  }
template<class TPixel >
virtual const char* itk::TensorDerivedMeasurementsFilter< TPixel >::GetClassName (  ) const [virtual]

Runtime information support.

template<class TPixel >
static Pointer itk::TensorDerivedMeasurementsFilter< TPixel >::New (  ) [static]

Method for creation through the object factory.

template<class TPixel >
virtual void itk::TensorDerivedMeasurementsFilter< TPixel >::SetMeasure ( Measure  _arg ) [virtual]

Member Data Documentation

template<class TPixel >
Measure itk::TensorDerivedMeasurementsFilter< TPixel >::m_Measure [protected]

Definition at line 70 of file itkTensorDerivedMeasurementsFilter.h.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Defines