clear all
xm=100;
ym=100;
sink.x=0.5*xm;  %location of sink on x-axis
sink.y=0.5*ym;  %location of sink on y-axis
n=100  %nodes+
P=0.1;  %probability of cluster heads
Eo=0.5;%initial energy
%
Echeck=Eo;
%
ETX=50*0.000000001;  %tx energy
ERX=50*0.000000001;  %rx energy
Efs=10*0.000000000001;  %free space loss
Emp=0.0013*0.000000000001;   %multipath loss
%Data Aggregation Energy
EDA=5*0.000000001;  %compression energy
a=1.5;   %fraction of energy enhancment of advance nodes
rmax=10000  %maximum number of rounds
do=sqrt(Efs/Emp);  %distance do is measured
Et=0;  %variable just use below 
m=0.5;
mo=0.4;
b=3;%fraction of extra energy in super nodes
normal=n*(1-m);
advance=n*m*(1-mo);
super=n*m*mo;






m7=0;
c=0.02;
mony=0;
A=0;
for i=1:1:n
    S1(i).xd=rand(1,1)*xm;
    S2(i).xd=S1(i).xd;
    S3(i).xd=S1(i).xd;
    S4(i).xd=S1(i).xd;
    XR4(i)=S4(i).xd;
    XR3(i)=S3(i).xd;
    XR2(i)=S2(i).xd;
    XR1(i)=S1(i).xd;
    S1(i).yd=rand(1,1)*ym;
    S2(i).yd=S1(i).yd;
    S3(i).yd=S1(i).yd;
    S4(i).yd=S1(i).yd;
    YR4(i)=S4(i).yd;
    S4(i).G=0;
    YR3(i)=S3(i).yd;
    S3(i).G=0;
    YR2(i)=S2(i).yd;
    YR1(i)=S1(i).yd;
    S1(i).G=0;
    S2(i).G=0;
    talha=rand*a;
    S1(i).E=Eo*(1+talha);
    S2(i).E=S1(i).E;
    S3(i).E=S1(i).E;
    S4(i).E=S3(i).E;
    E3(i)= S3(i).E;
    E4(i)= S4(i).E;
    Et=Et+E3(i);
    E(i)= S1(i).E;
    A=A+talha;   % for deec algo probabilty function
    
    
    S1(i).A=talha;
    S2(i).A=talha;
    if (S2(i).E>Echeck)
        mony=mony+1;
    end
    
    Et=Et+E3(i);

    %initially there are no cluster heads only nodes
    S1(i).type='N';
    S2(i).type='N';
    S3(i).type='N';
    S4(i).type='N';
end

