obraz = imread ('klocki2.png'); 
%obraz = imread ('klocki1.png'); 
wymiary=size(obraz);
if(wymiary(1)~=358 || wymiary(2)~=433)
   skalax=wymiary(1)/358;
   skalay=wymiary(2)/433;
   if(skalax<skalay&&skalax>0&&skalay>0)
       x=floor(wymiary(1)/skalay);
       y=floor(wymiary(2)/skalay);       
       obraz=imresize(obraz,[x y]);
      %  x=385-x;
       % if(x>0)
       % obraz=[obraz, uint8(zeros(x,y,3))];
        %end
        %y=433-y;
        %if(y>0)
       %     for k=1:385
       %     obraz(k)=[obraz(k,:,:), uint8(zeros(y,3))];
        %    end
  % end
   elseif(skalax>skalay&&skalax>0&&skalay>0)
       y=ceil(wymiary(2)/skalax);
       x=ceil(wymiary(1)/skalax);
       obraz=imresize(obraz,[x y]);
       % x=385-x;
       % if(x>0)
       % obraz=[obraz, uint8(zeros(x,y,3))];
       % end
       % y=433-y;
        %if(y>0)
         %   for k=1:385
          %      obraz(k)=[obraz(k,:,:), uint8(zeros(y,3))];
         % 
          %  end
 %       end
   %else
    %     x=385-x;
     %   if(x>0)
      %  obraz=[obraz, uint8(zeros(x,y,3))];
       % end
        %y=433-y;
        %if(y>0)
         %    for k=1:385
          %       obraz(k)=[obraz(k,:,:), uint8(zeros(y,3))];
          %
           % end
            
    %    end
   end
end

gray=rgb2gray(obraz); 
BW = im2bw(obraz);
BW1 = edge(gray,'canny');
BW2 = edge(gray,'sobel');
SE = strel('octagon',3);
% BW3 = imdilate(BW,SE);
% BW4 = imerode(BW3,SE);
BW5 = -1 *(BW -1);
BW6 = imerode(BW5,SE);
%imshow(BW6 + BW4)
BW_filled2 = imfill(BW6,'holes');
%BW7 = imerode(BW_filled2,SE);
BW7 = imdilate(BW_filled2,SE);
%boundaries = bwboundaries(BW_filled2);
boundaries = bwboundaries(BW7);
%imshow(obraz)
%hold on
rozmiar=size(boundaries);
%for k=1:rozmiar(1)
%b = boundaries{k};
%   plot(b(:,2),b(:,1),'g','LineWidth',1);
%end
%hold off
%dystans=1;
for k=1:rozmiar(1)
   b = boundaries{k};
   rozmiar2=size(b);
   nachylenie=zeros(rozmiar2(1)-1,4);
   for l=1:(rozmiar2(1)-1)
     nachylenie(l,:)=[(b(l,1) - b(l+1,1)),(b(l,2) - b(l+1,2)), b(l,1), b(l,2)];
   end
end
k1=1;
k2=1;
rozmiar=size(boundaries);
for k=1:rozmiar(1)
   b = boundaries{k};
   rozmiar2=size(b);
   pochodne=zeros(rozmiar2(1)-1,4);
   for l=1:(rozmiar2(1)-1)
       %disp(l);
       %disp(round((rozmiar2(1)-4)/10)-1)
     if(nachylenie(l,1)~=0)  
     pochodnex(k1,:)=[nachylenie(l,3), nachylenie(l,4), (nachylenie(l,2)/nachylenie(l,1))];%, 1];
     k1=k1+1;
     else
     pochodney(k2,:)=[nachylenie(l,3), nachylenie(l,4), (nachylenie(l,1)/nachylenie(l,2))];%, 0];  
     k2=k2+1;
     end
     %disp('koniec')
   end
end

wspolczynnikix=unique(pochodnex(:,3));
lambda1=max(wspolczynnikix);
lambda2=min(wspolczynnikix);
lambda3=mean(wspolczynnikix);
k1=1;
k2=1;
k3=1;
for k=1:size(pochodnex,1)
 if(pochodnex(k,3)==lambda1)
    gora(k1,:)=[pochodnex(k,1) pochodnex(k,2)];
    k1=k1+1;
 elseif(pochodnex(k,3)==lambda2)
    dol(k2,:)=[pochodnex(k,1) pochodnex(k,2)];
    k2=k2+1;
 else
    pion(k3,:)=[pochodnex(k,1) pochodnex(k,2)];
    k3=k3+1;
 end
