Slicer  4.8
Slicer is a multi-platform, free and open source software package for visualization and medical image computing
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SegmentEditorSmoothingEffect.py
Go to the documentation of this file.
1 import os
2 import vtk, qt, ctk, slicer
3 import logging
4 from SegmentEditorEffects import *
5 
7  """ SmoothingEffect is an Effect that smoothes a selected segment
8  """
9 
10  def __init__(self, scriptedEffect):
11  scriptedEffect.name = 'Smoothing'
12  AbstractScriptedSegmentEditorEffect.__init__(self, scriptedEffect)
13 
14  def clone(self):
15  import qSlicerSegmentationsEditorEffectsPythonQt as effects
16  clonedEffect = effects.qSlicerSegmentEditorScriptedEffect(None)
17  clonedEffect.setPythonSource(__file__.replace('\\','/'))
18  return clonedEffect
19 
20  def icon(self):
21  iconPath = os.path.join(os.path.dirname(__file__), 'Resources/Icons/Smoothing.png')
22  if os.path.exists(iconPath):
23  return qt.QIcon(iconPath)
24  return qt.QIcon()
25 
26  def helpText(self):
27  return """<html>Make segment boundaries smoother<br> by removing extrusions and filling small holes. Available methods:<p>
28 <ul style="margin: 0">
29 <li><b>Median:</b> removes small details while keeps smooth contours mostly unchanged. Applied to selected segment only.</li>
30 <li><b>Opening:</b> removes extrusions smaller than the specified kernel size. Applied to selected segment only.</li>
31 <li><b>Closing:</b> fills sharp corners and holes smaller than the specified kernel size. Applied to selected segment only.</li>
32 <li><b>Gaussian:</b> smoothes all contours, tends to shrink the segment. Applied to selected segment only.</li>
33 <li><b>Joint smoothing:</b> smoothes multiple segments at once, preserving watertight interface between them. Masking settings are bypassed.
34 If segments overlap, segment higher in the segments table will have priority. <b>Applied to all visible segments.</b></li>
35 </ul><p></html>"""
36 
37  def setupOptionsFrame(self):
38 
39  self.methodSelectorComboBox = qt.QComboBox()
40  self.methodSelectorComboBox.addItem("Median", MEDIAN)
41  self.methodSelectorComboBox.addItem("Opening (remove extrusions)", MORPHOLOGICAL_OPENING)
42  self.methodSelectorComboBox.addItem("Closing (fill holes)", MORPHOLOGICAL_CLOSING)
43  self.methodSelectorComboBox.addItem("Gaussian", GAUSSIAN)
44  self.methodSelectorComboBox.addItem("Joint smoothing", JOINT_TAUBIN)
45  self.scriptedEffect.addLabeledOptionsWidget("Smoothing method:", self.methodSelectorComboBox)
46 
47  self.kernelSizeMmSpinBox = slicer.qMRMLSpinBox()
48  self.kernelSizeMmSpinBox.setMRMLScene(slicer.mrmlScene)
49  self.kernelSizeMmSpinBox.setToolTip("Diameter of the neighborhood that will be considered around each voxel. Higher value makes smoothing stronger (more details are suppressed).")
50  self.kernelSizeMmSpinBox.quantity = "length"
51  self.kernelSizeMmSpinBox.unitAwareProperties &= ~slicer.qMRMLSpinBox.MinimumValue # disable setting deafult minimum value (it would be a large negative value)
52  self.kernelSizeMmSpinBox.minimum = 0.0
53  self.kernelSizeMmSpinBox.value = 3.0
54  self.kernelSizeMmSpinBox.singleStep = 1.0
55 
56  self.kernelSizePixel = qt.QLabel()
57  self.kernelSizePixel.setToolTip("Diameter of the neighborhood in pixels. Computed from the segment's spacing and the specified kernel size.")
58 
59  kernelSizeFrame = qt.QHBoxLayout()
60  kernelSizeFrame.addWidget(self.kernelSizeMmSpinBox)
61  kernelSizeFrame.addWidget(self.kernelSizePixel)
62  self.kernelSizeMmLabel = self.scriptedEffect.addLabeledOptionsWidget("Kernel size:", kernelSizeFrame)
63 
64  self.gaussianStandardDeviationMmSpinBox = slicer.qMRMLSpinBox()
65  self.gaussianStandardDeviationMmSpinBox.setMRMLScene(slicer.mrmlScene)
66  self.gaussianStandardDeviationMmSpinBox.setToolTip("Standard deviation of the Gaussian smoothing filter coefficients. Higher value makes smoothing stronger (more details are suppressed).")
67  self.gaussianStandardDeviationMmSpinBox.quantity = "length"
68  self.gaussianStandardDeviationMmSpinBox.value = 3.0
69  self.gaussianStandardDeviationMmSpinBox.singleStep = 1.0
70  self.gaussianStandardDeviationMmLabel = self.scriptedEffect.addLabeledOptionsWidget("Standard deviation:", self.gaussianStandardDeviationMmSpinBox)
71 
72  self.jointTaubinSmoothingFactorSlider = ctk.ctkSliderWidget()
73  self.jointTaubinSmoothingFactorSlider.setToolTip("Higher value means stronger smoothing.")
74  self.jointTaubinSmoothingFactorSlider.minimum = 0.01
75  self.jointTaubinSmoothingFactorSlider.maximum = 1.0
76  self.jointTaubinSmoothingFactorSlider.value = 0.5
77  self.jointTaubinSmoothingFactorSlider.singleStep = 0.01
78  self.jointTaubinSmoothingFactorSlider.pageStep = 0.1
79  self.jointTaubinSmoothingFactorLabel = self.scriptedEffect.addLabeledOptionsWidget("Smoothing factor:", self.jointTaubinSmoothingFactorSlider)
80 
81  self.applyButton = qt.QPushButton("Apply")
82  self.applyButton.objectName = self.__class__.__name__ + 'Apply'
83  self.applyButton.setToolTip("Apply smoothing to selected segment")
84  self.scriptedEffect.addOptionsWidget(self.applyButton)
85 
86  self.methodSelectorComboBox.connect("currentIndexChanged(int)", self.updateMRMLFromGUI)
87  self.kernelSizeMmSpinBox.connect("valueChanged(double)", self.updateMRMLFromGUI)
88  self.gaussianStandardDeviationMmSpinBox.connect("valueChanged(double)", self.updateMRMLFromGUI)
89  self.jointTaubinSmoothingFactorSlider.connect("valueChanged(double)", self.updateMRMLFromGUI)
90  self.applyButton.connect('clicked()', self.onApply)
91 
92  def createCursor(self, widget):
93  # Turn off effect-specific cursor for this effect
94  return slicer.util.mainWindow().cursor
95 
96  def setMRMLDefaults(self):
97  self.scriptedEffect.setParameterDefault("SmoothingMethod", MEDIAN)
98  self.scriptedEffect.setParameterDefault("KernelSizeMm", 3)
99  self.scriptedEffect.setParameterDefault("GaussianStandardDeviationMm", 3)
100  self.scriptedEffect.setParameterDefault("JointTaubinSmoothingFactor", 0.5)
101 
103  methodIndex = self.methodSelectorComboBox.currentIndex
104  smoothingMethod = self.methodSelectorComboBox.itemData(methodIndex)
105  morphologicalMethod = (smoothingMethod==MEDIAN or smoothingMethod==MORPHOLOGICAL_OPENING or smoothingMethod==MORPHOLOGICAL_CLOSING)
106  self.kernelSizeMmLabel.setVisible(morphologicalMethod)
107  self.kernelSizeMmSpinBox.setVisible(morphologicalMethod)
108  self.kernelSizePixel.setVisible(morphologicalMethod)
109  self.gaussianStandardDeviationMmLabel.setVisible(smoothingMethod==GAUSSIAN)
110  self.gaussianStandardDeviationMmSpinBox.setVisible(smoothingMethod==GAUSSIAN)
111  self.jointTaubinSmoothingFactorLabel.setVisible(smoothingMethod==JOINT_TAUBIN)
112  self.jointTaubinSmoothingFactorSlider.setVisible(smoothingMethod==JOINT_TAUBIN)
113 
115  selectedSegmentLabelmapSpacing = [1.0, 1.0, 1.0]
116  selectedSegmentLabelmap = self.scriptedEffect.selectedSegmentLabelmap()
117  if selectedSegmentLabelmap:
118  selectedSegmentLabelmapSpacing = selectedSegmentLabelmap.GetSpacing()
119 
120  # size rounded to nearest odd number. If kernel size is even then image gets shifted.
121  kernelSizeMm = self.scriptedEffect.doubleParameter("KernelSizeMm")
122  kernelSizePixel = [int(round((kernelSizeMm / selectedSegmentLabelmapSpacing[componentIndex]+1)/2)*2-1) for componentIndex in range(3)]
123  return kernelSizePixel
124 
125  def updateGUIFromMRML(self):
126  methodIndex = self.methodSelectorComboBox.findData(self.scriptedEffect.parameter("SmoothingMethod"))
127  wasBlocked = self.methodSelectorComboBox.blockSignals(True)
128  self.methodSelectorComboBox.setCurrentIndex(methodIndex)
129  self.methodSelectorComboBox.blockSignals(wasBlocked)
130 
131  wasBlocked = self.kernelSizeMmSpinBox.blockSignals(True)
132  self.kernelSizeMmSpinBox.value = self.scriptedEffect.doubleParameter("KernelSizeMm")
133  self.kernelSizeMmSpinBox.blockSignals(wasBlocked)
134  kernelSizePixel = self.getKernelSizePixel()
135  self.kernelSizePixel.text = "{0}x{1}x{2} pixels".format(kernelSizePixel[0], kernelSizePixel[1], kernelSizePixel[2])
136 
137  wasBlocked = self.gaussianStandardDeviationMmSpinBox.blockSignals(True)
138  self.gaussianStandardDeviationMmSpinBox.value = self.scriptedEffect.doubleParameter("GaussianStandardDeviationMm")
139  self.gaussianStandardDeviationMmSpinBox.blockSignals(wasBlocked)
140 
141  wasBlocked = self.jointTaubinSmoothingFactorSlider.blockSignals(True)
142  self.jointTaubinSmoothingFactorSlider.value = self.scriptedEffect.doubleParameter("JointTaubinSmoothingFactor")
143  self.jointTaubinSmoothingFactorSlider.blockSignals(wasBlocked)
144 
146 
147  selectedSegmentLabelmap = self.scriptedEffect.selectedSegmentLabelmap()
148  if selectedSegmentLabelmap:
149  import math
150  selectedSegmentLabelmapSpacing = selectedSegmentLabelmap.GetSpacing()
151  singleStep = min(selectedSegmentLabelmapSpacing)
152  # round to power of 10 (for example: 0.2 -> 0.1; 0.09 -> 0.01) to show "nice" values
153  singleStep = pow(10,math.floor(math.log(singleStep)/math.log(10)))
154  self.kernelSizeMmSpinBox.singleStep = singleStep
155  self.gaussianStandardDeviationMmSpinBox.singleStep = singleStep
156 
157  def updateMRMLFromGUI(self):
158  methodIndex = self.methodSelectorComboBox.currentIndex
159  smoothingMethod = self.methodSelectorComboBox.itemData(methodIndex)
160  self.scriptedEffect.setParameter("SmoothingMethod", smoothingMethod)
161  self.scriptedEffect.setParameter("KernelSizeMm", self.kernelSizeMmSpinBox.value)
162  self.scriptedEffect.setParameter("GaussianStandardDeviationMm", self.gaussianStandardDeviationMmSpinBox.value)
163  self.scriptedEffect.setParameter("JointTaubinSmoothingFactor", self.jointTaubinSmoothingFactorSlider.value)
164 
166 
167  #
168  # Effect specific methods (the above ones are the API methods to override)
169  #
170 
171  def onApply(self):
172  try:
173  # This can be a long operation - indicate it to the user
174  qt.QApplication.setOverrideCursor(qt.Qt.WaitCursor)
175 
176  self.scriptedEffect.saveStateForUndo()
177 
178  smoothingMethod = self.scriptedEffect.parameter("SmoothingMethod")
179  if smoothingMethod == JOINT_TAUBIN:
181  else:
182  self.smoothSelectedSegment()
183  finally:
184  qt.QApplication.restoreOverrideCursor()
185 
187  try:
188 
189  # Get master volume image data
190  import vtkSegmentationCorePython
191 
192  # Get modifier labelmap
193  modifierLabelmap = self.scriptedEffect.defaultModifierLabelmap()
194  selectedSegmentLabelmap = self.scriptedEffect.selectedSegmentLabelmap()
195 
196  smoothingMethod = self.scriptedEffect.parameter("SmoothingMethod")
197 
198  if smoothingMethod == GAUSSIAN:
199  maxValue = 255
200 
201  thresh = vtk.vtkImageThreshold()
202  thresh.SetInputData(selectedSegmentLabelmap)
203  thresh.ThresholdByLower(0)
204  thresh.SetInValue(0)
205  thresh.SetOutValue(maxValue)
206  thresh.SetOutputScalarType(vtk.VTK_UNSIGNED_CHAR)
207 
208  standardDeviationMm = self.scriptedEffect.doubleParameter("GaussianStandardDeviationMm")
209  gaussianFilter = vtk.vtkImageGaussianSmooth()
210  gaussianFilter.SetInputConnection(thresh.GetOutputPort())
211  gaussianFilter.SetStandardDeviation(standardDeviationMm)
212  gaussianFilter.SetRadiusFactor(4)
213 
214  thresh2 = vtk.vtkImageThreshold()
215  thresh2.SetInputConnection(gaussianFilter.GetOutputPort())
216  thresh2.ThresholdByUpper(maxValue/2)
217  thresh2.SetInValue(1)
218  thresh2.SetOutValue(0)
219  thresh2.SetOutputScalarType(selectedSegmentLabelmap.GetScalarType())
220  thresh2.Update()
221  modifierLabelmap.DeepCopy(thresh2.GetOutput())
222 
223  else:
224  # size rounded to nearest odd number. If kernel size is even then image gets shifted.
225  kernelSizePixel = self.getKernelSizePixel()
226 
227  if smoothingMethod == MEDIAN:
228  # Median filter does not require a particular label value
229  smoothingFilter = vtk.vtkImageMedian3D()
230  smoothingFilter.SetInputData(selectedSegmentLabelmap)
231 
232  else:
233  # We need to know exactly the value of the segment voxels, apply threshold to make force the selected label value
234  labelValue = 1
235  backgroundValue = 0
236  thresh = vtk.vtkImageThreshold()
237  thresh.SetInputData(selectedSegmentLabelmap)
238  thresh.ThresholdByLower(0)
239  thresh.SetInValue(backgroundValue)
240  thresh.SetOutValue(labelValue)
241  thresh.SetOutputScalarType(selectedSegmentLabelmap.GetScalarType())
242 
243  smoothingFilter = vtk.vtkImageOpenClose3D()
244  smoothingFilter.SetInputConnection(thresh.GetOutputPort())
245  if smoothingMethod == MORPHOLOGICAL_OPENING:
246  smoothingFilter.SetOpenValue(labelValue)
247  smoothingFilter.SetCloseValue(backgroundValue)
248  else: # must be smoothingMethod == MORPHOLOGICAL_CLOSING:
249  smoothingFilter.SetOpenValue(backgroundValue)
250  smoothingFilter.SetCloseValue(labelValue)
251 
252  smoothingFilter.SetKernelSize(kernelSizePixel[0],kernelSizePixel[1],kernelSizePixel[2])
253  smoothingFilter.Update()
254  modifierLabelmap.DeepCopy(smoothingFilter.GetOutput())
255 
256  except IndexError:
257  logging.error('apply: Failed to apply smoothing')
258 
259  # Apply changes
260  self.scriptedEffect.modifySelectedSegmentByLabelmap(modifierLabelmap, slicer.qSlicerSegmentEditorAbstractEffect.ModificationModeSet)
261 
263  import vtkSegmentationCorePython as vtkSegmentationCore
264 
265  # Generate merged labelmap of all visible segments
266  segmentationNode = self.scriptedEffect.parameterSetNode().GetSegmentationNode()
267  visibleSegmentIds = vtk.vtkStringArray()
268  segmentationNode.GetDisplayNode().GetVisibleSegmentIDs(visibleSegmentIds)
269  if visibleSegmentIds.GetNumberOfValues() == 0:
270  logging.info("Smoothing operation skipped: there are no visible segments")
271  return
272 
273  mergedImage = vtkSegmentationCore.vtkOrientedImageData()
274  if not segmentationNode.GenerateMergedLabelmapForAllSegments(mergedImage,
275  vtkSegmentationCore.vtkSegmentation.EXTENT_UNION_OF_SEGMENTS_PADDED,
276  None, visibleSegmentIds):
277  logging.error('Failed to apply smoothing: cannot get list of visible segments')
278  return
279 
280  segmentLabelValues = [] # list of [segmentId, labelValue]
281  for i in range(visibleSegmentIds.GetNumberOfValues()):
282  segmentId = visibleSegmentIds.GetValue(i)
283  segmentLabelValues.append([segmentId, i+1])
284 
285  # Perform smoothing in voxel space
286  ici = vtk.vtkImageChangeInformation()
287  ici.SetInputData(mergedImage)
288  ici.SetOutputSpacing(1, 1, 1)
289  ici.SetOutputOrigin(0, 0, 0)
290 
291  # Convert labelmap to combined polydata
292  convertToPolyData = vtk.vtkDiscreteMarchingCubes()
293  convertToPolyData.SetInputConnection(ici.GetOutputPort())
294  convertToPolyData.SetNumberOfContours(len(segmentLabelValues))
295  contourIndex = 0
296  for segmentId, labelValue in segmentLabelValues:
297  convertToPolyData.SetValue(contourIndex, labelValue)
298  contourIndex += 1
299 
300  # Low-pass filtering using Taubin's method
301  smoothingFactor = self.scriptedEffect.doubleParameter("JointTaubinSmoothingFactor")
302  smoothingIterations = 100 # according to VTK documentation 10-20 iterations could be enough but we use a higher value to reduce chance of shrinking
303  passBand = pow(10.0, -4.0*smoothingFactor) # gives a nice range of 1-0.0001 from a user input of 0-1
304  smoother = vtk.vtkWindowedSincPolyDataFilter()
305  smoother.SetInputConnection(convertToPolyData.GetOutputPort())
306  smoother.SetNumberOfIterations(smoothingIterations)
307  smoother.BoundarySmoothingOff()
308  smoother.FeatureEdgeSmoothingOff()
309  smoother.SetFeatureAngle(90.0)
310  smoother.SetPassBand(passBand)
311  smoother.NonManifoldSmoothingOn()
312  smoother.NormalizeCoordinatesOn()
313 
314  # Extract a label
315  threshold = vtk.vtkThreshold()
316  threshold.SetInputConnection(smoother.GetOutputPort())
317 
318  # Convert to polydata
319  geometryFilter = vtk.vtkGeometryFilter()
320  geometryFilter.SetInputConnection(threshold.GetOutputPort())
321 
322  # Convert polydata to stencil
323  polyDataToImageStencil = vtk.vtkPolyDataToImageStencil()
324  polyDataToImageStencil.SetInputConnection(geometryFilter.GetOutputPort())
325  polyDataToImageStencil.SetOutputSpacing(1,1,1)
326  polyDataToImageStencil.SetOutputOrigin(0,0,0)
327  polyDataToImageStencil.SetOutputWholeExtent(mergedImage.GetExtent())
328 
329  # Convert stencil to image
330  stencil = vtk.vtkImageStencil()
331  emptyBinaryLabelMap = vtk.vtkImageData()
332  emptyBinaryLabelMap.SetExtent(mergedImage.GetExtent())
333  emptyBinaryLabelMap.AllocateScalars(vtk.VTK_UNSIGNED_CHAR, 1)
334  vtkSegmentationCore.vtkOrientedImageDataResample.FillImage(emptyBinaryLabelMap, 0)
335  stencil.SetInputData(emptyBinaryLabelMap)
336  stencil.SetStencilConnection(polyDataToImageStencil.GetOutputPort())
337  stencil.ReverseStencilOn()
338  stencil.SetBackgroundValue(1) # General foreground value is 1 (background value because of reverse stencil)
339 
340  imageToWorldMatrix = vtk.vtkMatrix4x4()
341  mergedImage.GetImageToWorldMatrix(imageToWorldMatrix)
342 
343  for segmentId, labelValue in segmentLabelValues:
344  threshold.ThresholdBetween(labelValue, labelValue)
345  stencil.Update()
346  smoothedBinaryLabelMap = vtkSegmentationCore.vtkOrientedImageData()
347  smoothedBinaryLabelMap.ShallowCopy(stencil.GetOutput())
348  smoothedBinaryLabelMap.SetImageToWorldMatrix(imageToWorldMatrix)
349  # Write results to segments directly, bypassing masking
350  slicer.vtkSlicerSegmentationsModuleLogic.SetBinaryLabelmapToSegment(smoothedBinaryLabelMap,
351  segmentationNode, segmentId, slicer.vtkSlicerSegmentationsModuleLogic.MODE_REPLACE, smoothedBinaryLabelMap.GetExtent())
352 
353 MEDIAN = 'MEDIAN'
354 GAUSSIAN = 'GAUSSIAN'
355 MORPHOLOGICAL_OPENING = 'MORPHOLOGICAL_OPENING'
356 MORPHOLOGICAL_CLOSING = 'MORPHOLOGICAL_CLOSING'
357 JOINT_TAUBIN = 'JOINT_TAUBIN'