diff --git a/Modules/DiffusionImaging/FiberTracking/Algorithms/mitkTractAnalyzer.cpp b/Modules/DiffusionImaging/FiberTracking/Algorithms/mitkTractAnalyzer.cpp
deleted file mode 100644
index 29d329292b..0000000000
--- a/Modules/DiffusionImaging/FiberTracking/Algorithms/mitkTractAnalyzer.cpp
+++ /dev/null
@@ -1,177 +0,0 @@
-/*===================================================================
-
-The Medical Imaging Interaction Toolkit (MITK)
-
-Copyright (c) German Cancer Research Center,
-Division of Medical and Biological Informatics.
-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.
-
-===================================================================*/
-
-#ifndef __mitkTractAnalyzer_cpp
-#define __mitkTractAnalyzer_cpp
-
-
-
-#include <mitkLogMacros.h>
-#include <mitkTractAnalyzer.h>
-
-#include <itkImageRegionIteratorWithIndex.h>
-
-
-#include <itkShortestPathCostFunctionTbss.h>
-#include <itkShortestPathImageFilter.h>
-
-#include <iostream>
-#include <fstream>
-
-
-using namespace std;
-
-namespace mitk {
-
-  TractAnalyzer::TractAnalyzer()
-  {
-
-  }
-
-
-  void TractAnalyzer::BuildGraph(itk::Index<3> startPoint, itk::Index<3> endPoint)
-  {
-
-    typedef itk::ShortestPathImageFilter<FloatImageType, CharImageType> ShortestPathFilterType;
-    typedef itk::ShortestPathCostFunctionTbss<FloatImageType> CostFunctionType;
-
-
-    FloatImageType::Pointer meanSkeleton;
-
-    mitk::CastToItkImage(m_InputImage, meanSkeleton);
-
-  // Only use the mitk image
-
-
-
-    if(meanSkeleton)
-    {
-      CostFunctionType::Pointer costFunction = CostFunctionType::New();
-      costFunction->SetImage(meanSkeleton);
-      costFunction->SetStartIndex(startPoint);
-      costFunction->SetEndIndex(endPoint);
-      costFunction->SetThreshold(m_Threshold);
-
-      ShortestPathFilterType::Pointer pathFinder = ShortestPathFilterType::New();
-      pathFinder->SetCostFunction(costFunction);
-      pathFinder->SetFullNeighborsMode(true);
-      //pathFinder->SetCalcMode(ShortestPathFilterType::A_STAR);
-      pathFinder->SetInput(meanSkeleton);
-      pathFinder->SetStartIndex(startPoint);
-      pathFinder->SetEndIndex(endPoint);
-      pathFinder->Update();
-
-      m_Path = pathFinder->GetVectorPath();
-
-
-      m_RoiImg = pathFinder->GetOutput();
-    }
-
-
-
-  }
-
-  void TractAnalyzer::MeasureRoi()
-  {
-
-    // Output two types
-    ProjectionsImageType::SizeType size = m_Projections->GetLargestPossibleRegion().GetSize();
-
-    std::ofstream file(m_FileName.c_str());
-
-    std::vector<std::string> individuals;
-    for(int i=0; i<m_Groups.size(); i++)
-    {
-      std::pair<std::string, int> group = m_Groups[i];
-      for(int j=0; j<group.second; j++)
-      {
-        individuals.push_back(group.first);
-      }
-    }
-
-
-    // Position ID
-
-    file << "position ";
-
-    for(int i=0; i<m_Roi.size(); i++)
-    {
-      file << "pos" << i << " ";
-    }
-
-    file << "\n";
-
-    // For every subect
-    for(int j=0; j<size[3]; j++)
-    {
-
-      file << individuals[j] << " ";
-
-      // For every point on the current ROI
-      for(int k=0; k<m_Roi.size(); k++)
-      {
-        itk::Index<3> ix = m_Roi[k];
-        itk::Index<4> ix4;
-        ix4[0] = ix[0];
-        ix4[1] = ix[1];
-        ix4[2] = ix[2];
-        ix4[3] = j;
-
-        float f = m_Projections->GetPixel(ix4);
-
-        file << f << " ";
-
-      }
-
-      file << "\n";
-
-    }
-
-    file.close();
-
-
-
-    // Write the long format output
-    std::ofstream fileLong(m_FileNameLong.c_str());
-
-    fileLong << "ID " << "group " << "position " << "value\n";
-
-    for(int i=0; i<m_Roi.size(); i++)
-    {
-      for (int j=0; j<individuals.size(); j++)
-      {
-        itk::Index<3> ix = m_Roi[i];
-        itk::Index<4> ix4;
-        ix4[0] = ix[0];
-        ix4[1] = ix[1];
-        ix4[2] = ix[2];
-        ix4[3] = j;
-
-        float f = m_Projections->GetPixel(ix4);
-
-        fileLong << "ID" << j << " " << individuals[j] << " " << "pos" << i << " "<< f << "\n";
-      }
-    }
-
-    fileLong.close();
-
-
-  }
-
-
-
-}
-#endif
diff --git a/Modules/DiffusionImaging/FiberTracking/files.cmake b/Modules/DiffusionImaging/FiberTracking/files.cmake
index c793736569..1792d5ee6d 100644
--- a/Modules/DiffusionImaging/FiberTracking/files.cmake
+++ b/Modules/DiffusionImaging/FiberTracking/files.cmake
@@ -1,103 +1,101 @@
 set(CPP_FILES
 
 # DataStructures -> FiberBundleX
   IODataStructures/FiberBundleX/mitkFiberBundleX.cpp
   IODataStructures/FiberBundleX/mitkFiberBundleXWriter.cpp
   IODataStructures/FiberBundleX/mitkFiberBundleXReader.cpp
   IODataStructures/FiberBundleX/mitkFiberBundleXIOFactory.cpp
   IODataStructures/FiberBundleX/mitkFiberBundleXWriterFactory.cpp
   IODataStructures/FiberBundleX/mitkFiberBundleXSerializer.cpp
   IODataStructures/FiberBundleX/mitkFiberBundleXThreadMonitor.cpp
 
   # DataStructures -> PlanarFigureComposite
   IODataStructures/PlanarFigureComposite/mitkPlanarFigureComposite.cpp
 
   # DataStructures
   IODataStructures/mitkFiberTrackingObjectFactory.cpp
 
   # Rendering
   Rendering/mitkFiberBundleXMapper2D.cpp
   Rendering/mitkFiberBundleXMapper3D.cpp
   Rendering/mitkFiberBundleXThreadMonitorMapper3D.cpp
   #Rendering/mitkPlanarFigureMapper3D.cpp
 
   # Interactions
   Interactions/mitkFiberBundleInteractor.cpp
 
-  # Algorithms
-  Algorithms/mitkTractAnalyzer.cpp
 
   # Tractography
   Algorithms/GibbsTracking/mitkParticleGrid.cpp
   Algorithms/GibbsTracking/mitkMetropolisHastingsSampler.cpp
   Algorithms/GibbsTracking/mitkEnergyComputer.cpp
   Algorithms/GibbsTracking/mitkGibbsEnergyComputer.cpp
   Algorithms/GibbsTracking/mitkFiberBuilder.cpp
   Algorithms/GibbsTracking/mitkSphereInterpolator.cpp
 )
 
 set(H_FILES
   # Rendering
   Rendering/mitkFiberBundleXMapper3D.h
   Rendering/mitkFiberBundleXMapper2D.h
   Rendering/mitkFiberBundleXThreadMonitorMapper3D.h
   #Rendering/mitkPlanarFigureMapper3D.h
 
   # DataStructures -> FiberBundleX
   IODataStructures/FiberBundleX/mitkFiberBundleX.h
   IODataStructures/FiberBundleX/mitkFiberBundleXWriter.h
   IODataStructures/FiberBundleX/mitkFiberBundleXReader.h
   IODataStructures/FiberBundleX/mitkFiberBundleXIOFactory.h
   IODataStructures/FiberBundleX/mitkFiberBundleXWriterFactory.h
   IODataStructures/FiberBundleX/mitkFiberBundleXSerializer.h
   IODataStructures/FiberBundleX/mitkFiberBundleXThreadMonitor.h
 
   IODataStructures/mitkFiberTrackingObjectFactory.h
 
   # Algorithms
   Algorithms/itkTractDensityImageFilter.h
   Algorithms/itkTractsToFiberEndingsImageFilter.h
   Algorithms/itkTractsToRgbaImageFilter.h
   Algorithms/itkElectrostaticRepulsionDiffusionGradientReductionFilter.h
   Algorithms/itkFibersFromPlanarFiguresFilter.h
   Algorithms/itkTractsToDWIImageFilter.h
   Algorithms/itkTractsToVectorImageFilter.h
   Algorithms/itkKspaceImageFilter.h
   Algorithms/itkDftImageFilter.h
   Algorithms/itkAddArtifactsToDwiImageFilter.h
   Algorithms/itkFieldmapGeneratorFilter.h
 
   # (old) Tractography
   Algorithms/itkGibbsTrackingFilter.h
   Algorithms/itkStochasticTractographyFilter.h
   Algorithms/itkStreamlineTrackingFilter.h
   Algorithms/GibbsTracking/mitkParticle.h
   Algorithms/GibbsTracking/mitkParticleGrid.h
   Algorithms/GibbsTracking/mitkMetropolisHastingsSampler.h
   Algorithms/GibbsTracking/mitkSimpSamp.h
   Algorithms/GibbsTracking/mitkEnergyComputer.h
   Algorithms/GibbsTracking/mitkGibbsEnergyComputer.h
   Algorithms/GibbsTracking/mitkSphereInterpolator.h
   Algorithms/GibbsTracking/mitkFiberBuilder.h
 
   # Signal Models
   SignalModels/mitkDiffusionSignalModel.h
   SignalModels/mitkTensorModel.h
   SignalModels/mitkBallModel.h
   SignalModels/mitkDotModel.h
   SignalModels/mitkAstroStickModel.h
   SignalModels/mitkStickModel.h
   SignalModels/mitkDiffusionNoiseModel.h
   SignalModels/mitkRicianNoiseModel.h
   SignalModels/mitkKspaceArtifact.h
   SignalModels/mitkGibbsRingingArtifact.h
 )
 
 set(RESOURCE_FILES
   # Binary directory resources
   FiberTrackingLUTBaryCoords.bin
   FiberTrackingLUTIndices.bin
 
   # Shaders
   Shaders/mitkShaderFiberClipping.xml
 )
