/normxcorr/trunk

function [grid_x,grid_y]=grid_generator(FileNameBase,PathNameBase);

% Code to generate the DIC analysis grid

% Completely rewritten by Chris

% Programmed first by Dan and Rob 

% 

% Last revision: 12/27/06

% The grid_generator function will help you create grids of markers. The

% dialog has different options allowing you to create a marker grid which is rectangular,

% circular, a line or two rectangels of a shape or contains only of two

% markers. After choosing one of the shapes you will be asked for the base

% image which is typically your first image. After opening that image you

% will be asked to click at the sites of interest and the markers will be

% plotted on top of your image. You can choose if you want to keep these

% markers or if you want to try again.

% It has to be noted that you can

% always generate your own marker positions. Therefore the marker position

% in pixel has to be saved as a text based format where the x-position is

% saved as grid_x.dat and the y-position saved as grid_y.dat.

%

% Prompt user for base image

if exist('FileNameBase')==0

[FileNameBase,PathNameBase] = uigetfile( ...

    {'*.bmp;*.tif;*.jpg;*.TIF;*.BMP;*.JPG','Image files (*.bmp,*.tif,*.jpg)';'*.*',  'All Files (*.*)'}, ...

    'Open base image for grid creation');

end

curdir = pwd;

cd(PathNameBase);

im_grid = imread(FileNameBase);

[grid_x,grid_y,FileNameBase,PathNameBase] = gridtypeselection(FileNameBase, PathNameBase, im_grid);

close all

cd(curdir);

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

%

% Decide which type of grid you want to create

function [grid_x,grid_y,FileNameBase,PathNameBase] = gridtypeselection(FileNameBase, PathNameBase, im_grid);

hold off

imshow(im_grid,'InitialMagnification', 100);

gridselection = menu(sprintf('Which type of grid do you want to use'),...

    'Rectangular','Circular','Two Markers','Line','Two Rectangles of Markers','Cancel');

if gridselection==1

    [grid_x,grid_y,FileNameBase,PathNameBase] = rect_grid(FileNameBase, PathNameBase, im_grid);

    return

end

if gridselection==2

    [grid_x,grid_y,FileNameBase,PathNameBase] = circ_grid(FileNameBase, PathNameBase, im_grid);

    return

end

if gridselection==3

    [grid_x,grid_y,FileNameBase,PathNameBase] = twop_grid(FileNameBase, PathNameBase, im_grid);

    return

end

if gridselection==4

    [grid_x,grid_y,FileNameBase,PathNameBase] = line_grid(FileNameBase, PathNameBase, im_grid);

    return

end

if gridselection==5

    [grid_x,grid_y,FileNameBase,PathNameBase] = tworect_grid(FileNameBase, PathNameBase, im_grid);

    return

end

if gridselection==6

    return;

end

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

%

% Define two rectangles and add them to one marker array

function [grid_x,grid_y,FileNameBase,PathNameBase] = tworect_grid(FileNameBase, PathNameBase, im_grid);

[grid_x1,grid_y1,FileNameBase,PathNameBase] = rect_grid(FileNameBase, PathNameBase, im_grid);

imshow(im_grid,'InitialMagnification', 100);

[grid_x2,grid_y2,FileNameBase,PathNameBase] = rect_grid(FileNameBase, PathNameBase, im_grid);

grid_x1=reshape(grid_x1,[],1);

grid_x2=reshape(grid_x2,[],1);

grid_y1=reshape(grid_y1,[],1);

grid_y2=reshape(grid_y2,[],1);

grid_x=[grid_x1; grid_x2];

grid_y=[grid_y1; grid_y2];

imshow(im_grid,'InitialMagnification', 100);

hold on

plot(grid_x,grid_y,'.')

title(['Selected grid has ',num2str(length(grid_x)), ' rasterpoints'])    % plot a title onto the image

% Accept the chosen markers, try again or give up 

confirmcircselection = menu(sprintf('Do you want to use these markers?'),...

    'Yes','No, try again','Go back to grid-type selection');