d1=0.765*xm/2;  %distance between cluster head and base station
K=sqrt(0.5*n*do/pi)*xm/d1^2; %optimal no. of cluster heads
d2=xm/sqrt(2*pi*K);  %distance between cluster members and cluster head
Er=4000*(2*n*ETX+n*EDA+K*Emp*d1^4+n*Efs*d2^2);  %energy desipated in a round
S1(n+1).xd=sink.x; %sink is a n+1 node, x-axis postion of a node
S1(n+1).yd=sink.y; %sink is a n+1 node, y-axis postion of a node
S2(n+1).xd=sink.x; %sink is a n+1 node, x-axis postion of a node
S2(n+1).yd=sink.y; %sink is a n+1 node, y-axis postion of a node
S3(n+1).xd=sink.x; %sink is a n+1 node, x-axis postion of a node
S3(n+1).yd=sink.y; %sink is a n+1 node, y-axis postion of a node
S4(n+1).xd=sink.x; %sink is a n+1 node, x-axis postion of a node
S4(n+1).yd=sink.y; %sink is a n+1 node, y-axis postion of a node
countCHs1=0;  %variable, counts the cluster head
countCHs2=0;  %variable, counts the cluster head
countCHs3=0;  %variable, counts the cluster head
countCHs4=0;  %variable, counts the cluster head
cluster1=1;  %cluster is initialized as 1
cluster2=1;  %cluster is initialized as 1
cluster3=1;  %cluster is initialized as 1
cluster4=1;  %cluster is initialized as 1
flag_first_dead1=0; %flag tells the first node dead
flag_first_dead2=0; %flag tells the first node dead
flag_first_dead3=0; %flag tells the first node dead
flag_first_dead4=0; %flag tells the first node dead
flag_teenth_dead1=0;  %flag tells the 10th node dead
flag_teenth_dead2=0;  %flag tells the 10th node dead
flag_teenth_dead3=0;  %flag tells the 10th node dead
flag_teenth_dead4=0;  %flag tells the 10th node dead
flag_all_dead1=0;  %flag tells all nodes dead
flag_all_dead2=0;  %flag tells all nodes dead
flag_all_dead3=0;  %flag tells all nodes dead
flag_all_dead4=0;  %flag tells all nodes dead
dead1=0;  %dead nodes count initialized to 0
dead2=0;  %dead nodes count initialized to 0
dead3=0;  %dead nodes count initialized to 0
dead4=0;  %dead nodes count initialized to 0
first_dead1=0;
first_dead2=0;
first_dead3=0;
first_dead4=0;
teenth_dead1=0;
teenth_dead2=0;
teenth_dead3=0;
teenth_dead4=0;
all_dead1=0;
all_dead2=0;
all_dead3=0;
all_dead4=0;
allive1=n;
allive2=n;
allive3=n;
allive4=n;
%counter for bit transmitted to Bases Station and to Cluster Heads
packets_TO_BS1=0;
packets_TO_BS2=0;
packets_TO_BS3=0;
packets_TO_BS4=0;
packets_TO_CH1=0;
packets_TO_CH2=0;
packets_TO_CH3=0;
packets_TO_CH4=0;



%Start of DEEC Algo




for r=0:1:rmax     
    r
  if(mod(r, round(1/P) )==0)
    for i=1:1:n
        S1(i).G=0;
        S1(i).cl=0;
    end
  end
Ea=Et*(1-r/rmax)/n;
dead1=0;
for i=1:1:n
   
    if (S1(i).E<=0)
        
        dead1=dead1+1; 
        if (dead1==1)
           if(flag_first_dead1==0)
              first_dead1=r;
              flag_first_dead1=1;
           end
        end
        if(dead1==0.1*n)
           if(flag_teenth_dead1==0)
              teenth_dead1=r;
              flag_teenth_dead1=1;
           end
        end
        if(dead1==n)
           if(flag_all_dead1==0)
              all_dead1=r;
              flag_all_dead1=1;
           end
        end
    end
    if S1(i).E>0
        S1(i).type='N';
    end
end
STATISTICS.DEAD1(r+1)=dead1;
STATISTICS.ALLIVE1(r+1)=allive1-dead1;
countCHs1=0;
cluster1=1;
for i=1:1:n
 if Ea>0
    
     
   p(i)=P*n*(1+S1(i).A)*E(i)/(n+A)*(Ea); 
   %p(i)=P*n*S(i).E*E(i)/(Et*Ea);
 
 
 if(S1(i).E>0)
   temp_rand=rand;     
   if ( (S1(i).G)<=0)  
        if(temp_rand<= (p(i)/(1-p(i)*mod(r,round(1/p(i))))))
            countCHs1=countCHs1+1;
            packets_TO_BS1=packets_TO_BS1+1;
            PACKETS_TO_BS1(r+1)=packets_TO_BS1;
             S1(i).type='C';
            S1(i).G=round(1/p(i))-1;
            C1(cluster1).xd=S1(i).xd;
            C1(cluster1).yd=S1(i).yd;
           distance=sqrt( (S1(i).xd-(S1(n+1).xd) )^2 + (S1(i).yd-(S1(n+1).yd) )^2 );
            C1(cluster1).distance=distance;
            C1(cluster1).id=i;
            X1(cluster1)=S1(i).xd;
            Y1(cluster1)=S1(i).yd;
            cluster1=cluster1+1;
           distance;
            if (distance>do)
                S1(i).E=S1(i).E- ( (ETX+EDA)*(4000) + Emp*4000*( distance*distance*distance*distance )); 
            end
            if (distance<=do)
                S1(i).E=S1(i).E- ( (ETX+EDA)*(4000)  + Efs*4000*( distance * distance )); 
            end
        end     
    
   end
   
 end 
 end
