Medical Imaging Interaction Toolkit  2016.11.0
Medical Imaging Interaction Toolkit
mitkFastMarchingTool.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,
6 Division of Medical and Biological Informatics.
7 All rights reserved.
8 
9 This software is distributed WITHOUT ANY WARRANTY; without
10 even the implied warranty of MERCHANTABILITY or FITNESS FOR
11 A PARTICULAR PURPOSE.
12 
13 See LICENSE.txt or http://www.mitk.org for details.
14 
15 ===================================================================*/
16 
17 #include "mitkFastMarchingTool.h"
18 #include "mitkToolManager.h"
19 
20 #include "mitkBaseRenderer.h"
21 #include "mitkInteractionConst.h"
22 #include "mitkRenderingManager.h"
23 
24 #include "itkOrImageFilter.h"
25 #include "mitkImageTimeSelector.h"
26 
27 // us
28 #include <usGetModuleContext.h>
29 #include <usModule.h>
30 #include <usModuleContext.h>
31 #include <usModuleResource.h>
32 
33 namespace mitk
34 {
36 }
37 
39  : FeedbackContourTool("PressMoveReleaseAndPointSetting"),
40  m_NeedUpdate(true),
41  m_CurrentTimeStep(0),
42  m_PositionEvent(0),
43  m_LowerThreshold(0),
44  m_UpperThreshold(200),
45  m_StoppingValue(100),
46  m_Sigma(1.0),
47  m_Alpha(-0.5),
48  m_Beta(3.0)
49 {
50 }
51 
53 {
54  if (this->m_SmoothFilter.IsNotNull())
55  this->m_SmoothFilter->RemoveAllObservers();
56 
57  if (this->m_SigmoidFilter.IsNotNull())
58  this->m_SigmoidFilter->RemoveAllObservers();
59 
60  if (this->m_GradientMagnitudeFilter.IsNotNull())
61  this->m_GradientMagnitudeFilter->RemoveAllObservers();
62 
63  if (this->m_FastMarchingFilter.IsNotNull())
64  this->m_FastMarchingFilter->RemoveAllObservers();
65 }
66 
68 {
69  CONNECT_FUNCTION("ShiftSecondaryButtonPressed", OnAddPoint);
70  CONNECT_FUNCTION("ShiftPrimaryButtonPressed", OnAddPoint);
71  CONNECT_FUNCTION("DeletePoint", OnDelete);
72 }
73 
74 const char **mitk::FastMarchingTool::GetXPM() const
75 {
76  return NULL; // mitkFastMarchingTool_xpm;
77 }
78 
80 {
82  us::ModuleResource resource = module->GetResource("FastMarching_48x48.png");
83  return resource;
84 }
85 
87 {
89  us::ModuleResource resource = module->GetResource("FastMarching_Cursor_32x32.png");
90  return resource;
91 }
92 
94 {
95  return "2D Fast Marching";
96 }
97 
99 {
100  m_ReferenceImageSliceAsITK = InternalImageType::New();
101 
102  m_ReferenceImageSlice = GetAffectedReferenceSlice(m_PositionEvent);
103  CastToItkImage(m_ReferenceImageSlice, m_ReferenceImageSliceAsITK);
104 
105  m_ProgressCommand = mitk::ToolCommand::New();
106 
107  m_SmoothFilter = SmoothingFilterType::New();
108  m_SmoothFilter->SetInput(m_ReferenceImageSliceAsITK);
109  m_SmoothFilter->SetTimeStep(0.05);
110  m_SmoothFilter->SetNumberOfIterations(2);
111  m_SmoothFilter->SetConductanceParameter(9.0);
112 
113  m_GradientMagnitudeFilter = GradientFilterType::New();
114  m_GradientMagnitudeFilter->SetSigma(m_Sigma);
115 
116  m_SigmoidFilter = SigmoidFilterType::New();
117  m_SigmoidFilter->SetAlpha(m_Alpha);
118  m_SigmoidFilter->SetBeta(m_Beta);
119  m_SigmoidFilter->SetOutputMinimum(0.0);
120  m_SigmoidFilter->SetOutputMaximum(1.0);
121 
122  m_FastMarchingFilter = FastMarchingFilterType::New();
123  m_FastMarchingFilter->SetStoppingValue(m_StoppingValue);
124 
125  m_ThresholdFilter = ThresholdingFilterType::New();
126  m_ThresholdFilter->SetLowerThreshold(m_LowerThreshold);
127  m_ThresholdFilter->SetUpperThreshold(m_UpperThreshold);
128  m_ThresholdFilter->SetOutsideValue(0);
129  m_ThresholdFilter->SetInsideValue(1.0);
130 
131  m_SeedContainer = NodeContainer::New();
132  m_SeedContainer->Initialize();
133  m_FastMarchingFilter->SetTrialPoints(m_SeedContainer);
134 
135  if (this->m_SmoothFilter.IsNotNull())
136  this->m_SmoothFilter->RemoveAllObservers();
137 
138  if (this->m_SigmoidFilter.IsNotNull())
139  this->m_SigmoidFilter->RemoveAllObservers();
140 
141  if (this->m_GradientMagnitudeFilter.IsNotNull())
142  this->m_GradientMagnitudeFilter->RemoveAllObservers();
143 
144  if (this->m_FastMarchingFilter.IsNotNull())
145  this->m_FastMarchingFilter->RemoveAllObservers();
146 
147  m_SmoothFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
148  m_GradientMagnitudeFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
149  m_SigmoidFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
150  m_FastMarchingFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
151 
152  m_SmoothFilter->SetInput(m_ReferenceImageSliceAsITK);
153  m_GradientMagnitudeFilter->SetInput(m_SmoothFilter->GetOutput());
154  m_SigmoidFilter->SetInput(m_GradientMagnitudeFilter->GetOutput());
155  m_FastMarchingFilter->SetInput(m_SigmoidFilter->GetOutput());
156  m_ThresholdFilter->SetInput(m_FastMarchingFilter->GetOutput());
157  m_ReferenceImageSliceAsITK = InternalImageType::New();
158 }
159 
161 {
162  if (m_UpperThreshold != value)
163  {
164  m_UpperThreshold = value / 10.0;
165  m_ThresholdFilter->SetUpperThreshold(m_UpperThreshold);
166  m_NeedUpdate = true;
167  }
168 }
169 
171 {
172  if (m_LowerThreshold != value)
173  {
174  m_LowerThreshold = value / 10.0;
175  m_ThresholdFilter->SetLowerThreshold(m_LowerThreshold);
176  m_NeedUpdate = true;
177  }
178 }
179 
181 {
182  if (m_Beta != value)
183  {
184  m_Beta = value;
185  m_SigmoidFilter->SetBeta(m_Beta);
186  m_NeedUpdate = true;
187  }
188 }
189 
191 {
192  if (m_Sigma != value)
193  {
194  m_Sigma = value;
195  m_GradientMagnitudeFilter->SetSigma(m_Sigma);
196  m_NeedUpdate = true;
197  }
198 }
199 
201 {
202  if (m_Alpha != value)
203  {
204  m_Alpha = value;
205  m_SigmoidFilter->SetAlpha(m_Alpha);
206  m_NeedUpdate = true;
207  }
208 }
209 
211 {
212  if (m_StoppingValue != value)
213  {
214  m_StoppingValue = value;
215  m_FastMarchingFilter->SetStoppingValue(m_StoppingValue);
216  m_NeedUpdate = true;
217  }
218 }
219 
221 {
222  Superclass::Activated();
223 
224  m_ResultImageNode = mitk::DataNode::New();
225  m_ResultImageNode->SetName("FastMarching_Preview");
226  m_ResultImageNode->SetBoolProperty("helper object", true);
227  m_ResultImageNode->SetColor(0.0, 1.0, 0.0);
228  m_ResultImageNode->SetVisibility(true);
229  m_ToolManager->GetDataStorage()->Add(this->m_ResultImageNode, m_ToolManager->GetReferenceData(0));
230 
231  m_SeedsAsPointSet = mitk::PointSet::New();
232  m_SeedsAsPointSetNode = mitk::DataNode::New();
233  m_SeedsAsPointSetNode->SetData(m_SeedsAsPointSet);
234  m_SeedsAsPointSetNode->SetName("Seeds_Preview");
235  m_SeedsAsPointSetNode->SetBoolProperty("helper object", true);
236  m_SeedsAsPointSetNode->SetColor(0.0, 1.0, 0.0);
237  m_SeedsAsPointSetNode->SetVisibility(true);
238  m_ToolManager->GetDataStorage()->Add(this->m_SeedsAsPointSetNode, m_ToolManager->GetReferenceData(0));
239 
240  this->Initialize();
241 }
242 
244 {
245  m_ToolManager->GetDataStorage()->Remove(this->m_ResultImageNode);
246  m_ToolManager->GetDataStorage()->Remove(this->m_SeedsAsPointSetNode);
247  this->ClearSeeds();
248  m_ResultImageNode = NULL;
249  m_SeedsAsPointSetNode = NULL;
251 
252  Superclass::Deactivated();
253 }
254 
256 {
257  m_ReferenceImage = dynamic_cast<mitk::Image *>(m_ToolManager->GetReferenceData(0)->GetData());
258  if (m_ReferenceImage->GetTimeGeometry()->CountTimeSteps() > 1)
259  {
261  timeSelector->SetInput(m_ReferenceImage);
262  timeSelector->SetTimeNr(m_CurrentTimeStep);
263  timeSelector->UpdateLargestPossibleRegion();
264  m_ReferenceImage = timeSelector->GetOutput();
265  }
266  m_NeedUpdate = true;
267 }
268 
270 {
271  // combine preview image with current working segmentation
272  if (dynamic_cast<mitk::Image *>(m_ResultImageNode->GetData()))
273  {
274  // logical or combination of preview and segmentation slice
275 
276  mitk::Image::Pointer workingImageSlice;
277  mitk::Image::Pointer workingImage = dynamic_cast<mitk::Image *>(this->m_ToolManager->GetWorkingData(0)->GetData());
278  workingImageSlice = GetAffectedImageSliceAs2DImage(m_WorkingPlane, workingImage, m_CurrentTimeStep);
279 
280  mitk::Image::Pointer segmentationResult = mitk::Image::New();
281 
282  bool isDeprecatedUnsignedCharSegmentation =
283  (workingImage->GetPixelType().GetComponentType() == itk::ImageIOBase::UCHAR);
284 
285  if (isDeprecatedUnsignedCharSegmentation)
286  {
287  typedef itk::Image<unsigned char, 2> OutputUCharImageType;
288  OutputUCharImageType::Pointer workingImageSliceInITK = OutputUCharImageType::New();
289 
290  CastToItkImage(workingImageSlice, workingImageSliceInITK);
291 
292  typedef itk::OrImageFilter<OutputImageType, OutputUCharImageType, OutputUCharImageType> OrImageFilterType;
294 
295  orFilter->SetInput1(m_ThresholdFilter->GetOutput());
296  orFilter->SetInput2(workingImageSliceInITK);
297  orFilter->Update();
298 
299  mitk::CastToMitkImage(orFilter->GetOutput(), segmentationResult);
300  }
301  else
302  {
303  OutputImageType::Pointer workingImageSliceInITK = OutputImageType::New();
304 
305  CastToItkImage(workingImageSlice, workingImageSliceInITK);
306 
307  typedef itk::OrImageFilter<OutputImageType, OutputImageType> OrImageFilterType;
309 
310  orFilter->SetInput(0, m_ThresholdFilter->GetOutput());
311  orFilter->SetInput(1, workingImageSliceInITK);
312  orFilter->Update();
313 
314  mitk::CastToMitkImage(orFilter->GetOutput(), segmentationResult);
315  }
316 
317  segmentationResult->GetGeometry()->SetOrigin(workingImageSlice->GetGeometry()->GetOrigin());
318  segmentationResult->GetGeometry()->SetIndexToWorldTransform(
319  workingImageSlice->GetGeometry()->GetIndexToWorldTransform());
320 
321  // write to segmentation volume and hide preview image
322  // again, current time step is not considered
323  this->WriteBackSegmentationResult(m_WorkingPlane, segmentationResult, m_CurrentTimeStep);
324  this->m_ResultImageNode->SetVisibility(false);
325 
326  this->ClearSeeds();
327  }
328 
330  m_ToolManager->ActivateTool(-1);
331 }
332 
334 {
335  // Add a new seed point for FastMarching algorithm
336  mitk::InteractionPositionEvent *positionEvent = dynamic_cast<mitk::InteractionPositionEvent *>(interactionEvent);
337  // const PositionEvent* p = dynamic_cast<const PositionEvent*>(stateEvent->GetEvent());
338  if (positionEvent == NULL)
339  return;
340 
341  if (m_PositionEvent.IsNotNull())
342  m_PositionEvent = NULL;
343 
344  m_PositionEvent =
345  InteractionPositionEvent::New(positionEvent->GetSender(), positionEvent->GetPointerPositionOnScreen());
346 
347  // if click was on another renderwindow or slice then reset pipeline and preview
348  if ((m_LastEventSender != m_PositionEvent->GetSender()) ||
349  (m_LastEventSlice != m_PositionEvent->GetSender()->GetSlice()))
350  {
351  this->BuildITKPipeline();
352  this->ClearSeeds();
353  }
354 
355  m_LastEventSender = m_PositionEvent->GetSender();
356  m_LastEventSlice = m_LastEventSender->GetSlice();
357  m_WorkingPlane = positionEvent->GetSender()->GetCurrentWorldPlaneGeometry()->Clone();
358 
359  mitk::Point3D clickInIndex;
360 
361  m_ReferenceImageSlice->GetGeometry()->WorldToIndex(m_PositionEvent->GetPositionInWorld(), clickInIndex);
362  itk::Index<2> seedPosition;
363  seedPosition[0] = clickInIndex[0];
364  seedPosition[1] = clickInIndex[1];
365 
366  NodeType node;
367  const double seedValue = 0.0;
368  node.SetValue(seedValue);
369  node.SetIndex(seedPosition);
370  this->m_SeedContainer->InsertElement(this->m_SeedContainer->Size(), node);
371  m_FastMarchingFilter->Modified();
372 
373  m_SeedsAsPointSet->InsertPoint(m_SeedsAsPointSet->GetSize(), m_PositionEvent->GetPositionInWorld());
374 
376 
377  m_NeedUpdate = true;
378 
379  this->Update();
380 
381  m_ReadyMessage.Send();
382 }
383 
385 {
386  // delete last seed point
387  if (!(this->m_SeedContainer->empty()))
388  {
389  // delete last element of seeds container
390  this->m_SeedContainer->pop_back();
391  m_FastMarchingFilter->Modified();
392 
393  // delete last point in pointset - somehow ugly
394  m_SeedsAsPointSet->GetPointSet()->GetPoints()->DeleteIndex(m_SeedsAsPointSet->GetSize() - 1);
395 
397 
398  m_NeedUpdate = true;
399 
400  this->Update();
401  }
402 }
403 
405 {
406  const unsigned int progress_steps = 20;
407 
408  // update FastMarching pipeline and show result
409  if (m_NeedUpdate)
410  {
411  m_ProgressCommand->AddStepsToDo(progress_steps);
412  CurrentlyBusy.Send(true);
413  try
414  {
415  m_ThresholdFilter->Update();
416  }
417  catch (itk::ExceptionObject &excep)
418  {
419  MITK_ERROR << "Exception caught: " << excep.GetDescription();
420 
421  // progress by max step count, will force
422  m_ProgressCommand->SetProgress(progress_steps);
423  CurrentlyBusy.Send(false);
424 
425  std::string msg = excep.GetDescription();
426  ErrorMessage.Send(msg);
427 
428  return;
429  }
430  m_ProgressCommand->SetProgress(progress_steps);
431  CurrentlyBusy.Send(false);
432 
433  // make output visible
435  CastToMitkImage(m_ThresholdFilter->GetOutput(), result);
436  result->GetGeometry()->SetOrigin(m_ReferenceImageSlice->GetGeometry()->GetOrigin());
437  result->GetGeometry()->SetIndexToWorldTransform(m_ReferenceImageSlice->GetGeometry()->GetIndexToWorldTransform());
438  m_ResultImageNode->SetData(result);
439  m_ResultImageNode->SetVisibility(true);
441  }
442 }
443 
445 {
446  // clear seeds for FastMarching as well as the PointSet for visualization
447  if (this->m_SeedContainer.IsNotNull())
448  this->m_SeedContainer->Initialize();
449 
450  if (this->m_SeedsAsPointSet.IsNotNull())
451  {
452  this->m_SeedsAsPointSet = mitk::PointSet::New();
453  this->m_SeedsAsPointSetNode->SetData(this->m_SeedsAsPointSet);
454  m_SeedsAsPointSetNode->SetName("Seeds_Preview");
455  m_SeedsAsPointSetNode->SetBoolProperty("helper object", true);
456  m_SeedsAsPointSetNode->SetColor(0.0, 1.0, 0.0);
457  m_SeedsAsPointSetNode->SetVisibility(true);
458  }
459 
460  if (this->m_FastMarchingFilter.IsNotNull())
461  m_FastMarchingFilter->Modified();
462 
463  this->m_NeedUpdate = true;
464 }
465 
467 {
468  // clear all seeds and preview empty result
469  this->ClearSeeds();
470 
471  m_ResultImageNode->SetVisibility(false);
472 
474 }
475 
477 {
478  if (m_CurrentTimeStep != t)
479  {
480  m_CurrentTimeStep = t;
481 
482  this->Initialize();
483  }
484 }
virtual us::ModuleResource GetCursorIconResource() const override
Returns the path of a cursor icon.
Super class for all position events.
itk::SmartPointer< Self > Pointer
BaseRenderer * GetSender() const
#define MITK_ERROR
Definition: mitkLogMacros.h:24
static Pointer New()
virtual void OnDelete(StateMachineAction *, InteractionEvent *interactionEvent)
Delete action of StateMachine pattern.
void Update()
Updates the itk pipeline and shows the result of FastMarching.
static void Update(vtkPolyData *)
Definition: mitkSurface.cpp:35
#define MITKSEGMENTATION_EXPORT
DataCollection - Class to facilitate loading/accessing structured data.
#define MITK_TOOL_MACRO(EXPORT_SPEC, CLASS_NAME, DESCRIPTION)
virtual const PlaneGeometry * GetCurrentWorldPlaneGeometry()
Get the current 2D-worldgeometry (m_CurrentWorldPlaneGeometry) used for 2D-rendering.
Constants for most interaction classes, due to the generic StateMachines.
Pointer Clone() const
static Pointer New(BaseRenderer *_arga, const Point2D &_argb)
virtual void SetCurrentTimeStep(int t)
Set the working time step.
static Pointer New()
us::ModuleResource GetIconResource() const override
Returns the tool button icon of the tool wrapped by a usModuleResource.
void Reset()
Reset all relevant inputs of the itk pipeline.
virtual void OnAddPoint(StateMachineAction *, InteractionEvent *interactionEvent)
Add point action of StateMachine pattern.
void SetLowerThreshold(double)
Set parameter used in Threshold filter.
void SetUpperThreshold(double)
Set parameter used in Threshold filter.
static Pointer New()
static RenderingManager * GetInstance()
Represents an action, that is executed after a certain event (in statemachine-mechanism) TODO: implem...
Module * GetModule() const
virtual const char * GetName() const override
Returns the name of this tool. Make it short!
Image class for storing images.
Definition: mitkImage.h:76
void ConnectActionsAndFunctions() override
static Pointer New()
void SetSigma(double)
Set parameter used in Gradient Magnitude filter.
void SetAlpha(double)
Set parameter used in Fast Marching filter.
virtual void ConfirmSegmentation()
Adds the feedback image to the current working image.
Base class for tools that use a contour for feedback.
FastMarchingFilterType::NodeType NodeType
void CastToMitkImage(const itk::SmartPointer< ItkOutputImageType > &itkimage, itk::SmartPointer< mitk::Image > &mitkoutputimage)
Cast an itk::Image (with a specific type) to an mitk::Image.
Definition: mitkImageCast.h:78
void MITKCORE_EXPORT CastToItkImage(const mitk::Image *mitkImage, itk::SmartPointer< ItkOutputImageType > &itkOutputImage)
Cast an mitk::Image to an itk::Image with a specific type.
void ClearSeeds()
Clear all seed points.
#define CONNECT_FUNCTION(a, f)
void SetStoppingValue(double)
Set parameter used in Fast Marching filter.
ModuleResource GetResource(const std::string &path) const
Definition: usModule.cpp:267
static ModuleContext * GetModuleContext()
Returns the module context of the calling module.
FastMarching semgentation tool.
void SetBeta(double)
Set parameter used in Fast Marching filter.
void RequestUpdateAll(RequestType type=REQUEST_UPDATE_ALL)
static Pointer New()
virtual void Activated() override
Called when the tool gets activated.
virtual const char ** GetXPM() const override
Returns an icon in the XPM format.
virtual void Deactivated() override
Called when the tool gets deactivated.
static itkEventMacro(BoundingShapeInteractionEvent, itk::AnyEvent) class MITKBOUNDINGSHAPE_EXPORT BoundingShapeInteractor Pointer New()
Basic interaction methods for mitk::GeometryData.