diff --git a/Modules/DiffusionImaging/Quantification/Algorithms/mitkTractAnalyzer.cpp b/Modules/DiffusionImaging/Quantification/Algorithms/mitkTractAnalyzer.cpp
new file mode 100644
index 0000000000..bfa25499af
--- /dev/null
+++ b/Modules/DiffusionImaging/Quantification/Algorithms/mitkTractAnalyzer.cpp
@@ -0,0 +1,212 @@
+/*===================================================================
+
+The Medical Imaging Interaction Toolkit (MITK)
+
+Copyright (c) German Cancer Research Center,
+Division of Medical and Biological Informatics.
+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.
+
+===================================================================*/
+
+#ifndef __mitkTractAnalyzer_cpp
+#define __mitkTractAnalyzer_cpp
+
+
+#include <mitkProgressBar.h>
+#include <mitkLogMacros.h>
+#include <mitkTractAnalyzer.h>
+
+#include <itkImageRegionIteratorWithIndex.h>
+
+
+#include <itkShortestPathCostFunctionTbss.h>
+#include <itkShortestPathImageFilter.h>
+
+#include <iostream>
+#include <fstream>
+
+
+#include <QMessageBox>
+
+
+
+using namespace std;
+
+namespace mitk {
+
+  TractAnalyzer::TractAnalyzer()
+  {
+
+  }
+
+
+  void TractAnalyzer::BuildGraph()
+  {
+
+    int n = 0;
+    if(m_PointSetNode.IsNotNull())
+    {
+      n = m_PointSetNode->GetSize();
+      if(n==0)
+      {
+        QMessageBox msgBox;
+        msgBox.setText("No points have been set yet.");
+        msgBox.exec();
+      }
+    }
+    else{
+      QMessageBox msgBox;
+      msgBox.setText("No points have been set yet.");
+      msgBox.exec();
+    }
+
+    std::string pathDescription = "";
+    std::vector< itk::Index<3> > totalPath;
+
+    if(n>0)
+    {
+      for(int i=0; i<n-1; ++i)
+      {
+
+        mitk::ProgressBar::GetInstance()->Progress();
+
+        mitk::Point3D p = m_PointSetNode->GetPoint(i);
+        mitk::Point3D p2 = m_PointSetNode->GetPoint(i+1);
+
+
+        itk::Index<3> startPoint;
+        itk::Index<3> endPoint;
+
+        m_InputImage->GetGeometry()->WorldToIndex(p,startPoint);
+        m_InputImage->GetGeometry()->WorldToIndex(p2,endPoint);
+
+        MITK_INFO << "create roi";
+
+        std::vector< itk::Index<3> > path = CreateSegment(startPoint, endPoint);
+
+        for(std::vector< itk::Index<3> >::iterator it = path.begin();
+            it != path.end(); it++)
+        {
+          itk::Index<3> ix = *it;
+
+          if (!(ix==endPoint))
+          {
+            mitk::ProgressBar::GetInstance()->Progress();
+
+            totalPath.push_back(ix);
+            std::stringstream ss;
+            ss << ix[0] << " " << ix[1] << " " << ix[2] << "\n";
+            pathDescription += ss.str();
+          }
+          else
+          {
+            // Only when dealing with the last segment the last point should be added. This one will not occur
+            // as the first point of the next roi segment.
+            if(i == (n-2))
+            {
+              totalPath.push_back(endPoint);
+              std::stringstream ss;
+              ss << endPoint[0] << " " << endPoint[1] << " " << endPoint[2] << "\n";
+              pathDescription += ss.str();
+            }
+
+          }
+
+        }
+
+      }
+
+
+      // save pathDescription to m_PathDescription
+      m_PathDescription = pathDescription;
+
+      FloatImageType::Pointer itkImg = FloatImageType::New();
+      mitk::CastToItkImage(m_InputImage, itkImg);
+
+      CharImageType::Pointer roiImg = CharImageType::New();
+      roiImg->SetRegions(itkImg->GetLargestPossibleRegion().GetSize());
+      roiImg->SetOrigin(itkImg->GetOrigin());
+      roiImg->SetSpacing(itkImg->GetSpacing());
+      roiImg->SetDirection(itkImg->GetDirection());
+      roiImg->Allocate();
+      roiImg->FillBuffer(0);
+
+
+      std::vector< itk::Index<3> > roi;
+
+      std::vector< itk::Index<3> >::iterator it;
+      for(it = totalPath.begin();
+          it != totalPath.end();
+          it++)
+      {
+        itk::Index<3> ix = *it;
+        roiImg->SetPixel(ix, 1);
+        roi.push_back(ix);
+      }
+
+
+      m_TbssRoi = mitk::TbssRoiImage::New();
+
+      m_TbssRoi->SetRoi(roi);
+
+      m_TbssRoi->SetImage(roiImg);
+
+      m_TbssRoi->InitializeFromImage();
+
+
+
+    }
+
+
+
+  }
+
+
+
+  std::vector< itk::Index<3> > TractAnalyzer::CreateSegment(itk::Index<3> startPoint, itk::Index<3> endPoint)
+  {
+
+      typedef itk::ShortestPathImageFilter<FloatImageType, CharImageType> ShortestPathFilterType;
+      typedef itk::ShortestPathCostFunctionTbss<FloatImageType> CostFunctionType;
+
+
+      FloatImageType::Pointer meanSkeleton;
+
+      mitk::CastToItkImage(m_InputImage, meanSkeleton);
+
+      // Only use the mitk image
+
+
+
+      if(meanSkeleton)
+      {
+        CostFunctionType::Pointer costFunction = CostFunctionType::New();
+        costFunction->SetImage(meanSkeleton);
+        costFunction->SetStartIndex(startPoint);
+        costFunction->SetEndIndex(endPoint);
+        costFunction->SetThreshold(m_Threshold);
+
+        ShortestPathFilterType::Pointer pathFinder = ShortestPathFilterType::New();
+        pathFinder->SetCostFunction(costFunction);
+        pathFinder->SetFullNeighborsMode(true);
+        //pathFinder->SetCalcMode(ShortestPathFilterType::A_STAR);
+        pathFinder->SetInput(meanSkeleton);
+        pathFinder->SetStartIndex(startPoint);
+        pathFinder->SetEndIndex(endPoint);
+        pathFinder->Update();
+
+        return pathFinder->GetVectorPath();
+
+
+      }
+  }
+
+
+}
+#endif
diff --git a/Modules/DiffusionImaging/FiberTracking/Algorithms/mitkTractAnalyzer.h b/Modules/DiffusionImaging/Quantification/Algorithms/mitkTractAnalyzer.h
similarity index 56%
rename from Modules/DiffusionImaging/FiberTracking/Algorithms/mitkTractAnalyzer.h
rename to Modules/DiffusionImaging/Quantification/Algorithms/mitkTractAnalyzer.h
index 8eac2743e3..09e543dc6b 100644
--- a/Modules/DiffusionImaging/FiberTracking/Algorithms/mitkTractAnalyzer.h
+++ b/Modules/DiffusionImaging/Quantification/Algorithms/mitkTractAnalyzer.h
@@ -1,133 +1,125 @@
 /*===================================================================
 
 The Medical Imaging Interaction Toolkit (MITK)
 
 Copyright (c) German Cancer Research Center,
 Division of Medical and Biological Informatics.
 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.
 
 ===================================================================*/
 
-/*===================================================================
-
-This file is based heavily on a corresponding ITK filter.
 
-===================================================================*/
 #ifndef __mitkTractAnalyzer_h_
 #define __mitkTractAnalyzer_h_
 
