Medical Imaging Interaction Toolkit  2018.4.99-b585543d
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 (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 
13 #include "mitkFastMarchingTool.h"
14 #include "mitkToolManager.h"
15 
16 #include "mitkBaseRenderer.h"
17 #include "mitkInteractionConst.h"
18 #include "mitkRenderingManager.h"
19 
20 #include "itkOrImageFilter.h"
21 #include "mitkImageTimeSelector.h"
22 
23 // us
24 #include <usGetModuleContext.h>
25 #include <usModule.h>
26 #include <usModuleContext.h>
27 #include <usModuleResource.h>
28 
29 namespace mitk
30 {
32 }
33 
35  : FeedbackContourTool("FastMarchingTool"),
36  m_NeedUpdate(true),
37  m_CurrentTimeStep(0),
38  m_PositionEvent(nullptr),
39  m_LowerThreshold(0),
40  m_UpperThreshold(200),
41  m_StoppingValue(100),
42  m_Sigma(1.0),
43  m_Alpha(-0.5),
44  m_Beta(3.0)
45 {
46 }
47 
49 {
50  if (this->m_SmoothFilter.IsNotNull())
51  this->m_SmoothFilter->RemoveAllObservers();
52 
53  if (this->m_SigmoidFilter.IsNotNull())
54  this->m_SigmoidFilter->RemoveAllObservers();
55 
56  if (this->m_GradientMagnitudeFilter.IsNotNull())
57  this->m_GradientMagnitudeFilter->RemoveAllObservers();
58 
59  if (this->m_FastMarchingFilter.IsNotNull())
60  this->m_FastMarchingFilter->RemoveAllObservers();
61 }
62 
64 {
65  CONNECT_FUNCTION("ShiftSecondaryButtonPressed", OnAddPoint);
66  CONNECT_FUNCTION("ShiftPrimaryButtonPressed", OnAddPoint);
67  CONNECT_FUNCTION("DeletePoint", OnDelete);
68 }
69 
70 const char **mitk::FastMarchingTool::GetXPM() const
71 {
72  return nullptr; // mitkFastMarchingTool_xpm;
73 }
74 
76 {
78  us::ModuleResource resource = module->GetResource("FastMarching_48x48.png");
79  return resource;
80 }
81 
83 {
85  us::ModuleResource resource = module->GetResource("FastMarching_Cursor_32x32.png");
86  return resource;
87 }
88 
90 {
91  return "2D Fast Marching";
92 }
93 
95 {
96  m_ReferenceImageSliceAsITK = InternalImageType::New();
97 
100 
102 
103  m_SmoothFilter = SmoothingFilterType::New();
105  m_SmoothFilter->SetTimeStep(0.05);
106  m_SmoothFilter->SetNumberOfIterations(2);
107  m_SmoothFilter->SetConductanceParameter(9.0);
108 
109  m_GradientMagnitudeFilter = GradientFilterType::New();
111 
112  m_SigmoidFilter = SigmoidFilterType::New();
113  m_SigmoidFilter->SetAlpha(m_Alpha);
114  m_SigmoidFilter->SetBeta(m_Beta);
115  m_SigmoidFilter->SetOutputMinimum(0.0);
116  m_SigmoidFilter->SetOutputMaximum(1.0);
117 
118  m_FastMarchingFilter = FastMarchingFilterType::New();
119  m_FastMarchingFilter->SetStoppingValue(m_StoppingValue);
120 
121  m_ThresholdFilter = ThresholdingFilterType::New();
122  m_ThresholdFilter->SetLowerThreshold(m_LowerThreshold);
123  m_ThresholdFilter->SetUpperThreshold(m_UpperThreshold);
124  m_ThresholdFilter->SetOutsideValue(0);
125  m_ThresholdFilter->SetInsideValue(1.0);
126 
127  m_SeedContainer = NodeContainer::New();
128  m_SeedContainer->Initialize();
129  m_FastMarchingFilter->SetTrialPoints(m_SeedContainer);
130 
131  if (this->m_SmoothFilter.IsNotNull())
132  this->m_SmoothFilter->RemoveAllObservers();
133 
134  if (this->m_SigmoidFilter.IsNotNull())
135  this->m_SigmoidFilter->RemoveAllObservers();
136 
137  if (this->m_GradientMagnitudeFilter.IsNotNull())
138  this->m_GradientMagnitudeFilter->RemoveAllObservers();
139 
140  if (this->m_FastMarchingFilter.IsNotNull())
141  this->m_FastMarchingFilter->RemoveAllObservers();
142 
143  m_SmoothFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
144  m_GradientMagnitudeFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
145  m_SigmoidFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
146  m_FastMarchingFilter->AddObserver(itk::ProgressEvent(), m_ProgressCommand);
147 
149  m_GradientMagnitudeFilter->SetInput(m_SmoothFilter->GetOutput());
150  m_SigmoidFilter->SetInput(m_GradientMagnitudeFilter->GetOutput());
151  m_FastMarchingFilter->SetInput(m_SigmoidFilter->GetOutput());
152  m_ThresholdFilter->SetInput(m_FastMarchingFilter->GetOutput());
153  m_ReferenceImageSliceAsITK = InternalImageType::New();
154 }
155 
157 {
158  if (m_UpperThreshold != value)
159  {
160  m_UpperThreshold = value / 10.0;
161  m_ThresholdFilter->SetUpperThreshold(m_UpperThreshold);
162  m_NeedUpdate = true;
163  }
164 }
165 
167 {
168  if (m_LowerThreshold != value)
169  {
170  m_LowerThreshold = value / 10.0;
171  m_ThresholdFilter->SetLowerThreshold(m_LowerThreshold);
172  m_NeedUpdate = true;
173  }
174 }
175 
177 {
178  if (m_Beta != value)
179  {
180  m_Beta = value;
181  m_SigmoidFilter->SetBeta(m_Beta);
182  m_NeedUpdate = true;
183  }
184 }
185 
187 {
188  if (m_Sigma != value)
189  {
190  m_Sigma = value;
192  m_NeedUpdate = true;
193  }
194 }
195 
197 {
198  if (m_Alpha != value)
199  {
200  m_Alpha = value;
201  m_SigmoidFilter->SetAlpha(m_Alpha);
202  m_NeedUpdate = true;
203  }
204 }
205 
207 {
208  if (m_StoppingValue != value)
209  {
210  m_StoppingValue = value;
211  m_FastMarchingFilter->SetStoppingValue(m_StoppingValue);
212  m_NeedUpdate = true;
213  }
214 }
215 
217 {
218  Superclass::Activated();
219 
221  m_ResultImageNode->SetName("FastMarching_Preview");
222  m_ResultImageNode->SetBoolProperty("helper object", true);
223  m_ResultImageNode->SetColor(0.0, 1.0, 0.0);
224  m_ResultImageNode->SetVisibility(true);
226 
230  m_SeedsAsPointSetNode->SetName("Seeds_Preview");
231  m_SeedsAsPointSetNode->SetBoolProperty("helper object", true);
232  m_SeedsAsPointSetNode->SetColor(0.0, 1.0, 0.0);
233  m_SeedsAsPointSetNode->SetVisibility(true);
235 
236  this->Initialize();
237 }
238 
240 {
243  this->ClearSeeds();
244  m_ResultImageNode = nullptr;
245  m_SeedsAsPointSetNode = nullptr;
247 
248  Superclass::Deactivated();
249 }
250 
252 {
253  m_ReferenceImage = dynamic_cast<mitk::Image *>(m_ToolManager->GetReferenceData(0)->GetData());
254  if (m_ReferenceImage->GetTimeGeometry()->CountTimeSteps() > 1)
255  {
257  timeSelector->SetInput(m_ReferenceImage);
258  timeSelector->SetTimeNr(m_CurrentTimeStep);
259  timeSelector->UpdateLargestPossibleRegion();
260  m_ReferenceImage = timeSelector->GetOutput();
261  }
262  m_NeedUpdate = true;
263 }
264 
266 {
267  // combine preview image with current working segmentation
268  if (dynamic_cast<mitk::Image *>(m_ResultImageNode->GetData()))
269  {
270  // logical or combination of preview and segmentation slice
271 
272  mitk::Image::Pointer workingImageSlice;
273  mitk::Image::Pointer workingImage = dynamic_cast<mitk::Image *>(this->m_ToolManager->GetWorkingData(0)->GetData());
274  workingImageSlice = GetAffectedImageSliceAs2DImage(m_WorkingPlane, workingImage, m_CurrentTimeStep);
275 
276  mitk::Image::Pointer segmentationResult = mitk::Image::New();
277 
278  bool isDeprecatedUnsignedCharSegmentation =
279  (workingImage->GetPixelType().GetComponentType() == itk::ImageIOBase::UCHAR);
280 
281  if (isDeprecatedUnsignedCharSegmentation)
282  {
283  typedef itk::Image<unsigned char, 2> OutputUCharImageType;
284  OutputUCharImageType::Pointer workingImageSliceInITK = OutputUCharImageType::New();
285 
286  CastToItkImage(workingImageSlice, workingImageSliceInITK);
287 
288  typedef itk::OrImageFilter<OutputImageType, OutputUCharImageType, OutputUCharImageType> OrImageFilterType;
289  OrImageFilterType::Pointer orFilter = OrImageFilterType::New();
290 
291  orFilter->SetInput1(m_ThresholdFilter->GetOutput());
292  orFilter->SetInput2(workingImageSliceInITK);
293  orFilter->Update();
294 
295  mitk::CastToMitkImage(orFilter->GetOutput(), segmentationResult);
296  }
297  else
298  {
299  OutputImageType::Pointer workingImageSliceInITK = OutputImageType::New();
300 
301  CastToItkImage(workingImageSlice, workingImageSliceInITK);
302 
303  typedef itk::OrImageFilter<OutputImageType, OutputImageType> OrImageFilterType;
304  OrImageFilterType::Pointer orFilter = OrImageFilterType::New();
305 
306  orFilter->SetInput(0, m_ThresholdFilter->GetOutput());
307  orFilter->SetInput(1, workingImageSliceInITK);
308  orFilter->Update();
309 
310  mitk::CastToMitkImage(orFilter->GetOutput(), segmentationResult);
311  }
312 
313  segmentationResult->GetGeometry()->SetOrigin(workingImageSlice->GetGeometry()->GetOrigin());
314  segmentationResult->GetGeometry()->SetIndexToWorldTransform(
315  workingImageSlice->GetGeometry()->GetIndexToWorldTransform());
316 
317  // write to segmentation volume and hide preview image
318  // again, current time step is not considered
319  this->WriteBackSegmentationResult(m_WorkingPlane, segmentationResult, m_CurrentTimeStep);
320  this->m_ResultImageNode->SetVisibility(false);
321 
322  this->ClearSeeds();
323  }
324 
327 }
328 
330 {
331  // Add a new seed point for FastMarching algorithm
332  auto *positionEvent = dynamic_cast<mitk::InteractionPositionEvent *>(interactionEvent);
333  // const PositionEvent* p = dynamic_cast<const PositionEvent*>(stateEvent->GetEvent());
334  if (positionEvent == nullptr)
335  return;
336 
337  if (m_PositionEvent.IsNotNull())
338  m_PositionEvent = nullptr;
339 
341  InteractionPositionEvent::New(positionEvent->GetSender(), positionEvent->GetPointerPositionOnScreen());
342 
343  // if click was on another renderwindow or slice then reset pipeline and preview
344  if ((m_LastEventSender != m_PositionEvent->GetSender()) ||
345  (m_LastEventSlice != m_PositionEvent->GetSender()->GetSlice()))
346  {
347  this->BuildITKPipeline();
348  this->ClearSeeds();
349  }
350 
351  m_LastEventSender = m_PositionEvent->GetSender();
353  m_WorkingPlane = positionEvent->GetSender()->GetCurrentWorldPlaneGeometry()->Clone();
354 
355  mitk::Point3D clickInIndex;
356 
357  m_ReferenceImageSlice->GetGeometry()->WorldToIndex(m_PositionEvent->GetPositionInWorld(), clickInIndex);
358  itk::Index<2> seedPosition;
359  seedPosition[0] = clickInIndex[0];
360  seedPosition[1] = clickInIndex[1];
361 
362  NodeType node;
363  const double seedValue = 0.0;
364  node.SetValue(seedValue);
365  node.SetIndex(seedPosition);
366  this->m_SeedContainer->InsertElement(this->m_SeedContainer->Size(), node);
367  m_FastMarchingFilter->Modified();
368 
369  m_SeedsAsPointSet->InsertPoint(m_SeedsAsPointSet->GetSize(), m_PositionEvent->GetPositionInWorld());
370 
372 
373  m_NeedUpdate = true;
374 
375  this->Update();
376 
377  m_ReadyMessage.Send();
378 }
379 
381 {
382  // delete last seed point
383  if (!(this->m_SeedContainer->empty()))
384  {
385  // delete last element of seeds container
386  this->m_SeedContainer->pop_back();
387  m_FastMarchingFilter->Modified();
388 
389  // delete last point in pointset - somehow ugly
390  m_SeedsAsPointSet->GetPointSet()->GetPoints()->DeleteIndex(m_SeedsAsPointSet->GetSize() - 1);
391 
393 
394  m_NeedUpdate = true;
395 
396  this->Update();
397  }
398 }
399 
401 {
402  const unsigned int progress_steps = 20;
403 
404  // update FastMarching pipeline and show result
405  if (m_NeedUpdate)
406  {
407  m_ProgressCommand->AddStepsToDo(progress_steps);
408  CurrentlyBusy.Send(true);
409  try
410  {
411  m_ThresholdFilter->Update();
412  }
413  catch (itk::ExceptionObject &excep)
414  {
415  MITK_ERROR << "Exception caught: " << excep.GetDescription();
416 
417  // progress by max step count, will force
418  m_ProgressCommand->SetProgress(progress_steps);
419  CurrentlyBusy.Send(false);
420 
421  std::string msg = excep.GetDescription();
422  ErrorMessage.Send(msg);
423 
424  return;
425  }
426  m_ProgressCommand->SetProgress(progress_steps);
427  CurrentlyBusy.Send(false);
428 
429  // make output visible
431  CastToMitkImage(m_ThresholdFilter->GetOutput(), result);
432  result->GetGeometry()->SetOrigin(m_ReferenceImageSlice->GetGeometry()->GetOrigin());
433  result->GetGeometry()->SetIndexToWorldTransform(m_ReferenceImageSlice->GetGeometry()->GetIndexToWorldTransform());
434  m_ResultImageNode->SetData(result);
435  m_ResultImageNode->SetVisibility(true);
437  }
438 }
439 
441 {
442  // clear seeds for FastMarching as well as the PointSet for visualization
443  if (this->m_SeedContainer.IsNotNull())
444  this->m_SeedContainer->Initialize();
445 
446  if (this->m_SeedsAsPointSet.IsNotNull())
447  {
449  this->m_SeedsAsPointSetNode->SetData(this->m_SeedsAsPointSet);
450  m_SeedsAsPointSetNode->SetName("Seeds_Preview");
451  m_SeedsAsPointSetNode->SetBoolProperty("helper object", true);
452  m_SeedsAsPointSetNode->SetColor(0.0, 1.0, 0.0);
453  m_SeedsAsPointSetNode->SetVisibility(true);
454  }
455 
456  if (this->m_FastMarchingFilter.IsNotNull())
457  m_FastMarchingFilter->Modified();
458 
459  this->m_NeedUpdate = true;
460 }
461 
463 {
464  // clear all seeds and preview empty result
465  this->ClearSeeds();
466 
467  m_ResultImageNode->SetVisibility(false);
468 
470 }
471 
473 {
474  if (m_CurrentTimeStep != t)
475  {
476  m_CurrentTimeStep = t;
477 
478  this->Initialize();
479  }
480 }
unsigned int m_LastEventSlice
us::ModuleResource GetCursorIconResource() const override
Returns the path of a cursor icon.
Super class for all position events.
mitk::InteractionPositionEvent::Pointer m_PositionEvent
Image::Pointer GetAffectedImageSliceAs2DImage(const InteractionPositionEvent *positionEvent, const Image *image, unsigned int component=0)
Extract the slice of an image that the user just scribbles on. The given component denotes the vector...
virtual unsigned int GetSlice() const
DataStorage * GetDataStorage()
Message1< bool > CurrentlyBusy
To send whether the tool is busy (to be shown by some GUI)
Definition: mitkTool.h:104
#define MITK_ERROR
Definition: mitkLogMacros.h:20
SigmoidFilterType::Pointer m_SigmoidFilter
static Pointer New()
virtual void OnDelete(StateMachineAction *, InteractionEvent *interactionEvent)
Delete action of StateMachine pattern.
mitk::ToolCommand::Pointer m_ProgressCommand
void Update()
Updates the itk pipeline and shows the result of FastMarching.
void Send(T t)
Definition: mitkMessage.h:602
#define MITKSEGMENTATION_EXPORT
DataCollection - Class to facilitate loading/accessing structured data.
ThresholdingFilterType::Pointer m_ThresholdFilter
Message1< std::string > ErrorMessage
To send error messages (to be shown by some GUI)
Definition: mitkTool.h:99
Constants for most interaction classes, due to the generic StateMachines.
BaseRenderer * m_LastEventSender
Image::Pointer GetAffectedReferenceSlice(const InteractionPositionEvent *)
Extract the slice of the currently selected reference image that the user just scribbles on...
virtual void Add(DataNode *node, const DataStorage::SetOfObjects *parents=nullptr)=0
Adds a DataNode containing a data object to its internal storage.
static Pointer New(BaseRenderer *_arga, const Point2D &_argb)
bool ActivateTool(int id)
virtual void SetCurrentTimeStep(int t)
Set the working time step.
static Pointer New()
mitk::PointSet::Pointer m_SeedsAsPointSet
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.
MITK_TOOL_MACRO(MITKSEGMENTATION_EXPORT, LiveWireTool2D, "LiveWire tool")
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...
const char * GetName() const override
Returns the name of this tool. Make it short!
Module * GetModule() const
Image class for storing images.
Definition: mitkImage.h:72
ModuleResource GetResource(const std::string &path) const
Definition: usModule.cpp:267
void WriteBackSegmentationResult(const InteractionPositionEvent *, Image *)
void ConnectActionsAndFunctions() override
InternalImageType::Pointer m_ReferenceImageSliceAsITK
SmoothingFilterType::Pointer m_SmoothFilter
static Pointer New()
Image::Pointer m_ReferenceImageSlice
void SetSigma(double)
Set parameter used in Gradient Magnitude filter.
virtual void Remove(const DataNode *node)=0
Removes node from the DataStorage.
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.
mitk::DataNode::Pointer m_ResultImageNode
FastMarchingFilterType::NodeType NodeType
GradientFilterType::Pointer m_GradientMagnitudeFilter
mitk::DataNode::Pointer m_SeedsAsPointSetNode
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:74
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.
FastMarchingFilterType::Pointer m_FastMarchingFilter
DataVectorType GetReferenceData()
ToolManager * m_ToolManager
Definition: mitkTool.h:228
NodeContainer::Pointer m_SeedContainer
#define CONNECT_FUNCTION(a, f)
void SetStoppingValue(double)
Set parameter used in Fast Marching filter.
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()
DataVectorType GetWorkingData()
void Activated() override
Called when the tool gets activated.
const char ** GetXPM() const override
Returns an icon in the XPM format.
void Deactivated() override
Called when the tool gets deactivated.