diff --git a/Modules/Core/src/Algorithms/mitkExtractSliceFilter.cpp b/Modules/Core/src/Algorithms/mitkExtractSliceFilter.cpp
index ec1081234d..6daf8f416d 100644
--- a/Modules/Core/src/Algorithms/mitkExtractSliceFilter.cpp
+++ b/Modules/Core/src/Algorithms/mitkExtractSliceFilter.cpp
@@ -1,504 +1,504 @@
 /*============================================================================
 
 The Medical Imaging Interaction Toolkit (MITK)
 
 Copyright (c) German Cancer Research Center (DKFZ)
 All rights reserved.
 
 Use of this source code is governed by a 3-clause BSD license that can be
 found in the LICENSE file.
 
 ============================================================================*/
 
 #include "mitkExtractSliceFilter.h"
 
 #include <mitkAbstractTransformGeometry.h>
 #include <mitkPlaneClipping.h>
 
 #include <vtkGeneralTransform.h>
 #include <vtkImageChangeInformation.h>
 #include <vtkImageData.h>
 #include <vtkImageExtractComponents.h>
 #include <vtkLinearTransform.h>
 
 mitk::ExtractSliceFilter::ExtractSliceFilter(vtkImageReslice *reslicer): m_XMin(0), m_XMax(0), m_YMin(0), m_YMax(0)
 {
   if (reslicer == nullptr)
   {
     m_Reslicer = vtkSmartPointer<vtkImageReslice>::New();
   }
   else
   {
     m_Reslicer = reslicer;
   }
 
   m_TimeStep = 0;
   m_Reslicer->ReleaseDataFlagOn();
   m_InterpolationMode = ExtractSliceFilter::RESLICE_NEAREST;
   m_ResliceTransform = nullptr;
   m_InPlaneResampleExtentByGeometry = false;
   m_OutPutSpacing = new mitk::ScalarType[2];
   m_OutputDimension = 2;
   m_ZSpacing = 1.0;
   m_ZMin = 0;
   m_ZMax = 0;
   m_VtkOutputRequested = false;
   m_BackgroundLevel = -32768.0;
   m_Component = 0;
 }
 
 mitk::ExtractSliceFilter::~ExtractSliceFilter()
 {
   m_ResliceTransform = nullptr;
   m_WorldGeometry = nullptr;
   delete[] m_OutPutSpacing;
 }
 
 void mitk::ExtractSliceFilter::GenerateOutputInformation()
 {
   Image::ConstPointer input = this->GetInput();
 
   if (input.IsNull())
   {
     return;
   }
 
   if (nullptr == m_WorldGeometry)
   {
     return;
   }
 
   Vector3D right, bottom;
   double widthInMM, heightInMM;
   Vector2D extent;
 
   // set the geometry from current worldgeometry for the resultimage
   // this is needed that the image has the correct mitk geometry
   // the sliceGeometry is the Geometry of the result slice
   PlaneGeometry::Pointer sliceGeometry = m_WorldGeometry->Clone();
 
   sliceGeometry->GetIndexToWorldTransform()->SetMatrix(m_WorldGeometry->GetIndexToWorldTransform()->GetMatrix());
 
   // the origin of the worldGeometry is transformed to center based coordinates to be an imageGeometry
   Point3D sliceOrigin = sliceGeometry->GetOrigin();
 
   auto abstractGeometry =
     dynamic_cast<const AbstractTransformGeometry *>(m_WorldGeometry.GetPointer());
 
   if (abstractGeometry != nullptr)
   {
     extent[0] = abstractGeometry->GetParametricExtent(0);
     extent[1] = abstractGeometry->GetParametricExtent(1);
 
     widthInMM = abstractGeometry->GetParametricExtentInMM(0);
     heightInMM = abstractGeometry->GetParametricExtentInMM(1);
 
     right = abstractGeometry->GetPlane()->GetAxisVector(0);
     bottom = abstractGeometry->GetPlane()->GetAxisVector(1);
   }
   else
   {
     // if the worldGeomatry is a PlaneGeometry everything is straight forward
     right = m_WorldGeometry->GetAxisVector(0);
     bottom = m_WorldGeometry->GetAxisVector(1);
 
     if (m_InPlaneResampleExtentByGeometry)
     {
       // Resampling grid corresponds to the current world geometry. This
       // means that the spacing of the output 2D image depends on the
       // currently selected world geometry, and *not* on the image itself.
       extent[0] = m_WorldGeometry->GetExtent(0);
       extent[1] = m_WorldGeometry->GetExtent(1);
     }
     else
     {
       const TimeGeometry *inputTimeGeometry = input->GetTimeGeometry();
       if ((inputTimeGeometry == nullptr) || (inputTimeGeometry->CountTimeSteps() <= 0))
       {
         itkWarningMacro(<< "Error reading input image TimeGeometry.");
         return;
       }
 
       // Resampling grid corresponds to the input geometry. This means that
       // the spacing of the output 2D image is directly derived from the
       // associated input image, regardless of the currently selected world
       // geometry.
       Vector3D rightInIndex, bottomInIndex;
       inputTimeGeometry->GetGeometryForTimeStep(m_TimeStep)->WorldToIndex(right, rightInIndex);
       inputTimeGeometry->GetGeometryForTimeStep(m_TimeStep)->WorldToIndex(bottom, bottomInIndex);
       extent[0] = rightInIndex.GetNorm();
       extent[1] = bottomInIndex.GetNorm();
     }
 
     // Get the extent of the current world geometry and calculate resampling
     // spacing therefrom.
     widthInMM = m_WorldGeometry->GetExtentInMM(0);
     heightInMM = m_WorldGeometry->GetExtentInMM(1);
   }
 
   right.Normalize();
   bottom.Normalize();
 
   m_OutPutSpacing[0] = widthInMM / extent[0];
   m_OutPutSpacing[1] = heightInMM / extent[1];
 
   /*========== BEGIN setup extent of the slice ==========*/
   int xMin, xMax, yMin, yMax;
 
   xMin = yMin = 0;
   xMax = static_cast<int>(extent[0]);
   yMax = static_cast<int>(extent[1]);
 
   if (m_WorldGeometry->GetReferenceGeometry())
   {
     double sliceBounds[6];
     for (auto &sliceBound : sliceBounds)
     {
       sliceBound = 0.0;
     }
 
     if (this->GetClippedPlaneBounds(m_WorldGeometry->GetReferenceGeometry(), m_WorldGeometry, sliceBounds))
     {
       // Calculate output extent (integer values)
       xMin = static_cast<int>(sliceBounds[0] / m_OutPutSpacing[0] + 0.5);
       xMax = static_cast<int>(sliceBounds[1] / m_OutPutSpacing[0] + 0.5);
       yMin = static_cast<int>(sliceBounds[2] / m_OutPutSpacing[1] + 0.5);
       yMax = static_cast<int>(sliceBounds[3] / m_OutPutSpacing[1] + 0.5);
     } // ELSE we use the default values
   }
 
 
   sliceOrigin += right * (m_OutPutSpacing[0] * 0.5);
   sliceOrigin += bottom * (m_OutPutSpacing[1] * 0.5);
 
   // a worldGeometry is no imageGeometry, thus it is manually set to true
   sliceGeometry->ImageGeometryOn();
 
   /*At this point we have to adjust the geometry because the origin isn't correct.
   The wrong origin is related to the rotation of the current world geometry plane.
   This causes errors on transferring world to index coordinates. We just shift the
   origin in each direction about the amount of the expanding (needed while rotating
   the plane).
   */
   Vector3D axis0 = sliceGeometry->GetAxisVector(0);
   Vector3D axis1 = sliceGeometry->GetAxisVector(1);
   axis0.Normalize();
   axis1.Normalize();
 
   // adapt the origin. Note that for orthogonal planes the minima are '0' and thus the origin stays the same.
   sliceOrigin += (axis0 * (xMin * m_OutPutSpacing[0])) + (axis1 * (yMin * m_OutPutSpacing[1]));
 
   sliceGeometry->SetOrigin(sliceOrigin);
 
   /*the bounds as well as the extent of the worldGeometry are not adapted correctly during crosshair rotation.
   This is only a quick fix and has to be evaluated.
   The new bounds are set via the max values of the calculated slice extent. It will look like [ 0, x, 0, y, 0, 1].
   */
   mitk::BoundingBox::BoundsArrayType boundsCopy;
   boundsCopy[0] = boundsCopy[2] = boundsCopy[4] = 0;
   boundsCopy[5] = 1;
   boundsCopy[1] = std::max(xMax - xMin, 1);
   boundsCopy[3] = std::max(yMax - yMin, 1);
   sliceGeometry->SetBounds(boundsCopy);
 
   sliceGeometry->Modified();
 
   Image::Pointer output = this->GetOutput();
-  output->Initialize(input->GetPixelType(), 2, *sliceGeometry);
+  output->Initialize(input->GetPixelType(), 1, *sliceGeometry);
 
   m_XMin = xMin;
   m_XMax = xMax;
   m_YMin = yMin;
   m_YMax = yMax;
   m_Right = right;
   m_Bottom = bottom;
 }
 
 void mitk::ExtractSliceFilter::GenerateInputRequestedRegion()
 {
   // As we want all pixel information fo the image in our plane, the requested region
   // is set to the largest possible region in the image.
   // This is needed because an oblique plane has a larger extent then the image
   // and the in pipeline it is checked via PropagateResquestedRegion(). But the
   // extent of the slice is actually fitting because it is oblique within the image.
   ImageToImageFilter::InputImagePointer input = this->GetInput();
   input->SetRequestedRegionToLargestPossibleRegion();
 }
 
 mitk::ScalarType *mitk::ExtractSliceFilter::GetOutputSpacing()
 {
   return m_OutPutSpacing;
 }
 
 void mitk::ExtractSliceFilter::GenerateData()
 {
   mitk::Image *input = this->GetInput();
 
   if (!input)
   {
     MITK_ERROR << "mitk::ExtractSliceFilter: No input image available. Please set the input!" << std::endl;
     itkExceptionMacro("mitk::ExtractSliceFilter: No input image available. Please set the input!");
     return;
   }
 
   if (!m_WorldGeometry)
   {
     MITK_ERROR << "mitk::ExtractSliceFilter: No Geometry for reslicing available." << std::endl;
     itkExceptionMacro("mitk::ExtractSliceFilter: No Geometry for reslicing available.");
     return;
   }
 
   const TimeGeometry *inputTimeGeometry = this->GetInput()->GetTimeGeometry();
   if ((inputTimeGeometry == nullptr) || (inputTimeGeometry->CountTimeSteps() <= 0))
   {
     itkWarningMacro(<< "Error reading input image TimeGeometry.");
     return;
   }
 
   // is it a valid timeStep?
   if (inputTimeGeometry->IsValidTimeStep(m_TimeStep) == false)
   {
     itkWarningMacro(<< "This is not a valid timestep: " << m_TimeStep);
     return;
   }
 
   // check if there is something to display.
   if (!input->IsVolumeSet(m_TimeStep))
   {
     itkWarningMacro(<< "No volume data existent at given timestep " << m_TimeStep);
     return;
   }
 
   /*================#BEGIN setup vtkImageReslice properties================*/
   Point3D origin;
   Vector3D normal;
 
   const auto *planeGeometry = dynamic_cast<const PlaneGeometry *>(m_WorldGeometry.GetPointer());
   // Code for curved planes, mostly taken 1:1 from imageVtkMapper2D and not tested yet.
   // Do we have an AbstractTransformGeometry?
   // This is the case for AbstractTransformGeometry's (e.g. a ThinPlateSplineCurvedGeometry )
   const auto *abstractGeometry =
     dynamic_cast<const AbstractTransformGeometry *>(m_WorldGeometry.GetPointer());
 
   if (abstractGeometry != nullptr)
   {
     m_ResliceTransform = abstractGeometry;
 
     origin = abstractGeometry->GetPlane()->GetOrigin();
 
     normal = abstractGeometry->GetPlane()->GetNormal();
     normal.Normalize();
 
     // Use a combination of the InputGeometry *and* the possible non-rigid
     // AbstractTransformGeometry for reslicing the 3D Image
     vtkSmartPointer<vtkGeneralTransform> composedResliceTransform = vtkSmartPointer<vtkGeneralTransform>::New();
     composedResliceTransform->Identity();
     composedResliceTransform->Concatenate(
       inputTimeGeometry->GetGeometryForTimeStep(m_TimeStep)->GetVtkTransform()->GetLinearInverse());
     composedResliceTransform->Concatenate(abstractGeometry->GetVtkAbstractTransform());
 
     m_Reslicer->SetResliceTransform(composedResliceTransform);
 
     // Set background level to BLACK instead of translucent, to avoid
     // boundary artifacts (see PlaneGeometryDataVtkMapper3D)
     // Note: Backgroundlevel was hardcoded before to -1023
     m_Reslicer->SetBackgroundLevel(m_BackgroundLevel);
   }
   else
   {
     if (planeGeometry != nullptr)
     {
       // if the worldGeomatry is a PlaneGeometry everything is straight forward
 
       origin = planeGeometry->GetOrigin();
       normal = planeGeometry->GetNormal();
       normal.Normalize();
 
       /*
       * Transform the origin to center based coordinates.
       * Note:
       * This is needed besause vtk's origin is center based too (!!!) ( see 'The VTK book' page 88 )
       * and the worldGeometry surrouding the image is no imageGeometry. So the worldGeometry
       * has its origin at the corner of the voxel and needs to be transformed.
       */
       origin += m_Right * (m_OutPutSpacing[0] * 0.5);
       origin += m_Bottom * (m_OutPutSpacing[1] * 0.5);
 
       // set the tranform for reslicing.
       // Use inverse transform of the input geometry for reslicing the 3D image.
       // This is needed if the image volume already transformed
       if (m_ResliceTransform.IsNotNull())
         m_Reslicer->SetResliceTransform(m_ResliceTransform->GetVtkTransform()->GetLinearInverse());
 
       // Set background level to TRANSLUCENT (see PlaneGeometryDataVtkMapper3D),
       // else the background of the image turns out gray
       // Note: Backgroundlevel was hardcoded to -32768
       m_Reslicer->SetBackgroundLevel(m_BackgroundLevel);
     }
     else
     {
       itkExceptionMacro("mitk::ExtractSliceFilter: No fitting geometry for reslice axis!");
       return;
     }
   }
 
   if (m_ResliceTransform.IsNotNull())
   {
     // if the resliceTransform is set the reslice axis are recalculated.
     // Thus the geometry information is not fitting. Therefor a unitSpacingFilter
     // is used to set up a global spacing of 1 and compensate the transform.
     vtkSmartPointer<vtkImageChangeInformation> unitSpacingImageFilter =
       vtkSmartPointer<vtkImageChangeInformation>::New();
     unitSpacingImageFilter->ReleaseDataFlagOn();
 
     unitSpacingImageFilter->SetOutputSpacing(1.0, 1.0, 1.0);
     unitSpacingImageFilter->SetInputData(input->GetVtkImageData(m_TimeStep));
 
     m_Reslicer->SetInputConnection(unitSpacingImageFilter->GetOutputPort());
   }
   else
   {
     // if no transform is set the image can be used directly
     m_Reslicer->SetInputData(input->GetVtkImageData(m_TimeStep));
   }
 
   /*setup the plane where vktImageReslice extracts the slice*/
 
   // ResliceAxesOrigin is the anchor point of the plane
   double originInVtk[3];
   itk2vtk(origin, originInVtk);
   m_Reslicer->SetResliceAxesOrigin(originInVtk);
 
   // the cosines define the plane: x and y are the direction vectors, n is the planes normal
   // this specifies a matrix 3x3
   //  x1 y1 n1
   //  x2 y2 n2
   //  x3 y3 n3
   double cosines[9];
 
   vnl2vtk(m_Right.GetVnlVector(), cosines); // x
 
   vnl2vtk(m_Bottom.GetVnlVector(), cosines + 3); // y
 
   vnl2vtk(normal.GetVnlVector(), cosines + 6); // n
 
   m_Reslicer->SetResliceAxesDirectionCosines(cosines);
 
   // we only have one slice, not a volume
   m_Reslicer->SetOutputDimensionality(m_OutputDimension);
 
   // set the interpolation mode for slicing
   switch (this->m_InterpolationMode)
   {
     case RESLICE_NEAREST:
       m_Reslicer->SetInterpolationModeToNearestNeighbor();
       break;
     case RESLICE_LINEAR:
       m_Reslicer->SetInterpolationModeToLinear();
       break;
     case RESLICE_CUBIC:
       m_Reslicer->SetInterpolationModeToCubic();
       break;
     default:
       // the default interpolation used by mitk
       m_Reslicer->SetInterpolationModeToNearestNeighbor();
   }
 
   /*========== BEGIN setup extent of the slice ==========*/
   // Set the output extents! First included pixel index and last included pixel index
   // xMax and yMax are one after the last pixel. so they have to be decremented by 1.
   // In case we have a 2D image, xMax or yMax might be 0. in this case, do not decrement, but take 0.
 
   m_Reslicer->SetOutputExtent(m_XMin, std::max(0, m_XMax - 1), m_YMin, std::max(0, m_YMax - 1), m_ZMin, m_ZMax);
   /*========== END setup extent of the slice ==========*/
 
   m_Reslicer->SetOutputOrigin(0.0, 0.0, 0.0);
 
   m_Reslicer->SetOutputSpacing(m_OutPutSpacing[0], m_OutPutSpacing[1], m_ZSpacing);
 
   // TODO check the following lines, they are responsible whether vtk error outputs appear or not
   m_Reslicer->UpdateWholeExtent(); // this produces a bad allocation error for 2D images
   // m_Reslicer->GetOutput()->UpdateInformation();
   // m_Reslicer->GetOutput()->SetUpdateExtentToWholeExtent();
 
   // start the pipeline
   m_Reslicer->Update();
   /*================ #END setup vtkImageReslice properties================*/
 
   if (m_VtkOutputRequested)
   {
     // no conversion to mitk
     // no mitk geometry will be set, as the output is vtkImageData only!!!
     // no image component will be extracted, as the caller might need the whole multi-component image as vtk output
     return;
   }
   else
   {
     auto reslicedImage = vtkSmartPointer<vtkImageData>::New();
     reslicedImage = m_Reslicer->GetOutput();
 
     if (nullptr == reslicedImage)
     {
       itkWarningMacro(<< "Reslicer returned empty image");
       return;
     }
 
     /*================ #BEGIN Extract component from image slice ================*/
     int numberOfScalarComponent = reslicedImage->GetNumberOfScalarComponents();
     if (numberOfScalarComponent > 1 && static_cast<unsigned int>(numberOfScalarComponent) >= m_Component)
     {
       // image has more than one component, extract the correct component information with the given 'component' parameter
       auto vectorComponentExtractor = vtkSmartPointer<vtkImageExtractComponents>::New();
       vectorComponentExtractor->SetInputData(reslicedImage);
       vectorComponentExtractor->SetComponents(m_Component);
       vectorComponentExtractor->Update();
 
       reslicedImage = vectorComponentExtractor->GetOutput();
     }
     /*================ #END Extract component from image slice ================*/
 
     /*================ #BEGIN Convert the slice to an mitk::Image ================*/
     mitk::Image::Pointer resultImage = GetOutput();
 
     /*Temporary store the geometry that is already correct (set in GeneratOutputInformation())
      but will be reset due to initialize.*/
     mitk::BaseGeometry::Pointer resultGeometry = resultImage->GetGeometry();
 
     // initialize resultimage with the specs of the vtkImageData object returned from vtkImageReslice
     if (reslicedImage->GetDataDimension() == 1)
     {
       // If original image was 2D, the slice might have an y extent of 0.
       // Still i want to ensure here that Image is 2D
       resultImage->Initialize(reslicedImage, 1, -1, -1, 1);
     }
     else
     {
       resultImage->Initialize(reslicedImage);
     }
 
     // transfer the voxel data
     resultImage->SetVolume(reslicedImage->GetScalarPointer());
 
     /*================ #END Convert the slice to an mitk::Image ================*/
 
     resultImage->SetGeometry(resultGeometry);
   }
 }
 
 bool mitk::ExtractSliceFilter::GetClippedPlaneBounds(double bounds[6])
 {
   if (!m_WorldGeometry || !this->GetInput())
     return false;
 
   return this->GetClippedPlaneBounds(
     m_WorldGeometry->GetReferenceGeometry(), m_WorldGeometry, bounds);
 }
 
 bool mitk::ExtractSliceFilter::GetClippedPlaneBounds(const BaseGeometry *boundingGeometry,
                                                      const PlaneGeometry *planeGeometry,
                                                      double *bounds)
 {
   bool b = mitk::PlaneClipping::CalculateClippedPlaneBounds(boundingGeometry, planeGeometry, bounds);
 
   return b;
 }