-#include "FiberTrackingExports.h"
+#include "QuantificationExports.h"
 #include <itkImage.h>
 #include <itkVectorImage.h>
 #include "mitkImage.h"
 #include "mitkImageCast.h"
 
-//#include <mitkTbssImage.h>
+#include <mitkTbssRoiImage.h>
+#include <mitkPointSet.h>
 
 
 namespace mitk{
 /** \class TractAnalyzer
  */
 
 
-class FiberTracking_EXPORT TractAnalyzer
+class Quantification_EXPORT TractAnalyzer
 {
 
 public:
 
 
   TractAnalyzer();
   ~TractAnalyzer() {};
 
   typedef itk::Image<char,3> CharImageType;
   typedef itk::Image<float,3> FloatImageType;
   typedef itk::Image<float,4> ProjectionsImageType;
   typedef itk::VectorImage<float, 3> VectorImageType;
 
+  void BuildGraph();
 
 
-  void SetProjections(ProjectionsImageType::Pointer projections)
-  {
-    m_Projections = projections;
-  }
-
-  void BuildGraph(itk::Index<3> startPoint, itk::Index<3> endPoint);
-
+ /*
   std::vector< itk::Index<3> > GetPath()
   {
     return m_Path;
   }
 
-  void SetFileName(std::string fname)
-  {
-    m_FileName = fname;
-  }
-
-  void SetFileNameLong(std::string fname)
-  {
-    m_FileNameLong = fname;
-  }
-
-  void SetRoi(std::vector< itk::Index<3> > roi)
-  {
-    m_Roi = roi;
-  }
 
   CharImageType::Pointer GetRoiImage()
   {
     return m_RoiImg;
   }
 
-  void SetGroups(std::vector< std::pair<std::string, int> > groups)
+
+
+*/
+
+  mitk::TbssRoiImage::Pointer GetRoiImage()
   {
-    m_Groups = groups;
+    return m_TbssRoi;
   }
 
-  void MeasureRoi();
 
   void SetInputImage(mitk::Image::Pointer inputImage)
   {
     m_InputImage = inputImage;
   }
 
+  void SetPointSet(mitk::PointSet::Pointer pointSet)
+  {
+    m_PointSetNode = pointSet;
+  }
+
 
   void SetThreshold(double threshold)
   {
     m_Threshold = threshold;
   }
 
+  std::string GetPathDescription()
+  {
+    return m_PathDescription;
+  }
 
 protected:
 
 
-  CharImageType::Pointer m_RoiImg;
-  ProjectionsImageType::Pointer m_Projections;
+
+  std::vector< itk::Index<3> > CreateSegment(itk::Index<3> startPoint, itk::Index<3> endPoint);
+
+
+  mitk::TbssRoiImage::Pointer m_TbssRoi;
+
+  //CharImageType::Pointer m_RoiImg;
 
   mitk::Image::Pointer m_InputImage;
 
   double m_Threshold;
 
-  std::vector< itk::Index<3> > m_Path;
+  //std::vector< itk::Index<3> > m_Path;
 
-  std::string m_FileName;
+  mitk::PointSet::Pointer m_PointSetNode;
 
-  std::string m_FileNameLong; // For the regression analysis 'friendly' file
+  std::string m_PathDescription;
 
-  std::vector< std::pair<std::string, int> > m_Groups;
 
-  std::vector< itk::Index<3> > m_Roi;
 
 private:
 
 };
 
 }
 
 #endif //__itkTractAnalyzer_h_
-
diff --git a/Modules/DiffusionImaging/Quantification/files.cmake b/Modules/DiffusionImaging/Quantification/files.cmake
index 817dd133ab..5f136e1209 100644
--- a/Modules/DiffusionImaging/Quantification/files.cmake
+++ b/Modules/DiffusionImaging/Quantification/files.cmake
@@ -1,35 +1,38 @@
 set(CPP_FILES
 
   # DataStructures -> Tbss
   IODataStructures/TbssImages/mitkTbssImageSource.cpp
   IODataStructures/TbssImages/mitkTbssRoiImageSource.cpp
   IODataStructures/TbssImages/mitkNrrdTbssImageReader.cpp
   IODataStructures/TbssImages/mitkNrrdTbssImageIOFactory.cpp
   IODataStructures/TbssImages/mitkNrrdTbssRoiImageReader.cpp
   IODataStructures/TbssImages/mitkNrrdTbssRoiImageIOFactory.cpp
   IODataStructures/TbssImages/mitkTbssImage.cpp
   IODataStructures/TbssImages/mitkTbssRoiImage.cpp
   IODataStructures/TbssImages/mitkNrrdTbssImageWriter.cpp
   IODataStructures/TbssImages/mitkNrrdTbssImageWriterFactory.cpp
   IODataStructures/TbssImages/mitkNrrdTbssRoiImageWriter.cpp
   IODataStructures/TbssImages/mitkNrrdTbssRoiImageWriterFactory.cpp
   IODataStructures/TbssImages/mitkTbssImporter.cpp
 
   # DataStructures
   IODataStructures/mitkQuantificationObjectFactory.cpp
 
   # Rendering
   Rendering/mitkTbssImageMapper.cpp
+
+  #Algorithms
+  Algorithms/mitkTractAnalyzer.cpp
 )
 
 set(H_FILES
   # DataStructures
   IODataStructures/TbssImages/mitkTbssImporter.h
 
   # Rendering
   Rendering/mitkTbssImageMapper.h
 
   Algorithms/itkRegularizedIVIMLocalVariationImageFilter.h
   Algorithms/itkRegularizedIVIMReconstructionFilter.h
   Algorithms/itkRegularizedIVIMReconstructionSingleIteration.h
 )
