use Tk;

$cell = 10; # Pixels per cell
$size = 50; # Number of cells along each edge of the simulation box

$base_drain = 0.02;  # Minimal drainage per time step
$max_drain  = 0.05;  # Range of variable drainage

$max_fluid = 1.0;  # Fluid content of a fresh cell
$auto_pop  = 0.0;  # if fluid content less than this, bubble will burst ("pop")
$pert_pop  = 0.45; # if less than this AND neighouring bubble just popped: pop

$drain_delay = 500; # Milliseconds between updates of graphics: drain mode
$pop_delay   = 100; # ... and pop mode

$idx    =  0; # Time step counter
@data   = (); # The simulation box
@popped = (); # Bubbles popped in the current avalance

$bub = 0;     # Count the total number of bubbles burst in current avalanche

# ==================================================

$mw = MainWindow->new;
$mw->configure( -width=>500, -height=>500 );
$mw->resizable( 0, 0 );
$canvas = $mw->Canvas( -cursor=>"crosshair", -background=>"white",
                       -width=>500, -height=>500 );
$canvas->pack( -side=>"right" );
$canvas->Tk::bind( "<Button-1>", \&savePS );

# =====

# Set up the simulation box.
# Each cell has fluid content, individual drain rate, and an associated widget
for my $y ( 0..$size ) {
  for my $x ( 0..$size ) {
    $data[$x][$y] = { fluid => $max_fluid, drain => new_drain() };
    $data[$x][$y]{widget} = $canvas->createRectangle( $cell*$x, $cell*$y,
						      $cell*($x+1),
						      $cell*($y+1),
						      -fill=> 'green',
						      -outline=> 'green' );
  }
}

# Smooth the drain rates over nearest neighbours - note: weighted average
for my $y ( 0..$size ) {
  for my $x ( 0..$size ) {
    $data[$x][$y]{drain} = avg_drain( $x, $y );
  }
}

dispatch();

MainLoop;

# ==================================================

# Create lists of neighbouring cells
sub nearest_nb {
  my ( $x, $y ) = @_;
  return ( [$x-1,$y], [($x+1)%$size,$y], [$x,$y-1], [$x,($y+1)%$size] );
}

sub extended_nb {
  my ( $x, $y ) = @_;
  return ( [$x-1,$y], [($x+1)%$size,$y], [$x,$y-1], [$x,($y+1)%$size],
	   [$x-1,$y-1], [$x-1,($y+1)%$size],
	   [($x+1)%$size,$y-1], [($x+1)%$size,($y+1)%$size],
	   [$x-2,$y], [($x+2)%$size,$y], [$x,$y-2], [$x,($y+2)%$size] );
}

sub new_drain {
  return $base_drain + rand $max_drain;
}

# Weighted average over nearest neighbours
sub avg_drain {
  my ( $x, $y ) = @_;

  my @nb = nearest_nb( $x, $y );

  my $d = 2*$data[$x][$y]{drain};
  for my $nb ( @nb ) { $d += $data[$$nb[0]][$$nb[1]]{drain}; }

  return $d/6.0;
}

# ==================================================

sub dispatch {
  ++$idx;

  redraw();

  # When there are popped bubbles, finish popping before continuing to drain
  if( scalar @popped ) {
    $bub += scalar @popped;

    push @popped, 'null'; # Separate different generations of popped bubbles

    do_pop();
    $mw->after( $pop_delay, \&dispatch );

  } else {
    # if( $bub > 0 ) { print "$bub\n"; $bub = 0; } # print size of avalanche

    do_drain();
    $mw->after( $drain_delay, \&dispatch );
  }
}

# =====

# Drain each cell and collect any cell with less fluid than the threshold
sub do_drain {
  for my $y ( 0..$size ) {
    for my $x ( 0..$size ) {
      $data[$x][$y]{fluid} -= $data[$x][$y]{drain};

      if( $data[$x][$y]{fluid} <= $auto_pop ) {
	push @popped, [ $x, $y ];
      }
    }
  }
}

# =====

sub do_pop {
  while( my $p = shift @popped ) {

    if( $p eq 'null' ) { last; }

    # Make a fresh bubble for any bubble that burst in the previous step
    my ( $x, $y ) = @$p;
    $data[$x][$y]{fluid} = $max_fluid;
    $data[$x][$y]{drain} = new_drain();

    $data[$x][$y]{drain} = avg_drain( $x, $y );

    # Neighbour list, includes nearest, next-nearest, and 3rd-level neighbours
    my @nb = extended_nb( $x, $y );
    for my $nb ( @nb ) {
      my ( $xx, $yy ) = @$nb;

      # Magic value to signal "popped" - this is nasty!
      if( $data[$xx][$yy]{fluid} <= $pert_pop &&
	  $data[$xx][$yy]{fluid} > -5 ){

	$data[$xx][$yy]{fluid} = -5;

	push @popped, [ $xx, $yy ];
      }
    }
  }
}

# =====

sub redraw {
  my ( $color, $tmp ) = ( '', 0 );

  for my $y ( 0..$size ) {
    for my $x ( 0..$size ) {
      if( $data[$x][$y]{fluid} <= 0 ) {
	$color = "red";

      } elsif( $data[$x][$y]{fluid} > $pert_pop ) {

	$tmp = ($data[$x][$y]{fluid}-$pert_pop)/(1 - $pert_pop);
	$color = hsvToRgbStr( 60*$tmp + 60, 1, 1 );

      } elsif( $data[$x][$y]{fluid} < $pert_pop ) {

	$tmp = ($pert_pop - $data[$x][$y]{fluid})/$pert_pop;
	$color = hsvToRgbStr( 60 - 15*$tmp, 1, 1 );

#	$color = 'yellow'; # Override - keep color fixed below $pert_pop
      }

      $canvas->itemconfigure( $data[$x][$y]{widget}, -fill => $color );
    }
  }
}

# Save current state to postscript file on mouse-click anywhere
sub savePS {
  my ( $gray, $color ) = ( $idx . '_gray.eps', $idx . '_color.eps' );
  $canvas->postscript( -file=>$gray, -colormode=>'gray' );
  $canvas->postscript( -file=>$color, -colormode=>'color' );
}

# ==================================================

sub hsvToRgbStr {
  my ( $r, $g, $b ) = hsvToRgb( @_ );
  return sprintf( "#%02x%02x%02x", $r, $g, $b );
}

sub hsvToRgb {
  my ( $h, $s, $v ) = @_;

  $v *= 255;

  if( $s == 0 ) { return ($v, $v, $v); } # achromatic (grey)

  my $i = int( $h/60 );  # sector 0 to 5
  my $f = ($h/60) - $i;  # fractional part of h/60

  my $p = $v * ( 1 - $s );
  my $q = $v * ( 1 - $s * $f );
  my $t = $v * ( 1 - $s * ( 1 - $f ) );

  if(    $i==0 ) { return ( $v, $t, $p ); }
  elsif( $i==1 ) { return ( $q, $v, $p ); }
  elsif( $i==2 ) { return ( $p, $v, $t ); }
  elsif( $i==3 ) { return ( $p, $q, $v ); }
  elsif( $i==4 ) { return ( $t, $p, $v ); }
  else           { return ( $v, $p, $q ); }
}