
clear all
clf

x = imread("3394-3.tif");

figure(2)
title("'Original' Image")
imshow(x)

% opening and closing by reconstruction

disk = strel('disk', 10);
opened_image = imopen(x, disk);

eroded_image = imerode(x, disk);
reconstructed_eroded_image = imreconstruct(eroded_image, x);

closed_image = imclose(opened_image, disk);

dilated_image = imdilate(reconstructed_eroded_image, disk);
reconstructed_dilated_image = imreconstruct(imcomplement(dilated_image), imcomplement(reconstructed_eroded_image));

final = reconstructed_dilated_image;

figure(3)
title("Relative Height")
imshow(final)

fibrofinal = final;
finalchange = process(final, 3.2); 
[row, col] = size(finalchange);


figure(4)
title("Binary Gradient - m0")
imshow(finalchange)

figure(5)
imshow(finalchange)
title("Binary Gradient - Overlayed Area")
[centers,radii] = imfindcircles(finalchange,[30 50],'ObjectPolarity','bright', ...
    'Sensitivity',0.97)
viscircles(centers, radii,'Color','b')

[o, p] = size(radii);

averageradii = sum(radii,1)/o;
cellarea = pi*(averageradii)^2;

% Cell Counter

cellpix = 0;
for a = 1:row
    for b = 1:col
        if (finalchange(a,b)== 255 || finalchange(a,b) == 256)
            cellpix = cellpix + 1;
            
        end
        
    end
end


m0coverage = cellpix/(row*col)
m0cellcount = round(cellpix/cellarea)

% Real stuff here

fibrochange = process(fibrofinal, 4);

figure(5)
imshow(fibrochange)
title('total image')

cellpix = 0;
for a = 1:row
    for b = 1:col
        if (fibrochange(a,b) > 252)
            cellpix = cellpix + 1;
            
        end
        
    end
end

totalcoverage = cellpix/(row*col)
fibroblastcoverage = totalcoverage - m0coverage
m0_3t3_ratio = 1/(fibroblastcoverage/m0coverage)
