www.pudn.com > testConvexHull.rar > testConvexHull.cpp
// testConvexHull.cpp : Defines the entry point for the console application.
//
#include "stdafx.h"
//
// This program uses the Graham scan algorithm to calculate the convex
// hull of a batch of points
//
#include
#include
#include
#include
#include
#include
std::ostream &operator <<(std::ostream &s, const std::pair &point )
{
s << "("
<< point.first
<< ","
<< point.second
<< ")";
return s;
}
class GrahamScan
{
public :
GrahamScan( size_t n, int xmin, int xmax, int ymin, int ymax )
: N( n )
, x_range( xmin, xmax )
, y_range( ymin, ymax )
{
//
// In this constructor I generate the N random points asked for
// by the caller. Their values are randomly assigned in the x/y
// ranges specified as arguments
//
srand( static_cast( time(NULL) ) );
for ( size_t i = 0 ; i < N ; i++ ) {
int x = ( rand() % ( x_range.second - x_range.first + 1 ) ) + x_range.first;
int y = ( rand() % ( y_range.second - y_range.first + 1 ) ) + y_range.first;
raw_points.push_back( std::make_pair( x, y ) );
}
}
//
// The initial array of points is stored in vectgor raw_points. I first
// sort it, which gives me the far left and far right points of the hull.
// These are special values, and they are stored off separately in the left
// and right members.
//
// I then go through the list of raw_points, and one by one determine whether
// each point is above or below the line formed by the right and left points.
// If it is above, the point is moved into the upper_partition_points sequence. If it
// is below, the point is moved into the lower_partition_points sequence. So the output
// of this routine is the left and right points, and the sorted points that are in the
// upper and lower partitions.
//
void partition_points()
{
//
// Step one in partitioning the points is to sort the raw data
//
std::sort( raw_points.begin(), raw_points.end() );
//
// The the far left and far right points, remove them from the
// sorted sequence and store them in special members
//
left = raw_points.front();
raw_points.erase( raw_points.begin() );
right = raw_points.back();
raw_points.pop_back();
//
// Now put the remaining points in one of the two output sequences
//
for ( size_t i = 0 ; i < raw_points.size() ; i++ )
{
int dir = direction( left, right, raw_points[ i ] );
if ( dir < 0 )
upper_partition_points.push_back( raw_points[ i ] );
else
lower_partition_points.push_back( raw_points[ i ] );
}
}
//
// Building the hull consists of two procedures: building the lower and
// then the upper hull. The two procedures are nearly identical - the main
// difference between the two is the test for convexity. When building the upper
// hull, our rull is that the middle point must always be *above* the line formed
// by its two closest neighbors. When building the lower hull, the rule is that point
// must be *below* its two closest neighbors. We pass this information to the
// building routine as the last parameter, which is either -1 or 1.
//
void build_hull( std::ofstream &f )
{
build_half_hull( f, lower_partition_points, lower_hull, 1 );
build_half_hull( f, upper_partition_points, upper_hull, -1 );
}
//
// This is the method that builds either the upper or the lower half convex
// hull. It takes as its input a copy of the input array, which will be the
// sorted list of points in one of the two halfs. It produces as output a list
// of the points in the corresponding convex hull.
//
// The factor should be 1 for the lower hull, and -1 for the upper hull.
//
void build_half_hull( std::ostream &f,
std::vector< std::pair > input,
std::vector< std::pair > &output,
int factor )
{
//
// The hull will always start with the left point, and end with the right
// point. According, we start by adding the left point as the first point
// in the output sequence, and make sure the right point is the last point
// in the input sequence.
//
input.push_back( right );
output.push_back( left );
//
// The construction loop runs until the input is exhausted
//
while ( input.size() != 0 ) {
//
// Repeatedly add the leftmost point to the null, then test to see
// if a convexity violation has occured. If it has, fix things up
// by removing the next-to-last point in the output suqeence until
// convexity is restored.
//
output.push_back( input.front() );
plot_hull( f, "adding a new point" );
input.erase( input.begin() );
while ( output.size() >= 3 ) {
size_t end = output.size() - 1;
if ( factor * direction( output[ end - 2 ],
output[ end ],
output[ end - 1 ] ) <= 0 ) {
output.erase( output.begin() + end - 1 );
plot_hull( f, "backtracking" );
}
else
break;
}
}
}
//
// In this program we frequently want to look at three consecutive
// points, p0, p1, and p2, and determine whether p2 has taken a turn
// to the left or a turn to the right.
//
// We can do this by by translating the points so that p1 is at the origin,
// then taking the cross product of p0 and p2. The result will be positive,
// negative, or 0, meaning respectively that p2 has turned right, left, or
// is on a straight line.
//
static int direction( std::pair p0,
std::pair p1,
std::pair p2 )
{
return ( (p0.first - p1.first ) * (p2.second - p1.second ) )
- ( (p2.first - p1.first ) * (p0.second - p1.second ) );
}
void log_raw_points( std::ostream &f )
{
f << "Creating raw points:\n";
for ( size_t i = 0 ; i < N ; i++ )
f << raw_points[ i ] << " ";
f << "\n";
}
void log_partitioned_points( std::ostream &f )
{
f << "Partitioned set:\n"
<< "Left : " << left << "\n"
<< "Right : " << right << "\n"
<< "Lower partition: ";
for ( size_t i = 0 ; i < lower_partition_points.size() ; i++ )
f << lower_partition_points[ i ];
f << "\n";
f << "Upper partition: ";
for ( i = 0 ; i < upper_partition_points.size() ; i++ )
f << upper_partition_points[ i ];
f << "\n";
}
void log_hull( std::ostream & f )
{
f << "Lower hull: ";
for ( size_t i = 0 ; i < lower_hull.size() ; i++ )
f << lower_hull[ i ];
f << "\n";
f << "Upper hull: ";
for ( i = 0 ; i < upper_hull.size() ; i++ )
f << upper_hull[ i ];
f << "\n";
f << "Convex hull: ";
for ( i = 0 ; i < lower_hull.size() ; i++ )
f << lower_hull[ i ] << " ";
for ( std::vector< std::pair >::reverse_iterator ii = upper_hull.rbegin() + 1 ;
ii != upper_hull.rend();
ii ++ )
f << *ii << " ";
f << "\n";
}
void plot_raw_points( std::ostream &f )
{
f << "set xrange ["
<< x_range.first
<< ":"
<< x_range.second
<< "]\n";
f << "set yrange ["
<< y_range.first
<< ":"
<< y_range.second
<< "]\n";
f << "unset mouse\n";
f << "set title 'The set of raw points in the set' font 'Arial,12'\n";
f << "set style line 1 pointtype 7 linecolor rgb 'red'\n";
f << "set style line 2 pointtype 7 linecolor rgb 'green'\n";
f << "set style line 3 pointtype 7 linecolor rgb 'black'\n";
f << "plot '-' ls 1 with points notitle\n";
for ( size_t i = 0 ; i < N ; i++ )
f << raw_points[ i ].first << " " << raw_points[ i ].second << "\n";
f << "e\n";
f << "pause -1 'Hit OK to move to the next state'\n";
}
void plot_partitioned_points( std::ostream &f )
{
f << "set title 'The points partitioned into an upper and lower hull' font 'Arial,12'\n";
f << "plot '-' ls 1 with points notitle, "
<< "'-' ls 2 with points notitle, "
<< "'-' ls 3 with linespoints notitle\n";
for ( size_t i = 0 ; i < lower_partition_points.size() ; i++ )
f << lower_partition_points[ i ].first
<< " "
<< lower_partition_points[ i ].second
<< "\n";
f << "e\n";
for ( i = 0 ; i < upper_partition_points.size() ; i++ )
f << upper_partition_points[ i ].first
<< " "
<< upper_partition_points[ i ].second
<< "\n";
f << "e\n";
f << left.first << " " << left.second << "\n";
f << right.first << " " << right.second << "\n";
f << "e\n";
f << "pause -1 'Hit OK to move to the next state'\n";
}
void plot_hull( std::ostream &f, std::string text )
{
f << "set title 'The hull in state: "
<< text
<< "' font 'Arial,12'\n";
f << "plot '-' ls 1 with points notitle, ";
if ( lower_hull.size() )
f << "'-' ls 3 with linespoints notitle, ";
if ( upper_hull.size() )
f << "'-' ls 3 with linespoints notitle, ";
f << "'-' ls 2 with points notitle\n";
for ( size_t i = 0 ; i < lower_partition_points.size() ; i++ )
f << lower_partition_points[ i ].first
<< " "
<< lower_partition_points[ i ].second
<< "\n";
f << right.first << " " << right.second << "\n";
f << "e\n";
if ( lower_hull.size() ) {
for ( size_t i = 0 ; i < lower_hull.size() ; i++ )
f << lower_hull[ i ].first
<< " "
<< lower_hull[ i ].second
<< "\n";
f << "e\n";
}
if ( upper_hull.size() ) {
for ( std::vector< std::pair >::reverse_iterator ii = upper_hull.rbegin();
ii != upper_hull.rend();
ii++ )
f << ii->first
<< " "
<< ii->second
<< "\n";
f << "e\n";
}
for ( i = 0 ; i < upper_partition_points.size() ; i++ )
f << upper_partition_points[ i ].first
<< " "
<< upper_partition_points[ i ].second
<< "\n";
f << "e\n";
f << "pause -1 'Hit OK to move to the next state'\n";
}
private :
//
// These values determine the range of numbers generated to
// provide the input data. The values are all passed in as part
// of the constructor
//
const size_t N;
const std::pair x_range;
const std::pair y_range;
// The raw data points generated by the constructor
std::vector< std::pair > raw_points;
//
// These values are used to represent the partitioned set. A special
// leftmost and rightmost value, and the sorted set of upper and lower
// partitioned points that lie inside those two points.
//
std::pair left;
std::pair right;
std::vector< std::pair > upper_partition_points;
std::vector< std::pair > lower_partition_points;
//
// After the convex hull is created, the lower hull and upper hull
// are stored in these sorted sequences. There is a bit of duplication
// between the two, because both sets include the leftmost and rightmost point.
//
std::vector< std::pair > lower_hull;
std::vector< std::pair > upper_hull;
};
int main(int argc, char* argv[])
{
std::ofstream gnuplot_file( "gnuplot.cmd" );
const int N = 10;
GrahamScan g( N, 0, 100, 0, 100 );
g.log_raw_points( std::cout );
g.plot_raw_points( gnuplot_file );
g.partition_points();
g.log_partitioned_points( std::cout );
g.plot_partitioned_points( gnuplot_file );
//
// Okay, time to start building the hull
//
g.build_hull( gnuplot_file );
g.log_hull( std::cout );
g.plot_hull( gnuplot_file, "complete" );
return 0;
}