function [validx, validy, displx, disply]=RTCorrCode(grid_x,grid_y,Firstimagename,ImageFolder)

% Real time Correlation Code
%
% Written by Chris
% OPTIMIZE parameter is controlling which optimizations should be used. 
% <3 - The original version, no optimizations
%  3 - Most of computations are performed on NVidia card
%  4 - Optimize FFT sizes for better performance (affects precision)
%  5 - Load images in C code
OPTIMIZE = 4;
CORRSIZE = 15;	
PRECISION = 1000;

warning off Images:initSize:adjustingMag

if OPTIMIZE > 2
    hwid = normxcorr_hw();

    if hwid > 0
    else
	OPTIMIZE = 0;
    end
else
    hwid = 0;
end

if OPTIMIZE > 2
    RTselection = menu(sprintf('End processing by end.txt or by last image?'),...
	'Stop with end.txt','Stop with image check','Exit');

    if RTselection==1
    end

    if RTselection==2
    end

    if RTselection==3
        if hwid > 0
	    normxcorr_hw(hwid);
	end
	return
    end
end

% Filename

if exist('Firstimagename')==0
    [Firstimagename ImageFolder]=uigetfile('*.tif','Open First Image');
end


if ~isempty(Firstimagename)
% Get the number of image name
letters=isletter(Firstimagename);
Pointposition=findstr(Firstimagename,'.');
Firstimagenamesize=size(Firstimagename);
counter=Pointposition-1;
counterpos=1;
letterstest=0;

while letterstest==0
    letterstest=letters(counter);
    if letterstest==1
        break
    end
    Numberpos(counterpos)=counter;
    counter=counter-1;
    counterpos=counterpos+1;
    if counter==0
        break
    end
end

Filename_first = Firstimagename(1:min(Numberpos)-1);
Firstfilenumber=Firstimagename(min(Numberpos):max(Numberpos));
Lastname_first = Firstimagename(max(Numberpos)+1:Firstimagenamesize(1,2));
Firstfilenumbersize=size(Firstfilenumber);
onemore=10^(Firstfilenumbersize(1,2));
filenamelist(1,:)=Firstimagename;


h=figure;
if exist('grid_x')==0
    fpstest=1;
    Filelist=[Firstimagename;Firstimagename];
    while fpstest==1
        [grid_x,grid_y]=grid_generator(Firstimagename,ImageFolder);
	if OPTIMIZE > 2
	    ncp = prod(size(grid_x));
	    err = normxcorr_hw(hwid, 1, ncp, CORRSIZE, PRECISION, OPTIMIZE);
	else
	    err = 0;
	end
	if err == 0
	    if OPTIMIZE > 2
		base_points_x=single(grid_x);
		base_points_y=single(grid_y);
		
		normxcorr_hw(hwid, 3, base_points_x, base_points_y);

		if OPTIMIZE > 4
		    normxcorr_hw(hwid, 4, strcat(imagedir, Firstimagename));
		else
		    base = uint8(mean(double(imread([ImageFolder, Firstimagename])),3));
		    normxcorr_hw(hwid, 4, base);
		end	
		
		normxcorr_hw(hwid, 12, base_points_x, base_points_y);
	    end
	    [processingtime]=fpstestfunc(hwid,OPTIMIZE,grid_x,grid_y,Filelist,ImageFolder);
    	    fpstest = menu(sprintf(['Processing the selected grid will allow ' , num2str(1/processingtime),' frames per second' ]),'Try again','Use the grid');
    	    if fpstest==1
        	clear grid_x; clear grid_y;
    	    end
	else
	    ASelection = menu(sprintf('CUDA initialization failed?'),...
	'Retry another grid','Continue in software mode','Exit');
	    if Aselection==1
		fptest = 1;
	    end

	    if Aselection==2
		fptest = 0;
		normxcorr_hw(hwid);
		OPTIMIZE = 0;
	    end

	    if Aselection==3
		return
	    end
	end
    end
