source: ETALON/CLIO/Analysis/Data19July/KK/PositionDirection.m @ 723

function [ChiV,Corr,ExterProfile]=PositionDirection(prof2,prof)
[ChiV1,Corr1,ExterProfile1]=Position(prof2,prof);
[ChiV2,Corr2,ExterProfile2]=Position(fliplr(prof2),prof);
if ChiV1<ChiV2
    ChiV=ChiV1;
    Corr=Corr1;
    ExterProfile=ExterProfile1;
else
    ChiV=ChiV2;
    Corr=Corr2;
    ExterProfile=ExterProfile2;
end
end


function [ChiV,Corr,ExterProfile]=Position(prof2,prof)
c1=find(cumsum(prof2)>=(sum(prof2)/2),1);
c2=find(cumsum(prof)>=(sum(prof)/2),1);
fi=(1+sqrt(5))/2;
  T=100;
a=c2-T; b=c2+T;
while(abs(b-a)>1)
x1=b-(b-a)/fi;
x2=a+(b-a)/fi;
y1=Chi(circshift(prof2,[ 1 round(x1)-c1]),prof);
y2=Chi(circshift(prof2,[ 1 round(x2)-c1]),prof);
if(y1<=y2)
b=x2; 
else
a=x1; 
end
end
y1=Chi(circshift(prof2,[ 1 round(a)-c1]),prof);
y2=Chi(circshift(prof2,[ 1 round(b)-c1]),prof);
if ((y1-y2)>0)  rez=b; ChiV=y2;
else rez=a;ChiV=y1
end
ExterProfile=circshift(prof2,[ 1 round(rez)-c1]);
Corr=sum(prof.*ExterProfile)/(sqrt(sum((prof).^2))*sqrt(sum((ExterProfile).^2)));
end

function Chi=Chi(funk1,funk2)
Chi=sum((funk1-funk2).^2./abs(funk1+funk2+1e-10)./(length(funk1)));
end
% %one side
% c1=centerOFweight(abs(funk1));%find center of weight of function
% c2=centerOFweight(funk2b);%base funk
% c1=c1-(T+1);
% for i=1:(2*T+1)
% c1=c1+1;
%     delta=c2-c1;
%     N1=zeros(size(funk1));
% if delta>0
%     N1(delta:end)=funk1(1:length(funk1)-delta+1);   
% end
% if c1==c2
%     N1=funk1;
% end
% if delta<0
%     N1(1:end+delta+1)=funk1(-delta:end);   
%  end
% 

% [chi1(i),cor1(i)]=ChiCor(N1,funk2b);
% end

% [chi1V,chi1P]=min(chi1);
% [cor1V,cor1P]=max(cor1);
% 
% %Other side
% funk1fliped=flip(funk1);
% c1=centerOFweight(abs(funk1fliped));%find center of weight of function
% c2=centerOFweight(funk2b);
% c1=c1-(T+1);
% for i=1:(2*T+1)
% c1=c1+1;
%     delta=c2-c1;
%     N1=zeros(size(funk1));
% if delta>0
%     N1(delta:length(N1))=funk1fliped(1:length(funk1)-delta+1);   
% end
% if c1==c2
%     N1=funk1;
% end
% if delta<0
%     N1(1:length(N1)+delta+1)=funk1fliped(-delta:length(funk1));   
%  end
% 
% [chi2(i),cor2(i)]=ChiCor(N1,funk2b);
% end
% [chi2V,chi2P]=min(chi2);
% [cor2V,cor2P]=max(cor2);
% 
% %compare both sides
% 
% 
% %by minimum of chi2
% 
% if chi1V<chi2V
%     c1=centerOFweight(abs(funk1))-(T+1)+chi1P;
%     c2=centerOFweight(funk2b);
%      delta=c2-c1;
%     N1=zeros(size(funk1));
% if delta>0
%     N1(delta:length(N1))=funk1(1:length(funk1)-delta+1);   
% end
% if c1==c2
%     N1=funk1;
% end
% if delta<0
%     N1(1:length(N1)+delta+1)=funk1(-delta:length(funk1));   
% end
% ChiVal1=chi1V;
% CorVal1=cor1V;
% else
%     
%     c1=centerOFweight(abs(funk1fliped))+chi2P-(T+1);
%     c2=centerOFweight(funk2b);
%         delta=c2-c1;
%     N1=zeros(size(funk1));
% if delta>0
%     N1(delta:length(N1))=funk1fliped(1:length(funk1)-delta+1);   
% end
% if c1==c2
%     N1=funk1;
% end
% if delta<0
%     N1(1:length(N1)+delta+1)=funk1fliped(-delta:length(funk1));   
% end
%   ChiVal1=chi2V;
%   CorVal1=cor2V;
% end
% 
% % %by maximum of corr
% % 
% % if cor1V>cor2V
% %     c1=centerOFweight(abs(funk1))-201+cor1P;
% %     c2=centerOFweight(funk2b);
% %      delta=c2-c1;
% %     N1c=zeros(size(funk1));
% % if delta>0
% %     N1c(delta:length(N1c))=funk1(1:length(funk1)-delta+1);   
% % end
% % if c1==c2
% %     N1c=funk1;
% % end
% % if delta<0
% %     N1c(1:length(N1c)+delta+1)=funk1(-delta:length(funk1));   
% % end
% % ChiVal2=chi1V;
% % CorVal2=cor1V;
% % else
% %     
% %     c1=centerOFweight(abs(funk1fliped))+cor2P-201;
% %     c2=centerOFweight(funk2b);
% %         delta=c2-c1;
% %     N1c=zeros(size(funk1));
% % if delta>0
% %     N1c(delta:length(N1c))=funk1fliped(1:length(funk1)-delta+1);   
% % end
% % if c1==c2
% %     N1c=funk1;
% % end
% % if delta<0
% %     N1c(1:length(N1c)+delta+1)=funk1fliped(-delta:length(funk1));   
% % end
% %   ChiVal2=chi2V;
% %   CorVal2=cor2V;
% % end
% 
% ExterProfile=N1;
% %ExterProfileCor=N1c;