Page MenuHomePhabricator

interpolationsegmentation.diff

Authored By
engelm
Mar 14 2013, 2:33 PM
Size
233 KB
Referenced Files
None
Subscribers
None

interpolationsegmentation.diff

diff --git a/Core/Code/Algorithms/mitkAutoCropImageFilter.cpp b/Core/Code/Algorithms/mitkAutoCropImageFilter.cpp
new file mode 100644
index 0000000..3fdb301
--- /dev/null
+++ b/Core/Code/Algorithms/mitkAutoCropImageFilter.cpp
@@ -0,0 +1,371 @@
+/*===================================================================
+
+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.
+
+===================================================================*/
+
+#include "mitkAutoCropImageFilter.h"
+
+#include "mitkImageCast.h"
+#include "mitkImageAccessByItk.h"
+#include "mitkGeometry3D.h"
+#include "mitkStatusBar.h"
+
+#include "mitkPlaneGeometry.h"
+
+#include <itkImageRegionConstIterator.h>
+#include <itkRegionOfInterestImageFilter.h>
+
+
+mitk::AutoCropImageFilter::AutoCropImageFilter()
+: m_BackgroundValue(0),
+ m_MarginFactor(1.0),
+ m_TimeSelector(NULL),
+ m_OverrideCroppingRegion(false)
+{
+
+}
+
+
+mitk::AutoCropImageFilter::~AutoCropImageFilter()
+{
+
+}
+
+
+template < typename TPixel, unsigned int VImageDimension>
+void mitk::AutoCropImageFilter::ITKCrop3DImage( itk::Image< TPixel, VImageDimension >* inputItkImage, unsigned int timestep)
+{
+
+ if (inputItkImage == NULL)
+ {
+ mitk::StatusBar::GetInstance()->DisplayErrorText ("An internal error occurred. Can't convert Image. Please report to bugs@mitk.org");
+ MITK_ERROR << "image is NULL...returning" << std::endl;
+ return;
+ }
+
+ typedef itk::Image< TPixel, VImageDimension > InternalImageType;
+ typedef typename InternalImageType::Pointer InternalImagePointer;
+
+ typedef itk::RegionOfInterestImageFilter < InternalImageType, InternalImageType > ROIFilterType;
+ typedef typename itk::RegionOfInterestImageFilter < InternalImageType, InternalImageType >::Pointer ROIFilterPointer;
+
+ InternalImagePointer outputItk = InternalImageType::New();
+
+ ROIFilterPointer roiFilter = ROIFilterType::New();
+ roiFilter->SetInput(0,inputItkImage);
+ roiFilter->SetRegionOfInterest(this->GetCroppingRegion());
+ roiFilter->Update();
+ outputItk = roiFilter->GetOutput();
+ outputItk->DisconnectPipeline();
+
+ mitk::Image::Pointer newMitkImage = mitk::Image::New();
+ mitk::CastToMitkImage( outputItk, newMitkImage );
+ MITK_INFO << "Crop-Output dimension: " << (newMitkImage->GetDimension() == 3) << " Filter-Output dimension: "<<this->GetOutput()->GetDimension()<< " Timestep: " << timestep;
+
+// const mitk::ChannelDescriptor desc = newMitkImage->GetChannelDescriptor(0);
+// unsigned char* image3D = desc.GetData();
+// this->GetOutput()->SetVolume( (void*) &image3D , timestep );
+
+ this->GetOutput()->SetVolume( newMitkImage->GetData(), timestep);
+// this->SetOutput(newMitkImage);
+}
+
+void mitk::AutoCropImageFilter::GenerateOutputInformation()
+{
+ mitk::Image::Pointer input = const_cast<mitk::Image*> (this->GetInput());
+ mitk::Image::Pointer output = this->GetOutput();
+
+ if(input->GetDimension() <= 2)
+ {
+ MITK_ERROR << "Only 3D any 4D images are supported." << std::endl;
+ return;
+ }
+
+ ComputeNewImageBounds();
+
+ if ((output->IsInitialized()) && (output->GetPipelineMTime() <= m_TimeOfHeaderInitialization.GetMTime()))
+ return;
+
+ itkDebugMacro(<<"GenerateOutputInformation()");
+
+ // PART I: initialize input requested region. We do this already here (and not
+ // later when GenerateInputRequestedRegion() is called), because we
+ // also need the information to setup the output.
+
+ // pre-initialize input-requested-region to largest-possible-region
+ // and correct time-region; spatial part will be cropped by
+ // bounding-box of bounding-object below
+ m_InputRequestedRegion = input->GetLargestPossibleRegion();
+
+ // build region out of index and size calculated in ComputeNewImageBounds()
+
+ mitk::SlicedData::IndexType index;
+ index[0] = m_RegionIndex[0];
+ index[1] = m_RegionIndex[1];
+ index[2] = m_RegionIndex[2];
+ index[3] = m_InputRequestedRegion.GetIndex()[3];
+ index[4] = m_InputRequestedRegion.GetIndex()[4];
+
+ mitk::SlicedData::SizeType size;
+ size[0] = m_RegionSize[0];
+ size[1] = m_RegionSize[1];
+ size[2] = m_RegionSize[2];
+ size[3] = m_InputRequestedRegion.GetSize()[3];
+ size[4] = m_InputRequestedRegion.GetSize()[4];
+
+ mitk::SlicedData::RegionType cropRegion(index, size);
+
+ // crop input-requested-region with cropping region computed from the image data
+ if(m_InputRequestedRegion.Crop(cropRegion)==false)
+ {
+ // crop not possible => do nothing: set time size to 0.
+ size.Fill(0);
+ m_InputRequestedRegion.SetSize(size);
+ return;
+ }
+
+ // set input-requested-region, because we access it later in
+ // GenerateInputRequestedRegion (there we just set the time)
+ input->SetRequestedRegion(&m_InputRequestedRegion);
+
+
+ // PART II: initialize output image
+ unsigned int dimension = input->GetDimension();
+ unsigned int *dimensions = new unsigned int [dimension];
+ itk2vtk(m_InputRequestedRegion.GetSize(), dimensions);
+ if(dimension>3)
+ memcpy(dimensions+3, input->GetDimensions()+3, (dimension-3)*sizeof(unsigned int));
+
+ // create basic slicedGeometry that will be initialized below
+ output->Initialize(mitk::PixelType( GetOutputPixelType() ), dimension, dimensions);
+ delete [] dimensions;
+
+ //clone the IndexToWorldTransform from the input, otherwise we will overwrite it, when adjusting the origin of the output image!!
+ itk::ScalableAffineTransform< mitk::ScalarType,3 >::Pointer cloneTransform = itk::ScalableAffineTransform< mitk::ScalarType,3 >::New();
+ cloneTransform->Compose(input->GetGeometry()->GetIndexToWorldTransform());
+ output->GetGeometry()->SetIndexToWorldTransform( cloneTransform.GetPointer() );
+
+ // Position the output Image to match the corresponding region of the input image
+ mitk::SlicedGeometry3D* slicedGeometry = output->GetSlicedGeometry();
+ mitk::SlicedGeometry3D::Pointer inputGeometry = input->GetSlicedGeometry();
+ const mitk::SlicedData::IndexType& start = m_InputRequestedRegion.GetIndex();
+ mitk::Point3D origin; vtk2itk(start, origin);
+ input->GetSlicedGeometry()->IndexToWorld(origin, origin);
+ slicedGeometry->SetOrigin(origin);
+
+ // get the PlaneGeometry for the first slice of the original image
+ mitk::PlaneGeometry::Pointer plane = dynamic_cast<mitk::PlaneGeometry*>( inputGeometry->GetGeometry2D( 0 )->Clone().GetPointer() );
+ assert( plane );
+
+ // re-initialize the plane according to the new requirements:
+ // dimensions of the cropped image
+ // right- and down-vector as well as spacing do not change, so use the ones from
+ // input image
+ ScalarType dimX = output->GetDimensions()[0];
+ ScalarType dimY = output->GetDimensions()[1];
+ mitk::Vector3D right = plane->GetAxisVector(0);
+ mitk::Vector3D down = plane->GetAxisVector(1);
+ mitk::Vector3D spacing = plane->GetSpacing();
+ plane->InitializeStandardPlane( dimX, dimY, right, down, &spacing );
+ // set the new origin on the PlaneGeometry as well
+ plane->SetOrigin(origin);
+
+ // re-initialize the slicedGeometry with the correct planeGeometry
+ // in order to get a fully initialized SlicedGeometry3D
+ slicedGeometry->InitializeEvenlySpaced( plane, inputGeometry->GetSpacing()[2], output->GetSlicedGeometry()->GetSlices() );
+
+
+ mitk::TimeSlicedGeometry* timeSlicedGeometry = output->GetTimeSlicedGeometry();
+ timeSlicedGeometry->InitializeEvenlyTimed(slicedGeometry, output->GetDimension(3));
+ timeSlicedGeometry->CopyTimes(input->GetTimeSlicedGeometry());
+
+ m_TimeOfHeaderInitialization.Modified();
+
+ output->SetPropertyList(input->GetPropertyList()->Clone());
+}
+
+void mitk::AutoCropImageFilter::GenerateData()
+{
+ mitk::Image::ConstPointer input = this->GetInput();
+ mitk::Image::Pointer output = this->GetOutput();
+
+ if(input.IsNull())
+ return;
+
+ if(input->GetDimension() <= 2)
+ {
+ MITK_ERROR << "Only 3D and 4D images supported";
+ return;
+ }
+
+ if((output->IsInitialized()==false) )
+ return;
+
+ if( m_TimeSelector.IsNull() ) m_TimeSelector = mitk::ImageTimeSelector::New();
+
+ m_TimeSelector->SetInput(input);
+
+ mitk::SlicedData::RegionType outputRegion = input->GetRequestedRegion();
+
+ int tstart = outputRegion.GetIndex(3);
+ int tmax = tstart + outputRegion.GetSize(3);
+
+ for( int timestep=tstart;timestep<tmax;++timestep )
+ {
+ m_TimeSelector->SetTimeNr(timestep);
+ m_TimeSelector->UpdateLargestPossibleRegion();
+
+ AccessFixedDimensionByItk_1( m_TimeSelector->GetOutput(), ITKCrop3DImage, 3, timestep );
+ }
+
+ // this->GetOutput()->Update(); // Not sure if this is necessary...
+
+ m_TimeOfHeaderInitialization.Modified();
+}
+
+void mitk::AutoCropImageFilter::ComputeNewImageBounds()
+{
+ mitk::Image::ConstPointer inputMitk = this->GetInput();
+
+ if (m_OverrideCroppingRegion)
+ {
+ for (unsigned int i=0; i<3; ++i)
+ {
+ m_RegionIndex[i] = m_CroppingRegion.GetIndex()[i];
+ m_RegionSize[i] = m_CroppingRegion.GetSize()[i];
+
+ if (m_RegionIndex[i] >= inputMitk->GetDimension(i))
+ {
+ itkExceptionMacro("Cropping index is not inside the image. "
+ << std::endl << "Index:"
+ << std::endl << m_CroppingRegion.GetIndex()
+ << std::endl << "Size:"
+ << std::endl << m_CroppingRegion.GetSize());
+ }
+
+ if (m_RegionIndex[i] + m_RegionSize[i] >= inputMitk->GetDimension(i))
+ {
+ m_RegionSize[i] = inputMitk->GetDimension(i) - m_RegionIndex[i];
+ }
+ }
+
+ for (unsigned int i=0; i<3; ++i)
+ {
+ m_RegionIndex[i] = m_CroppingRegion.GetIndex()[i];
+ m_RegionSize[i] = m_CroppingRegion.GetSize()[i];
+
+ }
+ }
+ else
+ {
+ // Check if a 3D or 4D image is present
+ unsigned int timeSteps = 1;
+ if (inputMitk->GetDimension() == 4 )
+ timeSteps = inputMitk->GetDimension(3);
+
+ ImageType::IndexType minima,maxima;
+
+ if (inputMitk->GetDimension() == 4)
+ {
+ // initialize with time step 0
+ m_TimeSelector = mitk::ImageTimeSelector::New();
+ m_TimeSelector->SetInput( inputMitk );
+ m_TimeSelector->SetTimeNr( 0 );
+ m_TimeSelector->UpdateLargestPossibleRegion();
+ inputMitk = m_TimeSelector->GetOutput();
+ }
+
+ ImagePointer inputItk = ImageType::New();
+ mitk::CastToItkImage( inputMitk , inputItk );
+
+ // it is assumed that all volumes in a time series have the same 3D dimensions
+ ImageType::RegionType origRegion = inputItk->GetLargestPossibleRegion();
+
+ // Initialize min and max on the first (or only) time step
+ maxima = inputItk->GetLargestPossibleRegion().GetIndex();
+ minima[0] = inputItk->GetLargestPossibleRegion().GetSize()[0];
+ minima[1] = inputItk->GetLargestPossibleRegion().GetSize()[1];
+ minima[2] = inputItk->GetLargestPossibleRegion().GetSize()[2];
+
+ typedef itk::ImageRegionConstIterator< ImageType > ConstIteratorType;
+
+ for(unsigned int idx = 0; idx < timeSteps; ++idx)
+ {
+ // if 4D image, update time step and itk image
+ if( idx > 0)
+ {
+ m_TimeSelector->SetTimeNr( idx );
+ m_TimeSelector->UpdateLargestPossibleRegion();
+ inputMitk = m_TimeSelector->GetOutput();
+ mitk::CastToItkImage( inputMitk , inputItk );
+ }
+
+ ConstIteratorType inIt( inputItk, origRegion );
+
+ for ( inIt.GoToBegin(); !inIt.IsAtEnd(); ++inIt)
+ {
+ float pix_val = inIt.Get();
+ if ( fabs(pix_val - m_BackgroundValue) > mitk::eps )
+ {
+ for (int i=0; i < 3; i++)
+ {
+ minima[i] = vnl_math_min((int)minima[i],(int)(inIt.GetIndex()[i]));
+ maxima[i] = vnl_math_max((int)maxima[i],(int)(inIt.GetIndex()[i]));
+ }
+ }
+ }
+ }
+
+ typedef ImageType::RegionType::SizeType::SizeValueType SizeValueType;
+
+ m_RegionSize[0] = (SizeValueType)(m_MarginFactor * (maxima[0] - minima[0] + 1 ));
+ m_RegionSize[1] = (SizeValueType)(m_MarginFactor * (maxima[1] - minima[1] + 1 ));
+ m_RegionSize[2] = (SizeValueType)(m_MarginFactor * (maxima[2] - minima[2] + 1 ));
+ m_RegionIndex = minima;
+
+ m_RegionIndex[0] -= (m_RegionSize[0] - maxima[0] + minima[0] - 1 )/2;
+ m_RegionIndex[1] -= (m_RegionSize[1] - maxima[1] + minima[1] - 1 )/2;
+ m_RegionIndex[2] -= (m_RegionSize[2] - maxima[2] + minima[2] - 1 )/2;
+
+ ImageType::RegionType cropRegion(m_RegionIndex,m_RegionSize);
+ origRegion.Crop(cropRegion);
+
+ m_RegionSize[0] = origRegion.GetSize()[0];
+ m_RegionSize[1] = origRegion.GetSize()[1];
+ m_RegionSize[2] = origRegion.GetSize()[2];
+
+ m_RegionIndex[0] = origRegion.GetIndex()[0];
+ m_RegionIndex[1] = origRegion.GetIndex()[1];
+ m_RegionIndex[2] = origRegion.GetIndex()[2];
+
+ m_CroppingRegion = origRegion;
+ }
+}
+
+
+void mitk::AutoCropImageFilter::GenerateInputRequestedRegion()
+{
+
+}
+
+const mitk::PixelType mitk::AutoCropImageFilter::GetOutputPixelType()
+{
+ return this->GetInput()->GetPixelType();
+}
+
+void mitk::AutoCropImageFilter::SetCroppingRegion(RegionType overrideRegion)
+{
+ m_CroppingRegion = overrideRegion;
+ m_OverrideCroppingRegion = true;
+}
diff --git a/Core/Code/Algorithms/mitkAutoCropImageFilter.h b/Core/Code/Algorithms/mitkAutoCropImageFilter.h
new file mode 100644
index 0000000..269645d
--- /dev/null
+++ b/Core/Code/Algorithms/mitkAutoCropImageFilter.h
@@ -0,0 +1,130 @@
+/*===================================================================
+
+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 __mitkAutoCropImageFilter_h_
+#define __mitkAutoCropImageFilter_h_
+
+#include "mitkCommon.h"
+#include "MitkExports.h"
+#include "mitkSubImageSelector.h"
+#include "mitkImageTimeSelector.h"
+
+#include <itkImageRegion.h>
+#include <itkImage.h>
+
+namespace mitk {
+
+/**
+ *
+ * @brief Shrink the image borders to a minimum considering a background color.
+ *
+ * This filter determines the smallest bounding box of all pixels different
+ * from the background, and returns an output image which has been cropped to this size.
+ * The box calculated this way is not the smallest possible box, but the box with the
+ * smallest sides perpendicular to the world coordinate system.
+ *
+ * The filter works on 3D and 4D image data. For the 4D case, the smallest box is
+ * calculated with side lengths as the maximum of single side lengths from all time steps.
+ *
+ * 2D images are not supported, and will never be.
+ *
+ * It is also possible to set the region to be cropped manually using the
+ * SetCroppingRegion() method.
+ *
+ * A margin can be set to enlarge the cropped region with a constant factor in all
+ * directions around the smallest possible.
+ *
+ *
+ * @ingroup Process
+ *
+ * @author Thomas Boettger; revised by Tobias Schwarz and Daniel Stein; Division of Medical
+ * and Biological Informatics
+ *
+ */
+
+class MITK_CORE_EXPORT AutoCropImageFilter : public SubImageSelector
+{
+public:
+
+ typedef itk::ImageRegion<3> RegionType;
+
+ mitkClassMacro(AutoCropImageFilter, SubImageSelector);
+
+ itkNewMacro(Self);
+
+ itkGetConstMacro(BackgroundValue,float);
+ itkSetMacro(BackgroundValue,float);
+
+ itkGetConstMacro(MarginFactor,float);
+ itkSetMacro(MarginFactor,float);
+
+ itkGetMacro(CroppingRegion, RegionType);
+
+ // Use this method to manually set a region
+ void SetCroppingRegion(RegionType overrideRegion);
+
+ virtual const PixelType GetOutputPixelType();
+
+protected:
+
+ // default constructor
+ AutoCropImageFilter();
+
+ // default destructor
+ virtual ~AutoCropImageFilter();
+
+ // This method calculates the actual smallest box
+ void ComputeNewImageBounds();
+
+ // Crops the image using the itk::RegionOfInterestImageFilter and creates the new output image
+ template < typename TPixel, unsigned int VImageDimension>
+ void ITKCrop3DImage( itk::Image< TPixel, VImageDimension >* inputItkImage, unsigned int timestep );
+
+ // Here, the output image is initialized by the input and the newly calculated region
+ virtual void GenerateOutputInformation();
+
+ // Purposely not implemented
+ virtual void GenerateInputRequestedRegion();
+
+ // Crops the image on all time steps
+ virtual void GenerateData();
+
+ float m_BackgroundValue;
+
+ RegionType m_CroppingRegion;
+
+ float m_MarginFactor;
+
+ typedef itk::Image<float,3> ImageType;
+ typedef ImageType::Pointer ImagePointer;
+
+ RegionType::SizeType m_RegionSize;
+ RegionType::IndexType m_RegionIndex;
+
+ mitk::ImageTimeSelector::Pointer m_TimeSelector;
+
+ mitk::SlicedData::RegionType m_InputRequestedRegion;
+ itk::TimeStamp m_TimeOfHeaderInitialization;
+
+ bool m_OverrideCroppingRegion;
+
+};
+
+} // namespace mitk
+
+#endif
+
+
diff --git a/Core/Code/Algorithms/mitkImageToSurfaceFilter.cpp b/Core/Code/Algorithms/mitkImageToSurfaceFilter.cpp
new file mode 100644
index 0000000..6b2633a
--- /dev/null
+++ b/Core/Code/Algorithms/mitkImageToSurfaceFilter.cpp
@@ -0,0 +1,234 @@
+/*===================================================================
+
+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.
+
+===================================================================*/
+
+#include <mitkImageToSurfaceFilter.h>
+#include "mitkException.h"
+#include <vtkImageData.h>
+#include <vtkDecimatePro.h>
+#include <vtkImageChangeInformation.h>
+#include <vtkLinearTransform.h>
+#include <vtkMath.h>
+#include <vtkMatrix4x4.h>
+#include <vtkQuadricDecimation.h>
+
+#include <vtkSmartPointer.h>
+#include <vtkPolyDataNormals.h>
+#include <vtkCleanPolyData.h>
+
+#include "mitkProgressBar.h"
+
+mitk::ImageToSurfaceFilter::ImageToSurfaceFilter():
+ m_Smooth(false),
+ m_Decimate( NoDecimation),
+ m_Threshold(1.0),
+ m_TargetReduction(0.95f),
+ m_SmoothIteration(50),
+ m_SmoothRelaxation(0.1)
+{
+}
+
+mitk::ImageToSurfaceFilter::~ImageToSurfaceFilter()
+{
+}
+
+void mitk::ImageToSurfaceFilter::CreateSurface(int time, vtkImageData *vtkimage, mitk::Surface * surface, const ScalarType threshold)
+{
+ vtkImageChangeInformation *indexCoordinatesImageFilter = vtkImageChangeInformation::New();
+ indexCoordinatesImageFilter->SetInput(vtkimage);
+ indexCoordinatesImageFilter->SetOutputOrigin(0.0,0.0,0.0);
+
+ //MarchingCube -->create Surface
+ vtkMarchingCubes *skinExtractor = vtkMarchingCubes::New();
+ skinExtractor->ComputeScalarsOff();
+ skinExtractor->SetInput(indexCoordinatesImageFilter->GetOutput());//RC++
+ indexCoordinatesImageFilter->Delete();
+ skinExtractor->SetValue(0, threshold);
+
+ vtkPolyData *polydata;
+ polydata = skinExtractor->GetOutput();
+ polydata->Register(NULL);//RC++
+ skinExtractor->Delete();
+
+ if (m_Smooth)
+ {
+ vtkSmoothPolyDataFilter *smoother = vtkSmoothPolyDataFilter::New();
+ //read poly1 (poly1 can be the original polygon, or the decimated polygon)
+ smoother->SetInput(polydata);//RC++
+ smoother->SetNumberOfIterations( m_SmoothIteration );
+ smoother->SetRelaxationFactor( m_SmoothRelaxation );
+ smoother->SetFeatureAngle( 60 );
+ smoother->FeatureEdgeSmoothingOff();
+ smoother->BoundarySmoothingOff();
+ smoother->SetConvergence( 0 );
+
+ polydata->Delete();//RC--
+ polydata = smoother->GetOutput();
+ polydata->Register(NULL);//RC++
+ smoother->Delete();
+ }
+ ProgressBar::GetInstance()->Progress();
+
+ //decimate = to reduce number of polygons
+ if(m_Decimate==DecimatePro)
+ {
+ vtkDecimatePro *decimate = vtkDecimatePro::New();
+ decimate->SplittingOff();
+ decimate->SetErrorIsAbsolute(5);
+ decimate->SetFeatureAngle(30);
+ decimate->PreserveTopologyOn();
+ decimate->BoundaryVertexDeletionOff();
+ decimate->SetDegree(10); //std-value is 25!
+
+ decimate->SetInput(polydata);//RC++
+ decimate->SetTargetReduction(m_TargetReduction);
+ decimate->SetMaximumError(0.002);
+
+ polydata->Delete();//RC--
+ polydata = decimate->GetOutput();
+ polydata->Register(NULL);//RC++
+ decimate->Delete();
+ }
+ else if (m_Decimate==QuadricDecimation)
+ {
+ vtkQuadricDecimation* decimate = vtkQuadricDecimation::New();
+ decimate->SetTargetReduction(m_TargetReduction);
+
+ decimate->SetInput(polydata);
+ polydata->Delete();
+ polydata = decimate->GetOutput();
+ polydata->Register(NULL);
+ decimate->Delete();
+ }
+
+ polydata->Update();
+ ProgressBar::GetInstance()->Progress();
+
+ polydata->SetSource(NULL);
+
+ if(polydata->GetNumberOfPoints() > 0)
+ {
+ mitk::Vector3D spacing = GetInput()->GetGeometry(time)->GetSpacing();
+
+ vtkPoints * points = polydata->GetPoints();
+ vtkMatrix4x4 *vtkmatrix = vtkMatrix4x4::New();
+ GetInput()->GetGeometry(time)->GetVtkTransform()->GetMatrix(vtkmatrix);
+ double (*matrix)[4] = vtkmatrix->Element;
+
+ unsigned int i,j;
+ for(i=0;i<3;++i)
+ for(j=0;j<3;++j)
+ matrix[i][j]/=spacing[j];
+
+ unsigned int n = points->GetNumberOfPoints();
+ vtkFloatingPointType point[3];
+
+ for (i = 0; i < n; i++)
+ {
+ points->GetPoint(i, point);
+ mitkVtkLinearTransformPoint(matrix,point,point);
+ points->SetPoint(i, point);
+ }
+ vtkmatrix->Delete();
+ }
+ ProgressBar::GetInstance()->Progress();
+
+ // determine point_data normals for the poly data points.
+ vtkSmartPointer<vtkPolyDataNormals> normalsGenerator = vtkSmartPointer<vtkPolyDataNormals>::New();
+ normalsGenerator->SetInput( polydata );
+
+ vtkSmartPointer<vtkCleanPolyData> cleanPolyDataFilter = vtkSmartPointer<vtkCleanPolyData>::New();
+ cleanPolyDataFilter->SetInput(normalsGenerator->GetOutput());
+ cleanPolyDataFilter->PieceInvariantOff();
+ cleanPolyDataFilter->ConvertLinesToPointsOff();
+ cleanPolyDataFilter->ConvertPolysToLinesOff();
+ cleanPolyDataFilter->ConvertStripsToPolysOff();
+ cleanPolyDataFilter->PointMergingOn();
+ cleanPolyDataFilter->Update();
+
+ surface->SetVtkPolyData(cleanPolyDataFilter->GetOutput(), time);
+ polydata->UnRegister(NULL);
+}
+
+
+void mitk::ImageToSurfaceFilter::GenerateData()
+{
+ mitk::Surface *surface = this->GetOutput();
+ mitk::Image * image = (mitk::Image*)GetInput();
+ if(image == NULL || !image->IsInitialized())
+ mitkThrow() << "No input image set, please set an valid input image!";
+
+ mitk::Image::RegionType outputRegion = image->GetRequestedRegion();
+
+ int tstart=outputRegion.GetIndex(3);
+ int tmax=tstart+outputRegion.GetSize(3); //GetSize()==1 - will aber 0 haben, wenn nicht zeitaufgeloest
+
+ if ((tmax-tstart) > 0)
+ {
+ ProgressBar::GetInstance()->AddStepsToDo( 4 * (tmax - tstart) );
+ }
+
+
+ int t;
+ for( t=tstart; t < tmax; ++t)
+ {
+ vtkImageData *vtkimagedata = image->GetVtkImageData(t);
+ CreateSurface(t,vtkimagedata,surface,m_Threshold);
+ ProgressBar::GetInstance()->Progress();
+ }
+}
+
+void mitk::ImageToSurfaceFilter::SetSmoothIteration(int smoothIteration)
+{
+ m_SmoothIteration = smoothIteration;
+}
+
+void mitk::ImageToSurfaceFilter::SetSmoothRelaxation(float smoothRelaxation)
+{
+ m_SmoothRelaxation = smoothRelaxation;
+}
+
+void mitk::ImageToSurfaceFilter::SetInput(const mitk::Image *image)
+{
+ // Process object is not const-correct so the const_cast is required here
+ this->ProcessObject::SetNthInput(0, const_cast< mitk::Image * >( image ) );
+}
+
+
+const mitk::Image *mitk::ImageToSurfaceFilter::GetInput(void)
+{
+ if (this->GetNumberOfInputs() < 1)
+ {
+ return 0;
+ }
+
+ return static_cast<const mitk::Image * >
+ ( this->ProcessObject::GetInput(0) );
+}
+
+
+void mitk::ImageToSurfaceFilter::GenerateOutputInformation()
+{
+ mitk::Image::ConstPointer inputImage =(mitk::Image*) this->GetInput();
+
+ //mitk::Image *inputImage = (mitk::Image*)this->GetImage();
+ mitk::Surface::Pointer output = this->GetOutput();
+
+ itkDebugMacro(<<"GenerateOutputInformation()");
+
+ if(inputImage.IsNull()) return;
+
+ //Set Data
+}
diff --git a/Core/Code/Algorithms/mitkImageToSurfaceFilter.h b/Core/Code/Algorithms/mitkImageToSurfaceFilter.h
new file mode 100644
index 0000000..be49699
--- /dev/null
+++ b/Core/Code/Algorithms/mitkImageToSurfaceFilter.h
@@ -0,0 +1,241 @@
+/*===================================================================
+
+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 _MITKIMAGETOSURFACEFILTER_h__
+#define _MITKIMAGETOSURFACEFILTER_h__
+
+#include <vtkPolyData.h>
+#include "MITKExports.h"
+#include <mitkCommon.h>
+#include <mitkSurfaceSource.h>
+#include <mitkSurface.h>
+
+#include <mitkImage.h>
+#include <vtkImageData.h>
+
+#include <vtkSmoothPolyDataFilter.h>
+#include <vtkMarchingCubes.h>
+
+
+namespace mitk {
+ /**
+ * @brief Converts pixel data to surface data by using a threshold
+ * The mitkImageToSurfaceFilter is used to create a new surface out of an mitk image. The filter
+ * uses a threshold to define the surface. It is based on the vtkMarchingCube algorithm. By default
+ * a vtkPolyData surface based on an input threshold for the input image will be created. Optional
+ * it is possible to reduce the number of triangles/polygones [SetDecimate(mitk::ImageToSurfaceFilter::DecimatePro) and SetTargetReduction (float _arg)]
+ * or smooth the surface-data [SetSmooth(true), SetSmoothIteration(int smoothIteration) and SetSmoothRelaxation(float smoothRelaxation)].
+ *
+ * The resulting vtk-surface has the same size as the input image. The surface
+ * can be generally smoothed by vtkDecimatePro reduce complexity of triangles
+ * and vtkSmoothPolyDataFilter to relax the mesh. Both are enabled by default
+ * and connected in the common way of pipelining in ITK. It's also possible
+ * to create time sliced surfaces.
+ *
+ * @ingroup ImageFilters
+ * @ingroup Process
+ */
+
+ class MITK_CORE_EXPORT ImageToSurfaceFilter : public SurfaceSource
+ {
+ public:
+
+ /*
+ * To decide whether a reduction of polygons in the created surface shall be
+ * done or not by using the vtkDecimatePro Filter. Till vtk 4.x an vtkDecimateFilter existed,
+ * but was patented. So since vtk 5.x it was replaced by the (worser?) vtkDecimateProFilter
+ * Maybe another Filter will come soon.
+ */
+ enum DecimationType {NoDecimation,DecimatePro,QuadricDecimation};
+
+ mitkClassMacro(ImageToSurfaceFilter, SurfaceSource);
+ itkNewMacro(Self);
+
+ /**
+ * For each image time slice a surface will be created. This method is
+ * called by Update().
+ */
+
+ virtual void GenerateData();
+
+ /**
+ * Initializes the output information ( i.e. the geometry information ) of
+ * the output of the filter
+ */
+ virtual void GenerateOutputInformation();
+
+ /**
+ * Returns a const reference to the input image (e.g. the original input image that ist used to create the surface)
+ */
+ const mitk::Image *GetInput(void);
+
+ /**
+ * Set the source image to create a surface for this filter class. As input every mitk
+ * 3D or 3D+t image can be used.
+ */
+ virtual void SetInput(const mitk::Image *image);
+
+ /**
+ * Set the number of iterations that is used to smooth the surface. Used is the vtkSmoothPolydataFilter that uses the laplacian filter. The higher the number of iterations that stronger the smooth-result
+ *
+ * @param smoothIteration As smoothIteration default in that case 50 was choosen. The VTK documentation recommends small relaxation factors and large numbers of iterations.
+ */
+ void SetSmoothIteration(int smoothIteration);
+
+ /**
+ * Set number of relaxation. Specify the relaxation factor for Laplacian
+ * smoothing. The VTK documentation recommends small relaxation factors
+ * and large numbers of iterations.
+ *
+ * @param smoothRelaxation As smoothRelaxation default in that case 0.1 was choosen. The VTK documentation recommends small relaxation factors and large numbers of iterations.
+ */
+ void SetSmoothRelaxation(float smoothRelaxation);
+
+
+ /**
+ * Threshold that is used to create the surface. All pixel in the input image that are higher than that
+ * value will be considered in the surface. The threshold referees to
+ * vtkMarchingCube. Default value is 1. See also SetThreshold (ScalarType _arg)
+ */
+ itkSetMacro(Threshold, ScalarType);
+
+ /**
+ * Get Threshold from vtkMarchingCube. Threshold can be manipulated by
+ * inherited classes.
+ */
+ itkGetConstMacro(Threshold, ScalarType);
+
+ /**
+ * Enables vtkSmoothPolyDataFilter. With Laplacian smoothing this filter
+ * will relax the surface. You can control the Filter by manipulating the
+ * number of iterations and the relaxing factor.
+ * */
+ itkSetMacro(Smooth,bool);
+
+ /*
+ * Enable/Disable surface smoothing.
+ */
+ itkBooleanMacro(Smooth);
+
+ /*
+ * Returns if surface smoothing is enabled
+ */
+ itkGetConstMacro(Smooth,bool);
+
+ /**
+ * Get the state of decimation mode to reduce triangle in the
+ * surface represantation. Modes can only be NoDecimation or DecimatePro
+ * (till vtk 4.x also Decimate)
+ * */
+ itkGetConstMacro(Decimate,DecimationType);
+
+ /**
+ * Enable the decimation filter to reduce the number of triangles in the
+ * mesh and produce a good approximation to the original image. The filter
+ * has support for vtk-5 and earlier versions. More detailed information
+ * check the vtkDecimatePro and vtkDecimate.
+ * */
+ itkSetMacro(Decimate,DecimationType);
+
+ /**
+ * Set desired TargetReduction of triangles in the range from 0.0 to 1.0.
+ * The destroyed triangles are in relation with the size of data. For example 0.9
+ * will reduce the data set to 10%.
+ *
+ * @param Set a TargetReduction float-value from 0.0 to 1.0
+ * */
+ itkSetMacro(TargetReduction, float);
+
+ /**
+ * Returns the reduction factor for the VtkDecimatePro Decimation Filter as a float value
+ */
+ itkGetConstMacro(TargetReduction, float);
+
+ /**
+ * Transforms a point by a 4x4 matrix
+ */
+ template <class T1, class T2, class T3>
+ inline void mitkVtkLinearTransformPoint(T1 matrix[4][4], T2 in[3], T3 out[3])
+ {
+ T3 x = matrix[0][0]*in[0]+matrix[0][1]*in[1]+matrix[0][2]*in[2]+matrix[0][3];
+ T3 y = matrix[1][0]*in[0]+matrix[1][1]*in[1]+matrix[1][2]*in[2]+matrix[1][3];
+ T3 z = matrix[2][0]*in[0]+matrix[2][1]*in[1]+matrix[2][2]*in[2]+matrix[2][3];
+ out[0] = x;
+ out[1] = y;
+ out[2] = z;
+ }
+
+ protected:
+
+
+ ImageToSurfaceFilter();
+
+ /**
+ * Destructor
+ * */
+ virtual ~ImageToSurfaceFilter();
+
+ /**
+ * With the given threshold vtkMarchingCube creates the surface. By default nothing a
+ * vtkPolyData surface based on a threshold of the input image will be created. Optional
+ * it is possible to reduce the number of triangles/polygones [SetDecimate(mitk::ImageToSurfaceFilter::DecimatePro) and SetTargetReduction (float _arg)]
+ * or smooth the data [SetSmooth(true), SetSmoothIteration(int smoothIteration) and SetSmoothRelaxation(float smoothRelaxation)].
+ *
+ * @param time selected slice or "0" for single
+ * @param *vtkimage input image
+ * @param *surface output
+ * @param threshold can be different from SetThreshold()
+ */
+ void CreateSurface(int time, vtkImageData *vtkimage, mitk::Surface * surface, const ScalarType threshold);
+
+ /**
+ * Flag whether the created surface shall be smoothed or not (default is "false"). SetSmooth (bool _arg)
+ * */
+ bool m_Smooth;
+
+ /**
+ * Decimation mode, default mode is "NoDecimation". See also SetDecimate (DecimationType _arg)
+ * */
+ DecimationType m_Decimate;
+
+ /**
+ * Threshold that is used to create the surface. All pixel in the input image that are higher than that
+ * value will be considered in the surface. Default value is 1. See also SetThreshold (ScalarType _arg)
+ * */
+ ScalarType m_Threshold;
+
+ /**
+ * The Reduction factor of the Decimation Filter for the created surface. See also SetTargetReduction (float _arg)
+ * */
+ float m_TargetReduction;
+
+ /**
+ * The Iteration value for the Smooth Filter of the created surface. See also SetSmoothIteration (int smoothIteration)
+ * */
+ int m_SmoothIteration;
+
+ /**
+ * The Relaxation value for the Smooth Filter of the created surface. See also SetSmoothRelaxation (float smoothRelaxation)
+ * */
+ float m_SmoothRelaxation;
+
+ };
+
+} // namespace mitk
+
+#endif //_MITKIMAGETOSURFACEFILTER_h__
+
+
diff --git a/Core/Code/Algorithms/mitkSurfaceToImageFilter.cpp b/Core/Code/Algorithms/mitkSurfaceToImageFilter.cpp
new file mode 100644
index 0000000..acf52c4
--- /dev/null
+++ b/Core/Code/Algorithms/mitkSurfaceToImageFilter.cpp
@@ -0,0 +1,209 @@
+/*===================================================================
+
+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.
+
+===================================================================*/
+
+#include "mitkSurfaceToImageFilter.h"
+#include "mitkTimeHelper.h"
+
+#include <vtkPolyData.h>
+#include <vtkPolyDataToImageStencil.h>
+#include <vtkImageStencil.h>
+#include <vtkImageData.h>
+#include <vtkPolyData.h>
+#include <vtkLinearTransform.h>
+#include <vtkTriangleFilter.h>
+#include <vtkLinearExtrusionFilter.h>
+#include <vtkDataSetTriangleFilter.h>
+#include <vtkImageThreshold.h>
+#include <vtkImageMathematics.h>
+#include <vtkPolyDataNormals.h>
+#include <vtkTransformPolyDataFilter.h>
+#include <vtkTransform.h>
+
+
+mitk::SurfaceToImageFilter::SurfaceToImageFilter()
+: m_MakeOutputBinary( false ),
+ m_BackgroundValue( -10000 )
+{
+}
+
+mitk::SurfaceToImageFilter::~SurfaceToImageFilter()
+{
+}
+
+void mitk::SurfaceToImageFilter::GenerateInputRequestedRegion()
+{
+ mitk::Image* output = this->GetOutput();
+ if((output->IsInitialized()==false) )
+ return;
+
+ GenerateTimeInInputRegion(output, const_cast< mitk::Image * > ( this->GetImage() ));
+}
+
+void mitk::SurfaceToImageFilter::GenerateOutputInformation()
+{
+ mitk::Image *inputImage = (mitk::Image*)this->GetImage();
+ mitk::Image::Pointer output = this->GetOutput();
+
+ itkDebugMacro(<<"GenerateOutputInformation()");
+
+ if((inputImage == NULL) ||
+ (inputImage->IsInitialized() == false) ||
+ (inputImage->GetTimeSlicedGeometry() == NULL)) return;
+
+ if (m_MakeOutputBinary)
+ output->Initialize(mitk::MakeScalarPixelType<unsigned char>() , *inputImage->GetTimeSlicedGeometry());
+ else
+ output->Initialize(inputImage->GetPixelType(), *inputImage->GetTimeSlicedGeometry());
+
+ output->SetPropertyList(inputImage->GetPropertyList()->Clone());
+}
+
+void mitk::SurfaceToImageFilter::GenerateData()
+{
+ mitk::Image::ConstPointer inputImage = this->GetImage();
+ mitk::Image::Pointer output = this->GetOutput();
+
+ if(inputImage.IsNull())
+ return;
+
+ if(output->IsInitialized()==false )
+ return;
+
+ mitk::Image::RegionType outputRegion = output->GetRequestedRegion();
+
+ int tstart=outputRegion.GetIndex(3);
+ int tmax=tstart+outputRegion.GetSize(3);
+
+ if ( tmax > 0)
+ {
+ int t;
+ for(t=tstart;t<tmax;++t)
+ {
+ Stencil3DImage( t );
+ }
+ }
+ else
+ {
+ Stencil3DImage( 0 );
+ }
+}
+
+void mitk::SurfaceToImageFilter::Stencil3DImage(int time)
+{
+ const mitk::TimeSlicedGeometry *surfaceTimeGeometry = GetInput()->GetTimeSlicedGeometry();
+ const mitk::TimeSlicedGeometry *imageTimeGeometry = GetImage()->GetTimeSlicedGeometry();
+
+ // Convert time step from image time-frame to surface time-frame
+ int surfaceTimeStep = surfaceTimeGeometry->TimeStepToTimeStep( imageTimeGeometry, time );
+
+ vtkPolyData * polydata = ( (mitk::Surface*)GetInput() )->GetVtkPolyData( surfaceTimeStep );
+
+ vtkTransformPolyDataFilter * move=vtkTransformPolyDataFilter::New();
+ move->SetInput(polydata);
+ move->ReleaseDataFlagOn();
+
+ vtkTransform *transform=vtkTransform::New();
+ Geometry3D* geometry = surfaceTimeGeometry->GetGeometry3D( surfaceTimeStep );
+ geometry->TransferItkToVtkTransform();
+ transform->PostMultiply();
+ transform->Concatenate(geometry->GetVtkTransform()->GetMatrix());
+ // take image geometry into account. vtk-Image information will be changed to unit spacing and zero origin below.
+ Geometry3D* imageGeometry = imageTimeGeometry->GetGeometry3D(time);
+ imageGeometry->TransferItkToVtkTransform();
+ transform->Concatenate(imageGeometry->GetVtkTransform()->GetLinearInverse());
+ move->SetTransform(transform);
+ transform->Delete();
+
+ vtkPolyDataNormals * normalsFilter = vtkPolyDataNormals::New();
+ normalsFilter->SetFeatureAngle(50);
+ normalsFilter->SetConsistency(1);
+ normalsFilter->SetSplitting(1);
+ normalsFilter->SetFlipNormals(0);
+ normalsFilter->ReleaseDataFlagOn();
+
+ normalsFilter->SetInput( move->GetOutput() );
+ move->Delete();
+
+ vtkPolyDataToImageStencil * surfaceConverter = vtkPolyDataToImageStencil::New();
+ surfaceConverter->SetTolerance( 0.0 );
+ surfaceConverter->ReleaseDataFlagOn();
+
+ surfaceConverter->SetInput( normalsFilter->GetOutput() );
+ normalsFilter->Delete();
+
+ mitk::Image::Pointer binaryImage = mitk::Image::New();
+
+ if (m_MakeOutputBinary)
+ {
+ binaryImage->Initialize(mitk::MakeScalarPixelType<unsigned char>(), *this->GetImage()->GetTimeSlicedGeometry());
+
+ unsigned int size = sizeof(unsigned char);
+ for (unsigned int i = 0; i < binaryImage->GetDimension(); ++i)
+ size *= binaryImage->GetDimension(i);
+
+ memset(binaryImage->GetData(), 1, size);
+ }
+
+ vtkImageData *image = m_MakeOutputBinary
+ ? binaryImage->GetVtkImageData(time)
+ : const_cast<mitk::Image *>(this->GetImage())->GetVtkImageData(time);
+
+ // Create stencil and use numerical minimum of pixel type as background value
+ vtkImageStencil *stencil = vtkImageStencil::New();
+ stencil->SetInput(image);
+ stencil->ReverseStencilOff();
+ stencil->ReleaseDataFlagOn();
+ stencil->SetStencil(surfaceConverter->GetOutput());
+ surfaceConverter->Delete();
+
+ stencil->SetBackgroundValue(m_MakeOutputBinary ? 0 : m_BackgroundValue);
+ stencil->Update();
+
+ mitk::Image::Pointer output = this->GetOutput();
+ output->SetVolume( stencil->GetOutput()->GetScalarPointer(), time );
+ MITK_INFO << "stencil ref count: " << stencil->GetReferenceCount() << std::endl;
+
+ stencil->Delete();
+}
+
+
+const mitk::Surface *mitk::SurfaceToImageFilter::GetInput(void)
+{
+ if (this->GetNumberOfInputs() < 1)
+ {
+ return 0;
+ }
+
+ return static_cast<const mitk::Surface * >
+ ( this->ProcessObject::GetInput(0) );
+}
+
+void mitk::SurfaceToImageFilter::SetInput(const mitk::Surface *input)
+{
+ // Process object is not const-correct so the const_cast is required here
+ this->ProcessObject::SetNthInput(0,
+ const_cast< mitk::Surface * >( input ) );
+}
+
+void mitk::SurfaceToImageFilter::SetImage(const mitk::Image *source)
+{
+ this->ProcessObject::SetNthInput( 1, const_cast< mitk::Image * >( source ) );
+}
+
+const mitk::Image *mitk::SurfaceToImageFilter::GetImage(void)
+{
+ return static_cast< const mitk::Image * >(this->ProcessObject::GetInput(1));
+}
diff --git a/Core/Code/Algorithms/mitkSurfaceToImageFilter.h b/Core/Code/Algorithms/mitkSurfaceToImageFilter.h
new file mode 100644
index 0000000..41b4192
--- /dev/null
+++ b/Core/Code/Algorithms/mitkSurfaceToImageFilter.h
@@ -0,0 +1,96 @@
+/*===================================================================
+
+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 _mitkSurfaceToImageFilter_h__
+#define _mitkSurfaceToImageFilter_h__
+
+#include "mitkCommon.h"
+#include "MitkExports.h"
+#include "mitkImageSource.h"
+#include "mitkSurface.h"
+//#include "mitkImage.h"
+
+class vtkPolyData;
+
+namespace mitk {
+
+//class Mesh;
+//class VectorOfContourLines;
+
+/**
+ *
+ * @brief Converts surface data to pixel data. Requires a surface and an
+ * image, which header information defines the output image.
+ *
+ * The resulting image has the same dimension, size, and Geometry3D
+ * as the input image. The image is cut using a vtkStencil.
+ * The user can decide if he wants to keep the original values or create a
+ * binary image by setting MakeBinaryOutputOn (default is \a false). If
+ * set to \a true all voxels inside the surface are set to one and all
+ * outside voxel are set to zero.
+ *
+ * NOTE: Since the reference input image is passed to the vtkStencil in
+ * any case, the image needs to be initialized with pixel values greater than
+ * the numerical minimum of the used pixel type (e.g. at least -127 for
+ * unsigned char images, etc.) to produce a correct binary image
+ * representation of the surface in MakeOutputBinary mode.
+ *
+ * @ingroup SurfaceFilters
+ * @ingroup Process
+ */
+class MITK_CORE_EXPORT SurfaceToImageFilter : public ImageSource
+{
+public:
+ mitkClassMacro(SurfaceToImageFilter, ImageSource);
+ itkNewMacro(Self);
+
+ itkSetMacro(MakeOutputBinary, bool);
+ itkGetMacro(MakeOutputBinary, bool);
+ itkBooleanMacro(MakeOutputBinary);
+
+ itkGetConstMacro(BackgroundValue,float);
+ itkSetMacro(BackgroundValue,float);
+
+ virtual void GenerateInputRequestedRegion();
+
+ virtual void GenerateOutputInformation();
+
+ virtual void GenerateData();
+
+ const mitk::Surface *GetInput(void);
+
+ virtual void SetInput(const mitk::Surface *surface);
+
+ void SetImage(const mitk::Image *source);
+
+ const mitk::Image *GetImage(void);
+
+protected:
+ SurfaceToImageFilter();
+
+ virtual ~SurfaceToImageFilter();
+
+ void Stencil3DImage(int time = 0);
+
+ bool m_MakeOutputBinary;
+
+ float m_BackgroundValue;
+
+};
+
+} // namespace mitk
+
+#endif /* MITKCOONSPATCHFILTER_H_HEADER_INCLUDED_C10B22CD */
diff --git a/Core/Code/Algorithms/mitkTimeHelper.h b/Core/Code/Algorithms/mitkTimeHelper.h
new file mode 100644
index 0000000..cb930d8
--- /dev/null
+++ b/Core/Code/Algorithms/mitkTimeHelper.h
@@ -0,0 +1,81 @@
+/*===================================================================
+
+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 MITKTIMEHELPER_H_HEADER_INCLUDED_C1C2FCD2
+#define MITKTIMEHELPER_H_HEADER_INCLUDED_C1C2FCD2
+
+namespace mitk
+{
+
+//## @brief convert the start- and end-index-time of output-region in
+//## start- and end-index-time of input-region via millisecond-time
+template <class TOutputRegion, class TInputRegion>
+void ITK_EXPORT GenerateTimeInInputRegion(const mitk::TimeSlicedGeometry *outputTimeGeometry, const TOutputRegion& outputRegion, const mitk::TimeSlicedGeometry *inputTimeGeometry, TInputRegion& inputRegion)
+{
+ assert(outputTimeGeometry!=NULL);
+ assert(inputTimeGeometry!=NULL);
+
+ // convert the start-index-time of output in start-index-time of input via millisecond-time
+ ScalarType timeInMS = outputTimeGeometry->TimeStepToMS(outputRegion.GetIndex(3));
+ int timestep = inputTimeGeometry->MSToTimeStep( timeInMS );
+ if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( inputTimeGeometry->IsValidTime( timestep ) ) )
+ inputRegion.SetIndex( 3, timestep );
+ else
+ inputRegion.SetIndex( 3, 0 );
+ // convert the end-index-time of output in end-index-time of input via millisecond-time
+ timeInMS = outputTimeGeometry->TimeStepToMS(outputRegion.GetIndex(3)+outputRegion.GetSize(3)-1);
+ timestep = inputTimeGeometry->MSToTimeStep( timeInMS );
+ if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( outputTimeGeometry->IsValidTime( timestep ) ) )
+ inputRegion.SetSize( 3, timestep - inputRegion.GetIndex(3) + 1 );
+ else
+ inputRegion.SetSize( 3, 1 );
+}
+
+//##Documentation
+//## @brief convert the start- and end-index-time of output in
+//## start- and end-index-time of input1 and input2 via millisecond-time
+template <class TOutputData, class TInputData>
+void ITK_EXPORT GenerateTimeInInputRegion(const TOutputData* output, TInputData* input)
+{
+ assert(output!=NULL);
+ assert(input!=NULL);
+
+ const typename TOutputData::RegionType& outputRegion = output->GetRequestedRegion();
+ typename TInputData::RegionType inputRegion;
+
+ if(outputRegion.GetSize(3)<1)
+ {
+ typename TInputData::RegionType::SizeType inputsize;
+ inputsize.Fill(0);
+ inputRegion.SetSize(inputsize);
+ input->SetRequestedRegion( &inputRegion );
+ }
+
+ // convert the start-index-time of output in start-index-time of input via millisecond-time
+ inputRegion = input->GetRequestedRegion();
+ GenerateTimeInInputRegion(output->GetTimeSlicedGeometry(), outputRegion, input->GetTimeSlicedGeometry(), inputRegion);
+ input->SetRequestedRegion( &inputRegion );
+}
+
+} // end namespace mitk
+
+//#ifndef ITK_MANUAL_INSTANTIATION
+//#include "mitkTimeHelper.txx"
+#include "MitkExports.h"
+//#endif
+
+#endif // MITKTIMEHELPER_H_HEADER_INCLUDED_C1C2FCD2
diff --git a/Core/Code/files.cmake b/Core/Code/files.cmake
index 5ca0ce8..8eed017 100644
--- a/Core/Code/files.cmake
+++ b/Core/Code/files.cmake
@@ -9,6 +9,7 @@ set(H_FILES
Algorithms/itkTotalVariationDenoisingImageFilter.txx
Algorithms/itkTotalVariationSingleIterationImageFilter.h
Algorithms/itkTotalVariationSingleIterationImageFilter.txx
+ Algorithms/mitkAutoCropImageFilter.h
Algorithms/mitkBilateralFilter.h
Algorithms/mitkBilateralFilter.cpp
Algorithms/mitkInstantiateAccessFunctions.h
@@ -73,6 +74,7 @@ set(H_FILES
)
set(CPP_FILES
+ Algorithms/mitkAutoCropImageFilter.cpp
Algorithms/mitkBaseDataSource.cpp
Algorithms/mitkBaseProcess.cpp
Algorithms/mitkDataNodeSource.cpp
@@ -83,11 +85,13 @@ set(CPP_FILES
Algorithms/mitkImageSource.cpp
Algorithms/mitkImageTimeSelector.cpp
Algorithms/mitkImageToImageFilter.cpp
+ Algorithms/mitkImageToSurfaceFilter.cpp
Algorithms/mitkPointSetSource.cpp
Algorithms/mitkPointSetToPointSetFilter.cpp
Algorithms/mitkRGBToRGBACastImageFilter.cpp
Algorithms/mitkSubImageSelector.cpp
Algorithms/mitkSurfaceSource.cpp
+ Algorithms/mitkSurfaceToImageFilter.cpp
Algorithms/mitkSurfaceToSurfaceFilter.cpp
Algorithms/mitkUIDGenerator.cpp
Algorithms/mitkVolumeCalculator.cpp
diff --git a/Modules/CMakeLists.txt b/Modules/CMakeLists.txt
index 9005b35..94024e2 100644
--- a/Modules/CMakeLists.txt
+++ b/Modules/CMakeLists.txt
@@ -8,6 +8,7 @@ set(module_dirs
ImageStatistics
LegacyAdaptors
IpPicSupport
+ SurfaceInterpolation
MitkExt
SceneSerialization
Segmentation
diff --git a/Modules/MitkExt/Algorithms/mitkAutoCropImageFilter.cpp b/Modules/MitkExt/Algorithms/mitkAutoCropImageFilter.cpp
deleted file mode 100644
index 3fdb301..0000000
--- a/Modules/MitkExt/Algorithms/mitkAutoCropImageFilter.cpp
+++ /dev/null
@@ -1,371 +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.
-
-===================================================================*/
-
-#include "mitkAutoCropImageFilter.h"
-
-#include "mitkImageCast.h"
-#include "mitkImageAccessByItk.h"
-#include "mitkGeometry3D.h"
-#include "mitkStatusBar.h"
-
-#include "mitkPlaneGeometry.h"
-
-#include <itkImageRegionConstIterator.h>
-#include <itkRegionOfInterestImageFilter.h>
-
-
-mitk::AutoCropImageFilter::AutoCropImageFilter()
-: m_BackgroundValue(0),
- m_MarginFactor(1.0),
- m_TimeSelector(NULL),
- m_OverrideCroppingRegion(false)
-{
-
-}
-
-
-mitk::AutoCropImageFilter::~AutoCropImageFilter()
-{
-
-}
-
-
-template < typename TPixel, unsigned int VImageDimension>
-void mitk::AutoCropImageFilter::ITKCrop3DImage( itk::Image< TPixel, VImageDimension >* inputItkImage, unsigned int timestep)
-{
-
- if (inputItkImage == NULL)
- {
- mitk::StatusBar::GetInstance()->DisplayErrorText ("An internal error occurred. Can't convert Image. Please report to bugs@mitk.org");
- MITK_ERROR << "image is NULL...returning" << std::endl;
- return;
- }
-
- typedef itk::Image< TPixel, VImageDimension > InternalImageType;
- typedef typename InternalImageType::Pointer InternalImagePointer;
-
- typedef itk::RegionOfInterestImageFilter < InternalImageType, InternalImageType > ROIFilterType;
- typedef typename itk::RegionOfInterestImageFilter < InternalImageType, InternalImageType >::Pointer ROIFilterPointer;
-
- InternalImagePointer outputItk = InternalImageType::New();
-
- ROIFilterPointer roiFilter = ROIFilterType::New();
- roiFilter->SetInput(0,inputItkImage);
- roiFilter->SetRegionOfInterest(this->GetCroppingRegion());
- roiFilter->Update();
- outputItk = roiFilter->GetOutput();
- outputItk->DisconnectPipeline();
-
- mitk::Image::Pointer newMitkImage = mitk::Image::New();
- mitk::CastToMitkImage( outputItk, newMitkImage );
- MITK_INFO << "Crop-Output dimension: " << (newMitkImage->GetDimension() == 3) << " Filter-Output dimension: "<<this->GetOutput()->GetDimension()<< " Timestep: " << timestep;
-
-// const mitk::ChannelDescriptor desc = newMitkImage->GetChannelDescriptor(0);
-// unsigned char* image3D = desc.GetData();
-// this->GetOutput()->SetVolume( (void*) &image3D , timestep );
-
- this->GetOutput()->SetVolume( newMitkImage->GetData(), timestep);
-// this->SetOutput(newMitkImage);
-}
-
-void mitk::AutoCropImageFilter::GenerateOutputInformation()
-{
- mitk::Image::Pointer input = const_cast<mitk::Image*> (this->GetInput());
- mitk::Image::Pointer output = this->GetOutput();
-
- if(input->GetDimension() <= 2)
- {
- MITK_ERROR << "Only 3D any 4D images are supported." << std::endl;
- return;
- }
-
- ComputeNewImageBounds();
-
- if ((output->IsInitialized()) && (output->GetPipelineMTime() <= m_TimeOfHeaderInitialization.GetMTime()))
- return;
-
- itkDebugMacro(<<"GenerateOutputInformation()");
-
- // PART I: initialize input requested region. We do this already here (and not
- // later when GenerateInputRequestedRegion() is called), because we
- // also need the information to setup the output.
-
- // pre-initialize input-requested-region to largest-possible-region
- // and correct time-region; spatial part will be cropped by
- // bounding-box of bounding-object below
- m_InputRequestedRegion = input->GetLargestPossibleRegion();
-
- // build region out of index and size calculated in ComputeNewImageBounds()
-
- mitk::SlicedData::IndexType index;
- index[0] = m_RegionIndex[0];
- index[1] = m_RegionIndex[1];
- index[2] = m_RegionIndex[2];
- index[3] = m_InputRequestedRegion.GetIndex()[3];
- index[4] = m_InputRequestedRegion.GetIndex()[4];
-
- mitk::SlicedData::SizeType size;
- size[0] = m_RegionSize[0];
- size[1] = m_RegionSize[1];
- size[2] = m_RegionSize[2];
- size[3] = m_InputRequestedRegion.GetSize()[3];
- size[4] = m_InputRequestedRegion.GetSize()[4];
-
- mitk::SlicedData::RegionType cropRegion(index, size);
-
- // crop input-requested-region with cropping region computed from the image data
- if(m_InputRequestedRegion.Crop(cropRegion)==false)
- {
- // crop not possible => do nothing: set time size to 0.
- size.Fill(0);
- m_InputRequestedRegion.SetSize(size);
- return;
- }
-
- // set input-requested-region, because we access it later in
- // GenerateInputRequestedRegion (there we just set the time)
- input->SetRequestedRegion(&m_InputRequestedRegion);
-
-
- // PART II: initialize output image
- unsigned int dimension = input->GetDimension();
- unsigned int *dimensions = new unsigned int [dimension];
- itk2vtk(m_InputRequestedRegion.GetSize(), dimensions);
- if(dimension>3)
- memcpy(dimensions+3, input->GetDimensions()+3, (dimension-3)*sizeof(unsigned int));
-
- // create basic slicedGeometry that will be initialized below
- output->Initialize(mitk::PixelType( GetOutputPixelType() ), dimension, dimensions);
- delete [] dimensions;
-
- //clone the IndexToWorldTransform from the input, otherwise we will overwrite it, when adjusting the origin of the output image!!
- itk::ScalableAffineTransform< mitk::ScalarType,3 >::Pointer cloneTransform = itk::ScalableAffineTransform< mitk::ScalarType,3 >::New();
- cloneTransform->Compose(input->GetGeometry()->GetIndexToWorldTransform());
- output->GetGeometry()->SetIndexToWorldTransform( cloneTransform.GetPointer() );
-
- // Position the output Image to match the corresponding region of the input image
- mitk::SlicedGeometry3D* slicedGeometry = output->GetSlicedGeometry();
- mitk::SlicedGeometry3D::Pointer inputGeometry = input->GetSlicedGeometry();
- const mitk::SlicedData::IndexType& start = m_InputRequestedRegion.GetIndex();
- mitk::Point3D origin; vtk2itk(start, origin);
- input->GetSlicedGeometry()->IndexToWorld(origin, origin);
- slicedGeometry->SetOrigin(origin);
-
- // get the PlaneGeometry for the first slice of the original image
- mitk::PlaneGeometry::Pointer plane = dynamic_cast<mitk::PlaneGeometry*>( inputGeometry->GetGeometry2D( 0 )->Clone().GetPointer() );
- assert( plane );
-
- // re-initialize the plane according to the new requirements:
- // dimensions of the cropped image
- // right- and down-vector as well as spacing do not change, so use the ones from
- // input image
- ScalarType dimX = output->GetDimensions()[0];
- ScalarType dimY = output->GetDimensions()[1];
- mitk::Vector3D right = plane->GetAxisVector(0);
- mitk::Vector3D down = plane->GetAxisVector(1);
- mitk::Vector3D spacing = plane->GetSpacing();
- plane->InitializeStandardPlane( dimX, dimY, right, down, &spacing );
- // set the new origin on the PlaneGeometry as well
- plane->SetOrigin(origin);
-
- // re-initialize the slicedGeometry with the correct planeGeometry
- // in order to get a fully initialized SlicedGeometry3D
- slicedGeometry->InitializeEvenlySpaced( plane, inputGeometry->GetSpacing()[2], output->GetSlicedGeometry()->GetSlices() );
-
-
- mitk::TimeSlicedGeometry* timeSlicedGeometry = output->GetTimeSlicedGeometry();
- timeSlicedGeometry->InitializeEvenlyTimed(slicedGeometry, output->GetDimension(3));
- timeSlicedGeometry->CopyTimes(input->GetTimeSlicedGeometry());
-
- m_TimeOfHeaderInitialization.Modified();
-
- output->SetPropertyList(input->GetPropertyList()->Clone());
-}
-
-void mitk::AutoCropImageFilter::GenerateData()
-{
- mitk::Image::ConstPointer input = this->GetInput();
- mitk::Image::Pointer output = this->GetOutput();
-
- if(input.IsNull())
- return;
-
- if(input->GetDimension() <= 2)
- {
- MITK_ERROR << "Only 3D and 4D images supported";
- return;
- }
-
- if((output->IsInitialized()==false) )
- return;
-
- if( m_TimeSelector.IsNull() ) m_TimeSelector = mitk::ImageTimeSelector::New();
-
- m_TimeSelector->SetInput(input);
-
- mitk::SlicedData::RegionType outputRegion = input->GetRequestedRegion();
-
- int tstart = outputRegion.GetIndex(3);
- int tmax = tstart + outputRegion.GetSize(3);
-
- for( int timestep=tstart;timestep<tmax;++timestep )
- {
- m_TimeSelector->SetTimeNr(timestep);
- m_TimeSelector->UpdateLargestPossibleRegion();
-
- AccessFixedDimensionByItk_1( m_TimeSelector->GetOutput(), ITKCrop3DImage, 3, timestep );
- }
-
- // this->GetOutput()->Update(); // Not sure if this is necessary...
-
- m_TimeOfHeaderInitialization.Modified();
-}
-
-void mitk::AutoCropImageFilter::ComputeNewImageBounds()
-{
- mitk::Image::ConstPointer inputMitk = this->GetInput();
-
- if (m_OverrideCroppingRegion)
- {
- for (unsigned int i=0; i<3; ++i)
- {
- m_RegionIndex[i] = m_CroppingRegion.GetIndex()[i];
- m_RegionSize[i] = m_CroppingRegion.GetSize()[i];
-
- if (m_RegionIndex[i] >= inputMitk->GetDimension(i))
- {
- itkExceptionMacro("Cropping index is not inside the image. "
- << std::endl << "Index:"
- << std::endl << m_CroppingRegion.GetIndex()
- << std::endl << "Size:"
- << std::endl << m_CroppingRegion.GetSize());
- }
-
- if (m_RegionIndex[i] + m_RegionSize[i] >= inputMitk->GetDimension(i))
- {
- m_RegionSize[i] = inputMitk->GetDimension(i) - m_RegionIndex[i];
- }
- }
-
- for (unsigned int i=0; i<3; ++i)
- {
- m_RegionIndex[i] = m_CroppingRegion.GetIndex()[i];
- m_RegionSize[i] = m_CroppingRegion.GetSize()[i];
-
- }
- }
- else
- {
- // Check if a 3D or 4D image is present
- unsigned int timeSteps = 1;
- if (inputMitk->GetDimension() == 4 )
- timeSteps = inputMitk->GetDimension(3);
-
- ImageType::IndexType minima,maxima;
-
- if (inputMitk->GetDimension() == 4)
- {
- // initialize with time step 0
- m_TimeSelector = mitk::ImageTimeSelector::New();
- m_TimeSelector->SetInput( inputMitk );
- m_TimeSelector->SetTimeNr( 0 );
- m_TimeSelector->UpdateLargestPossibleRegion();
- inputMitk = m_TimeSelector->GetOutput();
- }
-
- ImagePointer inputItk = ImageType::New();
- mitk::CastToItkImage( inputMitk , inputItk );
-
- // it is assumed that all volumes in a time series have the same 3D dimensions
- ImageType::RegionType origRegion = inputItk->GetLargestPossibleRegion();
-
- // Initialize min and max on the first (or only) time step
- maxima = inputItk->GetLargestPossibleRegion().GetIndex();
- minima[0] = inputItk->GetLargestPossibleRegion().GetSize()[0];
- minima[1] = inputItk->GetLargestPossibleRegion().GetSize()[1];
- minima[2] = inputItk->GetLargestPossibleRegion().GetSize()[2];
-
- typedef itk::ImageRegionConstIterator< ImageType > ConstIteratorType;
-
- for(unsigned int idx = 0; idx < timeSteps; ++idx)
- {
- // if 4D image, update time step and itk image
- if( idx > 0)
- {
- m_TimeSelector->SetTimeNr( idx );
- m_TimeSelector->UpdateLargestPossibleRegion();
- inputMitk = m_TimeSelector->GetOutput();
- mitk::CastToItkImage( inputMitk , inputItk );
- }
-
- ConstIteratorType inIt( inputItk, origRegion );
-
- for ( inIt.GoToBegin(); !inIt.IsAtEnd(); ++inIt)
- {
- float pix_val = inIt.Get();
- if ( fabs(pix_val - m_BackgroundValue) > mitk::eps )
- {
- for (int i=0; i < 3; i++)
- {
- minima[i] = vnl_math_min((int)minima[i],(int)(inIt.GetIndex()[i]));
- maxima[i] = vnl_math_max((int)maxima[i],(int)(inIt.GetIndex()[i]));
- }
- }
- }
- }
-
- typedef ImageType::RegionType::SizeType::SizeValueType SizeValueType;
-
- m_RegionSize[0] = (SizeValueType)(m_MarginFactor * (maxima[0] - minima[0] + 1 ));
- m_RegionSize[1] = (SizeValueType)(m_MarginFactor * (maxima[1] - minima[1] + 1 ));
- m_RegionSize[2] = (SizeValueType)(m_MarginFactor * (maxima[2] - minima[2] + 1 ));
- m_RegionIndex = minima;
-
- m_RegionIndex[0] -= (m_RegionSize[0] - maxima[0] + minima[0] - 1 )/2;
- m_RegionIndex[1] -= (m_RegionSize[1] - maxima[1] + minima[1] - 1 )/2;
- m_RegionIndex[2] -= (m_RegionSize[2] - maxima[2] + minima[2] - 1 )/2;
-
- ImageType::RegionType cropRegion(m_RegionIndex,m_RegionSize);
- origRegion.Crop(cropRegion);
-
- m_RegionSize[0] = origRegion.GetSize()[0];
- m_RegionSize[1] = origRegion.GetSize()[1];
- m_RegionSize[2] = origRegion.GetSize()[2];
-
- m_RegionIndex[0] = origRegion.GetIndex()[0];
- m_RegionIndex[1] = origRegion.GetIndex()[1];
- m_RegionIndex[2] = origRegion.GetIndex()[2];
-
- m_CroppingRegion = origRegion;
- }
-}
-
-
-void mitk::AutoCropImageFilter::GenerateInputRequestedRegion()
-{
-
-}
-
-const mitk::PixelType mitk::AutoCropImageFilter::GetOutputPixelType()
-{
- return this->GetInput()->GetPixelType();
-}
-
-void mitk::AutoCropImageFilter::SetCroppingRegion(RegionType overrideRegion)
-{
- m_CroppingRegion = overrideRegion;
- m_OverrideCroppingRegion = true;
-}
diff --git a/Modules/MitkExt/Algorithms/mitkAutoCropImageFilter.h b/Modules/MitkExt/Algorithms/mitkAutoCropImageFilter.h
deleted file mode 100644
index 3cc44d5..0000000
--- a/Modules/MitkExt/Algorithms/mitkAutoCropImageFilter.h
+++ /dev/null
@@ -1,130 +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 __mitkAutoCropImageFilter_h_
-#define __mitkAutoCropImageFilter_h_
-
-#include "mitkCommon.h"
-#include "MitkExtExports.h"
-#include "mitkSubImageSelector.h"
-#include "mitkImageTimeSelector.h"
-
-#include <itkImageRegion.h>
-#include <itkImage.h>
-
-namespace mitk {
-
-/**
- *
- * @brief Shrink the image borders to a minimum considering a background color.
- *
- * This filter determines the smallest bounding box of all pixels different
- * from the background, and returns an output image which has been cropped to this size.
- * The box calculated this way is not the smallest possible box, but the box with the
- * smallest sides perpendicular to the world coordinate system.
- *
- * The filter works on 3D and 4D image data. For the 4D case, the smallest box is
- * calculated with side lengths as the maximum of single side lengths from all time steps.
- *
- * 2D images are not supported, and will never be.
- *
- * It is also possible to set the region to be cropped manually using the
- * SetCroppingRegion() method.
- *
- * A margin can be set to enlarge the cropped region with a constant factor in all
- * directions around the smallest possible.
- *
- *
- * @ingroup Process
- *
- * @author Thomas Boettger; revised by Tobias Schwarz and Daniel Stein; Division of Medical
- * and Biological Informatics
- *
- */
-
-class MitkExt_EXPORT AutoCropImageFilter : public SubImageSelector
-{
-public:
-
- typedef itk::ImageRegion<3> RegionType;
-
- mitkClassMacro(AutoCropImageFilter, SubImageSelector);
-
- itkNewMacro(Self);
-
- itkGetConstMacro(BackgroundValue,float);
- itkSetMacro(BackgroundValue,float);
-
- itkGetConstMacro(MarginFactor,float);
- itkSetMacro(MarginFactor,float);
-
- itkGetMacro(CroppingRegion, RegionType);
-
- // Use this method to manually set a region
- void SetCroppingRegion(RegionType overrideRegion);
-
- virtual const PixelType GetOutputPixelType();
-
-protected:
-
- // default constructor
- AutoCropImageFilter();
-
- // default destructor
- virtual ~AutoCropImageFilter();
-
- // This method calculates the actual smallest box
- void ComputeNewImageBounds();
-
- // Crops the image using the itk::RegionOfInterestImageFilter and creates the new output image
- template < typename TPixel, unsigned int VImageDimension>
- void ITKCrop3DImage( itk::Image< TPixel, VImageDimension >* inputItkImage, unsigned int timestep );
-
- // Here, the output image is initialized by the input and the newly calculated region
- virtual void GenerateOutputInformation();
-
- // Purposely not implemented
- virtual void GenerateInputRequestedRegion();
-
- // Crops the image on all time steps
- virtual void GenerateData();
-
- float m_BackgroundValue;
-
- RegionType m_CroppingRegion;
-
- float m_MarginFactor;
-
- typedef itk::Image<float,3> ImageType;
- typedef ImageType::Pointer ImagePointer;
-
- RegionType::SizeType m_RegionSize;
- RegionType::IndexType m_RegionIndex;
-
- mitk::ImageTimeSelector::Pointer m_TimeSelector;
-
- mitk::SlicedData::RegionType m_InputRequestedRegion;
- itk::TimeStamp m_TimeOfHeaderInitialization;
-
- bool m_OverrideCroppingRegion;
-
-};
-
-} // namespace mitk
-
-#endif
-
-
diff --git a/Modules/MitkExt/Algorithms/mitkImageToSurfaceFilter.cpp b/Modules/MitkExt/Algorithms/mitkImageToSurfaceFilter.cpp
deleted file mode 100644
index 6b2633a..0000000
--- a/Modules/MitkExt/Algorithms/mitkImageToSurfaceFilter.cpp
+++ /dev/null
@@ -1,234 +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.
-
-===================================================================*/
-
-#include <mitkImageToSurfaceFilter.h>
-#include "mitkException.h"
-#include <vtkImageData.h>
-#include <vtkDecimatePro.h>
-#include <vtkImageChangeInformation.h>
-#include <vtkLinearTransform.h>
-#include <vtkMath.h>
-#include <vtkMatrix4x4.h>
-#include <vtkQuadricDecimation.h>
-
-#include <vtkSmartPointer.h>
-#include <vtkPolyDataNormals.h>
-#include <vtkCleanPolyData.h>
-
-#include "mitkProgressBar.h"
-
-mitk::ImageToSurfaceFilter::ImageToSurfaceFilter():
- m_Smooth(false),
- m_Decimate( NoDecimation),
- m_Threshold(1.0),
- m_TargetReduction(0.95f),
- m_SmoothIteration(50),
- m_SmoothRelaxation(0.1)
-{
-}
-
-mitk::ImageToSurfaceFilter::~ImageToSurfaceFilter()
-{
-}
-
-void mitk::ImageToSurfaceFilter::CreateSurface(int time, vtkImageData *vtkimage, mitk::Surface * surface, const ScalarType threshold)
-{
- vtkImageChangeInformation *indexCoordinatesImageFilter = vtkImageChangeInformation::New();
- indexCoordinatesImageFilter->SetInput(vtkimage);
- indexCoordinatesImageFilter->SetOutputOrigin(0.0,0.0,0.0);
-
- //MarchingCube -->create Surface
- vtkMarchingCubes *skinExtractor = vtkMarchingCubes::New();
- skinExtractor->ComputeScalarsOff();
- skinExtractor->SetInput(indexCoordinatesImageFilter->GetOutput());//RC++
- indexCoordinatesImageFilter->Delete();
- skinExtractor->SetValue(0, threshold);
-
- vtkPolyData *polydata;
- polydata = skinExtractor->GetOutput();
- polydata->Register(NULL);//RC++
- skinExtractor->Delete();
-
- if (m_Smooth)
- {
- vtkSmoothPolyDataFilter *smoother = vtkSmoothPolyDataFilter::New();
- //read poly1 (poly1 can be the original polygon, or the decimated polygon)
- smoother->SetInput(polydata);//RC++
- smoother->SetNumberOfIterations( m_SmoothIteration );
- smoother->SetRelaxationFactor( m_SmoothRelaxation );
- smoother->SetFeatureAngle( 60 );
- smoother->FeatureEdgeSmoothingOff();
- smoother->BoundarySmoothingOff();
- smoother->SetConvergence( 0 );
-
- polydata->Delete();//RC--
- polydata = smoother->GetOutput();
- polydata->Register(NULL);//RC++
- smoother->Delete();
- }
- ProgressBar::GetInstance()->Progress();
-
- //decimate = to reduce number of polygons
- if(m_Decimate==DecimatePro)
- {
- vtkDecimatePro *decimate = vtkDecimatePro::New();
- decimate->SplittingOff();
- decimate->SetErrorIsAbsolute(5);
- decimate->SetFeatureAngle(30);
- decimate->PreserveTopologyOn();
- decimate->BoundaryVertexDeletionOff();
- decimate->SetDegree(10); //std-value is 25!
-
- decimate->SetInput(polydata);//RC++
- decimate->SetTargetReduction(m_TargetReduction);
- decimate->SetMaximumError(0.002);
-
- polydata->Delete();//RC--
- polydata = decimate->GetOutput();
- polydata->Register(NULL);//RC++
- decimate->Delete();
- }
- else if (m_Decimate==QuadricDecimation)
- {
- vtkQuadricDecimation* decimate = vtkQuadricDecimation::New();
- decimate->SetTargetReduction(m_TargetReduction);
-
- decimate->SetInput(polydata);
- polydata->Delete();
- polydata = decimate->GetOutput();
- polydata->Register(NULL);
- decimate->Delete();
- }
-
- polydata->Update();
- ProgressBar::GetInstance()->Progress();
-
- polydata->SetSource(NULL);
-
- if(polydata->GetNumberOfPoints() > 0)
- {
- mitk::Vector3D spacing = GetInput()->GetGeometry(time)->GetSpacing();
-
- vtkPoints * points = polydata->GetPoints();
- vtkMatrix4x4 *vtkmatrix = vtkMatrix4x4::New();
- GetInput()->GetGeometry(time)->GetVtkTransform()->GetMatrix(vtkmatrix);
- double (*matrix)[4] = vtkmatrix->Element;
-
- unsigned int i,j;
- for(i=0;i<3;++i)
- for(j=0;j<3;++j)
- matrix[i][j]/=spacing[j];
-
- unsigned int n = points->GetNumberOfPoints();
- vtkFloatingPointType point[3];
-
- for (i = 0; i < n; i++)
- {
- points->GetPoint(i, point);
- mitkVtkLinearTransformPoint(matrix,point,point);
- points->SetPoint(i, point);
- }
- vtkmatrix->Delete();
- }
- ProgressBar::GetInstance()->Progress();
-
- // determine point_data normals for the poly data points.
- vtkSmartPointer<vtkPolyDataNormals> normalsGenerator = vtkSmartPointer<vtkPolyDataNormals>::New();
- normalsGenerator->SetInput( polydata );
-
- vtkSmartPointer<vtkCleanPolyData> cleanPolyDataFilter = vtkSmartPointer<vtkCleanPolyData>::New();
- cleanPolyDataFilter->SetInput(normalsGenerator->GetOutput());
- cleanPolyDataFilter->PieceInvariantOff();
- cleanPolyDataFilter->ConvertLinesToPointsOff();
- cleanPolyDataFilter->ConvertPolysToLinesOff();
- cleanPolyDataFilter->ConvertStripsToPolysOff();
- cleanPolyDataFilter->PointMergingOn();
- cleanPolyDataFilter->Update();
-
- surface->SetVtkPolyData(cleanPolyDataFilter->GetOutput(), time);
- polydata->UnRegister(NULL);
-}
-
-
-void mitk::ImageToSurfaceFilter::GenerateData()
-{
- mitk::Surface *surface = this->GetOutput();
- mitk::Image * image = (mitk::Image*)GetInput();
- if(image == NULL || !image->IsInitialized())
- mitkThrow() << "No input image set, please set an valid input image!";
-
- mitk::Image::RegionType outputRegion = image->GetRequestedRegion();
-
- int tstart=outputRegion.GetIndex(3);
- int tmax=tstart+outputRegion.GetSize(3); //GetSize()==1 - will aber 0 haben, wenn nicht zeitaufgeloest
-
- if ((tmax-tstart) > 0)
- {
- ProgressBar::GetInstance()->AddStepsToDo( 4 * (tmax - tstart) );
- }
-
-
- int t;
- for( t=tstart; t < tmax; ++t)
- {
- vtkImageData *vtkimagedata = image->GetVtkImageData(t);
- CreateSurface(t,vtkimagedata,surface,m_Threshold);
- ProgressBar::GetInstance()->Progress();
- }
-}
-
-void mitk::ImageToSurfaceFilter::SetSmoothIteration(int smoothIteration)
-{
- m_SmoothIteration = smoothIteration;
-}
-
-void mitk::ImageToSurfaceFilter::SetSmoothRelaxation(float smoothRelaxation)
-{
- m_SmoothRelaxation = smoothRelaxation;
-}
-
-void mitk::ImageToSurfaceFilter::SetInput(const mitk::Image *image)
-{
- // Process object is not const-correct so the const_cast is required here
- this->ProcessObject::SetNthInput(0, const_cast< mitk::Image * >( image ) );
-}
-
-
-const mitk::Image *mitk::ImageToSurfaceFilter::GetInput(void)
-{
- if (this->GetNumberOfInputs() < 1)
- {
- return 0;
- }
-
- return static_cast<const mitk::Image * >
- ( this->ProcessObject::GetInput(0) );
-}
-
-
-void mitk::ImageToSurfaceFilter::GenerateOutputInformation()
-{
- mitk::Image::ConstPointer inputImage =(mitk::Image*) this->GetInput();
-
- //mitk::Image *inputImage = (mitk::Image*)this->GetImage();
- mitk::Surface::Pointer output = this->GetOutput();
-
- itkDebugMacro(<<"GenerateOutputInformation()");
-
- if(inputImage.IsNull()) return;
-
- //Set Data
-}
diff --git a/Modules/MitkExt/Algorithms/mitkImageToSurfaceFilter.h b/Modules/MitkExt/Algorithms/mitkImageToSurfaceFilter.h
deleted file mode 100644
index 484ebc5..0000000
--- a/Modules/MitkExt/Algorithms/mitkImageToSurfaceFilter.h
+++ /dev/null
@@ -1,241 +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 _MITKIMAGETOSURFACEFILTER_h__
-#define _MITKIMAGETOSURFACEFILTER_h__
-
-#include <vtkPolyData.h>
-#include "MitkExtExports.h"
-#include <mitkCommon.h>
-#include <mitkSurfaceSource.h>
-#include <mitkSurface.h>
-
-#include <mitkImage.h>
-#include <vtkImageData.h>
-
-#include <vtkSmoothPolyDataFilter.h>
-#include <vtkMarchingCubes.h>
-
-
-namespace mitk {
- /**
- * @brief Converts pixel data to surface data by using a threshold
- * The mitkImageToSurfaceFilter is used to create a new surface out of an mitk image. The filter
- * uses a threshold to define the surface. It is based on the vtkMarchingCube algorithm. By default
- * a vtkPolyData surface based on an input threshold for the input image will be created. Optional
- * it is possible to reduce the number of triangles/polygones [SetDecimate(mitk::ImageToSurfaceFilter::DecimatePro) and SetTargetReduction (float _arg)]
- * or smooth the surface-data [SetSmooth(true), SetSmoothIteration(int smoothIteration) and SetSmoothRelaxation(float smoothRelaxation)].
- *
- * The resulting vtk-surface has the same size as the input image. The surface
- * can be generally smoothed by vtkDecimatePro reduce complexity of triangles
- * and vtkSmoothPolyDataFilter to relax the mesh. Both are enabled by default
- * and connected in the common way of pipelining in ITK. It's also possible
- * to create time sliced surfaces.
- *
- * @ingroup ImageFilters
- * @ingroup Process
- */
-
- class MitkExt_EXPORT ImageToSurfaceFilter : public SurfaceSource
- {
- public:
-
- /*
- * To decide whether a reduction of polygons in the created surface shall be
- * done or not by using the vtkDecimatePro Filter. Till vtk 4.x an vtkDecimateFilter existed,
- * but was patented. So since vtk 5.x it was replaced by the (worser?) vtkDecimateProFilter
- * Maybe another Filter will come soon.
- */
- enum DecimationType {NoDecimation,DecimatePro,QuadricDecimation};
-
- mitkClassMacro(ImageToSurfaceFilter, SurfaceSource);
- itkNewMacro(Self);
-
- /**
- * For each image time slice a surface will be created. This method is
- * called by Update().
- */
-
- virtual void GenerateData();
-
- /**
- * Initializes the output information ( i.e. the geometry information ) of
- * the output of the filter
- */
- virtual void GenerateOutputInformation();
-
- /**
- * Returns a const reference to the input image (e.g. the original input image that ist used to create the surface)
- */
- const mitk::Image *GetInput(void);
-
- /**
- * Set the source image to create a surface for this filter class. As input every mitk
- * 3D or 3D+t image can be used.
- */
- virtual void SetInput(const mitk::Image *image);
-
- /**
- * Set the number of iterations that is used to smooth the surface. Used is the vtkSmoothPolydataFilter that uses the laplacian filter. The higher the number of iterations that stronger the smooth-result
- *
- * @param smoothIteration As smoothIteration default in that case 50 was choosen. The VTK documentation recommends small relaxation factors and large numbers of iterations.
- */
- void SetSmoothIteration(int smoothIteration);
-
- /**
- * Set number of relaxation. Specify the relaxation factor for Laplacian
- * smoothing. The VTK documentation recommends small relaxation factors
- * and large numbers of iterations.
- *
- * @param smoothRelaxation As smoothRelaxation default in that case 0.1 was choosen. The VTK documentation recommends small relaxation factors and large numbers of iterations.
- */
- void SetSmoothRelaxation(float smoothRelaxation);
-
-
- /**
- * Threshold that is used to create the surface. All pixel in the input image that are higher than that
- * value will be considered in the surface. The threshold referees to
- * vtkMarchingCube. Default value is 1. See also SetThreshold (ScalarType _arg)
- */
- itkSetMacro(Threshold, ScalarType);
-
- /**
- * Get Threshold from vtkMarchingCube. Threshold can be manipulated by
- * inherited classes.
- */
- itkGetConstMacro(Threshold, ScalarType);
-
- /**
- * Enables vtkSmoothPolyDataFilter. With Laplacian smoothing this filter
- * will relax the surface. You can control the Filter by manipulating the
- * number of iterations and the relaxing factor.
- * */
- itkSetMacro(Smooth,bool);
-
- /*
- * Enable/Disable surface smoothing.
- */
- itkBooleanMacro(Smooth);
-
- /*
- * Returns if surface smoothing is enabled
- */
- itkGetConstMacro(Smooth,bool);
-
- /**
- * Get the state of decimation mode to reduce triangle in the
- * surface represantation. Modes can only be NoDecimation or DecimatePro
- * (till vtk 4.x also Decimate)
- * */
- itkGetConstMacro(Decimate,DecimationType);
-
- /**
- * Enable the decimation filter to reduce the number of triangles in the
- * mesh and produce a good approximation to the original image. The filter
- * has support for vtk-5 and earlier versions. More detailed information
- * check the vtkDecimatePro and vtkDecimate.
- * */
- itkSetMacro(Decimate,DecimationType);
-
- /**
- * Set desired TargetReduction of triangles in the range from 0.0 to 1.0.
- * The destroyed triangles are in relation with the size of data. For example 0.9
- * will reduce the data set to 10%.
- *
- * @param Set a TargetReduction float-value from 0.0 to 1.0
- * */
- itkSetMacro(TargetReduction, float);
-
- /**
- * Returns the reduction factor for the VtkDecimatePro Decimation Filter as a float value
- */
- itkGetConstMacro(TargetReduction, float);
-
- /**
- * Transforms a point by a 4x4 matrix
- */
- template <class T1, class T2, class T3>
- inline void mitkVtkLinearTransformPoint(T1 matrix[4][4], T2 in[3], T3 out[3])
- {
- T3 x = matrix[0][0]*in[0]+matrix[0][1]*in[1]+matrix[0][2]*in[2]+matrix[0][3];
- T3 y = matrix[1][0]*in[0]+matrix[1][1]*in[1]+matrix[1][2]*in[2]+matrix[1][3];
- T3 z = matrix[2][0]*in[0]+matrix[2][1]*in[1]+matrix[2][2]*in[2]+matrix[2][3];
- out[0] = x;
- out[1] = y;
- out[2] = z;
- }
-
- protected:
-
-
- ImageToSurfaceFilter();
-
- /**
- * Destructor
- * */
- virtual ~ImageToSurfaceFilter();
-
- /**
- * With the given threshold vtkMarchingCube creates the surface. By default nothing a
- * vtkPolyData surface based on a threshold of the input image will be created. Optional
- * it is possible to reduce the number of triangles/polygones [SetDecimate(mitk::ImageToSurfaceFilter::DecimatePro) and SetTargetReduction (float _arg)]
- * or smooth the data [SetSmooth(true), SetSmoothIteration(int smoothIteration) and SetSmoothRelaxation(float smoothRelaxation)].
- *
- * @param time selected slice or "0" for single
- * @param *vtkimage input image
- * @param *surface output
- * @param threshold can be different from SetThreshold()
- */
- void CreateSurface(int time, vtkImageData *vtkimage, mitk::Surface * surface, const ScalarType threshold);
-
- /**
- * Flag whether the created surface shall be smoothed or not (default is "false"). SetSmooth (bool _arg)
- * */
- bool m_Smooth;
-
- /**
- * Decimation mode, default mode is "NoDecimation". See also SetDecimate (DecimationType _arg)
- * */
- DecimationType m_Decimate;
-
- /**
- * Threshold that is used to create the surface. All pixel in the input image that are higher than that
- * value will be considered in the surface. Default value is 1. See also SetThreshold (ScalarType _arg)
- * */
- ScalarType m_Threshold;
-
- /**
- * The Reduction factor of the Decimation Filter for the created surface. See also SetTargetReduction (float _arg)
- * */
- float m_TargetReduction;
-
- /**
- * The Iteration value for the Smooth Filter of the created surface. See also SetSmoothIteration (int smoothIteration)
- * */
- int m_SmoothIteration;
-
- /**
- * The Relaxation value for the Smooth Filter of the created surface. See also SetSmoothRelaxation (float smoothRelaxation)
- * */
- float m_SmoothRelaxation;
-
- };
-
-} // namespace mitk
-
-#endif //_MITKIMAGETOSURFACEFILTER_h__
-
-
diff --git a/Modules/MitkExt/Algorithms/mitkSurfaceToImageFilter.cpp b/Modules/MitkExt/Algorithms/mitkSurfaceToImageFilter.cpp
deleted file mode 100644
index acf52c4..0000000
--- a/Modules/MitkExt/Algorithms/mitkSurfaceToImageFilter.cpp
+++ /dev/null
@@ -1,209 +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.
-
-===================================================================*/
-
-#include "mitkSurfaceToImageFilter.h"
-#include "mitkTimeHelper.h"
-
-#include <vtkPolyData.h>
-#include <vtkPolyDataToImageStencil.h>
-#include <vtkImageStencil.h>
-#include <vtkImageData.h>
-#include <vtkPolyData.h>
-#include <vtkLinearTransform.h>
-#include <vtkTriangleFilter.h>
-#include <vtkLinearExtrusionFilter.h>
-#include <vtkDataSetTriangleFilter.h>
-#include <vtkImageThreshold.h>
-#include <vtkImageMathematics.h>
-#include <vtkPolyDataNormals.h>
-#include <vtkTransformPolyDataFilter.h>
-#include <vtkTransform.h>
-
-
-mitk::SurfaceToImageFilter::SurfaceToImageFilter()
-: m_MakeOutputBinary( false ),
- m_BackgroundValue( -10000 )
-{
-}
-
-mitk::SurfaceToImageFilter::~SurfaceToImageFilter()
-{
-}
-
-void mitk::SurfaceToImageFilter::GenerateInputRequestedRegion()
-{
- mitk::Image* output = this->GetOutput();
- if((output->IsInitialized()==false) )
- return;
-
- GenerateTimeInInputRegion(output, const_cast< mitk::Image * > ( this->GetImage() ));
-}
-
-void mitk::SurfaceToImageFilter::GenerateOutputInformation()
-{
- mitk::Image *inputImage = (mitk::Image*)this->GetImage();
- mitk::Image::Pointer output = this->GetOutput();
-
- itkDebugMacro(<<"GenerateOutputInformation()");
-
- if((inputImage == NULL) ||
- (inputImage->IsInitialized() == false) ||
- (inputImage->GetTimeSlicedGeometry() == NULL)) return;
-
- if (m_MakeOutputBinary)
- output->Initialize(mitk::MakeScalarPixelType<unsigned char>() , *inputImage->GetTimeSlicedGeometry());
- else
- output->Initialize(inputImage->GetPixelType(), *inputImage->GetTimeSlicedGeometry());
-
- output->SetPropertyList(inputImage->GetPropertyList()->Clone());
-}
-
-void mitk::SurfaceToImageFilter::GenerateData()
-{
- mitk::Image::ConstPointer inputImage = this->GetImage();
- mitk::Image::Pointer output = this->GetOutput();
-
- if(inputImage.IsNull())
- return;
-
- if(output->IsInitialized()==false )
- return;
-
- mitk::Image::RegionType outputRegion = output->GetRequestedRegion();
-
- int tstart=outputRegion.GetIndex(3);
- int tmax=tstart+outputRegion.GetSize(3);
-
- if ( tmax > 0)
- {
- int t;
- for(t=tstart;t<tmax;++t)
- {
- Stencil3DImage( t );
- }
- }
- else
- {
- Stencil3DImage( 0 );
- }
-}
-
-void mitk::SurfaceToImageFilter::Stencil3DImage(int time)
-{
- const mitk::TimeSlicedGeometry *surfaceTimeGeometry = GetInput()->GetTimeSlicedGeometry();
- const mitk::TimeSlicedGeometry *imageTimeGeometry = GetImage()->GetTimeSlicedGeometry();
-
- // Convert time step from image time-frame to surface time-frame
- int surfaceTimeStep = surfaceTimeGeometry->TimeStepToTimeStep( imageTimeGeometry, time );
-
- vtkPolyData * polydata = ( (mitk::Surface*)GetInput() )->GetVtkPolyData( surfaceTimeStep );
-
- vtkTransformPolyDataFilter * move=vtkTransformPolyDataFilter::New();
- move->SetInput(polydata);
- move->ReleaseDataFlagOn();
-
- vtkTransform *transform=vtkTransform::New();
- Geometry3D* geometry = surfaceTimeGeometry->GetGeometry3D( surfaceTimeStep );
- geometry->TransferItkToVtkTransform();
- transform->PostMultiply();
- transform->Concatenate(geometry->GetVtkTransform()->GetMatrix());
- // take image geometry into account. vtk-Image information will be changed to unit spacing and zero origin below.
- Geometry3D* imageGeometry = imageTimeGeometry->GetGeometry3D(time);
- imageGeometry->TransferItkToVtkTransform();
- transform->Concatenate(imageGeometry->GetVtkTransform()->GetLinearInverse());
- move->SetTransform(transform);
- transform->Delete();
-
- vtkPolyDataNormals * normalsFilter = vtkPolyDataNormals::New();
- normalsFilter->SetFeatureAngle(50);
- normalsFilter->SetConsistency(1);
- normalsFilter->SetSplitting(1);
- normalsFilter->SetFlipNormals(0);
- normalsFilter->ReleaseDataFlagOn();
-
- normalsFilter->SetInput( move->GetOutput() );
- move->Delete();
-
- vtkPolyDataToImageStencil * surfaceConverter = vtkPolyDataToImageStencil::New();
- surfaceConverter->SetTolerance( 0.0 );
- surfaceConverter->ReleaseDataFlagOn();
-
- surfaceConverter->SetInput( normalsFilter->GetOutput() );
- normalsFilter->Delete();
-
- mitk::Image::Pointer binaryImage = mitk::Image::New();
-
- if (m_MakeOutputBinary)
- {
- binaryImage->Initialize(mitk::MakeScalarPixelType<unsigned char>(), *this->GetImage()->GetTimeSlicedGeometry());
-
- unsigned int size = sizeof(unsigned char);
- for (unsigned int i = 0; i < binaryImage->GetDimension(); ++i)
- size *= binaryImage->GetDimension(i);
-
- memset(binaryImage->GetData(), 1, size);
- }
-
- vtkImageData *image = m_MakeOutputBinary
- ? binaryImage->GetVtkImageData(time)
- : const_cast<mitk::Image *>(this->GetImage())->GetVtkImageData(time);
-
- // Create stencil and use numerical minimum of pixel type as background value
- vtkImageStencil *stencil = vtkImageStencil::New();
- stencil->SetInput(image);
- stencil->ReverseStencilOff();
- stencil->ReleaseDataFlagOn();
- stencil->SetStencil(surfaceConverter->GetOutput());
- surfaceConverter->Delete();
-
- stencil->SetBackgroundValue(m_MakeOutputBinary ? 0 : m_BackgroundValue);
- stencil->Update();
-
- mitk::Image::Pointer output = this->GetOutput();
- output->SetVolume( stencil->GetOutput()->GetScalarPointer(), time );
- MITK_INFO << "stencil ref count: " << stencil->GetReferenceCount() << std::endl;
-
- stencil->Delete();
-}
-
-
-const mitk::Surface *mitk::SurfaceToImageFilter::GetInput(void)
-{
- if (this->GetNumberOfInputs() < 1)
- {
- return 0;
- }
-
- return static_cast<const mitk::Surface * >
- ( this->ProcessObject::GetInput(0) );
-}
-
-void mitk::SurfaceToImageFilter::SetInput(const mitk::Surface *input)
-{
- // Process object is not const-correct so the const_cast is required here
- this->ProcessObject::SetNthInput(0,
- const_cast< mitk::Surface * >( input ) );
-}
-
-void mitk::SurfaceToImageFilter::SetImage(const mitk::Image *source)
-{
- this->ProcessObject::SetNthInput( 1, const_cast< mitk::Image * >( source ) );
-}
-
-const mitk::Image *mitk::SurfaceToImageFilter::GetImage(void)
-{
- return static_cast< const mitk::Image * >(this->ProcessObject::GetInput(1));
-}
diff --git a/Modules/MitkExt/Algorithms/mitkSurfaceToImageFilter.h b/Modules/MitkExt/Algorithms/mitkSurfaceToImageFilter.h
deleted file mode 100644
index 591ac89..0000000
--- a/Modules/MitkExt/Algorithms/mitkSurfaceToImageFilter.h
+++ /dev/null
@@ -1,96 +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 _mitkSurfaceToImageFilter_h__
-#define _mitkSurfaceToImageFilter_h__
-
-#include "mitkCommon.h"
-#include "MitkExtExports.h"
-#include "mitkImageSource.h"
-#include "mitkSurface.h"
-//#include "mitkImage.h"
-
-class vtkPolyData;
-
-namespace mitk {
-
-//class Mesh;
-//class VectorOfContourLines;
-
-/**
- *
- * @brief Converts surface data to pixel data. Requires a surface and an
- * image, which header information defines the output image.
- *
- * The resulting image has the same dimension, size, and Geometry3D
- * as the input image. The image is cut using a vtkStencil.
- * The user can decide if he wants to keep the original values or create a
- * binary image by setting MakeBinaryOutputOn (default is \a false). If
- * set to \a true all voxels inside the surface are set to one and all
- * outside voxel are set to zero.
- *
- * NOTE: Since the reference input image is passed to the vtkStencil in
- * any case, the image needs to be initialized with pixel values greater than
- * the numerical minimum of the used pixel type (e.g. at least -127 for
- * unsigned char images, etc.) to produce a correct binary image
- * representation of the surface in MakeOutputBinary mode.
- *
- * @ingroup SurfaceFilters
- * @ingroup Process
- */
-class MitkExt_EXPORT SurfaceToImageFilter : public ImageSource
-{
-public:
- mitkClassMacro(SurfaceToImageFilter, ImageSource);
- itkNewMacro(Self);
-
- itkSetMacro(MakeOutputBinary, bool);
- itkGetMacro(MakeOutputBinary, bool);
- itkBooleanMacro(MakeOutputBinary);
-
- itkGetConstMacro(BackgroundValue,float);
- itkSetMacro(BackgroundValue,float);
-
- virtual void GenerateInputRequestedRegion();
-
- virtual void GenerateOutputInformation();
-
- virtual void GenerateData();
-
- const mitk::Surface *GetInput(void);
-
- virtual void SetInput(const mitk::Surface *surface);
-
- void SetImage(const mitk::Image *source);
-
- const mitk::Image *GetImage(void);
-
-protected:
- SurfaceToImageFilter();
-
- virtual ~SurfaceToImageFilter();
-
- void Stencil3DImage(int time = 0);
-
- bool m_MakeOutputBinary;
-
- float m_BackgroundValue;
-
-};
-
-} // namespace mitk
-
-#endif /* MITKCOONSPATCHFILTER_H_HEADER_INCLUDED_C10B22CD */
diff --git a/Modules/MitkExt/Algorithms/mitkTimeHelper.h b/Modules/MitkExt/Algorithms/mitkTimeHelper.h
deleted file mode 100644
index 5313018..0000000
--- a/Modules/MitkExt/Algorithms/mitkTimeHelper.h
+++ /dev/null
@@ -1,81 +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 MITKTIMEHELPER_H_HEADER_INCLUDED_C1C2FCD2
-#define MITKTIMEHELPER_H_HEADER_INCLUDED_C1C2FCD2
-
-namespace mitk
-{
-
-//## @brief convert the start- and end-index-time of output-region in
-//## start- and end-index-time of input-region via millisecond-time
-template <class TOutputRegion, class TInputRegion>
-void ITK_EXPORT GenerateTimeInInputRegion(const mitk::TimeSlicedGeometry *outputTimeGeometry, const TOutputRegion& outputRegion, const mitk::TimeSlicedGeometry *inputTimeGeometry, TInputRegion& inputRegion)
-{
- assert(outputTimeGeometry!=NULL);
- assert(inputTimeGeometry!=NULL);
-
- // convert the start-index-time of output in start-index-time of input via millisecond-time
- ScalarType timeInMS = outputTimeGeometry->TimeStepToMS(outputRegion.GetIndex(3));
- int timestep = inputTimeGeometry->MSToTimeStep( timeInMS );
- if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( inputTimeGeometry->IsValidTime( timestep ) ) )
- inputRegion.SetIndex( 3, timestep );
- else
- inputRegion.SetIndex( 3, 0 );
- // convert the end-index-time of output in end-index-time of input via millisecond-time
- timeInMS = outputTimeGeometry->TimeStepToMS(outputRegion.GetIndex(3)+outputRegion.GetSize(3)-1);
- timestep = inputTimeGeometry->MSToTimeStep( timeInMS );
- if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( outputTimeGeometry->IsValidTime( timestep ) ) )
- inputRegion.SetSize( 3, timestep - inputRegion.GetIndex(3) + 1 );
- else
- inputRegion.SetSize( 3, 1 );
-}
-
-//##Documentation
-//## @brief convert the start- and end-index-time of output in
-//## start- and end-index-time of input1 and input2 via millisecond-time
-template <class TOutputData, class TInputData>
-void ITK_EXPORT GenerateTimeInInputRegion(const TOutputData* output, TInputData* input)
-{
- assert(output!=NULL);
- assert(input!=NULL);
-
- const typename TOutputData::RegionType& outputRegion = output->GetRequestedRegion();
- typename TInputData::RegionType inputRegion;
-
- if(outputRegion.GetSize(3)<1)
- {
- typename TInputData::RegionType::SizeType inputsize;
- inputsize.Fill(0);
- inputRegion.SetSize(inputsize);
- input->SetRequestedRegion( &inputRegion );
- }
-
- // convert the start-index-time of output in start-index-time of input via millisecond-time
- inputRegion = input->GetRequestedRegion();
- GenerateTimeInInputRegion(output->GetTimeSlicedGeometry(), outputRegion, input->GetTimeSlicedGeometry(), inputRegion);
- input->SetRequestedRegion( &inputRegion );
-}
-
-} // end namespace mitk
-
-//#ifndef ITK_MANUAL_INSTANTIATION
-//#include "mitkTimeHelper.txx"
-#include "MitkExtExports.h"
-//#endif
-
-#endif // MITKTIMEHELPER_H_HEADER_INCLUDED_C1C2FCD2
diff --git a/Modules/MitkExt/files.cmake b/Modules/MitkExt/files.cmake
index 4061c11..6092fac 100644
--- a/Modules/MitkExt/files.cmake
+++ b/Modules/MitkExt/files.cmake
@@ -4,7 +4,6 @@ set(CPP_FILES
Algorithms/vtkPointSetSlicer.cxx
Algorithms/mitkCoreExtObjectFactory.cpp
Algorithms/mitkAngleCorrectByPointFilter.cpp
- Algorithms/mitkAutoCropImageFilter.cpp
Algorithms/mitkBoundingObjectCutter.cpp
Algorithms/mitkCylindricToCartesianFilter.cpp
Algorithms/mitkDopplerToStrainRateFilter.cpp
@@ -12,7 +11,6 @@ set(CPP_FILES
Algorithms/mitkGeometryDataSource.cpp
Algorithms/mitkHeightFieldSurfaceClipImageFilter.cpp
Algorithms/mitkImageToLookupTableFilter.cpp
- Algorithms/mitkImageToSurfaceFilter.cpp
Algorithms/mitkInterpolateLinesFilter.cpp
Algorithms/mitkLabeledImageToSurfaceFilter.cpp
Algorithms/mitkLabeledImageVolumeCalculator.cpp
@@ -33,7 +31,6 @@ set(CPP_FILES
Algorithms/mitkProbeFilter.cpp
Algorithms/mitkSimpleHistogram.cpp
Algorithms/mitkSimpleUnstructuredGridHistogram.cpp
- Algorithms/mitkSurfaceToImageFilter.cpp
Algorithms/mitkUnstructuredGridHistogram.cpp
Algorithms/mitkUnstructuredGridSource.cpp
Algorithms/mitkVolumeVisualizationImagePreprocessor.cpp
diff --git a/Modules/Segmentation/Algorithms/mitkComputeContourSetNormalsFilter.cpp b/Modules/Segmentation/Algorithms/mitkComputeContourSetNormalsFilter.cpp
deleted file mode 100644
index 9797416..0000000
--- a/Modules/Segmentation/Algorithms/mitkComputeContourSetNormalsFilter.cpp
+++ /dev/null
@@ -1,305 +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.
-
-===================================================================*/
-
-#include "mitkComputeContourSetNormalsFilter.h"
-
-
-mitk::ComputeContourSetNormalsFilter::ComputeContourSetNormalsFilter()
-{
- m_MaxSpacing = 5;
- this->m_UseProgressBar = false;
- this->m_ProgressStepSize = 1;
-}
-
-mitk::ComputeContourSetNormalsFilter::~ComputeContourSetNormalsFilter()
-{
-}
-
-void mitk::ComputeContourSetNormalsFilter::GenerateData()
-{
- unsigned int numberOfInputs = this->GetNumberOfInputs();
- this->CreateOutputsForAllInputs(numberOfInputs);
-
- //Iterating over each input
- for(unsigned int i = 0; i < numberOfInputs; i++)
- {
- //Getting the inputs polydata and polygons
- Surface* currentSurface = const_cast<Surface*>( this->GetInput(i) );
- vtkPolyData* polyData = currentSurface->GetVtkPolyData();
-
- vtkSmartPointer<vtkCellArray> existingPolys = polyData->GetPolys();
-
- vtkSmartPointer<vtkPoints> existingPoints = polyData->GetPoints();
-
- existingPolys->InitTraversal();
-
- vtkIdType* cell (NULL);
- vtkIdType cellSize (0);
-
- //The array that contains all the vertex normals of the current polygon
- vtkSmartPointer<vtkDoubleArray> normals = vtkSmartPointer<vtkDoubleArray>::New();
- normals->SetNumberOfComponents(3);
- normals->SetNumberOfTuples(polyData->GetNumberOfPoints());
-
- //If the current contour is an inner contour then the direction is -1
- //A contour lies inside another one if the pixel values in the direction of the normal is 1
- m_NegativeNormalCounter = 0;
- m_PositiveNormalCounter = 0;
-
- //Iterating over each polygon
- for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
- {
- if(cellSize < 3)continue;
-
- //First we calculate the current polygon's normal
- double polygonNormal[3] = {0.0};
-
- double p1[3];
- double p2[3];
-
- double v1[3];
- double v2[3];
-
- existingPoints->GetPoint(cell[0], p1);
- unsigned int index = cellSize*0.5;
- existingPoints->GetPoint(cell[index], p2);
-
- v1[0] = p2[0]-p1[0];
- v1[1] = p2[1]-p1[1];
- v1[2] = p2[2]-p1[2];
-
- for (unsigned int k = 2; k < cellSize; k++)
- {
- index = cellSize*0.25;
- existingPoints->GetPoint(cell[index], p1);
- index = cellSize*0.75;
- existingPoints->GetPoint(cell[index], p2);
-
- v2[0] = p2[0]-p1[0];
- v2[1] = p2[1]-p1[1];
- v2[2] = p2[2]-p1[2];
-
- vtkMath::Cross(v1,v2,polygonNormal);
- if (vtkMath::Norm(polygonNormal) != 0)
- break;
- }
-
- vtkMath::Normalize(polygonNormal);
-
- //Now we start computing the normal for each vertex
-
- double vertexNormalTemp[3];
- existingPoints->GetPoint(cell[0], p1);
- existingPoints->GetPoint(cell[1], p2);
-
- v1[0] = p2[0]-p1[0];
- v1[1] = p2[1]-p1[1];
- v1[2] = p2[2]-p1[2];
-
- vtkMath::Cross(v1,polygonNormal,vertexNormalTemp);
-
- vtkMath::Normalize(vertexNormalTemp);
-
- double vertexNormal[3];
- for (unsigned j = 0; j < cellSize-2; j++)
- {
- existingPoints->GetPoint(cell[j+1], p1);
- existingPoints->GetPoint(cell[j+2], p2);
-
- v1[0] = p2[0]-p1[0];
- v1[1] = p2[1]-p1[1];
- v1[2] = p2[2]-p1[2];
-
- vtkMath::Cross(v1,polygonNormal,vertexNormal);
-
- vtkMath::Normalize(vertexNormal);
-
- double finalNormal[3];
-
- finalNormal[0] = (vertexNormal[0] + vertexNormalTemp[0])*0.5;
- finalNormal[1] = (vertexNormal[1] + vertexNormalTemp[1])*0.5;
- finalNormal[2] = (vertexNormal[2] + vertexNormalTemp[2])*0.5;
-
- //Here we determine the direction of the normal
- if (j == 0 && m_SegmentationBinaryImage)
- {
- Point3D worldCoord;
- worldCoord[0] = p1[0]+finalNormal[0]*m_MaxSpacing;
- worldCoord[1] = p1[1]+finalNormal[1]*m_MaxSpacing;
- worldCoord[2] = p1[2]+finalNormal[2]*m_MaxSpacing;
-
- double val = m_SegmentationBinaryImage->GetPixelValueByWorldCoordinate(worldCoord);
- if (val == 1.0)
- {
- ++m_PositiveNormalCounter;
- }
- else
- {
- ++m_NegativeNormalCounter;
- }
- }
-
- vertexNormalTemp[0] = vertexNormal[0];
- vertexNormalTemp[1] = vertexNormal[1];
- vertexNormalTemp[2] = vertexNormal[2];
-
- vtkIdType id = cell[j+1];
- normals->SetTuple(id,finalNormal);
- }
-
- existingPoints->GetPoint(cell[0], p1);
- existingPoints->GetPoint(cell[1], p2);
-
- v1[0] = p2[0]-p1[0];
- v1[1] = p2[1]-p1[1];
- v1[2] = p2[2]-p1[2];
-
- vtkMath::Cross(v1,polygonNormal,vertexNormal);
-
- vtkMath::Normalize(vertexNormal);
-
- vertexNormal[0] = (vertexNormal[0] + vertexNormalTemp[0])*0.5;
- vertexNormal[1] = (vertexNormal[1] + vertexNormalTemp[1])*0.5;
- vertexNormal[2] = (vertexNormal[2] + vertexNormalTemp[2])*0.5;
-
- vtkIdType id = cell[0];
- normals->SetTuple(id,vertexNormal);
- id = cell[cellSize-1];
- normals->SetTuple(id,vertexNormal);
-
- int normalDirection(-1);
-
- if(m_NegativeNormalCounter < m_PositiveNormalCounter)
- {
- normalDirection = 1;
- }
-
- for(unsigned int n = 0; n < normals->GetNumberOfTuples(); n++)
- {
- double normal[3];
- normals->GetTuple(n,normal);
- normal[0] = normalDirection*normal[0];
- normal[1] = normalDirection*normal[1];
- normal[2] = normalDirection*normal[2];
- }
-
-
- }//end for all cells
-
- Surface::Pointer surface = this->GetOutput(i);
- surface->GetVtkPolyData()->GetCellData()->SetNormals(normals);
- }//end for all inputs
-
- //Setting progressbar
- if (this->m_UseProgressBar)
- mitk::ProgressBar::GetInstance()->Progress(this->m_ProgressStepSize);
-}
-
-
-mitk::Surface::Pointer mitk::ComputeContourSetNormalsFilter::GetNormalsAsSurface()
-{
- //Just for debugging:
- vtkSmartPointer<vtkPolyData> newPolyData = vtkSmartPointer<vtkPolyData>::New();
- vtkSmartPointer<vtkCellArray> newLines = vtkSmartPointer<vtkCellArray>::New();
- vtkSmartPointer<vtkPoints> newPoints = vtkSmartPointer<vtkPoints>::New();
- unsigned int idCounter (0);
- //Debug end
-
- for (unsigned int i = 0; i < this->GetNumberOfOutputs(); i++)
- {
- Surface* currentSurface = const_cast<Surface*>( this->GetOutput(i) );
- vtkPolyData* polyData = currentSurface->GetVtkPolyData();
-
- vtkSmartPointer<vtkDoubleArray> currentCellNormals = vtkDoubleArray::SafeDownCast(polyData->GetCellData()->GetNormals());
-
- vtkSmartPointer<vtkCellArray> existingPolys = polyData->GetPolys();
-
- vtkSmartPointer<vtkPoints> existingPoints = polyData->GetPoints();
-
- existingPolys->InitTraversal();
-
- vtkIdType* cell (NULL);
- vtkIdType cellSize (0);
-
- for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
- {
- for ( unsigned int j = 0; j < cellSize; j++ )
- {
- double currentNormal[3];
- currentCellNormals->GetTuple(cell[j], currentNormal);
- vtkSmartPointer<vtkLine> line = vtkSmartPointer<vtkLine>::New();
- line->GetPointIds()->SetNumberOfIds(2);
- double newPoint[3];
- double p0[3];
- existingPoints->GetPoint(cell[j], p0);
- newPoint[0] = p0[0] + currentNormal[0];
- newPoint[1] = p0[1] + currentNormal[1];
- newPoint[2] = p0[2] + currentNormal[2];
-
- line->GetPointIds()->SetId(0, idCounter);
- newPoints->InsertPoint(idCounter, p0);
- idCounter++;
-
- line->GetPointIds()->SetId(1, idCounter);
- newPoints->InsertPoint(idCounter, newPoint);
- idCounter++;
-
- newLines->InsertNextCell(line);
- }//end for all points
- }//end for all cells
- }//end for all outputs
-
- newPolyData->SetPoints(newPoints);
- newPolyData->SetLines(newLines);
- newPolyData->BuildCells();
-
-
- mitk::Surface::Pointer surface = mitk::Surface::New();
- surface->SetVtkPolyData(newPolyData);
-
- return surface;
-
-}
-
-void mitk::ComputeContourSetNormalsFilter::SetMaxSpacing(double maxSpacing)
-{
- m_MaxSpacing = maxSpacing;
-}
-
-void mitk::ComputeContourSetNormalsFilter::GenerateOutputInformation()
-{
- Superclass::GenerateOutputInformation();
-}
-
-void mitk::ComputeContourSetNormalsFilter::Reset()
-{
- for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
- {
- this->PopBackInput();
- }
- this->SetNumberOfInputs(0);
- this->SetNumberOfOutputs(0);
-}
-
-void mitk::ComputeContourSetNormalsFilter::SetUseProgressBar(bool status)
-{
- this->m_UseProgressBar = status;
-}
-
-void mitk::ComputeContourSetNormalsFilter::SetProgressStepSize(unsigned int stepSize)
-{
- this->m_ProgressStepSize = stepSize;
-}
diff --git a/Modules/Segmentation/Algorithms/mitkComputeContourSetNormalsFilter.h b/Modules/Segmentation/Algorithms/mitkComputeContourSetNormalsFilter.h
deleted file mode 100644
index a0c454b..0000000
--- a/Modules/Segmentation/Algorithms/mitkComputeContourSetNormalsFilter.h
+++ /dev/null
@@ -1,108 +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 mitkComputeContourSetNormalsFilter_h_Included
-#define mitkComputeContourSetNormalsFilter_h_Included
-
-#include "SegmentationExports.h"
-#include "mitkSurfaceToSurfaceFilter.h"
-#include "mitkProgressBar.h"
-#include "mitkSurface.h"
-
-#include "vtkCellArray.h"
-#include "vtkPolyData.h"
-#include "vtkSmartPointer.h"
-#include "vtkDoubleArray.h"
-#include "vtkMath.h"
-#include "vtkCellData.h"
-#include "vtkLine.h"
-
-#include "mitkImage.h"
-
-namespace mitk {
-
- /**
- \brief Filter to compute the normales for contours based on vtkPolygons
-
-
-
- This filter takes a number of extracted contours and computes the normals for each
- contour edge point. The normals can be accessed by calling:
-
- filter->GetOutput(i)->GetVtkPolyData()->GetCellData()->GetNormals();
-
- See also the method GetNormalsAsSurface()
-
- Note: If a segmentation binary image is provided this filter assures that the computed normals
- do not point into the segmentation image
-
- $Author: fetzer$
-*/
-class Segmentation_EXPORT ComputeContourSetNormalsFilter : public SurfaceToSurfaceFilter
-{
-public:
-
- mitkClassMacro(ComputeContourSetNormalsFilter,SurfaceToSurfaceFilter);
- itkNewMacro(Self);
-
- itkSetMacro(SegmentationBinaryImage, mitk::Image::Pointer);
-
- /*
- \brief Returns the computed normals as a surface
- */
- mitk::Surface::Pointer GetNormalsAsSurface();
-
- //Resets the filter, i.e. removes all inputs and outputs
- void Reset();
-
- void SetMaxSpacing(double);
-
- /**
- \brief Set whether the mitkProgressBar should be used
-
- \a Parameter true for using the progress bar, false otherwise
- */
- void SetUseProgressBar(bool);
-
- /**
- \brief Set the stepsize which the progress bar should proceed
-
- \a Parameter The stepsize for progressing
- */
- void SetProgressStepSize(unsigned int stepSize);
-
-protected:
- ComputeContourSetNormalsFilter();
- virtual ~ComputeContourSetNormalsFilter();
- virtual void GenerateData();
- virtual void GenerateOutputInformation();
-
-private:
-
- //The segmentation out of which the contours were extracted. Can be used to determine the direction of the normals
- mitk::Image::Pointer m_SegmentationBinaryImage;
- double m_MaxSpacing;
-
- unsigned int m_NegativeNormalCounter;
- unsigned int m_PositiveNormalCounter;
-
- bool m_UseProgressBar;
- unsigned int m_ProgressStepSize;
-
-};//class
-
-}//namespace
-#endif
diff --git a/Modules/Segmentation/Algorithms/mitkCreateDistanceImageFromSurfaceFilter.cpp b/Modules/Segmentation/Algorithms/mitkCreateDistanceImageFromSurfaceFilter.cpp
deleted file mode 100644
index bd29ae0..0000000
--- a/Modules/Segmentation/Algorithms/mitkCreateDistanceImageFromSurfaceFilter.cpp
+++ /dev/null
@@ -1,509 +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.
-
-===================================================================*/
-
-#include "mitkCreateDistanceImageFromSurfaceFilter.h"
-
-
-mitk::CreateDistanceImageFromSurfaceFilter::CreateDistanceImageFromSurfaceFilter()
-{
- m_DistanceImageVolume = 50000;
- this->m_UseProgressBar = false;
- this->m_ProgressStepSize = 5;
-}
-
-
-mitk::CreateDistanceImageFromSurfaceFilter::~CreateDistanceImageFromSurfaceFilter()
-{
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::GenerateData()
-{
- //First of all we have to build the equation-system from the existing contour-edge-points
- this->CreateSolutionMatrixAndFunctionValues();
-
- //Then we solve the equation-system via QR - decomposition. The interpolation weights are obtained in that way
- vnl_qr<double> solver (m_SolutionMatrix);
- m_Weights = solver.solve(m_FunctionValues);
-
- //Setting progressbar
- if (this->m_UseProgressBar)
- mitk::ProgressBar::GetInstance()->Progress(2);
-
- //The last step is to create the distance map with the interpolated distance function
- this->CreateDistanceImage();
- m_Centers.clear();
- m_FunctionValues.clear();
- m_Normals.clear();
- m_Weights.clear();
- m_SolutionMatrix.clear();
-
- //Setting progressbar
- if (this->m_UseProgressBar)
- mitk::ProgressBar::GetInstance()->Progress(3);
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::CreateSolutionMatrixAndFunctionValues()
-{
- unsigned int numberOfInputs = this->GetNumberOfInputs();
-
- if (numberOfInputs == 0)
- {
- MITK_ERROR << "mitk::CreateDistanceImageFromSurfaceFilter: No input available. Please set an input!" << std::endl;
- itkExceptionMacro("mitk::CreateDistanceImageFromSurfaceFilter: No input available. Please set an input!");
- return;
- }
-
- //First of all we have to extract the nomals and the surface points.
- //Duplicated points can be eliminated
-
- Surface* currentSurface;
- vtkSmartPointer<vtkPolyData> polyData;
- vtkSmartPointer<vtkDoubleArray> currentCellNormals;
- vtkSmartPointer<vtkCellArray> existingPolys;
- vtkSmartPointer<vtkPoints> existingPoints;
-
- double p[3];
- PointType currentPoint;
- PointType normal;
-
- for (unsigned int i = 0; i < numberOfInputs; i++)
- {
- currentSurface = const_cast<Surface*>( this->GetInput(i) );
- polyData = currentSurface->GetVtkPolyData();
-
- if (polyData->GetNumberOfPolys() == 0)
- {
- MITK_INFO << "mitk::CreateDistanceImageFromSurfaceFilter: No input-polygons available. Please be sure the input surface consists of polygons!" << std::endl;
- }
-
- currentCellNormals = vtkDoubleArray::SafeDownCast(polyData->GetCellData()->GetNormals());
-
- existingPolys = polyData->GetPolys();
-
- existingPoints = polyData->GetPoints();
-
- existingPolys->InitTraversal();
-
- vtkIdType* cell (NULL);
- vtkIdType cellSize (0);
-
- for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
- {
- for ( unsigned int j = 0; j < cellSize; j++ )
- {
- existingPoints->GetPoint(cell[j], p);
-
- currentPoint.copy_in(p);
-
- int count = std::count(m_Centers.begin() ,m_Centers.end(),currentPoint);
-
- if (count == 0)
- {
- double currentNormal[3];
- currentCellNormals->GetTuple(cell[j], currentNormal);
-
- normal.copy_in(currentNormal);
-
- m_Normals.push_back(normal);
-
- m_Centers.push_back(currentPoint);
- }
-
- }//end for all points
- }//end for all cells
- }//end for all outputs
-
- //For we can now calculate the exact size of the centers we initialize the data structures
- unsigned int numberOfCenters = m_Centers.size();
- m_Centers.reserve(numberOfCenters*3);
-
- m_FunctionValues.set_size(numberOfCenters*3);
- m_FunctionValues.fill(0);
-
- //Create inner points
- for (unsigned int i = 0; i < numberOfCenters; i++)
- {
- currentPoint = m_Centers.at(i);
- normal = m_Normals.at(i);
-
- currentPoint[0] = currentPoint[0] - normal[0];
- currentPoint[1] = currentPoint[1] - normal[1];
- currentPoint[2] = currentPoint[2] - normal[2];
-
- m_Centers.push_back(currentPoint);
-
- m_FunctionValues.put(numberOfCenters+i, -1);
- }
-
- //Create outer points
- for (unsigned int i = 0; i < numberOfCenters; i++)
- {
- currentPoint = m_Centers.at(i);
- normal = m_Normals.at(i);
-
- currentPoint[0] = currentPoint[0] + normal[0];
- currentPoint[1] = currentPoint[1] + normal[1];
- currentPoint[2] = currentPoint[2] + normal[2];
-
- m_Centers.push_back(currentPoint);
-
- m_FunctionValues.put(numberOfCenters*2+i, 1);
- }
-
- //Now we have created all centers and all function values. Next step is to create the solution matrix
- numberOfCenters = m_Centers.size();
- m_SolutionMatrix.set_size(numberOfCenters, numberOfCenters);
-
- m_Weights.set_size(numberOfCenters);
-
- PointType p1;
- PointType p2;
- double norm;
-
- for (unsigned int i = 0; i < numberOfCenters; i++)
- {
- for (unsigned int j = 0; j < numberOfCenters; j++)
- {
- //Calculate the RBF value. Currently using Phi(r) = r with r is the euclidian distance between two points
- p1 = m_Centers.at(i);
- p2 = m_Centers.at(j);
- p1 = p1 - p2;
- norm = p1.two_norm();
- m_SolutionMatrix(i,j) = norm;
-
- }
- }
-
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::CreateDistanceImage()
-{
- typedef itk::Image<double, 3> DistanceImageType;
- typedef itk::ImageRegionIteratorWithIndex<DistanceImageType> ImageIterator;
- typedef itk::NeighborhoodIterator<DistanceImageType> NeighborhoodImageIterator;
-
- DistanceImageType::Pointer distanceImg = DistanceImageType::New();
-
- //Determin the bounding box of the delineated contours
- double xmin = m_Centers.at(0)[0];
- double ymin = m_Centers.at(0)[1];
- double zmin = m_Centers.at(0)[2];
- double xmax = m_Centers.at(0)[0];
- double ymax = m_Centers.at(0)[1];
- double zmax = m_Centers.at(0)[2];
-
- for (unsigned int i = 1; i < m_Centers.size(); i++)
- {
- if (xmin > m_Centers.at(i)[0])
- {
- xmin = m_Centers.at(i)[0];
- }
- if (ymin > m_Centers.at(i)[1])
- {
- ymin = m_Centers.at(i)[1];
- }
- if (zmin > m_Centers.at(i)[2])
- {
- zmin = m_Centers.at(i)[2];
- }
- if (xmax < m_Centers.at(i)[0])
- {
- xmax = m_Centers.at(i)[0];
- }
- if (ymax < m_Centers.at(i)[1])
- {
- ymax = m_Centers.at(i)[1];
- }
- if (zmax < m_Centers.at(i)[2])
- {
- zmax = m_Centers.at(i)[2];
- }
- }
-
- Vector3D extentMM;
- extentMM[0] = xmax - xmin + 5;
- extentMM[1] = ymax - ymin + 5;
- extentMM[2] = zmax - zmin + 5;
-
- //Shifting the distance image's offest to achieve an exact distance calculation
- xmin = xmin - 5;
- ymin = ymin - 5;
- zmin = zmin - 5;
-
- /*
- Now create an empty distance image. The create image will always have the same size, independent from
- the original image (e.g. always consists of 500000 pixels) and will have an isotropic spacing.
- The spacing is calculated like the following:
- The image's volume = 500000 Pixels = extentX*spacing*extentY*spacing*extentZ*spacing
- So the spacing is: spacing = ( 500000 / extentX*extentY*extentZ )^(1/3)
- */
-
- double basis = (extentMM[0]*extentMM[1]*extentMM[2]) / m_DistanceImageVolume;
- double exponent = 1.0/3.0;
- double distImgSpacing = pow(basis, exponent);
- int tempSpacing = (distImgSpacing+0.05)*10;
- m_DistanceImageSpacing = (double)tempSpacing/10.0;
-
- unsigned int numberOfXPixel = extentMM[0] / m_DistanceImageSpacing;
- unsigned int numberOfYPixel = extentMM[1] / m_DistanceImageSpacing;
- unsigned int numberOfZPixel = extentMM[2] / m_DistanceImageSpacing;
-
- DistanceImageType::SizeType size;
-
- //Increase the distance image's size a little bit to achieve an exact distance calculation
- size[0] = numberOfXPixel + 5;
- size[1] = numberOfYPixel + 5;
- size[2] = numberOfZPixel + 5;
-
- DistanceImageType::IndexType start;
- start[0] = 0;
- start[1] = 0;
- start[2] = 0;
-
- DistanceImageType::RegionType lpRegion;
-
- lpRegion.SetSize(size);
- lpRegion.SetIndex(start);
-
- distanceImg->SetRegions( lpRegion );
- distanceImg->SetSpacing( m_DistanceImageSpacing );
- distanceImg->Allocate();
-
- //First of all the image is initialized with the value 10 for each pixel
- distanceImg->FillBuffer(10);
-
- /*
- Now we must caculate the distance for each pixel. But instead of calculating the distance value
- for all of the image's pixels we proceed similar to the region growing algorithm:
-
- 1. Take the first pixel from the narrowband_point_list and calculate the distance for each neighbor (6er)
- 2. If the current index's distance value is below a certain threshold push it into the list
- 3. Next iteration take the next index from the list and start with 1. again
-
- This is done until the narrowband_point_list is empty.
- */
- std::queue<DistanceImageType::IndexType> narrowbandPoints;
- PointType currentPoint = m_Centers.at(0);
- double distance = this->CalculateDistanceValue(currentPoint);
-
- DistanceImageType::IndexType currentIndex;
- currentIndex[0] = ( currentPoint[0]-xmin ) / m_DistanceImageSpacing;
- currentIndex[1] = ( currentPoint[1]-ymin ) / m_DistanceImageSpacing;
- currentIndex[2] = ( currentPoint[2]-zmin ) / m_DistanceImageSpacing;
-
- narrowbandPoints.push(currentIndex);
- distanceImg->SetPixel(currentIndex, distance);
-
- NeighborhoodImageIterator::RadiusType radius;
- radius.Fill(1);
- NeighborhoodImageIterator nIt(radius, distanceImg, distanceImg->GetLargestPossibleRegion());
- unsigned int relativeNbIdx[] = {4, 10, 12, 14, 16, 22};
-
- bool isInBounds = false;
-
- while ( !narrowbandPoints.empty() )
- {
-
- nIt.SetLocation(narrowbandPoints.front());
- narrowbandPoints.pop();
-
- for (int i = 0; i < 6; i++)
- {
- nIt.GetPixel(relativeNbIdx[i], isInBounds);
- if( isInBounds && nIt.GetPixel(relativeNbIdx[i]) == 10)
- {
- currentIndex = nIt.GetIndex(relativeNbIdx[i]);
-
- currentPoint[0] = currentIndex[0]*m_DistanceImageSpacing + xmin;
- currentPoint[1] = currentIndex[1]*m_DistanceImageSpacing + ymin;
- currentPoint[2] = currentIndex[2]*m_DistanceImageSpacing + zmin;
-
- distance = this->CalculateDistanceValue(currentPoint);
- if ( abs(distance) <= m_DistanceImageSpacing*2 )
- {
- nIt.SetPixel(relativeNbIdx[i], distance);
- narrowbandPoints.push(currentIndex);
- }
- }
- }
- }
-
- ImageIterator imgRegionIterator (distanceImg, distanceImg->GetLargestPossibleRegion());
- imgRegionIterator.GoToBegin();
-
- double prevPixelVal = 1;
-
- unsigned int _size[3] = { (unsigned int)(size[0] - 1), (unsigned int)(size[1] - 1), (unsigned int)(size[2] - 1) };
-
- //Set every pixel inside the surface to -10 except the edge point (so that the received surface is closed)
- while (!imgRegionIterator.IsAtEnd()) {
-
- if ( imgRegionIterator.Get() == 10 && prevPixelVal < 0 )
- {
-
- while (imgRegionIterator.Get() == 10)
- {
- if (imgRegionIterator.GetIndex()[0] == _size[0] || imgRegionIterator.GetIndex()[1] == _size[1] || imgRegionIterator.GetIndex()[2] == _size[2]
- || imgRegionIterator.GetIndex()[0] == 0U || imgRegionIterator.GetIndex()[1] == 0U || imgRegionIterator.GetIndex()[2] == 0U )
- {
- imgRegionIterator.Set(10);
- prevPixelVal = 10;
- ++imgRegionIterator;
- break;
- }
- else
- {
- imgRegionIterator.Set(-10);
- ++imgRegionIterator;
- prevPixelVal = -10;
- }
-
- }
-
- }
- else if (imgRegionIterator.GetIndex()[0] == _size[0] || imgRegionIterator.GetIndex()[1] == _size[1] || imgRegionIterator.GetIndex()[2] == _size[2]
- || imgRegionIterator.GetIndex()[0] == 0U || imgRegionIterator.GetIndex()[1] == 0U || imgRegionIterator.GetIndex()[2] == 0U)
- {
- imgRegionIterator.Set(10);
- prevPixelVal = 10;
- ++imgRegionIterator;
- }
- else {
- prevPixelVal = imgRegionIterator.Get();
- ++imgRegionIterator;
- }
-
- }
-
- Image::Pointer resultImage = this->GetOutput();
-
- Point3D origin;
- origin[0] = xmin;
- origin[1] = ymin;
- origin[2] = zmin;
-
- CastToMitkImage(distanceImg, resultImage);
- resultImage->GetGeometry()->SetOrigin(origin);
-}
-
-
-double mitk::CreateDistanceImageFromSurfaceFilter::CalculateDistanceValue(PointType p)
-{
- double distanceValue (0);
- PointType p1;
- PointType p2;
- double norm;
-
- for (unsigned int i = 0; i < m_Centers.size(); i++)
- {
- p1 = m_Centers.at(i);
- p2 = p-p1;
- norm = p2.two_norm();
- distanceValue = distanceValue + norm*m_Weights.get(i);
- }
- return distanceValue;
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::GenerateOutputInformation()
-{
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::PrintEquationSystem()
-{
- std::ofstream esfile;
- esfile.open("C:/Users/fetzer/Desktop/equationSystem/es.txt");
- esfile<<"Nummber of rows: "<<m_SolutionMatrix.rows()<<" ****** Number of columns: "<<m_SolutionMatrix.columns()<<endl;
- esfile<<"[ ";
- for (unsigned int i = 0; i < m_SolutionMatrix.rows(); i++)
- {
- for (unsigned int j = 0; j < m_SolutionMatrix.columns(); j++)
- {
- esfile<<m_SolutionMatrix(i,j)<<" ";
- }
- esfile<<";"<<endl;
- }
- esfile<<" ]";
- esfile.close();
-
- std::ofstream centersFile;
- centersFile.open("C:/Users/fetzer/Desktop/equationSystem/centers.txt");
- for (unsigned int i = 0; i < m_Centers.size(); i++)
- {
- centersFile<<m_Centers.at(i)<<";"<<endl;
- }
- centersFile.close();
-
-}
-
-
-
-void mitk::CreateDistanceImageFromSurfaceFilter::SetInput( const mitk::Surface* surface )
-{
- this->SetInput( 0, const_cast<mitk::Surface*>( surface ) );
-}
-
-
-void mitk::CreateDistanceImageFromSurfaceFilter::SetInput( unsigned int idx, const mitk::Surface* surface )
-{
- if ( this->GetInput(idx) != surface )
- {
- this->SetNthInput( idx, const_cast<mitk::Surface*>( surface ) );
- this->Modified();
- }
-}
-
-
-const mitk::Surface* mitk::CreateDistanceImageFromSurfaceFilter::GetInput()
-{
- if (this->GetNumberOfInputs() < 1)
- return NULL;
-
- return static_cast<const mitk::Surface*>(this->ProcessObject::GetInput(0));
-}
-
-
-const mitk::Surface* mitk::CreateDistanceImageFromSurfaceFilter::GetInput( unsigned int idx)
-{
- if (this->GetNumberOfInputs() < 1)
- return NULL;
-
- return static_cast<const mitk::Surface*>(this->ProcessObject::GetInput(idx));
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::RemoveInputs(mitk::Surface* input)
-{
- this->RemoveInput(input);
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::Reset()
-{
- for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
- {
- this->PopBackInput();
- }
- this->SetNumberOfInputs(0);
- this->SetNumberOfOutputs(1);
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::SetUseProgressBar(bool status)
-{
- this->m_UseProgressBar = status;
-}
-
-void mitk::CreateDistanceImageFromSurfaceFilter::SetProgressStepSize(unsigned int stepSize)
-{
- this->m_ProgressStepSize = stepSize;
-}
diff --git a/Modules/Segmentation/Algorithms/mitkCreateDistanceImageFromSurfaceFilter.h b/Modules/Segmentation/Algorithms/mitkCreateDistanceImageFromSurfaceFilter.h
deleted file mode 100644
index a822411..0000000
--- a/Modules/Segmentation/Algorithms/mitkCreateDistanceImageFromSurfaceFilter.h
+++ /dev/null
@@ -1,161 +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 mitkCreateDistanceImageFromSurfaceFilter_h_Included
-#define mitkCreateDistanceImageFromSurfaceFilter_h_Included
-
-#include "SegmentationExports.h"
-
-#include "mitkImageSource.h"
-#include "mitkSurface.h"
-#include "mitkProgressBar.h"
-
-#include "vtkSmartPointer.h"
-#include "vtkDoubleArray.h"
-#include "vtkCellArray.h"
-#include "vtkCellData.h"
-#include "vtkPolyData.h"
-
-#include "vnl/vnl_matrix.h"
-#include "vnl/vnl_vector.h"
-#include "vnl/vnl_vector_fixed.h"
-#include "vnl/algo/vnl_qr.h"
-
-#include "itkImage.h"
-#include "itkImageRegionIteratorWithIndex.h"
-#include "itkNeighborhoodIterator.h"
-
-#include <queue>
-
-namespace mitk {
-
- /**
- \brief This filter interpolates the 3D surface for a segmented area. The basis for the interpolation
- are the edge-points of contours that are drawn into an image.
-
- The interpolation itself is performed via Radial Basis Function Interpolation.
-
- ATTENTION:
- This filter needs beside the edge points of the delineated contours additionally the normals for each
- edge point.
-
- \sa mitkSurfaceInterpolationController
-
- Based on the contour edge points and their normal this filter calculates a distance function with the following
- properties:
- - Putting a point into the distance function that lies inside the considered surface gives a negativ scalar value
- - Putting a point into the distance function that lies outside the considered surface gives a positive scalar value
- - Putting a point into the distance function that lies exactly on the considered surface gives the value zero
-
- With this interpolated distance function a distance image will be created. The desired surface can then be extract e.g.
- with the marching cubes algorithm. (Within the distance image the surface goes exactly where the pixelvalues are zero)
-
- Note that the obtained distance image has always an isotropig spacing. The size (in this case volume) of the image can be
- adjusted by calling SetDistanceImageVolume(unsigned int volume) which specifies the number ob pixels enclosed by the image.
-
- \ingroup Process
-
- $Author: fetzer$
- */
- class Segmentation_EXPORT CreateDistanceImageFromSurfaceFilter : public ImageSource
- {
-
- public:
-
- typedef vnl_vector_fixed<double,3> PointType;
-
- typedef std::vector< PointType > NormalList;
- typedef std::vector< PointType > CenterList;
-
- typedef vnl_matrix<double> SolutionMatrix;
- typedef vnl_vector<double> FunctionValues;
- typedef vnl_vector<double> InterpolationWeights;
-
- typedef std::vector<Surface::Pointer> SurfaceList;
-
- mitkClassMacro(CreateDistanceImageFromSurfaceFilter,ImageSource);
- itkNewMacro(Self);
-
- //Methods copied from mitkSurfaceToSurfaceFilter
- virtual void SetInput( const mitk::Surface* surface );
-
- virtual void SetInput( unsigned int idx, const mitk::Surface* surface );
-
- virtual const mitk::Surface* GetInput();
-
- virtual const mitk::Surface* GetInput( unsigned int idx );
-
- virtual void RemoveInputs(mitk::Surface* input);
-
-
- /**
- \brief Set the size of the output distance image. The size is specified by the image's volume
- (i.e. in this case how many pixels are enclosed by the image)
- If non is set, the volume will be 500000 pixels.
- */
- itkSetMacro(DistanceImageVolume, unsigned int);
-
- void PrintEquationSystem();
-
- //Resets the filter, i.e. removes all inputs and outputs
- void Reset();
-
- /**
- \brief Set whether the mitkProgressBar should be used
-
- \a Parameter true for using the progress bar, false otherwise
- */
- void SetUseProgressBar(bool);
-
- /**
- \brief Set the stepsize which the progress bar should proceed
-
- \a Parameter The stepsize for progressing
- */
- void SetProgressStepSize(unsigned int stepSize);
-
- protected:
- CreateDistanceImageFromSurfaceFilter();
- virtual ~CreateDistanceImageFromSurfaceFilter();
- virtual void GenerateData();
- virtual void GenerateOutputInformation();
-
-
- private:
-
- void CreateSolutionMatrixAndFunctionValues();
- double CalculateDistanceValue(PointType p);
-
- void CreateDistanceImage ();
-
- //Datastructures for the interpolation
- CenterList m_Centers;
- NormalList m_Normals;
- FunctionValues m_FunctionValues;
- InterpolationWeights m_Weights;
- SolutionMatrix m_SolutionMatrix;
- double m_DistanceImageSpacing;
-
- unsigned int m_DistanceImageVolume;
-
- bool m_UseProgressBar;
- unsigned int m_ProgressStepSize;
-};
-
-}//namespace
-
-
-#endif
diff --git a/Modules/Segmentation/Algorithms/mitkReduceContourSetFilter.cpp b/Modules/Segmentation/Algorithms/mitkReduceContourSetFilter.cpp
deleted file mode 100644
index 8332034..0000000
--- a/Modules/Segmentation/Algorithms/mitkReduceContourSetFilter.cpp
+++ /dev/null
@@ -1,492 +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.
-
-===================================================================*/
-
-#include "mitkReduceContourSetFilter.h"
-
-mitk::ReduceContourSetFilter::ReduceContourSetFilter()
-{
- m_MaxSegmentLenght = 0;
- m_StepSize = 10;
- m_Tolerance = -1;
- m_ReductionType = DOUGLAS_PEUCKER;
- m_MaxSpacing = -1;
- m_MinSpacing = -1;
- this->m_UseProgressBar = false;
- this->m_ProgressStepSize = 1;
- m_NumberOfPointsAfterReduction = 0;
-}
-
-mitk::ReduceContourSetFilter::~ReduceContourSetFilter()
-{
-}
-
-void mitk::ReduceContourSetFilter::GenerateData()
-{
- unsigned int numberOfInputs = this->GetNumberOfInputs();
- unsigned int numberOfOutputs (0);
-
- vtkSmartPointer<vtkPolyData> newPolyData;
- vtkSmartPointer<vtkCellArray> newPolygons;
- vtkSmartPointer<vtkPoints> newPoints;
-
- //For the purpose of evaluation
-// unsigned int numberOfPointsBefore (0);
- m_NumberOfPointsAfterReduction=0;
-
- for(unsigned int i = 0; i < numberOfInputs; i++)
- {
- mitk::Surface* currentSurface = const_cast<mitk::Surface*>( this->GetInput(i) );
- vtkSmartPointer<vtkPolyData> polyData = currentSurface->GetVtkPolyData();
-
- newPolyData = vtkPolyData::New();
- newPolygons = vtkCellArray::New();
- newPoints = vtkPoints::New();
-
- vtkSmartPointer<vtkCellArray> existingPolys = polyData->GetPolys();
-
- vtkSmartPointer<vtkPoints> existingPoints = polyData->GetPoints();
-
- existingPolys->InitTraversal();
-
- vtkIdType* cell (NULL);
- vtkIdType cellSize (0);
-
- for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
- {
- bool incorporatePolygon = this->CheckForIntersection(cell,cellSize,existingPoints, /*numberOfIntersections, intersectionPoints, */i);
- if ( !incorporatePolygon ) continue;
-
- vtkSmartPointer<vtkPolygon> newPolygon = vtkSmartPointer<vtkPolygon>::New();
-
- if(m_ReductionType == NTH_POINT)
- {
- this->ReduceNumberOfPointsByNthPoint(cellSize, cell, existingPoints, newPolygon, newPoints);
- if (newPolygon->GetPointIds()->GetNumberOfIds() != 0)
- {
- newPolygons->InsertNextCell(newPolygon);
- }
- }
- else if (m_ReductionType == DOUGLAS_PEUCKER)
- {
- this->ReduceNumberOfPointsByDouglasPeucker(cellSize, cell, existingPoints, newPolygon, newPoints);
- if (newPolygon->GetPointIds()->GetNumberOfIds() > 3)
- {
- newPolygons->InsertNextCell(newPolygon);
- }
- }
-
- //Again for evaluation
-// numberOfPointsBefore += cellSize;
- m_NumberOfPointsAfterReduction += newPolygon->GetPointIds()->GetNumberOfIds();
-
- }
-
- if (newPolygons->GetNumberOfCells() != 0)
- {
- newPolyData->SetPolys(newPolygons);
- newPolyData->SetPoints(newPoints);
- newPolyData->BuildLinks();
-
- Surface::Pointer surface = this->GetOutput(numberOfOutputs);
- surface->SetVtkPolyData(newPolyData);
- numberOfOutputs++;
- }
-
- }
-// MITK_INFO<<"Points before: "<<numberOfPointsBefore<<" ##### Points after: "<<numberOfPointsAfter;
- this->SetNumberOfOutputs(numberOfOutputs);
-
- //Setting progressbar
- if (this->m_UseProgressBar)
- mitk::ProgressBar::GetInstance()->Progress(this->m_ProgressStepSize);
-}
-
-void mitk::ReduceContourSetFilter::ReduceNumberOfPointsByNthPoint (vtkIdType cellSize, vtkIdType* cell, vtkPoints* points, vtkPolygon* reducedPolygon, vtkPoints* reducedPoints)
-{
-
- unsigned int newNumberOfPoints (0);
- unsigned int mod = cellSize%m_StepSize;
-
- if(mod == 0)
- {
- newNumberOfPoints = cellSize/m_StepSize;
- }
- else
- {
- newNumberOfPoints = ( (cellSize-mod)/m_StepSize )+1;
- }
-
- if (newNumberOfPoints <= 3)
- {
- return;
- }
- reducedPolygon->GetPointIds()->SetNumberOfIds(newNumberOfPoints);
- reducedPolygon->GetPoints()->SetNumberOfPoints(newNumberOfPoints);
-
- for (unsigned int i = 0; i < cellSize; i++)
- {
- if (i%m_StepSize == 0)
- {
- double point[3];
- points->GetPoint(cell[i], point);
- vtkIdType id = reducedPoints->InsertNextPoint(point);
- reducedPolygon->GetPointIds()->SetId(i/m_StepSize, id);
- }
-
- }
- vtkIdType id = cell[0];
- double point[3];
- points->GetPoint(id, point);
- id = reducedPoints->InsertNextPoint(point);
- reducedPolygon->GetPointIds()->SetId(newNumberOfPoints-1, id);
-
-}
-
-void mitk::ReduceContourSetFilter::ReduceNumberOfPointsByDouglasPeucker(vtkIdType cellSize, vtkIdType* cell, vtkPoints* points,
- vtkPolygon* reducedPolygon, vtkPoints* reducedPoints)
-{
- //If the cell is too small to obtain a reduced polygon with the given stepsize return
- if (cellSize <= m_StepSize*3)return;
-
- /*
- What we do now is (see the Douglas Peucker Algorithm):
-
- 1. Divide the current contour in two line segments (start - middle; middle - end), put them into the stack
-
- 2. Fetch first line segment and create the following vectors:
-
- - v1 = (start;end)
- - v2 = (start;currentPoint) -> for each point of the current line segment!
-
- 3. Calculate the distance from the currentPoint to v1:
-
- a. Determine the length of the orthogonal projection of v2 to v1 by:
-
- l = v2 * (normalized v1)
-
- b. There a three possibilities for the distance then:
-
- d = sqrt(lenght(v2)^2 - l^2) if l > 0 and l < length(v1)
- d = lenght(v2-v1) if l > 0 and l > lenght(v1)
- d = length(v2) if l < 0 because v2 is then pointing in a different direction than v1
-
- 4. Memorize the point with the biggest distance and create two new line segments with it at the end of the iteration
- and put it into the stack
-
- 5. If the distance value D <= m_Tolerance, then add the start and end index and the corresponding points to the reduced ones
- */
-
- //First of all set tolerance if none is specified
- if(m_Tolerance < 0)
- {
- if(m_MaxSpacing > 0)
- {
- m_Tolerance = m_MinSpacing;
- }
- else
- {
- m_Tolerance = 1.5;
- }
- }
-
- std::stack<LineSegment> lineSegments;
-
- //1. Divide in line segments
-
- LineSegment ls2;
- ls2.StartIndex = cell[cellSize/2];
- ls2.EndIndex = cell[cellSize-1];
- lineSegments.push(ls2);
-
- LineSegment ls1;
- ls1.StartIndex = cell[0];
- ls1.EndIndex = cell[cellSize/2];
- lineSegments.push(ls1);
-
- LineSegment currentSegment;
- double v1[3];
- double v2[3];
- double tempV[3];
- double lenghtV1;
-
- double currentMaxDistance (0);
- vtkIdType currentMaxDistanceIndex (0);
-
- double l;
- double d;
-
- vtkIdType pointId (0);
- //Add the start index to the reduced points. From now on just the end indices will be added
- pointId = reducedPoints->InsertNextPoint(points->GetPoint(cell[0]));
- reducedPolygon->GetPointIds()->InsertNextId(pointId);
-
- while (!lineSegments.empty())
- {
- currentSegment = lineSegments.top();
- lineSegments.pop();
-
- //2. Create vectors
-
- points->GetPoint(currentSegment.EndIndex, tempV);
- points->GetPoint(currentSegment.StartIndex, v1);
-
- v1[0] = tempV[0]-v1[0];
- v1[1] = tempV[1]-v1[1];
- v1[2] = tempV[2]-v1[2];
-
- lenghtV1 = vtkMath::Norm(v1);
-
- vtkMath::Normalize(v1);
- int range = currentSegment.EndIndex - currentSegment.StartIndex;
- for (int i = 1; i < abs(range); ++i)
- {
- points->GetPoint(currentSegment.StartIndex+i, tempV);
- points->GetPoint(currentSegment.StartIndex, v2);
-
- v2[0] = tempV[0]-v2[0];
- v2[1] = tempV[1]-v2[1];
- v2[2] = tempV[2]-v2[2];
-
- //3. Calculate the distance
-
- l = vtkMath::Dot(v2, v1);
-
- d = vtkMath::Norm(v2);
-
- if (l > 0 && l < lenghtV1)
- {
- d = sqrt((d*d-l*l));
- }
- else if (l > 0 && l > lenghtV1)
- {
- tempV[0] = lenghtV1*v1[0] - v2[0];
- tempV[1] = lenghtV1*v1[1] - v2[1];
- tempV[2] = lenghtV1*v1[2] - v2[2];
- d = vtkMath::Norm(tempV);
- }
-
- //4. Memorize maximum distance
- if (d > currentMaxDistance)
- {
- currentMaxDistance = d;
- currentMaxDistanceIndex = currentSegment.StartIndex+i;
- }
-
- }
-
- //4. & 5.
- if (currentMaxDistance <= m_Tolerance)
- {
-
- //double temp[3];
- int segmentLenght = currentSegment.EndIndex - currentSegment.StartIndex;
- if (segmentLenght > (int)m_MaxSegmentLenght)
- {
- m_MaxSegmentLenght = (unsigned int)segmentLenght;
- }
-
-// MITK_INFO<<"Lenght: "<<abs(segmentLenght);
- if (abs(segmentLenght) > 25)
- {
- unsigned int newLenght(segmentLenght);
- while (newLenght > 25)
- {
- newLenght = newLenght*0.5;
- }
- unsigned int divisions = abs(segmentLenght)/newLenght;
-// MITK_INFO<<"Divisions: "<<divisions;
-
- for (unsigned int i = 1; i<=divisions; ++i)
- {
-// MITK_INFO<<"Inserting MIDDLE: "<<(currentSegment.StartIndex + newLenght*i);
- pointId = reducedPoints->InsertNextPoint(points->GetPoint(currentSegment.StartIndex + newLenght*i));
- reducedPolygon->GetPointIds()->InsertNextId(pointId);
- }
- }
-// MITK_INFO<<"Inserting END: "<<currentSegment.EndIndex;
- pointId = reducedPoints->InsertNextPoint(points->GetPoint(currentSegment.EndIndex));
- reducedPolygon->GetPointIds()->InsertNextId(pointId);
- }
- else
- {
- ls2.StartIndex = currentMaxDistanceIndex;
- ls2.EndIndex = currentSegment.EndIndex;
- lineSegments.push(ls2);
-
- ls1.StartIndex = currentSegment.StartIndex;
- ls1.EndIndex = currentMaxDistanceIndex;
- lineSegments.push(ls1);
- }
- currentMaxDistance = 0;
-
- }
-
-}
-
-bool mitk::ReduceContourSetFilter::CheckForIntersection (vtkIdType* currentCell, vtkIdType currentCellSize, vtkPoints* currentPoints,/* vtkIdType numberOfIntersections, vtkIdType* intersectionPoints,*/ unsigned int currentInputIndex)
-{
- /*
- If we check the current cell for intersections then we have to consider three possibilies:
- 1. There is another cell among all the other input surfaces which intersects the current polygon:
- - That means we have to save the intersection points because these points should not be eliminated
- 2. There current polygon exists just because of an intersection of another polygon with the current plane defined by the current polygon
- - That means the current polygon should not be incorporated and all of its points should be eliminated
- 3. There is no intersection
- - That mean we can just reduce the current polygons points without considering any intersections
- */
-
- for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
- {
- //Don't check for intersection with the polygon itself
- if (i == currentInputIndex) continue;
-
- //Get the next polydata to check for intersection
- vtkSmartPointer<vtkPolyData> poly = const_cast<Surface*>( this->GetInput(i) )->GetVtkPolyData();
- vtkSmartPointer<vtkCellArray> polygonArray = poly->GetPolys();
- polygonArray->InitTraversal();
- vtkIdType anotherInputPolygonSize (0);
- vtkIdType* anotherInputPolygonIDs(NULL);
-
- /*
- The procedure is:
- - Create the equation of the plane, defined by the points of next input
- - Calculate the distance of each point of the current polygon to the plane
- - If the maximum distance is not bigger than 1.5 of the maximum spacing AND the minimal distance is not bigger
- than 0.5 of the minimum spacing then the current contour is an intersection contour
- */
-
- for( polygonArray->InitTraversal(); polygonArray->GetNextCell(anotherInputPolygonSize, anotherInputPolygonIDs);)
- {
- //Choosing three plane points to calculate the plane vectors
- double p1[3];
- double p2[3];
- double p3[3];
-
- //The plane vectors
- double v1[3];
- double v2[3] = { 0 };
- //The plane normal
- double normal[3];
-
- //Create first Vector
- poly->GetPoint(anotherInputPolygonIDs[0], p1);
- poly->GetPoint(anotherInputPolygonIDs[1], p2);
-
- v1[0] = p2[0]-p1[0];
- v1[1] = p2[1]-p1[1];
- v1[2] = p2[2]-p1[2];
-
- //Find 3rd point for 2nd vector (The angle between the two plane vectors should be bigger than 30 degrees)
-
- double maxDistance (0);
- double minDistance (10000);
-
- for (unsigned int j = 2; j < anotherInputPolygonSize; j++)
- {
- poly->GetPoint(anotherInputPolygonIDs[j], p3);
-
- v2[0] = p3[0]-p1[0];
- v2[1] = p3[1]-p1[1];
- v2[2] = p3[2]-p1[2];
-
- //Calculate the angle between the two vector for the current point
- double dotV1V2 = vtkMath::Dot(v1,v2);
- double absV1 = sqrt(vtkMath::Dot(v1,v1));
- double absV2 = sqrt(vtkMath::Dot(v2,v2));
- double cosV1V2 = dotV1V2/(absV1*absV2);
-
- double arccos = acos(cosV1V2);
- double degree = vtkMath::DegreesFromRadians(arccos);
-
- //If angle is bigger than 30 degrees break
- if (degree > 30) break;
-
- }//for (to find 3rd point)
-
- //Calculate normal of the plane by taking the cross product of the two vectors
- vtkMath::Cross(v1,v2,normal);
- vtkMath::Normalize(normal);
-
- //Determine position of the plane
- double lambda = vtkMath::Dot(normal, p1);
-
- /*
- Calculate the distance to the plane for each point of the current polygon
- If the distance is zero then save the currentPoint as intersection point
- */
- for (unsigned int k = 0; k < currentCellSize; k++)
- {
- double currentPoint[3];
- currentPoints->GetPoint(currentCell[k], currentPoint);
-
- double tempPoint[3];
- tempPoint[0] = normal[0]*currentPoint[0];
- tempPoint[1] = normal[1]*currentPoint[1];
- tempPoint[2] = normal[2]*currentPoint[2];
-
- double temp = tempPoint[0]+tempPoint[1]+tempPoint[2]-lambda;
- double distance = fabs(temp);
-
- if (distance > maxDistance)
- {
- maxDistance = distance;
- }
- if (distance < minDistance)
- {
- minDistance = distance;
- }
- }//for (to calculate distance and intersections with currentPolygon)
-
- if (maxDistance < 1.5*m_MaxSpacing && minDistance < 0.5*m_MinSpacing)
- {
- return false;
- }
-
- //Because we are considering the plane defined by the acual input polygon only one iteration is sufficient
- //We do not need to consider each cell of the plane
- break;
- }//for (to traverse through all cells of actualInputPolyData)
-
- }//for (to iterate through all inputs)
-
- return true;
-}
-
-void mitk::ReduceContourSetFilter::GenerateOutputInformation()
-{
- Superclass::GenerateOutputInformation();
-}
-
-void mitk::ReduceContourSetFilter::Reset()
-{
- for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
- {
- this->PopBackInput();
- }
- this->SetNumberOfInputs(0);
- this->SetNumberOfOutputs(0);
- m_NumberOfPointsAfterReduction = 0;
-}
-
-void mitk::ReduceContourSetFilter::SetUseProgressBar(bool status)
-{
- this->m_UseProgressBar = status;
-}
-
-void mitk::ReduceContourSetFilter::SetProgressStepSize(unsigned int stepSize)
-{
- this->m_ProgressStepSize = stepSize;
-}
diff --git a/Modules/Segmentation/Algorithms/mitkReduceContourSetFilter.h b/Modules/Segmentation/Algorithms/mitkReduceContourSetFilter.h
deleted file mode 100644
index 2f18782..0000000
--- a/Modules/Segmentation/Algorithms/mitkReduceContourSetFilter.h
+++ /dev/null
@@ -1,128 +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 mitkReduceContourSetFilter_h_Included
-#define mitkReduceContourSetFilter_h_Included
-
-#include "SegmentationExports.h"
-#include "mitkSurfaceToSurfaceFilter.h"
-#include "mitkSurface.h"
-#include "mitkProgressBar.h"
-
-#include "vtkSmartPointer.h"
-#include "vtkPolyData.h"
-#include "vtkPolygon.h"
-#include "vtkPoints.h"
-#include "vtkCellArray.h"
-#include "vtkMath.h"
-
-#include <stack>
-
-namespace mitk {
-
-/**
- \brief A filter that reduces the number of points of contours represented by a mitk::Surface
-
- The type of the reduction can be set via SetReductionType. The two ways provided by this filter is the
-
- \li NTH_POINT Algorithm which reduces the contours according to a certain stepsize
- \li DOUGLAS_PEUCKER Algorithm which incorpates an error tolerance into the reduction.
-
- Stepsize and error tolerance can be set via SetStepSize and SetTolerance.
-
- Additional if more than one input contour is provided this filter tries detect contours which occur just because
- of an intersection. These intersection contours are eliminated. In oder to ensure a correct elimination the min and max
- spacing of the original image must be provided.
-
- The output is a mitk::Surface.
-
- $Author: fetzer$
-*/
-
- class Segmentation_EXPORT ReduceContourSetFilter : public SurfaceToSurfaceFilter
- {
-
- public:
-
- enum Reduction_Type
- {
- NTH_POINT, DOUGLAS_PEUCKER
- };
-
- struct LineSegment
- {
- unsigned int StartIndex;
- unsigned int EndIndex;
- };
-
- mitkClassMacro(ReduceContourSetFilter,SurfaceToSurfaceFilter);
- itkNewMacro(Self);
-
- itkSetMacro(MinSpacing, double);
- itkSetMacro(MaxSpacing, double);
- itkSetMacro(ReductionType, Reduction_Type);
- itkSetMacro(StepSize, unsigned int);
- itkSetMacro(Tolerance, double);
-
- itkGetMacro(NumberOfPointsAfterReduction, unsigned int);
-
- //Resets the filter, i.e. removes all inputs and outputs
- void Reset();
-
- /**
- \brief Set whether the mitkProgressBar should be used
-
- \a Parameter true for using the progress bar, false otherwise
- */
- void SetUseProgressBar(bool);
-
- /**
- \brief Set the stepsize which the progress bar should proceed
-
- \a Parameter The stepsize for progressing
- */
- void SetProgressStepSize(unsigned int stepSize);
-
- protected:
- ReduceContourSetFilter();
- virtual ~ReduceContourSetFilter();
- virtual void GenerateData();
- virtual void GenerateOutputInformation();
-
- private:
- void ReduceNumberOfPointsByNthPoint (vtkIdType cellSize, vtkIdType* cell, vtkPoints* points, vtkPolygon* reducedPolygon, vtkPoints* reducedPoints);
-
- void ReduceNumberOfPointsByDouglasPeucker (vtkIdType cellSize, vtkIdType* cell, vtkPoints* points, vtkPolygon* reducedPolygon, vtkPoints* reducedPoints);
-
- bool CheckForIntersection (vtkIdType* currentCell, vtkIdType currentCellSize, vtkPoints* currentPoints, /*vtkIdType numberOfIntersections, vtkIdType* intersectionPoints,*/ unsigned int currentInputIndex);
-
- double m_MinSpacing;
- double m_MaxSpacing;
-
- Reduction_Type m_ReductionType;
- unsigned int m_StepSize;
- double m_Tolerance;
- unsigned int m_MaxSegmentLenght;
-
- bool m_UseProgressBar;
- unsigned int m_ProgressStepSize;
-
- unsigned int m_NumberOfPointsAfterReduction;
-
- };//class
-
-}//namespace
-#endif
diff --git a/Modules/Segmentation/CMakeLists.txt b/Modules/Segmentation/CMakeLists.txt
index 133b39a..cc67575 100644
--- a/Modules/Segmentation/CMakeLists.txt
+++ b/Modules/Segmentation/CMakeLists.txt
@@ -3,8 +3,8 @@
#configure_file(${PROJECT_SOURCE_DIR}/CMake/ToolGUIExtensionITKFactory.cpp.in $#{PROJECT_BINARY_DIR}/ToolGUIExtensionITKFactory.cpp.in COPYONLY)
MITK_CREATE_MODULE( Segmentation
- INCLUDE_DIRS Algorithms Controllers DataManagement Interactions IO Rendering
- DEPENDS Mitk ipSegmentation mitkIpFunc MitkExt
+ INCLUDE_DIRS Algorithms DataManagement Interactions IO Rendering
+ DEPENDS Mitk ipSegmentation mitkIpFunc SurfaceInterpolation MitkExt
)
add_subdirectory(Testing)
\ No newline at end of file
diff --git a/Modules/Segmentation/Controllers/mitkSurfaceInterpolationController.cpp b/Modules/Segmentation/Controllers/mitkSurfaceInterpolationController.cpp
deleted file mode 100644
index 770f81f..0000000
--- a/Modules/Segmentation/Controllers/mitkSurfaceInterpolationController.cpp
+++ /dev/null
@@ -1,266 +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.
-
-===================================================================*/
-
-#include "mitkSurfaceInterpolationController.h"
-#include "mitkMemoryUtilities.h"
-
-mitk::SurfaceInterpolationController::SurfaceInterpolationController()
- :m_SelectedSegmentation(0)
-{
- m_ReduceFilter = ReduceContourSetFilter::New();
- m_NormalsFilter = ComputeContourSetNormalsFilter::New();
- m_InterpolateSurfaceFilter = CreateDistanceImageFromSurfaceFilter::New();
-
- m_ReduceFilter->SetUseProgressBar(false);
- m_NormalsFilter->SetUseProgressBar(false);
- m_InterpolateSurfaceFilter->SetUseProgressBar(false);
-
- m_Contours = Surface::New();
-
- m_PolyData = vtkSmartPointer<vtkPolyData>::New();
- m_PolyData->SetPoints(vtkPoints::New());
-
- m_InterpolationResult = 0;
- m_CurrentNumberOfReducedContours = 0;
-}
-
-mitk::SurfaceInterpolationController::~SurfaceInterpolationController()
-{
-
- for (ContourListMap::iterator it = m_MapOfContourLists.begin(); it != m_MapOfContourLists.end(); it++)
- {
- for (unsigned int j = 0; j < m_MapOfContourLists[(*it).first].size(); ++j)
- {
- delete(m_MapOfContourLists[(*it).first].at(j).position);
- }
- m_MapOfContourLists.erase(it);
- }
-}
-
-mitk::SurfaceInterpolationController* mitk::SurfaceInterpolationController::GetInstance()
-{
- static mitk::SurfaceInterpolationController* m_Instance;
-
- if ( m_Instance == 0)
- {
- m_Instance = new SurfaceInterpolationController();
- }
-
- return m_Instance;
-}
-
-void mitk::SurfaceInterpolationController::AddNewContour (mitk::Surface::Pointer newContour ,RestorePlanePositionOperation* op)
-{
- AffineTransform3D::Pointer transform = AffineTransform3D::New();
- transform = op->GetTransform();
-
- mitk::Vector3D direction = op->GetDirectionVector();
- int pos (-1);
-
- for (unsigned int i = 0; i < m_MapOfContourLists[m_SelectedSegmentation].size(); i++)
- {
- itk::Matrix<float> diffM = transform->GetMatrix()-m_MapOfContourLists[m_SelectedSegmentation].at(i).position->GetTransform()->GetMatrix();
- bool isSameMatrix(true);
- for (unsigned int j = 0; j < 3; j++)
- {
- if (fabs(diffM[j][0]) > 0.0001 && fabs(diffM[j][1]) > 0.0001 && fabs(diffM[j][2]) > 0.0001)
- {
- isSameMatrix = false;
- break;
- }
- }
- itk::Vector<float> diffV = m_MapOfContourLists[m_SelectedSegmentation].at(i).position->GetTransform()->GetOffset()-transform->GetOffset();
- if ( isSameMatrix && m_MapOfContourLists[m_SelectedSegmentation].at(i).position->GetPos() == op->GetPos() && (fabs(diffV[0]) < 0.0001 && fabs(diffV[1]) < 0.0001 && fabs(diffV[2]) < 0.0001) )
- {
- pos = i;
- break;
- }
-
- }
-
- if (pos == -1)
- {
- //MITK_INFO<<"New Contour";
- mitk::RestorePlanePositionOperation* newOp = new mitk::RestorePlanePositionOperation (OpRESTOREPLANEPOSITION, op->GetWidth(),
- op->GetHeight(), op->GetSpacing(), op->GetPos(), direction, transform);
- ContourPositionPair newData;
- newData.contour = newContour;
- newData.position = newOp;
-
- m_ReduceFilter->SetInput(m_MapOfContourLists[m_SelectedSegmentation].size(), newContour);
- m_MapOfContourLists[m_SelectedSegmentation].push_back(newData);
- }
- else
- {
- //MITK_INFO<<"Modified Contour";
- m_MapOfContourLists[m_SelectedSegmentation].at(pos).contour = newContour;
- m_ReduceFilter->SetInput(pos, newContour);
- }
-
- m_ReduceFilter->Update();
- m_CurrentNumberOfReducedContours = m_ReduceFilter->GetNumberOfOutputs();
-
- for (unsigned int i = 0; i < m_CurrentNumberOfReducedContours; i++)
- {
- m_NormalsFilter->SetInput(i, m_ReduceFilter->GetOutput(i));
- m_InterpolateSurfaceFilter->SetInput(i, m_NormalsFilter->GetOutput(i));
- }
-
- this->Modified();
-}
-
-void mitk::SurfaceInterpolationController::Interpolate()
-{
- if (m_CurrentNumberOfReducedContours< 2)
- return;
-
- //Setting up progress bar
- /*
- * Removed due to bug 12441. ProgressBar messes around with Qt event queue which is fatal for segmentation
- */
- //mitk::ProgressBar::GetInstance()->AddStepsToDo(8);
-
- m_InterpolateSurfaceFilter->Update();
-
- Image::Pointer distanceImage = m_InterpolateSurfaceFilter->GetOutput();
-
- vtkSmartPointer<vtkMarchingCubes> mcFilter = vtkSmartPointer<vtkMarchingCubes>::New();
- mcFilter->SetInput(distanceImage->GetVtkImageData());
- mcFilter->SetValue(0,0);
- mcFilter->Update();
-
- m_InterpolationResult = 0;
- m_InterpolationResult = mitk::Surface::New();
- m_InterpolationResult->SetVtkPolyData(mcFilter->GetOutput());
- m_InterpolationResult->GetGeometry()->SetOrigin(distanceImage->GetGeometry()->GetOrigin());
-
- vtkSmartPointer<vtkAppendPolyData> polyDataAppender = vtkSmartPointer<vtkAppendPolyData>::New();
- for (unsigned int i = 0; i < m_ReduceFilter->GetNumberOfOutputs(); i++)
- {
- polyDataAppender->AddInput(m_ReduceFilter->GetOutput(i)->GetVtkPolyData());
- }
- polyDataAppender->Update();
- m_Contours->SetVtkPolyData(polyDataAppender->GetOutput());
-
- //Last progress step
- /*
- * Removed due to bug 12441. ProgressBar messes around with Qt event queue which is fatal for segmentation
- */
- //mitk::ProgressBar::GetInstance()->Progress(8);
-
- m_InterpolationResult->DisconnectPipeline();
-}
-
-mitk::Surface::Pointer mitk::SurfaceInterpolationController::GetInterpolationResult()
-{
- return m_InterpolationResult;
-}
-
-mitk::Surface* mitk::SurfaceInterpolationController::GetContoursAsSurface()
-{
- return m_Contours;
-}
-
-void mitk::SurfaceInterpolationController::SetDataStorage(DataStorage &ds)
-{
- m_DataStorage = &ds;
-}
-
-void mitk::SurfaceInterpolationController::SetMinSpacing(double minSpacing)
-{
- m_ReduceFilter->SetMinSpacing(minSpacing);
-}
-
-void mitk::SurfaceInterpolationController::SetMaxSpacing(double maxSpacing)
-{
- m_ReduceFilter->SetMaxSpacing(maxSpacing);
- m_NormalsFilter->SetMaxSpacing(maxSpacing);
-}
-
-void mitk::SurfaceInterpolationController::SetDistanceImageVolume(unsigned int distImgVolume)
-{
- m_InterpolateSurfaceFilter->SetDistanceImageVolume(distImgVolume);
-}
-
-void mitk::SurfaceInterpolationController::SetSegmentationImage(Image* workingImage)
-{
- m_NormalsFilter->SetSegmentationBinaryImage(workingImage);
-}
-
-mitk::Image* mitk::SurfaceInterpolationController::GetImage()
-{
- return m_InterpolateSurfaceFilter->GetOutput();
-}
-
-double mitk::SurfaceInterpolationController::EstimatePortionOfNeededMemory()
-{
- double numberOfPointsAfterReduction = m_ReduceFilter->GetNumberOfPointsAfterReduction()*3;
- double sizeOfPoints = pow(numberOfPointsAfterReduction,2)*sizeof(double);
- double totalMem = mitk::MemoryUtilities::GetTotalSizeOfPhysicalRam();
- double percentage = sizeOfPoints/totalMem;
- return percentage;
-}
-
-void mitk::SurfaceInterpolationController::SetCurrentSegmentationInterpolationList(mitk::Image* segmentation)
-{
- if (segmentation == m_SelectedSegmentation)
- return;
-
- if (segmentation == 0)
- return;
-
- ContourListMap::iterator it = m_MapOfContourLists.find(segmentation);
-
- m_SelectedSegmentation = segmentation;
-
- m_ReduceFilter->Reset();
- m_NormalsFilter->Reset();
- m_InterpolateSurfaceFilter->Reset();
-
- if (it == m_MapOfContourLists.end())
- {
- ContourPositionPairList newList;
- m_MapOfContourLists.insert(std::pair<mitk::Image*, ContourPositionPairList>(segmentation, newList));
- m_InterpolationResult = 0;
- m_CurrentNumberOfReducedContours = 0;
- }
- else
- {
- for (unsigned int i = 0; i < m_MapOfContourLists[m_SelectedSegmentation].size(); i++)
- {
- m_ReduceFilter->SetInput(i, m_MapOfContourLists[m_SelectedSegmentation].at(i).contour);
- }
-
- m_ReduceFilter->Update();
-
- m_CurrentNumberOfReducedContours = m_ReduceFilter->GetNumberOfOutputs();
-
- for (unsigned int i = 0; i < m_CurrentNumberOfReducedContours; i++)
- {
- m_NormalsFilter->SetInput(i, m_ReduceFilter->GetOutput(i));
- m_InterpolateSurfaceFilter->SetInput(i, m_NormalsFilter->GetOutput(i));
- }
- }
- Modified();
-}
-
-void mitk::SurfaceInterpolationController::RemoveSegmentationFromContourList(mitk::Image *segmentation)
-{
- if (segmentation != 0)
- {
- m_MapOfContourLists.erase(segmentation);
- }
-}
diff --git a/Modules/Segmentation/Controllers/mitkSurfaceInterpolationController.h b/Modules/Segmentation/Controllers/mitkSurfaceInterpolationController.h
deleted file mode 100644
index 6483939..0000000
--- a/Modules/Segmentation/Controllers/mitkSurfaceInterpolationController.h
+++ /dev/null
@@ -1,167 +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 mitkSurfaceInterpolationController_h_Included
-#define mitkSurfaceInterpolationController_h_Included
-
-#include "mitkCommon.h"
-#include "SegmentationExports.h"
-#include "mitkRestorePlanePositionOperation.h"
-#include "mitkSurface.h"
-#include "mitkInteractionConst.h"
-#include "mitkColorProperty.h"
-#include "mitkProperties.h"
-
-#include "mitkCreateDistanceImageFromSurfaceFilter.h"
-#include "mitkReduceContourSetFilter.h"
-#include "mitkComputeContourSetNormalsFilter.h"
-
-#include "mitkDataNode.h"
-#include "mitkDataStorage.h"
-#include "mitkWeakPointer.h"
-
-#include "vtkPolygon.h"
-#include "vtkPoints.h"
-#include "vtkCellArray.h"
-#include "vtkPolyData.h"
-#include "vtkSmartPointer.h"
-#include "vtkAppendPolyData.h"
-
-#include "vtkMarchingCubes.h"
-#include "vtkImageData.h"
-#include "mitkVtkRepresentationProperty.h"
-#include "vtkProperty.h"
-
-#include "mitkProgressBar.h"
-
-namespace mitk
-{
-
- class Segmentation_EXPORT SurfaceInterpolationController : public itk::Object
- {
-
- public:
-
- mitkClassMacro(SurfaceInterpolationController, itk::Object)
- itkNewMacro(Self)
-
- static SurfaceInterpolationController* GetInstance();
-
- /**
- * Adds a new extracted contour to the list
- */
- void AddNewContour(Surface::Pointer newContour, RestorePlanePositionOperation *op);
-
- /**
- * Interpolates the 3D surface from the given extracted contours
- */
- void Interpolate ();
-
- mitk::Surface::Pointer GetInterpolationResult();
-
- /**
- * Sets the minimum spacing of the current selected segmentation
- * This is needed since the contour points we reduced before they are used to interpolate the surface
- */
- void SetMinSpacing(double minSpacing);
-
- /**
- * Sets the minimum spacing of the current selected segmentation
- * This is needed since the contour points we reduced before they are used to interpolate the surface
- */
- void SetMaxSpacing(double maxSpacing);
-
- /**
- * Sets the volume i.e. the number of pixels that the distance image should have
- * By evaluation we found out that 50.000 pixel delivers a good result
- */
- void SetDistanceImageVolume(unsigned int distImageVolume);
-
- /**
- * Sets the current segmentation which is used by the interpolation
- * This is needed because the calculation of the normals needs to now wheather a normal points inside a segmentation or not
- */
- void SetSegmentationImage(Image* workingImage);
-
- Surface* GetContoursAsSurface();
-
- void SetDataStorage(DataStorage &ds);
-
- /**
- * Sets the current list of contourpoints which is used for the surface interpolation
- * @param segmentation The current selected segmentation
- */
- void SetCurrentSegmentationInterpolationList(mitk::Image* segmentation);
-
- /**
- * Removes the segmentation and all its contours from the list
- * @param segmentation The segmentation to be removed
- */
- void RemoveSegmentationFromContourList(mitk::Image* segmentation);
-
- mitk::Image* GetImage();
-
- /**
- * Estimates the memory which is needed to build up the equationsystem for the interpolation.
- * \returns The percentage of the real memory which will be used by the interpolation
- */
- double EstimatePortionOfNeededMemory();
-
- protected:
-
- SurfaceInterpolationController();
-
- ~SurfaceInterpolationController();
-
- private:
-
- struct ContourPositionPair {
- Surface::Pointer contour;
- RestorePlanePositionOperation* position;
- };
-
- typedef std::vector<ContourPositionPair> ContourPositionPairList;
- typedef std::map<mitk::Image* , ContourPositionPairList> ContourListMap;
-
- ContourPositionPairList::iterator m_Iterator;
-
- ReduceContourSetFilter::Pointer m_ReduceFilter;
- ComputeContourSetNormalsFilter::Pointer m_NormalsFilter;
- CreateDistanceImageFromSurfaceFilter::Pointer m_InterpolateSurfaceFilter;
-
- double m_MinSpacing;
- double m_MaxSpacing;
-
- const Image* m_WorkingImage;
-
- Surface::Pointer m_Contours;
-
- vtkSmartPointer<vtkPolyData> m_PolyData;
-
- unsigned int m_DistImageVolume;
-
- mitk::WeakPointer<mitk::DataStorage> m_DataStorage;
-
- ContourListMap m_MapOfContourLists;
-
- mitk::Surface::Pointer m_InterpolationResult;
-
- unsigned int m_CurrentNumberOfReducedContours;
-
- mitk::Image* m_SelectedSegmentation;
- };
-}
-#endif
diff --git a/Modules/Segmentation/files.cmake b/Modules/Segmentation/files.cmake
index 5d11450..6caf1c5 100644
--- a/Modules/Segmentation/files.cmake
+++ b/Modules/Segmentation/files.cmake
@@ -1,16 +1,13 @@
set(CPP_FILES
Algorithms/mitkCalculateSegmentationVolume.cpp
-Algorithms/mitkComputeContourSetNormalsFilter.cpp
Algorithms/mitkContourSetToPointSetFilter.cpp
Algorithms/mitkContourUtils.cpp
Algorithms/mitkCorrectorAlgorithm.cpp
-Algorithms/mitkCreateDistanceImageFromSurfaceFilter.cpp
Algorithms/mitkDiffImageApplier.cpp
Algorithms/mitkImageToContourFilter.cpp
Algorithms/mitkManualSegmentationToSurfaceFilter.cpp
Algorithms/mitkOverwriteDirectedPlaneImageFilter.cpp
Algorithms/mitkOverwriteSliceImageFilter.cpp
-Algorithms/mitkReduceContourSetFilter.cpp
Algorithms/mitkSegmentationObjectFactory.cpp
Algorithms/mitkSegmentationSink.cpp
Algorithms/mitkShapeBasedInterpolationAlgorithm.cpp
@@ -24,7 +21,6 @@ Algorithms/mitkContourModelToPointSetFilter.cpp
Algorithms/mitkContourModelToSurfaceFilter.cpp
Algorithms/mitkContourModelSubDivisionFilter.cpp
Controllers/mitkSegmentationInterpolationController.cpp
-Controllers/mitkSurfaceInterpolationController.cpp
# DataManagement/mitkApplyDiffImageOperation.cpp
DataManagement/mitkContour.cpp
DataManagement/mitkContourSet.cpp
diff --git a/Modules/SurfaceInterpolation/CMakeLists.txt b/Modules/SurfaceInterpolation/CMakeLists.txt
new file mode 100644
index 0000000..0aba049
--- /dev/null
+++ b/Modules/SurfaceInterpolation/CMakeLists.txt
@@ -0,0 +1,3 @@
+MITK_CREATE_MODULE( SurfaceInterpolation
+ DEPENDS Mitk ImageExtraction
+)
diff --git a/Modules/SurfaceInterpolation/files.cmake b/Modules/SurfaceInterpolation/files.cmake
new file mode 100644
index 0000000..247737c
--- /dev/null
+++ b/Modules/SurfaceInterpolation/files.cmake
@@ -0,0 +1,6 @@
+set(CPP_FILES
+mitkComputeContourSetNormalsFilter.cpp
+mitkCreateDistanceImageFromSurfaceFilter.cpp
+mitkReduceContourSetFilter.cpp
+mitkSurfaceInterpolationController.cpp
+)
diff --git a/Modules/SurfaceInterpolation/mitkComputeContourSetNormalsFilter.cpp b/Modules/SurfaceInterpolation/mitkComputeContourSetNormalsFilter.cpp
new file mode 100644
index 0000000..9797416
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkComputeContourSetNormalsFilter.cpp
@@ -0,0 +1,305 @@
+/*===================================================================
+
+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.
+
+===================================================================*/
+
+#include "mitkComputeContourSetNormalsFilter.h"
+
+
+mitk::ComputeContourSetNormalsFilter::ComputeContourSetNormalsFilter()
+{
+ m_MaxSpacing = 5;
+ this->m_UseProgressBar = false;
+ this->m_ProgressStepSize = 1;
+}
+
+mitk::ComputeContourSetNormalsFilter::~ComputeContourSetNormalsFilter()
+{
+}
+
+void mitk::ComputeContourSetNormalsFilter::GenerateData()
+{
+ unsigned int numberOfInputs = this->GetNumberOfInputs();
+ this->CreateOutputsForAllInputs(numberOfInputs);
+
+ //Iterating over each input
+ for(unsigned int i = 0; i < numberOfInputs; i++)
+ {
+ //Getting the inputs polydata and polygons
+ Surface* currentSurface = const_cast<Surface*>( this->GetInput(i) );
+ vtkPolyData* polyData = currentSurface->GetVtkPolyData();
+
+ vtkSmartPointer<vtkCellArray> existingPolys = polyData->GetPolys();
+
+ vtkSmartPointer<vtkPoints> existingPoints = polyData->GetPoints();
+
+ existingPolys->InitTraversal();
+
+ vtkIdType* cell (NULL);
+ vtkIdType cellSize (0);
+
+ //The array that contains all the vertex normals of the current polygon
+ vtkSmartPointer<vtkDoubleArray> normals = vtkSmartPointer<vtkDoubleArray>::New();
+ normals->SetNumberOfComponents(3);
+ normals->SetNumberOfTuples(polyData->GetNumberOfPoints());
+
+ //If the current contour is an inner contour then the direction is -1
+ //A contour lies inside another one if the pixel values in the direction of the normal is 1
+ m_NegativeNormalCounter = 0;
+ m_PositiveNormalCounter = 0;
+
+ //Iterating over each polygon
+ for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
+ {
+ if(cellSize < 3)continue;
+
+ //First we calculate the current polygon's normal
+ double polygonNormal[3] = {0.0};
+
+ double p1[3];
+ double p2[3];
+
+ double v1[3];
+ double v2[3];
+
+ existingPoints->GetPoint(cell[0], p1);
+ unsigned int index = cellSize*0.5;
+ existingPoints->GetPoint(cell[index], p2);
+
+ v1[0] = p2[0]-p1[0];
+ v1[1] = p2[1]-p1[1];
+ v1[2] = p2[2]-p1[2];
+
+ for (unsigned int k = 2; k < cellSize; k++)
+ {
+ index = cellSize*0.25;
+ existingPoints->GetPoint(cell[index], p1);
+ index = cellSize*0.75;
+ existingPoints->GetPoint(cell[index], p2);
+
+ v2[0] = p2[0]-p1[0];
+ v2[1] = p2[1]-p1[1];
+ v2[2] = p2[2]-p1[2];
+
+ vtkMath::Cross(v1,v2,polygonNormal);
+ if (vtkMath::Norm(polygonNormal) != 0)
+ break;
+ }
+
+ vtkMath::Normalize(polygonNormal);
+
+ //Now we start computing the normal for each vertex
+
+ double vertexNormalTemp[3];
+ existingPoints->GetPoint(cell[0], p1);
+ existingPoints->GetPoint(cell[1], p2);
+
+ v1[0] = p2[0]-p1[0];
+ v1[1] = p2[1]-p1[1];
+ v1[2] = p2[2]-p1[2];
+
+ vtkMath::Cross(v1,polygonNormal,vertexNormalTemp);
+
+ vtkMath::Normalize(vertexNormalTemp);
+
+ double vertexNormal[3];
+ for (unsigned j = 0; j < cellSize-2; j++)
+ {
+ existingPoints->GetPoint(cell[j+1], p1);
+ existingPoints->GetPoint(cell[j+2], p2);
+
+ v1[0] = p2[0]-p1[0];
+ v1[1] = p2[1]-p1[1];
+ v1[2] = p2[2]-p1[2];
+
+ vtkMath::Cross(v1,polygonNormal,vertexNormal);
+
+ vtkMath::Normalize(vertexNormal);
+
+ double finalNormal[3];
+
+ finalNormal[0] = (vertexNormal[0] + vertexNormalTemp[0])*0.5;
+ finalNormal[1] = (vertexNormal[1] + vertexNormalTemp[1])*0.5;
+ finalNormal[2] = (vertexNormal[2] + vertexNormalTemp[2])*0.5;
+
+ //Here we determine the direction of the normal
+ if (j == 0 && m_SegmentationBinaryImage)
+ {
+ Point3D worldCoord;
+ worldCoord[0] = p1[0]+finalNormal[0]*m_MaxSpacing;
+ worldCoord[1] = p1[1]+finalNormal[1]*m_MaxSpacing;
+ worldCoord[2] = p1[2]+finalNormal[2]*m_MaxSpacing;
+
+ double val = m_SegmentationBinaryImage->GetPixelValueByWorldCoordinate(worldCoord);
+ if (val == 1.0)
+ {
+ ++m_PositiveNormalCounter;
+ }
+ else
+ {
+ ++m_NegativeNormalCounter;
+ }
+ }
+
+ vertexNormalTemp[0] = vertexNormal[0];
+ vertexNormalTemp[1] = vertexNormal[1];
+ vertexNormalTemp[2] = vertexNormal[2];
+
+ vtkIdType id = cell[j+1];
+ normals->SetTuple(id,finalNormal);
+ }
+
+ existingPoints->GetPoint(cell[0], p1);
+ existingPoints->GetPoint(cell[1], p2);
+
+ v1[0] = p2[0]-p1[0];
+ v1[1] = p2[1]-p1[1];
+ v1[2] = p2[2]-p1[2];
+
+ vtkMath::Cross(v1,polygonNormal,vertexNormal);
+
+ vtkMath::Normalize(vertexNormal);
+
+ vertexNormal[0] = (vertexNormal[0] + vertexNormalTemp[0])*0.5;
+ vertexNormal[1] = (vertexNormal[1] + vertexNormalTemp[1])*0.5;
+ vertexNormal[2] = (vertexNormal[2] + vertexNormalTemp[2])*0.5;
+
+ vtkIdType id = cell[0];
+ normals->SetTuple(id,vertexNormal);
+ id = cell[cellSize-1];
+ normals->SetTuple(id,vertexNormal);
+
+ int normalDirection(-1);
+
+ if(m_NegativeNormalCounter < m_PositiveNormalCounter)
+ {
+ normalDirection = 1;
+ }
+
+ for(unsigned int n = 0; n < normals->GetNumberOfTuples(); n++)
+ {
+ double normal[3];
+ normals->GetTuple(n,normal);
+ normal[0] = normalDirection*normal[0];
+ normal[1] = normalDirection*normal[1];
+ normal[2] = normalDirection*normal[2];
+ }
+
+
+ }//end for all cells
+
+ Surface::Pointer surface = this->GetOutput(i);
+ surface->GetVtkPolyData()->GetCellData()->SetNormals(normals);
+ }//end for all inputs
+
+ //Setting progressbar
+ if (this->m_UseProgressBar)
+ mitk::ProgressBar::GetInstance()->Progress(this->m_ProgressStepSize);
+}
+
+
+mitk::Surface::Pointer mitk::ComputeContourSetNormalsFilter::GetNormalsAsSurface()
+{
+ //Just for debugging:
+ vtkSmartPointer<vtkPolyData> newPolyData = vtkSmartPointer<vtkPolyData>::New();
+ vtkSmartPointer<vtkCellArray> newLines = vtkSmartPointer<vtkCellArray>::New();
+ vtkSmartPointer<vtkPoints> newPoints = vtkSmartPointer<vtkPoints>::New();
+ unsigned int idCounter (0);
+ //Debug end
+
+ for (unsigned int i = 0; i < this->GetNumberOfOutputs(); i++)
+ {
+ Surface* currentSurface = const_cast<Surface*>( this->GetOutput(i) );
+ vtkPolyData* polyData = currentSurface->GetVtkPolyData();
+
+ vtkSmartPointer<vtkDoubleArray> currentCellNormals = vtkDoubleArray::SafeDownCast(polyData->GetCellData()->GetNormals());
+
+ vtkSmartPointer<vtkCellArray> existingPolys = polyData->GetPolys();
+
+ vtkSmartPointer<vtkPoints> existingPoints = polyData->GetPoints();
+
+ existingPolys->InitTraversal();
+
+ vtkIdType* cell (NULL);
+ vtkIdType cellSize (0);
+
+ for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
+ {
+ for ( unsigned int j = 0; j < cellSize; j++ )
+ {
+ double currentNormal[3];
+ currentCellNormals->GetTuple(cell[j], currentNormal);
+ vtkSmartPointer<vtkLine> line = vtkSmartPointer<vtkLine>::New();
+ line->GetPointIds()->SetNumberOfIds(2);
+ double newPoint[3];
+ double p0[3];
+ existingPoints->GetPoint(cell[j], p0);
+ newPoint[0] = p0[0] + currentNormal[0];
+ newPoint[1] = p0[1] + currentNormal[1];
+ newPoint[2] = p0[2] + currentNormal[2];
+
+ line->GetPointIds()->SetId(0, idCounter);
+ newPoints->InsertPoint(idCounter, p0);
+ idCounter++;
+
+ line->GetPointIds()->SetId(1, idCounter);
+ newPoints->InsertPoint(idCounter, newPoint);
+ idCounter++;
+
+ newLines->InsertNextCell(line);
+ }//end for all points
+ }//end for all cells
+ }//end for all outputs
+
+ newPolyData->SetPoints(newPoints);
+ newPolyData->SetLines(newLines);
+ newPolyData->BuildCells();
+
+
+ mitk::Surface::Pointer surface = mitk::Surface::New();
+ surface->SetVtkPolyData(newPolyData);
+
+ return surface;
+
+}
+
+void mitk::ComputeContourSetNormalsFilter::SetMaxSpacing(double maxSpacing)
+{
+ m_MaxSpacing = maxSpacing;
+}
+
+void mitk::ComputeContourSetNormalsFilter::GenerateOutputInformation()
+{
+ Superclass::GenerateOutputInformation();
+}
+
+void mitk::ComputeContourSetNormalsFilter::Reset()
+{
+ for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
+ {
+ this->PopBackInput();
+ }
+ this->SetNumberOfInputs(0);
+ this->SetNumberOfOutputs(0);
+}
+
+void mitk::ComputeContourSetNormalsFilter::SetUseProgressBar(bool status)
+{
+ this->m_UseProgressBar = status;
+}
+
+void mitk::ComputeContourSetNormalsFilter::SetProgressStepSize(unsigned int stepSize)
+{
+ this->m_ProgressStepSize = stepSize;
+}
diff --git a/Modules/SurfaceInterpolation/mitkComputeContourSetNormalsFilter.h b/Modules/SurfaceInterpolation/mitkComputeContourSetNormalsFilter.h
new file mode 100644
index 0000000..42daf6e
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkComputeContourSetNormalsFilter.h
@@ -0,0 +1,108 @@
+/*===================================================================
+
+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 mitkComputeContourSetNormalsFilter_h_Included
+#define mitkComputeContourSetNormalsFilter_h_Included
+
+#include "SurfaceInterpolationExports.h"
+#include "mitkSurfaceToSurfaceFilter.h"
+#include "mitkProgressBar.h"
+#include "mitkSurface.h"
+
+#include "vtkCellArray.h"
+#include "vtkPolyData.h"
+#include "vtkSmartPointer.h"
+#include "vtkDoubleArray.h"
+#include "vtkMath.h"
+#include "vtkCellData.h"
+#include "vtkLine.h"
+
+#include "mitkImage.h"
+
+namespace mitk {
+
+ /**
+ \brief Filter to compute the normales for contours based on vtkPolygons
+
+
+
+ This filter takes a number of extracted contours and computes the normals for each
+ contour edge point. The normals can be accessed by calling:
+
+ filter->GetOutput(i)->GetVtkPolyData()->GetCellData()->GetNormals();
+
+ See also the method GetNormalsAsSurface()
+
+ Note: If a segmentation binary image is provided this filter assures that the computed normals
+ do not point into the segmentation image
+
+ $Author: fetzer$
+*/
+class SurfaceInterpolation_EXPORT ComputeContourSetNormalsFilter : public SurfaceToSurfaceFilter
+{
+public:
+
+ mitkClassMacro(ComputeContourSetNormalsFilter,SurfaceToSurfaceFilter);
+ itkNewMacro(Self);
+
+ itkSetMacro(SegmentationBinaryImage, mitk::Image::Pointer);
+
+ /*
+ \brief Returns the computed normals as a surface
+ */
+ mitk::Surface::Pointer GetNormalsAsSurface();
+
+ //Resets the filter, i.e. removes all inputs and outputs
+ void Reset();
+
+ void SetMaxSpacing(double);
+
+ /**
+ \brief Set whether the mitkProgressBar should be used
+
+ \a Parameter true for using the progress bar, false otherwise
+ */
+ void SetUseProgressBar(bool);
+
+ /**
+ \brief Set the stepsize which the progress bar should proceed
+
+ \a Parameter The stepsize for progressing
+ */
+ void SetProgressStepSize(unsigned int stepSize);
+
+protected:
+ ComputeContourSetNormalsFilter();
+ virtual ~ComputeContourSetNormalsFilter();
+ virtual void GenerateData();
+ virtual void GenerateOutputInformation();
+
+private:
+
+ //The segmentation out of which the contours were extracted. Can be used to determine the direction of the normals
+ mitk::Image::Pointer m_SegmentationBinaryImage;
+ double m_MaxSpacing;
+
+ unsigned int m_NegativeNormalCounter;
+ unsigned int m_PositiveNormalCounter;
+
+ bool m_UseProgressBar;
+ unsigned int m_ProgressStepSize;
+
+};//class
+
+}//namespace
+#endif
diff --git a/Modules/SurfaceInterpolation/mitkCreateDistanceImageFromSurfaceFilter.cpp b/Modules/SurfaceInterpolation/mitkCreateDistanceImageFromSurfaceFilter.cpp
new file mode 100644
index 0000000..bd29ae0
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkCreateDistanceImageFromSurfaceFilter.cpp
@@ -0,0 +1,509 @@
+/*===================================================================
+
+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.
+
+===================================================================*/
+
+#include "mitkCreateDistanceImageFromSurfaceFilter.h"
+
+
+mitk::CreateDistanceImageFromSurfaceFilter::CreateDistanceImageFromSurfaceFilter()
+{
+ m_DistanceImageVolume = 50000;
+ this->m_UseProgressBar = false;
+ this->m_ProgressStepSize = 5;
+}
+
+
+mitk::CreateDistanceImageFromSurfaceFilter::~CreateDistanceImageFromSurfaceFilter()
+{
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::GenerateData()
+{
+ //First of all we have to build the equation-system from the existing contour-edge-points
+ this->CreateSolutionMatrixAndFunctionValues();
+
+ //Then we solve the equation-system via QR - decomposition. The interpolation weights are obtained in that way
+ vnl_qr<double> solver (m_SolutionMatrix);
+ m_Weights = solver.solve(m_FunctionValues);
+
+ //Setting progressbar
+ if (this->m_UseProgressBar)
+ mitk::ProgressBar::GetInstance()->Progress(2);
+
+ //The last step is to create the distance map with the interpolated distance function
+ this->CreateDistanceImage();
+ m_Centers.clear();
+ m_FunctionValues.clear();
+ m_Normals.clear();
+ m_Weights.clear();
+ m_SolutionMatrix.clear();
+
+ //Setting progressbar
+ if (this->m_UseProgressBar)
+ mitk::ProgressBar::GetInstance()->Progress(3);
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::CreateSolutionMatrixAndFunctionValues()
+{
+ unsigned int numberOfInputs = this->GetNumberOfInputs();
+
+ if (numberOfInputs == 0)
+ {
+ MITK_ERROR << "mitk::CreateDistanceImageFromSurfaceFilter: No input available. Please set an input!" << std::endl;
+ itkExceptionMacro("mitk::CreateDistanceImageFromSurfaceFilter: No input available. Please set an input!");
+ return;
+ }
+
+ //First of all we have to extract the nomals and the surface points.
+ //Duplicated points can be eliminated
+
+ Surface* currentSurface;
+ vtkSmartPointer<vtkPolyData> polyData;
+ vtkSmartPointer<vtkDoubleArray> currentCellNormals;
+ vtkSmartPointer<vtkCellArray> existingPolys;
+ vtkSmartPointer<vtkPoints> existingPoints;
+
+ double p[3];
+ PointType currentPoint;
+ PointType normal;
+
+ for (unsigned int i = 0; i < numberOfInputs; i++)
+ {
+ currentSurface = const_cast<Surface*>( this->GetInput(i) );
+ polyData = currentSurface->GetVtkPolyData();
+
+ if (polyData->GetNumberOfPolys() == 0)
+ {
+ MITK_INFO << "mitk::CreateDistanceImageFromSurfaceFilter: No input-polygons available. Please be sure the input surface consists of polygons!" << std::endl;
+ }
+
+ currentCellNormals = vtkDoubleArray::SafeDownCast(polyData->GetCellData()->GetNormals());
+
+ existingPolys = polyData->GetPolys();
+
+ existingPoints = polyData->GetPoints();
+
+ existingPolys->InitTraversal();
+
+ vtkIdType* cell (NULL);
+ vtkIdType cellSize (0);
+
+ for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
+ {
+ for ( unsigned int j = 0; j < cellSize; j++ )
+ {
+ existingPoints->GetPoint(cell[j], p);
+
+ currentPoint.copy_in(p);
+
+ int count = std::count(m_Centers.begin() ,m_Centers.end(),currentPoint);
+
+ if (count == 0)
+ {
+ double currentNormal[3];
+ currentCellNormals->GetTuple(cell[j], currentNormal);
+
+ normal.copy_in(currentNormal);
+
+ m_Normals.push_back(normal);
+
+ m_Centers.push_back(currentPoint);
+ }
+
+ }//end for all points
+ }//end for all cells
+ }//end for all outputs
+
+ //For we can now calculate the exact size of the centers we initialize the data structures
+ unsigned int numberOfCenters = m_Centers.size();
+ m_Centers.reserve(numberOfCenters*3);
+
+ m_FunctionValues.set_size(numberOfCenters*3);
+ m_FunctionValues.fill(0);
+
+ //Create inner points
+ for (unsigned int i = 0; i < numberOfCenters; i++)
+ {
+ currentPoint = m_Centers.at(i);
+ normal = m_Normals.at(i);
+
+ currentPoint[0] = currentPoint[0] - normal[0];
+ currentPoint[1] = currentPoint[1] - normal[1];
+ currentPoint[2] = currentPoint[2] - normal[2];
+
+ m_Centers.push_back(currentPoint);
+
+ m_FunctionValues.put(numberOfCenters+i, -1);
+ }
+
+ //Create outer points
+ for (unsigned int i = 0; i < numberOfCenters; i++)
+ {
+ currentPoint = m_Centers.at(i);
+ normal = m_Normals.at(i);
+
+ currentPoint[0] = currentPoint[0] + normal[0];
+ currentPoint[1] = currentPoint[1] + normal[1];
+ currentPoint[2] = currentPoint[2] + normal[2];
+
+ m_Centers.push_back(currentPoint);
+
+ m_FunctionValues.put(numberOfCenters*2+i, 1);
+ }
+
+ //Now we have created all centers and all function values. Next step is to create the solution matrix
+ numberOfCenters = m_Centers.size();
+ m_SolutionMatrix.set_size(numberOfCenters, numberOfCenters);
+
+ m_Weights.set_size(numberOfCenters);
+
+ PointType p1;
+ PointType p2;
+ double norm;
+
+ for (unsigned int i = 0; i < numberOfCenters; i++)
+ {
+ for (unsigned int j = 0; j < numberOfCenters; j++)
+ {
+ //Calculate the RBF value. Currently using Phi(r) = r with r is the euclidian distance between two points
+ p1 = m_Centers.at(i);
+ p2 = m_Centers.at(j);
+ p1 = p1 - p2;
+ norm = p1.two_norm();
+ m_SolutionMatrix(i,j) = norm;
+
+ }
+ }
+
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::CreateDistanceImage()
+{
+ typedef itk::Image<double, 3> DistanceImageType;
+ typedef itk::ImageRegionIteratorWithIndex<DistanceImageType> ImageIterator;
+ typedef itk::NeighborhoodIterator<DistanceImageType> NeighborhoodImageIterator;
+
+ DistanceImageType::Pointer distanceImg = DistanceImageType::New();
+
+ //Determin the bounding box of the delineated contours
+ double xmin = m_Centers.at(0)[0];
+ double ymin = m_Centers.at(0)[1];
+ double zmin = m_Centers.at(0)[2];
+ double xmax = m_Centers.at(0)[0];
+ double ymax = m_Centers.at(0)[1];
+ double zmax = m_Centers.at(0)[2];
+
+ for (unsigned int i = 1; i < m_Centers.size(); i++)
+ {
+ if (xmin > m_Centers.at(i)[0])
+ {
+ xmin = m_Centers.at(i)[0];
+ }
+ if (ymin > m_Centers.at(i)[1])
+ {
+ ymin = m_Centers.at(i)[1];
+ }
+ if (zmin > m_Centers.at(i)[2])
+ {
+ zmin = m_Centers.at(i)[2];
+ }
+ if (xmax < m_Centers.at(i)[0])
+ {
+ xmax = m_Centers.at(i)[0];
+ }
+ if (ymax < m_Centers.at(i)[1])
+ {
+ ymax = m_Centers.at(i)[1];
+ }
+ if (zmax < m_Centers.at(i)[2])
+ {
+ zmax = m_Centers.at(i)[2];
+ }
+ }
+
+ Vector3D extentMM;
+ extentMM[0] = xmax - xmin + 5;
+ extentMM[1] = ymax - ymin + 5;
+ extentMM[2] = zmax - zmin + 5;
+
+ //Shifting the distance image's offest to achieve an exact distance calculation
+ xmin = xmin - 5;
+ ymin = ymin - 5;
+ zmin = zmin - 5;
+
+ /*
+ Now create an empty distance image. The create image will always have the same size, independent from
+ the original image (e.g. always consists of 500000 pixels) and will have an isotropic spacing.
+ The spacing is calculated like the following:
+ The image's volume = 500000 Pixels = extentX*spacing*extentY*spacing*extentZ*spacing
+ So the spacing is: spacing = ( 500000 / extentX*extentY*extentZ )^(1/3)
+ */
+
+ double basis = (extentMM[0]*extentMM[1]*extentMM[2]) / m_DistanceImageVolume;
+ double exponent = 1.0/3.0;
+ double distImgSpacing = pow(basis, exponent);
+ int tempSpacing = (distImgSpacing+0.05)*10;
+ m_DistanceImageSpacing = (double)tempSpacing/10.0;
+
+ unsigned int numberOfXPixel = extentMM[0] / m_DistanceImageSpacing;
+ unsigned int numberOfYPixel = extentMM[1] / m_DistanceImageSpacing;
+ unsigned int numberOfZPixel = extentMM[2] / m_DistanceImageSpacing;
+
+ DistanceImageType::SizeType size;
+
+ //Increase the distance image's size a little bit to achieve an exact distance calculation
+ size[0] = numberOfXPixel + 5;
+ size[1] = numberOfYPixel + 5;
+ size[2] = numberOfZPixel + 5;
+
+ DistanceImageType::IndexType start;
+ start[0] = 0;
+ start[1] = 0;
+ start[2] = 0;
+
+ DistanceImageType::RegionType lpRegion;
+
+ lpRegion.SetSize(size);
+ lpRegion.SetIndex(start);
+
+ distanceImg->SetRegions( lpRegion );
+ distanceImg->SetSpacing( m_DistanceImageSpacing );
+ distanceImg->Allocate();
+
+ //First of all the image is initialized with the value 10 for each pixel
+ distanceImg->FillBuffer(10);
+
+ /*
+ Now we must caculate the distance for each pixel. But instead of calculating the distance value
+ for all of the image's pixels we proceed similar to the region growing algorithm:
+
+ 1. Take the first pixel from the narrowband_point_list and calculate the distance for each neighbor (6er)
+ 2. If the current index's distance value is below a certain threshold push it into the list
+ 3. Next iteration take the next index from the list and start with 1. again
+
+ This is done until the narrowband_point_list is empty.
+ */
+ std::queue<DistanceImageType::IndexType> narrowbandPoints;
+ PointType currentPoint = m_Centers.at(0);
+ double distance = this->CalculateDistanceValue(currentPoint);
+
+ DistanceImageType::IndexType currentIndex;
+ currentIndex[0] = ( currentPoint[0]-xmin ) / m_DistanceImageSpacing;
+ currentIndex[1] = ( currentPoint[1]-ymin ) / m_DistanceImageSpacing;
+ currentIndex[2] = ( currentPoint[2]-zmin ) / m_DistanceImageSpacing;
+
+ narrowbandPoints.push(currentIndex);
+ distanceImg->SetPixel(currentIndex, distance);
+
+ NeighborhoodImageIterator::RadiusType radius;
+ radius.Fill(1);
+ NeighborhoodImageIterator nIt(radius, distanceImg, distanceImg->GetLargestPossibleRegion());
+ unsigned int relativeNbIdx[] = {4, 10, 12, 14, 16, 22};
+
+ bool isInBounds = false;
+
+ while ( !narrowbandPoints.empty() )
+ {
+
+ nIt.SetLocation(narrowbandPoints.front());
+ narrowbandPoints.pop();
+
+ for (int i = 0; i < 6; i++)
+ {
+ nIt.GetPixel(relativeNbIdx[i], isInBounds);
+ if( isInBounds && nIt.GetPixel(relativeNbIdx[i]) == 10)
+ {
+ currentIndex = nIt.GetIndex(relativeNbIdx[i]);
+
+ currentPoint[0] = currentIndex[0]*m_DistanceImageSpacing + xmin;
+ currentPoint[1] = currentIndex[1]*m_DistanceImageSpacing + ymin;
+ currentPoint[2] = currentIndex[2]*m_DistanceImageSpacing + zmin;
+
+ distance = this->CalculateDistanceValue(currentPoint);
+ if ( abs(distance) <= m_DistanceImageSpacing*2 )
+ {
+ nIt.SetPixel(relativeNbIdx[i], distance);
+ narrowbandPoints.push(currentIndex);
+ }
+ }
+ }
+ }
+
+ ImageIterator imgRegionIterator (distanceImg, distanceImg->GetLargestPossibleRegion());
+ imgRegionIterator.GoToBegin();
+
+ double prevPixelVal = 1;
+
+ unsigned int _size[3] = { (unsigned int)(size[0] - 1), (unsigned int)(size[1] - 1), (unsigned int)(size[2] - 1) };
+
+ //Set every pixel inside the surface to -10 except the edge point (so that the received surface is closed)
+ while (!imgRegionIterator.IsAtEnd()) {
+
+ if ( imgRegionIterator.Get() == 10 && prevPixelVal < 0 )
+ {
+
+ while (imgRegionIterator.Get() == 10)
+ {
+ if (imgRegionIterator.GetIndex()[0] == _size[0] || imgRegionIterator.GetIndex()[1] == _size[1] || imgRegionIterator.GetIndex()[2] == _size[2]
+ || imgRegionIterator.GetIndex()[0] == 0U || imgRegionIterator.GetIndex()[1] == 0U || imgRegionIterator.GetIndex()[2] == 0U )
+ {
+ imgRegionIterator.Set(10);
+ prevPixelVal = 10;
+ ++imgRegionIterator;
+ break;
+ }
+ else
+ {
+ imgRegionIterator.Set(-10);
+ ++imgRegionIterator;
+ prevPixelVal = -10;
+ }
+
+ }
+
+ }
+ else if (imgRegionIterator.GetIndex()[0] == _size[0] || imgRegionIterator.GetIndex()[1] == _size[1] || imgRegionIterator.GetIndex()[2] == _size[2]
+ || imgRegionIterator.GetIndex()[0] == 0U || imgRegionIterator.GetIndex()[1] == 0U || imgRegionIterator.GetIndex()[2] == 0U)
+ {
+ imgRegionIterator.Set(10);
+ prevPixelVal = 10;
+ ++imgRegionIterator;
+ }
+ else {
+ prevPixelVal = imgRegionIterator.Get();
+ ++imgRegionIterator;
+ }
+
+ }
+
+ Image::Pointer resultImage = this->GetOutput();
+
+ Point3D origin;
+ origin[0] = xmin;
+ origin[1] = ymin;
+ origin[2] = zmin;
+
+ CastToMitkImage(distanceImg, resultImage);
+ resultImage->GetGeometry()->SetOrigin(origin);
+}
+
+
+double mitk::CreateDistanceImageFromSurfaceFilter::CalculateDistanceValue(PointType p)
+{
+ double distanceValue (0);
+ PointType p1;
+ PointType p2;
+ double norm;
+
+ for (unsigned int i = 0; i < m_Centers.size(); i++)
+ {
+ p1 = m_Centers.at(i);
+ p2 = p-p1;
+ norm = p2.two_norm();
+ distanceValue = distanceValue + norm*m_Weights.get(i);
+ }
+ return distanceValue;
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::GenerateOutputInformation()
+{
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::PrintEquationSystem()
+{
+ std::ofstream esfile;
+ esfile.open("C:/Users/fetzer/Desktop/equationSystem/es.txt");
+ esfile<<"Nummber of rows: "<<m_SolutionMatrix.rows()<<" ****** Number of columns: "<<m_SolutionMatrix.columns()<<endl;
+ esfile<<"[ ";
+ for (unsigned int i = 0; i < m_SolutionMatrix.rows(); i++)
+ {
+ for (unsigned int j = 0; j < m_SolutionMatrix.columns(); j++)
+ {
+ esfile<<m_SolutionMatrix(i,j)<<" ";
+ }
+ esfile<<";"<<endl;
+ }
+ esfile<<" ]";
+ esfile.close();
+
+ std::ofstream centersFile;
+ centersFile.open("C:/Users/fetzer/Desktop/equationSystem/centers.txt");
+ for (unsigned int i = 0; i < m_Centers.size(); i++)
+ {
+ centersFile<<m_Centers.at(i)<<";"<<endl;
+ }
+ centersFile.close();
+
+}
+
+
+
+void mitk::CreateDistanceImageFromSurfaceFilter::SetInput( const mitk::Surface* surface )
+{
+ this->SetInput( 0, const_cast<mitk::Surface*>( surface ) );
+}
+
+
+void mitk::CreateDistanceImageFromSurfaceFilter::SetInput( unsigned int idx, const mitk::Surface* surface )
+{
+ if ( this->GetInput(idx) != surface )
+ {
+ this->SetNthInput( idx, const_cast<mitk::Surface*>( surface ) );
+ this->Modified();
+ }
+}
+
+
+const mitk::Surface* mitk::CreateDistanceImageFromSurfaceFilter::GetInput()
+{
+ if (this->GetNumberOfInputs() < 1)
+ return NULL;
+
+ return static_cast<const mitk::Surface*>(this->ProcessObject::GetInput(0));
+}
+
+
+const mitk::Surface* mitk::CreateDistanceImageFromSurfaceFilter::GetInput( unsigned int idx)
+{
+ if (this->GetNumberOfInputs() < 1)
+ return NULL;
+
+ return static_cast<const mitk::Surface*>(this->ProcessObject::GetInput(idx));
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::RemoveInputs(mitk::Surface* input)
+{
+ this->RemoveInput(input);
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::Reset()
+{
+ for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
+ {
+ this->PopBackInput();
+ }
+ this->SetNumberOfInputs(0);
+ this->SetNumberOfOutputs(1);
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::SetUseProgressBar(bool status)
+{
+ this->m_UseProgressBar = status;
+}
+
+void mitk::CreateDistanceImageFromSurfaceFilter::SetProgressStepSize(unsigned int stepSize)
+{
+ this->m_ProgressStepSize = stepSize;
+}
diff --git a/Modules/SurfaceInterpolation/mitkCreateDistanceImageFromSurfaceFilter.h b/Modules/SurfaceInterpolation/mitkCreateDistanceImageFromSurfaceFilter.h
new file mode 100644
index 0000000..04c39b4
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkCreateDistanceImageFromSurfaceFilter.h
@@ -0,0 +1,161 @@
+/*===================================================================
+
+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 mitkCreateDistanceImageFromSurfaceFilter_h_Included
+#define mitkCreateDistanceImageFromSurfaceFilter_h_Included
+
+#include "SurfaceInterpolationExports.h"
+
+#include "mitkImageSource.h"
+#include "mitkSurface.h"
+#include "mitkProgressBar.h"
+
+#include "vtkSmartPointer.h"
+#include "vtkDoubleArray.h"
+#include "vtkCellArray.h"
+#include "vtkCellData.h"
+#include "vtkPolyData.h"
+
+#include "vnl/vnl_matrix.h"
+#include "vnl/vnl_vector.h"
+#include "vnl/vnl_vector_fixed.h"
+#include "vnl/algo/vnl_qr.h"
+
+#include "itkImage.h"
+#include "itkImageRegionIteratorWithIndex.h"
+#include "itkNeighborhoodIterator.h"
+
+#include <queue>
+
+namespace mitk {
+
+ /**
+ \brief This filter interpolates the 3D surface for a segmented area. The basis for the interpolation
+ are the edge-points of contours that are drawn into an image.
+
+ The interpolation itself is performed via Radial Basis Function Interpolation.
+
+ ATTENTION:
+ This filter needs beside the edge points of the delineated contours additionally the normals for each
+ edge point.
+
+ \sa mitkSurfaceInterpolationController
+
+ Based on the contour edge points and their normal this filter calculates a distance function with the following
+ properties:
+ - Putting a point into the distance function that lies inside the considered surface gives a negativ scalar value
+ - Putting a point into the distance function that lies outside the considered surface gives a positive scalar value
+ - Putting a point into the distance function that lies exactly on the considered surface gives the value zero
+
+ With this interpolated distance function a distance image will be created. The desired surface can then be extract e.g.
+ with the marching cubes algorithm. (Within the distance image the surface goes exactly where the pixelvalues are zero)
+
+ Note that the obtained distance image has always an isotropig spacing. The size (in this case volume) of the image can be
+ adjusted by calling SetDistanceImageVolume(unsigned int volume) which specifies the number ob pixels enclosed by the image.
+
+ \ingroup Process
+
+ $Author: fetzer$
+ */
+ class SurfaceInterpolation_EXPORT CreateDistanceImageFromSurfaceFilter : public ImageSource
+ {
+
+ public:
+
+ typedef vnl_vector_fixed<double,3> PointType;
+
+ typedef std::vector< PointType > NormalList;
+ typedef std::vector< PointType > CenterList;
+
+ typedef vnl_matrix<double> SolutionMatrix;
+ typedef vnl_vector<double> FunctionValues;
+ typedef vnl_vector<double> InterpolationWeights;
+
+ typedef std::vector<Surface::Pointer> SurfaceList;
+
+ mitkClassMacro(CreateDistanceImageFromSurfaceFilter,ImageSource);
+ itkNewMacro(Self);
+
+ //Methods copied from mitkSurfaceToSurfaceFilter
+ virtual void SetInput( const mitk::Surface* surface );
+
+ virtual void SetInput( unsigned int idx, const mitk::Surface* surface );
+
+ virtual const mitk::Surface* GetInput();
+
+ virtual const mitk::Surface* GetInput( unsigned int idx );
+
+ virtual void RemoveInputs(mitk::Surface* input);
+
+
+ /**
+ \brief Set the size of the output distance image. The size is specified by the image's volume
+ (i.e. in this case how many pixels are enclosed by the image)
+ If non is set, the volume will be 500000 pixels.
+ */
+ itkSetMacro(DistanceImageVolume, unsigned int);
+
+ void PrintEquationSystem();
+
+ //Resets the filter, i.e. removes all inputs and outputs
+ void Reset();
+
+ /**
+ \brief Set whether the mitkProgressBar should be used
+
+ \a Parameter true for using the progress bar, false otherwise
+ */
+ void SetUseProgressBar(bool);
+
+ /**
+ \brief Set the stepsize which the progress bar should proceed
+
+ \a Parameter The stepsize for progressing
+ */
+ void SetProgressStepSize(unsigned int stepSize);
+
+ protected:
+ CreateDistanceImageFromSurfaceFilter();
+ virtual ~CreateDistanceImageFromSurfaceFilter();
+ virtual void GenerateData();
+ virtual void GenerateOutputInformation();
+
+
+ private:
+
+ void CreateSolutionMatrixAndFunctionValues();
+ double CalculateDistanceValue(PointType p);
+
+ void CreateDistanceImage ();
+
+ //Datastructures for the interpolation
+ CenterList m_Centers;
+ NormalList m_Normals;
+ FunctionValues m_FunctionValues;
+ InterpolationWeights m_Weights;
+ SolutionMatrix m_SolutionMatrix;
+ double m_DistanceImageSpacing;
+
+ unsigned int m_DistanceImageVolume;
+
+ bool m_UseProgressBar;
+ unsigned int m_ProgressStepSize;
+};
+
+}//namespace
+
+
+#endif
diff --git a/Modules/SurfaceInterpolation/mitkReduceContourSetFilter.cpp b/Modules/SurfaceInterpolation/mitkReduceContourSetFilter.cpp
new file mode 100644
index 0000000..8332034
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkReduceContourSetFilter.cpp
@@ -0,0 +1,492 @@
+/*===================================================================
+
+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.
+
+===================================================================*/
+
+#include "mitkReduceContourSetFilter.h"
+
+mitk::ReduceContourSetFilter::ReduceContourSetFilter()
+{
+ m_MaxSegmentLenght = 0;
+ m_StepSize = 10;
+ m_Tolerance = -1;
+ m_ReductionType = DOUGLAS_PEUCKER;
+ m_MaxSpacing = -1;
+ m_MinSpacing = -1;
+ this->m_UseProgressBar = false;
+ this->m_ProgressStepSize = 1;
+ m_NumberOfPointsAfterReduction = 0;
+}
+
+mitk::ReduceContourSetFilter::~ReduceContourSetFilter()
+{
+}
+
+void mitk::ReduceContourSetFilter::GenerateData()
+{
+ unsigned int numberOfInputs = this->GetNumberOfInputs();
+ unsigned int numberOfOutputs (0);
+
+ vtkSmartPointer<vtkPolyData> newPolyData;
+ vtkSmartPointer<vtkCellArray> newPolygons;
+ vtkSmartPointer<vtkPoints> newPoints;
+
+ //For the purpose of evaluation
+// unsigned int numberOfPointsBefore (0);
+ m_NumberOfPointsAfterReduction=0;
+
+ for(unsigned int i = 0; i < numberOfInputs; i++)
+ {
+ mitk::Surface* currentSurface = const_cast<mitk::Surface*>( this->GetInput(i) );
+ vtkSmartPointer<vtkPolyData> polyData = currentSurface->GetVtkPolyData();
+
+ newPolyData = vtkPolyData::New();
+ newPolygons = vtkCellArray::New();
+ newPoints = vtkPoints::New();
+
+ vtkSmartPointer<vtkCellArray> existingPolys = polyData->GetPolys();
+
+ vtkSmartPointer<vtkPoints> existingPoints = polyData->GetPoints();
+
+ existingPolys->InitTraversal();
+
+ vtkIdType* cell (NULL);
+ vtkIdType cellSize (0);
+
+ for( existingPolys->InitTraversal(); existingPolys->GetNextCell(cellSize, cell);)
+ {
+ bool incorporatePolygon = this->CheckForIntersection(cell,cellSize,existingPoints, /*numberOfIntersections, intersectionPoints, */i);
+ if ( !incorporatePolygon ) continue;
+
+ vtkSmartPointer<vtkPolygon> newPolygon = vtkSmartPointer<vtkPolygon>::New();
+
+ if(m_ReductionType == NTH_POINT)
+ {
+ this->ReduceNumberOfPointsByNthPoint(cellSize, cell, existingPoints, newPolygon, newPoints);
+ if (newPolygon->GetPointIds()->GetNumberOfIds() != 0)
+ {
+ newPolygons->InsertNextCell(newPolygon);
+ }
+ }
+ else if (m_ReductionType == DOUGLAS_PEUCKER)
+ {
+ this->ReduceNumberOfPointsByDouglasPeucker(cellSize, cell, existingPoints, newPolygon, newPoints);
+ if (newPolygon->GetPointIds()->GetNumberOfIds() > 3)
+ {
+ newPolygons->InsertNextCell(newPolygon);
+ }
+ }
+
+ //Again for evaluation
+// numberOfPointsBefore += cellSize;
+ m_NumberOfPointsAfterReduction += newPolygon->GetPointIds()->GetNumberOfIds();
+
+ }
+
+ if (newPolygons->GetNumberOfCells() != 0)
+ {
+ newPolyData->SetPolys(newPolygons);
+ newPolyData->SetPoints(newPoints);
+ newPolyData->BuildLinks();
+
+ Surface::Pointer surface = this->GetOutput(numberOfOutputs);
+ surface->SetVtkPolyData(newPolyData);
+ numberOfOutputs++;
+ }
+
+ }
+// MITK_INFO<<"Points before: "<<numberOfPointsBefore<<" ##### Points after: "<<numberOfPointsAfter;
+ this->SetNumberOfOutputs(numberOfOutputs);
+
+ //Setting progressbar
+ if (this->m_UseProgressBar)
+ mitk::ProgressBar::GetInstance()->Progress(this->m_ProgressStepSize);
+}
+
+void mitk::ReduceContourSetFilter::ReduceNumberOfPointsByNthPoint (vtkIdType cellSize, vtkIdType* cell, vtkPoints* points, vtkPolygon* reducedPolygon, vtkPoints* reducedPoints)
+{
+
+ unsigned int newNumberOfPoints (0);
+ unsigned int mod = cellSize%m_StepSize;
+
+ if(mod == 0)
+ {
+ newNumberOfPoints = cellSize/m_StepSize;
+ }
+ else
+ {
+ newNumberOfPoints = ( (cellSize-mod)/m_StepSize )+1;
+ }
+
+ if (newNumberOfPoints <= 3)
+ {
+ return;
+ }
+ reducedPolygon->GetPointIds()->SetNumberOfIds(newNumberOfPoints);
+ reducedPolygon->GetPoints()->SetNumberOfPoints(newNumberOfPoints);
+
+ for (unsigned int i = 0; i < cellSize; i++)
+ {
+ if (i%m_StepSize == 0)
+ {
+ double point[3];
+ points->GetPoint(cell[i], point);
+ vtkIdType id = reducedPoints->InsertNextPoint(point);
+ reducedPolygon->GetPointIds()->SetId(i/m_StepSize, id);
+ }
+
+ }
+ vtkIdType id = cell[0];
+ double point[3];
+ points->GetPoint(id, point);
+ id = reducedPoints->InsertNextPoint(point);
+ reducedPolygon->GetPointIds()->SetId(newNumberOfPoints-1, id);
+
+}
+
+void mitk::ReduceContourSetFilter::ReduceNumberOfPointsByDouglasPeucker(vtkIdType cellSize, vtkIdType* cell, vtkPoints* points,
+ vtkPolygon* reducedPolygon, vtkPoints* reducedPoints)
+{
+ //If the cell is too small to obtain a reduced polygon with the given stepsize return
+ if (cellSize <= m_StepSize*3)return;
+
+ /*
+ What we do now is (see the Douglas Peucker Algorithm):
+
+ 1. Divide the current contour in two line segments (start - middle; middle - end), put them into the stack
+
+ 2. Fetch first line segment and create the following vectors:
+
+ - v1 = (start;end)
+ - v2 = (start;currentPoint) -> for each point of the current line segment!
+
+ 3. Calculate the distance from the currentPoint to v1:
+
+ a. Determine the length of the orthogonal projection of v2 to v1 by:
+
+ l = v2 * (normalized v1)
+
+ b. There a three possibilities for the distance then:
+
+ d = sqrt(lenght(v2)^2 - l^2) if l > 0 and l < length(v1)
+ d = lenght(v2-v1) if l > 0 and l > lenght(v1)
+ d = length(v2) if l < 0 because v2 is then pointing in a different direction than v1
+
+ 4. Memorize the point with the biggest distance and create two new line segments with it at the end of the iteration
+ and put it into the stack
+
+ 5. If the distance value D <= m_Tolerance, then add the start and end index and the corresponding points to the reduced ones
+ */
+
+ //First of all set tolerance if none is specified
+ if(m_Tolerance < 0)
+ {
+ if(m_MaxSpacing > 0)
+ {
+ m_Tolerance = m_MinSpacing;
+ }
+ else
+ {
+ m_Tolerance = 1.5;
+ }
+ }
+
+ std::stack<LineSegment> lineSegments;
+
+ //1. Divide in line segments
+
+ LineSegment ls2;
+ ls2.StartIndex = cell[cellSize/2];
+ ls2.EndIndex = cell[cellSize-1];
+ lineSegments.push(ls2);
+
+ LineSegment ls1;
+ ls1.StartIndex = cell[0];
+ ls1.EndIndex = cell[cellSize/2];
+ lineSegments.push(ls1);
+
+ LineSegment currentSegment;
+ double v1[3];
+ double v2[3];
+ double tempV[3];
+ double lenghtV1;
+
+ double currentMaxDistance (0);
+ vtkIdType currentMaxDistanceIndex (0);
+
+ double l;
+ double d;
+
+ vtkIdType pointId (0);
+ //Add the start index to the reduced points. From now on just the end indices will be added
+ pointId = reducedPoints->InsertNextPoint(points->GetPoint(cell[0]));
+ reducedPolygon->GetPointIds()->InsertNextId(pointId);
+
+ while (!lineSegments.empty())
+ {
+ currentSegment = lineSegments.top();
+ lineSegments.pop();
+
+ //2. Create vectors
+
+ points->GetPoint(currentSegment.EndIndex, tempV);
+ points->GetPoint(currentSegment.StartIndex, v1);
+
+ v1[0] = tempV[0]-v1[0];
+ v1[1] = tempV[1]-v1[1];
+ v1[2] = tempV[2]-v1[2];
+
+ lenghtV1 = vtkMath::Norm(v1);
+
+ vtkMath::Normalize(v1);
+ int range = currentSegment.EndIndex - currentSegment.StartIndex;
+ for (int i = 1; i < abs(range); ++i)
+ {
+ points->GetPoint(currentSegment.StartIndex+i, tempV);
+ points->GetPoint(currentSegment.StartIndex, v2);
+
+ v2[0] = tempV[0]-v2[0];
+ v2[1] = tempV[1]-v2[1];
+ v2[2] = tempV[2]-v2[2];
+
+ //3. Calculate the distance
+
+ l = vtkMath::Dot(v2, v1);
+
+ d = vtkMath::Norm(v2);
+
+ if (l > 0 && l < lenghtV1)
+ {
+ d = sqrt((d*d-l*l));
+ }
+ else if (l > 0 && l > lenghtV1)
+ {
+ tempV[0] = lenghtV1*v1[0] - v2[0];
+ tempV[1] = lenghtV1*v1[1] - v2[1];
+ tempV[2] = lenghtV1*v1[2] - v2[2];
+ d = vtkMath::Norm(tempV);
+ }
+
+ //4. Memorize maximum distance
+ if (d > currentMaxDistance)
+ {
+ currentMaxDistance = d;
+ currentMaxDistanceIndex = currentSegment.StartIndex+i;
+ }
+
+ }
+
+ //4. & 5.
+ if (currentMaxDistance <= m_Tolerance)
+ {
+
+ //double temp[3];
+ int segmentLenght = currentSegment.EndIndex - currentSegment.StartIndex;
+ if (segmentLenght > (int)m_MaxSegmentLenght)
+ {
+ m_MaxSegmentLenght = (unsigned int)segmentLenght;
+ }
+
+// MITK_INFO<<"Lenght: "<<abs(segmentLenght);
+ if (abs(segmentLenght) > 25)
+ {
+ unsigned int newLenght(segmentLenght);
+ while (newLenght > 25)
+ {
+ newLenght = newLenght*0.5;
+ }
+ unsigned int divisions = abs(segmentLenght)/newLenght;
+// MITK_INFO<<"Divisions: "<<divisions;
+
+ for (unsigned int i = 1; i<=divisions; ++i)
+ {
+// MITK_INFO<<"Inserting MIDDLE: "<<(currentSegment.StartIndex + newLenght*i);
+ pointId = reducedPoints->InsertNextPoint(points->GetPoint(currentSegment.StartIndex + newLenght*i));
+ reducedPolygon->GetPointIds()->InsertNextId(pointId);
+ }
+ }
+// MITK_INFO<<"Inserting END: "<<currentSegment.EndIndex;
+ pointId = reducedPoints->InsertNextPoint(points->GetPoint(currentSegment.EndIndex));
+ reducedPolygon->GetPointIds()->InsertNextId(pointId);
+ }
+ else
+ {
+ ls2.StartIndex = currentMaxDistanceIndex;
+ ls2.EndIndex = currentSegment.EndIndex;
+ lineSegments.push(ls2);
+
+ ls1.StartIndex = currentSegment.StartIndex;
+ ls1.EndIndex = currentMaxDistanceIndex;
+ lineSegments.push(ls1);
+ }
+ currentMaxDistance = 0;
+
+ }
+
+}
+
+bool mitk::ReduceContourSetFilter::CheckForIntersection (vtkIdType* currentCell, vtkIdType currentCellSize, vtkPoints* currentPoints,/* vtkIdType numberOfIntersections, vtkIdType* intersectionPoints,*/ unsigned int currentInputIndex)
+{
+ /*
+ If we check the current cell for intersections then we have to consider three possibilies:
+ 1. There is another cell among all the other input surfaces which intersects the current polygon:
+ - That means we have to save the intersection points because these points should not be eliminated
+ 2. There current polygon exists just because of an intersection of another polygon with the current plane defined by the current polygon
+ - That means the current polygon should not be incorporated and all of its points should be eliminated
+ 3. There is no intersection
+ - That mean we can just reduce the current polygons points without considering any intersections
+ */
+
+ for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
+ {
+ //Don't check for intersection with the polygon itself
+ if (i == currentInputIndex) continue;
+
+ //Get the next polydata to check for intersection
+ vtkSmartPointer<vtkPolyData> poly = const_cast<Surface*>( this->GetInput(i) )->GetVtkPolyData();
+ vtkSmartPointer<vtkCellArray> polygonArray = poly->GetPolys();
+ polygonArray->InitTraversal();
+ vtkIdType anotherInputPolygonSize (0);
+ vtkIdType* anotherInputPolygonIDs(NULL);
+
+ /*
+ The procedure is:
+ - Create the equation of the plane, defined by the points of next input
+ - Calculate the distance of each point of the current polygon to the plane
+ - If the maximum distance is not bigger than 1.5 of the maximum spacing AND the minimal distance is not bigger
+ than 0.5 of the minimum spacing then the current contour is an intersection contour
+ */
+
+ for( polygonArray->InitTraversal(); polygonArray->GetNextCell(anotherInputPolygonSize, anotherInputPolygonIDs);)
+ {
+ //Choosing three plane points to calculate the plane vectors
+ double p1[3];
+ double p2[3];
+ double p3[3];
+
+ //The plane vectors
+ double v1[3];
+ double v2[3] = { 0 };
+ //The plane normal
+ double normal[3];
+
+ //Create first Vector
+ poly->GetPoint(anotherInputPolygonIDs[0], p1);
+ poly->GetPoint(anotherInputPolygonIDs[1], p2);
+
+ v1[0] = p2[0]-p1[0];
+ v1[1] = p2[1]-p1[1];
+ v1[2] = p2[2]-p1[2];
+
+ //Find 3rd point for 2nd vector (The angle between the two plane vectors should be bigger than 30 degrees)
+
+ double maxDistance (0);
+ double minDistance (10000);
+
+ for (unsigned int j = 2; j < anotherInputPolygonSize; j++)
+ {
+ poly->GetPoint(anotherInputPolygonIDs[j], p3);
+
+ v2[0] = p3[0]-p1[0];
+ v2[1] = p3[1]-p1[1];
+ v2[2] = p3[2]-p1[2];
+
+ //Calculate the angle between the two vector for the current point
+ double dotV1V2 = vtkMath::Dot(v1,v2);
+ double absV1 = sqrt(vtkMath::Dot(v1,v1));
+ double absV2 = sqrt(vtkMath::Dot(v2,v2));
+ double cosV1V2 = dotV1V2/(absV1*absV2);
+
+ double arccos = acos(cosV1V2);
+ double degree = vtkMath::DegreesFromRadians(arccos);
+
+ //If angle is bigger than 30 degrees break
+ if (degree > 30) break;
+
+ }//for (to find 3rd point)
+
+ //Calculate normal of the plane by taking the cross product of the two vectors
+ vtkMath::Cross(v1,v2,normal);
+ vtkMath::Normalize(normal);
+
+ //Determine position of the plane
+ double lambda = vtkMath::Dot(normal, p1);
+
+ /*
+ Calculate the distance to the plane for each point of the current polygon
+ If the distance is zero then save the currentPoint as intersection point
+ */
+ for (unsigned int k = 0; k < currentCellSize; k++)
+ {
+ double currentPoint[3];
+ currentPoints->GetPoint(currentCell[k], currentPoint);
+
+ double tempPoint[3];
+ tempPoint[0] = normal[0]*currentPoint[0];
+ tempPoint[1] = normal[1]*currentPoint[1];
+ tempPoint[2] = normal[2]*currentPoint[2];
+
+ double temp = tempPoint[0]+tempPoint[1]+tempPoint[2]-lambda;
+ double distance = fabs(temp);
+
+ if (distance > maxDistance)
+ {
+ maxDistance = distance;
+ }
+ if (distance < minDistance)
+ {
+ minDistance = distance;
+ }
+ }//for (to calculate distance and intersections with currentPolygon)
+
+ if (maxDistance < 1.5*m_MaxSpacing && minDistance < 0.5*m_MinSpacing)
+ {
+ return false;
+ }
+
+ //Because we are considering the plane defined by the acual input polygon only one iteration is sufficient
+ //We do not need to consider each cell of the plane
+ break;
+ }//for (to traverse through all cells of actualInputPolyData)
+
+ }//for (to iterate through all inputs)
+
+ return true;
+}
+
+void mitk::ReduceContourSetFilter::GenerateOutputInformation()
+{
+ Superclass::GenerateOutputInformation();
+}
+
+void mitk::ReduceContourSetFilter::Reset()
+{
+ for (unsigned int i = 0; i < this->GetNumberOfInputs(); i++)
+ {
+ this->PopBackInput();
+ }
+ this->SetNumberOfInputs(0);
+ this->SetNumberOfOutputs(0);
+ m_NumberOfPointsAfterReduction = 0;
+}
+
+void mitk::ReduceContourSetFilter::SetUseProgressBar(bool status)
+{
+ this->m_UseProgressBar = status;
+}
+
+void mitk::ReduceContourSetFilter::SetProgressStepSize(unsigned int stepSize)
+{
+ this->m_ProgressStepSize = stepSize;
+}
diff --git a/Modules/SurfaceInterpolation/mitkReduceContourSetFilter.h b/Modules/SurfaceInterpolation/mitkReduceContourSetFilter.h
new file mode 100644
index 0000000..8ad1b22
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkReduceContourSetFilter.h
@@ -0,0 +1,128 @@
+/*===================================================================
+
+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 mitkReduceContourSetFilter_h_Included
+#define mitkReduceContourSetFilter_h_Included
+
+#include "SurfaceInterpolationExports.h"
+#include "mitkSurfaceToSurfaceFilter.h"
+#include "mitkSurface.h"
+#include "mitkProgressBar.h"
+
+#include "vtkSmartPointer.h"
+#include "vtkPolyData.h"
+#include "vtkPolygon.h"
+#include "vtkPoints.h"
+#include "vtkCellArray.h"
+#include "vtkMath.h"
+
+#include <stack>
+
+namespace mitk {
+
+/**
+ \brief A filter that reduces the number of points of contours represented by a mitk::Surface
+
+ The type of the reduction can be set via SetReductionType. The two ways provided by this filter is the
+
+ \li NTH_POINT Algorithm which reduces the contours according to a certain stepsize
+ \li DOUGLAS_PEUCKER Algorithm which incorpates an error tolerance into the reduction.
+
+ Stepsize and error tolerance can be set via SetStepSize and SetTolerance.
+
+ Additional if more than one input contour is provided this filter tries detect contours which occur just because
+ of an intersection. These intersection contours are eliminated. In oder to ensure a correct elimination the min and max
+ spacing of the original image must be provided.
+
+ The output is a mitk::Surface.
+
+ $Author: fetzer$
+*/
+
+ class SurfaceInterpolation_EXPORT ReduceContourSetFilter : public SurfaceToSurfaceFilter
+ {
+
+ public:
+
+ enum Reduction_Type
+ {
+ NTH_POINT, DOUGLAS_PEUCKER
+ };
+
+ struct LineSegment
+ {
+ unsigned int StartIndex;
+ unsigned int EndIndex;
+ };
+
+ mitkClassMacro(ReduceContourSetFilter,SurfaceToSurfaceFilter);
+ itkNewMacro(Self);
+
+ itkSetMacro(MinSpacing, double);
+ itkSetMacro(MaxSpacing, double);
+ itkSetMacro(ReductionType, Reduction_Type);
+ itkSetMacro(StepSize, unsigned int);
+ itkSetMacro(Tolerance, double);
+
+ itkGetMacro(NumberOfPointsAfterReduction, unsigned int);
+
+ //Resets the filter, i.e. removes all inputs and outputs
+ void Reset();
+
+ /**
+ \brief Set whether the mitkProgressBar should be used
+
+ \a Parameter true for using the progress bar, false otherwise
+ */
+ void SetUseProgressBar(bool);
+
+ /**
+ \brief Set the stepsize which the progress bar should proceed
+
+ \a Parameter The stepsize for progressing
+ */
+ void SetProgressStepSize(unsigned int stepSize);
+
+ protected:
+ ReduceContourSetFilter();
+ virtual ~ReduceContourSetFilter();
+ virtual void GenerateData();
+ virtual void GenerateOutputInformation();
+
+ private:
+ void ReduceNumberOfPointsByNthPoint (vtkIdType cellSize, vtkIdType* cell, vtkPoints* points, vtkPolygon* reducedPolygon, vtkPoints* reducedPoints);
+
+ void ReduceNumberOfPointsByDouglasPeucker (vtkIdType cellSize, vtkIdType* cell, vtkPoints* points, vtkPolygon* reducedPolygon, vtkPoints* reducedPoints);
+
+ bool CheckForIntersection (vtkIdType* currentCell, vtkIdType currentCellSize, vtkPoints* currentPoints, /*vtkIdType numberOfIntersections, vtkIdType* intersectionPoints,*/ unsigned int currentInputIndex);
+
+ double m_MinSpacing;
+ double m_MaxSpacing;
+
+ Reduction_Type m_ReductionType;
+ unsigned int m_StepSize;
+ double m_Tolerance;
+ unsigned int m_MaxSegmentLenght;
+
+ bool m_UseProgressBar;
+ unsigned int m_ProgressStepSize;
+
+ unsigned int m_NumberOfPointsAfterReduction;
+
+ };//class
+
+}//namespace
+#endif
diff --git a/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.cpp b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.cpp
new file mode 100644
index 0000000..770f81f
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.cpp
@@ -0,0 +1,266 @@
+/*===================================================================
+
+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.
+
+===================================================================*/
+
+#include "mitkSurfaceInterpolationController.h"
+#include "mitkMemoryUtilities.h"
+
+mitk::SurfaceInterpolationController::SurfaceInterpolationController()
+ :m_SelectedSegmentation(0)
+{
+ m_ReduceFilter = ReduceContourSetFilter::New();
+ m_NormalsFilter = ComputeContourSetNormalsFilter::New();
+ m_InterpolateSurfaceFilter = CreateDistanceImageFromSurfaceFilter::New();
+
+ m_ReduceFilter->SetUseProgressBar(false);
+ m_NormalsFilter->SetUseProgressBar(false);
+ m_InterpolateSurfaceFilter->SetUseProgressBar(false);
+
+ m_Contours = Surface::New();
+
+ m_PolyData = vtkSmartPointer<vtkPolyData>::New();
+ m_PolyData->SetPoints(vtkPoints::New());
+
+ m_InterpolationResult = 0;
+ m_CurrentNumberOfReducedContours = 0;
+}
+
+mitk::SurfaceInterpolationController::~SurfaceInterpolationController()
+{
+
+ for (ContourListMap::iterator it = m_MapOfContourLists.begin(); it != m_MapOfContourLists.end(); it++)
+ {
+ for (unsigned int j = 0; j < m_MapOfContourLists[(*it).first].size(); ++j)
+ {
+ delete(m_MapOfContourLists[(*it).first].at(j).position);
+ }
+ m_MapOfContourLists.erase(it);
+ }
+}
+
+mitk::SurfaceInterpolationController* mitk::SurfaceInterpolationController::GetInstance()
+{
+ static mitk::SurfaceInterpolationController* m_Instance;
+
+ if ( m_Instance == 0)
+ {
+ m_Instance = new SurfaceInterpolationController();
+ }
+
+ return m_Instance;
+}
+
+void mitk::SurfaceInterpolationController::AddNewContour (mitk::Surface::Pointer newContour ,RestorePlanePositionOperation* op)
+{
+ AffineTransform3D::Pointer transform = AffineTransform3D::New();
+ transform = op->GetTransform();
+
+ mitk::Vector3D direction = op->GetDirectionVector();
+ int pos (-1);
+
+ for (unsigned int i = 0; i < m_MapOfContourLists[m_SelectedSegmentation].size(); i++)
+ {
+ itk::Matrix<float> diffM = transform->GetMatrix()-m_MapOfContourLists[m_SelectedSegmentation].at(i).position->GetTransform()->GetMatrix();
+ bool isSameMatrix(true);
+ for (unsigned int j = 0; j < 3; j++)
+ {
+ if (fabs(diffM[j][0]) > 0.0001 && fabs(diffM[j][1]) > 0.0001 && fabs(diffM[j][2]) > 0.0001)
+ {
+ isSameMatrix = false;
+ break;
+ }
+ }
+ itk::Vector<float> diffV = m_MapOfContourLists[m_SelectedSegmentation].at(i).position->GetTransform()->GetOffset()-transform->GetOffset();
+ if ( isSameMatrix && m_MapOfContourLists[m_SelectedSegmentation].at(i).position->GetPos() == op->GetPos() && (fabs(diffV[0]) < 0.0001 && fabs(diffV[1]) < 0.0001 && fabs(diffV[2]) < 0.0001) )
+ {
+ pos = i;
+ break;
+ }
+
+ }
+
+ if (pos == -1)
+ {
+ //MITK_INFO<<"New Contour";
+ mitk::RestorePlanePositionOperation* newOp = new mitk::RestorePlanePositionOperation (OpRESTOREPLANEPOSITION, op->GetWidth(),
+ op->GetHeight(), op->GetSpacing(), op->GetPos(), direction, transform);
+ ContourPositionPair newData;
+ newData.contour = newContour;
+ newData.position = newOp;
+
+ m_ReduceFilter->SetInput(m_MapOfContourLists[m_SelectedSegmentation].size(), newContour);
+ m_MapOfContourLists[m_SelectedSegmentation].push_back(newData);
+ }
+ else
+ {
+ //MITK_INFO<<"Modified Contour";
+ m_MapOfContourLists[m_SelectedSegmentation].at(pos).contour = newContour;
+ m_ReduceFilter->SetInput(pos, newContour);
+ }
+
+ m_ReduceFilter->Update();
+ m_CurrentNumberOfReducedContours = m_ReduceFilter->GetNumberOfOutputs();
+
+ for (unsigned int i = 0; i < m_CurrentNumberOfReducedContours; i++)
+ {
+ m_NormalsFilter->SetInput(i, m_ReduceFilter->GetOutput(i));
+ m_InterpolateSurfaceFilter->SetInput(i, m_NormalsFilter->GetOutput(i));
+ }
+
+ this->Modified();
+}
+
+void mitk::SurfaceInterpolationController::Interpolate()
+{
+ if (m_CurrentNumberOfReducedContours< 2)
+ return;
+
+ //Setting up progress bar
+ /*
+ * Removed due to bug 12441. ProgressBar messes around with Qt event queue which is fatal for segmentation
+ */
+ //mitk::ProgressBar::GetInstance()->AddStepsToDo(8);
+
+ m_InterpolateSurfaceFilter->Update();
+
+ Image::Pointer distanceImage = m_InterpolateSurfaceFilter->GetOutput();
+
+ vtkSmartPointer<vtkMarchingCubes> mcFilter = vtkSmartPointer<vtkMarchingCubes>::New();
+ mcFilter->SetInput(distanceImage->GetVtkImageData());
+ mcFilter->SetValue(0,0);
+ mcFilter->Update();
+
+ m_InterpolationResult = 0;
+ m_InterpolationResult = mitk::Surface::New();
+ m_InterpolationResult->SetVtkPolyData(mcFilter->GetOutput());
+ m_InterpolationResult->GetGeometry()->SetOrigin(distanceImage->GetGeometry()->GetOrigin());
+
+ vtkSmartPointer<vtkAppendPolyData> polyDataAppender = vtkSmartPointer<vtkAppendPolyData>::New();
+ for (unsigned int i = 0; i < m_ReduceFilter->GetNumberOfOutputs(); i++)
+ {
+ polyDataAppender->AddInput(m_ReduceFilter->GetOutput(i)->GetVtkPolyData());
+ }
+ polyDataAppender->Update();
+ m_Contours->SetVtkPolyData(polyDataAppender->GetOutput());
+
+ //Last progress step
+ /*
+ * Removed due to bug 12441. ProgressBar messes around with Qt event queue which is fatal for segmentation
+ */
+ //mitk::ProgressBar::GetInstance()->Progress(8);
+
+ m_InterpolationResult->DisconnectPipeline();
+}
+
+mitk::Surface::Pointer mitk::SurfaceInterpolationController::GetInterpolationResult()
+{
+ return m_InterpolationResult;
+}
+
+mitk::Surface* mitk::SurfaceInterpolationController::GetContoursAsSurface()
+{
+ return m_Contours;
+}
+
+void mitk::SurfaceInterpolationController::SetDataStorage(DataStorage &ds)
+{
+ m_DataStorage = &ds;
+}
+
+void mitk::SurfaceInterpolationController::SetMinSpacing(double minSpacing)
+{
+ m_ReduceFilter->SetMinSpacing(minSpacing);
+}
+
+void mitk::SurfaceInterpolationController::SetMaxSpacing(double maxSpacing)
+{
+ m_ReduceFilter->SetMaxSpacing(maxSpacing);
+ m_NormalsFilter->SetMaxSpacing(maxSpacing);
+}
+
+void mitk::SurfaceInterpolationController::SetDistanceImageVolume(unsigned int distImgVolume)
+{
+ m_InterpolateSurfaceFilter->SetDistanceImageVolume(distImgVolume);
+}
+
+void mitk::SurfaceInterpolationController::SetSegmentationImage(Image* workingImage)
+{
+ m_NormalsFilter->SetSegmentationBinaryImage(workingImage);
+}
+
+mitk::Image* mitk::SurfaceInterpolationController::GetImage()
+{
+ return m_InterpolateSurfaceFilter->GetOutput();
+}
+
+double mitk::SurfaceInterpolationController::EstimatePortionOfNeededMemory()
+{
+ double numberOfPointsAfterReduction = m_ReduceFilter->GetNumberOfPointsAfterReduction()*3;
+ double sizeOfPoints = pow(numberOfPointsAfterReduction,2)*sizeof(double);
+ double totalMem = mitk::MemoryUtilities::GetTotalSizeOfPhysicalRam();
+ double percentage = sizeOfPoints/totalMem;
+ return percentage;
+}
+
+void mitk::SurfaceInterpolationController::SetCurrentSegmentationInterpolationList(mitk::Image* segmentation)
+{
+ if (segmentation == m_SelectedSegmentation)
+ return;
+
+ if (segmentation == 0)
+ return;
+
+ ContourListMap::iterator it = m_MapOfContourLists.find(segmentation);
+
+ m_SelectedSegmentation = segmentation;
+
+ m_ReduceFilter->Reset();
+ m_NormalsFilter->Reset();
+ m_InterpolateSurfaceFilter->Reset();
+
+ if (it == m_MapOfContourLists.end())
+ {
+ ContourPositionPairList newList;
+ m_MapOfContourLists.insert(std::pair<mitk::Image*, ContourPositionPairList>(segmentation, newList));
+ m_InterpolationResult = 0;
+ m_CurrentNumberOfReducedContours = 0;
+ }
+ else
+ {
+ for (unsigned int i = 0; i < m_MapOfContourLists[m_SelectedSegmentation].size(); i++)
+ {
+ m_ReduceFilter->SetInput(i, m_MapOfContourLists[m_SelectedSegmentation].at(i).contour);
+ }
+
+ m_ReduceFilter->Update();
+
+ m_CurrentNumberOfReducedContours = m_ReduceFilter->GetNumberOfOutputs();
+
+ for (unsigned int i = 0; i < m_CurrentNumberOfReducedContours; i++)
+ {
+ m_NormalsFilter->SetInput(i, m_ReduceFilter->GetOutput(i));
+ m_InterpolateSurfaceFilter->SetInput(i, m_NormalsFilter->GetOutput(i));
+ }
+ }
+ Modified();
+}
+
+void mitk::SurfaceInterpolationController::RemoveSegmentationFromContourList(mitk::Image *segmentation)
+{
+ if (segmentation != 0)
+ {
+ m_MapOfContourLists.erase(segmentation);
+ }
+}
diff --git a/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.h b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.h
new file mode 100644
index 0000000..4207f56
--- /dev/null
+++ b/Modules/SurfaceInterpolation/mitkSurfaceInterpolationController.h
@@ -0,0 +1,167 @@
+/*===================================================================
+
+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 mitkSurfaceInterpolationController_h_Included
+#define mitkSurfaceInterpolationController_h_Included
+
+#include "mitkCommon.h"
+#include "SurfaceInterpolationExports.h"
+#include "mitkRestorePlanePositionOperation.h"
+#include "mitkSurface.h"
+#include "mitkInteractionConst.h"
+#include "mitkColorProperty.h"
+#include "mitkProperties.h"
+
+#include "mitkCreateDistanceImageFromSurfaceFilter.h"
+#include "mitkReduceContourSetFilter.h"
+#include "mitkComputeContourSetNormalsFilter.h"
+
+#include "mitkDataNode.h"
+#include "mitkDataStorage.h"
+#include "mitkWeakPointer.h"
+
+#include "vtkPolygon.h"
+#include "vtkPoints.h"
+#include "vtkCellArray.h"
+#include "vtkPolyData.h"
+#include "vtkSmartPointer.h"
+#include "vtkAppendPolyData.h"
+
+#include "vtkMarchingCubes.h"
+#include "vtkImageData.h"
+#include "mitkVtkRepresentationProperty.h"
+#include "vtkProperty.h"
+
+#include "mitkProgressBar.h"
+
+namespace mitk
+{
+
+ class SurfaceInterpolation_EXPORT SurfaceInterpolationController : public itk::Object
+ {
+
+ public:
+
+ mitkClassMacro(SurfaceInterpolationController, itk::Object)
+ itkNewMacro(Self)
+
+ static SurfaceInterpolationController* GetInstance();
+
+ /**
+ * Adds a new extracted contour to the list
+ */
+ void AddNewContour(Surface::Pointer newContour, RestorePlanePositionOperation *op);
+
+ /**
+ * Interpolates the 3D surface from the given extracted contours
+ */
+ void Interpolate ();
+
+ mitk::Surface::Pointer GetInterpolationResult();
+
+ /**
+ * Sets the minimum spacing of the current selected segmentation
+ * This is needed since the contour points we reduced before they are used to interpolate the surface
+ */
+ void SetMinSpacing(double minSpacing);
+
+ /**
+ * Sets the minimum spacing of the current selected segmentation
+ * This is needed since the contour points we reduced before they are used to interpolate the surface
+ */
+ void SetMaxSpacing(double maxSpacing);
+
+ /**
+ * Sets the volume i.e. the number of pixels that the distance image should have
+ * By evaluation we found out that 50.000 pixel delivers a good result
+ */
+ void SetDistanceImageVolume(unsigned int distImageVolume);
+
+ /**
+ * Sets the current segmentation which is used by the interpolation
+ * This is needed because the calculation of the normals needs to now wheather a normal points inside a segmentation or not
+ */
+ void SetSegmentationImage(Image* workingImage);
+
+ Surface* GetContoursAsSurface();
+
+ void SetDataStorage(DataStorage &ds);
+
+ /**
+ * Sets the current list of contourpoints which is used for the surface interpolation
+ * @param segmentation The current selected segmentation
+ */
+ void SetCurrentSegmentationInterpolationList(mitk::Image* segmentation);
+
+ /**
+ * Removes the segmentation and all its contours from the list
+ * @param segmentation The segmentation to be removed
+ */
+ void RemoveSegmentationFromContourList(mitk::Image* segmentation);
+
+ mitk::Image* GetImage();
+
+ /**
+ * Estimates the memory which is needed to build up the equationsystem for the interpolation.
+ * \returns The percentage of the real memory which will be used by the interpolation
+ */
+ double EstimatePortionOfNeededMemory();
+
+ protected:
+
+ SurfaceInterpolationController();
+
+ ~SurfaceInterpolationController();
+
+ private:
+
+ struct ContourPositionPair {
+ Surface::Pointer contour;
+ RestorePlanePositionOperation* position;
+ };
+
+ typedef std::vector<ContourPositionPair> ContourPositionPairList;
+ typedef std::map<mitk::Image* , ContourPositionPairList> ContourListMap;
+
+ ContourPositionPairList::iterator m_Iterator;
+
+ ReduceContourSetFilter::Pointer m_ReduceFilter;
+ ComputeContourSetNormalsFilter::Pointer m_NormalsFilter;
+ CreateDistanceImageFromSurfaceFilter::Pointer m_InterpolateSurfaceFilter;
+
+ double m_MinSpacing;
+ double m_MaxSpacing;
+
+ const Image* m_WorkingImage;
+
+ Surface::Pointer m_Contours;
+
+ vtkSmartPointer<vtkPolyData> m_PolyData;
+
+ unsigned int m_DistImageVolume;
+
+ mitk::WeakPointer<mitk::DataStorage> m_DataStorage;
+
+ ContourListMap m_MapOfContourLists;
+
+ mitk::Surface::Pointer m_InterpolationResult;
+
+ unsigned int m_CurrentNumberOfReducedContours;
+
+ mitk::Image* m_SelectedSegmentation;
+ };
+}
+#endif

File Metadata

Mime Type
text/plain
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
940
Default Alt Text
interpolationsegmentation.diff (233 KB)

Event Timeline