function out=tree(n,p,A,b,cc,drawall)
% Growing a fractal tree. Self-demonstrating, or tree(n,p,A,b,cc,drawall).
% Calls the modifications pplot and ffill of plot and fill.
%
% Given an initial polygon p=[p_1,...p_n] in the plane, a collec-
% tion of m 2x2 matrices A_i and 2x1 vectors b_i, this procedure 
% iterates the map   pnew_i = A_i * p + b_i   n-times and overlays
% the resulting m^n polygons.
% The original purpose was to sketch the attractors for such
% multivalued affine maps; A should be a contraction (the abso-
% lute value of its largest eigenvalue should be less than 1), 
% and the collection of new polygons pnew_1, ... pnew_m should
% closely resemble the outline of the original polygon.
% This procedure was motivated by the article "Fractals in Linear
% Algebra" by J.Walsh in the College Math Journal, Sept.1996.
%
% n is the number of iterations to be carried out (maximal 6)
% p is a 2xn matrix of the coordinates of the initial polygon
% A is a 2*2m matrix, holding m 2x2 matrices   A_1, ... A_m.
% b is a 2*m matrix,  containing m 2*1 vectors b_1, ... b_m.
% cc is a plot/color-style passed through to the plot command.
% drawall determines whether all stages are drawn or only the 
%         last iteration(in the case of drawall==0)
% If no arguments or only the number of iterations n is supplied,
% the procedure will use a default set of points and matrices.
%
% Author: Matt Kawski, http://math.la.asu.edu/~kawski Sept.1996
%
if nargin<=4, cc='g-';   end      % set default color to green
if nargin<=5, drawall=0; end      % default draws only last iteration
if nargin<=1                      % selfdemonstrating: example data
   if nargin==0, n=3;             % default 3 iterations
      else, n=min([n,4]);         % but not more than 4 for selfdemo
      end;
   p=[-.15,-0.10,-1.14,-1.69,-1.74,-1.42,-0.89,-1.51,-0.98,.03,.70,...
        1.39,1.67,1.14,.82,1.2,2.02,2.16,1.63,1.09,.56,0.40,.41,.4,-.15;...
         -2.45,.07,.25,.81,1.54,2.01,2.23,2.8,3.48,3.83,3.8,3.74,3.39,...
        2.56,2.24,1.86,1.54,0.92,.28,.1,.07,-.78,-1.18,-2.27,-2.44];
   A=[.25,-.43,.25,.43,.29,.17,.27,-0.2,.6,0;...
     .43,.25,-.43,.25,-.17,.29,.20,.27,0,.8];
   b=[0,.3,.7,-.2,.02;.6,.7,2.7,2.3,0];
   end;
n=min([n,6]);                      % number of iterations is maximal 6
sa=size(A);
sb=size(b);                           
nn=min([sa(1,2),sb(1,2)]);         % number of "children" to be drawn
hold on
if n==0, ffill(p,cc)               % after last iteration draw result
   else
   if drawall, ffill(p,cc), end    % draw all iterations if drawall=1
   while nn>0
   tree(n-1,A(:,2*nn-1:2*nn)*p+b(:,nn)*ones(size(p(1,:))),A,b,cc,drawall)
   nn=nn-1;                        %counter for level of iteration
   end
end;