% The script reads Tracking.mat and outputs 
%  (1) mean speed of tracers vs. time as a PNG file
%  (2) time-averaged mean speed, displayed in the Command Window
%
% Usage:
%  (1) Copy the script to the same directory as Tracking.mat
%  (2) Fill the parameters in "Set parameters" section
%  (3) Execute the script by pressing F5
%
% Written by Kun-Ta Wu, June 2019

%% Set parameters
clc; close all; clear all;

unitT = ;               % The time period between frames (sec)
pixel = ;               % Image pixel size (um)

xmin = ; xmax = ;       % Plotting range in x axis
ymin = ; ymax = ;       % Plotting range in y axis

% Average window
averageTimeMin = ;      % minimum bound (sec)
averageTimeMax = ;      % maximum bound (sec)

%% Load Data

load('Tracking.mat');
speed = sqrt(u.^2+v.^2); % Calculate tracer speeds

clearvars -except timeIndices t u v x y speed unitT pixel xmin xmax ...
                  ymin ymax averageTimeMin averageTimeMax;

%% Group data according to frames

rxT    = cellfun(@(a) x(a),timeIndices,'UniformOutput',false); 
ryT    = cellfun(@(a) y(a),timeIndices,'UniformOutput',false);
vxT    = cellfun(@(a) u(a),timeIndices,'UniformOutput',false);
vyT    = cellfun(@(a) v(a),timeIndices,'UniformOutput',false);
speedT = cellfun(@(a) speed(a),timeIndices,'UniformOutput',false);

%% Plot mean speed vs time

% calculate meanSpeedT
meanSpeedT = squeeze(            ...
                     cell2mat(   ...
                              cellfun(@(a) mean(a), speedT,'UniformOutput',false)));
time = squeeze((1:length(meanSpeedT))');

% Plot mean speed vs time
scatter(time*unitT,meanSpeedT*pixel/unitT,10,'filled','MarkerFaceColor',[0 0 1]);

% Label plot axes
box on;
xlabel('Time (sec)','FontSize',15,'Color',[0 0 0]);
ylabel('Mean Speed (\mum/s)','FontSize',15,'Color',[0,0,0]);
set(gca,'FontSize',15);

% Define plotting ranges in both axes
xlim([xmin xmax]);
ylim([ymin ymax]);

% Export figure
saveas(gcf,'Mean Speed vs Time.png','png');

%% Calculate time-averaged mean speeds

% Convert the averaging window from seconds to frame numbers
minFrame = floor(averageTimeMin/unitT);
maxFrame = floor(averageTimeMax/unitT);

selectedMeanSpeeds = meanSpeedT(minFrame:maxFrame); % Select mean speeds within the given window

% Calculate the mean speed and its standard deviation
timeAveragedMeanSpeed = mean(selectedMeanSpeeds);
standardDeviation     = std(selectedMeanSpeeds);

% Display the result in the Command Window
clc;
display(sprintf('\n Time-Averaged Mean Speed from %d to %d sec = %f %c %f %cm/s\n',...
        averageTimeMin,averageTimeMax,timeAveragedMeanSpeed,char(177),standardDeviation,char(956)));
