diff --git a/Core/Code/Algorithms/mitkExtractSliceFilter.cpp b/Core/Code/Algorithms/mitkExtractSliceFilter.cpp index c48e194ad6..1eefe28de0 100644 --- a/Core/Code/Algorithms/mitkExtractSliceFilter.cpp +++ b/Core/Code/Algorithms/mitkExtractSliceFilter.cpp @@ -1,552 +1,563 @@ /*========================================================================= Program: Medical Imaging & Interaction Toolkit Language: C++ Date: $Date$ Version: $Revision: $ Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. See MITKCopyright.txt or http://www.mitk.org/copyright.html for details. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notices for more information. =========================================================================*/ #include "mitkExtractSliceFilter.h" #include #include #include #include #include #include mitk::ExtractSliceFilter::ExtractSliceFilter(vtkImageReslice* reslicer ){ if(reslicer == NULL){ m_Reslicer = vtkSmartPointer::New(); } else { m_Reslicer = reslicer; } m_TimeStep = 0; m_Reslicer->ReleaseDataFlagOn(); m_InterpolationMode = ExtractSliceFilter::RESLICE_NEAREST; m_ResliceTransform = NULL; m_InPlaneResampleExtentByGeometry = false; m_OutPutSpacing = new mitk::ScalarType[2]; m_OutputDimension = 2; m_ZSpacing = 1.0; m_ZMin = 0; m_ZMax = 0; m_VtkOutputRequested = false; - + } mitk::ExtractSliceFilter::~ExtractSliceFilter(){ m_ResliceTransform = NULL; m_WorldGeometry = NULL; delete [] m_OutPutSpacing; } void mitk::ExtractSliceFilter::GenerateOutputInformation(){ //TODO try figure out how to set the specs of the slice before it is actually extracted /*Image::Pointer output = this->GetOutput(); Image::ConstPointer input = this->GetInput(); if (input.IsNull()) return; unsigned int dimensions[2]; dimensions[0] = m_WorldGeometry->GetExtent(0); dimensions[1] = m_WorldGeometry->GetExtent(1); output->Initialize(input->GetPixelType(), 2, dimensions, 1);*/ } void mitk::ExtractSliceFilter::GenerateInputRequestedRegion(){ //As we want all pixel information fo the image in our plane, the requested region //is set to the largest possible region in the image. //This is needed because an oblique plane has a larger extent then the image //and the in pipeline it is checked via PropagateResquestedRegion(). But the //extent of the slice is actually fitting because it is oblique within the image. ImageToImageFilter::InputImagePointer input = const_cast< ImageToImageFilter::InputImageType* > ( this->GetInput() ); input->SetRequestedRegionToLargestPossibleRegion(); } mitk::ScalarType* mitk::ExtractSliceFilter::GetOutputSpacing(){ return m_OutPutSpacing; } void mitk::ExtractSliceFilter::GenerateData(){ - + mitk::Image *input = const_cast< mitk::Image * >( this->GetInput() ); - + if (!input) { MITK_ERROR << "mitk::ExtractSliceFilter: No input image available. Please set the input!" << std::endl; - itkExceptionMacro("mitk::ExtractSliceFilter: No input image available. Please set the input!"); - return; - } + itkExceptionMacro("mitk::ExtractSliceFilter: No input image available. Please set the input!"); + return; + } if(!m_WorldGeometry) { MITK_ERROR << "mitk::ExtractSliceFilter: No Geometry for reslicing available." << std::endl; - itkExceptionMacro("mitk::ExtractSliceFilter: No Geometry for reslicing available."); - return; - } - + itkExceptionMacro("mitk::ExtractSliceFilter: No Geometry for reslicing available."); + return; + } + const TimeSlicedGeometry *inputTimeGeometry = this->GetInput()->GetTimeSlicedGeometry(); if ( ( inputTimeGeometry == NULL ) || ( inputTimeGeometry->GetTimeSteps() == 0 ) ) { itkWarningMacro(<<"Error reading input image TimeSlicedGeometry."); return; } // is it a valid timeStep? if ( inputTimeGeometry->IsValidTime( m_TimeStep ) == false ) { itkWarningMacro(<<"This is not a valid timestep: "<< m_TimeStep ); return; } - // check if there is something to display. - if ( ! input->IsVolumeSet( m_TimeStep ) ) - { + // check if there is something to display. + if ( ! input->IsVolumeSet( m_TimeStep ) ) + { itkWarningMacro(<<"No volume data existent at given timestep "<< m_TimeStep ); return; - } + } - -/*================#BEGIN setup vtkImageRslice properties================*/ - + /*================#BEGIN setup vtkImageRslice properties================*/ + + Point3D origin; Vector3D right, bottom, normal; double widthInMM, heightInMM; Vector2D extent; - - + + const PlaneGeometry* planeGeometry = dynamic_cast(m_WorldGeometry); - + if ( planeGeometry != NULL ) { //if the worldGeomatry is a PlaneGeometry everthing is straight forward - + origin = planeGeometry->GetOrigin(); right = planeGeometry->GetAxisVector( 0 ); bottom = planeGeometry->GetAxisVector( 1 ); normal = planeGeometry->GetNormal(); if ( m_InPlaneResampleExtentByGeometry ) { // Resampling grid corresponds to the current world geometry. This // means that the spacing of the output 2D image depends on the // currently selected world geometry, and *not* on the image itself. extent[0] = m_WorldGeometry->GetExtent( 0 ); extent[1] = m_WorldGeometry->GetExtent( 1 ); } else { // Resampling grid corresponds to the input geometry. This means that // the spacing of the output 2D image is directly derived from the // associated input image, regardless of the currently selected world // geometry. Vector3D rightInIndex, bottomInIndex; inputTimeGeometry->GetGeometry3D( m_TimeStep )->WorldToIndex( right, rightInIndex ); inputTimeGeometry->GetGeometry3D( m_TimeStep )->WorldToIndex( bottom, bottomInIndex ); extent[0] = rightInIndex.GetNorm(); extent[1] = bottomInIndex.GetNorm(); } // Get the extent of the current world geometry and calculate resampling // spacing therefrom. widthInMM = m_WorldGeometry->GetExtentInMM( 0 ); heightInMM = m_WorldGeometry->GetExtentInMM( 1 ); - + m_OutPutSpacing[0] = widthInMM / extent[0]; m_OutPutSpacing[1] = heightInMM / extent[1]; right.Normalize(); bottom.Normalize(); normal.Normalize(); - + /* - * Transform the origin to center based coordinates. - * Note: - * This is needed besause vtk's origin is center based too (!!!) ( see 'The VTK book' page 88 ) - * and the worldGeometry surrouding the image is no imageGeometry. So the worldGeometry - * has its origin at the corner of the voxel and needs to be transformed. - */ + * Transform the origin to center based coordinates. + * Note: + * This is needed besause vtk's origin is center based too (!!!) ( see 'The VTK book' page 88 ) + * and the worldGeometry surrouding the image is no imageGeometry. So the worldGeometry + * has its origin at the corner of the voxel and needs to be transformed. + */ origin += right * ( m_OutPutSpacing[0] * 0.5 ); origin += bottom * ( m_OutPutSpacing[1] * 0.5 ); - - + + //set the tranform for reslicing. // Use inverse transform of the input geometry for reslicing the 3D image. // This is needed if the image volume already transformed if(m_ResliceTransform != NULL) m_Reslicer->SetResliceTransform(m_ResliceTransform->GetVtkTransform()->GetLinearInverse()); // Set background level to TRANSLUCENT (see Geometry2DDataVtkMapper3D), // else the background of the image turns out gray m_Reslicer->SetBackgroundLevel( -32768 ); } else{ //Code for curved planes, mostly taken 1:1 from imageVtkMapper2D and not tested yet. // Do we have an AbstractTransformGeometry? // This is the case for AbstractTransformGeometry's (e.g. a ThinPlateSplineCurvedGeometry ) const mitk::AbstractTransformGeometry* abstractGeometry = dynamic_cast< const AbstractTransformGeometry * >(m_WorldGeometry); if(abstractGeometry != NULL) { m_ResliceTransform = abstractGeometry; extent[0] = abstractGeometry->GetParametricExtent(0); extent[1] = abstractGeometry->GetParametricExtent(1); widthInMM = abstractGeometry->GetParametricExtentInMM(0); heightInMM = abstractGeometry->GetParametricExtentInMM(1); m_OutPutSpacing[0] = widthInMM / extent[0]; m_OutPutSpacing[1] = heightInMM / extent[1]; origin = abstractGeometry->GetPlane()->GetOrigin(); right = abstractGeometry->GetPlane()->GetAxisVector(0); right.Normalize(); bottom = abstractGeometry->GetPlane()->GetAxisVector(1); bottom.Normalize(); normal = abstractGeometry->GetPlane()->GetNormal(); normal.Normalize(); // Use a combination of the InputGeometry *and* the possible non-rigid // AbstractTransformGeometry for reslicing the 3D Image vtkSmartPointer composedResliceTransform = vtkGeneralTransform::New(); composedResliceTransform->Identity(); composedResliceTransform->Concatenate( inputTimeGeometry->GetGeometry3D( m_TimeStep )->GetVtkTransform()->GetLinearInverse() ); composedResliceTransform->Concatenate( abstractGeometry->GetVtkAbstractTransform() ); m_Reslicer->SetResliceTransform( composedResliceTransform ); composedResliceTransform->UnRegister( NULL ); // decrease RC // Set background level to BLACK instead of translucent, to avoid // boundary artifacts (see Geometry2DDataVtkMapper3D) m_Reslicer->SetBackgroundLevel( -1023 ); } else { itkExceptionMacro("mitk::ExtractSliceFilter: No fitting geometry for reslice axis!"); return; } } if(m_ResliceTransform != NULL){ //if the resliceTransform is set the reslice axis are recalculated. //Thus the geometry information is not fitting. Therefor a unitSpacingFilter //is used to set up a global spacing of 1 and compensate the transform. vtkSmartPointer unitSpacingImageFilter = vtkSmartPointer::New() ; unitSpacingImageFilter->ReleaseDataFlagOn(); unitSpacingImageFilter->SetOutputSpacing( 1.0, 1.0, 1.0 ); unitSpacingImageFilter->SetInput( input->GetVtkImageData(m_TimeStep) ); m_Reslicer->SetInput(unitSpacingImageFilter->GetOutput() ); } else { //if no tranform is set the image can be used directly m_Reslicer->SetInput(input->GetVtkImageData(m_TimeStep)); } /*setup the plane where vktImageReslice extracts the slice*/ //ResliceAxesOrigin is the ancor point of the plane double originInVtk[3]; itk2vtk(origin, originInVtk); m_Reslicer->SetResliceAxesOrigin(originInVtk); //the cosines define the plane: x and y are the direction vectors, n is the planes normal //this specifies a matrix 3x3 // x1 y1 n1 // x2 y2 n2 // x3 y3 n3 double cosines[9]; vnl2vtk(right.GetVnlVector(), cosines);//x vnl2vtk(bottom.GetVnlVector(), cosines + 3);//y vnl2vtk(normal.GetVnlVector(), cosines + 6);//n m_Reslicer->SetResliceAxesDirectionCosines(cosines); //we only have one slice, not a volume m_Reslicer->SetOutputDimensionality(m_OutputDimension); //set the interpolation mode for slicing switch(this->m_InterpolationMode){ case RESLICE_NEAREST: m_Reslicer->SetInterpolationModeToNearestNeighbor(); break; case RESLICE_LINEAR: m_Reslicer->SetInterpolationModeToLinear(); break; case RESLICE_CUBIC: m_Reslicer->SetInterpolationModeToCubic(); default: //the default interpolation used by mitk m_Reslicer->SetInterpolationModeToNearestNeighbor(); } - + /*========== BEGIN setup extent of the slice ==========*/ int xMin, xMax, yMin, yMax; bool boundsInitialized = false; if(m_WorldGeometry->GetReferenceGeometry()){ vtkFloatingPointType sliceBounds[6]; - + for ( int i = 0; i < 6; ++i ) { sliceBounds[i] = 0.0; } this->CalculateClippedPlaneBounds( m_WorldGeometry->GetReferenceGeometry(), planeGeometry, sliceBounds ); - + // Calculate output extent (integer values) xMin = static_cast< int >( sliceBounds[0] / m_OutPutSpacing[0] + 0.5 ); xMax = static_cast< int >( sliceBounds[1] / m_OutPutSpacing[0] + 0.5 ); yMin = static_cast< int >( sliceBounds[2] / m_OutPutSpacing[1] + 0.5 ); yMax = static_cast< int >( sliceBounds[3] / m_OutPutSpacing[1] + 0.5 ); - + }else { // If no reference geometry is available, we also don't know about the // maximum plane size; xMin = yMin = 0; xMax = static_cast< int >( extent[0]); yMax = static_cast< int >( extent[1]); } - - + + m_Reslicer->SetOutputExtent(xMin, xMax-1, yMin, yMax-1, m_ZMin, m_ZMax ); /*========== END setup extent of the slice ==========*/ - + m_Reslicer->SetOutputOrigin( 0.0, 0.0, 0.0 ); m_Reslicer->SetOutputSpacing( m_OutPutSpacing[0], m_OutPutSpacing[1], m_ZSpacing ); - + //TODO check the following lines, they are responsible wether vtk error outputs appear or not m_Reslicer->UpdateWholeExtent(); //this produces a bad allocation error for 2D images //m_Reslicer->GetOutput()->UpdateInformation(); //m_Reslicer->GetOutput()->SetUpdateExtentToWholeExtent(); - + //start the pipeline m_Reslicer->Update(); -/*================ #END setup vtkImageRslice properties================*/ + /*================ #END setup vtkImageRslice properties================*/ if(m_VtkOutputRequested){ return; //no converting to mitk //no mitk geometry will be set, as the output is vtkImageData only!!! } else { /*================ #BEGIN Get the slice from vtkImageReslice and convert it to mit::Image================*/ vtkImageData* reslicedImage; reslicedImage = m_Reslicer->GetOutput(); if(!reslicedImage) { itkWarningMacro(<<"Reslicer returned empty image"); return; } mitk::Image::Pointer resultImage = this->GetOutput(); //initialize resultimage with the specs of the vtkImageData object returned from vtkImageReslice resultImage->Initialize(reslicedImage); //transfer the voxel data resultImage->SetVolume(reslicedImage->GetScalarPointer()); //set the geometry from current worldgeometry for the reusultimage //this is needed that the image has the correct mitk geometry //the originalGeometry is the Geometry of the result slice AffineGeometryFrame3D::Pointer originalGeometryAGF = m_WorldGeometry->Clone(); Geometry2D::Pointer originalGeometry = dynamic_cast( originalGeometryAGF.GetPointer() ); //the origin of the worldGeometry is transformed to center based coordinates to be an imageGeometry Point3D sliceOrigin = originalGeometry->GetOrigin(); sliceOrigin += right * ( m_OutPutSpacing[0] * 0.5 ); sliceOrigin += bottom * ( m_OutPutSpacing[1] * 0.5 ); //a worldGeometry is no imageGeometry, thus it is manually set to true originalGeometry->ImageGeometryOn(); originalGeometry->SetOrigin(sliceOrigin); resultImage->SetGeometry( originalGeometry ); //resultImage->GetGeometry()->TransferItkToVtkTransform(); resultImage->GetGeometry()->Modified(); /*================ #END Get the slice from vtkImageReslice and convert it to mitk Image================*/ } } -bool mitk::ExtractSliceFilter::GetBounds(vtkFloatingPointType bounds[6]){ - +bool mitk::ExtractSliceFilter::GetClippedPlaneBounds(vtkFloatingPointType bounds[6]){ + + if(!m_WorldGeometry || !this->GetInput()) + return false; + return this->CalculateClippedPlaneBounds(m_WorldGeometry->GetReferenceGeometry(), dynamic_cast< const PlaneGeometry * >( m_WorldGeometry ), bounds); } + +bool mitk::ExtractSliceFilter::GetClippedPlaneBounds( const Geometry3D *boundingGeometry, + const PlaneGeometry *planeGeometry, vtkFloatingPointType *bounds ) +{ + return this->CalculateClippedPlaneBounds(boundingGeometry, planeGeometry, bounds); +} + + bool mitk::ExtractSliceFilter ::CalculateClippedPlaneBounds( const Geometry3D *boundingGeometry, - const PlaneGeometry *planeGeometry, vtkFloatingPointType *bounds ) + const PlaneGeometry *planeGeometry, vtkFloatingPointType *bounds ) { // Clip the plane with the bounding geometry. To do so, the corner points // of the bounding box are transformed by the inverse transformation // matrix, and the transformed bounding box edges derived therefrom are // clipped with the plane z=0. The resulting min/max values are taken as // bounds for the image reslicer. const BoundingBox *boundingBox = boundingGeometry->GetBoundingBox(); BoundingBox::PointType bbMin = boundingBox->GetMinimum(); BoundingBox::PointType bbMax = boundingBox->GetMaximum(); BoundingBox::PointType bbCenter = boundingBox->GetCenter(); vtkPoints *points = vtkPoints::New(); if(boundingGeometry->GetImageGeometry()) { points->InsertPoint( 0, bbMin[0]-0.5, bbMin[1]-0.5, bbMin[2]-0.5 ); points->InsertPoint( 1, bbMin[0]-0.5, bbMin[1]-0.5, bbMax[2]-0.5 ); points->InsertPoint( 2, bbMin[0]-0.5, bbMax[1]-0.5, bbMax[2]-0.5 ); points->InsertPoint( 3, bbMin[0]-0.5, bbMax[1]-0.5, bbMin[2]-0.5 ); points->InsertPoint( 4, bbMax[0]-0.5, bbMin[1]-0.5, bbMin[2]-0.5 ); points->InsertPoint( 5, bbMax[0]-0.5, bbMin[1]-0.5, bbMax[2]-0.5 ); points->InsertPoint( 6, bbMax[0]-0.5, bbMax[1]-0.5, bbMax[2]-0.5 ); points->InsertPoint( 7, bbMax[0]-0.5, bbMax[1]-0.5, bbMin[2]-0.5 ); } else { points->InsertPoint( 0, bbMin[0], bbMin[1], bbMin[2] ); points->InsertPoint( 1, bbMin[0], bbMin[1], bbMax[2] ); points->InsertPoint( 2, bbMin[0], bbMax[1], bbMax[2] ); points->InsertPoint( 3, bbMin[0], bbMax[1], bbMin[2] ); points->InsertPoint( 4, bbMax[0], bbMin[1], bbMin[2] ); points->InsertPoint( 5, bbMax[0], bbMin[1], bbMax[2] ); points->InsertPoint( 6, bbMax[0], bbMax[1], bbMax[2] ); points->InsertPoint( 7, bbMax[0], bbMax[1], bbMin[2] ); } vtkPoints *newPoints = vtkPoints::New(); vtkTransform *transform = vtkTransform::New(); transform->Identity(); transform->Concatenate( planeGeometry->GetVtkTransform()->GetLinearInverse() ); transform->Concatenate( boundingGeometry->GetVtkTransform() ); transform->TransformPoints( points, newPoints ); transform->Delete(); bounds[0] = bounds[2] = 10000000.0; bounds[1] = bounds[3] = -10000000.0; bounds[4] = bounds[5] = 0.0; this->LineIntersectZero( newPoints, 0, 1, bounds ); this->LineIntersectZero( newPoints, 1, 2, bounds ); this->LineIntersectZero( newPoints, 2, 3, bounds ); this->LineIntersectZero( newPoints, 3, 0, bounds ); this->LineIntersectZero( newPoints, 0, 4, bounds ); this->LineIntersectZero( newPoints, 1, 5, bounds ); this->LineIntersectZero( newPoints, 2, 6, bounds ); this->LineIntersectZero( newPoints, 3, 7, bounds ); this->LineIntersectZero( newPoints, 4, 5, bounds ); this->LineIntersectZero( newPoints, 5, 6, bounds ); this->LineIntersectZero( newPoints, 6, 7, bounds ); this->LineIntersectZero( newPoints, 7, 4, bounds ); // clean up vtk data points->Delete(); newPoints->Delete(); if ( (bounds[0] > 9999999.0) || (bounds[2] > 9999999.0) || (bounds[1] < -9999999.0) || (bounds[3] < -9999999.0) ) { return false; } else { // The resulting bounds must be adjusted by the plane spacing, since we // we have so far dealt with index coordinates const float *planeSpacing = planeGeometry->GetFloatSpacing(); bounds[0] *= planeSpacing[0]; bounds[1] *= planeSpacing[0]; bounds[2] *= planeSpacing[1]; bounds[3] *= planeSpacing[1]; bounds[4] *= planeSpacing[2]; bounds[5] *= planeSpacing[2]; return true; } } bool mitk::ExtractSliceFilter ::LineIntersectZero( vtkPoints *points, int p1, int p2, - vtkFloatingPointType *bounds ) + vtkFloatingPointType *bounds ) { vtkFloatingPointType point1[3]; vtkFloatingPointType point2[3]; points->GetPoint( p1, point1 ); points->GetPoint( p2, point2 ); if ( (point1[2] * point2[2] <= 0.0) && (point1[2] != point2[2]) ) { double x, y; x = ( point1[0] * point2[2] - point1[2] * point2[0] ) / ( point2[2] - point1[2] ); y = ( point1[1] * point2[2] - point1[2] * point2[1] ) / ( point2[2] - point1[2] ); if ( x < bounds[0] ) { bounds[0] = x; } if ( x > bounds[1] ) { bounds[1] = x; } if ( y < bounds[2] ) { bounds[2] = y; } if ( y > bounds[3] ) { bounds[3] = y; } bounds[4] = bounds[5] = 0.0; return true; } return false; } \ No newline at end of file diff --git a/Core/Code/Algorithms/mitkExtractSliceFilter.h b/Core/Code/Algorithms/mitkExtractSliceFilter.h index 44d4bdd7a8..62aad42035 100644 --- a/Core/Code/Algorithms/mitkExtractSliceFilter.h +++ b/Core/Code/Algorithms/mitkExtractSliceFilter.h @@ -1,162 +1,185 @@ /*========================================================================= - + Program: Medical Imaging & Interaction Toolkit Language: C++ Date: $Date$ Version: $Revision: $ - + Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. See MITKCopyright.txt or http://www.mitk.org/copyright.html for details. - + This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notices for more information. - + =========================================================================*/ #ifndef mitkExtractSliceFilter_h_Included #define mitkExtractSliceFilter_h_Included #include "MitkExports.h" #include "mitkImageToImageFilter.h" #include #include #include #include #include #include #include namespace mitk { /** - \brief ExtractSliceFilter extracts a 2D abitrary oriented slice from a 3D volume. - - The filter can reslice in all orthogonal planes such as sagittal, coronal and transversal, - and is also able to reslice a abitrary oriented oblique plane. - Curved planes are specified via an AbstractTransformGeometry as the input worldgeometry. - - The convinient workflow is: - 1. Set an image as input. - 2. Set the worldGeometry2D. This defines a grid where the slice is being extracted - 3. And then start the pipeline. - - There are a few more properties that can be set to modify the behavior of the slicing. - The properties are: - - interpolation mode either Nearestneighbor, Linear or Cubic. - - a transform this is a convinient way to adapt the reslice axis for the case - that the image is transformed e.g. rotated. - - time step the time step in a timesliced volume. - - resample by geometry wether the resampling grid corresponds to the specs of the - worldgeometry or is directly derived from the input image - - By default the properties are set to: - - interpolation mode Nearestneighbor. - - a transform NULL (No transform is set). - - time step 0. - - resample by geometry false (Corresponds to input image). + \brief ExtractSliceFilter extracts a 2D abitrary oriented slice from a 3D volume. + + The filter can reslice in all orthogonal planes such as sagittal, coronal and transversal, + and is also able to reslice a abitrary oriented oblique plane. + Curved planes are specified via an AbstractTransformGeometry as the input worldgeometry. + + The convinient workflow is: + 1. Set an image as input. + 2. Set the worldGeometry2D. This defines a grid where the slice is being extracted + 3. And then start the pipeline. + + There are a few more properties that can be set to modify the behavior of the slicing. + The properties are: + - interpolation mode either Nearestneighbor, Linear or Cubic. + - a transform this is a convinient way to adapt the reslice axis for the case + that the image is transformed e.g. rotated. + - time step the time step in a timesliced volume. + - resample by geometry wether the resampling grid corresponds to the specs of the + worldgeometry or is directly derived from the input image + + By default the properties are set to: + - interpolation mode Nearestneighbor. + - a transform NULL (No transform is set). + - time step 0. + - resample by geometry false (Corresponds to input image). */ class MITK_CORE_EXPORT ExtractSliceFilter : public ImageToImageFilter { public: mitkClassMacro(ExtractSliceFilter, ImageToImageFilter); itkNewMacro(ExtractSliceFilter); mitkNewMacro1Param(Self, vtkImageReslice*); - + /** \brief Set the axis where to reslice at.*/ void SetWorldGeometry(const Geometry2D* geometry ){ this->m_WorldGeometry = geometry; } /** \brief Set the time step in the 4D volume */ void SetTimeStep( unsigned int timestep){ this->m_TimeStep = timestep; } unsigned int GetTimeStep(){ return this->m_TimeStep; } - /** - \brief Set a transform for the reslice axes. - - This transform is needed if the image volume itself is transformed. (Effects the reslice axis) + /** \brief Set a transform for the reslice axes. + * This transform is needed if the image volume itself is transformed. (Effects the reslice axis) */ void SetResliceTransformByGeometry(const Geometry3D* transform){ this->m_ResliceTransform = transform; } /** \brief Resampling grid corresponds to: false->image true->worldgeometry*/ void SetInPlaneResampleExtentByGeometry(bool inPlaneResampleExtentByGeometry){ this->m_InPlaneResampleExtentByGeometry = inPlaneResampleExtentByGeometry; } /** \brief Sets the output dimension of the slice*/ void SetOutputDimensionality(unsigned int dimension){ this->m_OutputDimension = dimension; } /** \brief Set the spacing in z direction manually. - Required if the outputDimension is > 2. + * Required if the outputDimension is > 2. */ void SetOutputSpacingZDirection(double zSpacing){ this->m_ZSpacing = zSpacing; } /** \brief Set the extent in pixel for direction z manualy. - Required if the output dimension is > 2. + Required if the output dimension is > 2. */ void SetOutputExtentZDirection(int zMin, int zMax) { this->m_ZMin = zMin; this->m_ZMax = zMax; } + /** \brief Get the bounding box of the slice [xMin, xMax, yMin, yMax, zMin, zMax] + * The method uses the input of the filter to calculate the bounds. + * It is recommended to use + * GetClippedPlaneBounds(const Geometry3D*, const PlaneGeometry*, vtkFloatingPointType*) + * if you are not sure about the input. + */ + bool GetClippedPlaneBounds(double bounds[6]); + /** \brief Get the bounding box of the slice [xMin, xMax, yMin, yMax, zMin, zMax]*/ - bool GetBounds(double bounds[6]); + bool GetClippedPlaneBounds( const Geometry3D *boundingGeometry, + const PlaneGeometry *planeGeometry, vtkFloatingPointType *bounds ); /** \brief Get the spacing of the slice. returns mitk::ScalarType[2] */ mitk::ScalarType* GetOutputSpacing(); - /** \brief Get Output as vtkImageData. */ + /** \brief Get Output as vtkImageData. + * Note: + * SetVtkOutputRequest(true) has to be called at least once before + * GetVtkOutput(). Otherwise the output is empty for the first update step. + */ vtkImageData* GetVtkOutput(){ m_VtkOutputRequested = true; return m_Reslicer->GetOutput(); } - /** */ + /** Set VtkOutPutRequest to suppress the convertion of the image. + * It is suggested to use this with GetVtkOutput(). + * Note: + * SetVtkOutputRequest(true) has to be called at least once before + * GetVtkOutput(). Otherwise the output is empty for the first update step. + */ void SetVtkOutputRequest(bool isRequested){ m_VtkOutputRequested = isRequested; } /** \brief Get the reslices axis matrix. - Note: the axis are recalculated when calling SetResliceTransformByGeometry. + * Note: the axis are recalculated when calling SetResliceTransformByGeometry. */ vtkMatrix4x4* GetResliceAxes(){ return this->m_Reslicer->GetResliceAxes(); } enum ResliceInterpolation { RESLICE_NEAREST, RESLICE_LINEAR, RESLICE_CUBIC }; - + void SetInterPolationMode( ExtractSliceFilter::ResliceInterpolation interpolation){ this->m_InterpolationMode = interpolation; } protected: ExtractSliceFilter(vtkImageReslice* reslicer = NULL); virtual ~ExtractSliceFilter(); virtual void GenerateData(); virtual void GenerateOutputInformation(); virtual void GenerateInputRequestedRegion(); const Geometry2D* m_WorldGeometry; vtkSmartPointer m_Reslicer; unsigned int m_TimeStep; unsigned int m_OutputDimension; double m_ZSpacing; int m_ZMin; int m_ZMax; ResliceInterpolation m_InterpolationMode; const Geometry3D* m_ResliceTransform; bool m_InPlaneResampleExtentByGeometry;//Resampling grid corresponds to: false->image true->worldgeometry mitk::ScalarType* m_OutPutSpacing; - + bool m_VtkOutputRequested; - bool CalculateClippedPlaneBounds( const Geometry3D *boundingGeometry, const PlaneGeometry *planeGeometry, vtkFloatingPointType *bounds ); - bool LineIntersectZero( vtkPoints *points, int p1, int p2, vtkFloatingPointType *bounds ); + /** \brief Internal helper method for intersection testing used only in CalculateClippedPlaneBounds() */ + bool LineIntersectZero( vtkPoints *points, int p1, int p2, + vtkFloatingPointType *bounds ); + + /** \brief Calculate the bounding box of the resliced image. This is necessary for + * arbitrarily rotated planes in an image volume. A rotated plane (e.g. in swivel mode) + * will have a new bounding box, which needs to be calculated. */ + bool CalculateClippedPlaneBounds( const Geometry3D *boundingGeometry, + const PlaneGeometry *planeGeometry, vtkFloatingPointType *bounds ); }; } #endif // mitkExtractSliceFilter_h_Included \ No newline at end of file diff --git a/Core/Code/Rendering/mitkImageVtkMapper2D.h b/Core/Code/Rendering/mitkImageVtkMapper2D.h index 16b7d7980e..96ae68a40d 100644 --- a/Core/Code/Rendering/mitkImageVtkMapper2D.h +++ b/Core/Code/Rendering/mitkImageVtkMapper2D.h @@ -1,269 +1,260 @@ /*========================================================================= Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. See MITKCopyright.txt or http://www.mitk.org/copyright.html for details. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notices for more information. =========================================================================*/ #ifndef MITKIMAGEVTKMAPPER2D_H_HEADER_INCLUDED_C10E906E #define MITKIMAGEVTKMAPPER2D_H_HEADER_INCLUDED_C10E906E //MITK #include //MITK Rendering #include "mitkBaseRenderer.h" #include "mitkVtkMapper2D.h" #include "mitkExtractSliceFilter.h" //VTK #include class vtkActor; class vtkPolyDataMapper; class vtkPlaneSource; class vtkImageData; class vtkLookupTable; class vtkImageReslice; class vtkImageChangeInformation; class vtkPoints; class vtkMitkThickSlicesFilter; class vtkPolyData; class vtkMitkApplyLevelWindowToRGBFilter; namespace mitk { /** \brief Mapper to resample and display 2D slices of a 3D image. * * The following image gives a brief overview of the mapping and the involved parts. * * \image html imageVtkMapper2Darchitecture.png * * First, the image is resliced by means of vtkImageReslice. The volume image * serves as input to the mapper in addition to spatial placement of the slice and a few other * properties such as thick slices. This code was already present in the old version * (mitkImageMapperGL2D). * * Next, the obtained slice (m_ReslicedImage) is used to create a texture * (m_Texture) and a plane onto which the texture is rendered (m_Plane). For * mapping purposes, a vtkPolyDataMapper (m_Mapper) is utilized. Orthographic * projection is applied to create the effect of a 2D image. The mapper and the * texture are assigned to the actor (m_Actor) which is passed to the VTK rendering * pipeline via the method GetVtkProp(). * * In order to transform the textured plane to the correct position in space, the * same transformation as used for reslicing is applied to both the camera and the * vtkActor. All important steps are explained in more detail below. The resulting * 2D image (by reslicing the underlying 3D input image appropriately) can either * be directly rendered in a 2D view or just be calculated to be used later by another * rendering entity, e.g. in texture mapping in a 3D view. * * Properties that can be set for images and influence the imageMapper2D are: * * - \b "opacity": (FloatProperty) Opacity of the image * - \b "color": (ColorProperty) Color of the image * - \b "use color": (BoolProperty) Use the color of the image or not * - \b "binary": (BoolProperty) is the image a binary image or not * - \b "outline binary": (BoolProperty) show outline of the image or not * - \b "texture interpolation": (BoolProperty) texture interpolation of the image * - \b "reslice interpolation": (VtkResliceInterpolationProperty) reslice interpolation of the image * - \b "in plane resample extent by geometry": (BoolProperty) Do it or not * - \b "bounding box": (BoolProperty) Is the Bounding Box of the image shown or not * - \b "layer": (IntProperty) Layer of the image * - \b "volume annotation color": (ColorProperty) color of the volume annotation, TODO has to be reimplemented * - \b "volume annotation unit": (StringProperty) annotation unit as string (does not implicit convert the unit!) unit is ml or cm3, TODO has to be reimplemented * The default properties are: * - \b "opacity", mitk::FloatProperty::New(0.3f), renderer, overwrite ) * - \b "color", ColorProperty::New(1.0,0.0,0.0), renderer, overwrite ) * - \b "use color", mitk::BoolProperty::New( true ), renderer, overwrite ) * - \b "binary", mitk::BoolProperty::New( true ), renderer, overwrite ) * - \b "outline binary", mitk::BoolProperty::New( false ), renderer, overwrite ) * - \b "texture interpolation", mitk::BoolProperty::New( mitk::DataNodeFactory::m_TextureInterpolationActive ) ) * - \b "reslice interpolation", mitk::VtkResliceInterpolationProperty::New() ) * - \b "in plane resample extent by geometry", mitk::BoolProperty::New( false ) ) * - \b "bounding box", mitk::BoolProperty::New( false ) ) * - \b "layer", mitk::IntProperty::New(10), renderer, overwrite) * If the modality-property is set for an image, the mapper uses modality-specific default properties, * e.g. color maps, if they are defined. * \ingroup Mapper */ class MITK_CORE_EXPORT ImageVtkMapper2D : public VtkMapper2D { public: /** Standard class typedefs. */ mitkClassMacro( ImageVtkMapper2D,VtkMapper2D ); /** Method for creation through the object factory. */ itkNewMacro(Self); /** \brief Get the Image to map */ const mitk::Image *GetInput(void); /** \brief Checks whether this mapper needs to update itself and generate * data. */ virtual void Update(mitk::BaseRenderer * renderer); //### methods of MITK-VTK rendering pipeline virtual vtkProp* GetVtkProp(mitk::BaseRenderer* renderer); virtual void MitkRenderOverlay(BaseRenderer* renderer); virtual void MitkRenderOpaqueGeometry(BaseRenderer* renderer); virtual void MitkRenderTranslucentGeometry(BaseRenderer* renderer); virtual void MitkRenderVolumetricGeometry(BaseRenderer* renderer); //### end of methods of MITK-VTK rendering pipeline /** \brief Internal class holding the mapper, actor, etc. for each of the 3 2D render windows */ /** * To render transveral, coronal, and sagittal, the mapper is called three times. * For performance reasons, the corresponding data for each view is saved in the * internal helper class LocalStorage. This allows rendering n views with just * 1 mitkMapper using n vtkMapper. * */ class MITK_CORE_EXPORT LocalStorage : public mitk::Mapper::BaseLocalStorage { public: /** \brief Actor of a 2D render window. */ vtkSmartPointer m_Actor; /** \brief Mapper of a 2D render window. */ vtkSmartPointer m_Mapper; /** \brief Current slice of a 2D render window. */ vtkSmartPointer m_ReslicedImage; /** \brief Plane on which the slice is rendered as texture. */ vtkSmartPointer m_Plane; /** \brief The texture which is used to render the current slice. */ vtkSmartPointer m_Texture; /** \brief The lookuptable for colors and level window */ vtkSmartPointer m_LookupTable; /** \brief The actual reslicer (one per renderer) */ mitk::ExtractSliceFilter::Pointer m_Reslicer; /** \brief Filter for thick slices */ vtkSmartPointer m_TSFilter; /** \brief PolyData object containg all lines/points needed for outlining the contour. This container is used to save a computed contour for the next rendering execution. For instance, if you zoom or pann, there is no need to recompute the contour. */ vtkSmartPointer m_OutlinePolyData; /** \brief Timestamp of last update of stored data. */ itk::TimeStamp m_LastUpdateTime; /** \brief mmPerPixel relation between pixel and mm. (World spacing).*/ mitk::ScalarType* m_mmPerPixel; /** \brief This filter is used to apply the level window to RBG(A) images. */ vtkMitkApplyLevelWindowToRGBFilter* m_LevelWindowToRGBFilterObject; /** \brief Default constructor of the local storage. */ LocalStorage(); /** \brief Default deconstructor of the local storage. */ ~LocalStorage() { } }; /** \brief The LocalStorageHandler holds all (three) LocalStorages for the three 2D render windows. */ mitk::Mapper::LocalStorageHandler m_LSH; /** \brief Get the LocalStorage corresponding to the current renderer. */ LocalStorage* GetLocalStorage(mitk::BaseRenderer* renderer); /** \brief Set the default properties for general image rendering. */ static void SetDefaultProperties(mitk::DataNode* node, mitk::BaseRenderer* renderer = NULL, bool overwrite = false); protected: /** \brief Transforms the actor to the actual position in 3D. * \param renderer The current renderer corresponding to the render window. */ void TransformActor(mitk::BaseRenderer* renderer); /** \brief Generates a plane according to the size of the resliced image in milimeters. * * \image html texturedPlane.png * * In VTK a vtkPlaneSource is defined through three points. The origin and two * points defining the axes of the plane (see VTK documentation). The origin is * set to (xMin; yMin; Z), where xMin and yMin are the minimal bounds of the * resliced image in space. Z is relevant for blending and the layer property. * The center of the plane (C) is also the center of the view plane (cf. the image above). * * \note For the standard MITK view with three 2D render windows showing three * different slices, three such planes are generated. All these planes are generated * in the XY-plane (even if they depict a YZ-slice of the volume). * */ void GeneratePlane(mitk::BaseRenderer* renderer, vtkFloatingPointType planeBounds[6]); /** \brief Generates a vtkPolyData object containing the outline of a given binary slice. \param binarySlice - The binary image slice. (Volumes are not supported.) \param mmPerPixel - Spacing of the binary image slice. Hence it's 2D, only in x/y-direction. \note This code has been taken from the deprecated library iil. */ vtkSmartPointer CreateOutlinePolyData(mitk::BaseRenderer* renderer); /** Default constructor */ ImageVtkMapper2D(); /** Default deconstructor */ virtual ~ImageVtkMapper2D(); /** \brief Does the actual resampling, without rendering the image yet. * All the data is generated inside this method. The vtkProp (or Actor) * is filled with content (i.e. the resliced image). * * After generation, a 4x4 transformation matrix(t) of the current slice is obtained * from the vtkResliceImage object via GetReslicesAxis(). This matrix is * applied to each textured plane (actor->SetUserTransform(t)) to transform everything * to the actual 3D position (cf. the following image). * * \image html cameraPositioning3D.png * */ virtual void GenerateDataForRenderer(mitk::BaseRenderer *renderer); - /** \brief Internal helper method for intersection testing used only in CalculateClippedPlaneBounds() */ - bool LineIntersectZero( vtkPoints *points, int p1, int p2, - vtkFloatingPointType *bounds ); - - /** \brief Calculate the bounding box of the resliced image. This is necessary for - arbitrarily rotated planes in an image volume. A rotated plane (e.g. in swivel mode) - will have a new bounding box, which needs to be calculated. */ - bool CalculateClippedPlaneBounds( const Geometry3D *boundingGeometry, - const PlaneGeometry *planeGeometry, vtkFloatingPointType *bounds ); /** \brief This method uses the vtkCamera clipping range and the layer property * to calcualte the depth of the object (e.g. image or contour). The depth is used * to keep the correct order for the final VTK rendering.*/ float CalculateLayerDepth(mitk::BaseRenderer* renderer); /** \brief This method applies a level window on RBG(A) images. * It should only be called for internally for RGB(A) images. */ void ApplyRBGALevelWindow( mitk::BaseRenderer* renderer ); /** \brief This method applies (or modifies) the lookuptable for all types of images. */ void ApplyLookuptable( mitk::BaseRenderer* renderer ); /** \brief This method applies a color transfer function, if no LookuptableProperty is set. Internally, a vtkColorTransferFunction is used. This is usefull for coloring continous images (e.g. float) */ void ApplyColorTransferFunction(mitk::BaseRenderer* renderer); /** \brief Set the color of the image/polydata */ void ApplyColor( mitk::BaseRenderer* renderer ); /** \brief Set the opacity of the actor. */ void ApplyOpacity( mitk::BaseRenderer* renderer ); }; } // namespace mitk #endif /* MITKIMAGEVTKMAPPER2D_H_HEADER_INCLUDED_C10E906E */