% Copyright: Wing Kam Liu, Adrian Kopacz, Eduard Karpov,Harold Park,
%            David Farrell, Mechanical Engineering,
%            Northwestern University, 2004 


function md2d_sim(hx,hy,sL,nIt,mI,tR,fN,fps);

% input
% positions (x,y)
% space dimension (sL)
% number of iterations (nIt)
% number of mirror images (mI)
% trace option (tR)
% movie output, frames per second (fps)

% output
% none

% set up vector plot
[nP,n]=size(hx);
sVx=zeros(nP,(1+2*mI));
sVy=zeros(nP,(1+2*mI));

sVxx=[];
sVyy=[];

for it=1:nP
    sVx(it,:)=hx(it,1):sL:((2*mI*sL)+hx(it,1));
    sVy(it,:)=hy(it,1):sL:((2*mI*sL)+hy(it,1));
end

sT=0;
for i=1:(1+2*mI)
    for j=1:(1+2*mI)
        sT=sT+1;
        sVxx=[sVxx; sVx(:,i)];
        sVyy=[sVyy; sVy(:,j)];         
    end
end

% set up plot
set(gcf,'color',[1 1 1]);
h=plot(sVxx,sVyy,'.k');
axis([0 (sL+2*mI*sL) 0 (sL+2*mI*sL)]); 
axis off;
set(h, 'EraseMode', 'xor', 'MarkerSize',10); 

% draw outer boundaries
line([0,0],[0,(sL+2*mI*sL)],'Color', [1,0,0],'LineWidth', 2);
line([(sL+2*mI*sL),0],[0,0],'Color', [1,0,0],'LineWidth', 2);
line([0,(sL+2*mI*sL)],[(sL+2*mI*sL),(sL+2*mI*sL)],'Color', [1,0,0],'LineWidth', 2);
line([(sL+2*mI*sL),(sL+2*mI*sL)],[(sL+2*mI*sL),0],'Color', [1,0,0],'LineWidth', 2);

% draw image boundaries
for it=1:(2*mI)
    
    % vertical line
    line([sL*it,sL*it],[0,(sL+2*mI*sL)],...
        'Color', [0,1,0],'LineWidth', 1);
    
    % horizontal line
    line([0,(sL+2*mI*sL)],[sL*it,sL*it],...
        'Color', [0,1,0],'LineWidth', 1);
end;


for k=1:nIt  
    
    % set up vector plot
    if (tR == 0)
        sVxx=[];
        sVyy=[];
    end

    for it=1:nP
        sVx(it,:)=hx(it,k):sL:((2*mI*sL)+hx(it,k));
        sVy(it,:)=hy(it,k):sL:((2*mI*sL)+hy(it,k));
    end
    
    sT=0;
    for i=1:(1+2*mI)
        for j=1:(1+2*mI)
            sT=sT+1;
            sVxx=[sVxx; sVx(:,i)];
            sVyy=[sVyy; sVy(:,j)];       
        end
    end
    
    % plot
    set(h,'XData',sVxx,'YData',sVyy); 
    drawnow; 
end;
