function papartis
% papartis Fraktalas, pagal M.F. Barnsley, L.P. Hurd,
% fractal Image Compression, A.K. Peters Ltd., 1993
% Papartis apibreziamas 24 skaiciais:
% w1: (x,y) ->(0.85*x+0.04*y,-0.04*x+0.85*y+0.40),
% w2: (x,y) ->(0.20*x-0.26*y,0.23*x+0.22*y+0.40),
% w3: (x,y) ->(-0.15*x+0.28*y,0.26*x+0.24*y+0.11),
% w4: (x,y) ->(0,-0.16*y),

N=64*16;
b=[1.6 1.6 0.44 0]; 
w{1}=[0.85  0.04; -0.04 0.85]; 
w{2}=[0.20 -0.26;  0.23 0.22];
w{3}=[-0.15 0.28; 0.28 0.24];
w{4}=[0 0; 0 0.16];

c='gbry';

clf; hold on; 
C=ones(1,N); x=zeros(2,N);
for l=1:16,
 disp(l)
 ii=randperm(N); ic=randperm(N);
 for k=2:N,
  i=mod(ii(k),4)+1;
  C(k)=c(mod(ic(k),4)+1);
  x(:,k)=w{i}*x(:,k-1)+[0;b(i)];
 end;
 for j=1:4,
  ii=find(C==c(j));
  plot(x(1,ii),x(2,ii),[c(j) '.'],'markersize',3); 
 end;
 C(1)=C(end); x(:,1)=x(:,end);
end;
axis equal; axis off;
