diff --git a/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.cpp b/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.cpp
index 58ca994..0b11005 100644
--- a/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.cpp
+++ b/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.cpp
@@ -1,459 +1,468 @@
 /*===================================================================
 
 The Medical Imaging Interaction Toolkit (MITK)
 
 Copyright (c) German Cancer Research Center.
 
 All rights reserved.
 
 This software is distributed WITHOUT ANY WARRANTY; without
 even the implied warranty of MERCHANTABILITY or FITNESS FOR
 A PARTICULAR PURPOSE.
 
 See LICENSE.txt or http://www.mitk.org for details.
 
 ===================================================================*/
 #include "itkTractParcellationFilter.h"
 
 // VTK
 #include <vtkPolyLine.h>
 #include <vtkCellArray.h>
 #include <vtkCellData.h>
 #include <boost/progress.hpp>
 #include <mitkDiffusionFunctionCollection.h>
 #include <itkTractDensityImageFilter.h>
 #include <itkImageRegionIterator.h>
 #include <mitkTractClusteringFilter.h>
 #include <mitkClusteringMetricEuclideanStd.h>
 #include <mitkClusteringMetricEuclideanMean.h>
 #include <mitkClusteringMetricLength.h>
 #include <itkImageRegionConstIterator.h>
 #include <itkTractsToVectorImageFilter.h>
 #include <mitkTractometry.h>
 
 namespace itk{
 
 template< class OutImageType, class InputImageType >
 TractParcellationFilter< OutImageType, InputImageType >::TractParcellationFilter()
   : m_UpsamplingFactor(1)
   , m_NumParcels(0)
   , m_NumCentroids(0)
   , m_StartClusterSize(5)
   , m_InputImage(nullptr)
 {
   this->SetNumberOfRequiredOutputs(2);
 }
 
 template< class OutImageType, class InputImageType >
 TractParcellationFilter< OutImageType, InputImageType >::~TractParcellationFilter()
 {
 }
 
 template< class OutImageType, class InputImageType >
 bool TractParcellationFilter< OutImageType, InputImageType >::Flip(vtkSmartPointer< vtkPolyData > polydata1, int i, vtkSmartPointer< vtkPolyData > ref_poly, int ref_i)
 {
   double d_direct = 0;
   double d_flipped = 0;
 
   vtkCell* cell1 = polydata1->GetCell(ref_i);
   if (ref_poly!=nullptr)
     cell1 = ref_poly->GetCell(ref_i);
   auto numPoints1 = cell1->GetNumberOfPoints();
   vtkPoints* points1 = cell1->GetPoints();
 
   std::vector<itk::Point<double, 3>> ref_points;
   for (int j=0; j<numPoints1; ++j)
   {
     double* p1 = points1->GetPoint(j);
     itk::Point<double, 3> itk_p;
     itk_p[0] = p1[0];
     itk_p[1] = p1[1];
     itk_p[2] = p1[2];
     ref_points.push_back(itk_p);
   }
 
   vtkCell* cell2 = polydata1->GetCell(i);
   vtkPoints* points2 = cell2->GetPoints();
 
   for (int j=0; j<numPoints1; ++j)
   {
     auto p1 = ref_points.at(j);
 
     double* p2 = points2->GetPoint(j);
     d_direct += (p1[0]-p2[0])*(p1[0]-p2[0]) + (p1[1]-p2[1])*(p1[1]-p2[1]) + (p1[2]-p2[2])*(p1[2]-p2[2]);
 
     double* p3 = points2->GetPoint(numPoints1-j-1);
     d_flipped += (p1[0]-p3[0])*(p1[0]-p3[0]) + (p1[1]-p3[1])*(p1[1]-p3[1]) + (p1[2]-p3[2])*(p1[2]-p3[2]);
   }
 
   if (d_direct>d_flipped)
     return true;
   return false;
 }
 
 template< class OutImageType, class InputImageType >
 mitk::FiberBundle::Pointer TractParcellationFilter< OutImageType, InputImageType >::GetWorkingFib()
 {
   mitk::FiberBundle::Pointer fib_static_resampled = m_InputTract->GetDeepCopy();
   mitk::Tractometry::ResampleIfNecessary(fib_static_resampled, m_NumParcels);
 
   // clustering
   std::vector< mitk::ClusteringMetric* > metrics;
   metrics.push_back({new mitk::ClusteringMetricEuclideanStd()});
 //  metrics.push_back({new mitk::ClusteringMetricEuclideanMean()});
 //  metrics.push_back({new mitk::ClusteringMetricLength()});
 
   if (m_NumCentroids>0)
   {
     std::vector<mitk::FiberBundle::Pointer> centroids;
     std::shared_ptr< mitk::TractClusteringFilter > clusterer = std::make_shared<mitk::TractClusteringFilter>();
     int c=0;
     while (c<30 && (centroids.empty() || centroids.size()>static_cast<unsigned long>(m_NumCentroids)))
     {
       float cluster_size = m_StartClusterSize + m_StartClusterSize*0.2*c;
       float max_d = 0;
       int i=1;
       std::vector< float > distances;
       while (max_d < m_InputTract->GetGeometry()->GetDiagonalLength()/2)
       {
         distances.push_back(cluster_size*i);
         max_d = cluster_size*i;
         ++i;
       }
 
       clusterer->SetDistances(distances);
       clusterer->SetTractogram(fib_static_resampled);
       clusterer->SetMetrics(metrics);
       clusterer->SetMergeDuplicateThreshold(cluster_size);
       clusterer->SetDoResampling(false);
       clusterer->SetNumPoints(m_NumParcels);
       //  clusterer->SetMaxClusters(m_Controls->m_MaxCentroids->value());
       clusterer->SetMinClusterSize(1);
       clusterer->Update();
       centroids = clusterer->GetOutCentroids();
       ++c;
     }
 
     mitk::FiberBundle::Pointer centroid_bundle = mitk::FiberBundle::New();
     centroid_bundle = centroid_bundle->AddBundles(centroids);
     return centroid_bundle;
   }
 
   return fib_static_resampled;
 }
 
 template< class OutImageType, class InputImageType >
 std::vector< typename itk::Image<unsigned char, 3>::Pointer > TractParcellationFilter< OutImageType, InputImageType >::GetBinarySplit(typename OutImageType::Pointer inImage)
 {
   std::vector< typename itk::Image<unsigned char, 3>::Pointer > binary_maps;
 
   for (unsigned int i=0; i<m_NumParcels; ++i)
   {
     typename itk::Image<unsigned char, 3>::Pointer parcel_image = itk::Image<unsigned char, 3>::New();
     parcel_image->SetSpacing( inImage->GetSpacing() );
     parcel_image->SetOrigin( inImage->GetOrigin() );
     parcel_image->SetDirection( inImage->GetDirection() );
     parcel_image->SetRegions( inImage->GetLargestPossibleRegion() );
     parcel_image->Allocate();
     parcel_image->FillBuffer(0);
 
     binary_maps.push_back(parcel_image);
   }
 
   itk::ImageRegionIterator< itk::Image<unsigned char, 3> > p_it(inImage, inImage->GetLargestPossibleRegion());
   while(!p_it.IsAtEnd())
   {
     if (p_it.Get()>0)
     {
       binary_maps.at(p_it.Get()-1)->SetPixel(p_it.GetIndex(), 1);
     }
     ++p_it;
   }
 
   return binary_maps;
 }
 
 template< class OutImageType, class InputImageType >
 typename OutImageType::Pointer TractParcellationFilter< OutImageType, InputImageType >::PostprocessParcellation(typename OutImageType::Pointer inImage)
 {
   itk::ImageRegionConstIterator< OutImageType > in_it(inImage, inImage->GetLargestPossibleRegion());
 
   typename OutImageType::Pointer outImage = OutImageType::New();
   outImage->SetSpacing( inImage->GetSpacing() );
   outImage->SetOrigin( inImage->GetOrigin() );
   outImage->SetDirection( inImage->GetDirection() );
   outImage->SetRegions( inImage->GetLargestPossibleRegion() );
   outImage->Allocate();
   outImage->FillBuffer(0);
   itk::ImageRegionIterator< OutImageType > out_it(outImage, outImage->GetLargestPossibleRegion());
 
   MITK_INFO << "Postprocessing parcellation";
   while( !in_it.IsAtEnd() )
   {
     if (in_it.Get()>0)
     {
       std::vector< OutPixelType > vote; vote.resize(m_NumParcels + 1, 0);
 
       typename OutImageType::SizeType regionSize;
       regionSize[0] = 3;
       regionSize[1] = 3;
       regionSize[2] = 3;
 
       typename OutImageType::IndexType regionIndex = in_it.GetIndex();
       regionIndex[0] -= 1;
       regionIndex[1] -= 1;
       regionIndex[2] -= 1;
 
       typename OutImageType::RegionType region;
       region.SetSize(regionSize);
       region.SetIndex(regionIndex);
 
       itk::ImageRegionConstIterator<OutImageType> rit(inImage, region);
       while( !rit.IsAtEnd() )
       {
         if (rit.GetIndex()!=regionIndex)
           vote[rit.Get()]++;
         ++rit;
       }
 
       OutPixelType max = 0;
       unsigned int max_parcel = -1;
       for (unsigned int i=1; i<vote.size(); ++i)
       {
         if (vote[i]>max)
         {
           max = vote[i];
           max_parcel = i;
         }
       }
 
       out_it.Set(max_parcel);
     }
     ++out_it;
     ++in_it;
   }
   MITK_INFO << "DONE";
   return outImage;
 }
 
 template< class OutImageType, class InputImageType >
 void TractParcellationFilter< OutImageType, InputImageType >::StaticResampleParcelVoting(typename OutImageType::Pointer outImage)
 {
   typename itk::TractDensityImageFilter< OutImageType >::Pointer generator = itk::TractDensityImageFilter< OutImageType >::New();
   generator->SetFiberBundle(m_InputTract);
   generator->SetMode(TDI_MODE::BINARY);
   generator->SetInputImage(outImage);
   generator->SetUseImageGeometry(true);
   generator->Update();
   auto tdi = generator->GetOutput();
 
   if (m_NumParcels < 3)
     m_NumParcels = mitk::Tractometry::EstimateNumSamplingPoints(tdi, m_InputTract, 3);
 
   itk::TractsToVectorImageFilter<float>::Pointer fOdfFilter = itk::TractsToVectorImageFilter<float>::New();
   fOdfFilter->SetMaskImage(tdi);
   fOdfFilter->SetFiberBundle(m_InputTract);
   fOdfFilter->SetNormalizationMethod(itk::TractsToVectorImageFilter<float>::NormalizationMethods::SINGLE_VEC_NORM);
   fOdfFilter->SetMaxNumDirections(1);
+  fOdfFilter->SetOnlyUseMaskGeometry(true);
   fOdfFilter->Update();
   itk::Image< float, 4 >::Pointer dir_image = fOdfFilter->GetDirectionImage();
 
   m_WorkingTract = GetWorkingFib();
   vtkSmartPointer< vtkPolyData > polydata = m_WorkingTract->GetFiberPolyData();
 
   float maxd = m_InputTract->GetMeanFiberLength()/(0.5*m_NumParcels);
   itk::ImageRegionIterator< OutImageType > it(outImage, outImage->GetLargestPossibleRegion());
   itk::ImageRegionIterator< OutImageType > it_tdi(tdi, tdi->GetLargestPossibleRegion());
 
   vtkSmartPointer< vtkPolyData > reference_polydata = nullptr;
   if (m_ReferenceTract.IsNotNull())
     reference_polydata = m_ReferenceTract->GetFiberPolyData();
 
   unsigned long num_vox = 0;
   while( !it_tdi.IsAtEnd() )
   {
     if (it_tdi.Get()>0)
       ++num_vox;
     ++it_tdi;
   }
   it_tdi.GoToBegin();
 
   MITK_INFO << "Parcellating tract";
   boost::progress_display disp(num_vox);
   while( !it.IsAtEnd() )
   {
     if (it_tdi.Get()>0)
     {
       int final_seg_id = -1;
       int mult = 1;
 
       itk::Image< float, 4 >::IndexType idx4;
       idx4[0] = it_tdi.GetIndex()[0];
       idx4[1] = it_tdi.GetIndex()[1];
       idx4[2] = it_tdi.GetIndex()[2];
 
       vnl_vector_fixed<float, 3> ref_dir;
       idx4[3] = 0;
       ref_dir[0] = dir_image->GetPixel(idx4);
       idx4[3] = 1;
       ref_dir[1] = dir_image->GetPixel(idx4);
       idx4[3] = 2;
       ref_dir[2] = dir_image->GetPixel(idx4);
 
-      while(final_seg_id<0)
+      if (ref_dir.magnitude()>0.01)
       {
-        std::vector<float> seg_vote; seg_vote.resize(m_NumParcels, 0);
-        typename OutImageType::PointType image_point;
-        tdi->TransformIndexToPhysicalPoint(it.GetIndex(), image_point);
+        ref_dir.normalize();
 
-        for (unsigned int i=0; i<m_WorkingTract->GetNumFibers(); ++i)
+        while(final_seg_id<0 && mult<5)
         {
-          vtkCell* cell = polydata->GetCell(i);
-          auto numPoints = cell->GetNumberOfPoints();
-          vtkPoints* points = cell->GetPoints();
+          std::vector<float> seg_vote; seg_vote.resize(m_NumParcels, 0);
+          typename OutImageType::PointType image_point;
+          tdi->TransformIndexToPhysicalPoint(it.GetIndex(), image_point);
 
-          bool flip = Flip(polydata, i, reference_polydata);
-
-          float local_d = 99999999;
-          int local_closest_seg = -1;
-//          float weight = 1.0;
-
-          for (int j=0; j<numPoints; j++)
+          for (unsigned int i=0; i<m_WorkingTract->GetNumFibers(); ++i)
           {
-            itk::Point<float, 3> p;
-            int segment_id = -1;
-            if (flip)
-            {
-              segment_id = numPoints - j - 1;
-              p = mitk::imv::GetItkPoint(points->GetPoint(segment_id));
-            }
-            else
-            {
-              p = mitk::imv::GetItkPoint(points->GetPoint(j));
-              segment_id = j;
-            }
+            vtkCell* cell = polydata->GetCell(i);
+            auto numPoints = cell->GetNumberOfPoints();
+            vtkPoints* points = cell->GetPoints();
 
-            float d = std::fabs( (p[0]-image_point[0]) ) + std::fabs( (p[1]-image_point[1]) ) + std::fabs( (p[2]-image_point[2]) );
+            bool flip = Flip(polydata, i, reference_polydata);
 
+            float local_d = 99999999;
+            int local_closest_seg = -1;
+  //          float weight = 1.0;
 
-            itk::Point<float, 3> p2;
-            if (segment_id<numPoints-1)
-            {
-              p2 = mitk::imv::GetItkPoint(points->GetPoint(segment_id+1));
-            }
-            else
+            for (int j=0; j<numPoints; j++)
             {
-              p2 = mitk::imv::GetItkPoint(points->GetPoint(segment_id-1));
+              itk::Point<float, 3> p;
+              int segment_id = -1;
+              if (flip)
+              {
+                segment_id = numPoints - j - 1;
+                p = mitk::imv::GetItkPoint(points->GetPoint(segment_id));
+              }
+              else
+              {
+                p = mitk::imv::GetItkPoint(points->GetPoint(j));
+                segment_id = j;
+              }
+
+              float d = std::fabs( (p[0]-image_point[0]) ) + std::fabs( (p[1]-image_point[1]) ) + std::fabs( (p[2]-image_point[2]) );
+
+
+              itk::Point<float, 3> p2;
+              if (segment_id<numPoints-1)
+              {
+                p2 = mitk::imv::GetItkPoint(points->GetPoint(segment_id+1));
+              }
+              else
+              {
+                p2 = mitk::imv::GetItkPoint(points->GetPoint(segment_id-1));
+              }
+
+              vnl_vector_fixed<float, 3> dir;
+              dir[0] = p[0]-p2[0];
+              dir[1] = p[1]-p2[1];
+              dir[2] = p[2]-p2[2];
+
+              if (dir.magnitude()<0.0000001)
+                continue;
+              dir.normalize();
+
+              float a = std::fabs(dot_product(dir, ref_dir));
+              if (a<0.0000001)
+                a += 0.0000001;
+              d += (1.0/a - 1.0) * maxd;
+
+              if (d<local_d)
+              {
+                local_d = d;
+                local_closest_seg = j;
+  //              typename OutImageType::IndexType p_idx;
+  //              typename OutImageType::PointType mitk_p;
+  //              mitk_p[0] = p[0];
+  //              mitk_p[1] = p[2];
+  //              mitk_p[2] = p[1];
+  //              tdi->TransformPhysicalPointToIndex(mitk_p, p_idx);
+  //              weight = tdi->GetPixel(p_idx);
+              }
             }
 
-            vnl_vector_fixed<float, 3> dir;
-            dir[0] = p[0]-p2[0];
-            dir[1] = p[1]-p2[1];
-            dir[2] = p[2]-p2[2];
-            dir.normalize();
-
-            float a = std::fabs(dot_product(dir, ref_dir));
-            if (a<0.0000001)
-              a += 0.0000001;
-            d += (1.0/a - 1.0) * maxd;
+            if (local_d<maxd*mult)
+              seg_vote[local_closest_seg] += 1.0/(local_d);
+          }
 
-            if (d<local_d)
+          float max_count = 0;
+          for (unsigned int i=0; i<seg_vote.size(); ++i)
+          {
+            if (seg_vote.at(i)>max_count)
             {
-              local_d = d;
-              local_closest_seg = j;
-//              typename OutImageType::IndexType p_idx;
-//              typename OutImageType::PointType mitk_p;
-//              mitk_p[0] = p[0];
-//              mitk_p[1] = p[2];
-//              mitk_p[2] = p[1];
-//              tdi->TransformPhysicalPointToIndex(mitk_p, p_idx);
-//              weight = tdi->GetPixel(p_idx);
+              final_seg_id = i;
+              max_count = seg_vote.at(i);
             }
           }
 
-          if (local_d<maxd*mult)
-            seg_vote[local_closest_seg] += 1.0/(local_d);
-        }
+          if (final_seg_id>=0)
+            it.Set(final_seg_id + 1);
 
-        float max_count = 0;
-        for (unsigned int i=0; i<seg_vote.size(); ++i)
-        {
-          if (seg_vote.at(i)>max_count)
-          {
-            final_seg_id = i;
-            max_count = seg_vote.at(i);
-          }
+          ++mult;
         }
-
-        if (final_seg_id>=0)
-          it.Set(final_seg_id + 1);
-
-        ++mult;
       }
 
       ++disp;
     }
     ++it;
     ++it_tdi;
   }
   MITK_INFO << "DONE";
 }
 
 template< class OutImageType, class InputImageType >
 void TractParcellationFilter< OutImageType, InputImageType >::GenerateData()
 {
   // generate upsampled image
   mitk::BaseGeometry::Pointer geometry = m_InputTract->GetGeometry();
 
   // calculate new image parameters
   itk::Vector<double, 3> newSpacing;
   mitk::Point3D newOrigin;
   itk::Matrix<double, 3, 3> newDirection;
   ImageRegion<3> upsampledRegion;
   if (!m_InputImage.IsNull())
   {
     newSpacing = m_InputImage->GetSpacing()/m_UpsamplingFactor;
     upsampledRegion = m_InputImage->GetLargestPossibleRegion();
     newOrigin = m_InputImage->GetOrigin();
     typename OutImageType::RegionType::SizeType size = upsampledRegion.GetSize();
     size[0] *= m_UpsamplingFactor;
     size[1] *= m_UpsamplingFactor;
     size[2] *= m_UpsamplingFactor;
     upsampledRegion.SetSize(size);
     newDirection = m_InputImage->GetDirection();
   }
   else
   {
     newSpacing = geometry->GetSpacing()/m_UpsamplingFactor;
     newOrigin = geometry->GetOrigin();
     mitk::Geometry3D::BoundsArrayType bounds = geometry->GetBounds();
     newOrigin[0] += bounds.GetElement(0);
     newOrigin[1] += bounds.GetElement(2);
     newOrigin[2] += bounds.GetElement(4);
 
     for (int i=0; i<3; i++)
       for (int j=0; j<3; j++)
         newDirection[j][i] = geometry->GetMatrixColumn(i)[j];
     upsampledRegion.SetSize(0, geometry->GetExtent(0)*m_UpsamplingFactor);
     upsampledRegion.SetSize(1, geometry->GetExtent(1)*m_UpsamplingFactor);
     upsampledRegion.SetSize(2, geometry->GetExtent(2)*m_UpsamplingFactor);
   }
 
   // apply new image parameters
   typename OutImageType::Pointer outImage = OutImageType::New();
   outImage->SetSpacing( newSpacing );
   outImage->SetOrigin( newOrigin );
   outImage->SetDirection( newDirection );
   outImage->SetRegions( upsampledRegion );
   outImage->Allocate();
   outImage->FillBuffer(0);
   this->SetNthOutput(0, outImage);
 
   StaticResampleParcelVoting(outImage);
 
   auto outImage_pp = PostprocessParcellation(outImage);
   outImage_pp = PostprocessParcellation(outImage_pp);
   outImage_pp = PostprocessParcellation(outImage_pp);
   this->SetNthOutput(1, outImage_pp);
 }
 }
