clear

% Initial Triangle
%v=[0 0; 1 0; 0.5 sqrt(3)/2 ];
v=[0 0; 1 0; 0 1];
N=size(v,1);
z=v;
patch(v(:,1), v(:,2),[1 0 0]);
axis([0 1 0 1])
hold

% Linear Transformations
A=[0.5 0; 0 0.5];

% Translations
b1=[0 0];
b2=[0.5 0];
%b3=[0.25 sqrt(3)/4];
b3=[0 0.5];

for i=1:8
    M=size(z,1)/N;   % Number of triangles
    for j=1:M
        v=z(3*(j-1)+1:3*j,:);
        % Triangle 1
        w=A*v';
        w=w'+repmat(b1,3,1);
        if j==1
            z2=w;
        else
            z2=[z2;w];
        end
        patch(w(:,1), w(:,2),[1-0.1*i 0.1*i 0.1*i]);
        
        % Triangle 2
        w=A*v';
        w=w'+repmat(b2,3,1);
        z2=[z2;w];
        patch(w(:,1), w(:,2),[1-0.1*i 0.1*i 0.1*i]);
        
        % Triangle 3
        w=A*v';
        w=w'+repmat(b3,3,1);
        z2=[z2;w];
        patch(w(:,1), w(:,2),[1-0.1*i 0.1*i 0.1*i]);
    end
    z=z2;
end
    
    hold