Medical Imaging Interaction Toolkit  2018.4.99-f51274ea
Medical Imaging Interaction Toolkit
mitkGeometryClipImageFilter.cpp
Go to the documentation of this file.
1 /*============================================================================
2 
3 The Medical Imaging Interaction Toolkit (MITK)
4 
5 Copyright (c) German Cancer Research Center (DKFZ)
6 All rights reserved.
7 
8 Use of this source code is governed by a 3-clause BSD license that can be
9 found in the LICENSE file.
10 
11 ============================================================================*/
12 
14 #include "mitkImageTimeSelector.h"
15 #include "mitkProperties.h"
16 #include "mitkTimeHelper.h"
17 
18 #include "mitkImageToItk.h"
19 
20 #include "itkImageRegionConstIterator.h"
21 #include "itkImageRegionIteratorWithIndex.h"
22 
23 #include <limits>
24 
26  : m_ClippingGeometry(nullptr),
27  m_ClipPartAboveGeometry(true),
28  m_OutsideValue(0),
29  m_AutoOutsideValue(false),
30  m_LabelBothSides(false),
31  m_AutoOrientLabels(false),
32  m_AboveGeometryLabel(1),
33  m_BelowGeometryLabel(2)
34 {
35  this->SetNumberOfIndexedInputs(2);
36  this->SetNumberOfRequiredInputs(2);
40 }
41 
43 {
44 }
45 
47 {
48  m_TimeClippingGeometry = timeClippingGeometry;
49  SetClippingGeometry(timeClippingGeometry->GetGeometryForTimeStep(0));
50 }
51 
53 {
54  if (aClippingGeometry != m_ClippingGeometry.GetPointer())
55  {
56  m_ClippingGeometry = aClippingGeometry;
57  m_ClippingGeometryData->SetGeometry(const_cast<mitk::BaseGeometry *>(aClippingGeometry));
58  SetNthInput(1, m_ClippingGeometryData);
59  Modified();
60  }
61 }
62 
64 {
65  return m_ClippingGeometry;
66 }
67 
69 {
71 }
72 
74 {
75  Superclass::GenerateInputRequestedRegion();
76 
77  mitk::Image *output = this->GetOutput();
78  mitk::Image *input = this->GetInput();
79  if ((output->IsInitialized() == false) || (m_ClippingGeometry.IsNull()))
80  return;
81 
83 
84  GenerateTimeInInputRegion(output, input);
85 }
86 
88 {
89  mitk::Image::ConstPointer input = this->GetInput();
90  mitk::Image::Pointer output = this->GetOutput();
91 
92  if ((output->IsInitialized()) && (this->GetMTime() <= m_TimeOfHeaderInitialization.GetMTime()))
93  return;
94 
95  itkDebugMacro(<< "GenerateOutputInformation()");
96 
97  unsigned int i;
98  auto tmpDimensions = new unsigned int[input->GetDimension()];
99 
100  for (i = 0; i < input->GetDimension(); ++i)
101  tmpDimensions[i] = input->GetDimension(i);
102 
103  output->Initialize(input->GetPixelType(), input->GetDimension(), tmpDimensions, input->GetNumberOfChannels());
104 
105  delete[] tmpDimensions;
106 
107  output->SetGeometry(static_cast<mitk::BaseGeometry *>(input->GetGeometry()->Clone().GetPointer()));
108 
109  output->SetPropertyList(input->GetPropertyList()->Clone());
110 
111  m_TimeOfHeaderInitialization.Modified();
112 }
113 
114 template <typename TPixel, unsigned int VImageDimension>
115 void mitk::GeometryClipImageFilter::_InternalComputeClippedImage(itk::Image<TPixel, VImageDimension> *inputItkImage,
116  mitk::GeometryClipImageFilter *geometryClipper,
117  const mitk::PlaneGeometry *clippingPlaneGeometry)
118 {
119  typedef itk::Image<TPixel, VImageDimension> ItkInputImageType;
120  typedef itk::Image<TPixel, VImageDimension> ItkOutputImageType;
121 
122  typedef itk::ImageRegionConstIteratorWithIndex<ItkInputImageType> ItkInputImageIteratorType;
123  typedef itk::ImageRegionIteratorWithIndex<ItkOutputImageType> ItkOutputImageIteratorType;
124 
126  outputimagetoitk->SetInput(geometryClipper->m_OutputTimeSelector->GetOutput());
127  outputimagetoitk->Update();
128  typename ItkOutputImageType::Pointer outputItkImage = outputimagetoitk->GetOutput();
129 
130  // create the iterators
131  typename ItkInputImageType::RegionType inputRegionOfInterest = inputItkImage->GetLargestPossibleRegion();
132  ItkInputImageIteratorType inputIt(inputItkImage, inputRegionOfInterest);
133  ItkOutputImageIteratorType outputIt(outputItkImage, inputRegionOfInterest);
134 
135  typename ItkOutputImageType::PixelType outsideValue;
136  if (geometryClipper->m_AutoOutsideValue)
138  else
139  outsideValue = (typename ItkOutputImageType::PixelType)geometryClipper->m_OutsideValue;
140 
141  mitk::BaseGeometry *inputGeometry = geometryClipper->m_InputTimeSelector->GetOutput()->GetGeometry();
142  typedef itk::Index<VImageDimension> IndexType;
143  Point3D indexPt;
144  indexPt.Fill(0);
145  int i, dim = IndexType::GetIndexDimension();
146  Point3D pointInMM;
147  bool above = geometryClipper->m_ClipPartAboveGeometry;
148  bool labelBothSides = geometryClipper->GetLabelBothSides();
149 
150  if (geometryClipper->GetAutoOrientLabels())
151  {
152  Point3D leftMostPoint;
153  leftMostPoint.Fill(std::numeric_limits<float>::min() / 2.0);
154  if (clippingPlaneGeometry->IsAbove(pointInMM) != above)
155  {
156  // invert meaning of above --> left is always the "above" side
157  above = !above;
158  MITK_INFO << leftMostPoint << " is BELOW geometry. Inverting meaning of above" << std::endl;
159  }
160  else
161  MITK_INFO << leftMostPoint << " is above geometry" << std::endl;
162  }
163 
164  auto aboveLabel =
165  (typename ItkOutputImageType::PixelType)geometryClipper->GetAboveGeometryLabel();
166  auto belowLabel =
167  (typename ItkOutputImageType::PixelType)geometryClipper->GetBelowGeometryLabel();
168 
169  for (inputIt.GoToBegin(), outputIt.GoToBegin(); !inputIt.IsAtEnd(); ++inputIt, ++outputIt)
170  {
171  if ((typename ItkOutputImageType::PixelType)inputIt.Get() == outsideValue)
172  {
173  outputIt.Set(outsideValue);
174  }
175  else
176  {
177  for (i = 0; i < dim; ++i)
178  indexPt[i] = (mitk::ScalarType)inputIt.GetIndex()[i];
179  inputGeometry->IndexToWorld(indexPt, pointInMM);
180  if (clippingPlaneGeometry->IsAbove(pointInMM) == above)
181  {
182  if (labelBothSides)
183  outputIt.Set(aboveLabel);
184  else
185  outputIt.Set(outsideValue);
186  }
187  else
188  {
189  if (labelBothSides)
190  outputIt.Set(belowLabel);
191  else
192  outputIt.Set(inputIt.Get());
193  }
194  }
195  }
196 }
197 
198 #include "mitkImageAccessByItk.h"
199 
201 {
202  Image::ConstPointer input = this->GetInput();
203  Image::Pointer output = this->GetOutput();
204 
205  if ((output->IsInitialized() == false) || (m_ClippingGeometry.IsNull()))
206  return;
207 
208  const PlaneGeometry *clippingGeometryOfCurrentTimeStep = nullptr;
209 
210  if (m_TimeClippingGeometry.IsNull())
211  {
212  clippingGeometryOfCurrentTimeStep = dynamic_cast<const PlaneGeometry *>(m_ClippingGeometry.GetPointer());
213  }
214  else
215  {
216  clippingGeometryOfCurrentTimeStep =
217  dynamic_cast<const PlaneGeometry *>(m_TimeClippingGeometry->GetGeometryForTimeStep(0).GetPointer());
218  }
219 
220  if (clippingGeometryOfCurrentTimeStep == nullptr)
221  return;
222 
223  m_InputTimeSelector->SetInput(input);
224  m_OutputTimeSelector->SetInput(this->GetOutput());
225 
226  mitk::Image::RegionType outputRegion = output->GetRequestedRegion();
227  const mitk::TimeGeometry *outputTimeGeometry = output->GetTimeGeometry();
228  const mitk::TimeGeometry *inputTimeGeometry = input->GetTimeGeometry();
229  ScalarType timeInMS;
230 
231  int timestep = 0;
232  int tstart = outputRegion.GetIndex(3);
233  int tmax = tstart + outputRegion.GetSize(3);
234 
235  int t;
236  for (t = tstart; t < tmax; ++t)
237  {
238  timeInMS = outputTimeGeometry->TimeStepToTimePoint(t);
239  timestep = inputTimeGeometry->TimePointToTimeStep(timeInMS);
240 
241  m_InputTimeSelector->SetTimeNr(timestep);
242  m_InputTimeSelector->UpdateLargestPossibleRegion();
243  m_OutputTimeSelector->SetTimeNr(t);
244  m_OutputTimeSelector->UpdateLargestPossibleRegion();
245 
246  if (m_TimeClippingGeometry.IsNotNull())
247  {
248  timestep = m_TimeClippingGeometry->TimePointToTimeStep(timeInMS);
249  if (m_TimeClippingGeometry->IsValidTimeStep(timestep) == false)
250  continue;
251 
252  clippingGeometryOfCurrentTimeStep =
253  dynamic_cast<const PlaneGeometry *>(m_TimeClippingGeometry->GetGeometryForTimeStep(timestep).GetPointer());
254  }
255 
257  m_InputTimeSelector->GetOutput(), _InternalComputeClippedImage, this, clippingGeometryOfCurrentTimeStep);
258  }
259 
260  m_TimeOfHeaderInitialization.Modified();
261 }
void SetRequestedRegionToLargestPossibleRegion() override
mitk::ImageTimeSelector::Pointer m_InputTimeSelector
#define MITK_INFO
Definition: mitkLogMacros.h:18
virtual bool GetLabelBothSides() const
double ScalarType
const mitk::BaseGeometry * GetClippingGeometry() const
const mitk::TimeGeometry * GetClippingTimeGeometry() const
virtual bool GetAutoOrientLabels() const
mitk::ImageTimeSelector::Pointer m_OutputTimeSelector
void SetClippingGeometry(const mitk::BaseGeometry *aClippingGeometry)
mitk::GeometryData::Pointer m_ClippingGeometryData
virtual bool IsAbove(const Point3D &pt3d_mm, bool considerBoundingBox=false) const
Calculates, whether a point is below or above the plane. There are two different calculation methods...
itk::ImageRegion< RegionDimension > RegionType
virtual TimeStepType TimePointToTimeStep(TimePointType timePoint) const =0
Converts a time point to the corresponding time step.
Image class for storing images.
Definition: mitkImage.h:72
virtual TimePointType TimeStepToTimePoint(TimeStepType timeStep) const =0
Converts a time step to a time point.
static Pointer New()
void GenerateTimeInInputRegion(const mitk::TimeGeometry *outputTimeGeometry, const TOutputRegion &outputRegion, const mitk::TimeGeometry *inputTimeGeometry, TInputRegion &inputRegion)
static T min(T x, T y)
Definition: svm.cpp:53
static Pointer New()
InputImageType * GetInput(void)
void _InternalComputeClippedImage(itk::Image< TPixel, VImageDimension > *itkImage, mitk::GeometryClipImageFilter *geometryClipper, const mitk::PlaneGeometry *clippingPlaneGeometry)
virtual ScalarType GetBelowGeometryLabel() const
virtual bool IsInitialized() const
Check whether the data has been initialized, i.e., at least the Geometry and other header data has be...
mitk::BaseGeometry::ConstPointer m_ClippingGeometry
virtual BaseGeometry::Pointer GetGeometryForTimeStep(TimeStepType timeStep) const =0
Returns the geometry which corresponds to the given time step.
OutputType * GetOutput()
Get the output data of this image source object.
Describes a two-dimensional, rectangular plane.
#define AccessByItk_2(mitkImage, itkImageTypeFunction, arg1, arg2)
mitk::TimeGeometry::ConstPointer m_TimeClippingGeometry
virtual ScalarType GetAboveGeometryLabel() const
BaseGeometry Describes the geometry of a data object.
Filter for clipping an image with a PlaneGeometry.
static Pointer New()