%%%% Justin Sutherland
%%%% 9/3/2019
%%%%
%%%% MATLAB script for implementing the adaptive thresholding image
%%%% processing algorithm with background removal/adjustment. Also - an
%%%% initial test of a k-medoid clustering algorithm as a possibile
%%%% canditate for a vein injection site algorithm.
​
%%% Start by pulling desired image into suitable data type
rawData1 = importdata('test_image_v2.png');
[~,name] = fileparts('test_image_v2.png');
image.(matlab.lang.makeValidName(name)) = rawData1;
%%% Convert to grayscale image
temp = image.test_image_v2(:, :, 1);
[B, L] = bwboundaries(temp);
cmap = gray(256);
cmap(9,:) = [1,0,0];
colormap(cmap); % Apply the custom colormap.
imagesc(temp)
title('Original Image')
​
%%% Now perform masking operation
%%% Create mask by thresholding the entire image into a binary image;
%%% If above threshold -> Black; below threshold -> White
%%% Adjust bias value until desired result occurs.
bias = 0.2;
level = graythresh(temp) - bias;
mask = imbinarize(temp,level);
figure
imshowpair(temp,mask,'montage')
title('Original Image versus Mask')
​
%%% Now create an image that hightlights the veins as effectively as
%%% possible. Use the mask to remove the remaining border that the adaptive
%%% thresholding sometimes is unable to remove effectively.
new_image = RemoveBackground(temp, mask);
%%% now perform adaptive filtering
level = adaptthresh(new_image, 0.67); %%% this 0.67 threshold is arbitrary and I am actively changing it from run to run...
seg_I = imbinarize(new_image,level);
seg_I = RemoveBackground(seg_I, mask);
figure
title('Image After Filtering')
imshow(seg_I, [])
​
​
%%% Now, let's do a k-medoid operation on this final image. The k-medoid
%%% operation should return some n-number of centroids.
%%% In order to do this, we need to make a 2D array that contains only
%%% values from the image greater than a certain grayscale threshold. Let's
%%% call this thresh.
%%% Also, this array needs have the following data: x and y coordinates.
datapoints = getPoints(seg_I, 0.5);
datapoints = transpose([datapoints(2, :); datapoints(1, :)]);
%%% Alright, with these datapoints, now perform the k-medoid operation.
​
%%% Literally stole this segment of code directly from Mathworks example on
%%% their website.
n_medoids = 1;
opts = statset('Display','final');
[idx,C] = kmedoids(datapoints,n_medoids,'Distance','cityblock',...
    'Replicates',5,'Options',opts);
figure;
plot(datapoints(:, 1), datapoints(:, 2), '.', 'LineWidth', 0.1)
hold on
plot(C(:,1),C(:,2),'kx',...
     'MarkerSize',15,'LineWidth',3)
legend('Cluster 1','Cluster 2','Centroids',...
       'Location','NW')
title 'Cluster Assignments and Centroids'
hold off
​
function [datapoints] = getPoints(image, thresh)
    %%% Function for converting image to set of 2D arrays
    i = 1;
    datapoints_x = [];
    datapoints_y = [];
    iterations = size(image);
    for j = 1:iterations(1)
        for k = 1:iterations(2)
            if image(j, k) < thresh
                datapoints_x(i) = -j;
                datapoints_y(i) = k;
                i = i + 1;
            end
        end
    end
    datapoints = [datapoints_x; datapoints_y];
end
​
function [new_image] = RemoveBackground(original, mask)
​
    image_size = size(original);
    if size(mask) ~= image_size
        'Cannot remove background of two different image sizes'
        return;
    end

    new_image = original;
    %%% iterate through each image, checking values vs mask. If the mask
    %%% value is > ~200 then make the respective pixel on the original
    %%% image
    for j = 1:(image_size(1))
        for k = 1:(image_size(2))
            if mask(j, k) < 0.5
                new_image(j, k) = 10;
            end
        end
    end
end