# Scikit-Image, Numpy, and Selecting Colors (python)

I'm trying to mask colored features from a photograph so that I can do some other processing on them. I've played with a few packages (scikit-image, mahotas, and openCV) and have settled on Scikit-image because it plays well with scikit-learn.

I'd like to accomplish two things, eventually: (1) get colored features so that I can compute things like length, and (2) remove colored features from the image while retaining the white spidering-looking veins (the lattice-looking stuff around the margin of the image will need to go, but that can be a later endeavour).

I'm new to image processing in Python, and I think I've mislead myself by using strange color spaces. Perhaps it would have been better to define an upper and lower boundary for each colour, and mask the image array? I'm definitely open to any flexible solution to the problem. Eventually I'll be applying this to hundreds of photos..

Original Photo

Code

from skimage import io, img_as_float, color, exposure

# Isolate paint marks

# Put image into LAB colour space
image_lab = color.rgb2lab(img)
img = exposure.rescale_intensity(img)

# Colours of interest
color_array = np.array([
[[[255, 255, 0.]]],    # Yellow stuff
[[[255, 190, 200.]]],    # Pink stuff
[[[255, 165, 0.]]],    # Orange stuff
[[[255, 0, 0.]]],    # Red stuff
])

# Loop through the color array and pick out the colored features
for i in range(0,color_array.ndim):
# Compute distance between the color of interest and the actual image
# http://scikit-image.org/docs/dev/api/skimage.color.html#skimage.color.deltaE_cmc
# "The usual values are kL=2, kC=1 for “acceptability”"
distance_color = color.deltaE_ciede2000(color_array[i], image_lab, kL=2, kC=1, kH=0.5)

# Normalise distance
distance_color = exposure.rescale_intensity(distance_color)

image_filtered = img.copy()
image_filtered[distance_color > 0.5] = 0

# Plot it up
print ("Filtered to: ", color_array[i])
f, (ax0, ax1, ax2) = plt.subplots(1, 3, figsize=(20, 10))
ax0.imshow(img)
ax1.imshow(distance_color, cmap='gray')
ax2.imshow(image_filtered)
plt.show()


Results

You can look at your picture in a different color space, e.g. the HSV space. In this space, each pixel has 3 components: hue, saturation and value. Hue defines the color's tone, saturation the color's strength and V the brightness of the pixel. For your particular problem, you can search for pixels that have e.g. a decent saturation (i.e. are reasonably colored and not just gray). The following code snippet uses OpenCV, but I believe you will be able to translate into scikit-image:

I = cv2.imread("/tmp/XApdz.jpg")  # Load the image

hsv = cv2.cvtColor(I, cv2.COLOR_BGR2HSV);  # Convert to HSV space
h, s, v = cv2.split(hsv);  # extract the single components

plt.figure(1, figsize=(10, 10))
plt.subplot(3,2,1); plt.imshow(h, cmap='gray'); plt.title('Hue')
plt.subplot(3,2,2); plt.imshow(s, cmap='gray'); plt.title('Saturation')
plt.subplot(3,2,3); plt.imshow(v, cmap='gray'); plt.title('Value')

mask = (s > 80).astype(np.uint8);  # Thresholding in the Saturation-channel

disk = cv2.getStructuringElement(cv2.MORPH_ELLIPSE,(3,3))  # remove very small regions
opened = cv2.morphologyEx(mask, cv2.MORPH_OPEN, disk);
plt.subplot(3,2,5); plt.imshow(opened); plt.title('Opened mask')

I2 = I.copy()
I2[opened.astype(bool), :] = 0; # Set the pixels to zero, where
plt.subplot(3,2,6); plt.imshow(I2)


Thanks to Maximilian Matthé, I have an answer. Below is his OpenCV code translated to Scikit-image. Note that I changed some parameters (e.g., reduce the size of the morphological disk) to replicate his results.

from skimage import color
from skimage.morphology import disk, opening, dilation

img_hsv = color.rgb2hsv(img) # Image into HSV colorspace
h = img_hsv[:,:,0] # Hue
s = img_hsv[:,:,1] # Saturation
v = img_hsv[:,:,2] # Value aka Lightness

plt.figure(1, figsize=(15, 15))

plt.subplot(4,2,1); plt.imshow(h, cmap='gray'); plt.title('Hue')
plt.subplot(4,2,2); plt.imshow(s, cmap='gray'); plt.title('Saturation')
plt.subplot(4,2,3); plt.imshow(v, cmap='gray'); plt.title('Value')
plt.tight_layout()

mask = (s > 0.35).astype(np.uint8);  # Thresholding in the Saturation-channel