% --- Executes on button press in analyse.
function analyse_Callback(hObject, eventdata, handles)
global electrode analyse_press
tic
go=1;

[disc,page]=size(electrode);
for i=1:page
    if electrode(i).voltage=='not'
        go=0;
        msgbox('Error: Voltages of all electrodes are not specified')
        break
    end
end
%above loop checks if all voltages have been entered
if page==0
    msgbox('Error: No electrodes drawn')
    go=0;
end

if go==1
    analyse_press=1; 
    msgbox('Press "OK", then wait till this msgbox disappears before proceeding ')
    ok(page)
    potanalyse()
end
toc


function ok(page)
global BC ID electrode
for i=1:page
    elem=electrode(i).elem;
    v=electrode(i).voltage;
    [melem,nelem]=size(elem);
    for s=1:melem
        ID(elem(s,2),elem(s,1))=1;
        BC(elem(s,2),elem(s,1))=v;
    end        
end
return


function potanalyse()
global mx my BC ID LOC col1 col2
global POT DX DY FXI FYI shape invisible
delx=1;    
dely=1;    
eps=8.854*1e-12;

beta=delx/dely;             
betasq=beta^2;              
gama=-2*(1+betasq);         

LOC=zeros(mx,my);
location=1;
for i=1:mx
    for j=1:my
        if ID(i,j)==1
            continue
        else
            LOC(i,j)=location;
            location=location+1;
        end
    end
end
location=location-1;
COF=zeros(location,location);
FV=zeros(location,1);

%ASSEMBLING CO-EFFICIENT MATRIX FOR CALCULATION OF POTENTIAL DISTRIBUTION%
for i=1:mx
    for j=1:my
        ptr=LOC(i,j);
        if ptr==0
            continue
        elseif (i==1)||(j==1)||(i==mx)||(j==my)
            if (i==1)&&(j>=1)&&(j<=my)
                COF(ptr,ptr)=COF(ptr,ptr)-1;
                COF(ptr,LOC(i+1,j))=1;
            end
            if (i==mx)&&(j>=1)&&(j<=my)
                COF(ptr,ptr)=COF(ptr,ptr)-1;
                COF(ptr,LOC(i-1,j))=1;
            end
            if (i>=1)&&(i<=mx)&&(j==1)
                COF(ptr,ptr)=COF(ptr,ptr)-1;
                COF(ptr,LOC(i,j+1))=1;
            end
            if (i>=1)&&(i<=mx)&&(j==my)
                COF(ptr,ptr)=COF(ptr,ptr)-1;
                COF(ptr,LOC(i,j-1))=1;
            end
        else
            COF(ptr,ptr)=gama;
            if ID(i,j+1)==1
                FV(ptr)=FV(ptr)-((BC(i,j+1))*1);
            else
                COF(ptr,LOC(i,j+1))=1;
            end
            
            if ID(i,j-1)==1
                FV(ptr)=FV(ptr)-((BC(i,j-1))*1);
            else
                COF(ptr,LOC(i,j-1))=1;
            end
        
            if ID(i+1,j)==1
                FV(ptr)=FV(ptr)-((BC(i+1,j))*betasq);
            else
                COF(ptr,LOC(i+1,j))=betasq;
            end
                      
            if ID(i-1,j)==1
                FV(ptr)=FV(ptr)-((BC(i-1,j))*betasq);
            else
                COF(ptr,LOC(i-1,j))=betasq;
            end
        end
    end
end
Q=COF\FV;

%ASSEMBLING POTENTIAL MATRIX%
POT=zeros(mx,my);
for i=1:mx
    for j=1:my
        ptr=LOC(i,j);
        if ptr==0
            POT(i,j)=BC(i,j);
        else
            POT(i,j)=Q(ptr);
        end
    end
end
clear COF FV Q
[DX,DY]=gradient(POT);
DX=-DX;
DY=-DY;

