%run ../../common.ipynb

from skimage.morphology import watershed

>>> from skimage.morphology import watershed
>>> from skimage.feature import peak_local_max
>>> from skimage.measure import label
>>>
>>> # Generate an initial image with two overlapping circles
>>> x, y = np.indices((80, 80))
>>> x1, y1, x2, y2 = 28, 28, 44, 52
>>> r1, r2 = 16, 20
>>> mask_circle1 = (x - x1)**2 + (y - y1)**2 < r1**2
>>> mask_circle2 = (x - x2)**2 + (y - y2)**2 < r2**2
>>> image = np.logical_or(mask_circle1, mask_circle2)
>>> # Now we want to separate the two objects in image
>>> # Generate the markers as local maxima of the distance
>>> # to the background
>>> from scipy import ndimage
>>> distance = ndimage.distance_transform_edt(image)
>>> local_maxi = peak_local_max(distance, indices=False, footprint=np.ones((3, 3)),labels=image)
>>> markers = label(local_maxi)
>>> labels_ws = watershed(-distance, markers, mask=image)

gimshow(labels_ws)

image = imread('../mp.tif')
local_max = peak_local_max(image, threshold_rel=0.2, min_distance=10, indices=False)
markers = label(local_max)
image_ws = watershed(-image, markers)

imshow(image_ws)
image_ws.dtype

# scale it to 3 colors * 8 bit
image_ws_rgb = image_ws / image_ws.max() * 2**24
image_ws_rgb = image_ws_rgb.astype(np.uint32)
# reshape magic, tuple concat: (1,2) + (3,) = (1,2,3)
image_ws_rgb = image_ws_rgb.view(np.uint8).reshape(image.shape + (4,))
# lets get rid of unused color dimension (alpha)
image_ws_rgb = image_ws_rgb[:,:,0:3]
imsave('watershed.png', image_ws_rgb)
Image('watershed.png')

imsave('watershed_gray.png', scale(image_ws))
Image('watershed_gray.png')