Slicer  4.8
Slicer is a multi-platform, free and open source software package for visualization and medical image computing
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'