diff --git a/Modules/Segmentation/Interactions/mitkRegionGrowingTool.cpp b/Modules/Segmentation/Interactions/mitkRegionGrowingTool.cpp
index 1c57d02b70..689b0dbded 100644
--- a/Modules/Segmentation/Interactions/mitkRegionGrowingTool.cpp
+++ b/Modules/Segmentation/Interactions/mitkRegionGrowingTool.cpp
@@ -1,640 +1,663 @@
 /*============================================================================
 
 The Medical Imaging Interaction Toolkit (MITK)
 
 Copyright (c) German Cancer Research Center (DKFZ)
 All rights reserved.
 
 Use of this source code is governed by a 3-clause BSD license that can be
 found in the LICENSE file.
 
 ============================================================================*/
 
 #include "mitkRegionGrowingTool.h"
 #include "mitkApplicationCursor.h"
 #include "mitkBaseRenderer.h"
 #include "mitkImageDataItem.h"
 #include "mitkImageToContourModelFilter.h"
 #include "mitkOverwriteSliceImageFilter.h"
 #include "mitkRegionGrowingTool.xpm"
 #include "mitkRenderingManager.h"
 #include "mitkToolManager.h"
 
 #include "mitkExtractDirectedPlaneImageFilterNew.h"
 #include "mitkLabelSetImage.h"
 #include "mitkOverwriteDirectedPlaneImageFilter.h"
 
 // us
 #include <usGetModuleContext.h>
 #include <usModule.h>
 #include <usModuleContext.h>
 #include <usModuleResource.h>
 
 // ITK
 #include "mitkITKImageImport.h"
 #include "mitkImageAccessByItk.h"
 #include <itkConnectedComponentImageFilter.h>
 #include <itkConnectedThresholdImageFilter.h>
 #include <itkImageRegionIteratorWithIndex.h>
 #include <itkNeighborhoodIterator.h>
 
 #include <itkImageDuplicator.h>
 
