declassified ITKwatershed

requirements.txt

-
 # numpy==1.20.1 # Python >=3.7
 numpy==1.19.5 # Python >=3.6
 
@@ -33,8 +32,8 @@ pyamg==4.0.0
 # mpi4py==3.0.3
 
 # Jupyter
-ipykernel==5.5.0 # Python >=3.5
-# ipython==7.21.0 # Python >=3.7
+#ipykernel==5.5.0 # Python >=3.5
+#ipython==7.21.0 # Python >=3.7
 ipython==7.16.1 # Python >=3.6
 
 # R (correlated random field)
@@ -43,9 +42,9 @@ rpy2==3.4.2
 #############
 # Graphical #
 #############
-PyQt5==5.15.4 # Python >=3.6
-PyQt5-sip==12.8.1 # Python >=3.5
-qimage2ndarray==1.8.3
+#PyQt5==5.15.4 # Python >=3.6
+#PyQt5-sip==12.8.1 # Python >=3.5
+#qimage2ndarray==1.8.3
 
 
 ##############

tools/label/ITKwatershed.py

@@ -16,132 +16,10 @@ You should have received a copy of the GNU General Public License along with
 this program.  If not, see <http://www.gnu.org/licenses/>.
 """
 
-from __future__ import print_function
-import numpy
-import tifffile
-import sys
-import os
 import SimpleITK
 import spam.label
 
-
-#if sys.version_info >= (3, 0):
-    #def xrange(*args, **kwargs):
-        #return iter(range(*args, **kwargs))
-
-"""
-Classes binding ITK morphological watershed algorithm implementation
-
-Note:
-Author is Olumide Okubadejo
-"""
-
-class _MWatershed(object):
-    def __init__(self, data, level=1, watershedLineOn=False, fullyConnected=True, fromMarkers=False, markerImage=None):
-        self._level = level
-        self._watershedLineOn = watershedLineOn
-        self._fullyConnected = fullyConnected
-        self._data = SimpleITK.GetImageFromArray(data)
-        self._fromMarkers = fromMarkers
-        if fromMarkers:
-            self._markerImage = SimpleITK.GetImageFromArray(markerImage)
-
-    #def _CalculateGradient(self):
-        #gradient = SimpleITK.GradientMagnitudeImageFilter()
-        #self._gradient = gradient.Execute(self._bdata)
-
-    def _Binarize(self):
-        threshold = SimpleITK.OtsuThresholdImageFilter()
-        self._bdata = threshold.Execute(self._data)
-        self._bdata = self._Rescale(self._bdata)
-        self._threshold = threshold.GetThreshold()
-
-    def _Rescale(self, data):
-        rescale = SimpleITK.RescaleIntensityImageFilter()
-        return rescale.Execute(data, 0, 65535)
-
-    def _DistanceMap(self):
-        distanceMap = SimpleITK.DanielssonDistanceMapImageFilter()
-        self._distance = distanceMap.Execute(self._bdata)
-        pass
-
-    def _InvertImage(self, data):
-        invert = SimpleITK.InvertIntensityImageFilter()
-        inverted = invert.Execute(data)
-        return inverted
-
-    def _MaskImage(self):
-        mask = SimpleITK.MaskImageFilter()
-        self._mask = mask.Execute(self._labelImage, self._bdata)
-
-    def _FillHoles(self):
-        fill = SimpleITK.BinaryFillholeImageFilter()
-        self._bdata = fill.Execute(self._bdata)
-
-    def _LabelOverlay(self):
-        overlay = SimpleITK.LabelOverlayImageFilter()
-        self._overlay = overlay.Execute(self._bdata, self._mask)
-
-    def _Watershed(self):
-        if (self._fromMarkers):
-            watershed = SimpleITK.MorphologicalWatershedFromMarkersImageFilter()
-        else:
-            watershed = SimpleITK.MorphologicalWatershedImageFilter()
-        watershed.SetFullyConnected(self._fullyConnected)
-        watershed.SetMarkWatershedLine(self._watershedLineOn)
-
-        if (self._fromMarkers):
-            self._labelImage = watershed.Execute(
-                self._distance, self._markerImage)
-        else:
-            watershed.SetLevel(self._level)
-            self._labelImage = watershed.Execute(self._distance)
-
-    #def GetThreshold(self):
-        #return self._threshold
-
-    #def GetLabeledImage(self):
-        #return SimpleITK.GetArrayFromImage(self._mask)
-
-    #def GetITKLabeledImage(self):
-        #return self._labelImage
-
-    #def GetColorMap(self):
-        #return SimpleITK.GetArrayFromImage(self._overlay)
-
-    #def GetBinaryImage(self):
-        #return SimpleITK.GetArrayFromImage(self._bdata)
-
-    #def GetDistanceMapImage(self):
-        #return SimpleITK.GetArrayFromImage(self._distance)
-
-    def GetMaskImage(self):
-        return SimpleITK.GetArrayFromImage(self._mask)
-
-    #def GetLabelParams(self):
-        #return self._labelParams.GetLabels()
-
-    #def GetLabelObject(self):
-        #return self._labelParams
-
-    #def GetRadiusArray(self):
-        #return self._labelParams.GetRadiusArray()
-
-    def run(self):
-        self._Binarize()
-        # self._bdata = self._InvertImage(self._bdata)
-        self._FillHoles()
-        self._DistanceMap()
-        self._distance = self._InvertImage(self._distance)
-        self._Watershed()
-        self._bdata = self._InvertImage(self._bdata)
-        self._MaskImage()
-        self._LabelOverlay()
-        #self._labelParams = _Labels(self._data, self._labelImage)
-        return SimpleITK.GetArrayFromImage(self._overlay)
-
-
-def watershed(binary, markers=None, verbose = False):
+def watershed(binary, markers=None, watershedLevel=1):
     """
     This function runs an ITK watershed on a binary image and returns a labelled image.
     This function uses an interpixel watershed.
