User Name Password Register
DaniWeb IT Discussion Community
All
What is DaniWeb IT Discussion Community?
You're currently browsing the Legacy and Other Languages section within the Software Development category of DaniWeb, a massive community of 455,974 software developers, web developers, Internet marketers, and tech gurus who are all enthusiastic about making contacts, networking, and learning from each other. In fact, there are 3,807 IT professionals currently interacting right now! Registration is free, only takes a minute and lets you enjoy all of the interactive features of the site.
Please support our Legacy and Other Languages advertiser: Programming Forums
Views: 3007 | Replies: 0
Reply
Join Date: Nov 2007
Posts: 1
Reputation: rick_miller82 is an unknown quantity at this point 
Rep Power: 0
Solved Threads: 0
rick_miller82 rick_miller82 is offline Offline
Newbie Poster

Help Matlab code using Genetic Algorithm for a digital images correlation

  #1  
Nov 27th, 2007
I am a masters student working on a project for one of my professors and i am coming across an issue. I am working on a 6 variable GA scheme for Digital image correlation. I have created a majority of the code but i want to be able to add in the ability to to subsize imaging. the code right now will only return hole numbers for the amount of displacement in the u and v direction. i would like to be able to have the ability to have fractions of displacement in terms of the pixels. i have the code that will interpolate the image to make the number of pixels smaller but i am having trouble inserting it into the code, keep coming up with errors. i would like someone who is more knowledgable than me with matlab to see if they could help. i know that the interpolation needs to go inside the while loop when it is looking at the the subset size. if anyone could help please let me know. here is the code for the DIC

clear; clc;
u_displacement = 5;
v_displacement = 6;
size = 100;
%%%%%%% generate artificial array a1 and a2 %%%%%%%%%%
IMAGE = imread('8a.bmp');
a1 = double(IMAGE(1:size,1:size));
for i = 1:size
for j = 1:size
newrow(i,j) = i + v_displacement + round((j-1)*0.1 + (i-1)*0.2);
newcol(i,j) = j + u_displacement + round((j-1)*0.3 + (i-1)*0.4);
end
end
for i = 1:size
for j = 1:size
a2(newrow(i,j),newcol(i,j)) = a1(i,j);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%% a1 and a2 generated %%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

tpointrow = 1;
tpointcol = 1;
num_lookup_table = 2000; % number of values for u and v generated in the look up tables
subset_size = 60;
ipop = 500; % initial population of u and v indices
iguesu_lower = 0; % initial guess for lower limit of u
iguesu_upper = 15; % initial guess for upper limit of u
iguesv_lower = 0; % initial guess for lower limit of v
iguesv_upper = 10; % initial guess for upper limit of u
generations = 100;
crossovertimes = 3;
crossoverrate = 0.5;
wpop = 100;
b = 0.5;
mutationrate = 0.5;

fact = subset_size;%ceil(subset_size/2); %%%factor to be added and subtracted to point of interest

%%%% generate initial table of values for u and v %%%%%%%
for i = 1num_lookup_table+1)
ulookup1(i) = floor(rand*(iguesu_upper-1.5))+iguesu_lower;
vlookup1(i) = floor(rand*(iguesv_upper-1.5))+iguesv_lower;
s11(i) = (round(rand*0.5*10))/10;
s21(i) = (round(rand*0.5*10))/10;
s31(i) = (round(rand*0.5*10))/10;
s41(i) = (round(rand*0.5*10))/10;
end
[ulookup indx_ulookup] = sort(ulookup1);
[vlookup indx_vlookup] = sort(vlookup1);
[s1lookup indx_s1lookup] = sort(s11);
[s2lookup indx_s2lookup] = sort(s21);
[s3lookup indx_s3lookup] = sort(s31);
[s4lookup indx_s4lookup] = sort(s41);

for k = 1:ipop
u(k) = ceil(rand*num_lookup_table)+1;
v(k) = ceil(rand*num_lookup_table)+1;
ps1(k) = ceil(rand*num_lookup_table)+1;
ps2(k) = ceil(rand*num_lookup_table)+1;
ps3(k) = ceil(rand*num_lookup_table)+1;
ps4(k) = ceil(rand*num_lookup_table)+1;
end

%%%generate initial pop%%%%%%%%%%
for k = 1:ipop
num = 0;
denum1 = 0;
denum2 = 0;
corr1(k)= 0;
for row = (tpointrow)tpointrow+fact)
for col = (tpointcol)tpointcol+fact)
DX = (col - tpointcol);
DY = (row - tpointrow);
dudx = s1lookup(ps1(k));
dudy = s2lookup(ps2(k));
dvdx = s3lookup(ps3(k));
dvdy = s4lookup(ps4(k));clc;
r = row+vlookup(v(k))+round(dudx*DX+dudy*DY);
c = col+ulookup(u(k))+round(dvdx*DX+dvdy*DY);
num = num + a1(row,col)*a2(r,c);
% a1(row,col)
% a2(r,c)
denum1 = denum1 + (a1(row,col))^2;
denum2 = denum2 + (a2(r,c))^2;
end
end
corr1(k) = -num/sqrt(double(denum1)*double(denum2));
end
[sortedcost1 indx] = sort(corr1);
newu = u;
newv = v;
news1 = ps1;
news2 = ps2;
news3 = ps3;
news4 = ps4;
g = 1;
flag = 0;
corr = sortedcost1(1:wpop);
bestu(1) = newu(indx(1));
bestv(1) = newv(indx(1));
bests1(1) = news1(indx(1));
bests2(1) = news2(indx(1));
bests3(1) = news3(indx(1));
bests4(1) = news4(indx(1));