diff --git a/Modules/FiberTracking/Algorithms/itkTractsToVectorImageFilter.cpp b/Modules/FiberTracking/Algorithms/itkTractsToVectorImageFilter.cpp
index 091edea..d104be4 100644
--- a/Modules/FiberTracking/Algorithms/itkTractsToVectorImageFilter.cpp
+++ b/Modules/FiberTracking/Algorithms/itkTractsToVectorImageFilter.cpp
@@ -1,396 +1,399 @@
 #include "itkTractsToVectorImageFilter.h"
 
 // VTK
 #include <vtkPolyLine.h>
 #include <vtkCellArray.h>
 #include <vtkCellData.h>
 
 // ITK
 #include <itkTimeProbe.h>
 #include <itkImageRegionIterator.h>
 
 // misc
 #define _USE_MATH_DEFINES
 #include <cmath>
 #include <boost/progress.hpp>
 #include <mitkDiffusionFunctionCollection.h>
 
 namespace itk{
 
 static inline bool CompareVectorLengths(const vnl_vector_fixed< double, 3 >& v1, const vnl_vector_fixed< double, 3 >& v2)
 {
   return (v1.magnitude()>v2.magnitude());
 }
 
 template< class PixelType >
 TractsToVectorImageFilter< PixelType >::TractsToVectorImageFilter():
   m_NormalizationMethod(GLOBAL_MAX),
   m_AngularThreshold(0.7),
   m_Epsilon(0.999),
   m_MaxNumDirections(3),
   m_SizeThreshold(0.3),
   m_OnlyUseMaskGeometry(false)
 {
   this->SetNumberOfRequiredOutputs(1);
 }
 
 
 template< class PixelType >
 TractsToVectorImageFilter< PixelType >::~TractsToVectorImageFilter()
 {
 }
 
 
 template< class PixelType >
 vnl_vector_fixed<double, 3> TractsToVectorImageFilter< PixelType >::GetVnlVector(double point[])
 {
   vnl_vector_fixed<double, 3> vnlVector;
   vnlVector[0] = point[0];
   vnlVector[1] = point[1];
   vnlVector[2] = point[2];
   return vnlVector;
 }
 
 template< class PixelType >
 void TractsToVectorImageFilter< PixelType >::GenerateData()
 {
   mitk::BaseGeometry::Pointer geometry = m_FiberBundle->GetGeometry();
 
   // calculate new image parameters
   itk::Vector<double> spacing3;
   itk::Point<double> origin3;
   itk::Matrix<double, 3, 3> direction3;
   ImageRegion<3> imageRegion3;
   if (!m_MaskImage.IsNull())
   {
     spacing3 = m_MaskImage->GetSpacing();
     imageRegion3 = m_MaskImage->GetLargestPossibleRegion();
     origin3 = m_MaskImage->GetOrigin();
     direction3 = m_MaskImage->GetDirection();
   }
   else
   {
     spacing3 = geometry->GetSpacing();
     origin3 = geometry->GetOrigin();
     mitk::BaseGeometry::BoundsArrayType bounds = geometry->GetBounds();
     origin3[0] += bounds.GetElement(0);
     origin3[1] += bounds.GetElement(2);
     origin3[2] += bounds.GetElement(4);
 
     for (int i=0; i<3; i++)
       for (int j=0; j<3; j++)
         direction3[j][i] = geometry->GetMatrixColumn(i)[j];
     imageRegion3.SetSize(0, geometry->GetExtent(0)+1);
     imageRegion3.SetSize(1, geometry->GetExtent(1)+1);
     imageRegion3.SetSize(2, geometry->GetExtent(2)+1);
 
 
     m_MaskImage = ItkUcharImgType::New();
     m_MaskImage->SetSpacing( spacing3 );
     m_MaskImage->SetOrigin( origin3 );
     m_MaskImage->SetDirection( direction3 );
     m_MaskImage->SetRegions( imageRegion3 );
     m_MaskImage->Allocate();
     m_MaskImage->FillBuffer(1);
   }
   OutputImageType::RegionType::SizeType outImageSize = imageRegion3.GetSize();
   m_OutImageSpacing = m_MaskImage->GetSpacing();
   m_ClusteredDirectionsContainer = ContainerType::New();
 
   // initialize num directions image
   m_NumDirectionsImage = ItkUcharImgType::New();
   m_NumDirectionsImage->SetSpacing( spacing3 );
   m_NumDirectionsImage->SetOrigin( origin3 );
   m_NumDirectionsImage->SetDirection( direction3 );
   m_NumDirectionsImage->SetRegions( imageRegion3 );
   m_NumDirectionsImage->Allocate();
   m_NumDirectionsImage->FillBuffer(0);
 
   itk::Vector<double, 4> spacing4;
   itk::Point<float, 4> origin4;
   itk::Matrix<double, 4, 4> direction4;
   itk::ImageRegion<4> imageRegion4;
 
   spacing4[0] = spacing3[0]; spacing4[1] = spacing3[1]; spacing4[2] = spacing3[2]; spacing4[3] = 1;
   origin4[0] = origin3[0]; origin4[1] = origin3[1]; origin4[2] = origin3[2]; origin3[3] = 0;
   for (int r=0; r<3; r++)
     for (int c=0; c<3; c++)
       direction4[r][c] = direction3[r][c];
   direction4[3][3] = 1;
   imageRegion4.SetSize(0, imageRegion3.GetSize()[0]);
   imageRegion4.SetSize(1, imageRegion3.GetSize()[1]);
   imageRegion4.SetSize(2, imageRegion3.GetSize()[2]);
   imageRegion4.SetSize(3, m_MaxNumDirections*3);
 
   m_DirectionImage = ItkDirectionImageType::New();
   m_DirectionImage->SetSpacing( spacing4 );
   m_DirectionImage->SetOrigin( origin4 );
   m_DirectionImage->SetDirection( direction4 );
   m_DirectionImage->SetRegions( imageRegion4 );
   m_DirectionImage->Allocate();
   m_DirectionImage->FillBuffer(0.0);
 
   // iterate over all fibers
   vtkSmartPointer<vtkPolyData> fiberPolyData = m_FiberBundle->GetFiberPolyData();
   int numFibers = m_FiberBundle->GetNumFibers();
   m_DirectionsContainer = ContainerType::New();
 
   VectorContainer< unsigned int, std::vector< double > >::Pointer peakLengths = VectorContainer< unsigned int, std::vector< double > >::New();
 
   MITK_INFO << "Generating directions from tractogram";
   boost::progress_display disp(numFibers);
   for( int i=0; i<numFibers; i++ )
   {
     ++disp;
     vtkCell* cell = fiberPolyData->GetCell(i);
     int numPoints = cell->GetNumberOfPoints();
     vtkPoints* points = cell->GetPoints();
     if (numPoints<2)
       continue;
 
     vnl_vector_fixed<double, 3> dir;
     vnl_vector<double> v;
 
     float fiberWeight = m_FiberBundle->GetFiberWeight(i);
 
     for( int j=0; j<numPoints-1; j++)
     {
       itk::Point<float, 3> startVertex = mitk::imv::GetItkPoint(points->GetPoint(j));
       itk::Index<3> startIndex;
       itk::ContinuousIndex<float, 3> startIndexCont;
       m_MaskImage->TransformPhysicalPointToIndex(startVertex, startIndex);
       m_MaskImage->TransformPhysicalPointToContinuousIndex(startVertex, startIndexCont);
 
       itk::Point<float, 3> endVertex = mitk::imv::GetItkPoint(points->GetPoint(j + 1));
       itk::Index<3> endIndex;
       itk::ContinuousIndex<float, 3> endIndexCont;
       m_MaskImage->TransformPhysicalPointToIndex(endVertex, endIndex);
       m_MaskImage->TransformPhysicalPointToContinuousIndex(endVertex, endIndexCont);
 
       dir[0] = endVertex[0]-startVertex[0];
       dir[1] = endVertex[1]-startVertex[1];
       dir[2] = endVertex[2]-startVertex[2];
       if (dir.is_zero())
         continue;
       dir.normalize();
 
       std::vector< std::pair< itk::Index<3>, double > > segments = mitk::imv::IntersectImage(spacing3, startIndex, endIndex, startIndexCont, endIndexCont);
       for (std::pair< itk::Index<3>, double > segment : segments)
       {
         if (!m_MaskImage->GetLargestPossibleRegion().IsInside(segment.first) || (!m_OnlyUseMaskGeometry && m_MaskImage->GetPixel(segment.first)==0))
           continue;
 
         // add direction to container
         unsigned int idx = segment.first[0] + outImageSize[0]*(segment.first[1] + outImageSize[1]*segment.first[2]);
         DirectionContainerType::Pointer dirCont;
         if (m_DirectionsContainer->IndexExists(idx))
         {
           peakLengths->ElementAt(idx).push_back(fiberWeight*segment.second);
 
           dirCont = m_DirectionsContainer->GetElement(idx);
           if (dirCont.IsNull())
           {
             dirCont = DirectionContainerType::New();
             dirCont->push_back(dir);
             m_DirectionsContainer->InsertElement(idx, dirCont);
           }
           else
             dirCont->push_back(dir);
         }
         else
         {
           dirCont = DirectionContainerType::New();
           dirCont->push_back(dir);
           m_DirectionsContainer->InsertElement(idx, dirCont);
 
           std::vector< double > lengths; lengths.push_back(fiberWeight*segment.second);
           peakLengths->InsertElement(idx, lengths);
         }
       }
 
     }
   }
 
   itk::ImageRegionIterator<ItkUcharImgType> dirIt(m_NumDirectionsImage, m_NumDirectionsImage->GetLargestPossibleRegion());
 
   MITK_INFO << "Clustering directions";
   float max_dir_mag = 0;
   boost::progress_display disp2(outImageSize[0]*outImageSize[1]*outImageSize[2]);
   while(!dirIt.IsAtEnd())
   {
     ++disp2;
     OutputImageType::IndexType idx3 = dirIt.GetIndex();
     int idx_lin = idx3[0]+(idx3[1]+idx3[2]*outImageSize[1])*outImageSize[0];
 
     itk::Index<4> idx4; idx4[0] = idx3[0]; idx4[1] = idx3[1]; idx4[2] = idx3[2];
 
     if (!m_DirectionsContainer->IndexExists(idx_lin))
     {
       ++dirIt;
       continue;
     }
     DirectionContainerType::Pointer dirCont = m_DirectionsContainer->GetElement(idx_lin);
     if (dirCont.IsNull() || dirCont->empty())
     {
       ++dirIt;
       continue;
     }
 
     DirectionContainerType::Pointer directions;
     if (m_MaxNumDirections>0)
     {
       directions = FastClustering(dirCont, peakLengths->GetElement(idx_lin));
       std::sort( directions->begin(), directions->end(), CompareVectorLengths );
     }
     else
       directions = dirCont;
 
     unsigned int numDir = directions->size();
     if (m_MaxNumDirections>0 && numDir>m_MaxNumDirections)
       numDir = m_MaxNumDirections;
 
     float voxel_max_mag = 0;
     for (unsigned int i=0; i<numDir; i++)
     {
       DirectionType dir = directions->at(i);
       float mag = dir.magnitude();
 
       if (mag>voxel_max_mag)
         voxel_max_mag = mag;
       if (mag>max_dir_mag)
         max_dir_mag = mag;
     }
 
     int count = 0;
     for (unsigned int i=0; i<numDir; i++)
     {
       DirectionType dir = directions->at(i);
       count++;
 
       float mag = dir.magnitude();
       if (m_NormalizationMethod==MAX_VEC_NORM && voxel_max_mag>mitk::eps)
         dir /= voxel_max_mag;
       else if (m_NormalizationMethod==SINGLE_VEC_NORM && mag>mitk::eps)
         dir.normalize();
 
       for (unsigned int j = 0; j<3; j++)
       {
         idx4[3] = i*3 + j;
         m_DirectionImage->SetPixel(idx4, dir[j]);
       }
     }
     dirIt.Set(count);
     ++dirIt;
   }
 
   if (m_NormalizationMethod==GLOBAL_MAX && max_dir_mag>0)
   {
     itk::ImageRegionIterator<ItkDirectionImageType> dirImgIt(m_DirectionImage, m_DirectionImage->GetLargestPossibleRegion());
     while(!dirImgIt.IsAtEnd())
     {
       dirImgIt.Set(dirImgIt.Get()/max_dir_mag);
       ++dirImgIt;
     }
   }
+
+  if (max_dir_mag < 0.00000001)
+    mitkThrow() << "Maximum peak size is 0.0!";
 }
 
 
 template< class PixelType >
 TractsToVectorImageFilter< PixelType >::DirectionContainerType::Pointer TractsToVectorImageFilter< PixelType >::FastClustering(DirectionContainerType::Pointer inDirs, std::vector< double > lengths)
 {
   DirectionContainerType::Pointer outDirs = DirectionContainerType::New();
   if (inDirs->size()<2)
   {
     if (inDirs->size()==1)
       inDirs->SetElement(0, inDirs->at(0)*lengths.at(0));
     return inDirs;
   }
 
   DirectionType oldMean, currentMean;
   std::vector< int > touched;
 
   // initialize
   touched.resize(inDirs->size(), 0);
   bool free = true;
   currentMean = inDirs->at(0);  // initialize first seed
   currentMean.normalize();
   double length = lengths.at(0);
   touched[0] = 1;
   std::vector< double > newLengths;
   bool meanChanged = false;
 
   double max = 0;
   while (free)
   {
     oldMean.fill(0.0);
 
     // start mean-shift clustering
     double angle = 0;
 
     while (fabs(dot_product(currentMean, oldMean))<0.99)
     {
       oldMean = currentMean;
       currentMean.fill(0.0);
       for (unsigned int i=0; i<inDirs->size(); i++)
       {
         angle = dot_product(oldMean, inDirs->at(i));
         if (angle>=m_AngularThreshold)
         {
           currentMean += inDirs->at(i);
           if (meanChanged)
             length += lengths.at(i);
           touched[i] = 1;
           meanChanged = true;
         }
         else if (-angle>=m_AngularThreshold)
         {
           currentMean -= inDirs->at(i);
           if (meanChanged)
             length += lengths.at(i);
           touched[i] = 1;
           meanChanged = true;
         }
       }
       if(!meanChanged)
         currentMean = oldMean;
       else
         currentMean.normalize();
     }
 
     // found stable mean
     outDirs->push_back(currentMean);
     newLengths.push_back(length);
     if (length>max)
       max = length;
 
     // find next unused seed
     free = false;
     for (unsigned int i=0; i<touched.size(); i++)
       if (touched[i]==0)
       {
         currentMean = inDirs->at(i);
         free = true;
         meanChanged = false;
         length = lengths.at(i);
         touched[i] = 1;
         break;
       }
   }
 
   if (inDirs->size()==outDirs->size())
   {
     if (max>0)
     {
       for (unsigned int i=0; i<outDirs->size(); i++)
         outDirs->SetElement(i, outDirs->at(i)*newLengths.at(i));
     }
     return outDirs;
   }
   else
     return FastClustering(outDirs, newLengths);
 }
 
 }
 
 
 
diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkFiberQuantificationView.cpp b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkFiberQuantificationView.cpp
index 32a7f27..4738129 100644
--- a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkFiberQuantificationView.cpp
+++ b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkFiberQuantificationView.cpp
@@ -1,543 +1,544 @@
 /*===================================================================
 
 The Medical Imaging Interaction Toolkit (MITK)
 
 Copyright (c) German Cancer Research Center.
 
 All rights reserved.
 
 This software is distributed WITHOUT ANY WARRANTY; without
 even the implied warranty of MERCHANTABILITY or FITNESS FOR
 A PARTICULAR PURPOSE.
 
 See LICENSE.txt or http://www.mitk.org for details.
 
 ===================================================================*/
 
 
 // Blueberry
 #include <berryISelectionService.h>
 #include <berryIWorkbenchWindow.h>
 
 // Qmitk
 #include "QmitkFiberQuantificationView.h"
 
 // Qt
 #include <QMessageBox>
 
 // MITK
 #include <mitkNodePredicateDataType.h>
 #include <mitkNodePredicateDimension.h>
 #include <mitkNodePredicateAnd.h>
 #include <mitkImageCast.h>
 #include <mitkPeakImage.h>
 #include <mitkLabelSetImage.h>
 #include <mitkDICOMSegmentationConstants.h>
 #include <mitkDICOMSegmentationPropertyHelper.cpp>
 #include <mitkDICOMQIPropertyHelper.h>
 
 // ITK
 #include <itkTractDensityImageFilter.h>
 #include <itkTractsToRgbaImageFilter.h>
 #include <itkTractsToVectorImageFilter.h>
 #include <itkTractsToFiberEndingsImageFilter.h>
 #include <itkSignedMaurerDistanceMapImageFilter.h>
 #include <itkBinaryThinningImageFilter.h>
 
 #include <mitkLexicalCast.h>
 
 const std::string QmitkFiberQuantificationView::VIEW_ID = "org.mitk.views.fiberquantification";
 using namespace mitk;
 
 QmitkFiberQuantificationView::QmitkFiberQuantificationView()
   : QmitkAbstractView()
   , m_Controls( 0 )
   , m_UpsamplingFactor(5)
   , m_Visible(false)
 {
 
 }
 
 // Destructor
 QmitkFiberQuantificationView::~QmitkFiberQuantificationView()
 {
 
 }
 
 void QmitkFiberQuantificationView::CreateQtPartControl( QWidget *parent )
 {
   // build up qt view, unless already done
   if ( !m_Controls )
   {
     // create GUI widgets from the Qt Designer's .ui file
     m_Controls = new Ui::QmitkFiberQuantificationViewControls;
     m_Controls->setupUi( parent );
 
     connect( m_Controls->m_ProcessFiberBundleButton, SIGNAL(clicked()), this, SLOT(ProcessSelectedBundles()) );
     connect( m_Controls->m_ExtractFiberPeaks, SIGNAL(clicked()), this, SLOT(CalculateFiberDirections()) );
 
     m_Controls->m_TractBox->SetDataStorage(this->GetDataStorage());
     mitk::TNodePredicateDataType<mitk::FiberBundle>::Pointer isFib = mitk::TNodePredicateDataType<mitk::FiberBundle>::New();
     m_Controls->m_TractBox->SetPredicate( isFib );
 
     m_Controls->m_ImageBox->SetDataStorage(this->GetDataStorage());
     m_Controls->m_ImageBox->SetZeroEntryText("--");
     mitk::TNodePredicateDataType<mitk::Image>::Pointer isImagePredicate = mitk::TNodePredicateDataType<mitk::Image>::New();
     mitk::NodePredicateDimension::Pointer is3D = mitk::NodePredicateDimension::New(3);
     m_Controls->m_ImageBox->SetPredicate( mitk::NodePredicateAnd::New(isImagePredicate, is3D) );
 
     connect( (QObject*)(m_Controls->m_TractBox), SIGNAL(currentIndexChanged(int)), this, SLOT(UpdateGui()));
     connect( (QObject*)(m_Controls->m_ImageBox), SIGNAL(currentIndexChanged(int)), this, SLOT(UpdateGui()));
   }
 }
 
 void QmitkFiberQuantificationView::Activated()
 {
 
 }
 
 void QmitkFiberQuantificationView::Deactivated()
 {
 
 }
 
 void QmitkFiberQuantificationView::Visible()
 {
   m_Visible = true;
 }
 
 void QmitkFiberQuantificationView::Hidden()
 {
   m_Visible = false;
 }
 
 void QmitkFiberQuantificationView::SetFocus()
 {
   m_Controls->m_ProcessFiberBundleButton->setFocus();
 }
 
 void QmitkFiberQuantificationView::CalculateFiberDirections()
 {
   typedef itk::Image<unsigned char, 3>                                            ItkUcharImgType;
 
   // load fiber bundle
   mitk::FiberBundle::Pointer inputTractogram = dynamic_cast<mitk::FiberBundle*>(m_SelectedFB.back()->GetData());
 
   itk::TractsToVectorImageFilter<float>::Pointer fOdfFilter = itk::TractsToVectorImageFilter<float>::New();
   if (m_SelectedImage.IsNotNull())
   {
     ItkUcharImgType::Pointer itkMaskImage = ItkUcharImgType::New();
     mitk::CastToItkImage<ItkUcharImgType>(m_SelectedImage, itkMaskImage);
     fOdfFilter->SetMaskImage(itkMaskImage);
   }
 
   // extract directions from fiber bundle
   fOdfFilter->SetFiberBundle(inputTractogram);
   fOdfFilter->SetAngularThreshold(cos(m_Controls->m_AngularThreshold->value()*itk::Math::pi/180));
   switch (m_Controls->m_FiberDirNormBox->currentIndex())
   {
   case 0:
     fOdfFilter->SetNormalizationMethod(itk::TractsToVectorImageFilter<float>::NormalizationMethods::GLOBAL_MAX);
     break;
   case 1:
     fOdfFilter->SetNormalizationMethod(itk::TractsToVectorImageFilter<float>::NormalizationMethods::SINGLE_VEC_NORM);
     break;
   case 2:
     fOdfFilter->SetNormalizationMethod(itk::TractsToVectorImageFilter<float>::NormalizationMethods::MAX_VEC_NORM);
     break;
   }
+  fOdfFilter->SetOnlyUseMaskGeometry(true);
   fOdfFilter->SetSizeThreshold(m_Controls->m_PeakThreshold->value());
   fOdfFilter->SetMaxNumDirections(m_Controls->m_MaxNumDirections->value());
   fOdfFilter->Update();
 
   QString name = m_SelectedFB.back()->GetName().c_str();
 
   if (m_Controls->m_NumDirectionsBox->isChecked())
   {
     mitk::Image::Pointer mitkImage = mitk::Image::New();
     mitkImage->InitializeByItk( fOdfFilter->GetNumDirectionsImage().GetPointer() );
     mitkImage->SetVolume( fOdfFilter->GetNumDirectionsImage()->GetBufferPointer() );
 
     mitk::DataNode::Pointer node = mitk::DataNode::New();
     node->SetData(mitkImage);
     node->SetName((name+"_NUM_DIRECTIONS").toStdString().c_str());
     GetDataStorage()->Add(node, m_SelectedFB.back());
   }
 
   Image::Pointer mitkImage = dynamic_cast<Image*>(PeakImage::New().GetPointer());
   mitk::CastToMitkImage(fOdfFilter->GetDirectionImage(), mitkImage);
   mitkImage->SetVolume(fOdfFilter->GetDirectionImage()->GetBufferPointer());
 
   mitk::DataNode::Pointer node = mitk::DataNode::New();
   node->SetData(mitkImage);
   node->SetName( (name+"_DIRECTIONS").toStdString().c_str());
   GetDataStorage()->Add(node, m_SelectedFB.back());
 }
 
 void QmitkFiberQuantificationView::UpdateGui()
 {
   m_SelectedFB.clear();
   if (m_Controls->m_TractBox->GetSelectedNode().IsNotNull())
     m_SelectedFB.push_back(m_Controls->m_TractBox->GetSelectedNode());
 
   m_SelectedImage = nullptr;
   if (m_Controls->m_ImageBox->GetSelectedNode().IsNotNull())
     m_SelectedImage = dynamic_cast<mitk::Image*>(m_Controls->m_ImageBox->GetSelectedNode()->GetData());
 
   m_Controls->m_ProcessFiberBundleButton->setEnabled(!m_SelectedFB.empty());
   m_Controls->m_ExtractFiberPeaks->setEnabled(!m_SelectedFB.empty());
 }
 
 void QmitkFiberQuantificationView::OnSelectionChanged(berry::IWorkbenchPart::Pointer /*part*/, const QList<mitk::DataNode::Pointer>& )
 {
   UpdateGui();
 }
 
 void QmitkFiberQuantificationView::ProcessSelectedBundles()
 {
   if ( m_SelectedFB.empty() ){
     QMessageBox::information( nullptr, "Warning", "No fibe bundle selected!");
     MITK_WARN("QmitkFiberQuantificationView") << "no fibe bundle selected";
     return;
   }
 
   int generationMethod = m_Controls->m_GenerationBox->currentIndex();
 
   for( unsigned int i=0; i<m_SelectedFB.size(); i++ )
   {
     mitk::DataNode::Pointer node = m_SelectedFB[i];
     if (node.IsNotNull() && dynamic_cast<mitk::FiberBundle*>(node->GetData()))
     {
       mitk::FiberBundle::Pointer fib = dynamic_cast<mitk::FiberBundle*>(node->GetData());
       QString name(node->GetName().c_str());
       DataNode::Pointer newNode = nullptr;
       switch(generationMethod){
       case 0:
         newNode = GenerateTractDensityImage(fib, TDI_MODE::DENSITY, true, node->GetName());
         name += "_TDI";
         break;
       case 1:
         newNode = GenerateTractDensityImage(fib, TDI_MODE::DENSITY, false, node->GetName());
         name += "_TDI";
         break;
       case 2:
         newNode = GenerateTractDensityImage(fib, TDI_MODE::BINARY, false, node->GetName());
         name += "_envelope";
         break;
       case 3:
         newNode = GenerateColorHeatmap(fib);
         break;
       case 4:
         newNode = GenerateFiberEndingsImage(fib);
         name += "_fiber_endings";
         break;
       case 5:
         newNode = GenerateFiberEndingsPointSet(fib);
         name += "_fiber_endings";
         break;
       case 6:
         newNode = GenerateDistanceMap(fib);
         name += "_distance_map";
         break;
       case 7:
         newNode = GenerateBinarySkeleton(fib);
         name += "_skeleton";
         break;
       case 8:
         newNode = GenerateTractDensityImage(fib, TDI_MODE::VISITATION_COUNT, true, node->GetName());
         name += "_visitations";
         break;
       }
       if (newNode.IsNotNull())
       {
         newNode->SetName(name.toStdString());
         GetDataStorage()->Add(newNode);
       }
     }
   }
 }
 
 // generate pointset displaying the fiber endings
 mitk::DataNode::Pointer QmitkFiberQuantificationView::GenerateFiberEndingsPointSet(mitk::FiberBundle::Pointer fib)
 {
   mitk::PointSet::Pointer pointSet = mitk::PointSet::New();
   vtkSmartPointer<vtkPolyData> fiberPolyData = fib->GetFiberPolyData();
 
   int count = 0;
   int numFibers = fib->GetNumFibers();
   for( int i=0; i<numFibers; i++ )
   {
     vtkCell* cell = fiberPolyData->GetCell(i);
     int numPoints = cell->GetNumberOfPoints();
     vtkPoints* points = cell->GetPoints();
 
     if (numPoints>0)
     {
       double* point = points->GetPoint(0);
       itk::Point<double,3> itkPoint = mitk::imv::GetItkPoint(point);
       pointSet->InsertPoint(count, itkPoint);
       count++;
     }
     if (numPoints>2)
     {
       double* point = points->GetPoint(numPoints-1);
       itk::Point<double,3> itkPoint = mitk::imv::GetItkPoint(point);
       pointSet->InsertPoint(count, itkPoint);
       count++;
     }
   }
 
   mitk::DataNode::Pointer node = mitk::DataNode::New();
   node->SetData( pointSet );
   return node;
 }
 
 // generate image displaying the fiber endings
 mitk::DataNode::Pointer QmitkFiberQuantificationView::GenerateFiberEndingsImage(mitk::FiberBundle::Pointer fib)
 {
   typedef unsigned int OutPixType;
   typedef itk::Image<OutPixType,3> OutImageType;
 
   typedef itk::TractsToFiberEndingsImageFilter< OutImageType > ImageGeneratorType;
   ImageGeneratorType::Pointer generator = ImageGeneratorType::New();
   generator->SetFiberBundle(fib);
   generator->SetUpsamplingFactor(m_Controls->m_UpsamplingSpinBox->value());
   if (m_SelectedImage.IsNotNull())
   {
     OutImageType::Pointer itkImage = OutImageType::New();
     CastToItkImage(m_SelectedImage, itkImage);
     generator->SetInputImage(itkImage);
     generator->SetUseImageGeometry(true);
   }
   generator->Update();
 
   // get output image
   OutImageType::Pointer outImg = generator->GetOutput();
   mitk::Image::Pointer img = mitk::Image::New();
   img->InitializeByItk(outImg.GetPointer());
   img->SetVolume(outImg->GetBufferPointer());
 
   // init data node
   mitk::DataNode::Pointer node = mitk::DataNode::New();
   node->SetData(img);
   return node;
 }
 
 // generate rgba heatmap from fiber bundle
 mitk::DataNode::Pointer QmitkFiberQuantificationView::GenerateColorHeatmap(mitk::FiberBundle::Pointer fib)
 {
   typedef itk::RGBAPixel<unsigned char> OutPixType;
   typedef itk::Image<OutPixType, 3> OutImageType;
   typedef itk::TractsToRgbaImageFilter< OutImageType > ImageGeneratorType;
   ImageGeneratorType::Pointer generator = ImageGeneratorType::New();
   generator->SetFiberBundle(fib);
   generator->SetUpsamplingFactor(m_Controls->m_UpsamplingSpinBox->value());
   if (m_SelectedImage.IsNotNull())
   {
     itk::Image<unsigned char, 3>::Pointer itkImage = itk::Image<unsigned char, 3>::New();
     CastToItkImage(m_SelectedImage, itkImage);
     generator->SetInputImage(itkImage);
     generator->SetUseImageGeometry(true);
   }
   generator->Update();
 
   // get output image
   typedef itk::Image<OutPixType,3> OutType;
   OutType::Pointer outImg = generator->GetOutput();
   mitk::Image::Pointer img = mitk::Image::New();
   img->InitializeByItk(outImg.GetPointer());
   img->SetVolume(outImg->GetBufferPointer());
 
   // init data node
   mitk::DataNode::Pointer node = mitk::DataNode::New();
   node->SetData(img);
   return node;
 }
 
 mitk::DataNode::Pointer QmitkFiberQuantificationView::GenerateBinarySkeleton(mitk::FiberBundle::Pointer fib)
 {
   typedef itk::Image<unsigned char, 3> UcharImageType;
 
   itk::TractDensityImageFilter< UcharImageType >::Pointer envelope_generator = itk::TractDensityImageFilter< UcharImageType >::New();
   envelope_generator->SetFiberBundle(fib);
   envelope_generator->SetMode(TDI_MODE::BINARY);
   envelope_generator->SetUpsamplingFactor(m_Controls->m_UpsamplingSpinBox->value());
   if (m_SelectedImage.IsNotNull())
   {
     UcharImageType::Pointer itkImage = UcharImageType::New();
     CastToItkImage(m_SelectedImage, itkImage);
     envelope_generator->SetInputImage(itkImage);
     envelope_generator->SetUseImageGeometry(true);
 
   }
   envelope_generator->Update();
 
   itk::BinaryThinningImageFilter<UcharImageType, UcharImageType>::Pointer map_generator = itk::BinaryThinningImageFilter<UcharImageType, UcharImageType>::New();
   map_generator->SetInput(envelope_generator->GetOutput());
   map_generator->Update();
 
   UcharImageType::Pointer outImg = map_generator->GetOutput();
   mitk::Image::Pointer img = mitk::Image::New();
   img->InitializeByItk(outImg.GetPointer());
   img->SetVolume(outImg->GetBufferPointer());
 
   mitk::DataNode::Pointer node = mitk::DataNode::New();
   node->SetData(img);
   return node;
 }
 
 mitk::DataNode::Pointer QmitkFiberQuantificationView::GenerateDistanceMap(mitk::FiberBundle::Pointer fib)
 {
   typedef itk::Image<unsigned char, 3> UcharImageType;
   typedef itk::Image<float, 3> FloatImageType;
 
   itk::TractDensityImageFilter< UcharImageType >::Pointer envelope_generator = itk::TractDensityImageFilter< UcharImageType >::New();
   envelope_generator->SetFiberBundle(fib);
   envelope_generator->SetMode(TDI_MODE::BINARY);
   envelope_generator->SetUpsamplingFactor(m_Controls->m_UpsamplingSpinBox->value());
   if (m_SelectedImage.IsNotNull())
   {
     UcharImageType::Pointer itkImage = UcharImageType::New();
     CastToItkImage(m_SelectedImage, itkImage);
     envelope_generator->SetInputImage(itkImage);
     envelope_generator->SetUseImageGeometry(true);
 
   }
   envelope_generator->Update();
 
   itk::SignedMaurerDistanceMapImageFilter<UcharImageType, FloatImageType>::Pointer map_generator = itk::SignedMaurerDistanceMapImageFilter<UcharImageType, FloatImageType>::New();
   map_generator->SetInput(envelope_generator->GetOutput());
   map_generator->SetUseImageSpacing(true);
   map_generator->SetSquaredDistance(false);
   map_generator->SetInsideIsPositive(true);
   map_generator->Update();
 
   FloatImageType::Pointer outImg = map_generator->GetOutput();
   mitk::Image::Pointer img = mitk::Image::New();
   img->InitializeByItk(outImg.GetPointer());
   img->SetVolume(outImg->GetBufferPointer());
 
   mitk::DataNode::Pointer node = mitk::DataNode::New();
   node->SetData(img);
   return node;
 }
 
 // generate tract density image from fiber bundle
 mitk::DataNode::Pointer QmitkFiberQuantificationView::GenerateTractDensityImage(mitk::FiberBundle::Pointer fib, TDI_MODE mode, bool absolute, std::string name)
 {
   mitk::DataNode::Pointer node = mitk::DataNode::New();
   if (mode==TDI_MODE::BINARY)
   {
     typedef unsigned char OutPixType;
     typedef itk::Image<OutPixType, 3> OutImageType;
 
     itk::TractDensityImageFilter< OutImageType >::Pointer generator = itk::TractDensityImageFilter< OutImageType >::New();
     generator->SetFiberBundle(fib);
     generator->SetMode(mode);
     generator->SetOutputAbsoluteValues(absolute);
     generator->SetUpsamplingFactor(m_Controls->m_UpsamplingSpinBox->value());
     if (m_SelectedImage.IsNotNull())
     {
       OutImageType::Pointer itkImage = OutImageType::New();
       CastToItkImage(m_SelectedImage, itkImage);
       generator->SetInputImage(itkImage);
       generator->SetUseImageGeometry(true);
 
     }
     generator->Update();
 
     // get output image
     typedef itk::Image<OutPixType,3> OutType;
     OutType::Pointer outImg = generator->GetOutput();
     mitk::Image::Pointer img = mitk::Image::New();
     img->InitializeByItk(outImg.GetPointer());
     img->SetVolume(outImg->GetBufferPointer());
 
 
     if (m_SelectedImage.IsNotNull())
     {
       mitk::LabelSetImage::Pointer multilabelImage = mitk::LabelSetImage::New();
       multilabelImage->InitializeByLabeledImage(img);
       multilabelImage->GetActiveLabelSet()->SetActiveLabel(1);
       mitk::Label::Pointer label = multilabelImage->GetActiveLabel();
       label->SetName("Tractogram");
 
       // Add Segmented Property Category Code Sequence tags (0062, 0003): Sequence defining the general category of this
       // segment.
       // (0008,0100) Code Value
       label->SetProperty(
             DICOMTagPathToPropertyName(DICOMSegmentationConstants::SEGMENT_CATEGORY_CODE_VALUE_PATH()).c_str(),
             TemporoSpatialStringProperty::New("T-D000A"));
 
       // (0008,0102) Coding Scheme Designator
       label->SetProperty(
             DICOMTagPathToPropertyName(DICOMSegmentationConstants::SEGMENT_CATEGORY_CODE_SCHEME_PATH()).c_str(),
             TemporoSpatialStringProperty::New("SRT"));
 
       // (0008,0104) Code Meaning
       label->SetProperty(
             DICOMTagPathToPropertyName(DICOMSegmentationConstants::SEGMENT_CATEGORY_CODE_MEANING_PATH()).c_str(),
             TemporoSpatialStringProperty::New("Anatomical Structure"));
       //------------------------------------------------------------
       // Add Segmented Property Type Code Sequence (0062, 000F): Sequence defining the specific property type of this
       // segment.
       // (0008,0100) Code Value
       label->SetProperty(
             DICOMTagPathToPropertyName(DICOMSegmentationConstants::SEGMENT_TYPE_CODE_VALUE_PATH()).c_str(),
             TemporoSpatialStringProperty::New("DUMMY"));
 
       // (0008,0102) Coding Scheme Designator
       label->SetProperty(
             DICOMTagPathToPropertyName(DICOMSegmentationConstants::SEGMENT_TYPE_CODE_SCHEME_PATH()).c_str(),
             TemporoSpatialStringProperty::New("SRT"));
 
       // (0008,0104) Code Meaning
       label->SetProperty(
             DICOMTagPathToPropertyName(DICOMSegmentationConstants::SEGMENT_TYPE_CODE_MEANING_PATH()).c_str(),
             TemporoSpatialStringProperty::New(name));
 
       //Error: is undeclared//      mitk::DICOMQIPropertyHandler::DeriveDICOMSourceProperties(m_SelectedImage, multilabelImage);
 
       // init data node
       node->SetData(multilabelImage);
     }
     else
     {
       // init data node
       node->SetData(img);
     }
 
   }
   else
   {
     typedef float OutPixType;
     typedef itk::Image<OutPixType, 3> OutImageType;
 
     itk::TractDensityImageFilter< OutImageType >::Pointer generator = itk::TractDensityImageFilter< OutImageType >::New();
     generator->SetFiberBundle(fib);
     generator->SetMode(mode);
     generator->SetOutputAbsoluteValues(absolute);
     generator->SetUpsamplingFactor(m_Controls->m_UpsamplingSpinBox->value());
     if (m_SelectedImage.IsNotNull())
     {
       OutImageType::Pointer itkImage = OutImageType::New();
       CastToItkImage(m_SelectedImage, itkImage);
       generator->SetInputImage(itkImage);
       generator->SetUseImageGeometry(true);
 
     }
     //generator->SetDoFiberResampling(false);
     generator->Update();
 
     // get output image
     typedef itk::Image<OutPixType,3> OutType;
     OutType::Pointer outImg = generator->GetOutput();
     mitk::Image::Pointer img = mitk::Image::New();
     img->InitializeByItk(outImg.GetPointer());
     img->SetVolume(outImg->GetBufferPointer());
 
     // init data node
     node->SetData(img);
   }
   return node;
 }