% Laboratory in Oceanography - Data and Methods
% SMAST, UMass Dartmouth
%
% For going through simple examples of image processing toolbox
%
% Written by Miles A. Sundermeyer, 4/28/09

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Read and write images
imfinfo('pout.tif')
I = imread('pout.tif'); 
imshow(I)
whos I

%%
imfinfo('StarWars.tif')
II = imread('StarWars.tif');
whos II
imshow(II(:,:,1:3));
imwrite (II(:,:,1:3), 'StarWars.png'); 

figure(1)
clf
subplot(2,2,1)
image(II(:,:,1:3))

subplot(2,2,2)
image(II(:,:,1))

subplot(2,2,3)
image(II(:,:,2))

subplot(2,2,4)
image(II(:,:,3))

%%
imtool(II(:,:,1:3))


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Read and display the grayscale image rice.png.
figure(1);
clf

%% Step 1: Read Image
I = imread('rice.png');
imshow(I)

%% Step 2: Use Morphological Opening to Estimate the Background
background = imopen(I,strel('disk',15));

%% Step 3: View the Background Approximation as a Surface
figure(2)
clf
surf(double(background(1:8:end,1:8:end))),zlim([0 255]);
colormap('jet')
set(gca,'ydir','reverse');

%%
imshow(background);

%% Step 4: Subtract the Background Image from the Original Image
I2 = I - background;
imshow(I2)

%% Step 5: Increase the Image Contrast
I3 = imadjust(I2);
imshow(I3);

figure(3)
I3too = (I2-min(I2(:)))/range(I2(:)) * 255;
imshow(I3too)




















%% 
I3too = uint8(double((I2-min(I2(:))))/double(range(I2(:))) * 255);
imshow(I3too)

%% Step 6: Threshold the Image
level = graythresh(I3);     % computes threshold to convert to b/w
bw = im2bw(I3,level);       % makes into binary image (logical 0's and 1's)
bw = bwareaopen(bw, 50);    % remove connected areas w/ fewer than P pixels
imshow(bw)

%%
bwtoo = bwareaopen(bw, 90);
imshow(bwtoo)

%%
figure(4)
clf
imhist(I3(:))
ax = axis;
hold on
plot(level*max(double(I3(:)))*[1 1],[ax(3:4)],'r');

%% Step 7: Identify Objects in the Image
if(strcmp(version('-release'),'2009a'))
  cc = bwconncomp(bw, 4)      % Matlab 2009

%% Step 8: Examine One Object
  grain = false(size(bw));
  grain(cc.PixelIdxList{50}) = true;
  imshow(grain);

%% Step 9: View All Objects
  labeled = labelmatrix(cc);
  whos labeled
  RGB_label = label2rgb(labeled, @spring, 'c', 'shuffle');
  imshow(RGB_label)
  
%% Step 10: Compute Area of Each Object
  graindata = regionprops(cc, 'basic')
  graindata(50).Area
  
%% Step 11: Compute Area-based Statistics
  grain_areas = [graindata.Area];
  [min_area, idx] = min(grain_areas)
  grain = false(size(bw));
  grain(cc.PixelIdxList{idx}) = true;
  imshow(grain);
end

%%
if(strcmp(version('-release'),'2007a'))
  [L,N] = bwlabel(bw,4);
  
%% Step 8: Examine One Object
  grain = false(size(bw));
  grain(L==50) = true;
  imshow(grain);

%% Step 9: View All Objects
labeled = uint8(L);
  whos labeled
  RGB_label = label2rgb(labeled, @spring, 'c', 'shuffle');
  imshow(RGB_label)
  
%% Step 10: Compute Area of Each Object
  graindata = regionprops(L, 'basic')
  graindata(50).Area
  
%% Step 11: Compute Area-based Statistics
  grain_areas = [graindata.Area];
  [min_area, idx] = min(grain_areas)
  grain = false(size(bw));
  grain(L==idx) = true;
  imshow(grain);
end

%% Step 12: Create Histogram of the Area
nbins = 20;
figure(5)
hist(grain_areas, nbins)
title('Histogram of Rice Grain Area');