% POSITION2TENSION : script 
%
% find correspondance between Vx,Vy (galvo tensions) and x,y (positions of the bead)
%
% Syntax: position2tension.m
%
% Inputs:
%    
% Outputs:
%    
% Other m-files required: 
% Subfunctions: 
% MAT-files required: 
%

% Author: Chardes.C, IBDM
% email address: claire.chardes@univ-amu.fr
% 2012 

%------------- BEGIN CODE --------------

% default folder
if ~ exist ('prevFolder', 'var')
    prevFolder ='E:\images\OT\SPIM\billes';
end

nMovies = str2double(cell2mat(inputdlg('How many calibration movies do you want to use ?', 'Calibration movies number', 1,{'1'})));

x=[];
y=[];
Vx = [];
Vy = [];
for i=1:nMovies

    nIm0 = str2double(cell2mat(inputdlg('first image number','First image',1,{'1'})));
    nIm1 = str2double(cell2mat(inputdlg('end image number','End image',1,{'1'})));

    % tension text file
    [fname, ffolder] = uigetfile('.txt', 'select text file', prevFolder);
    % extract tensions from file
    [Vxi, Vyi] = getVxVy(ffolder, fname);
    
    % calibration images
    [fname, ffolder] = uigetfile('.tif', 'select calibration tif', ffolder);
    info = imfinfo([ffolder,fname]);
    if (size(info,1)<numel(Vxi))
        Vxi = Vxi(end-size(info,1)+1:end);
        Vyi = Vyi(end-size(info,1)+1:end);
    end
    % localize bead positions on each images
    [xi, yi] = getXYPos(ffolder, fname, nIm0, nIm1);
    
    Vxi = Vxi(nIm0:nIm1)';
    Vyi = Vyi(nIm0:nIm1)';
    x = [x; xi];
    y = [y; yi];
    Vx = [Vx; Vxi];
    Vy = [Vy; Vyi];

end

ind = find (isnan(x)==0);
x = x(ind);
y = y(ind);
Vx = Vx(ind);
Vy = Vy(ind);

% we are looking for a function that pass through all the points
% If we have enough points, we can then make a meshgrid to visualize it

% Interpolate scattered data and find a function to convert tension to pos
fx = TriScatteredInterp(Vx,Vy,x);
fy = TriScatteredInterp(Vx,Vy,y);

% save the converting function
filename = [ffolder, 'convertFct.mat'];
save(filename, 'fx', 'fy');

% Obtain a function to Vx and Vy in function of x and y
fVx = TriScatteredInterp(x,y,Vx);
fVy = TriScatteredInterp(x,y,Vy);
filename = [ffolder, 'convertFctTension.mat'];
save(filename, 'fVx', 'fVy');

% Display
figure(1);
[Xinterp,Yinterp] = meshgrid(-0.4:0.05:0.4,-0.4:0.05:0.4);
Zinterp = fx(Xinterp,Yinterp);
contourf(gca,Xinterp,Yinterp,Zinterp)
pcolor(gca, Xinterp,Yinterp,Zinterp);
shading(gca,'interp')
col = colorbar;
set(get(col,'title'),'String','[x pos]','fontweight','bold')
set(gca,'tickdir','out','fontweight','bold')
xlabel('Vx (V)')
ylabel('Vy (V)')
title('X position vs tensions')
box('off')
shg
filename = [ffolder, 'XvsTensions.jpg'];
saveas(gcf, filename, 'jpg');

figure(2);
[Xinterp,Yinterp] = meshgrid(-0.4:0.05:0.4,-0.4:0.05:0.4);
Zinterp = fy(Xinterp,Yinterp);
contourf(gca,Xinterp,Yinterp,Zinterp)
pcolor(gca, Xinterp,Yinterp,Zinterp);
shading(gca,'interp')
col = colorbar;
set(get(col,'title'),'String','[y pos]','fontweight','bold')
set(gca,'tickdir','out','fontweight','bold')
xlabel('Vx (V)')
ylabel('Vy (V)')
title('Y position vs tensions')
box('off')
shg
filename = [ffolder, 'YvsTensions.jpg'];
saveas(gcf, filename, 'jpg');

figure(3);
ti = 0:0.5:512;
[qx,qy] = meshgrid(ti,ti);
qzx = fVx(qx,qy);
mesh(qx,qy,qzx);
hold on;
plot3(x,y,Vx,'o');
hold off

figure(4);
ti = 0:0.5:512;
[qx,qy] = meshgrid(ti,ti);
qzy = fVy(qx,qy);
mesh(qx,qy,qzy);
hold on;
plot3(x,y,Vy,'o');
hold off

% save a text file with x, y, Vx, Vy for optical tweezer positionning
% during experiment
f = fopen([ffolder, 'xy_VxVy_bis.txt'], 'a+','n' );
ind = ~isnan(qzx);
qzx = qzx(ind);
qzy = qzy(ind);
qx = qx(ind);
qy = qy(ind);
for i=1:numel(qx)
  
    fprintf(f, '%6.3f %6.3f %5.4f %5.4f\n', qx(i), 512-qy(i), qzx(i), qzy(i));
  
end
 fclose(f);

% Then, we give the list of tensions (Vx, Vy), and we calculate the values
% with fx and fy. We take the values from the text files, that we give as 
% inputs to the functions fx and fy : 
% qx = fx(xi, yi);
% qy = fy(xi,yi);


%------------- END OF CODE --------------