end
STATISTICS.COUNTCHS1(r+1)=countCHs1;

for i=1:1:n
   if ( S1(i).type=='N' && S1(i).E>0 )
     if(cluster1-1>=1)
       min_dis=sqrt( (S1(i).xd-S1(n+1).xd)^2 + (S1(i).yd-S1(n+1).yd)^2 );
       min_dis_cluster=0;
       for c=1:1:cluster1-1
           temp=min(min_dis,sqrt( (S1(i).xd-C1(c).xd)^2 + (S1(i).yd-C1(c).yd)^2 ) );
           if ( temp<min_dis )
               min_dis=temp;
               min_dis_cluster=c;
           end
       end
       
       if(min_dis_cluster~=0)    
            min_dis;
            if (min_dis>do)
                S1(i).E=S1(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S1(i).E=S1(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
      
            S1(C1(min_dis_cluster).id).E = S1(C1(min_dis_cluster).id).E- ( (ERX + EDA)*4000 ); 
            packets_TO_CH1=packets_TO_CH1+1;
       else 
            min_dis;
            if (min_dis>do)
                S1(i).E=S1(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S1(i).E=S1(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS1=packets_TO_BS1+1;
            
       end
        S1(i).min_dis=min_dis;
       S1(i).min_dis_cluster=min_dis_cluster;
   else
            min_dis=sqrt( (S1(i).xd-S1(n+1).xd)^2 + (S1(i).yd-S1(n+1).yd)^2 );
            if (min_dis>do)
                S1(i).E=S1(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S1(i).E=S1(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS1=packets_TO_BS1+1;
   end
  end
end
STATISTICS.PACKETS_TO_CH1(r+1)=packets_TO_CH1;
STATISTICS.PACKETS_TO_BS1(r+1)=packets_TO_BS1;
end
first_dead1
teenth_dead1
all_dead1
STATISTICS.DEAD1(r+1)
STATISTICS.ALLIVE1(r+1)
STATISTICS.PACKETS_TO_CH1(r+1)
STATISTICS.PACKETS_TO_BS1(r+1)
STATISTICS.COUNTCHS1(r+1)








%Start of DDEEC algo









m7=mony/100;
for r=0:1:rmax     
    r
  if(mod(r, round(1/P) )==0)
    for i=1:1:n
        S2(i).G=0;
        S2(i).cl=0;
    end
  end
 Ea=Et*(1-r/rmax)/n;
 dead2=0;
for i=1:1:n
   
    if (S2(i).E<=0)
        dead2=dead2+1; 
        if (dead2==1)
           if(flag_first_dead2==0)
              first_dead2=r;
              flag_first_dead2=1;
           end
        end
        if(dead2==0.1*n)
           if(flag_teenth_dead2==0)
              teenth_dead2=r;
              flag_teenth_dead2=1;
           end
        end
        if(dead2==n)
           if(flag_all_dead2==0)
              all_dead2=r;
              flag_all_dead2=1;
           end
        end
    end
    if S2(i).E>0
        S2(i).type='N';
    end
end

STATISTICS.DEAD2(r+1)=dead2;
STATISTICS.ALLIVE2(r+1)=allive2-dead2;
countCHs2=0;
cluster2=1;
for i=1:1:n
    
   
 if Ea>0
     
     if (S2(i).E<= (0.7*Eo))
         p(i)=c*(1+S2(i).A)*P*E(i)/((1+(S2(i).A*m7))*Ea);
     end
         
     if (S2(i).E<=Eo)
         p(i)=P*E(i)/((1+(S2(i).A*m7))*Ea);
     end
     if (S2(i).E>Eo)
         p(i)=(1+S2(i).A)*P*E(i)/((1+(S2(i).A*m7))*Ea);
     end
         
 
 if(S2(i).E>0)
   temp_rand=rand;     
   if ( (S2(i).G)<=0)  
        if(temp_rand<= (p(i)/(1-p(i)*mod(r,round(1/p(i))))))
            countCHs2=countCHs2+1;
            packets_TO_BS2=packets_TO_BS2+1;
            PACKETS_TO_BS2(r+1)=packets_TO_BS2;
             S2(i).type='C';
            S2(i).G=round(1/p(i))-1;
            C2(cluster2).xd=S2(i).xd;
            C2(cluster2).yd=S2(i).yd;
           distance=sqrt( (S2(i).xd-(S2(n+1).xd) )^2 + (S2(i).yd-(S2(n+1).yd) )^2 );
            C2(cluster2).distance=distance;
            C2(cluster2).id=i;
            X2(cluster2)=S2(i).xd;
            Y2(cluster2)=S2(i).yd;
            cluster2=cluster2+1;
           distance;
            if (distance>do)
                S2(i).E=S2(i).E- ( (ETX+EDA)*(4000) + Emp*4000*( distance*distance*distance*distance )); 
            end
            if (distance<=do)
                S2(i).E=S2(i).E- ( (ETX+EDA)*(4000)  + Efs*4000*( distance * distance )); 
            end
        end     
    
   end
   
 end 
 end
 end
STATISTICS.COUNTCHS2(r+1)=countCHs2;

for i=1:1:n
   if ( S2(i).type=='N' && S2(i).E>0 )
     if(cluster2-1>=1)
       min_dis=sqrt( (S2(i).xd-S2(n+1).xd)^2 + (S2(i).yd-S2(n+1).yd)^2 );
       min_dis_cluster=0;
       for c=1:1:cluster2-1
           temp=min(min_dis,sqrt( (S2(i).xd-C2(c).xd)^2 + (S2(i).yd-C2(c).yd)^2 ) );
           if ( temp<min_dis )
               min_dis=temp;
               min_dis_cluster=c;
           end
       end
       
       if(min_dis_cluster~=0)    
            min_dis;
            if (min_dis>do)
                S2(i).E=S2(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S2(i).E=S2(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
      
            S2(C2(min_dis_cluster).id).E = S2(C2(min_dis_cluster).id).E- ( (ERX + EDA)*4000 ); 
            packets_TO_CH2=packets_TO_CH2+1;
       else 
            min_dis;
            if (min_dis>do)
                S2(i).E=S2(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S2(i).E=S2(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS2=packets_TO_BS2+1;
            
       end
        S2(i).min_dis=min_dis;
       S2(i).min_dis_cluster=min_dis_cluster;
   else
            min_dis=sqrt( (S2(i).xd-S2(n+1).xd)^2 + (S2(i).yd-S2(n+1).yd)^2 );
            if (min_dis>do)
                S2(i).E=S2(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S2(i).E=S2(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS2=packets_TO_BS2+1;
   end
  end
end
STATISTICS.PACKETS_TO_CH2(r+1)=packets_TO_CH2;
STATISTICS.PACKETS_TO_BS2(r+1)=packets_TO_BS2;
end
first_dead2
teenth_dead2
all_dead2
STATISTICS.DEAD2(r+1)
STATISTICS.ALLIVE2(r+1)
STATISTICS.PACKETS_TO_CH2(r+1)
STATISTICS.PACKETS_TO_BS2(r+1)
STATISTICS.COUNTCHS2(r+1)













%Strat of EDEEC algo




for i=1:1:super
    
    
    S3(i).E=Eo*(1+b);
    
end
talha1=super+advance;
for i=super:1:talha1
    
    S3(i).E=Eo*(1+a);
    
end
for i=talha1:1:n
    
    S3(i).E=Eo;
    
end




for r=0:1:rmax     
    r
  if(mod(r, round(1/P) )==0)
    for i=1:1:n
        S3(i).G=0;
        S3(i).cl=0;
    end
  end
 Ea=Et*(1-r/rmax)/n;
 dead3=0;
for i=1:1:n
   
    if (S3(i).E<=0)
        dead3=dead3+1; 
        if (dead3==1)
           if(flag_first_dead3==0)
              first_dead3=r;
              flag_first_dead3=1;
           end
        end
        if(dead3==0.1*n)
           if(flag_teenth_dead3==0)
              teenth_dead3=r;
              flag_teenth_dead3=1;
           end
        end
        if(dead3==n)
           if(flag_all_dead3==0)
              all_dead3=r;
              flag_all_dead3=1;
           end
        end
    end
    if S3(i).E>0
        S3(i).type='N';
    end
end

STATISTICS.DEAD3(r+1)=dead3;
STATISTICS.ALLIVE3(r+1)=allive3-dead3;
countCHs3=0;
cluster3=1;
for i=1:1:n
    
 
 if Ea>0
     if (S3(i).E<= Eo)
         p(i)=P*E(i)/(1+m*(a+mo*b))*Ea;
     end
         
    else if (S3(i).E<=Eo*(1+a))
         p(i)=P*(1+a)*E(i)/(1+m*(a+mo*b))*Ea;
        end
    if (S3(i).E<=Eo*(1+b))
         p(i)=P*(1+b)*E(i)/(1+m*(a+mo*b))*Ea;
     end
     
 if(S3(i).E>0)
   temp_rand=rand;     
   if ( (S3(i).G)<=0)  
        if(temp_rand<= (p(i)/(1-(p(i)*(r*mod(1,p(i)))))))
            countCHs3=countCHs3+1;
            packets_TO_BS3=packets_TO_BS3+1;
            PACKETS_TO_BS3(r+1)=packets_TO_BS3;
             S3(i).type='C';
            S3(i).G=round(1/p(i))-1;
            C3(cluster3).xd=S3(i).xd;
            C3(cluster3).yd=S3(i).yd;
           distance=sqrt( (S3(i).xd-(S3(n+1).xd) )^2 + (S3(i).yd-(S3(n+1).yd) )^2 );
            C3(cluster3).distance=distance;
            C3(cluster3).id=i;
            X(cluster3)=S3(i).xd;
            Y(cluster3)=S3(i).yd;
            cluster3=cluster3+1;
           distance;
            if (distance>do)
                S3(i).E=S3(i).E- ( (ETX+EDA)*(4000) + Emp*4000*( distance*distance*distance*distance )); 
            end
            if (distance<=do)
                S3(i).E=S3(i).E- ( (ETX+EDA)*(4000)  + Efs*4000*( distance * distance )); 
            end
        end     
    
   end
   
 end 
 end
 end
STATISTICS.COUNTCHS3(r+1)=countCHs3;

for i=1:1:n
   if ( S3(i).type=='N' && S3(i).E>0 )
     if(cluster3-1>=1)
       min_dis=sqrt( (S3(i).xd-S3(n+1).xd)^2 + (S3(i).yd-S3(n+1).yd)^2 );
       min_dis_cluster=0;
       for c=1:1:cluster3-1
           temp=min(min_dis,sqrt( (S3(i).xd-C3(c).xd)^2 + (S3(i).yd-C3(c).yd)^2 ) );
           if ( temp<min_dis )
               min_dis=temp;
               min_dis_cluster=c;
           end
       end
      
       if(min_dis_cluster~=0)    
            min_dis;
            if (min_dis>do)
                S3(i).E=S3(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S3(i).E=S3(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
      
            S3(C3(min_dis_cluster).id).E = S3(C3(min_dis_cluster).id).E- ( (ERX + EDA)*4000 ); 
            packets_TO_CH3=packets_TO_CH3+1;
       else 
            min_dis;
            if (min_dis>do)
                S3(i).E=S3(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S3(i).E=S3(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS3=packets_TO_BS3+1;
            
       end
        S3(i).min_dis=min_dis;
       S3(i).min_dis_cluster=min_dis_cluster;
   else
            min_dis=sqrt( (S3(i).xd-S3(n+1).xd)^2 + (S3(i).yd-S3(n+1).yd)^2 );
            if (min_dis>do)
                S3(i).E=S3(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S3(i).E=S3(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS3=packets_TO_BS3+1;
   end
  end
end
STATISTICS.PACKETS_TO_CH3(r+1)=packets_TO_CH3;
STATISTICS.PACKETS_TO_BS3(r+1)=packets_TO_BS3;
end
first_dead3
teenth_dead3
all_dead3
STATISTICS.DEAD3(r+1)
STATISTICS.ALLIVE3(r+1)
STATISTICS.PACKETS_TO_CH3(r+1)
STATISTICS.PACKETS_TO_BS3(r+1)
STATISTICS.COUNTCHS3(r+1)







%Start of TDEEC algo














for i=1:1:super
  
    S4(i).E=Eo*(1+b);
    
    

end
talha1=super+advance;
for i=super:1:talha1
    
    S4(i).E=Eo*(1+a);
    
end
for i=talha1:1:n
    
    S4(i).E=Eo;
    
end





for r=0:1:rmax     
    r
  if(mod(r, round(1/P) )==0)
    for i=1:1:n
        S4(i).G=0;
        S4(i).cl=0;
    end
  end
 Ea=Et*(1-r/rmax)/n;
 dead4=0;
for i=1:1:n
   
    if (S4(i).E<=0)
        dead4=dead4+1; 
        if (dead4==1)
           if(flag_first_dead4==0)
              first_dead4=r;
              flag_first_dead4=1;
           end
        end
        if(dead4==0.1*n)
           if(flag_teenth_dead4==0)
              teenth_dead4=r;
              flag_teenth_dead4=1;
           end
        end
        if(dead4==n)
           if(flag_all_dead4==0)
              all_dead4=r;
              flag_all_dead4=1;
           end
        end
    end
    if S4(i).E>0
        S4(i).type='N';
    end
end

STATISTICS.DEAD4(r+1)=dead4;
STATISTICS.ALLIVE4(r+1)=allive4-dead4;
countCHs4=0;
cluster4=1;
for i=1:1:n
    
 
 if Ea>0
     if (S4(i).E<= Eo)
         p(i)=P*E(i)/(1+m*(a+mo*b))*Ea;
     end
         
     else if (S4(i).E<=Eo*(1+a))
         p(i)=P*(1+a)*E(i)/(1+m*(a+mo*b))*Ea;
         end
    if (S4(i).E<=Eo*(1+b))
         p(i)=P*(1+b)*E(i)/(1+m*(a+mo*b))*Ea;
     end
     

 if(S4(i).E>0)
   temp_rand=rand;     
   if ( (S4(i).G)<=0)  
        if((temp_rand<= (p(i)*E(i)*K/(1-p(i)*r*mod(1,p(i))*Ea))))
            countCHs4=countCHs4+1;
            packets_TO_BS4=packets_TO_BS4+1;
            PACKETS_TO_BS4(r+1)=packets_TO_BS4;
             S4(i).type='C';
            S4(i).G=round(1/p(i))-1;
            C4(cluster4).xd=S4(i).xd;
            C4(cluster4).yd=S4(i).yd;
           distance=sqrt( (S4(i).xd-(S4(n+1).xd) )^2 + (S4(i).yd-(S4(n+1).yd) )^2 );
            C4(cluster4).distance=distance;
            C4(cluster4).id=i;
            X4(cluster4)=S4(i).xd;
            Y4(cluster4)=S4(i).yd;
            cluster4=cluster4+1;
           distance;
            if (distance>do)
                S4(i).E=S4(i).E- ( (ETX+EDA)*(4000) + Emp*4000*( distance*distance*distance*distance )); 
            end
            if (distance<=do)
                S4(i).E=S4(i).E- ( (ETX+EDA)*(4000)  + Efs*4000*( distance * distance )); 
            end
        end     
    
   end
   
 end 
 end
 end
STATISTICS.COUNTCHS4(r+1)=countCHs4;


for i=1:1:n
   if ( S4(i).type=='N' && S4(i).E>0 )
     if(cluster4-1>=1)
       min_dis=sqrt( (S4(i).xd-S4(n+1).xd)^2 + (S4(i).yd-S4(n+1).yd)^2 );
       min_dis_cluster=0;
       for c=1:1:cluster4-1
           temp=min(min_dis,sqrt( (S4(i).xd-C4(c).xd)^2 + (S4(i).yd-C4(c).yd)^2 ) );
           if ( temp<min_dis )
               min_dis=temp;
               min_dis_cluster=c;
           end
       end
       
       
       if(min_dis_cluster~=0)    
            min_dis;
            if (min_dis>do)
                S4(i).E=S4(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S4(i).E=S4(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
      
            S4(C4(min_dis_cluster).id).E = S4(C4(min_dis_cluster).id).E- ( (ERX + EDA)*4000 ); 
            packets_TO_CH4=packets_TO_CH4+1;
       else 
            min_dis;
            if (min_dis>do)
                S4(i).E=S4(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S4(i).E=S4(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS4=packets_TO_BS4+1;
            
       end
        S4(i).min_dis=min_dis;
       S4(i).min_dis_cluster=min_dis_cluster;
   else
            min_dis=sqrt( (S4(i).xd-S4(n+1).xd)^2 + (S4(i).yd-S4(n+1).yd)^2 );
            if (min_dis>do)
                S4(i).E=S4(i).E- ( ETX*(4000) + Emp*4000*( min_dis * min_dis * min_dis * min_dis)); 
            end
            if (min_dis<=do)
                S4(i).E=S4(i).E- ( ETX*(4000) + Efs*4000*( min_dis * min_dis)); 
            end
            packets_TO_BS4=packets_TO_BS4+1;
   end
  end
end
STATISTICS.PACKETS_TO_CH4(r+1)=packets_TO_CH4;
STATISTICS.PACKETS_TO_BS4(r+1)=packets_TO_BS4;
end
first_dead4
teenth_dead4
all_dead4
STATISTICS.DEAD4(r+1)
STATISTICS.ALLIVE4(r+1)
STATISTICS.PACKETS_TO_CH4(r+1)
STATISTICS.PACKETS_TO_BS4(r+1)
STATISTICS.COUNTCHS4(r+1)


r=0:10000;
subplot(2,2,1);
plot(r,STATISTICS.DEAD1,':g',r,STATISTICS.DEAD2,'-r',r,STATISTICS.DEAD3,'--b',r,STATISTICS.DEAD4,'-.m');
legend('DEEC','DDEEC','EDEEC','TDEEC');
xlabel('x(rounds)');
ylabel('y(nodes dead)');
title('Nodes dead during rounds');
subplot(2,2,2);
plot(r,STATISTICS.ALLIVE1,':g',r,STATISTICS.ALLIVE2,'-r',r,STATISTICS.ALLIVE3,'--b',r,STATISTICS.ALLIVE4,'-.m');
legend('DEEC','DDEEC','EDEEC','TDEEC');
xlabel('x(rounds)');
ylabel('y(nodes alive)');
title('Nodes alive during rounds');
subplot(2,2,3);
plot(r,STATISTICS.PACKETS_TO_BS1,':g',r,STATISTICS.PACKETS_TO_BS2,'-r',r,STATISTICS.PACKETS_TO_BS3,'--b',r,STATISTICS.PACKETS_TO_BS4,'-.m');
legend('DEEC','DDEEC','EDEEC','TDEEC');
xlabel('x(rounds)');
ylabel('y(packets sent)');
title('Pckets sent to the base station');
subplot(2,2,4);
plot(r,STATISTICS.COUNTCHS1,':g',r,STATISTICS.COUNTCHS2,'-r',r,STATISTICS.COUNTCHS3,'--b',r,STATISTICS.COUNTCHS4,'-.m');
legend('DEEC','DDEEC','EDEEC','TDEEC');
xlabel('x(rounds)');
ylabel('y(cluster heads)');
title('Clusterheads');