function[temp]=SingleTrajectoryViewer(temp_fn, n, flag, k)

% "Olfactory behaviors assayed by computer tracking of Drosophila in a 
% four-quadrant olfactometer".
% 
% AUTHORS: 
% Chun-Chieh Lin, Olena Riabinina, Christopher J.Potter
% 
% The Solomon H. Snyder Department of Neuroscience, Center for Sensory 
% Biology, Johns Hopkins University School of Medicine, 
% Baltimore, MD 21205, USA. cpotter@jhmi.edu


% The output of the IFlySpy fly tracker groups tracked datapoints into
% trajectories that were generated by a fly. This code allows to plot out
% all datapoints, and then view trajectories one by one, by either clicking
% on an arbitrary datapoint, or selecting them automatically.

% temp_fn is the path to the data file

% n is the number or a range of numbers of trajectories that you want to
% view

% flag=0 means that trajectories will be displayed consequtively; flag=1
% means that all trajectories will be displayed together

% k sets the cutofff threshold: only trajectories that have >k datapoints
% will be displayed.

% example of use:
%
%  SingleTrajectoryViewer_250513('C:\Users\C1X-2.mat', 1:100, 0, 200)
%
% opens file C1X.mat, shown trajectories 1 to 100 one-by-one, skipping the
% ones that have less than 200 datapoints.
%
% If you would like to select trajectories by clicking on them, set the n
% large, i.e. n=1:100000, to make sure that the maximum of n is larger than
% the number of trajectories in the dataset.

% load data file that contains cleaned up and filtered data: 
% x, y coordinates of the datapoints,
% corresponding time t, and variable iTS, that is a cell array, each cell
% of which contains indexes of x,y and t that belong to the same
% trajectory.

global fn
fn=temp_fn;
load (fn)

% create the trajectory ID's to display
if nargin<2
    n=1:10;
end;

if nargin<3
    flag=0;
end;

if nargin<4
    k=100;
end;



%create array of all useful indexes, inds
inds=cell2mat(iTS);

% plot "good" datapoins in grey 
figure('WindowStyle', 'docked')

% Plot all datapoints as grey dots
plot(x,y,'.', 'MarkerEdgeColor', [0.3 0.3 0.3], 'MarkerSize',3)
axis equal
hold on

% Display filename in the title of the figure, and user instructions in the
% y-axis label
sy{1}='Press Enter to continue, "b" to go step back, "q" to quit';
sy{2}='or click on a point to highlight trajectory';
ylabel(sy)

st{1}='Start of trajectory is blue, end - red';
title(st,'Interpreter', 'none')


  

set(gca,'NextPlot', 'add')

i=n(1); j=1;

b=[];

% Go through the first n trajectories of the dataset
while i<=n(end)

    % inds contains indexes of the i-th trajectory
    inds=[iTS{i}];
    
    % coordinates of datapoints within a single i-th trajectory
    x_tr=x(inds);
    y_tr=y(inds);
    
    % check if this trajectory is longer than the cut-off k
    if length(x_tr)>k
        trajectory_length=length(x_tr)/30  
        
        % create a color map, to progres from blue atthe beginning ofthe
        % trajectory to red at the end
        cmap=buildcmap(1:length(x_tr));
        h_loc=[]; hold on;
        
        % plot out the trajectory
        for loc_ind=1:length(x_tr)
            h_loc(loc_ind)=plot(x_tr(loc_ind),y_tr(loc_ind),'.',...
                'MarkerEdgeColor', cmap(loc_ind,:));
        end;
    
    % plot first and last points of the trajectory as larger circles
    if ~isempty(x_tr)
    h_loc1=plot(x_tr(1), y_tr(1), 'ob');
    h_loc2=plot(x_tr(end), y_tr(end), 'or');
    end;
    
    
    if flag==0 
        % display trajectory number, and set the pointer to record a
        % coordinate, ic any is clicked in the figure
        st{2}=['Trajectory number ', num2str(i)];
        title(st,'Interpreter', 'none')
        [p, q, b]=myginput(1,'arrow')
        
        if b==98 % go back
            j=j-2;
        end;
        
        if b==113 % quit
            break
        end;
            
        if b==1 % a point was clicked on the figure. Show trajectory for that point

            % find a trajectory, closes to the point clicked
            [~,ind_f]=min(sqrt((x-p).^2+(y-q).^2));
            a=cellfun(@(x) find(x==ind_f), iTS,'UniformOutput', false);

            for l=1:size(a,1)
                if ~isempty(a{l})
                    ind_tr=l
                end;
            end;
            
            j=j; i=ind_tr;
        else j=j+1; i=n(j);
        end;
            
        % delete previous trajectory   
        delete(h_loc)
        if ~isempty(x_tr)
            delete(h_loc1);
            delete(h_loc2);
        end;
    end;
    
    elseif b==98
        j=j-1; i=n(j);
    else 
        j=j+1; i=n(j);
    end;

end;