end

if OPTIMIZE < 3
    [validx, validy, displx, disply] = RTCorrCode_orig(grid_x, grid_y, Firstimagename, ImageFolder)
    return
end


Firstfilenumber=str2num(Firstfilenumber);
u=1+onemore+Firstfilenumber;
ustr=num2str(u);
filenamelist(2,:)=[Filename_first ustr(2:Firstfilenumbersize(1,2)+1) Lastname_first];
numberofimages=2;

counter=1;

input_points_x=single(grid_x);
input_points_y=single(grid_y);
normxcorr_hw(hwid, 12, input_points_x, input_points_y);

fit_options = optimset('Display', 'off');

numberofmarkers=max(size(grid_x))*min(size(grid_x));
validx(:,1)=reshape(grid_x,[],1);
displx=zeros(numberofmarkers,1);
validy(:,1)=reshape(grid_y,[],1);
disply=zeros(numberofmarkers,1);
tic;

while exist([ImageFolder, 'end.txt'],'file') ==0;
    pause(0.01);

    if exist([ImageFolder, filenamelist((counter+1),:)],'file') ==2;
%        warning(['# Processed Images: ', num2str(numberofimages-1),'; # markers:',num2str(numberofmarkers), '; Processing Image: ',filenamelist(counter+1,:)]);

	if OPTIMIZE > 4
	    normxcorr_hw(hwid, 13, strcat(ImageFolder, filenamelist((counter+1),:)));
	    input_correl = normxcorr_hw(hwid, 14);
	else
	    input = uint8(mean(double(imread([ImageFolder, filenamelist((counter+1),:)])),3));
	    normxcorr_hw(hwid, 13, input);
	    input_correl = normxcorr_hw(hwid, 14);
	end

	input_correl_x=double(input_correl(:,1));
        input_correl_y=double(input_correl(:,2));
	
        validx(:,counter+1)=input_correl_x;                                                     % lets save the data
        savelinex=input_correl_x';
    	dlmwrite([ImageFolder, 'resultsimcorrx.txt'], savelinex , 'delimiter', '\t', '-append');       % Here we save the result from each image; if you are desperately want to run this function with e.g. matlab 6.5 then you should comment this line out. If you do that the data will be saved at the end of the correlation step - good luck ;-)

    	validy(:,counter+1)=input_correl_y;
	saveliney=input_correl_y';
    	dlmwrite([ImageFolder, 'resultsimcorry.txt'], saveliney , 'delimiter', '\t', '-append');

        subplot(2,2,1);
        imshow([ImageFolder, filenamelist(counter+1,:)]);
        hold on;
        plot(grid_x,grid_y,'g+');
        plot(input_correl_x,input_correl_y,'r+');
        hold off;
        drawnow;

        displx(:,counter+1)=validx(:,counter+1)-validx(:,1);
        disply(:,counter+1)=validy(:,counter+1)-validy(:,1);

        subplot(2,2,2);
        xdata=validx(:,counter+1);
        ydata=displx(:,counter+1);
        if counter==1
            x(1)=0;
            x(2)=0;
        end
        [x,resnormx,residual,exitflagx,output]  = lsqcurvefit(@linearfit, [x(1) x(2)], xdata, ydata, [], [], fit_options);
        plot(xdata,ydata,'.');
        hold on;
        ydatafit=x(1)*xdata+x(2);
        plot(xdata,ydatafit,'r');
        hold off;
        xlabel('x-pos [pixel]');
        ylabel('x-displ [pixel]');
        title('x displ. versus x pos. in [pixel]');

        slopex(counter,:)=[i x(1) x(2)];

        subplot(2,2,4);
        xdata=validy(:,counter+1);
        ydata=disply(:,counter+1);
        if counter==1
            y(1)=0;
            y(2)=0;
        end
        [y,resnormx,residual,exitflagx,output]  = lsqcurvefit(@linearfit, [y(1) y(2)], xdata, ydata, [], [], fit_options);
        plot(xdata,ydata,'.g');
        hold on;
        ydatafit=y(1)*xdata+y(2);
        plot(xdata,ydatafit,'r');
        hold off
        xlabel('y-pos [pixel]')
        ylabel('y-displ [pixel]')
        title('y displ. versus y pos. in [pixel]');

        slopey(counter,:)=[i y(1) y(2)];

        subplot(2,2,3);
        plot(slopex(:,2),'-b');
        hold on;
        plot(slopey(:,2),'-g');
        hold off;
        xlabel('Image # [ ]');
        ylabel('x- and y-strain [ ]');
        title('Strain in x and y direction versus Image #');

        counter=counter+1;

        u=1+u;
        ustr=num2str(u);
        filenamelist(counter+1,:)=[Filename_first ustr(2:Firstfilenumbersize(1,2)+1) Lastname_first];
        [numberofmarkers numberofimages]=size(validx);
        
        if RTselection==2
            if exist([ImageFolder, filenamelist((counter+1),:)],'file') ==0;
                save ([ImageFolder, 'validx.dat'], 'validx', '-ascii', '-tabs');
                save ([ImageFolder, 'validy.dat'], 'validy', '-ascii', '-tabs');
                %warning('Last image detected, RTCorrCode stopped');
		normxcorr_hw(hwid);
                return
            end
        end
        
        
        subplot(2,2,1),title(['# Processed Images: ', num2str(numberofimages-1),'; fps: ', num2str((numberofimages-1)/toc),'; # markers:',num2str(numberofmarkers), '; Waiting for Image: ',filenamelist(counter+1,:)]);

    end