+#include <limits>
+
 namespace mitk
 {
   MITK_TOOL_MACRO(MITKSEGMENTATION_EXPORT, RegionGrowingTool, "Region growing tool");
 }
 
 #define ROUND(a) ((a) > 0 ? (int)((a) + 0.5) : -(int)(0.5 - (a)))
 
 mitk::RegionGrowingTool::RegionGrowingTool()
   : FeedbackContourTool("PressMoveRelease"),
     m_SeedValue(0),
     m_ScreenYDifference(0),
     m_ScreenXDifference(0),
     m_MouseDistanceScaleFactor(0.5),
     m_PaintingPixelValue(0),
     m_FillFeedbackContour(true),
     m_ConnectedComponentValue(1)
 {
 }
 
 mitk::RegionGrowingTool::~RegionGrowingTool()
 {
 }
 
 void mitk::RegionGrowingTool::ConnectActionsAndFunctions()
 {
   CONNECT_FUNCTION("PrimaryButtonPressed", OnMousePressed);
   CONNECT_FUNCTION("Move", OnMouseMoved);
   CONNECT_FUNCTION("Release", OnMouseReleased);
 }
 
 const char **mitk::RegionGrowingTool::GetXPM() const
 {
   return mitkRegionGrowingTool_xpm;
 }
 
 us::ModuleResource mitk::RegionGrowingTool::GetIconResource() const
 {
   us::Module *module = us::GetModuleContext()->GetModule();
   us::ModuleResource resource = module->GetResource("RegionGrowing_48x48.png");
   return resource;
 }
 
 us::ModuleResource mitk::RegionGrowingTool::GetCursorIconResource() const
 {
   us::Module *module = us::GetModuleContext()->GetModule();
   us::ModuleResource resource = module->GetResource("RegionGrowing_Cursor_32x32.png");
   return resource;
 }
 
 const char *mitk::RegionGrowingTool::GetName() const
 {
   return "Region Growing";
 }
 
 void mitk::RegionGrowingTool::Activated()
 {
   Superclass::Activated();
 }
 
 void mitk::RegionGrowingTool::Deactivated()
 {
   Superclass::Deactivated();
 }
 
 // Get the average pixel value of square/cube with radius=neighborhood around index
 template <typename TPixel, unsigned int imageDimension>
