%% connect particles in time

%% User input

fn='track_*';  % Track files to process (nb: only 1 supported)

%% Load track data
dr=dir(fn);
Ntrk=length(dr);
if(Ntrk~=1)
	fprintf('More than one track file.\n');
	return;
end

if(~exist('tpx','var'))
	fprintf('Loading track data...');
	load(dr(1).name,'ppx','ppy','base','Npf','D','w');
	Np=fix(mean(Npf));
	if(Np~=mean(Npf))
		fprintf('Connect Error:: Missing Particles');
		return;
	end
	tpx=squeeze(mean(ppx(1:Np,2:end,:),2));  % only connect central particles in tetrames
	tpy=squeeze(mean(ppy(1:Np,2:end,:),2));
	fprintf('Done\n');
end

clear('ppx','ppy');

%%  Load background file

bgfn=['bg_' base '.mat'];
load(bgfn);

%% Calculated Parameters

sfn=['connect_' base '.mat'];

Np=max(Npf);

[MaxBalls Nf]=size(tpx);   % get Max balls and Number of frames

%% Define

cc=zeros(Np,Nf);
sxx=tpx;
syy=tpy;

%% Main Loop

fprintf(1,'Processing...');
Nc=length(sprintf('%d',Nf+1));
fmt=sprintf('[%%0%dd/%%0%dd]',Nc,Nc);
fprintf(1,fmt,0,Nf+1);

cc(:,1)=1:Np;
for nf=1:Nf-1           % loop over all frames
  rx1=repmat(tpx(:,nf),1,Np); 
	ry1=repmat(tpy(:,nf),1,Np); 
	rx2=repmat(tpx(:,nf+1),1,Np); 
	ry2=repmat(tpy(:,nf+1),1,Np); 
	r=(rx1.'-rx2).^2+(ry1.'-ry2).^2;
	
	[junk tmp]=min(r);
	cc(:,nf+1)=tmp(cc(:,nf));
	sxx(:,nf+1)=tpx(cc(:,nf+1),nf+1);
	syy(:,nf+1)=tpy(cc(:,nf+1),nf+1);
	if(rem(nf,100)==0)
		backspace(1,fmt,nf,Nf-1);
		fprintf(1,fmt,nf,Nf-1);
% 		plot(sxx(1,1:nf),syy(1,1:nf));
% 		axis('equal');
% 		drawnow;
	end
end
backspace(1,fmt,nf,Nf-1);
fprintf(1,fmt,nf,Nf-1);
fprintf(1,'Done.\n');
    

%% Save results

fprintf('Saving connection infomation...');
save(sfn,'cc','tpx','tpy','sxx','syy');
fprintf('Done.\n');








