10 #include <itkCastImageFilter.h>
11 #include <itkVTKImageExport.h>
12 #include <itkVTKImageImport.h>
13 #include <itkImageDuplicator.h>
14 #include <itkExceptionObject.h>
16 #include <vtkPoints.h>
17 #include <vtkImageStencil.h>
18 #include <vtkImageImport.h>
19 #include <vtkImageExport.h>
20 #include <vtkLassoStencilSource.h>
21 #include <vtkSmartPointer.h>
31 if ( planarFigure.IsNull() )
33 throw std::runtime_error(
"Error: planar figure empty!" );
35 if ( !planarFigure->IsClosed() )
37 throw std::runtime_error(
"Masking not possible for non-closed figures" );
40 const PlaneGeometry *planarFigurePlaneGeometry = planarFigure->GetPlaneGeometry();
41 if ( planarFigurePlaneGeometry ==
nullptr )
43 throw std::runtime_error(
"Planar-Figure not yet initialized!" );
47 dynamic_cast< const PlaneGeometry *
>( planarFigurePlaneGeometry );
48 if ( planarFigureGeometry ==
nullptr )
50 throw std::runtime_error(
"Non-planar planar figures not supported!" );
53 if (planarFigure != m_PlanarFigure)
56 m_PlanarFigure = planarFigure;
63 if (IsUpdateRequired())
65 this->CalculateMask();
67 return m_ReferenceImage;
70 template <
typename TPixel,
unsigned int VImageDimension >
71 void PlanarFigureMaskGenerator::InternalCalculateMaskFromPlanarFigure(
72 const itk::Image< TPixel, VImageDimension > *image,
unsigned int axis )
74 typedef itk::Image< TPixel, VImageDimension >
ImageType;
75 typedef itk::Image< unsigned short, 2 > MaskImage2DType;
78 maskImage->SetOrigin(image->GetOrigin());
79 maskImage->SetSpacing(image->GetSpacing());
80 maskImage->SetLargestPossibleRegion(image->GetLargestPossibleRegion());
81 maskImage->SetBufferedRegion(image->GetBufferedRegion());
82 maskImage->SetDirection(image->GetDirection());
83 maskImage->SetNumberOfComponentsPerPixel(image->GetNumberOfComponentsPerPixel());
84 maskImage->Allocate();
85 maskImage->FillBuffer(1);
96 if (m_PlanarFigure->GetPolyLinesSize() == 2)
97 planarFigureHolePolyline = m_PlanarFigure->GetPolyLine(1);
123 bool outOfBounds =
false;
125 typename PlanarFigure::PolyLineType::const_iterator it;
126 for ( it = planarFigurePolyline.begin();
127 it != planarFigurePolyline.end();
137 planarFigurePlaneGeometry->
Map( *it, point3D );
141 if ( !imageGeometry3D->
IsInside( point3D ) )
148 points->InsertNextPoint( point3D[i0], point3D[i1], 0 );
151 vtkSmartPointer<vtkPoints> holePoints =
nullptr;
153 if (!planarFigureHolePolyline.empty())
158 PlanarFigure::PolyLineType::const_iterator end = planarFigureHolePolyline.end();
160 for (it = planarFigureHolePolyline.begin(); it != end; ++it)
163 planarFigurePlaneGeometry->
Map(*it, point3D);
165 holePoints->InsertNextPoint(point3D[i0], point3D[i1], 0);
171 double bounds[6] = {0, 0, 0, 0, 0, 0};
172 points->GetBounds( bounds );
173 bool extent_x = (fabs(bounds[0] - bounds[1])) <
mitk::eps;
174 bool extent_y = (fabs(bounds[2] - bounds[3])) <
mitk::eps;
175 bool extent_z = (fabs(bounds[4] - bounds[5])) <
mitk::eps;
178 if ( m_PlanarFigure->IsClosed() &&
179 ((extent_x && extent_y) || (extent_x && extent_z) || (extent_y && extent_z)))
181 mitkThrow() <<
"Figure has a zero area and cannot be used for masking.";
186 throw std::runtime_error(
"Figure at least partially outside of image bounds!" );
191 lassoStencil->SetShapeToPolygon();
192 lassoStencil->SetPoints( points );
194 vtkSmartPointer<vtkLassoStencilSource> holeLassoStencil =
nullptr;
196 if (holePoints.GetPointer() !=
nullptr)
199 holeLassoStencil->SetShapeToPolygon();
200 holeLassoStencil->SetPoints(holePoints);
204 typedef itk::VTKImageImport< MaskImage2DType > ImageImportType;
205 typedef itk::VTKImageExport< MaskImage2DType > ImageExportType;
208 itkExporter->SetInput( maskImage );
212 this->ConnectPipelines( itkExporter, vtkImporter );
216 imageStencilFilter->SetInputConnection( vtkImporter->GetOutputPort() );
217 imageStencilFilter->SetStencilConnection(lassoStencil->GetOutputPort());
218 imageStencilFilter->ReverseStencilOff();
219 imageStencilFilter->SetBackgroundValue( 0 );
220 imageStencilFilter->Update();
222 vtkSmartPointer<vtkImageStencil> holeStencilFilter =
nullptr;
224 if (holeLassoStencil.GetPointer() !=
nullptr)
227 holeStencilFilter->SetInputConnection(imageStencilFilter->GetOutputPort());
228 holeStencilFilter->SetStencilConnection(holeLassoStencil->GetOutputPort());
229 holeStencilFilter->ReverseStencilOn();
230 holeStencilFilter->SetBackgroundValue(0);
231 holeStencilFilter->Update();
236 vtkExporter->SetInputConnection( holeStencilFilter.GetPointer() ==
nullptr
237 ? imageStencilFilter->GetOutputPort()
238 : holeStencilFilter->GetOutputPort());
239 vtkExporter->Update();
242 this->ConnectPipelines( vtkExporter, itkImporter );
243 itkImporter->Update();
245 typedef itk::ImageDuplicator< ImageImportType::OutputImageType > DuplicatorType;
247 duplicator->SetInputImage( itkImporter->GetOutput() );
248 duplicator->Update();
251 m_InternalITKImageMask2D = duplicator->GetOutput();
254 bool PlanarFigureMaskGenerator::GetPrincipalAxis(
255 const BaseGeometry *geometry,
Vector3D vector,
259 for (
unsigned int i = 0; i < 3; ++i )
261 Vector3D axisVector = geometry->GetAxisVector( i );
262 axisVector.Normalize();
264 if ( fabs( fabs( axisVector * vector ) - 1.0) <
mitk::eps )
274 void PlanarFigureMaskGenerator::CalculateMask()
281 if (m_PlanarFigure.IsNull())
288 MITK_WARN <<
"Multiple TimeSteps are not supported in PlanarFigureMaskGenerator (yet).";
291 const BaseGeometry *imageGeometry =
m_inputImage->GetGeometry();
292 if ( imageGeometry ==
nullptr )
294 throw std::runtime_error(
"Image geometry invalid!" );
302 imgTimeSel->UpdateLargestPossibleRegion();
303 m_InternalTimeSliceImage = imgTimeSel->GetOutput();
310 m_InternalITKImageMask2D =
nullptr;
311 const PlaneGeometry *planarFigurePlaneGeometry = m_PlanarFigure->GetPlaneGeometry();
312 const PlaneGeometry *planarFigureGeometry =
dynamic_cast< const PlaneGeometry *
>( planarFigurePlaneGeometry );
317 if ( !this->GetPrincipalAxis( imageGeometry,
318 planarFigureGeometry->GetNormal(), axis ) )
320 throw std::runtime_error(
"Non-aligned planar figures not supported!" );
322 m_PlanarFigureAxis = axis;
325 itk::Image< unsigned short, 3 >::IndexType index;
326 imageGeometry->WorldToIndex( planarFigureGeometry->GetOrigin(), index );
328 unsigned int slice = index[axis];
335 InternalCalculateMaskFromPlanarFigure,
348 m_ReferenceImage = inputImageSlice;
363 if (IsUpdateRequired())
365 this->CalculateMask();
369 m_InternalMaskUpdateTime = this->GetMTime();
373 mitk::Image::Pointer PlanarFigureMaskGenerator::extract2DImageSlice(
unsigned int axis,
unsigned int slice)
376 unsigned int dimension = m_InternalTimeSliceImage->GetDimension();
382 imageExtractor->SetInput( m_InternalTimeSliceImage );
383 imageExtractor->SetSliceDimension( axis );
384 imageExtractor->SetSliceIndex( slice );
385 imageExtractor->Update();
386 imageSlice = imageExtractor->GetOutput();
388 else if(dimension == 2)
390 imageSlice = m_InternalTimeSliceImage;
394 MITK_ERROR <<
"Unsupported image dimension. Dimension is: " << dimension <<
". Only 2D and 3D images are supported.";
400 bool PlanarFigureMaskGenerator::IsUpdateRequired()
const
402 unsigned long thisClassTimeStamp = this->GetMTime();
404 unsigned long planarFigureTimeStamp = m_PlanarFigure->GetMTime();
405 unsigned long inputImageTimeStamp =
m_inputImage->GetMTime();
407 if (thisClassTimeStamp > m_InternalMaskUpdateTime)
412 if (m_InternalMaskUpdateTime < planarFigureTimeStamp || m_InternalMaskUpdateTime < inputImageTimeStamp)
417 if (internalMaskTimeStamp > m_InternalMaskUpdateTime)
mitk::Image::Pointer m_inputImage
itk::SmartPointer< Self > Pointer
virtual bool Map(const mitk::Point3D &pt3d_mm, mitk::Point2D &pt2d_mm) const
Project a 3D point given in mm (pt3d_mm) onto the 2D geometry. The result is a 2D point in mm (pt2d_m...
DataCollection - Class to facilitate loading/accessing structured data.
Image::Pointer GrabItkImageMemory(itk::SmartPointer< ItkOutputImageType > &itkimage, mitk::Image *mitkImage=nullptr, const BaseGeometry *geometry=nullptr, bool update=true)
Grabs the memory of an itk::Image (with a specific type) and puts it into an mitk::Image.The memory is managed by the mitk::Image after calling this function. The itk::Image remains valid until the mitk::Image decides to free the memory.
map::core::discrete::Elements< 3 >::InternalImageType ImageType
Vector< ScalarType, 3 > Vector3D
mitk::Image::Pointer m_InternalMask
Point< ScalarType, 3 > Point3D
#define AccessFixedDimensionByItk_1(mitkImage, itkImageTypeFunction, dimension, arg1)
MITKCORE_EXPORT const ScalarType eps
Describes a two-dimensional, rectangular plane.
bool IsInside(const mitk::Point3D &p) const
Test whether the point p (world coordinates in mm) is inside the bounding box.
BaseGeometry Describes the geometry of a data object.
void WorldToIndex(const mitk::Point3D &pt_mm, mitk::Point3D &pt_units) const
Convert world coordinates (in mm) of a point to (continuous!) index coordinates.
static itkEventMacro(BoundingShapeInteractionEvent, itk::AnyEvent) class MITKBOUNDINGSHAPE_EXPORT BoundingShapeInteractor Pointer New()
Basic interaction methods for mitk::GeometryData.