if confirmcircselection==2

    close all

    hold off

    imshow(im_grid,'InitialMagnification', 100);

    tworect_grid(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==3

    close all

    gridtypeselection(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==1

    close all

    save grid_x.dat grid_x -ascii -tabs

    save grid_y.dat grid_y -ascii -tabs

end

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

%

% Define line and create markers

function [grid_x,grid_y,FileNameBase,PathNameBase] = line_grid(FileNameBase, PathNameBase, im_grid);

title(sprintf('Pick two points on the sample.') )

[x(1,1),y(1,1)]=ginput(1);

hold on

plot(x(1,1),y(1,1),'+g')

[x(2,1),y(2,1)]=ginput(1);

plot(x(2,1),y(2,1),'+g')

linelength=sqrt((x(2,1)-x(1,1))*(x(2,1)-x(1,1))+(y(2,1)-y(1,1))*(y(2,1)-y(1,1)));

lineslope=(y(2,1)-y(1,1))/(x(2,1)-x(1,1));

intersecty=y(1,1)-lineslope*x(1,1);

ycalc=zeros(2,1);

ycalc=lineslope*x+intersecty;

plot(x(:,1),ycalc(:,1),'-b')

prompt = {'Enter the number of intersections between markers on the line:'};

dlg_title = 'Input for grid creation';

num_lines= 1;

def     = {'30'};

answer = inputdlg(prompt,dlg_title,num_lines,def);

linediv = str2num(cell2mat(answer(1,1)));

linestep=((max(x)-min(x))/linediv);

grid_x(1:linediv+1)=min(x)+linestep*(1:linediv+1)-linestep;

grid_y=lineslope*grid_x+intersecty;

plot(grid_x,grid_y,'ob')

title(['Selected grid has ',num2str(linediv), ' rasterpoints'])    % plot a title onto the image

% Accept the chosen markers, try again or give up 

confirmcircselection = menu(sprintf('Do you want to use these markers?'),...

    'Yes','No, try again','Go back to grid-type selection');

if confirmcircselection==2

    close all

    hold off

    imshow(im_grid,'InitialMagnification', 100);

    twop_grid(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==3

    close all

    gridtypeselection(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==1

    save grid_x.dat grid_x -ascii -tabs

    save grid_y.dat grid_y -ascii -tabs

end

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

%

% Select two markers

function [grid_x,grid_y,FileNameBase,PathNameBase] = twop_grid(FileNameBase, PathNameBase, im_grid);

title(sprintf('Pick two points on the sample.') )

[x(1,1),y(1,1)]=ginput(1);

hold on

plot(x(1,1),y(1,1),'+g')

[x(2,1),y(2,1)]=ginput(1);

plot(x(2,1),y(2,1),'+g')

% Accept the chosen markers, try again or give up 

confirmcircselection = menu(sprintf('Do you want to use these two markers?'),...

    'Yes','No, try again','Go back to grid-type selection');

if confirmcircselection==2

    close all

    hold off

    imshow(im_grid,'InitialMagnification', 100);

    twop_grid(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==3

    close all

    gridtypeselection(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==1

    grid_x=x;

    grid_y=y;

    save grid_x.dat grid_x -ascii -tabs

    save grid_y.dat grid_y -ascii -tabs

end

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

%

% Select a circular area

function [grid_x,grid_y,FileNameBase,PathNameBase] = circ_grid(FileNameBase, PathNameBase, im_grid);

title(sprintf('Pick three points on the circle in clockwise order at the upper boundary of the sample.') )

[x(1,1),y(1,1)]=ginput(1);

hold on

plot(x(1,1),y(1,1),'+g')

[x(2,1),y(2,1)]=ginput(1);

plot(x(2,1),y(2,1),'+g')

[x(3,1),y(3,1)]=ginput(1);

plot(x(3,1),y(3,1),'+g')

xnew=x;

ynew=y;

% Calculate center between the 3 sorted points and the normal slope of the vectors

slope12=-1/((ynew(2,1)-ynew(1,1))/(xnew(2,1)-xnew(1,1)));

slope23=-1/((ynew(3,1)-ynew(2,1))/(xnew(3,1)-xnew(2,1)));

center12(1,1)=(xnew(2,1)-xnew(1,1))/2+xnew(1,1);

center12(1,2)=(ynew(2,1)-ynew(1,1))/2+ynew(1,1);

center23(1,1)=(xnew(3,1)-xnew(2,1))/2+xnew(2,1);

center23(1,2)=(ynew(3,1)-ynew(2,1))/2+ynew(2,1);

% plot(center12(1,1),center12(1,2),'+b')

% plot(center23(1,1),center23(1,2),'+b')

if slope12==slope23

    return

end

% Calculate the crossing point of the two vectors

achsenabschnitt1=center12(1,2)-center12(1,1)*slope12;

achsenabschnitt2=center23(1,2)-center23(1,1)*slope23;

xdata=min(x):max(x);

ydata1=achsenabschnitt1+slope12*xdata;

ydata2=achsenabschnitt2+slope23*xdata;

% plot(xdata,ydata1,'-b')

% plot(xdata,ydata2,'-b')

xcross=(achsenabschnitt2-achsenabschnitt1)/(slope12-slope23);

ycross=slope12*xcross+achsenabschnitt1;

plot(xcross,ycross,'or')

% Calculate radius and plot circle

R=sqrt((xcross-xnew(1,1))*(xcross-xnew(1,1))+(ycross-ynew(1,1))*(ycross-ynew(1,1)));

% ydata=ycross-sqrt(R*R-(xdata-xcross).*(xdata-xcross));

% plot(xdata,ydata,'-b')

% Calculate angle between vectors

xvector=[1;0];

x1vec(1,1)=xnew(1,1)-xcross;x1vec(2,1)=ynew(1,1)-ycross

x3vec(1,1)=xnew(3,1)-xcross;x3vec(2,1)=ynew(3,1)-ycross

alpha13=acos((dot(x1vec,x3vec))/(sqrt(x1vec'*x1vec)*sqrt(x3vec'*x3vec)))*180/pi;

alpha01=acos((dot(xvector,x1vec))/(sqrt(x1vec'*x1vec)*sqrt(xvector'*xvector)))*180/pi;

alpha03=acos((dot(xvector,x3vec))/(sqrt(xvector'*xvector)*sqrt(x3vec'*x3vec)))*180/pi;

totalangle=alpha13;

minangle=alpha01;

maxangle=alpha03;

angldiv=abs(round(totalangle))*10;

anglstep=(totalangle/angldiv);

anglall(1:angldiv+1)=maxangle+anglstep*(1:angldiv+1)-anglstep;

xcircle(1:angldiv+1)=xcross+R*cos(-anglall(1:angldiv+1)/180*pi);

ycircle(1:angldiv+1)=ycross+R*sin(-anglall(1:angldiv+1)/180*pi);

plot(xcircle,ycircle,'-b')

drawnow

title(['Segment of circle spreads over ',num2str(totalangle),'°'])

% Accept the chosen circle, try again or give up 

confirmcircselection = menu(sprintf('Do you want to use this circle as basis?'),...

    'Yes','No, try again','Go back to grid-type selection');

if confirmcircselection==2

    close all

    imshow(im_grid,'InitialMagnification', 100);

    circ_grid(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==3

    close all

    gridtypeselection(FileNameBase, PathNameBase, im_grid);

end

if confirmcircselection==1

    prompt = {'Enter the number of intersections between markers on the circle:'};

    dlg_title = 'Input for grid creation';

    num_lines= 1;

    def     = {'30'};

    answer = inputdlg(prompt,dlg_title,num_lines,def);

    angldiv = str2num(cell2mat(answer(1,1)));

    anglstep=(totalangle/angldiv);

    anglall(1:angldiv+1)=maxangle+anglstep*(1:angldiv+1)-anglstep;

    markerxpos(1:angldiv+1)=xcross+R*cos(-anglall(1:angldiv+1)/180*pi);

    markerypos(1:angldiv+1)=ycross+R*sin(-anglall(1:angldiv+1)/180*pi);

    plot(markerxpos,markerypos,'ob');

    % Pick the lower bound in the image

    title(sprintf('Pick three points lying on the circle in clockwise order. The first and last one define the width of the raster') )

    [x(4,1),y(4,1)]=ginput(1);

    hold on

    plot(x(1,1),y(1,1),'+r')

    lowboundx=x(4,1);

    lowboundy=y(4,1);

    R2=sqrt((xcross-lowboundx(1,1))*(xcross-lowboundx(1,1))+(ycross-lowboundy(1,1))*(ycross-lowboundy(1,1)));

    markerxposlb(1:angldiv+1)=xcross+R2*cos(-anglall(1:angldiv+1)/180*pi);

    markeryposlb(1:angldiv+1)=ycross+R2*sin(-anglall(1:angldiv+1)/180*pi);

    plot(markerxposlb,markeryposlb,'ob');

    prompt = {'Enter the number of intersections between the upper and lower bound:'};

    dlg_title = 'Input for grid creation';

    num_lines= 1;

    def     = {'5'};

    answer = inputdlg(prompt,dlg_title,num_lines,def);

    Rdiv = str2num(cell2mat(answer(1,1)));

    Rstep=((R-R2)/Rdiv);

    Rall(1:Rdiv+1)=R2+Rstep*(1:Rdiv+1)-Rstep;

    grid_x=ones(Rdiv+1,angldiv+1)*xcross;

    grid_y=ones(Rdiv+1,angldiv+1)*ycross;

    A=Rall;

    B=cos(-anglall(1:angldiv+1)/180*pi);

    C=A'*B;

    grid_x=grid_x+Rall'*cos(-anglall(1:angldiv+1)/180*pi);

    grid_y=grid_y+Rall'*sin(-anglall(1:angldiv+1)/180*pi);

    close all

    imshow(im_grid,'InitialMagnification', 100);

    hold on

    plot(grid_x,grid_y,'.b')    

    title(['Selected grid has ',num2str(angldiv*Rdiv), ' rasterpoints'])    % plot a title onto the image

    % Do you want to keep the grid?

    confirmselection = menu(sprintf('Do you want to use this grid?'),...

        'Yes','No, try again','Go back to grid-type selection');

    if confirmselection==1

        % Save settings and grid files in the image directory for visualization/plotting later

        %         save settings.dat xspacing yspacing xmin_new xmax_new ymin_new ymax_new -ascii -tabs

        save grid_x.dat grid_x -ascii -tabs

        save grid_y.dat grid_y -ascii -tabs

    end

    if confirmselection==2

        close all

        hold off

        imshow(im_grid,'InitialMagnification', 100);

        circ_grid(FileNameBase, PathNameBase, im_grid);

    end

    if confirmselection==3

        gridtypeselection(FileNameBase, PathNameBase, im_grid);

    end

end

return

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

%

function [grid_x,grid_y,FileNameBase,PathNameBase] = rect_grid(FileNameBase, PathNameBase, im_grid);

title(sprintf('Define the region of interest.  Pick (single click) a point in the LOWER LEFT region of the gage section.\n  Do the same for a point in the UPPER RIGHT portion of the gage section.'))

[x(1,1),y(1,1)]=ginput(1);

hold on

plot(x(1,1),y(1,1),'+b')

[x(2,1),y(2,1)]=ginput(1);

hold on

plot(x(2,1),y(2,1),'+b')

drawnow

xmin = min(x);

xmax = max(x);

ymin = min(y);

ymax = max(y);

lowerline=[xmin ymin; xmax ymin];

upperline=[xmin ymax; xmax ymax];

leftline=[xmin ymin; xmin ymax];

rightline=[xmax ymin; xmax ymax];

plot(lowerline(:,1),lowerline(:,2),'-b')

plot(upperline(:,1),upperline(:,2),'-b')

plot(leftline(:,1),leftline(:,2),'-b')

plot(rightline(:,1),rightline(:,2),'-b')

% closereq

cd(PathNameBase)

% Prompt user for grid spacing/resolution

prompt = {'Enter horizontal (x) resolution for image analysis [pixels]:', ...

        'Enter vertical (y) resolution for image analysis [pixels]:'};

dlg_title = 'Input for grid creation';

num_lines= 1;

def     = {'50','50'};

answer = inputdlg(prompt,dlg_title,num_lines,def);

xspacing = str2num(cell2mat(answer(1,1)));

yspacing = str2num(cell2mat(answer(2,1)));

% Round xmin,xmax and ymin,ymax "up" based on selected spacing

numXelem = ceil((xmax-xmin)/xspacing)-1;

numYelem = ceil((ymax-ymin)/yspacing)-1;

xmin_new = (xmax+xmin)/2-((numXelem/2)*xspacing);

xmax_new = (xmax+xmin)/2+((numXelem/2)*xspacing);

ymin_new = (ymax+ymin)/2-((numYelem/2)*yspacing);

ymax_new = (ymax+ymin)/2+((numYelem/2)*yspacing);

% Create the analysis grid and show user

[x,y] = meshgrid(xmin_new:xspacing:xmax_new,ymin_new:yspacing:ymax_new);

[rows columns] = size(x);

zdummy = 200.*ones(rows,columns);

imshow(FileNameBase)

title(['Selected grid has ',num2str(rows*columns), ' rasterpoints'])    % plot a title onto the image

hold on;

plot(x,y,'+b')

grid_x=x;

grid_y=y;

% Do you want to keep the grid?

confirmselection = menu(sprintf('Do you want to use this grid?'),...

    'Yes','No, try again','Go back to grid-type selection');

if confirmselection==1

    % Save settings and grid files in the image directory for visualization/plotting later

    save settings.dat xspacing yspacing xmin_new xmax_new ymin_new ymax_new -ascii -tabs

    save grid_x.dat x -ascii -tabs

    save grid_y.dat y -ascii -tabs

    close all

    hold off

end

if confirmselection==2

    close all

    hold off

    imshow(im_grid,'InitialMagnification', 100);

    rect_grid(FileNameBase, PathNameBase, im_grid);

end

if confirmselection==3

    close all

    hold off

    gridtypeselection(FileNameBase, PathNameBase, im_grid);

end