% Tracking pumice depth below water over time by turning images into binary images
% Created by Behnaz Hosseini, November 2015
% Adapted from Kristen Fauria's code tracking the dimensions of a density current, January 2015

clear all
close all

%% DEFINE PRELIMINARY VARIABLES OF INTEREST

pixel_scaling = 150; % pixels per centimeter
step = 100; % step size for images
dt = 60*step; % time interval
water_height = 130; % manually defined water level (rough estimate)
f = 0;

%% ACCESS FOLDER OF INTEREST

folder = uigetdir; % display folders inside current working directory
files = dir(fullfile(folder,'/*.jpg')); % extract file attributes in directory of interest
numfiles = length(files); % total number of image files
index = zeros(floor(numfiles/step),'double'); % create array of zeros

%% SET UP VIDEO OBJECT AND NAME VIDEO

writerObj = VideoWriter('FloatingPumice.avi');
open(writerObj);
set(figure, 'Position', [800, 1200, 800, 1300]); % set figure position for consistent axes during video

for i = 1:numfiles/step % loop through selected images
    I = imread(fullfile(folder,files(i*step).name)); % read in image
    
    disp(files(i*step).name) % display file names
    
    %% PROCESS IMAGE (I)
    
    % adjust skewed and misaligned images by rotation and translation
    s = files(i*step).name;
    snew = s(1:5);
    if str2num(snew) > 0 && str2num(snew) <= 4200
        Inew = imtranslate(I,[8,0]);
    elseif str2num(snew) > 4200 && str2num(snew) <= 5200
        Inew = imrotate(I,1);
        Inew = imtranslate(Inew,[2,-30]);
    else
        Inew = imrotate(I,1);
        Inew = imtranslate(Inew,[4,-30]);
    end
    
    % [J,newOrigin] = undistortImage(I,cameraParams); % lens distortion correction -- unnecessary on account of image cropping
    Inew = Inew(800:1400,925:1600,:); % crop image to pumice area of interest
    
    Inew_adj = imadjust(Inew,[0.1; 0.26],[0.99; 0.2]); % adjust image intensity values
    threshold = graythresh(Inew_adj); % define intensity threshold
    bw = im2bw(Inew_adj,threshold); % render image binary
    
    se = strel('disk', 5);
    bwnew = imclose(bw,se); % remove noise from image
    
    %% EXTRACT INFORMATION OF INTEREST
    
    % specify column of interest for series of images
    if str2num(snew) > 0 && str2num(snew) <= 4100
        col = 556;
    elseif str2num(snew) > 4100 && str2num(snew) <= 4900;
        col = 560;
    elseif str2num(snew) > 4900 && str2num(snew) <= 5000;
        col = 550;
    else str2num(snew) > 5000 && str2num(snew) <= 8000;
        col = 536;
    end
    
    column = bwnew(:,col); % specify column of interest
    position = find(column); % find first occurring nonzero in specified column (note: use ~ to find first '0')
    index(i) = position(1); % index the nonzero position
    
    depth(i)=(index(i)-water_height)/pixel_scaling; % define submerged depth using known water levels and depth of pumice from top of plot, and convert to centimeters
    
    %% CREATE VIDEOS AND PLOTS
    
    subplot(2,2,1), image(Inew); % original video
    title('Floating Pumice');
    subplot(2,2,2), image(bwnew); % binary video
    title('Floating Pumice in Binary');
    hold on
    plot(col,index(i),'ro'); % plot reference point tracking specified area of pumice below water
    hold on
    refline(0,130); % plot water level reference line
    %subplot(2,2,3), image(bwnew);
    %set(gca,'XLim',[300,600],'YLim',[200,500]);
    %hold on
    %plot(col,index(i),'ro');
    
    subplot(2,2,3),
    y = medfilt1(depth); % apply signal smoothing filter to plotted curve of depth vs. time
    yy = smooth(y); % smoothen data curve
    t(i) = i*dt; % time in seconds at each image
    plot(t,yy,'r');
    xlabel('Time (seconds)');
    ylabel('Submerged Depth (centimeters)');
    title('Submerged Depth of Havre Pumice as a Function of Time','Fontsize',10);
    axis([0,numfiles*60,1,2]);
    
    %{
    subplot(2,2,4), plot(y,'r');
    xlabel('Time (seconds)');
    ylabel('Saturation Level (%)');
    title('Saturation Level of Havre Pumice as a Function of Time');
    axis([0,numfiles*60,1,100]);
    %}
    
    colormap(bone(2.0));
    
    frame = getframe(gcf);
    writeVideo(writerObj,frame);
    
    f = f+1
    axis([0 numfiles*60 1 2])
    saveas(gcf,['pumice_',sprintf('%03d',f),'.png'])
end
close(writerObj);