@@ -154,40 +32,77 @@ def watershed(binary, markers=None, verbose = False):
         markers : 3D numpy array (optional, default = None)
             Not implemented yet, but try!
 
-	verbose : boolean (optional, default = False)
-	        True for printing the evolution of process
-	        False for not printing the evolution of process
+        watershedLevel : int, optional
+            Watershed level for merging maxima
 
     Returns
     --------
         labelled : 3D numpy array of ints
             3D array where each object is numbered
     """
+    import matplotlib.pyplot as plt
+    def plotme(im):
+        im = SimpleITK.GetArrayFromImage(im)
+        plt.imshow(im[im.shape[0]//2])
+        plt.show()
+
     binary = binary > 0
 
     # Let's convert it 8-bit
     binary = binary.astype('<u1') * 255
 
+    bdata = SimpleITK.GetImageFromArray(binary)
+
     if markers is not None:
-        if verbose:
-            print("\tITKwatershed.watershed(): Running watershed with your markers...", end='')
-        if markers.max() > 65535: markers = markers.astype('<u4')
-        else:                     markers = markers.astype('<u2')
+        markers = markers.astype('<u4') if markers.max() > 65535 else markers.astype('<u2')
+        markers = SimpleITK.GetImageFromArray(markers)
+
+    watershedLevel=1
+    watershedLineOn=False
+    fullyConnected=True
+
+    #thresholdFilt = SimpleITK.OtsuThresholdImageFilter()
+    #bdata = thresholdFilt.Execute(data)
+
+    #rescaleFilt = SimpleITK.RescaleIntensityImageFilter()
+    #rescaleFilt.SetOutputMinimum(0)
+    #rescaleFilt.SetOutputMaximum(65535)
+    #bdata = rescaleFilt.Execute(data)
 
-        mWatershed = _MWatershed(binary, markerImage=markers.astype(spam.label.labelType), fromMarkers=True,
-                                 level=1, watershedLineOn=False, fullyConnected=True)
+    #threshold = thresholdFilt.GetThreshold()
+
+    #fillFilt = SimpleITK.BinaryFillholeImageFilter()
+    #bdata = fillFilt.Execute(bdata)
+
+    invertFilt = SimpleITK.InvertIntensityImageFilter()
+    bdata = invertFilt.Execute(bdata)
+
+    distanceMapFilt = SimpleITK.DanielssonDistanceMapImageFilter()
+    distance = distanceMapFilt.Execute(bdata)
+
+    distance = invertFilt.Execute(distance)
+
+    if markers is None:
+        watershedFilt = SimpleITK.MorphologicalWatershedImageFilter()
     else:
-        if verbose:
-            print("\tITKwatershed.watershed(): Running watershed...", end='')
-        mWatershed = _MWatershed(
-            binary, level=1, watershedLineOn=False, fullyConnected=True)
-    labels = mWatershed.run()
-    if verbose:
-        print("done.")
-    if verbose:
-        print("\tITKwatershed.watershed(): Collecting labelled image...", end='')
-    lab = mWatershed.GetMaskImage().astype('<u4')
-    if verbose:
-        print("done.")
+        watershedFilt = SimpleITK.MorphologicalWatershedFromMarkersImageFilter()
+    watershedFilt.SetFullyConnected(fullyConnected)
+    watershedFilt.SetMarkWatershedLine(watershedLineOn)
+
+    if markers is None:
+        watershedFilt.SetLevel(watershedLevel)
+        labelImage = watershedFilt.Execute(distance)
+    else:
+        labelImage = watershedFilt.Execute(distance, markers)
+
+    bdata = invertFilt.Execute(bdata)
+
+    maskFilt = SimpleITK.MaskImageFilter()
+    mask = maskFilt.Execute(labelImage, bdata)
+
+    #overlayFilt = SimpleITK.LabelOverlayImageFilter()
+    #overlay = overlayFilt.Execute(bdata, mask)
+
+    lab = SimpleITK.GetArrayFromImage(mask).astype(spam.label.labelType)
 
     return lab

tools/tests/test_label.py

@@ -601,7 +601,7 @@ class TestFunctionLabel(unittest.TestCase):
         box = numpy.zeros((100, 100, 100), dtype='<f8')
         spam.kalisphera.makeSphere(box, [[50, 50, 50-20], [50, 50, 50+20]], [20, 20])
         box = box > 0.5
-        lab = spam.label.watershed(box, verbose=True)
+        lab = spam.label.watershed(box)
         labVolumes = spam.label.volumes(lab)
 
         setVoronoi = spam.label.setVoronoi(lab)