\ No newline at end of file
diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkTractbasedSpatialStatisticsView.cpp b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkTractbasedSpatialStatisticsView.cpp
index 37eeac6637..79bc0e16dc 100644
--- a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkTractbasedSpatialStatisticsView.cpp
+++ b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkTractbasedSpatialStatisticsView.cpp
@@ -1,1427 +1,1315 @@
 /*===================================================================
 
 The Medical Imaging Interaction Toolkit (MITK)
 
 Copyright (c) German Cancer Research Center,
 Division of Medical and Biological Informatics.
 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 "berryIWorkbenchWindow.h"
 #include "berryIWorkbenchPage.h"
 #include "berryISelectionService.h"
 #include "berryConstants.h"
 #include "berryPlatformUI.h"
 
 // Qmitk
 #include "QmitkTractbasedSpatialStatisticsView.h"
 #include "QmitkStdMultiWidget.h"
 
 #include "mitkDataNodeObject.h"
 #include <itkCastImageFilter.h>
 
 // Qt
 #include <QMessageBox>
 #include <qfiledialog.h>
 #include <QDirIterator>
 #include <QInputDialog>
 #include <QDir>
 #include <QStringList>
 #include <QMessageBox>
 #include <QSortFilterProxyModel>
 #include <QClipboard>
 
 #include <itkNiftiImageIO.h>
 #include <itkImageFileReader.h>
 
 #include <iostream>
 #include <fstream>
 #include <limits>
 
 
 #include <itkMultiplyImageFilter.h>
 #include <mitkTractAnalyzer.h>
 #include <mitkTbssImporter.h>
 #include <mitkProgressBar.h>
 #include <mitkPlanarCircle.h>
 #include <mitkPlanarFigureInteractor.h>
 
 
 #include <mitkVectorImageMapper2D.h>
 #include "vtkFloatArray.h"
 #include "vtkLinearTransform.h"
 #include "vtkPoints.h"
 #include "mitkSurface.h"
 #include <vtkGlyph3D.h>
 #include "vtkArrowSource.h"
 #include "vtkUnstructuredGrid.h"
 #include "vtkPointData.h"
 #include <vtkCellArray.h>
 #include <vtkPolyLine.h>
 
 #include <boost/graph/graph_traits.hpp>
 #include <boost/graph/adjacency_list.hpp>
 #include <boost/graph/dijkstra_shortest_paths.hpp>
 
 #include <QMessageBox>
 
 #include "mitkITKImageImport.h"
 // #include "mitkImageMapperGL2D.h"
 #include "mitkVolumeDataVtkMapper3D.h"
 #include "mitkImageAccessByItk.h"
 #include "mitkTensorImage.h"
 
 #include "itkDiffusionTensor3D.h"
 
 
 #define SEARCHSIGMA 10 // length in linear voxel dimens
 
 #define MAXSEARCHLENGTH (3*SEARCHSIGMA)
 
 const std::string QmitkTractbasedSpatialStatisticsView::VIEW_ID = "org.mitk.views.tractbasedspatialstatistics";
 
 using namespace berry;
 
 
 QmitkTractbasedSpatialStatisticsView::QmitkTractbasedSpatialStatisticsView()
 : QmitkFunctionality()
 , m_Controls( 0 )
 , m_MultiWidget( NULL )
 {
 
 }
 
 QmitkTractbasedSpatialStatisticsView::~QmitkTractbasedSpatialStatisticsView()
 {
 }
 
 
 void QmitkTractbasedSpatialStatisticsView::PerformChange()
 {
   m_Controls->m_RoiPlotWidget->ModifyPlot(m_Controls->m_Segments->value(), m_Controls->m_Average->isChecked());
 }
 
 void QmitkTractbasedSpatialStatisticsView::OnSelectionChanged(std::vector<mitk::DataNode*> nodes)
 {
   //datamanager selection changed
   if (!this->IsActivated())
     return;
 
 
   // Get DataManagerSelection
   if (!this->GetDataManagerSelection().empty())
   {
     mitk::DataNode::Pointer sourceImageNode = this->GetDataManagerSelection().front();
     mitk::Image::Pointer sourceImage = dynamic_cast<mitk::Image*>(sourceImageNode->GetData());
 
     if (!sourceImage)
     {
       m_Controls->m_TbssImageLabel->setText(
         QString( sourceImageNode->GetName().c_str() ) + " is no image"
         );
 
       return;
     }
 
     // set Text
     m_Controls->m_TbssImageLabel->setText(
       QString( sourceImageNode->GetName().c_str() ) + " (" +
       QString::number(sourceImage->GetDimension()) + "D)"
       );
 
 
   }
   else
   {
     m_Controls->m_TbssImageLabel->setText("Please select an image");
   }
 
 
   bool foundTbssRoi = false;
   bool foundTbss = false;
   bool found3dImage = false;
   bool found4dImage = false;
   bool foundFiberBundle = false;
   bool foundStartRoi = false;
   bool foundEndRoi = false;
 
   mitk::TbssRoiImage* roiImage;
   mitk::TbssImage* image;
   mitk::Image* img;
   mitk::FiberBundleX* fib;
   mitk::PlanarFigure* start;
   mitk::PlanarFigure* end;
 
   m_CurrentStartRoi = NULL;
   m_CurrentEndRoi = NULL;
 
 
 
   for ( int i=0; i<nodes.size(); i++ )
   {
 
     // only look at interesting types
     // check for valid data
     mitk::BaseData* nodeData = nodes[i]->GetData();
     if( nodeData )
     {
       if(QString("TbssRoiImage").compare(nodeData->GetNameOfClass())==0)
       {
         foundTbssRoi = true;
         roiImage = static_cast<mitk::TbssRoiImage*>(nodeData);
       }
       else if (QString("TbssImage").compare(nodeData->GetNameOfClass())==0)
       {
         foundTbss = true;
         image = static_cast<mitk::TbssImage*>(nodeData);
       }
       else if(QString("Image").compare(nodeData->GetNameOfClass())==0)
       {
         img = static_cast<mitk::Image*>(nodeData);
         if(img->GetDimension() == 3)
         {
           found3dImage = true;
         }
         else if(img->GetDimension() == 4)
         {
           found4dImage = true;
         }
       }
 
       else if (QString("FiberBundleX").compare(nodeData->GetNameOfClass())==0)
       {
         foundFiberBundle = true;
         fib = static_cast<mitk::FiberBundleX*>(nodeData);
         this->m_CurrentFiberNode = nodes[i];
       }
 
 
       if(QString("PlanarCircle").compare(nodeData->GetNameOfClass())==0)
       {
         if(!foundStartRoi)
         {
           start = dynamic_cast<mitk::PlanarFigure*>(nodeData);
           this->m_CurrentStartRoi = nodes[i];
           foundStartRoi =  true;
         }
         else
         {
           end = dynamic_cast<mitk::PlanarFigure*>(nodeData);
           this->m_CurrentEndRoi = nodes[i];
           foundEndRoi = true;
         }
       }
     }
   }
 
 
   this->m_Controls->m_CreateRoi->setEnabled(found3dImage);
   this->m_Controls->m_ImportFsl->setEnabled(found4dImage);
 
 
 
   if(foundTbss && foundTbssRoi)
   {
     this->Plot(image, roiImage);
   }
 
   if(found3dImage && foundFiberBundle && foundStartRoi && foundEndRoi)
   {
     this->PlotFiberBundle(fib, img, start, end);
   }
 
   else if(found3dImage == true && foundFiberBundle)
   {
     this->PlotFiberBundle(fib, img);
   }
 
   if(found3dImage)
   {
     this->InitPointsets();
   }
 
   this->m_Controls->m_Cut->setEnabled(foundFiberBundle && foundStartRoi && foundEndRoi);
   this->m_Controls->m_SegmentLabel->setEnabled(foundFiberBundle && foundStartRoi && foundEndRoi && found3dImage);
   this->m_Controls->m_Segments->setEnabled(foundFiberBundle && foundStartRoi && foundEndRoi && found3dImage);
   this->m_Controls->m_Average->setEnabled(foundFiberBundle && foundStartRoi && foundEndRoi && found3dImage);
 
 }
 
 void QmitkTractbasedSpatialStatisticsView::InitPointsets()
 {
   // Check if PointSetStart exsits, if not create it.
   m_P1 = this->GetDefaultDataStorage()->GetNamedNode("PointSetNode");
 
   if (m_PointSetNode)
   {
     //m_PointSetNode = dynamic_cast<mitk::PointSet*>(m_P1->GetData());
     return;
   }
 
   if ((!m_P1) || (!m_PointSetNode))
   {
     // create new ones
     m_PointSetNode = mitk::PointSet::New();
     m_P1 = mitk::DataNode::New();
     m_P1->SetData( m_PointSetNode );
     m_P1->SetProperty( "name", mitk::StringProperty::New( "PointSet" ) );
     m_P1->SetProperty( "opacity", mitk::FloatProperty::New( 1 ) );
     m_P1->SetProperty( "helper object", mitk::BoolProperty::New(true) ); // CHANGE if wanted
     m_P1->SetProperty( "pointsize", mitk::FloatProperty::New( 0.1 ) );
     m_P1->SetColor( 1.0, 0.0, 0.0 );
     this->GetDefaultDataStorage()->Add(m_P1);
     m_Controls->m_PointWidget->SetPointSetNode(m_P1);
     m_Controls->m_PointWidget->SetMultiWidget(GetActiveStdMultiWidget());
   }
 }
 
 void QmitkTractbasedSpatialStatisticsView::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::QmitkTractbasedSpatialStatisticsViewControls;
     m_Controls->setupUi( parent );
     this->CreateConnections();
   }
 
 
   m_IsInitialized = false;
 
   // Table for the FSL TBSS import
   m_GroupModel = new QmitkTbssTableModel();
   m_Controls->m_GroupInfo->setModel(m_GroupModel);
 
 
 }
 
 void QmitkTractbasedSpatialStatisticsView::Activated()
 {
   QmitkFunctionality::Activated();
 
   berry::ISelection::ConstPointer sel(
     this->GetSite()->GetWorkbenchWindow()->GetSelectionService()->GetSelection("org.mitk.views.datamanager"));
   //m_CurrentSelection = sel.Cast<const IStructuredSelection>();
   //m_SelListener.Cast<TbssSelListener>()->DoSelectionChanged(sel);
 }
 
 void QmitkTractbasedSpatialStatisticsView::Deactivated()
 {
   QmitkFunctionality::Deactivated();
 }
 
 void QmitkTractbasedSpatialStatisticsView::CreateConnections()
 {
   if ( m_Controls )
   {
     connect( (QObject*)(m_Controls->m_CreateRoi), SIGNAL(clicked()), this, SLOT(CreateRoi()) );
     connect( (QObject*)(m_Controls->m_ImportFsl), SIGNAL(clicked()), this, SLOT(TbssImport()) );
     connect( (QObject*)(m_Controls->m_AddGroup), SIGNAL(clicked()), this, SLOT(AddGroup()) );
     connect( (QObject*)(m_Controls->m_RemoveGroup), SIGNAL(clicked()), this, SLOT(RemoveGroup()) );
     connect( (QObject*)(m_Controls->m_Clipboard), SIGNAL(clicked()), this, SLOT(CopyToClipboard()) );
     connect( m_Controls->m_RoiPlotWidget->m_PlotPicker, SIGNAL(selected(const QwtDoublePoint&)), SLOT(Clicked(const QwtDoublePoint&) ) );
     connect( m_Controls->m_RoiPlotWidget->m_PlotPicker, SIGNAL(moved(const QwtDoublePoint&)), SLOT(Clicked(const QwtDoublePoint&) ) );
     connect( (QObject*)(m_Controls->m_Cut), SIGNAL(clicked()), this, SLOT(Cut()) );
     connect( (QObject*)(m_Controls->m_Average), SIGNAL(stateChanged(int)), this, SLOT(PerformChange()) );
     connect( (QObject*)(m_Controls->m_Segments), SIGNAL(valueChanged(int)), this, SLOT(PerformChange()) );
 
   }
 }
 
 
 void QmitkTractbasedSpatialStatisticsView::CopyToClipboard()
 {
 
 
 
   if(m_Controls->m_RoiPlotWidget->IsPlottingFiber())
   {
     // Working with fiber bundles
     std::vector <std::vector<double> > profiles = m_Controls->m_RoiPlotWidget->GetIndividualProfiles();
 
 
 
     QString clipboardText;
     for (std::vector<std::vector<double> >::iterator it = profiles.begin(); it
                                                            != profiles.end(); ++it)
     {
       for (std::vector<double>::iterator it2 = (*it).begin(); it2 !=
             (*it).end(); ++it2)
       {
         clipboardText.append(QString("%1 \t").arg(*it2));
       }
       clipboardText.append(QString("\n"));
     }
 
 
 
     if(m_Controls->m_Average->isChecked())
     {
       std::vector<double> averages = m_Controls->m_RoiPlotWidget->GetAverageProfile();
       clipboardText.append(QString("\nAverage\n"));
 
       for (std::vector<double>::iterator it2 = averages.begin(); it2 !=
             averages.end(); ++it2)
       {
         clipboardText.append(QString("%1 \t").arg(*it2));
 
       }
     }
 
     QApplication::clipboard()->setText(clipboardText, QClipboard::Clipboard);
 
   }
 
   else{
 
     // Working with TBSS Data
     if(m_Controls->m_Average->isChecked())
     {
       std::vector<std::vector<double> > vals = m_Controls->m_RoiPlotWidget->GetVals();
       QString clipboardText;
       for (std::vector<std::vector<double> >::iterator it = vals.begin(); it
                                                              != vals.end(); ++it)
       {
         for (std::vector<double>::iterator it2 = (*it).begin(); it2 !=
               (*it).end(); ++it2)
         {
           clipboardText.append(QString("%1 \t").arg(*it2));
 
           double d = *it2;
           std::cout << d <<std::endl;
         }
         clipboardText.append(QString("\n"));
       }
 
       QApplication::clipboard()->setText(clipboardText, QClipboard::Clipboard);
     }
     else
     {
       std::vector<std::vector<double> > vals = m_Controls->m_RoiPlotWidget->GetIndividualProfiles();
       QString clipboardText;
       for (std::vector<std::vector<double> >::iterator it = vals.begin(); it
                                                              != vals.end(); ++it)
       {
         for (std::vector<double>::iterator it2 = (*it).begin(); it2 !=
               (*it).end(); ++it2)
         {
           clipboardText.append(QString("%1 \t").arg(*it2));
 
           double d = *it2;
           std::cout << d <<std::endl;
         }
         clipboardText.append(QString("\n"));
       }
       QApplication::clipboard()->setText(clipboardText, QClipboard::Clipboard);
     }
 
   }
 
 
 }
 
 void QmitkTractbasedSpatialStatisticsView::RemoveGroup()
 {
 
   QTableView *temp = static_cast<QTableView*>(m_Controls->m_GroupInfo);
  // QSortFilterProxyModel *proxy = static_cast<QSortFilterProxyModel*>(temp->model());
   QItemSelectionModel *selectionModel = temp->selectionModel();
 
   QModelIndexList indices = selectionModel->selectedRows();
 
   QModelIndex index;
 
   foreach(index, indices)
   {
     int row = index.row();
     m_GroupModel->removeRows(row, 1, QModelIndex());
   }
 
 }
 
 std::string QmitkTractbasedSpatialStatisticsView::ReadFile(std::string whatfile)
 {
   std::string s = "Select a" + whatfile;
   QFileDialog* w = new QFileDialog(this->m_Controls->m_ImportFsl, QString(s.c_str()) );
   w->setFileMode(QFileDialog::ExistingFiles);
   w->setDirectory("/home");
 
   if(whatfile == "gradient image")
   {
     w->setNameFilter("Tbss gradient images (*.tgi)");
   }
 
   // RETRIEVE SELECTION
   if ( w->exec() != QDialog::Accepted )
   {
     return "";
     MITK_INFO << "Failed to load";
   }
 
   QStringList filenames = w->selectedFiles();
   if (filenames.size() > 0)
   {
     std::string retval = filenames.at(0).toStdString();
     return retval;
   }
   return "";
 }
 
 void QmitkTractbasedSpatialStatisticsView::AddGroup()
 {
   QString group("Group");
   int number = 0;
   QPair<QString, int> pair(group, number);
   QList< QPair<QString, int> >list = m_GroupModel->getList();
 
 
   if(!list.contains(pair))
   {
     m_GroupModel->insertRows(0, 1, QModelIndex());
 
     QModelIndex index = m_GroupModel->index(0, 0, QModelIndex());
     m_GroupModel->setData(index, group, Qt::EditRole);
     index = m_GroupModel->index(0, 1, QModelIndex());
     m_GroupModel->setData(index, number, Qt::EditRole);
 
   }
   else
   {
     //QMessageBox::information(this, "Duplicate name");
   }
 
 
 }
 
 
 void QmitkTractbasedSpatialStatisticsView::TbssImport()
 {
 
   // Read groups from the interface
   mitk::TbssImporter::Pointer importer = mitk::TbssImporter::New();
 
   QList< QPair<QString, int> >list = m_GroupModel->getList();
 
   if(list.size() == 0)
   {
     QMessageBox msgBox;
     msgBox.setText("No study group information has been set yet.");
     msgBox.exec();
     return;
   }
 
 
   std::vector < std::pair<std::string, int> > groups;
   for(int i=0; i<list.size(); i++)
   {
     QPair<QString, int> pair = list.at(i);
     std::string s = pair.first.toStdString();
     int n = pair.second;
 
     std::pair<std::string, int> p;
     p.first = s;
     p.second = n;
     groups.push_back(p);
   }
 
   importer->SetGroupInfo(groups);
 
   std::string minfo = m_Controls->m_MeasurementInfo->text().toStdString();
   importer->SetMeasurementInfo(minfo);
 
 
   std::string name = "";
 
   std::vector<mitk::DataNode*> nodes = this->GetDataManagerSelection();
 
   for ( int i=0; i<nodes.size(); i++ )
   {
     if(QString("Image").compare(nodes[i]->GetData()->GetNameOfClass())==0)
     {
       mitk::Image* img = static_cast<mitk::Image*>(nodes[i]->GetData());
       if(img->GetDimension() == 4)
       {
         importer->SetImportVolume(img);
         name = nodes[i]->GetName();
       }
     }
   }
 
 
   mitk::TbssImage::Pointer tbssImage;
 
   tbssImage = importer->Import();
   name += "_tbss";
   AddTbssToDataStorage(tbssImage, name);
 
 
 }
 
 
 void QmitkTractbasedSpatialStatisticsView::AddTbssToDataStorage(mitk::Image* image, std::string name)
 {
   mitk::LevelWindow levelwindow;
   levelwindow.SetAuto( image );
   mitk::LevelWindowProperty::Pointer levWinProp = mitk::LevelWindowProperty::New();
   levWinProp->SetLevelWindow( levelwindow );
 
   mitk::DataNode::Pointer result = mitk::DataNode::New();
   result->SetProperty( "name", mitk::StringProperty::New(name) );
   result->SetData( image );
   result->SetProperty( "levelwindow", levWinProp );
 
 
   // add new image to data storage and set as active to ease further processing
   GetDefaultDataStorage()->Add( result );
 
   // show the results
   mitk::RenderingManager::GetInstance()->RequestUpdateAll();
 }
 
 void QmitkTractbasedSpatialStatisticsView::Clicked(const QwtDoublePoint& pos)
 {
   int index = (int)pos.x();
 
   if(m_Roi.size() > 0 && m_CurrentGeometry != NULL && !m_Controls->m_RoiPlotWidget->IsPlottingFiber() )
   {
 
     index = std::min( (int)m_Roi.size()-1, std::max(0, index) );
     itk::Index<3> ix = m_Roi.at(index);
 
     mitk::Vector3D i;
     i[0] = ix[0];
     i[1] = ix[1];
     i[2] = ix[2];
 
     mitk::Vector3D w;
 
     m_CurrentGeometry->IndexToWorld(i, w);
 
     mitk::Point3D origin = m_CurrentGeometry->GetOrigin();
 
     mitk::Point3D p;
     p[0] = w[0] + origin[0];
     p[1] = w[1] + origin[1];
     p[2] = w[2] + origin[2];
 
     m_MultiWidget->MoveCrossToPosition(p);
     m_Controls->m_RoiPlotWidget->drawBar(index);
   }
 
   else if(m_Controls->m_RoiPlotWidget->IsPlottingFiber() )
   {
 
     mitk::Point3D point = m_Controls->m_RoiPlotWidget->GetPositionInWorld(index);
     m_MultiWidget->MoveCrossToPosition(point);
 
 
   }
 
 
 }
 
 void QmitkTractbasedSpatialStatisticsView::Cut()
 {
   mitk::BaseData* fibData = m_CurrentFiberNode->GetData();
   mitk::FiberBundleX* fib = static_cast<mitk::FiberBundleX*>(fibData);
 
   mitk::BaseData* startData = m_CurrentStartRoi->GetData();
   mitk::PlanarFigure* startRoi = static_cast<mitk::PlanarFigure*>(startData);
   mitk::PlaneGeometry* startGeometry2D = dynamic_cast<mitk::PlaneGeometry*>( const_cast<mitk::Geometry2D*>(startRoi->GetGeometry2D()) );
 
   mitk::BaseData* endData = m_CurrentEndRoi->GetData();
   mitk::PlanarFigure* endRoi = static_cast<mitk::PlanarFigure*>(endData);
   mitk::PlaneGeometry* endGeometry2D = dynamic_cast<mitk::PlaneGeometry*>( const_cast<mitk::Geometry2D*>(endRoi->GetGeometry2D()) );
 
   mitk::Point3D startCenter = startRoi->GetWorldControlPoint(0); //center Point of start roi
   mitk::Point3D endCenter = endRoi->GetWorldControlPoint(0); //center Point of end roi
 
   mitk::FiberBundleX::Pointer inStart = fib->ExtractFiberSubset(startRoi);
   mitk::FiberBundleX::Pointer inBoth = inStart->ExtractFiberSubset(endRoi);
 
   int num = inBoth->GetNumFibers();
 
   vtkSmartPointer<vtkPolyData> fiberPolyData = inBoth->GetFiberPolyData();
   vtkCellArray* lines = fiberPolyData->GetLines();
   lines->InitTraversal();
 
 
   // initialize new vtk polydata
   vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
   vtkSmartPointer<vtkPolyData> polyData = vtkSmartPointer<vtkPolyData>::New();
   vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
 
   int pointIndex=0;
 
 
   // find start and endpoint
   for( int fiberID( 0 ); fiberID < num; fiberID++ )
   {
     vtkIdType   numPointsInCell(0);
     vtkIdType*  pointsInCell(NULL);
     lines->GetNextCell ( numPointsInCell, pointsInCell );
 
     int startId = 0;
     int endId = 0;
 
     float minDistStart = std::numeric_limits<float>::max();
     float minDistEnd = std::numeric_limits<float>::max();
 
 
     vtkSmartPointer<vtkPolyLine> polyLine = vtkSmartPointer<vtkPolyLine>::New();
     int lineIndex=0;
 
 
     for( int pointInCellID( 0 ); pointInCellID < numPointsInCell ; pointInCellID++)
     {
       double *p = fiberPolyData->GetPoint( pointsInCell[ pointInCellID ] );
 
       mitk::Point3D point;
       point[0] = p[0];
       point[1] = p[1];
       point[2] = p[2];
 
       float distanceToStart = point.EuclideanDistanceTo(startCenter);
       float distanceToEnd = point.EuclideanDistanceTo(endCenter);
 
       if(distanceToStart < minDistStart)
       {
         minDistStart = distanceToStart;
         startId = pointInCellID;
       }
 
       if(distanceToEnd < minDistEnd)
       {
         minDistEnd = distanceToEnd;
         endId = pointInCellID;
       }
 
 
 
     }
 
 
 
     /* We found the start and end points of of the part that should be plottet for
        the current fiber. now we need to plot them. If the endId is smaller than the startId the plot order
        must be reversed*/
 
     if(startId < endId)
     {
 
       double *p = fiberPolyData->GetPoint( pointsInCell[ startId ] );
 
       mitk::Vector3D p0;
       p0[0] = p[0];
       p0[1] = p[1];
       p0[2] = p[2];
 
       p = fiberPolyData->GetPoint( pointsInCell[ startId+1 ] );
 
       mitk::Vector3D p1;
       p1[0] = p[0];
       p1[1] = p[1];
       p1[2] = p[2];
 
 
       // Check if p and p2 are both on the same side of the plane
       mitk::Vector3D normal = startGeometry2D->GetNormal();
 
       mitk::Point3D pStart;
       pStart[0] = p0[0];
       pStart[1] = p0[1];
       pStart[2] = p0[2];
 
       bool startOnPositive = startGeometry2D->IsAbove(pStart);
 
       mitk::Point3D pSecond;
       pSecond[0] = p1[0];
       pSecond[1] = p1[1];
       pSecond[2] = p1[2];
 
       bool secondOnPositive = startGeometry2D->IsAbove(pSecond);
 
 
       // Calculate intersection with the plane
 
       mitk::Vector3D onPlane;
       onPlane[0] = startCenter[0];
       onPlane[1] = startCenter[1];
       onPlane[2] = startCenter[2];
 
 
       if(! (secondOnPositive ^ startOnPositive) )
       {
         /* startId and startId+1 lie on the same side of the plane, so we need
            need startId-1 to calculate the intersection with the planar figure*/
         p = fiberPolyData->GetPoint( pointsInCell[ startId-1 ] );
         p1[0] = p[0];
         p1[1] = p[1];
         p1[2] = p[2];
       }
 
 
       double d = ( (onPlane-p0)*normal) / ( (p0-p1) * normal );
 
       mitk::Vector3D newPoint = (p0-p1);
 
       newPoint[0] = d*newPoint[0] + p0[0];
       newPoint[1] = d*newPoint[1] + p0[1];
       newPoint[2] = d*newPoint[2] + p0[2];
 
 
       double insertPoint[3];
       insertPoint[0] = newPoint[0];
       insertPoint[1] = newPoint[1];
       insertPoint[2] = newPoint[2];
 
 
       // First insert the intersection with the start roi
       points->InsertNextPoint(insertPoint);
       polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
       lineIndex++;
       pointIndex++;
 
       if(! (secondOnPositive ^ startOnPositive) )
       {
         /* StartId and startId+1 lie on the same side of the plane
            so startId is also part of the ROI*/
 
         double *start = fiberPolyData->GetPoint( pointsInCell[startId] );
         points->InsertNextPoint(start);
         polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
         lineIndex++;
         pointIndex++;
 
 
       }
 
       // Insert the rest up and to including endId-1
       for( int pointInCellID( startId+1 ); pointInCellID < endId ; pointInCellID++)
       {
         // create new polyline for new polydata
         double *p = fiberPolyData->GetPoint( pointsInCell[ pointInCellID ] );
         points->InsertNextPoint(p);
 
 
         // add point to line
         polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
         lineIndex++;
         pointIndex++;
 
       }
 
 
 
       /* endId must be included if endId and endId-1 lie on the same side of the
          plane defined by endRoi*/
 
 
       p = fiberPolyData->GetPoint( pointsInCell[ endId ] );
       p0[0] = p[0];      p0[1] = p[1];      p0[2] = p[2];
 
 
       p = fiberPolyData->GetPoint( pointsInCell[ endId-1 ] );
       p1[0] = p[0];      p1[1] = p[1];      p1[2] = p[2];
 
 
       mitk::Point3D pLast;
       pLast[0] = p0[0];      pLast[1] = p0[1];      pLast[2] = p0[2];
 
       mitk::Point3D pBeforeLast;
       pBeforeLast[0] = p1[0];      pBeforeLast[1] = p1[1];      pBeforeLast[2] = p1[2];
 
       bool lastOnPositive = endGeometry2D->IsAbove(pLast);
       bool secondLastOnPositive = endGeometry2D->IsAbove(pBeforeLast);
       normal = endGeometry2D->GetNormal();
 
 
       onPlane[0] = endCenter[0];
       onPlane[1] = endCenter[1];
       onPlane[2] = endCenter[2];
 
 
 
       if(! (lastOnPositive ^ secondLastOnPositive) )
       {
         /* endId and endId-1 lie on the same side of the plane, so we need
            need endId+1 to calculate the intersection with the planar figure.
            this should exist since we know that the fiber crosses the planar figure
            endId is part of the roi so can also be included here*/
         p = fiberPolyData->GetPoint( pointsInCell[ endId+1 ] );
         p1[0] = p[0];
         p1[1] = p[1];
         p1[2] = p[2];
 
 
         double *end = fiberPolyData->GetPoint( pointsInCell[endId] );
         points->InsertNextPoint(end);
         polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
         lineIndex++;
         pointIndex++;
       }
 
       d = ( (onPlane-p0)*normal) / ( (p0-p1) * normal );
 
       newPoint = (p0-p1);
 
       newPoint[0] = d*newPoint[0] + p0[0];
       newPoint[1] = d*newPoint[1] + p0[1];
       newPoint[2] = d*newPoint[2] + p0[2];
 
       insertPoint[0] = newPoint[0];
       insertPoint[1] = newPoint[1];
       insertPoint[2] = newPoint[2];
 
       //Insert the Last Point (intersection with the end roi)
 
       points->InsertNextPoint(insertPoint);
       polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
       lineIndex++;
       pointIndex++;
 
 
     }
 
     // Need to reverse walking order
     else{
       double *p = fiberPolyData->GetPoint( pointsInCell[ startId ] );
 
       mitk::Vector3D p0;
       p0[0] = p[0];
       p0[1] = p[1];
       p0[2] = p[2];
 
       p = fiberPolyData->GetPoint( pointsInCell[ startId-1 ] );
 
       mitk::Vector3D p1;
       p1[0] = p[0];
       p1[1] = p[1];
       p1[2] = p[2];
 
 
       // Check if p and p2 are both on the same side of the plane
       mitk::Vector3D normal = startGeometry2D->GetNormal();
 
       mitk::Point3D pStart;
       pStart[0] = p0[0];
       pStart[1] = p0[1];
       pStart[2] = p0[2];
 
       bool startOnPositive = startGeometry2D->IsAbove(pStart);
 
       mitk::Point3D pSecond;
       pSecond[0] = p1[0];
       pSecond[1] = p1[1];
       pSecond[2] = p1[2];
 
       bool secondOnPositive = startGeometry2D->IsAbove(pSecond);
 
 
       // Calculate intersection with the plane
 
       mitk::Vector3D onPlane;
       onPlane[0] = startCenter[0];
       onPlane[1] = startCenter[1];
       onPlane[2] = startCenter[2];
 
 
       if(! (secondOnPositive ^ startOnPositive) )
       {
         /* startId and startId-1 lie on the same side of the plane, so we need
            need startId+1 to calculate the intersection with the planar figure*/
         p = fiberPolyData->GetPoint( pointsInCell[ startId+1 ] );
         p1[0] = p[0];
         p1[1] = p[1];
         p1[2] = p[2];
       }
 
 
       double d = ( (onPlane-p0)*normal) / ( (p0-p1) * normal );
 
       mitk::Vector3D newPoint = (p0-p1);
 
       newPoint[0] = d*newPoint[0] + p0[0];
       newPoint[1] = d*newPoint[1] + p0[1];
       newPoint[2] = d*newPoint[2] + p0[2];
 
 
       double insertPoint[3];
       insertPoint[0] = newPoint[0];
       insertPoint[1] = newPoint[1];
       insertPoint[2] = newPoint[2];
 
 
       // First insert the intersection with the start roi
       points->InsertNextPoint(insertPoint);
       polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
       lineIndex++;
       pointIndex++;
 
       if(! (secondOnPositive ^ startOnPositive) )
       {
         /* startId and startId-1 lie on the same side of the plane
            so endId is also part of the ROI*/
 
         double *start = fiberPolyData->GetPoint( pointsInCell[startId] );
         points->InsertNextPoint(start);
         polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
         lineIndex++;
         pointIndex++;
 
 
       }
 
       // Insert the rest up and to including endId-1
       for( int pointInCellID( startId-1 ); pointInCellID > endId ; pointInCellID--)
       {
         // create new polyline for new polydata
         double *p = fiberPolyData->GetPoint( pointsInCell[ pointInCellID ] );
         points->InsertNextPoint(p);
 
 
         // add point to line
         polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
         lineIndex++;
         pointIndex++;
 
       }
 
 
 
       /* startId must be included if startId and startId+ lie on the same side of the
          plane defined by endRoi*/
 
 
       p = fiberPolyData->GetPoint( pointsInCell[ endId ] );
       p0[0] = p[0];
       p0[1] = p[1];
       p0[2] = p[2];
 
 
       p = fiberPolyData->GetPoint( pointsInCell[ endId+1 ] );
       p1[0] = p[0];
       p1[1] = p[1];
       p1[2] = p[2];
 
 
       mitk::Point3D pLast;
       pLast[0] = p0[0];
       pLast[1] = p0[1];
       pLast[2] = p0[2];
 
       bool lastOnPositive = endGeometry2D->IsAbove(pLast);
 
       mitk::Point3D pBeforeLast;
       pBeforeLast[0] = p1[0];
       pBeforeLast[1] = p1[1];
       pBeforeLast[2] = p1[2];
 
       bool secondLastOnPositive = endGeometry2D->IsAbove(pBeforeLast);
 
       onPlane[0] = endCenter[0];
       onPlane[1] = endCenter[1];
       onPlane[2] = endCenter[2];
 
 
 
       if(! (lastOnPositive ^ secondLastOnPositive) )
       {
         /* endId and endId+1 lie on the same side of the plane, so we need
            need endId-1 to calculate the intersection with the planar figure.
            this should exist since we know that the fiber crosses the planar figure*/
         p = fiberPolyData->GetPoint( pointsInCell[ endId-1 ] );
         p1[0] = p[0];
         p1[1] = p[1];
         p1[2] = p[2];
 
 
         /* endId and endId+1 lie on the same side of the plane
            so startId is also part of the ROI*/
 
         double *end = fiberPolyData->GetPoint( pointsInCell[endId] );
         points->InsertNextPoint(end);
         polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
         lineIndex++;
         pointIndex++;
       }
 
       d = ( (onPlane-p0)*normal) / ( (p0-p1) * normal );
 
       newPoint = (p0-p1);
 
       newPoint[0] = d*newPoint[0] + p0[0];
       newPoint[1] = d*newPoint[1] + p0[1];
       newPoint[2] = d*newPoint[2] + p0[2];
 
 
       insertPoint[0] = newPoint[0];
       insertPoint[1] = newPoint[1];
       insertPoint[2] = newPoint[2];
 
 
 
       //Insert the Last Point (intersection with the end roi)
       points->InsertNextPoint(insertPoint);
       polyLine->GetPointIds()->InsertId(lineIndex,pointIndex);
       lineIndex++;
       pointIndex++;
 
 
 
     }
 
 
     // add polyline to vtk cell array
      cells->InsertNextCell(polyLine);
 
   }
 
 
   // Add the points to the dataset
   polyData->SetPoints(points);
 
   // Add the lines to the dataset
   polyData->SetLines(cells);
 
   mitk::FiberBundleX::Pointer cutBundle = mitk::FiberBundleX::New(polyData);
 
 
   mitk::DataNode::Pointer cutNode = mitk::DataNode::New();
   cutNode->SetData(cutBundle);
   std::string name = "fiberCut";
   cutNode->SetName(name);
   GetDataStorage()->Add(cutNode);
 
 
 
 }
 
 
 void QmitkTractbasedSpatialStatisticsView::StdMultiWidgetAvailable (QmitkStdMultiWidget &stdMultiWidget)
 {
   m_MultiWidget = &stdMultiWidget;
 }
 
 
 void QmitkTractbasedSpatialStatisticsView::StdMultiWidgetNotAvailable()
 {
   m_MultiWidget = NULL;
 }
 
 /*
 void QmitkTractbasedSpatialStatisticsView::AdjustPlotMeasure(const QString & text)
 {
   berry::ISelection::ConstPointer sel(
     this->GetSite()->GetWorkbenchWindow()->GetSelectionService()->GetSelection("org.mitk.views.datamanager"));
   m_CurrentSelection = sel.Cast<const IStructuredSelection>();
   m_SelListener.Cast<TbssSelListener>()->DoSelectionChanged(sel);
 }*/
 
 
 
 void QmitkTractbasedSpatialStatisticsView::CreateRoi()
 {
 
   // It is important to load the MeanFASkeletonMask image in MITK to make sure that point selection and
   // pathfinding is done on the same image
   //string filename = m_TbssWorkspaceManager.GetInputDir().toStdString() + "/stats/" + m_TbssWorkspaceManager.GetMeanFASkeletonMask().toStdString();
 
   // Implement a way to obtain skeleton and skeletonFA without sml workspace
 
   bool ok;
   double threshold = QInputDialog::getDouble(m_Controls->m_CreateRoi, tr("Set an FA threshold"),
                                                       tr("Threshold:"), 0.2, 0.0, 1.0, 2, &ok);
 
   if(!ok)
     return;
 
   mitk::Image::Pointer image;
 
-  /*
-  for (IStructuredSelection::iterator i = m_CurrentSelection->Begin();
-    i != m_CurrentSelection->End(); ++i)
-  {
-    // extract datatree node
-    if (mitk::DataNodeObject::Pointer nodeObj = i->Cast<mitk::DataNodeObject>())
-    {
-      mitk::DataNode::Pointer node = nodeObj->GetDataNode();
-
-      if(QString("Image").compare(node->GetData()->GetNameOfClass())==0)
-      {
-        mitk::Image* img = static_cast<mitk::Image*>(node->GetData());
-        if(img->GetDimension() == 3)
-        {
-          image = img;
-        }
-      }
-    }
-  }
-  */
 
   std::vector<mitk::DataNode*> nodes = this->GetDataManagerSelection();
 
   for ( int i=0; i<nodes.size(); i++ )
   {
     if(QString("Image").compare(nodes[i]->GetData()->GetNameOfClass())==0)
     {
         mitk::Image* img = static_cast<mitk::Image*>(nodes[i]->GetData());
         if(img->GetDimension() == 3)
         {
           image = img;
         }
     }
   }
 
 
 
   if(image.IsNull())
   {
     return;
   }
 
   mitk::TractAnalyzer analyzer;
   analyzer.SetInputImage(image);
   analyzer.SetThreshold(threshold);
 
+  // Set Pointset to analyzer
+  analyzer.SetPointSet(m_PointSetNode);
 
+  // Run Analyzer
+  analyzer.BuildGraph();
 
-  int n = 0;
-  if(m_PointSetNode.IsNotNull())
-  {
-    n = m_PointSetNode->GetSize();
-    if(n==0)
-    {
-      QMessageBox msgBox;
-      msgBox.setText("No points have been set yet.");
-      msgBox.exec();
-    }
-  }
-  else{
-    QMessageBox msgBox;
-    msgBox.setText("No points have been set yet.");
-    msgBox.exec();
-  }
-
-  std::string pathDescription = "";
-  std::vector< itk::Index<3> > totalPath;
-
-  if(n>0)
-  {
-    for(int i=0; i<n-1; ++i)
-    {
-      mitk::Point3D p = m_PointSetNode->GetPoint(i);
-      mitk::Point3D p2 = m_PointSetNode->GetPoint(i+1);
-
-
-      itk::Index<3> StartPoint; mitk::ProgressBar::GetInstance()->Progress();
-
-      itk::Index<3> EndPoint;
-      image->GetGeometry()->WorldToIndex(p,StartPoint);
-      image->GetGeometry()->WorldToIndex(p2,EndPoint);
-      MITK_INFO << "create roi";
+  // Obtain tbss roi image from analyzer
+  mitk::TbssRoiImage::Pointer tbssRoi = analyzer.GetRoiImage();
 
+  tbssRoi->SetStructure(m_Controls->m_Structure->text().toStdString());
 
-      analyzer.BuildGraph(StartPoint, EndPoint);
-      std::vector< itk::Index<3> > path = analyzer.GetPath();
+  // get path description and set to interface
+  std::string pathDescription = analyzer.GetPathDescription();
+  m_Controls->m_PathTextEdit->setPlainText(QString(pathDescription.c_str()));
 
 
-      for(std::vector< itk::Index<3> >::iterator it = path.begin();
-          it != path.end(); it++)
-      {
-        itk::Index<3> ix = *it;
-
-        if (!(ix==EndPoint))
-        { mitk::ProgressBar::GetInstance()->Progress();
-
-          totalPath.push_back(ix);
-          std::stringstream ss;
-          ss << ix[0] << " " << ix[1] << " " << ix[2] << "\n";
-          pathDescription += ss.str();
-        }
-        else
-        {
-          // Only when dealing with the last segment the last point should be added. This one will not occur
-          // as the first point of the next roi segment.
-          if(i == (n-2))
-          {
-            totalPath.push_back(EndPoint);
-            std::stringstream ss;
-            ss << EndPoint[0] << " " << EndPoint[1] << " " << EndPoint[2] << "\n";
-            pathDescription += ss.str();
-          }
-
-        }
-      }
-
-
-
-    }
+  // Add roi image to datastorage
+  AddTbssToDataStorage(tbssRoi, m_Controls->m_RoiName->text().toStdString());
 
-    m_Controls->m_PathTextEdit->setPlainText(QString(pathDescription.c_str()));
 
 
-    FloatImageType::Pointer itkImg = FloatImageType::New();
-    mitk::CastToItkImage(image, itkImg);
-
-    CharImageType::Pointer roiImg = CharImageType::New();
-    roiImg->SetRegions(itkImg->GetLargestPossibleRegion().GetSize());
-    roiImg->SetOrigin(itkImg->GetOrigin());
-    roiImg->SetSpacing(itkImg->GetSpacing());
-    roiImg->SetDirection(itkImg->GetDirection());
-    roiImg->Allocate();
-    roiImg->FillBuffer(0);
-
-
-    std::vector< itk::Index<3> > roi;
-
-    std::vector< itk::Index<3> >::iterator it;
-    for(it = totalPath.begin();
-        it != totalPath.end();
-        it++)
-    {
-      itk::Index<3> ix = *it;
-      roiImg->SetPixel(ix, 1);
-      roi.push_back(ix);
-    }
-
-
-    mitk::TbssRoiImage::Pointer tbssRoi = mitk::TbssRoiImage::New();
-
-    tbssRoi->SetRoi(roi);
-    tbssRoi->SetImage(roiImg);
-    tbssRoi->SetStructure(m_Controls->m_Structure->text().toStdString());
-    tbssRoi->InitializeFromImage();
-
-    AddTbssToDataStorage(tbssRoi, m_Controls->m_RoiName->text().toStdString());
-
-
-  }
 
 }
 
 
 
 void QmitkTractbasedSpatialStatisticsView:: PlotFiberBundle(mitk::FiberBundleX *fib, mitk::Image* img,
                                                            mitk::PlanarFigure* startRoi, mitk::PlanarFigure* endRoi)
 {
   bool avg = m_Controls->m_Average->isChecked();
   int segments = m_Controls->m_Segments->value();
   m_Controls->m_RoiPlotWidget->PlotFiberBetweenRois(fib, img, startRoi ,endRoi, avg, segments);
   m_Controls->m_RoiPlotWidget->SetPlottingFiber(true);
   mitk::RenderingManager::GetInstance()->ForceImmediateUpdateAll();
 }
 
 
 void QmitkTractbasedSpatialStatisticsView::Plot(mitk::TbssImage* image, mitk::TbssRoiImage* roiImage)
 {
   if(m_Controls->m_TabWidget->currentWidget() == m_Controls->m_MeasureTAB)
   {
 
     std::vector< itk::Index<3> > roi = roiImage->GetRoi();
     m_Roi = roi;
     m_CurrentGeometry = image->GetGeometry();
 
 
     std::string resultfile = "";
     std::string structure = roiImage->GetStructure();
 
     m_Controls->m_RoiPlotWidget->SetGroups(image->GetGroupInfo());
     m_Controls->m_RoiPlotWidget->SetProjections(image->GetImage());
     m_Controls->m_RoiPlotWidget->SetRoi(roi);
     m_Controls->m_RoiPlotWidget->SetStructure(structure);
     m_Controls->m_RoiPlotWidget->SetMeasure( image->GetMeasurementInfo() );
     m_Controls->m_RoiPlotWidget->DrawProfiles(resultfile);
   }
 
   m_Controls->m_RoiPlotWidget->SetPlottingFiber(false);
 
 }
 
 
 
 
 void QmitkTractbasedSpatialStatisticsView::Masking()
 {
   //QString filename = m_Controls->m_WorkingDirectory->text();
   QString filename = "E:/Experiments/tbss";
   QString faFiles = filename + "/AxD";
   QString maskFiles = filename + "/bin_masks";
 
 
   QDirIterator faDirIt(faFiles, QDir::Files | QDir::NoSymLinks, QDirIterator::Subdirectories);
   QDirIterator maskDirIt(maskFiles, QDir::Files | QDir::NoSymLinks, QDirIterator::Subdirectories);
 
   std::vector<std::string> faFilenames;
   std::vector<std::string> maskFilenames;
   std::vector<std::string> outputFilenames;
 
   while(faDirIt.hasNext() && maskDirIt.hasNext())
   {
     faDirIt.next();
     maskDirIt.next();
     if((faDirIt.fileInfo().completeSuffix() == "nii"
        || faDirIt.fileInfo().completeSuffix() == "mhd"
        || faDirIt.fileInfo().completeSuffix() == "nii.gz")
        && (maskDirIt.fileInfo().completeSuffix() == "nii"
        || maskDirIt.fileInfo().completeSuffix() == "mhd"
        || maskDirIt.fileInfo().completeSuffix() == "nii.gz"))
     {
       faFilenames.push_back(faDirIt.filePath().toStdString());
       outputFilenames.push_back(faDirIt.fileName().toStdString());
       maskFilenames.push_back(maskDirIt.filePath().toStdString());
     }
   }
 
   std::vector<std::string>::iterator faIt = faFilenames.begin();
   std::vector<std::string>::iterator maskIt = maskFilenames.begin();
   std::vector<std::string>::iterator outputIt = outputFilenames.begin();
 
   // Now multiply all FA images with their corresponding masks
 
   QString outputDir = filename;
   while(faIt != faFilenames.end() && maskIt != maskFilenames.end() && outputIt != outputFilenames.end())
   {
     std::cout << "Mask " << *faIt << " with " << *maskIt << std::endl;
 
     typedef itk::MultiplyImageFilter<FloatImageType, CharImageType, FloatImageType> MultiplicationFilterType;
 
     FloatReaderType::Pointer floatReader = FloatReaderType::New();
     CharReaderType::Pointer charReader = CharReaderType::New();
 
     floatReader->SetFileName(*faIt);
     //floatReader->Update();
     //FloatImageType::Pointer faImage = floatReader->GetOutput();
 
     charReader->SetFileName(*maskIt);
     //charReader->Update();
     // CharImageType::Pointer maskImage = charReader->GetOutput();
 
     MultiplicationFilterType::Pointer multiplicationFilter = MultiplicationFilterType::New();
     multiplicationFilter->SetInput1(floatReader->GetOutput());
     multiplicationFilter->SetInput2(charReader->GetOutput());
     multiplicationFilter->Update();
 
     //FloatImageType::Pointer maskedImage = FloatImageType::New();
     //maskedImage = MultiplicationFilter->GetOutput();
 
     FloatWriterType::Pointer floatWriter = FloatWriterType::New();
     std::string s = faFiles.toStdString().append("/"+*outputIt);
     floatWriter->SetFileName(s.c_str());
     floatWriter->SetInput(multiplicationFilter->GetOutput());
     floatWriter->Update();
 
     ++faIt;
     ++maskIt;
     ++outputIt;
   }
 }