Image processing

% An RGB image is rows × columns × 3: one layer each for red, green and blue
I = testimage('shapes');            % built-in test image (uint8, 0–255)
size(I)
R = I(:, :, 1); G = I(:, :, 2); B = I(:, :, 3);
grey = rgb2gray(I);
subplot(2, 3, 1), imshow(I), title('RGB')
subplot(2, 3, 2), imshow(R), title('red')
subplot(2, 3, 3), imshow(G), title('green')
subplot(2, 3, 4), imshow(B), title('blue')
subplot(2, 3, 5), imshow(grey), title('rgb2gray')
subplot(2, 3, 6), imshow(255 - I), title('negative')

% Gaussian noise vs salt-and-pepper noise, and which filter suits which
rng(3)
I = rgb2gray(testimage('shapes'));
gn = imnoise(I, 'gaussian', 0, 0.01);
sp = imnoise(I, 'salt & pepper', 0.06);
subplot(2, 3, 1), imshow(gn), title('Gaussian noise')
subplot(2, 3, 2), imshow(imfilter(gn, fspecial('average', 5), 'replicate')), title('5×5 average')
subplot(2, 3, 3), imshow(imgaussfilt(gn, 1.5)), title('Gaussian blur σ = 1.5')
subplot(2, 3, 4), imshow(sp), title('salt & pepper')
subplot(2, 3, 5), imshow(imgaussfilt(sp, 1.5)), title('Gaussian: smears dots')
subplot(2, 3, 6), imshow(medfilt2(sp)), title('median 3×3: removes them')

I = rgb2gray(testimage('noisy'));
[Gmag, Gdir] = imgradient(I);
subplot(2, 3, 1), imshow(I), title('input (noisy)')
subplot(2, 3, 2), imshow(Gmag, []), title('gradient magnitude')
subplot(2, 3, 3), imshow(Gdir, []), title('gradient direction')
subplot(2, 3, 4), imshow(edge(I, 'sobel')), title('Sobel')
subplot(2, 3, 5), imshow(edge(I, 'canny')), title('Canny')
subplot(2, 3, 6), imshow(imoverlay(testimage('shapes'), edge(I, 'canny'), 'k')), title('Canny on the image')

% A dull, low-contrast picture (values squeezed between 90 and 160) and two ways to fix it
I = imgaussfilt(rgb2gray(testimage('noisy')), 1);
dull = uint8(90 + 0.3 * double(I));
stretched = imadjust(dull);
equalised = histeq(dull);
subplot(2, 3, 1), imshow(dull), title('low contrast')
subplot(2, 3, 2), imshow(stretched), title('imadjust (stretch)')
subplot(2, 3, 3), imshow(equalised), title('histeq (equalise)')
subplot(2, 3, 4), imhist(dull), title('histogram')
subplot(2, 3, 5), imhist(stretched)
subplot(2, 3, 6), imhist(equalised)

% Threshold, clean up, label and measure: the classic way to count cells, grains or particles
rng(7)
I = testimage('blobs');
BW = imbinarize(I);                       % Otsu's threshold
BW = imopen(BW, strel('disk', 2));        % remove specks
BW = imfill(BW, 'holes');
[L, n] = bwlabel(BW);
s = regionprops(L, 'Area', 'Centroid', 'EquivDiameter');
subplot(1, 2, 1), imshow(I), title('input')
subplot(1, 2, 2), imshow(label2rgb(L)), hold on
plot(s.Centroid(:, 1), s.Centroid(:, 2), 'k+', 'MarkerSize', 8), hold off
title(sprintf('%d objects, mean diameter %.1f px', n, mean(s.EquivDiameter)))

% Greek coins from Pompeii (public domain, scikit-image test data) — fetched live
I = imread('https://raw.githubusercontent.com/scikit-image/scikit-image/v0.22.0/skimage/data/coins.png');
BW = imbinarize(I);                                % Otsu's global threshold
BW = imfill(imopen(BW, strel('disk', 2)), 'holes');
BW = imclearborder(bwareaopen(BW, 800));           % drop pieces smaller than a coin, and the bright edge
s = regionprops(BW, 'Area', 'Centroid', 'BoundingBox');
subplot(1, 2, 1), imshow(I), title('coins.png')
subplot(1, 2, 2), imshow(imoverlay(I, bwperim(BW), 'r')), hold on
plot(s.Centroid(:, 1), s.Centroid(:, 2), 'y+', 'MarkerSize', 10), hold off
title(sprintf('%d of the 24 coins found', s.Count))
% The dark coin at the bottom left is about as bright as the background, so one global
% threshold misses it: a known limit of this simple method (try lowering the threshold).

% The standard test image for CT reconstruction, with noise, denoising and edges
rng(1)
P = phantom(256);
noisy = imnoise(P, 'gaussian', 0, 0.005);
smooth = imgaussfilt(noisy, 1.2);
subplot(2, 2, 1), imshow(P), title('phantom(256)')
subplot(2, 2, 2), imshow(noisy), title('with noise')
subplot(2, 2, 3), imshow(smooth), title('Gaussian filtered')
subplot(2, 2, 4), imshow(edge(smooth, 'canny')), title('Canny edges')

% Eileen Collins, NASA astronaut (public domain, scikit-image test data) — fetched live
I = imread('https://raw.githubusercontent.com/scikit-image/scikit-image/v0.22.0/skimage/data/astronaut.png');
small = imresize(I, 0.5);
hsv = rgb2hsv(small);
hsv(:, :, 2) = min(1, hsv(:, :, 2) * 1.6);          % more saturated colours
vivid = hsv2rgb(hsv);
subplot(2, 2, 1), imshow(small), title(sprintf('%d×%d', size(small, 2), size(small, 1)))
subplot(2, 2, 2), imshow(imrotate(small, 20, 'bilinear')), title('imrotate 20°')
subplot(2, 2, 3), imshow(imsharpen(small, 'Amount', 1.5)), title('imsharpen')
subplot(2, 2, 4), imshow(vivid), title('saturation × 1.6 (HSV)')