#!/usr/bin/perl
#

use strict;
use warnings;
use GD;
use Math::Trig qw( acos );

use constant EPSILON => 1e-10;
sub clockwise( $$$$$$ );
sub convex_hull_graham( @ );

my $pi = 4 * atan2( 1, 1 );

my ( $x, $y, $r );
my @points = map {
   my $r  = rand 200;
   my $th = rand 2 * $pi;
   ( $r * cos( $th ) + 300, $r * sin( $th ) + 300 );
} 1 .. 4;

my @hull = convex_hull_graham( @points );

my $graph = new GD::Image( 600, 600 );
my $white  = $graph->colorAllocate( 255, 255, 255 );
my $green  = $graph->colorAllocate( 0,   255, 0 );
my $blue   = $graph->colorAllocate( 0,   0,   255 );

push @hull, $hull[ 0 ], $hull[ 1 ];
for ( my $i = 0; $i < $#hull - 2; $i += 2 )
   {
   $graph->arc( $hull[ $i ], $hull[ $i + 1 ], 10, 10, 0, 360, $blue );
   $graph->fill( $hull[ $i ], $hull[ $i + 1 ], $blue );
   }
for ( my $i = 0; $i < $#hull - 2; $i += 2 )
   {
   $graph->line( $hull[ $i ],     $hull[ $i + 1 ], $hull[ $i + 2 ],
                 $hull[ $i + 3 ], $green
                 );
   }

open FILE, ">hull.png";
binmode FILE;
print FILE $graph->png;
close FILE;

sub convex_hull_graham( @ )
   {
   my @xy = @_;

   my $n = @xy / 2;
   my @i = map { 2 * $_ } 0 .. ( $n - 1 );
   my @x = map { $xy[ $_ ] } @i;
   my @y = map { $xy[ $_ + 1 ] } @i;

   my ( $ymin, $mini, $xmin, $i );

   for ( $ymin = $y[ 0 ], $mini = 0, $i = 1; $i < $n; ++$i )
      {
      if ( ( $y[ $i ] + EPSILON ) < $ymin )
         {
         $ymin = $y[ $i ];
         $mini = $i;
         }
      elsif ( abs( $y[ $i ] - $ymin ) < EPSILON )
         {
         $mini = $i if $x[ $i ] < $x[ $mini ];
         }
      }

   $xmin = $x[ $mini ];
   splice @x, $mini, 1;
   splice @y, $mini, 1;

   my @a = map { atan2( $y[ $_ ] - $ymin, $x[ $_ ] - $xmin ) } 0 .. $#
+x;

   my @j = map { $_->[ 0 ] } sort {
      return $a->[ 1 ] <=> $b->[ 1 ] || $x[ $a->[ 0 ] ] <=> $x[ $b->[ 
+0 ] ]
         || $y[ $a->[ 0 ] ] <=> $y[ $b->[ 0 ] ];
   } map { [ $_, $a[ $_ ] ] } 0 .. $#a;

   @x = @x[ @j ];
   @y = @y[ @j ];
   @a = @a[ @j ];

   unshift @x, $xmin;
   unshift @y, $ymin;
   unshift @a, 0;

   my @h = ( 0, 1 );
   my $cw;

   for ( $i = 2; $i < $n; ++$i )
      {
      while (
         clockwise(
                    $x[ $h[ $#h - 1 ] ], $y[ $h[ $#h - 1 ] ],
                    $x[ $h[ $#h ] ],     $y[ $h[ $#h ] ],
                    $x[ $i ],            $y[ $i ]
         ) > EPSILON and @h >= 2
         )
         {
         pop @h, scalar @h;
         }
      push @h, $i;
      scalar @h;
      }
   return map { ( $x[ $_ ], $y[ $_ ] ) } @h;
   }

sub clockwise( $$$$$$ )
   {
   my ( $x0, $y0, $x1, $y1, $x2, $y2 ) = @_;
   return ( $x2 - $x0 ) * ( $y1 - $y0 ) - ( $x1 - $x0 ) * ( $y2 - $y0 
+);
   }
[download]