trackcost(1) = sortedcost1(1);
if (corr(1) <=-0.99)
flag = 1;
end
while (g<=generations)&(flag~=1)
newu1 = newu; clear newu;
newv1 = newv; clear newv;
news11 = news1; clear news1;
news21 = news2; clear news2;
news31 = news3; clear news3;
news41 = news4; clear news4;
for i = 1wpop)
newu(i) = newu1(indx(i));
newv(i) = newv1(indx(i));
news1(i) = news11(indx(i));
news2(i) = news21(indx(i));
news3(i) = news31(indx(i));
news4(i) = news41(indx(i));
end
clear sortedcost indx;
for crosstimes = 1:crossovertimes
if(rand<=crossoverrate)
indextocrossoveru = ceil(rand*wpop);
indextocrossoverv = ceil(rand*wpop);
indextocrossovers1 = ceil(rand*wpop);
indextocrossovers2 = ceil(rand*wpop);
indextocrossovers3 = ceil(rand*wpop);
indextocrossovers4 = ceil(rand*wpop);
swap = newu(indextocrossoveru);
newu(indextocrossoveru) = newv(indextocrossoverv);
newv(indextocrossoverv) = swap;
swap = news1(indextocrossovers1);
news1(indextocrossovers1) = news2(indextocrossovers2);
news2(indextocrossovers2) = swap;
swap = news3(indextocrossovers3);
news3(indextocrossovers3) = news4(indextocrossovers4);
news4(indextocrossovers4) = swap;
swap = news2(indextocrossovers2);
news2(indextocrossovers2) = news4(indextocrossovers4);
news4(indextocrossovers4) = swap;
end
end
for i = ((wpop/2)+1):2wpop-1)
newu(i) = ceil((b*newu(i-(wpop/2))+(1-b)*newu(i-(wpop/2-1))));
newu(i+1) = ceil(((1-b)*newu(i-(wpop/2))+(b)*newu(i-(wpop/2-1))));
newv(i) = ceil((b*newv(i-(wpop/2))+(1-b)*newv(i-(wpop/2-1))));
newv(i+1) = ceil(((1-b)*newv(i-(wpop/2))+(b)*newv(i-(wpop/2-1))));
news1(i) = ceil((b*news1(i-(wpop/2))+(1-b)*news1(i-(wpop/2-1))));
news1(i+1) = ceil(((1-b)*news1(i-(wpop/2))+(b)*news1(i-(wpop/2-1))));
news2(i) = ceil((b*news2(i-(wpop/2))+(1-b)*news2(i-(wpop/2-1))));
news2(i+1) = ceil(((1-b)*news2(i-(wpop/2))+(b)*news2(i-(wpop/2-1))));
news3(i) = ceil((b*news3(i-(wpop/2))+(1-b)*news3(i-(wpop/2-1))));
news3(i+1) = ceil(((1-b)*news3(i-(wpop/2))+(b)*news3(i-(wpop/2-1))));
news4(i) = ceil((b*news4(i-(wpop/2))+(1-b)*news4(i-(wpop/2-1))));
news4(i+1) = ceil(((1-b)*news4(i-(wpop/2))+(b)*news4(i-(wpop/2-1))));
end
%%%%% mutate%%%%%%
if(rand<=mutationrate)
newu(ceil(rand*wpop)) = floor(rand*(num_lookup_table-1.5))+1;
newv(ceil(rand*wpop)) = floor(rand*(num_lookup_table-1.5))+1;
news1(ceil(rand*wpop)) = floor(rand*(num_lookup_table-1.5))+1;
news2(ceil(rand*wpop)) = floor(rand*(num_lookup_table-1.5))+1;
news3(ceil(rand*wpop)) = floor(rand*(num_lookup_table-1.5))+1;
news4(ceil(rand*wpop)) = floor(rand*(num_lookup_table-1.5))+1;
end
%%%%%%%%%%
%%%%%%%% objective function
num = 0;
denum1 = 0;
denum2 = 0;
% if(rand>0.5)
% sign = -1;
% else
sign = 1;
% end
for k = 1:wpop
corr(k) = 0;
for row = (tpointrow)tpointrow+fact)
for col = (tpointcol)tpointcol+fact)
DX = (col - tpointcol);
DY = (row - tpointrow);
dudx = s1lookup(news1(k));
dudy = s2lookup(news2(k));
dvdx = s3lookup(news3(k));
dvdy = s4lookup(news4(k));
r = round(row+vlookup(newv(k))+dudx*DX+dudy*DY);
c = round(col+ulookup(newu(k))+dvdx*DX+dvdy*DY);
num = num + a1(row,col)*a2(r,c);
denum1 = denum1 + (a1(row,col))^2;
denum2 = denum2 + (a2(r,c))^2;
end
end
corr(k) = -num/sqrt((denum1)*(denum2));
end
[sortedcost indx] = sort(corr);
trackcost(g+1) = sortedcost(1);
bestu(g+1) = newu(indx(1));
bestv(g+1) = newv(indx(1));
bests1(g+1) = news1(indx(1));
bests2(g+1) = news2(indx(1));
bests3(g+1) = news3(indx(1));
bests4(g+1) = news4(indx(1));
if (sortedcost(1) <=-0.99)
flag = 1;
end
g = g+1;
end