end 
poziom=[pochodney(:,1) pochodney(:,2)];
k1=1;
k2=1;
for k=1:(size(gora,1)-1)
    if(abs(gora(k,1)-gora(k+1,1))<2)
      prosta_gora(k1,k2,:)=[gora(k,1) gora(k,2)];
      k2=k2+1;
    else
      prosta_gora(k1,k2,:)=[gora(k,1) gora(k,2)];
      k1=k1+1; 
      k2=1;
    end
end
k1=1;
k2=1;
for k=1:(size(dol,1)-1)
    if(abs(dol(k,1)-dol(k+1,1))<2)
      prosta_dol(k1,k2,:)=[dol(k,1) dol(k,2)];
      k2=k2+1;
    else
      prosta_dol(k1,k2,:)=[dol(k,1) dol(k,2)];
      k1=k1+1; 
      k2=1;

    end
end
k1=1;
k2=1;
for k=1:(size(pion,1)-1)
    if(abs(pion(k,2)-pion(k+1,2))<2)
      prosta_pion(k1,k2,:)=[pion(k,1) pion(k,2)];
      k2=k2+1;
    else
      prosta_pion(k1,k2,:)=[pion(k,1) pion(k,2)];
      k1=k1+1; 
      k2=1;

    end
end

k1=1;
for k=1:size(prosta_dol,1)
    %if(prosta_dol(k,floor(size(prosta_dol,2)/3),1)>0)    
        if(prosta_dol(k,5,1)>0)
        for k2=1:size(prosta_dol,2)
            dol_fitowanie(k1,k2,:)=[prosta_dol(k,k2,1) prosta_dol(k,k2,2)];
       end
        k1=k1 + 1;  
    end
end

k1=1;
for k=1:size(prosta_gora,1)
    %if(prosta_gora(k,floor(size(prosta_gora,2)/3),1)>0)
    if(prosta_gora(k,5,1)>0)
    for k2=1:size(prosta_gora,2)
            gora_fitowanie(k1,k2,:)=[prosta_gora(k,k2,1) prosta_gora(k,k2,2)];
       end
        k1=k1 + 1;  
    end
end

k1=1;
for k=1:size(prosta_pion,1)
    %if(prosta_pion(k,floor(size(prosta_pion,2)/3),1)>0)    
    if(prosta_pion(k,5,1)>0)
    for k2=1:size(prosta_pion,2)
            pion_fitowanie(k1,k2,:)=[prosta_pion(k,k2,1) prosta_pion(k,k2,2)];
       end
        k1=k1 + 1;  
    end
end

for k=1:size(gora_fitowanie,1)
    for k2=1:size(gora_fitowanie,2)
        if(gora_fitowanie(k,k2,1)~=0 && gora_fitowanie(k,k2,2)~=0)
            dane_fit(k2,:)=[gora_fitowanie(k,k2,1) gora_fitowanie(k,k2,2)];
        end
    end
    wielomian_gora(k,:)=polyfit(dane_fit(:,1), dane_fit(:,2), 1);
    clear dane_fit;
end

for k=1:size(dol_fitowanie,1)
    for k2=1:size(dol_fitowanie,2)
        if(dol_fitowanie(k,k2,1)~=0 && dol_fitowanie(k,k2,2)~=0)
            dane_fit(k2,:)=[dol_fitowanie(k,k2,1) dol_fitowanie(k,k2,2)];
        end
    end
    wielomian_dol(k,:)=polyfit(dane_fit(:,1), dane_fit(:,2), 1); 
    clear dane_fit;
end

for k=1:size(pion_fitowanie,1)
    for k2=1:size(pion_fitowanie,2)
        
        if(pion_fitowanie(k,k2,1)~=0 && pion_fitowanie(k,k2,2)~=0)
            dane_fit(k2,:)=[pion_fitowanie(k,k2,1) pion_fitowanie(k,k2,2)];
        end
    end
    wielomian_pion(k,:)=polyfit(dane_fit(:,1), dane_fit(:,2), 0); 
    clear dane_fit;
end
%wielomian_pion=unique(wielomian_pion);
wielomian_pion=sort(wielomian_pion);
k2=1;
k=1;
while k < size(wielomian_pion,1)+1
    %disp(k)
    ilosc=1;
    srednia=wielomian_pion(k);
    for k1=k+1:size(wielomian_pion,1)
        if(abs(wielomian_pion(k)-wielomian_pion(k1))<10)
            srednia=srednia+wielomian_pion(k1);
            ilosc=ilosc+1;
        end
    end
    wielomian_pion_usredn(k2,:)=srednia/ilosc;
    k2=k2+1;
    k=k+ilosc;
end