%OBTAINING POINTS ONLY ON ELECTRODE BOUNDARIES
boundary=[0 0];
[disc,page]=size(shape);
for i=1:page
    elem=shape(i).elem;
    x=min(elem(:,1));
    y=min(elem(:,2));
    xmax=max(elem(:,1));
    ymax=max(elem(:,2));
    w=xmax-x;
    h=ymax-y;
    if shape(i).id==0           %ELLIPSE
        xc=(x+xmax)/2;
        yc=(y+ymax)/2;
        t=1;
        w=w/2;
        h=h/2;
        for s=x:xmax
            q=(sqrt(1-(((s-xc)/w)^2)))*h;
            a=yc+q;
            ell(t,1)=s;
            ell(t,2)=a;
            t=t+1;
        end
        for s=(xmax-1):-1:(x+1)
            q=(sqrt(1-(((s-xc)/w)^2)))*h;
            b=yc-q;
            ell(t,1)=s;
            ell(t,2)=b;
            t=t+1;
        end
        boundary=vertcat(boundary,ell);
        clear ell
    elseif shape(i).id==1      %RECTANGLE
        t=1;
        for s=x:xmax
           rect(t,1)=s;
           rect(t,2)=ymax;
           t=t+1;
        end
        for s=(ymax-1):-1:y
            rect(t,1)=xmax;
            rect(t,2)=s;
            t=t+1;
        end
        for s=(xmax-1):-1:x
            rect(t,1)=s;
            rect(t,2)=y;
            t=t+1;
        end
        for s=(y+1):(ymax-1)
            rect(t,1)=x;
            rect(t,2)=s;
            t=t+1;
        end
        boundary=vertcat(boundary,rect);
        clear rect 
    end
end
boundary(1,:)=[];
[disc,page]=size(invisible);
for i=1:page
    elem=invisible(i).elem;
    x=min(elem(:,1));
    y=min(elem(:,2));
    xmax=max(elem(:,1));
    ymax=max(elem(:,2));
    w=xmax-x;
    h=ymax-y;
    if invisible(i).id==0           %ELLIPSE
        xc=(x+xmax)/2;
        yc=(y+ymax)/2;
        t=1;
        w=w/2;
        h=h/2;
        for s=x:xmax
            q=(sqrt(1-(((s-xc)/w)^2)))*h;
            a=yc+q;
            ell(t,1)=s;
            ell(t,2)=a;
            t=t+1;
        end
        for s=(xmax-1):-1:(x+1)
            q=(sqrt(1-(((s-xc)/w)^2)))*h;
            b=yc-q;
            ell(t,1)=s;
            ell(t,2)=b;
            t=t+1;
        end
        boundary=vertcat(boundary,ell);
        clear ell
    elseif invisible(i).id==1      %RECTANGLE
        t=1;
        for s=x:xmax
           rect(t,1)=s;
           rect(t,2)=ymax;
           t=t+1;
        end
        for s=(ymax-1):-1:y
            rect(t,1)=xmax;
            rect(t,2)=s;
            t=t+1;
        end
        for s=(xmax-1):-1:x
            rect(t,1)=s;
            rect(t,2)=y;
            t=t+1;
        end
        for s=(y+1):(ymax-1)
            rect(t,1)=x;
            rect(t,2)=s;
            t=t+1;
        end
        boundary=vertcat(boundary,rect);
        clear rect 
    end
end
col1=boundary(:,1);
col2=boundary(:,2);

%ASSEMBLING FORCE MATRIX ONLY ON ELECTRODES%
FX=zeros(mx,my);
FY=zeros(mx,my);
for i=1:mx
    for j=1:my
        if LOC(i,j)==0
            FX(i,j)=(eps*DX(i,j)^2)/2;
            if LOC(i,j+1)==0
                n=-1;
            else n=1;
            end
            FX(i,j)=n*FX(i,j);
            FY(i,j)=(eps*DY(i,j)^2)/2;
            if LOC(i+1,j)==0
                n=-1;
            else n=1;
            end
            FY(i,j)=n*FY(i,j);
        end
    end
end

FXI=interp2(FX,col1,col2);
FYI=interp2(FY,col1,col2);

save analys_var POT FXI FYI DX DY LOC col1 col2 ID
clear POT LOC col1 col2 ID
return

