watershed does not like flat images

[jonwright]

... we are just running through the tutorial in here: https://www.silx.org/doc/pyFAI/dev/usage/tutorial/Detector/CCD_Calibration/CCD_calibration.html

... and we noticed that watershed does not like a flat image and produces lots peaks so it takes a really really long time to run when we have masked borders on our data.

Simple fix might be to ignore all pixels less than threshold?

import numpy as np, pylab as pl

try: #depends if the version of pyFAI you are using
    from pyFAI.watershed import InverseWatershed
except:
    from pyFAI.ext.watershed import InverseWatershed
    #Version of pyFAI newer than feb 2016

g = 100 * np.exp( -np.linspace(-3,3)**2 )
pk = np.outer( g, g )

bad = np.where( pk < 1000, 1000, pk )

LN = pl.matplotlib.colors.LogNorm
f,a = pl.subplots(2,2)
f.colorbar( a[0,0].imshow(pk, norm=LN()), ax=a[0,0])
f.colorbar( a[0,1].imshow(bad, norm=LN()), ax=a[0,1])

iw = InverseWatershed(pk, thres=2000)
iw.init()
iw.merge_singleton()
all_regions = set(iw.regions.values())
mini = 10 
regions = [i for i in all_regions if i.size>mini]

f.colorbar( a[1,0].imshow(iw.labels%19, cmap='tab20'), ax=a[1,0])


print("Number of region segmented: %s"%len(all_regions))
print("Number of large enough regions : %s"%len(regions))


iw2 = InverseWatershed(bad, thres=2000)
iw2.init()
iw2.merge_singleton()
all_regions2 = set(iw2.regions.values())
mini = 10 
regions = [i for i in all_regions2 if i.size>mini]

f.colorbar( a[1,1].imshow(iw2.labels%19, cmap='tab20'), ax=a[1,1] )

print("Number of region segmented: %s"%len(all_regions))
print("Number of large enough regions : %s"%len(regions))

pl.savefig('flat_regions_peaky.png')

[kif]

For a flat region, after the init, every-single pixel is considered as a peak and defines a region of size one. In your example, there are 1920 regions found.

After the merge_singleton, all those single-pixel regions are all merged together but the algorithm is in O(2) speed, so fairly slow with large flat regions.

I have to admit the "labels" are not properly updated... this deserves fixing.

[jonwright]

Given that scipy is a dependency already, then scipy.ndimage.watershed_ift might be worth testing? There is also something in skimage.

In ImageD11 we have ImageD11.cImageD11.localmaxlabel for assigning pixels to maxima and labeling, and also a sparse version of the same code.