-void mitk::RegionGrowingTool::GetNeighborhoodAverage(itk::Image<TPixel, imageDimension> *itkImage,
-                                                     itk::Index<imageDimension> index,
+void mitk::RegionGrowingTool::GetNeighborhoodAverage(const itk::Image<TPixel, imageDimension> *itkImage,
+                                                     const itk::Index<imageDimension>& index,
                                                      ScalarType *result,
                                                      unsigned int neighborhood)
 {
   // maybe assert that image dimension is only 2 or 3?
   auto neighborhoodInt = (int)neighborhood;
   TPixel averageValue(0);
   unsigned int numberOfPixels = (2 * neighborhood + 1) * (2 * neighborhood + 1);
   if (imageDimension == 3)
   {
     numberOfPixels *= (2 * neighborhood + 1);
   }
 
   MITK_DEBUG << "Getting neighborhood of " << numberOfPixels << " pixels around " << index;
 
   itk::Index<imageDimension> currentIndex;
 
   for (int i = (0 - neighborhoodInt); i <= neighborhoodInt; ++i)
   {
     currentIndex[0] = index[0] + i;
 
     for (int j = (0 - neighborhoodInt); j <= neighborhoodInt; ++j)
     {
       currentIndex[1] = index[1] + j;
 
       if (imageDimension == 3)
       {
         for (int k = (0 - neighborhoodInt); k <= neighborhoodInt; ++k)
         {
           currentIndex[2] = index[2] + k;
 
           if (itkImage->GetLargestPossibleRegion().IsInside(currentIndex))
           {
             averageValue += itkImage->GetPixel(currentIndex);
           }
           else
           {
             numberOfPixels -= 1;
           }
         }
       }
       else
       {
         if (itkImage->GetLargestPossibleRegion().IsInside(currentIndex))
         {
           averageValue += itkImage->GetPixel(currentIndex);
         }
         else
         {
           numberOfPixels -= 1;
         }
       }
     }
   }
 
   *result = (ScalarType)averageValue;
   *result /= numberOfPixels;
 }
 
 // Check whether index lies inside a segmentation
 template <typename TPixel, unsigned int imageDimension>
-void mitk::RegionGrowingTool::IsInsideSegmentation(itk::Image<TPixel, imageDimension> *itkImage,
-                                                   itk::Index<imageDimension> index,
+void mitk::RegionGrowingTool::IsInsideSegmentation(const itk::Image<TPixel, imageDimension> *itkImage,
+                                                   const itk::Index<imageDimension>& index,
                                                    bool *result)
 {
   if (itkImage->GetPixel(index) > 0)
   {
     *result = true;
   }
   else
   {
     *result = false;
   }
 }
 
 // Do the region growing (i.e. call an ITK filter that does it)
 template <typename TPixel, unsigned int imageDimension>
-void mitk::RegionGrowingTool::StartRegionGrowing(itk::Image<TPixel, imageDimension> *inputImage,
-                                                 itk::Index<imageDimension> seedIndex,
-                                                 std::array<ScalarType, 2> thresholds,
+void mitk::RegionGrowingTool::StartRegionGrowing(const itk::Image<TPixel, imageDimension> *inputImage,
+                                                 const itk::Index<imageDimension>& seedIndex,
+                                                 const std::array<ScalarType, 2>& thresholds,
                                                  mitk::Image::Pointer &outputImage)
 {
   MITK_DEBUG << "Starting region growing at index " << seedIndex << " with lower threshold " << thresholds[0]
              << " and upper threshold " << thresholds[1];
 
   typedef itk::Image<TPixel, imageDimension> InputImageType;
   typedef itk::Image<DefaultSegmentationDataType, imageDimension> OutputImageType;
 
   typedef itk::ConnectedThresholdImageFilter<InputImageType, OutputImageType> RegionGrowingFilterType;
   typename RegionGrowingFilterType::Pointer regionGrower = RegionGrowingFilterType::New();
 
   // perform region growing in desired segmented region
   regionGrower->SetInput(inputImage);
-  regionGrower->AddSeed(seedIndex);
+  regionGrower->SetSeed(seedIndex);
 
   regionGrower->SetLower(thresholds[0]);
   regionGrower->SetUpper(thresholds[1]);
 
   try
   {
     regionGrower->Update();
   }
   catch (...)
   {
     return; // Should we do something?
   }
 
   typename OutputImageType::Pointer resultImage = regionGrower->GetOutput();
 
   // Smooth result: Every pixel is replaced by the majority of the neighborhood
   typedef itk::NeighborhoodIterator<OutputImageType> NeighborhoodIteratorType;
   typedef itk::ImageRegionIterator<OutputImageType> ImageIteratorType;
 
   typename NeighborhoodIteratorType::RadiusType radius;
   radius.Fill(2); // for now, maybe make this something the user can adjust in the preferences?
 
   typedef itk::ImageDuplicator< OutputImageType > DuplicatorType;
   typename DuplicatorType::Pointer duplicator = DuplicatorType::New();
   duplicator->SetInputImage(resultImage);
   duplicator->Update();
 
   typename OutputImageType::Pointer resultDup = duplicator->GetOutput();
 
   NeighborhoodIteratorType neighborhoodIterator(radius, resultDup, resultDup->GetRequestedRegion());
   ImageIteratorType imageIterator(resultImage, resultImage->GetRequestedRegion());
 
   for (neighborhoodIterator.GoToBegin(), imageIterator.GoToBegin(); !neighborhoodIterator.IsAtEnd();
        ++neighborhoodIterator, ++imageIterator)
   {
     DefaultSegmentationDataType voteYes(0);
     DefaultSegmentationDataType voteNo(0);
 
     for (unsigned int i = 0; i < neighborhoodIterator.Size(); ++i)
     {
       if (neighborhoodIterator.GetPixel(i) > 0)
       {
         voteYes += 1;
       }
       else
       {
         voteNo += 1;
       }
     }
 
     if (voteYes > voteNo)
     {
       imageIterator.Set(1);
     }
     else
     {
       imageIterator.Set(0);
     }
   }
 
   if (resultImage.IsNull())
   {
     MITK_DEBUG << "Region growing result is empty.";
   }
 
   // Can potentially have multiple regions, use connected component image filter to label disjunct regions
   typedef itk::ConnectedComponentImageFilter<OutputImageType, OutputImageType> ConnectedComponentImageFilterType;
   typename ConnectedComponentImageFilterType::Pointer connectedComponentFilter =
     ConnectedComponentImageFilterType::New();
   connectedComponentFilter->SetInput(resultImage);
   connectedComponentFilter->Update();
   typename OutputImageType::Pointer resultImageCC = connectedComponentFilter->GetOutput();
   m_ConnectedComponentValue = resultImageCC->GetPixel(seedIndex);
 
   outputImage = mitk::GrabItkImageMemory(resultImageCC);
 }
 
+template <typename TPixel, unsigned int imageDimension>
+void mitk::RegionGrowingTool::CalculateInitialThresholds(const itk::Image<TPixel, imageDimension>*)
+{
+  LevelWindow levelWindow;
+  m_ToolManager->GetReferenceData(0)->GetLevelWindow(levelWindow);
+
+  m_ThresholdExtrema = { std::numeric_limits<TPixel>::lowest(), std::numeric_limits<TPixel>::max() };
+
+  const ScalarType lowerWindowBound = std::max(m_ThresholdExtrema[0], levelWindow.GetLowerWindowBound());
+  const ScalarType upperWindowBound = std::min(m_ThresholdExtrema[1], levelWindow.GetUpperWindowBound());
+
+  if (m_SeedValue < lowerWindowBound)
+  {
+    m_InitialThresholds = { m_ThresholdExtrema[0], lowerWindowBound };
+  }
+  else if (m_SeedValue > upperWindowBound)
+  {
+    m_InitialThresholds = { upperWindowBound, m_ThresholdExtrema[1] };
+  }
+  else
+  {
+    const ScalarType range = 0.1 * (upperWindowBound - lowerWindowBound); // 10% of the visible window
+
+    m_InitialThresholds[0] = std::min(std::max(lowerWindowBound, m_SeedValue - 0.5 * range), upperWindowBound - range);
+    m_InitialThresholds[1] = m_InitialThresholds[0] + range;
+  }
+}
+
 void mitk::RegionGrowingTool::OnMousePressed(StateMachineAction *, InteractionEvent *interactionEvent)
 {
   auto *positionEvent = dynamic_cast<mitk::InteractionPositionEvent *>(interactionEvent);
   if (!positionEvent)
     return;
 
   MITK_DEBUG << "OnMousePressed";
 
   m_LastEventSender = positionEvent->GetSender();
   m_LastEventSlice = m_LastEventSender->GetSlice();
   m_LastScreenPosition = positionEvent->GetPointerPositionOnScreen();
 
   // ReferenceSlice is from the underlying image, WorkingSlice from the active segmentation (can be empty)
   m_ReferenceSlice = FeedbackContourTool::GetAffectedReferenceSlice(positionEvent);
   m_WorkingSlice = FeedbackContourTool::GetAffectedWorkingSlice(positionEvent);
 
   if (m_WorkingSlice.IsNotNull()) // can't do anything without a working slice (i.e. a possibly empty segmentation)
   {
     MITK_DEBUG << "OnMousePressed: got working slice";
 
     // 2. Determine if the user clicked inside or outside of the segmentation/working slice (i.e. the whole volume)
     mitk::BaseGeometry::Pointer workingSliceGeometry;
     workingSliceGeometry = m_WorkingSlice->GetGeometry();
     workingSliceGeometry->WorldToIndex(positionEvent->GetPositionInWorld(), m_SeedPoint);
     itk::Index<2> indexInWorkingSlice2D;
     indexInWorkingSlice2D[0] = m_SeedPoint[0];
     indexInWorkingSlice2D[1] = m_SeedPoint[1];
 
     if (workingSliceGeometry->IsIndexInside(m_SeedPoint))
     {
       MITK_DEBUG << "OnMousePressed: point " << positionEvent->GetPositionInWorld() << " (index coordinates "
                  << m_SeedPoint << ") is inside working slice";
 
       // 3. determine the pixel value under the last click to determine what to do
       bool inside(true);
       AccessFixedDimensionByItk_2(m_WorkingSlice, IsInsideSegmentation, 2, indexInWorkingSlice2D, &inside);
       m_PaintingPixelValue = inside ? 0 : 1;
 
       if (inside)
       {
         MITK_DEBUG << "Clicked inside segmentation";
         // For now, we're doing nothing when the user clicks inside the segmentation. Behaviour can be implemented via
         // OnMousePressedInside()
         // When you do, be sure to remove the m_PaintingPixelValue check in OnMouseMoved() and OnMouseReleased()
         return;
       }
       else
       {
         MITK_DEBUG << "Clicked outside of segmentation";
         OnMousePressedOutside(nullptr, interactionEvent);
       }
     }
   }
 }
 
 // Use this to implement a behaviour for when the user clicks inside a segmentation (for example remove something)
 // Old IpPic code is kept as comment for reference
 void mitk::RegionGrowingTool::OnMousePressedInside()
 {
   //    mitk::InteractionPositionEvent* positionEvent = dynamic_cast<mitk::InteractionPositionEvent*>( interactionEvent
   //    );
   //    //const PositionEvent* positionEvent = dynamic_cast<const PositionEvent*>(stateEvent->GetEvent()); // checked in
   //    OnMousePressed
   //    // 3.1.1. Create a skeletonization of the segmentation and try to find a nice cut
   //    // apply the skeletonization-and-cut algorithm
   //    // generate contour to remove
   //    // set m_ReferenceSlice = nullptr so nothing will happen during mouse move
   //    // remember to fill the contour with 0 in mouserelease
   //    mitkIpPicDescriptor* segmentationHistory = ipMITKSegmentationCreateGrowerHistory( workingPicSlice,
   //    m_LastWorkingSeed, nullptr ); // free again
   //    if (segmentationHistory)
   //    {
   //        tCutResult cutContour = ipMITKSegmentationGetCutPoints( workingPicSlice, segmentationHistory,
   //        initialWorkingOffset ); // tCutResult is a ipSegmentation type
   //        mitkIpPicFree( segmentationHistory );
   //        if (cutContour.cutIt)
   //        {
   //            int timestep = positionEvent->GetSender()->GetTimeStep();
   //            // 3.1.2 copy point from float* to mitk::Contour
   //            ContourModel::Pointer contourInImageIndexCoordinates = ContourModel::New();
   //            contourInImageIndexCoordinates->Expand(timestep + 1);
   //            contourInImageIndexCoordinates->SetClosed(true, timestep);
   //            Point3D newPoint;
   //            for (int index = 0; index < cutContour.deleteSize; ++index)
   //            {
   //                newPoint[0] = cutContour.deleteCurve[ 2 * index + 0 ] - 0.5;//correction is needed because the
   //                output of the algorithm is center based
   //                newPoint[1] = cutContour.deleteCurve[ 2 * index + 1 ] - 0.5;//and we want our contour displayed
   //                corner based.
   //                newPoint[2] = 0.0;
 
   //                contourInImageIndexCoordinates->AddVertex( newPoint, timestep );
   //            }
 
   //            free(cutContour.traceline);
   //            free(cutContour.deleteCurve); // perhaps visualize this for fun?
   //            free(cutContour.onGradient);
 
   //            ContourModel::Pointer contourInWorldCoordinates = FeedbackContourTool::BackProjectContourFrom2DSlice(
   //            m_WorkingSlice->GetGeometry(), contourInImageIndexCoordinates, true ); // true: sub 0.5 for
   //            ipSegmentation correction
 
   //            FeedbackContourTool::SetFeedbackContour( contourInWorldCoordinates );
   //            FeedbackContourTool::SetFeedbackContourVisible(true);
   //            mitk::RenderingManager::GetInstance()->RequestUpdate( positionEvent->GetSender()->GetRenderWindow() );
   //            m_FillFeedbackContour = true;
   //        }
   //        else
   //        {
   //            m_FillFeedbackContour = false;
   //        }
 
   //    }
   //    else
   //    {
   //        m_FillFeedbackContour = false;
   //    }
 
   //    m_ReferenceSlice = nullptr;
 
   //    return true;
 }
 
 void mitk::RegionGrowingTool::OnMousePressedOutside(StateMachineAction *, InteractionEvent *interactionEvent)
 {
   auto *positionEvent = dynamic_cast<mitk::InteractionPositionEvent *>(interactionEvent);
 
   if (positionEvent)
   {
     // Get geometry and indices
     mitk::BaseGeometry::Pointer workingSliceGeometry;
     workingSliceGeometry = m_WorkingSlice->GetGeometry();
     itk::Index<2> indexInWorkingSlice2D;
     indexInWorkingSlice2D[0] = m_SeedPoint[0];
     indexInWorkingSlice2D[1] = m_SeedPoint[1];
 
     mitk::BaseGeometry::Pointer referenceSliceGeometry;
     referenceSliceGeometry =
       m_ReferenceSlice->GetGeometry();
     itk::Index<3> indexInReferenceSlice;
     itk::Index<2> indexInReferenceSlice2D;
     referenceSliceGeometry->WorldToIndex(positionEvent->GetPositionInWorld(), indexInReferenceSlice);
     indexInReferenceSlice2D[0] = indexInReferenceSlice[0];
     indexInReferenceSlice2D[1] = indexInReferenceSlice[1];
 
     // Get seed neighborhood
     ScalarType averageValue(0);
     AccessFixedDimensionByItk_3(m_ReferenceSlice, GetNeighborhoodAverage, 2, indexInReferenceSlice2D, &averageValue, 1);
     m_SeedValue = averageValue;
     MITK_DEBUG << "Seed value is " << m_SeedValue;
 
-    // Get level window settings
-    LevelWindow lw(0, 500); // default window 0 to 500, can we do something smarter here?
-    m_ToolManager->GetReferenceData(0)->GetLevelWindow(
-      lw); // will fill lw if levelwindow property is present, otherwise won't touch it.
-    ScalarType currentVisibleWindow = lw.GetWindow();
-    MITK_DEBUG << "Level window width is " << currentVisibleWindow;
-    m_InitialThresholds[0] = m_SeedValue - currentVisibleWindow / 20.0; // 20 is arbitrary (though works reasonably
-                                                                        // well), is there a better alternative (maybe
-                                                                        // option in preferences)?
-    m_InitialThresholds[1] = m_SeedValue + currentVisibleWindow / 20.0;
+    // Calculate initial thresholds
+    AccessFixedDimensionByItk(m_ReferenceSlice, CalculateInitialThresholds, 2);
     m_Thresholds[0] = m_InitialThresholds[0];
     m_Thresholds[1] = m_InitialThresholds[1];
 
     // Perform region growing
     mitk::Image::Pointer resultImage = mitk::Image::New();
     AccessFixedDimensionByItk_3(
       m_ReferenceSlice, StartRegionGrowing, 2, indexInWorkingSlice2D, m_Thresholds, resultImage);
     resultImage->SetGeometry(workingSliceGeometry);
 
     // Extract contour
     if (resultImage.IsNotNull() && m_ConnectedComponentValue >= 1)
     {
       float isoOffset = 0.33;
 
       mitk::ImageToContourModelFilter::Pointer contourExtractor = mitk::ImageToContourModelFilter::New();
       contourExtractor->SetInput(resultImage);
       contourExtractor->SetContourValue(m_ConnectedComponentValue - isoOffset);
       contourExtractor->Update();
       ContourModel::Pointer resultContour = ContourModel::New();
       resultContour = contourExtractor->GetOutput();
 
       // Show contour
       if (resultContour.IsNotNull())
       {
         ContourModel::Pointer resultContourWorld = FeedbackContourTool::BackProjectContourFrom2DSlice(
           workingSliceGeometry, FeedbackContourTool::ProjectContourTo2DSlice(m_WorkingSlice, resultContour));
 
         // this is not a beautiful solution, just one that works, check T22412 for details
         auto t = positionEvent->GetSender()->GetTimeStep();
 
         FeedbackContourTool::SetFeedbackContour(0 != t
           ? ContourModelUtils::MoveZerothContourTimeStep(resultContourWorld, t)
           : resultContourWorld);
 
         FeedbackContourTool::SetFeedbackContourVisible(true);
         mitk::RenderingManager::GetInstance()->RequestUpdate(m_LastEventSender->GetRenderWindow());
       }
     }
   }
 }
 
 void mitk::RegionGrowingTool::OnMouseMoved(StateMachineAction *, InteractionEvent *interactionEvent)
 {
   // Until OnMousePressedInside() implements a behaviour, we're just returning here whenever m_PaintingPixelValue is 0,
   // i.e. when the user clicked inside the segmentation
   if (m_PaintingPixelValue == 0)
   {
     return;
   }
 
   auto *positionEvent = dynamic_cast<mitk::InteractionPositionEvent *>(interactionEvent);
 
   if (m_ReferenceSlice.IsNotNull() && positionEvent)
   {
     // Get geometry and indices
     mitk::BaseGeometry::Pointer workingSliceGeometry;
     workingSliceGeometry = m_WorkingSlice->GetGeometry();
     itk::Index<2> indexInWorkingSlice2D;
     indexInWorkingSlice2D[0] = m_SeedPoint[0];
     indexInWorkingSlice2D[1] = m_SeedPoint[1];
 
     m_ScreenYDifference += positionEvent->GetPointerPositionOnScreen()[1] - m_LastScreenPosition[1];
     m_ScreenXDifference += positionEvent->GetPointerPositionOnScreen()[0] - m_LastScreenPosition[0];
     m_LastScreenPosition = positionEvent->GetPointerPositionOnScreen();
 
-    // Moving the mouse up and down adjusts the width of the threshold window, moving it left and right shifts the
-    // threshold window
-    m_Thresholds[0] = std::min<ScalarType>(
-      m_SeedValue, m_InitialThresholds[0] - (m_ScreenYDifference - m_ScreenXDifference) * m_MouseDistanceScaleFactor);
-    m_Thresholds[1] = std::max<ScalarType>(
-      m_SeedValue, m_InitialThresholds[1] + (m_ScreenYDifference + m_ScreenXDifference) * m_MouseDistanceScaleFactor);
-    MITK_DEBUG << "Screen difference X: " << m_ScreenXDifference;
+    // Moving the mouse up and down adjusts the width of the threshold window,
+    // moving it left and right shifts the threshold window
+    m_Thresholds[0] = std::min(m_SeedValue, m_InitialThresholds[0] - (m_ScreenYDifference - m_ScreenXDifference) * m_MouseDistanceScaleFactor);
+    m_Thresholds[1] = std::max(m_SeedValue, m_InitialThresholds[1] + (m_ScreenYDifference + m_ScreenXDifference) * m_MouseDistanceScaleFactor);
+
+    // Do not exceed the pixel type extrema of the reference slice, though
+    m_Thresholds[0] = std::max(m_ThresholdExtrema[0], m_Thresholds[0]);
+    m_Thresholds[1] = std::min(m_ThresholdExtrema[1], m_Thresholds[1]);
 
     // Perform region growing again and show the result
     mitk::Image::Pointer resultImage = mitk::Image::New();
     AccessFixedDimensionByItk_3(
       m_ReferenceSlice, StartRegionGrowing, 2, indexInWorkingSlice2D, m_Thresholds, resultImage);
     resultImage->SetGeometry(workingSliceGeometry);
 
     // Update the contour
     if (resultImage.IsNotNull() && m_ConnectedComponentValue >= 1)
     {
       float isoOffset = 0.33;
 
       mitk::ImageToContourModelFilter::Pointer contourExtractor = mitk::ImageToContourModelFilter::New();
       contourExtractor->SetInput(resultImage);
       contourExtractor->SetContourValue(m_ConnectedComponentValue - isoOffset);
       contourExtractor->Update();
       ContourModel::Pointer resultContour = ContourModel::New();
       resultContour = contourExtractor->GetOutput();
 
       // Show contour
       if (resultContour.IsNotNull())
       {
         ContourModel::Pointer resultContourWorld = FeedbackContourTool::BackProjectContourFrom2DSlice(
           workingSliceGeometry, FeedbackContourTool::ProjectContourTo2DSlice(m_WorkingSlice, resultContour));
 
         // this is not a beautiful solution, just one that works, check T22412 for details
         int timestep = positionEvent->GetSender()->GetTimeStep();
         if (0 != timestep)
         {
           int size = resultContourWorld->GetNumberOfVertices(0);
           auto resultContourTimeWorld = mitk::ContourModel::New();
           resultContourTimeWorld->Expand(timestep + 1);
           for (int loop = 0; loop < size; ++loop)
           {
             resultContourTimeWorld->AddVertex(resultContourWorld->GetVertexAt(loop, 0), timestep);
           }
           FeedbackContourTool::SetFeedbackContour(resultContourTimeWorld);
         }
         else
         {
           FeedbackContourTool::SetFeedbackContour(resultContourWorld);
         }
 
         FeedbackContourTool::SetFeedbackContourVisible(true);
         mitk::RenderingManager::GetInstance()->ForceImmediateUpdate(positionEvent->GetSender()->GetRenderWindow());
       }
     }
   }
 }
 
 void mitk::RegionGrowingTool::OnMouseReleased(StateMachineAction *, InteractionEvent *interactionEvent)
 {
   // Until OnMousePressedInside() implements a behaviour, we're just returning here whenever m_PaintingPixelValue is 0,
   // i.e. when the user clicked inside the segmentation
   if (m_PaintingPixelValue == 0)
   {
     return;
   }
 
   auto *positionEvent = dynamic_cast<mitk::InteractionPositionEvent *>(interactionEvent);
 
   if (m_WorkingSlice.IsNotNull() && m_FillFeedbackContour && positionEvent)
   {
     // Project contour into working slice
     ContourModel *feedbackContour(FeedbackContourTool::GetFeedbackContour());
 
     ContourModel::Pointer projectedContour;
 
     // this is not a beautiful solution, just one that works, check T22412 for details
     int timestep = positionEvent->GetSender()->GetTimeStep();
     if (0 != timestep)
     {
       int size = feedbackContour->GetNumberOfVertices(timestep);
       auto feedbackContourTime = mitk::ContourModel::New();
       feedbackContourTime->Expand(timestep + 1);
       for (int loop = 0; loop < size; ++loop)
       {
         feedbackContourTime->AddVertex(feedbackContour->GetVertexAt(loop, timestep), 0);
       }
 
       projectedContour =
         FeedbackContourTool::ProjectContourTo2DSlice(m_WorkingSlice, feedbackContourTime, false, false);
     }
     else
     {
       projectedContour =
         FeedbackContourTool::ProjectContourTo2DSlice(m_WorkingSlice, feedbackContour, false, false);
     }
 
     // If there is a projected contour, fill it
     if (projectedContour.IsNotNull())
     {
       // Get working data to pass to following method so we don't overwrite locked labels in a LabelSetImage
       mitk::DataNode *workingNode(m_ToolManager->GetWorkingData(0));
       mitk::LabelSetImage *labelImage = workingNode != nullptr
         ? dynamic_cast<mitk::LabelSetImage*>(workingNode->GetData())
         : nullptr;
 
       MITK_DEBUG << "Filling Segmentation";
 
       if (labelImage != nullptr)
       {
         // m_PaintingPixelValue only decides whether to paint or not
         // For LabelSetImages we want to paint with the active label value
         auto activeLabel = labelImage->GetActiveLabel(labelImage->GetActiveLayer())->GetValue();
         mitk::ContourModelUtils::FillContourInSlice(projectedContour,
                                                     0,
                                                     m_WorkingSlice,
                                                     labelImage,
                                                     m_PaintingPixelValue * activeLabel);
       }
       else
       {
         mitk::ContourModelUtils::FillContourInSlice(projectedContour,
                                                     0,
                                                     m_WorkingSlice,
                                                     m_WorkingSlice,
                                                     m_PaintingPixelValue);
       }
       this->WriteBackSegmentationResult(positionEvent, m_WorkingSlice);
       FeedbackContourTool::SetFeedbackContourVisible(false);
     }
 
     m_ScreenYDifference = 0;
     m_ScreenXDifference = 0;
   }
 }
diff --git a/Modules/Segmentation/Interactions/mitkRegionGrowingTool.h b/Modules/Segmentation/Interactions/mitkRegionGrowingTool.h
index 86f8f074f6..6ba9999718 100644
--- a/Modules/Segmentation/Interactions/mitkRegionGrowingTool.h
+++ b/Modules/Segmentation/Interactions/mitkRegionGrowingTool.h
@@ -1,153 +1,160 @@
 /*============================================================================
 
 The Medical Imaging Interaction Toolkit (MITK)
 
 Copyright (c) German Cancer Research Center (DKFZ)
 All rights reserved.
 
 Use of this source code is governed by a 3-clause BSD license that can be
 found in the LICENSE file.
 
 ============================================================================*/
 
 #ifndef mitkRegionGrowingTool_h_Included
 #define mitkRegionGrowingTool_h_Included
 
 #include "mitkFeedbackContourTool.h"
 #include <MitkSegmentationExports.h>
 #include <array>
 
 namespace us
 {
   class ModuleResource;
 }
 
 namespace mitk
 {
   /**
     \brief A slice based region growing tool.
 
     \sa FeedbackContourTool
 
     \ingroup Interaction
     \ingroup ToolManagerEtAl
 
     When the user presses the mouse button, RegionGrowingTool will use the gray values at that position
     to initialize a region growing algorithm (in the affected 2D slice).
 
     By moving the mouse up and down while the button is still pressed, the user can widen or narrow the threshold
     window, i.e. select more or less within the desired region.
     The current result of region growing will always be shown as a contour to the user.
 
     After releasing the button, the current result of the region growing algorithm will be written to the
     working image of this tool's ToolManager.
 
     If the first click is <i>inside</i> a segmentation, nothing will happen (other behaviour, for example removal of a
     region, can be implemented via OnMousePressedInside()).
 
     \warning Only to be instantiated by mitk::ToolManager.
 
     $Author$
   */
   class MITKSEGMENTATION_EXPORT RegionGrowingTool : public FeedbackContourTool
   {
   public:
     mitkClassMacro(RegionGrowingTool, FeedbackContourTool);
     itkFactorylessNewMacro(Self);
     itkCloneMacro(Self);
 
       const char **GetXPM() const override;
     us::ModuleResource GetCursorIconResource() const override;
     us::ModuleResource GetIconResource() const override;
 
     const char *GetName() const override;
 
   protected:
     RegionGrowingTool(); // purposely hidden
     ~RegionGrowingTool() override;
 
     void ConnectActionsAndFunctions() override;
 
     void Activated() override;
     void Deactivated() override;
 
     /**
      * @brief OnMousePressed is called when the user clicks.
      * Calls either OnMousePressedInside() or OnMousePressedOutside().
      */
     virtual void OnMousePressed(StateMachineAction *, InteractionEvent *interactionEvent);
 
     /**
      * @brief OnMousePressedInside can be used to implement behaviour when the user clicks inside a segmentation.
      */
     virtual void OnMousePressedInside();
 
     /**
      * @brief OnMousePressedOutside is called when the user clicks outside of the segmented area.
      * Grows a region.
      */
     virtual void OnMousePressedOutside(StateMachineAction *, InteractionEvent *interactionEvent);
 
     /**
      * @brief OnMouseMoved is called when the user moves the mouse with the left mouse button pressed.
      * Adjusts the thresholds.
      * Up: Increase upper threshold, decrease lower threshold.
      * Down: Decrease upper threshold, increase lower threshold.
      * Right: Increase both thresholds.
      * Left: Decrease both thresholds.
      */
     virtual void OnMouseMoved(StateMachineAction *, InteractionEvent *interactionEvent);
 
     /**
      * @brief OnMouseReleased converts the feedback contour to a segmentation.
      */
     virtual void OnMouseReleased(StateMachineAction *, InteractionEvent *interactionEvent);
 
     /**
      * @brief Template to calculate average pixel value around index using a square/cube with radius neighborhood.
      * Example: 1 = 3x3 pixels, 2 = 5x5 pixels, etc.
      */
     template <typename TPixel, unsigned int imageDimension>
-    void GetNeighborhoodAverage(itk::Image<TPixel, imageDimension> *itkImage,
-                                itk::Index<imageDimension> index,
+    void GetNeighborhoodAverage(const itk::Image<TPixel, imageDimension> *itkImage,
+                                const itk::Index<imageDimension>& index,
                                 ScalarType *result,
                                 unsigned int neighborhood = 1);
 
     /**
      * @brief Template to check whether index is inside already segmented area.
      */
     template <typename TPixel, unsigned int imageDimension>
-    void IsInsideSegmentation(itk::Image<TPixel, imageDimension> *itkImage,
-                              itk::Index<imageDimension> index,
+    void IsInsideSegmentation(const itk::Image<TPixel, imageDimension> *itkImage,
+                              const itk::Index<imageDimension>& index,
                               bool *result);
 
     /**
      * @brief Template that calls an ITK filter to do the region growing.
      */
     template <typename TPixel, unsigned int imageDimension>
-    void StartRegionGrowing(itk::Image<TPixel, imageDimension> *itkImage,
-                            itk::Index<imageDimension> seedPoint,
-                            std::array<ScalarType, 2> thresholds,
+    void StartRegionGrowing(const itk::Image<TPixel, imageDimension> *itkImage,
+                            const itk::Index<imageDimension>& seedPoint,
+                            const std::array<ScalarType, 2>& thresholds,
                             mitk::Image::Pointer &outputImage);
 
+    /**
+     * @brief Template to calculate the initial thresholds for region growing.
+     */
+    template <typename TPixel, unsigned int imageDimension>
+    void CalculateInitialThresholds(const itk::Image<TPixel, imageDimension>* itkImage);
+
     Image::Pointer m_ReferenceSlice;
     Image::Pointer m_WorkingSlice;
 
     ScalarType m_SeedValue;
     itk::Index<3> m_SeedPoint;
+    std::array<ScalarType, 2> m_ThresholdExtrema;
     std::array<ScalarType, 2> m_Thresholds;
     std::array<ScalarType, 2> m_InitialThresholds;
     Point2I m_LastScreenPosition;
     int m_ScreenYDifference;
     int m_ScreenXDifference;
 
   private:
     ScalarType m_MouseDistanceScaleFactor;
     int m_PaintingPixelValue;
     bool m_FillFeedbackContour;
     int m_ConnectedComponentValue;
   };
 
 } // namespace
 
 #endif