k2=1;
k=1;
while k < size(wielomian_gora,1)+1
    %disp(k)
    ilosc=1;
    srednia=[wielomian_gora(k,1) wielomian_gora(k,2)];
    for k1=k+1:size(wielomian_gora,1)
        if(abs(wielomian_gora(k,1)-wielomian_gora(k1,1))<0.05 && abs(wielomian_gora(k,2)-wielomian_gora(k1,2))<15)
            srednia=srednia+[wielomian_gora(k1,1) wielomian_gora(k1,2)];
            ilosc=ilosc+1;
        end
    end
    wielomian_gora_usredn(k2,:)=srednia/ilosc;
    k2=k2+1;
    k=k+ilosc;
end

k2=1;
k=1;
while k < size(wielomian_dol,1)+1
    %disp(k)
    ilosc=1;
    srednia=[wielomian_dol(k,1) wielomian_dol(k,2)];
    for k1=k+1:size(wielomian_dol,1)
        if(abs(wielomian_dol(k,1)-wielomian_dol(k1,1))<0.2 && abs(wielomian_dol(k,2)-wielomian_dol(k1,2))<15)
            srednia=srednia+[wielomian_dol(k1,1) wielomian_dol(k1,2)];
            ilosc=ilosc+1;
        end
    end
    wielomian_dol_usredn(k2,:)=srednia/ilosc;
    k2=k2+1;
    k=k+ilosc;
end

hold off
imshow(obraz)
hold on
for k=1:size(wielomian_dol_usredn,1)
    x=1:443;
   plot(wielomian_dol_usredn(k,1)*x+wielomian_dol_usredn(k,2),x,'-r');
end

for k=1:size(wielomian_gora_usredn,1)
    x=1:443;
   plot(wielomian_gora_usredn(k,1)*x+wielomian_gora_usredn(k,2),x,'-b');
end

for k=1:size(wielomian_pion_usredn,1)
    x=1:443;
   plot(wielomian_pion_usredn(k),x,'-g');
end


dol_usredn=sort(wielomian_dol_usredn(:,2));

for k1=1:size(dol_usredn,1)
    for k2=1:size(dol_usredn,1)
        if(wielomian_dol_usredn(k2,2)==dol_usredn(k1))
         wielomian_dol(k1,:)=[wielomian_dol_usredn(k2,1), wielomian_dol_usredn(k2,2)];
        end
    end 
end

gora_usredn=sort(wielomian_gora_usredn(:,2));

for k1=1:size(gora_usredn,1)
    for k2=1:size(gora_usredn,1)
        if(wielomian_gora_usredn(k2,2)==gora_usredn(k1))
         wielomian_gora(k1,:)=[wielomian_gora_usredn(k2,1), wielomian_gora_usredn(k2,2)];
        end
    end 
end

point_left_up=[wielomian_pion_usredn(1), (wielomian_pion_usredn(1)- wielomian_dol(1,2))/wielomian_dol(1,1)];
point_left_down=[wielomian_pion_usredn(1), (wielomian_pion_usredn(1)- wielomian_gora(1,2))/wielomian_gora(1,1)];
point_right_down=[wielomian_pion_usredn(size(wielomian_pion_usredn,1)), (wielomian_pion_usredn(size(wielomian_pion_usredn,1))- wielomian_dol(size(wielomian_dol,1),2))/wielomian_dol(size(wielomian_dol,1),1)];
%point_right_down=[wielomian_pion_usredn(size(wielomian_pion_usredn)), (wielomian_pion_usredn(size(wielomian_pion_usredn))- wielomian_dol(size(wielomian_gora,1),2))/wielomian_gora(size(wielomian_gora,1),1)];
point_right_up=[wielomian_pion_usredn(size(wielomian_pion_usredn,1)), (wielomian_pion_usredn(size(wielomian_pion_usredn,1))- wielomian_gora(size(wielomian_gora,1),2))/wielomian_gora(size(wielomian_gora,1),1)];

point_up=[(wielomian_dol(1,2)-wielomian_gora(size(wielomian_gora,1),2))/(wielomian_gora(size(wielomian_gora,1),1)-wielomian_dol(1,1))*wielomian_dol(1,1)+wielomian_dol(1,2), (wielomian_dol(1,2)-wielomian_gora(size(wielomian_gora,1),2))/(wielomian_gora(size(wielomian_gora,1),1)-wielomian_dol(1,1))];

