function S=sier2(A,p1,p2,b11,b22)
%Sierpinski triangle.  -JJL
if nargin==0,A=[0 0]';end
if nargin<2,p1=.5;end   %x-scaling
if nargin<3,p2=p1;end   %y-scaling
if nargin<4,b11=.5;end   %x-translation
if nargin<5,b22=b11;end   %xy-translation
[M,N]=size(A);
if M~=2,error('Incorrect data format'),end;
b1=repmat([b11,0]',1,N);
b2=repmat([0,b22]',1,N);
y1=[p1,0;0,p2]*A;
y2=y1+b1;
y3=y1+b2;
x=[y1(1,:),y2(1,:),y3(1,:)];y=[y1(2,:),y2(2,:),y3(2,:)];
plot(x,y,'.')
if nargout>0,S=[x;y];end