There some problem in the code of random packing of circles [closed]

this matlab code is not doing its job of drawing non overlapping circles according to the algo. i have taken this code from https://github.com/BluesBlues213/random-close-packing but there circles are always drawn in a square but I was trying this code to draw circles in a rectangle but overlapping is present in the output.please help where i am wrong.

function [NWP NWR] = package(pool)
mm=mean(pool);
rr=numel(pool);
pool(rr);
area1=(sum(pi*(pool.^2))-(pi*(pool(rr).^2)));
div=mm/(10^-5);
pool=pool/div;
area=sum(pi*(pool.^2))/0.9;
side=sqrt(area);
o=sqrt(area/6);
width=2*o;
length=3*o;
side1=2*sqrt(area/6);
side2=3*sqrt(area/6);
res=1400; 
% ST=side/res;
ST1=side1/res;
ST2=side2/res;
[x,y]=meshgrid(0:ST1:side1, 0:ST2:side2);
p0=[side1/2,side2/2];
r0=pool(1);
pool=pool(2:numel(pool));
Wp(1,:)=p0;
Wr=r0;
D0=zeros(size(x,1),size(y,1));
D1=zeros(size(x,1),size(y,1));
m_M=zeros(size(x,1),size(y,1));
R1=zeros(size(x,1),size(y,1));
i=1;
j=1;
MASK=zeros(size(x));
px=Wp(:,1);
py=Wp(:,2);
down=0.0*10^-4;
% up=side;
% t=linspace(down,up,res+1)
% SC=t(5)-t(4);
% [x,y]=meshgrid(down:SC:up , down:SC:up);
up1=side1;
up2=side2;
 t1=linspace(down,up1,res+1);
 t2=linspace(down,up2,res+1);
SC1=t1(5)-t1(4);
SC2=t2(5)-t2(4);
[x,y]=meshgrid(down:SC1:up1 , down:SC2:up2);
routter=Wr.*(1)
W=zeros(size(x));
for i=1:numel(Wr)
    tic;
    R = sqrt( (x-px(i)).^2 + (y-py(i)).^2 );
    ONES=zeros(size(x));
    ONES(R<routter(i))=inf;
    W=W+ONES;
    toc   
end
MASK11=W;
while i<=numel(pool)
        tic;
    WPT=vertcat(Wp,Wp);
    WR=vertcat(Wr,Wr);
    rnext=pool(i);
        avg_center=[mean(WPT(:,1)) mean(WPT(:,2))];
        d1=sqrt( (x-avg_center(1)).^2 + (y-avg_center(2)).^2 );
    d_t=d1;

      m_M=zeros(size(x));
    s=size(m_M)
    for n=1:numel(WR)
        w=zeros(size(x,1),size(y,1));
        w=sqrt( (x-WPT(n,1)).^2 + (y-WPT(n,2)).^2 ) <= ...
            WR(n)+rnext;
                      m_M=m_M+w;
    end
    m_M=m_M+MASK11;
    d_t(m_M>0)=1;
    [mx my]=find(d_t==min(min(d_t)));
    new_pt=[y(my(1),1) x(1,mx(1))];

    Wp=[Wp;new_pt];
    Wr=vertcat(Wr,rnext);

    i=i+1;
    display(['iteration #' num2str(i)])
    etoc(:,i)=toc;


    Wpc=Wp;
    while numel(Wp(:,1))>2
        mask=poly2mask(round(Wp(:,1)./ST1),round(Wp(:,2)./ST2),res+1,res+1);
        r=size(mask)
        MASK11=MASK11+mask;
        if i>2
            break
        end

    end
    toc
end
idx= ~(WPT(:,1)==0 & WPT(:,2)==0);
NWP=WPT(idx,:);
NWR=WR(idx);
NWP=NWP*div;
NWR=NWR*div;
% close all;
figure;
 axis equal; 
title('packing');
hold on
for i=1:numel(NWR);
    [H,areao,a,b]=circle(NWP(i,:),NWR(i));
end
 px=NWP(:,1);
 py=NWP(:,2);
 cx=mean(px);
 cy=mean(py);
dr=sqrt((px-cx).^2 + (py-cy).^2);
%  figure;
q=area1/areao
width=a
length=b
% rectangle('position', [0,0,side, side])
% plot(dr,NWR,'o');axis on;axis square;xlabel('distance to center');ylabel('radii size');
end