if (flag == 1)
fprintf('Correct displacements were predicted after %d generations\n',(g))
fprintf(' u = %d pixels v = %d pixels\ns1 = %5.3f s2 = %5.3f s3 = %5.3f s4 = %5.3f',ulookup(bestu(g)),vlookup(bestv(g)),s1lookup(bests1(g)),s2lookup(bests2(g)),s3lookup(bests3(g)),s4lookup(bests4(g)))
else
fprintf('Correct displacements could not be predicted with sufficient accuracy even after %d generations\n',(g-1))
end
figure (1); plot(trackcost,' - o');
title('movement of correlation');
for count = 1:numel(bestu)
toplotu(count) = ulookup(bestu(count));
toplotv(count) = vlookup(bestv(count));
toplots1(count) = s1lookup(bests1(count));
toplots2(count) = s2lookup(bests2(count));
toplots3(count) = s3lookup(bests3(count));
toplots4(count) = s4lookup(bests4(count));
end
figure (4); plot(toplotu,' - *');
title('movement of x-displacement');
figure (5); plot(toplotv,' - *');
title('movement of y-displacement');
figure (6); plot(toplots1,' - o');
title('movement of s1');
figure (7); plot(toplots2,' - o');
title('movement of s2');
figure (8); plot(toplots3,' - o');
title('movement of s3');
figure (9); plot(toplots4,' - o');
title('movement of s4');
figure(10)
imshow(uint8(a1));
figure(11);
imshow(uint8(a2));

and here is the code for the interpolation (it does alot more than i need but the part of it that does the interpolation to create, lets say a 20 x 20 into a 200 x 200 pixel, is what i need.

% COR.M
% Digital Image Correlation File
% % Uses a tif file - 'speckle1.tif'
%
clear;
%
speckle=imread('speckle1.tif');
mesh(speckle);
pause;
%imshow(speckle);pause;
%imhist(speckle);pause;
s=size(speckle)
% change int to double to work with mesh
% to change back use uint8
a=double(speckle);
s1=size(a)
mesh(a),rotate3d on;
pause;
%surf(a), %b = a > 128;
%imshow(b); pause;
c=1000+a;
mesh(c),rotate3d on;
axis([0 300 0 300 0 1400]);
pause;
% Looking at a 10x10 area
d=zeros(240,320);
roi=zeros(10,10);
start_x=100;
start_y=100;
for i=start_x:start_x+9
for j=start_y:start_y+9
d(i,j)=1500;
roi(i-start_x+1,j-start_y+1)=c(i,j);
end
end
mesh(roi);
axis([1 10 1 10 0 1400]);
pause;
e=c+d;
mesh(e),rotate3d on;
axis([0 300 0 300 0 3000]);
pause;
%
% This increases resolution of the region of interest (100x100)
f=zeros(100,100);
for i=0:10:90
for l=0:10:90
for j=1:10
for k=1:10
f(j+i,k+l)=c(start_x+i/10,start_y+l/10);
end
end
end
end
mesh(f);pause;
%
%Using Interpolation (Can use 'bilinear', 'bicubic' or 'nearest')
%
xa=1:1:10;
ya=1:1:10;
[x,y]=meshgrid(xa,ya);
z=roi;
mesh(z);pause;
xa1=1:0.1:10;
ya1=1:0.1:10;
[x1,y1]=meshgrid(xa1,ya1);
z1=interp2(x,y,z,x1,y1,'bilinear'); %default is 'bilinear'
mesh(z1);
pause;
AddThis Social Bookmark Button
Reply With Quote  
Reply

Only community members can participate in forum threads. You must register or log in to contribute.

DaniWeb Legacy and Other Languages Marketplace
Currently Active Users Viewing This Thread: 1 (0 members and 1 guests)

 

Thread Tools Display Modes

Similar Threads
Other Threads in the Legacy and Other Languages Forum

All times are GMT -4. The time now is 9:15 am.
Forum system based on vBulletin Copyright ©2000 - 2008, Jelsoft Enterprises Ltd.
©2003 - 2008 DaniWeb® LLC