%point_down=[(wielomian_dol(size(wielomian_dol,1),2)-wielomian_gora(1,2))/(wielomian_gora(size(1,1)-wielomian_dol(wielomian_dol(size(wielomian_dol,1),1))*wielomian_dol(wielomian_dol(size(wielomian_dol,1),1)+wielomian_dol(wielomian_dol(size(wielomian_dol,1),2), (wielomian_dol(wielomian_dol(size(wielomian_dol,1),2)-wielomian_gora(1,2))/(wielomian_gora(1,1)-wielomian_dol(wielomian_dol(size(wielomian_dol,1),1))];

point_down=[(wielomian_dol(size(wielomian_dol,1),2)- wielomian_gora(1,2) ) / ( wielomian_gora(1,1)- wielomian_dol(size(wielomian_dol,1),1) ) * wielomian_dol(size(wielomian_dol,1),1) + wielomian_dol(size(wielomian_dol,1),2) , ( wielomian_dol(size(wielomian_dol,1),2) - wielomian_gora(1,2) )/( wielomian_gora(1,1)- wielomian_dol(size(wielomian_dol,1),1) )];

plot(point_left_up(1),point_left_up(2),'om');
plot(point_left_down(1),point_left_down(2),'om');
plot(point_right_up(1),point_right_up(2),'om');
plot(point_right_down(1),point_right_down(2),'om');
plot(point_up(1),point_up(2),'om');
plot(point_down(1),point_down(2),'om');
hold off
volume1=[sqrt((point_right_down(1)-point_down(1))^2+(point_down(2)-point_right_down(2))^2),sqrt((point_down(1)-point_left_down(1))^2+(point_down(2)-point_left_down(2))^2), sqrt((point_right_up(1)-point_right_down(1))^2+(point_right_down(2)-point_right_up(2))^2) ];
volume2=[sqrt((point_up(1)-point_left_up(1))^2+(point_left_up(2)-point_up(2))^2),sqrt((point_right_up(1)-point_up(1))^2+(point_up(2)-point_right_up(2))^2), sqrt((point_left_up(1)-point_left_down(1))^2+(point_left_down(2)-point_left_up(2))^2) ];
objetosc=[(volume1(1)+volume2(1))/2,(volume1(2)+volume2(2))/2,(volume1(3)+volume2(3))/2];
figure(2)
%hold on


%plot3([0, 0, 0, 0, objetosc(1), objetosc(1), objetosc(1), objetosc(1)],[0, 0, objetosc(2), objetosc(2), objetosc(2), objetosc(2), 0, 0 ],[0, objetosc(3), objetosc(3), 0, 0, objetosc(3), objetosc(3), 0 ], 'or');
%hold off
%x=1;
%y=1;
%for k1=-20:5:objetosc(1)+20
%    for k2=-20:5:objetosc(2)+20
%        if((k1>0||k1<objetosc(1))&&(k2>0||k2<objetosc(2)))
%            Z(x,y)=objetosc(3);
%        else
%            Z(x,y)=0;
%        end
%        y=y+1;
%    end
%    x=x+1;
%end
%X = -20:5:objetosc(1)+20;
%Y = -20:5:objetosc(2)+20)
%
%[X Y Z]=meshgrid(-20:10:objetosc(1)+20,-20:10:objetosc(2)+20,-20:10:objetosc(3)+20);
x=-1:1:objetosc(1)+1;
y=-1:1:objetosc(2)+1;
%z=-20:5:objetosc(3)+20;

%,-20:10:objetosc(3)+20);
%Z=zeros(size(x,2),size(y,2));
for k1=1:size(x,2)
    for k2=1:size(y,2)
        %for k3=1:size(X,3)
    
        if(x(k1)>0 && x(k1)<objetosc(1) && y(k2)>0 && y(k2)<objetosc(2))            
%        if((X(k1,k2,k3)>0||X(k1,k2,k3)<objetosc(1))&&(Y(k1,k2,k3)>0||Y(k1,k2,k3)<objetosc(2))      
                Z(k2,k1)=objetosc(3);
        elseif((x(k1)==0||x(k1)==objetosc(1))&&(y(k2)==0||y(k2)==objetosc(2))) 
        z(k2,k1)=objetosc(3)/2;
        else    
                Z(k2,k1)=0;
        
            end
        
        end
    %end
    
    
end
[X Y]=meshgrid(x,y);
surf(X,Y,Z);





%objetosc(3)+20)


%[0, 0,0, 0, objetosc(1), objetosc(1), objetosc(1), objetosc(1) ; 0, 0, objetosc(2), objetosc(2), objetosc(2), objetosc(2), 0, 0 ; 0, objetosc(3), objetosc(3), 0, 0, objetosc(3), objetosc(3), 0 ]
%[0 0 0;
%0 0 volume(3);
%0  volume(2) 0;
%0 volume(2) volume(3);
%volume(1) 0 0;
%volume(1) 0 volume(3);
%volume(1) volume(2) 0;
%volume(1) volume(2) volume(3)])

%[volume(3) volume(3)]
%[volume(3) volume(3)]
%[volume(3) volume(3)]
%[0 volume(3)]
%[0 volume(3)]
%[0 volume(3)]
%[0 volume(3)]

%point_up=(wielomian_dol(1,1)-(wielomian_gora(size(wielomian_gora,1),1))/(wielomian_gora(size(wielomian_gora,1),2)-wielomian_dol(1,2))*wielomian_dol(1,2)+wielomian_dol(1,1),(wielomian_dol(1,1)-wielomian_gora(size(wielomian_gora),1))/(wielomian_gora(size(wielomian_gora,1),2)-wielomian_dol(1,2))];
%point_down=
%point_right_up=[wielomian_pion_usredn(1), (wielomian_pion_usredn(1)- wielomian_gora(1,2))/wielomian_gora(1,1)]

%point_left_down=(wielomian_pion_usredn(1), (wielomian_pion_usredn(1)- wielomian_dol(1,2))/wielomian_dol(1,1))


%for k=1:size(prosta_dol,1)
%    if(prosta_dol(k,floor(size(prosta_dol,2)/3),1)>0)    
%        for k2=1:size(prosta_dol,2)
%            if(prosta_dol(k,k2,1)~=0 && prosta_dol(k,k2,2)~=0)
%                dane_fit(k2,:)=[prosta_dol(k,k2,1) prosta_dol(k,k2,2)]
%            end
%        end
%        wielomian_dol(k,:)=polyfit(dane_fit(:,2), dane_fit(:,1), 1); 
%    end
%end

%k1=1;
%for k=1:size(prosta_pion,1)
%    if(prosta_pion(k,floor(size(prosta_pion,2)/3),1)>0)    
%        for k2=1:size(prosta_pion,2)
%            pion_fitowanie(k1,k2,:)=[prosta_pion(k,k2,1) prosta_pion(k,k2,2)];
%       end
%        k1=k1 + 1;  
%    end
%end
%   gora=pochodnex(:,3,:)==lambda1;
%dol=pochodnex(:,3,:)lambda2;

%k=size(pochodney)
%if(k(1)>0)
%lambda3=0;
%else
%lambda3=mean(wspolczynnikix);
%end 
%B1 =imrotate(BW7,atan(lambda1)*180/pi);
%B2 =imrotate(BW7,atan(lambda2)*180/pi);
%B3 =imrotate(BW7,atan(lambda3)*180/pi);
%B1 = atan(lambda1)*180/pi;
%B2 = atan(lambda2)*180/pi;
%B3 = atan(lambda3)*180/pi;
%wspolczynnikiy=unique(pochodney(:,3))
%k1=1;
%k2=1;
%k3=1;
%rozmiar3=size(pochodnex);
%   for l=1:rozmiar3(1)
%       if (pochodnex(l,3)==lambda1)
%           przesuniecia1(k1,:)=[pochodnex(l,1), pochodnex(l,2), (pochodnex(l,1)-lambda1*pochodne(l,2))];
%       k1=k1+1;
%       elseif(pochodnex(l,3)==lambda2)
%           przesuniecia2(k2,:)=[pochodnex(l,1), pochodnex(l,2), (pochodnex(l,1)-lambda2*pochodne(l,2))];
%       k2=k2+1;  
%       elseif(pochodne(l,3)>=lambda3-0.00000000000000001 || pochodne(l,3)<=lambda3+0.00000000000000001 )
%           przesuniecia3(k3,:)=[pochodnex(l,1), pochodnex(l,2), (pochodnex(l,1)-lambda3*pochodne(l,2))];
%       k3=k3+1;
%       end
%   end
%   
%       disp(l);
%przesuniecia1=(pochodnex(3)=lambda1);
%przesuniecia2=(pochodnex(3)=lambda2);
%przesuniecia3=(pochodnex(3)=lambda2);

%pion1=unique(przesuniecia1(:,3));
%pion2=unique(przesuniecia2(:,3));
%pion3=unique(przesuniecia3(:,3));


%theta = 0:359;
%[R,xp] = radon(BW7,theta);
%imagesc(theta, xp, R); colormap(hot);
%xlabel('\theta (degrees)'); ylabel('x\prime');
%title('R_{\theta} (x\prime)');
%colorbar

%polyfit(boundary1(:,2), boundary1(:,1), 1);