end

normxcorr_hw(hwid);

save ([ImageFolder, 'validx.dat'], 'validx', '-ascii', '-tabs');
save ([ImageFolder, 'validy.dat'], 'validy', '-ascii', '-tabs');
msgboxwicon=msgbox('end.txt file detected, RTCorrCode stopped','Processing stopped!');
%warning('end.txt file detected, RTCorrCode stopped');
end

%----------------------------------
%

function [processingtime]=fpstestfunc(hwid, OPTIMIZE, grid_x, grid_y, filenamelist, ImageFolder)

tic;

if hwid > 0
    if OPTIMIZE > 4
	normxcorr_hw(hwid, 13, strcat(ImageFolder, filenamelist(2,:)));
	input_correl = normxcorr_hw(hwid, 14);
    else
	input = uint8(mean(double(imread([ImageFolder, filenamelist(2,:)])),3));
	normxcorr_hw(hwid, 13, input)
	input_correl = normxcorr_hw(hwid, 14);
    end
    input_correl_x=double(input_correl(:,1));
    input_correl_y=double(input_correl(:,2));
else
    input_points_x=grid_x;
    base_points_x=grid_x;

    input_points_y=grid_y;
    base_points_y=grid_y;

    base = uint8(mean(double(imread([ImageFolder, filenamelist(1,:)])),3));            % read in the base image ( which is always  image number one. You might want to change that to improve correlation results in case the light conditions are changing during the experiment
    input = uint8(mean(double(imread([ImageFolder, filenamelist(2,:)])),3));       % read in the image which has to be correlated

    input_points_for(:,1)=reshape(input_points_x,[],1);         % we reshape the input points to one row of values since this is the shape cpcorr will accept
    input_points_for(:,2)=reshape(input_points_y,[],1);
    base_points_for(:,1)=reshape(base_points_x,[],1);
    base_points_for(:,2)=reshape(base_points_y,[],1);
    input_correl(:,:)=cpcorr(input_points_for, base_points_for, input, base);           % here we go and give all the markers and images to process to cpcorr.m which ic a function provided by the matlab image processing toolbox
    input_correl_x=input_correl(:,1);                                       % the results we get from cpcorr for the x-direction
    input_correl_y=input_correl(:,2);                                       % the results we get from cpcorr for the y-direction
end

processingtime=toc;