#!/usr/bin/perl -w

use 5.034;

my $N = 1000;
my (@x, @y, @z, @c);
for (1..$N){
    $x[$_]=rand();
    $y[$_]=rand();
    $z[$_]= 1000*rand();   # peso de cada punto: demanda, habitantes, etc.
    $c[$_]=0;  # al principio, no pertenecen
}   

my $cluster =0;    # counter
my ($uno,$otro);
my $iteracion = 0;
while (1){
    my $min = 9e99;
    for my $i (1..$N+$cluster){
	next unless ($c[$i]==0);
	for my $j ($i+1..$N+$cluster){
	    next unless ($c[$j]==0);
	    my $dist = sqrt ( ($x[$i]-$x[$j])**2 + ($y[$i]-$y[$j])**2 );
	    if ($dist < $min){ 
		$min = $dist;
		$uno = $i;
		$otro = $j;
	    }
	}
    }
    my $count=0;
    my $sum = 0;
    if ($uno<=$N && $otro<=$N){    ###                        CREACION
	$cluster++;
	$c[$uno]=$N+$cluster;
	$c[$otro]=$N+$cluster;
	$x[$N+$cluster]=($x[$uno]+$x[$otro])/2;
	$y[$N+$cluster]=($y[$uno]+$y[$otro])/2;
	$z[$N+$cluster] = $z[$uno] + $z[$otro];
	$c[$N+$cluster]=0;   # esta línea es dudosa pero el centroide NO pretenece a ningún cluster porque no es un punto real
    }elsif (($uno<=$N && $otro>$N) || ($uno>$N && $otro<=$N)){         ###                         ABSORCION
	if ($uno<=$N){
	    ($uno, $otro) = ($otro, $uno);
	}
	$c[$otro]=$uno;
	$x[$uno] = $y[$uno]= $z[$uno] = 0;
	for my $i (1..$N){
	    if ($c[$i]==$uno){
		$x[$uno] += $x[$i]*$z[$i];
		$y[$uno] += $y[$i]*$z[$i];
		$z[$uno] += $z[$i];
		$count++;
		$sum += $z[$i];
	    }
	}
	$x[$uno]/=$sum;
	$y[$uno]/=$sum;
    }elsif ($uno>$N && $otro>$N){     ###                                FUSION
	if ($uno>$otro){ 
	    ($uno,$otro) = ($otro,$uno);
	}
	$c[$otro] = -1;   # cluster borrado y no hace falta más: sus números irán salteados
#	$cluster--;
#	delete $c[$otro];
#	delete $x[$otro];
#	delete $y[$otro];
	$x[$uno]= $y[$uno]= $z[$uno] = 0;
	for my $i (1..$N){
	    if ($c[$i]==$otro){ $c[$i]=$uno; }
	    if ($c[$i]==$uno){
		$x[$uno] += $x[$i]*$z[$i];
		$y[$uno] += $y[$i]*$z[$i];
		$z[$uno] += $z[$i];
		$count++;
		$sum += $z[$i];		
	    }
	}
	$x[$uno]/=$sum;
	$y[$uno]/=$sum;
    }else{
	die "Error del programa";
    }
    ++$iteracion;
    print "=" x 80, "\n";
    say "Iteration $iteracion";
    for my $i ($N+1..$N+$cluster){
	next unless ($c[$i]==0);
	say "\t\tcluster $i, centroid x=$x[$i] y=$y[$i]   z=$z[$i]"; 
	print "points: ";
	for my $j (1..$N){
	    if ($c[$j]==$i){ print "$j "; }
	}
	print "\n";       # ,"-" x 70,"\n";
    }
    print "Outliers: ";
    for my $i (1..$N){
	print "$i " if ($c[$i]==0);
    }
    print "\n\n";
    my $counter = 0;
    for (1..$N+$cluster){
	if ($c[$_]==0) { $counter++; }
    }
    last if ($counter == 1);
}

__END__