fastjet 2.4.5
|
concrete implementation for finding closest pairs in 2D -- will use Chan's (hopefully efficient) shuffle based structures More...
#include <ClosestPair2D.hh>
Classes | |
class | Point |
class for representing all info needed about a point More... | |
class | Shuffle |
class that will take care of ordering of shuffles for us More... | |
class | triplet |
since sets of three objects will crop up repeatedly, useful to have a triplet class? More... | |
Public Member Functions | |
ClosestPair2D (const std::vector< Coord2D > &positions, const Coord2D &left_corner, const Coord2D &right_corner) | |
constructor from a vector of 2D positions -- number of objects after insertion and deletion must never exceed positions.size(); objects are given IDs that correspond to their index in the vector of positions | |
ClosestPair2D (const std::vector< Coord2D > &positions, const Coord2D &left_corner, const Coord2D &right_corner, const unsigned int max_size) | |
constructor which allows structure to grow beyond positions.size(), up to max_size | |
void | closest_pair (unsigned int &ID1, unsigned int &ID2, double &distance2) const |
provides the IDs of the closest pair as well as the distance between them | |
void | remove (unsigned int ID) |
removes the entry labelled by ID from the object; | |
unsigned int | insert (const Coord2D &) |
inserts the position into the closest pair structure and returns the ID that has been allocated for the object. | |
virtual unsigned int | replace (unsigned int ID1, unsigned int ID2, const Coord2D &position) |
removes ID1 and ID2 and inserts position, returning the ID corresponding to position... | |
virtual void | replace_many (const std::vector< unsigned int > &IDs_to_remove, const std::vector< Coord2D > &new_positions, std::vector< unsigned int > &new_IDs) |
replaces IDs_to_remove with points at the new_positions indicating the IDs allocated to the new points in new_IDs | |
void | print_tree_depths (std::ostream &outdev) const |
unsigned int | size () |
Private Types | |
typedef SearchTree< Shuffle > | Tree |
typedef Tree::circulator | circulator |
typedef Tree::const_circulator | const_circulator |
Private Member Functions | |
void | _initialize (const std::vector< Coord2D > &positions, const Coord2D &left_corner, const Coord2D &right_corner, const unsigned int max_size) |
void | _add_label (Point *point, unsigned int review_flag) |
add a label to a point as to the nature of review needed (includes adding it to list of points needing review) [doesn't affect other labels already set for the point] | |
void | _set_label (Point *point, unsigned int review_flag) |
sets the label for the point to be exclusively this review flag (and adds it to list of points needing review if not already there) | |
void | _deal_with_points_to_review () |
for all entries of the _points_under_review[] vector, carry out the actions indicated by its review flag; the points are then removed from _points_under_review[] and their flags set to zero | |
void | _remove_from_search_tree (Point *point_to_remove) |
carry out the search-tree related operations of point removal | |
void | _insert_into_search_tree (Point *new_point) |
carry out the search-tree related operations of point insertion | |
void | _point2shuffle (Point &, Shuffle &, unsigned int shift) |
takes a point and creates a shuffle with the given shift | |
int | _ID (const Point *) const |
returns the ID for the specified point... | |
Private Attributes | |
triplet< std::auto_ptr< Tree > > | _trees |
std::auto_ptr< MinHeap > | _heap |
std::vector< Point > | _points |
std::stack< Point * > | _available_points |
std::vector< Point * > | _points_under_review |
points that are "under review" in some way | |
Coord2D | _left_corner |
pieces needed for converting coordinates to integer | |
double | _range |
triplet< unsigned int > | _shifts |
triplet< unsigned int > | _rel_shifts |
unsigned int | _cp_search_range |
Static Private Attributes | |
static const unsigned int | _nshift = 3 |
static const unsigned int | _remove_heap_entry = 1 |
static const unsigned int | _review_heap_entry = 2 |
static const unsigned int | _review_neighbour = 4 |
concrete implementation for finding closest pairs in 2D -- will use Chan's (hopefully efficient) shuffle based structures
Definition at line 46 of file ClosestPair2D.hh.
typedef Tree::circulator fastjet::ClosestPair2D::circulator [private] |
Definition at line 128 of file ClosestPair2D.hh.
typedef Tree::const_circulator fastjet::ClosestPair2D::const_circulator [private] |
Definition at line 129 of file ClosestPair2D.hh.
typedef SearchTree<Shuffle> fastjet::ClosestPair2D::Tree [private] |
Definition at line 127 of file ClosestPair2D.hh.
fastjet::ClosestPair2D::ClosestPair2D | ( | const std::vector< Coord2D > & | positions, |
const Coord2D & | left_corner, | ||
const Coord2D & | right_corner | ||
) | [inline] |
constructor from a vector of 2D positions -- number of objects after insertion and deletion must never exceed positions.size(); objects are given IDs that correspond to their index in the vector of positions
Definition at line 52 of file ClosestPair2D.hh.
{ _initialize(positions, left_corner, right_corner, positions.size()); };
fastjet::ClosestPair2D::ClosestPair2D | ( | const std::vector< Coord2D > & | positions, |
const Coord2D & | left_corner, | ||
const Coord2D & | right_corner, | ||
const unsigned int | max_size | ||
) | [inline] |
constructor which allows structure to grow beyond positions.size(), up to max_size
Definition at line 59 of file ClosestPair2D.hh.
{ _initialize(positions, left_corner, right_corner, max_size); };
void fastjet::ClosestPair2D::_add_label | ( | Point * | point, |
unsigned int | review_flag | ||
) | [inline, private] |
add a label to a point as to the nature of review needed (includes adding it to list of points needing review) [doesn't affect other labels already set for the point]
Definition at line 192 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::review_flag.
{ // if it's not already under review, then put it on the list of // points needing review if (point->review_flag == 0) _points_under_review.push_back(point); // OR the point's current flag with the requested review flag point->review_flag |= review_flag; }
void fastjet::ClosestPair2D::_deal_with_points_to_review | ( | ) | [private] |
for all entries of the _points_under_review[] vector, carry out the actions indicated by its review flag; the points are then removed from _points_under_review[] and their flags set to zero
Definition at line 296 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::circ, fastjet::ClosestPair2D::Point::distance2(), fastjet::d0::inline_maths::min(), fastjet::ClosestPair2D::Point::neighbour, fastjet::ClosestPair2D::Point::neighbour_dist2, and fastjet::ClosestPair2D::Point::review_flag.
{ // the range in which we carry out searches for new neighbours on // the search tree unsigned int CP_range = min(_cp_search_range, size()-1); // now deal with the points that are "under review" in some way // (have lost their neighbour, or need their heap entry updating) while(_points_under_review.size() > 0) { // get the point to be considered Point * this_point = _points_under_review.back(); // remove it from the list _points_under_review.pop_back(); if (this_point->review_flag & _remove_heap_entry) { // make sure no other flags are on (it wouldn't be consistent?) assert(!(this_point->review_flag ^ _remove_heap_entry)); _heap->remove(_ID(this_point)); } // check to see if the _review_neighbour flag is on else { if (this_point->review_flag & _review_neighbour) { this_point->neighbour_dist2 = numeric_limits<double>::max(); // among all three shifts for (unsigned int ishift = 0; ishift < _nshift; ishift++) { circulator other = this_point->circ[ishift]; // among points within CP_range for (unsigned i=0; i < CP_range; i++) { ++other; double dist2 = this_point->distance2(*other->point); if (dist2 < this_point->neighbour_dist2) { this_point->neighbour_dist2 = dist2; this_point->neighbour = other->point; } } } } // for any non-zero review flag we'll have to update the heap _heap->update(_ID(this_point), this_point->neighbour_dist2); } // "delabel" the point this_point->review_flag = 0; } }
int fastjet::ClosestPair2D::_ID | ( | const Point * | point | ) | const [inline, private] |
returns the ID for the specified point...
Definition at line 220 of file ClosestPair2D.hh.
{ return point - &(_points[0]); }
void fastjet::ClosestPair2D::_initialize | ( | const std::vector< Coord2D > & | positions, |
const Coord2D & | left_corner, | ||
const Coord2D & | right_corner, | ||
const unsigned int | max_size | ||
) | [private] |
Definition at line 79 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::circ, fastjet::ClosestPair2D::Point::distance2(), fastjet::d0::inline_maths::min(), fastjet::ClosestPair2D::Point::neighbour, fastjet::ClosestPair2D::Point::neighbour_dist2, fastjet::twopow31, fastjet::Coord2D::x, and fastjet::Coord2D::y.
{ unsigned int n_positions = positions.size(); assert(max_size >= n_positions); //_points(positions.size()) // allow the points array to grow to the following size _points.resize(max_size); // currently unused points are immediately made available on the // stack for (unsigned int i = n_positions; i < max_size; i++) { _available_points.push(&(_points[i])); } _left_corner = left_corner; _range = max((right_corner.x - left_corner.x), (right_corner.y - left_corner.y)); // initialise the coordinates for the points and create the zero-shifted // shuffle array vector<Shuffle> shuffles(n_positions); for (unsigned int i = 0; i < n_positions; i++) { // set up the points _points[i].coord = positions[i]; _points[i].neighbour_dist2 = numeric_limits<double>::max(); _points[i].review_flag = 0; // create shuffle with 0 shift. _point2shuffle(_points[i], shuffles[i], 0); } // establish what our shifts will be for (unsigned ishift = 0; ishift < _nshift; ishift++) { // make sure we use double-precision for calculating the shifts // since otherwise we will hit integer overflow. _shifts[ishift] = static_cast<unsigned int>(((twopow31*1.0)*ishift)/_nshift); if (ishift == 0) {_rel_shifts[ishift] = 0;} else {_rel_shifts[ishift] = _shifts[ishift] - _shifts[ishift-1];} } //_shifts[0] = 0; //_shifts[1] = static_cast<unsigned int>((twopow31*1.0)/3.0); //_shifts[2] = static_cast<unsigned int>((twopow31*2.0)/3.0); //_rel_shifts[0] = 0; //_rel_shifts[1] = _shifts[1]; //_rel_shifts[2] = _shifts[2]-_shifts[1]; // and how we will search... //_cp_search_range = 49; _cp_search_range = 30; _points_under_review.reserve(_nshift * _cp_search_range); // now initialise the three trees for (unsigned int ishift = 0; ishift < _nshift; ishift++) { // shift the shuffles if need be. if (ishift > 0) { unsigned rel_shift = _rel_shifts[ishift]; for (unsigned int i = 0; i < shuffles.size(); i++) { shuffles[i] += rel_shift; } } // sort the shuffles sort(shuffles.begin(), shuffles.end()); // and create the search tree _trees[ishift] = auto_ptr<Tree>(new Tree(shuffles, max_size)); // now we look for the closest-pair candidates on this tree circulator circ = _trees[ishift]->somewhere(), start=circ; // the actual range in which we search unsigned int CP_range = min(_cp_search_range, n_positions-1); do { Point * this_point = circ->point; //cout << _ID(this_point) << " "; this_point->circ[ishift] = circ; // run over all points within _cp_search_range of this_point on tree circulator other = circ; for (unsigned i=0; i < CP_range; i++) { ++other; double dist2 = this_point->distance2(*other->point); if (dist2 < this_point->neighbour_dist2) { this_point->neighbour_dist2 = dist2; this_point->neighbour = other->point; } } } while (++circ != start); //cout << endl<<endl; } // now initialise the heap object... vector<double> mindists2(n_positions); for (unsigned int i = 0; i < n_positions; i++) { mindists2[i] = _points[i].neighbour_dist2;} _heap = auto_ptr<MinHeap>(new MinHeap(mindists2, max_size)); }
void fastjet::ClosestPair2D::_insert_into_search_tree | ( | Point * | new_point | ) | [private] |
carry out the search-tree related operations of point insertion
Definition at line 431 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::circ, fastjet::ClosestPair2D::Point::distance2(), fastjet::d0::inline_maths::min(), fastjet::ClosestPair2D::Point::neighbour, and fastjet::ClosestPair2D::Point::neighbour_dist2.
{ // this point will have to have it's heap entry reviewed... _set_label(new_point, _review_heap_entry); // set the current distance to "infinity" new_point->neighbour_dist2 = numeric_limits<double>::max(); // establish how far we will be searching; unsigned int CP_range = min(_cp_search_range, size()-1); for (unsigned ishift = 0; ishift < _nshift; ishift++) { // create the shuffle Shuffle new_shuffle; _point2shuffle(*new_point, new_shuffle, _shifts[ishift]); // insert it into the tree circulator new_circ = _trees[ishift]->insert(new_shuffle); new_point->circ[ishift] = new_circ; // now get hold of the right and left edges of the region we will be // looking at (cf CCN28-43) circulator right_edge = new_circ; right_edge++; circulator left_edge = new_circ; for (unsigned int i = 0; i < CP_range; i++) {left_edge--;} // now do { Point * left_point = left_edge->point; Point * right_point = right_edge->point; // see if the new point is closer to the left-edge than the latter's // current neighbour double new_dist2 = left_point->distance2(*new_point); if (new_dist2 < left_point->neighbour_dist2) { left_point->neighbour_dist2 = new_dist2; left_point->neighbour = new_point; _add_label(left_point, _review_heap_entry); } // see if the right-point is closer to the new point than it's current // neighbour new_dist2 = new_point->distance2(*right_point); if (new_dist2 < new_point->neighbour_dist2) { new_point->neighbour_dist2 = new_dist2; new_point->neighbour = right_point; } // if the right-edge point was the left-edge's neighbour, then // then it's just gone off-radar and the left-point will need to // have its neighbour recalculated [actually, this is overdoing // it a little, since right point may be an less "distant" // (circulator distance) in one of the other shifts -- but not // sure how to deal with this...] if (left_point->neighbour == right_point) { _add_label(left_point, _review_neighbour); } // shift the left and right edges until left edge hits new_circ right_edge++; } while (++left_edge != new_circ); } }
void fastjet::ClosestPair2D::_point2shuffle | ( | Point & | point, |
Shuffle & | shuffle, | ||
unsigned int | shift | ||
) | [private] |
takes a point and creates a shuffle with the given shift
takes a point and sets a shuffle with the given shift...
Definition at line 47 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::coord, fastjet::ClosestPair2D::Shuffle::point, fastjet::twopow31, fastjet::ClosestPair2D::Shuffle::x, fastjet::Coord2D::x, fastjet::ClosestPair2D::Shuffle::y, and fastjet::Coord2D::y.
{ Coord2D renorm_point = (point.coord - _left_corner)/_range; // make sure the point is sensible //cerr << point.coord.x <<" "<<point.coord.y<<endl; assert(renorm_point.x >=0); assert(renorm_point.x <=1); assert(renorm_point.y >=0); assert(renorm_point.y <=1); shuffle.x = static_cast<unsigned int>(twopow31 * renorm_point.x) + shift; shuffle.y = static_cast<unsigned int>(twopow31 * renorm_point.y) + shift; shuffle.point = &point; }
void fastjet::ClosestPair2D::_remove_from_search_tree | ( | Point * | point_to_remove | ) | [private] |
carry out the search-tree related operations of point removal
Definition at line 230 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::circ, fastjet::ClosestPair2D::Point::distance2(), fastjet::d0::inline_maths::min(), fastjet::ClosestPair2D::Point::neighbour, fastjet::ClosestPair2D::Point::neighbour_dist2, and fastjet::SearchTree< T >::circulator::next().
{ // add this point to the list of "available" points (this is // relevant also from the point of view of determining the range // over which we circulate). _available_points.push(point_to_remove); // label it so that it goes from the heap _set_label(point_to_remove, _remove_heap_entry); // establish the range over which we search (a) for points that have // acquired a new neighbour and (b) for points which had ID as their // neighbour; unsigned int CP_range = min(_cp_search_range, size()-1); // then, for each shift for (unsigned int ishift = 0; ishift < _nshift; ishift++) { //cout << " ishift = " << ishift <<endl; // get the circulator for the point we'll remove and its successor circulator removed_circ = point_to_remove->circ[ishift]; circulator right_end = removed_circ.next(); // remove the point _trees[ishift]->remove(removed_circ); // next find the point CP_range points to the left circulator left_end = right_end, orig_right_end = right_end; for (unsigned int i = 0; i < CP_range; i++) {left_end--;} if (size()-1 < _cp_search_range) { // we have a smaller range now than before -- but when seeing who // could have had ID as a neighbour, we still need to use the old // range for seeing how far back we search (but new separation between // points). [cf CCN28-42] left_end--; right_end--; } // and then for each left-end point: establish if the removed // point was its neighbour [in which case a new neighbour must be // found], otherwise see if the right-end point is a closer neighbour do { Point * left_point = left_end->point; //cout << " visited " << setw(3)<<_ID(left_point)<<" (its neighbour was "<< setw(3)<< _ID(left_point->neighbour)<<")"<<endl; if (left_point->neighbour == point_to_remove) { // we'll deal with it later... _add_label(left_point, _review_neighbour); } else { // check to see if right point has become its closest neighbour double dist2 = left_point->distance2(*right_end->point); if (dist2 < left_point->neighbour_dist2) { left_point->neighbour = right_end->point; left_point->neighbour_dist2 = dist2; // NB: (LESSER) REVIEW NEEDED HERE TOO... _add_label(left_point, _review_heap_entry); } } ++right_end; } while (++left_end != orig_right_end); } // ishift... }
void fastjet::ClosestPair2D::_set_label | ( | Point * | point, |
unsigned int | review_flag | ||
) | [inline, private] |
sets the label for the point to be exclusively this review flag (and adds it to list of points needing review if not already there)
Definition at line 203 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::review_flag.
{ // if it's not already under review, then put it on the list of // points needing review if (point->review_flag == 0) _points_under_review.push_back(point); // SET the point's current flag to the requested review flag point->review_flag = review_flag; }
void fastjet::ClosestPair2D::closest_pair | ( | unsigned int & | ID1, |
unsigned int & | ID2, | ||
double & | distance2 | ||
) | const [virtual] |
provides the IDs of the closest pair as well as the distance between them
Implements fastjet::ClosestPair2DBase.
Definition at line 182 of file ClosestPair2D.cc.
Referenced by fastjet::ClusterSequence::_CP2DChan_cluster(), and fastjet::ClusterSequence::_CP2DChan_limited_cluster().
unsigned int fastjet::ClosestPair2D::insert | ( | const Coord2D & | new_coord | ) | [virtual] |
inserts the position into the closest pair structure and returns the ID that has been allocated for the object.
Implements fastjet::ClosestPair2DBase.
Definition at line 346 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::coord.
{ // get hold of a point assert(_available_points.size() > 0); Point * new_point = _available_points.top(); _available_points.pop(); // set the point's coordinate new_point->coord = new_coord; // now find it's neighbour in the search tree _insert_into_search_tree(new_point); // sort out other points that may have been affected by this, // and/or for which the heap needs to be updated _deal_with_points_to_review(); // return _ID(new_point); }
void fastjet::ClosestPair2D::print_tree_depths | ( | std::ostream & | outdev | ) | const [inline] |
Definition at line 89 of file ClosestPair2D.hh.
void fastjet::ClosestPair2D::remove | ( | unsigned int | ID | ) | [virtual] |
removes the entry labelled by ID from the object;
Implements fastjet::ClosestPair2DBase.
Definition at line 214 of file ClosestPair2D.cc.
{ //cout << "While removing " << ID <<endl; Point * point_to_remove = & (_points[ID]); // remove this point from the search tree _remove_from_search_tree(point_to_remove); // the above statement labels certain points as needing "review" -- // deal with them... _deal_with_points_to_review(); }
unsigned int fastjet::ClosestPair2D::replace | ( | unsigned int | ID1, |
unsigned int | ID2, | ||
const Coord2D & | position | ||
) | [virtual] |
removes ID1 and ID2 and inserts position, returning the ID corresponding to position...
Reimplemented from fastjet::ClosestPair2DBase.
Definition at line 368 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::coord.
Referenced by fastjet::ClusterSequence::_CP2DChan_cluster().
{ // deletion from tree... Point * point_to_remove = & (_points[ID1]); _remove_from_search_tree(point_to_remove); point_to_remove = & (_points[ID2]); _remove_from_search_tree(point_to_remove); // insertion into tree // get hold of a point Point * new_point = _available_points.top(); _available_points.pop(); // set the point's coordinate new_point->coord = position; // now find it's neighbour in the search tree _insert_into_search_tree(new_point); // the above statement labels certain points as needing "review" -- // deal with them... _deal_with_points_to_review(); // return _ID(new_point); }
void fastjet::ClosestPair2D::replace_many | ( | const std::vector< unsigned int > & | IDs_to_remove, |
const std::vector< Coord2D > & | new_positions, | ||
std::vector< unsigned int > & | new_IDs | ||
) | [virtual] |
replaces IDs_to_remove with points at the new_positions indicating the IDs allocated to the new points in new_IDs
Reimplemented from fastjet::ClosestPair2DBase.
Definition at line 400 of file ClosestPair2D.cc.
References fastjet::ClosestPair2D::Point::coord.
Referenced by fastjet::ClusterSequence::_CP2DChan_limited_cluster().
{ // deletion from tree... for (unsigned int i = 0; i < IDs_to_remove.size(); i++) { _remove_from_search_tree(& (_points[IDs_to_remove[i]])); } // insertion into tree new_IDs.resize(0); for (unsigned int i = 0; i < new_positions.size(); i++) { Point * new_point = _available_points.top(); _available_points.pop(); // set the point's coordinate new_point->coord = new_positions[i]; // now find it's neighbour in the search tree _insert_into_search_tree(new_point); // record the ID new_IDs.push_back(_ID(new_point)); } // the above statement labels certain points as needing "review" -- // deal with them... _deal_with_points_to_review(); }
unsigned int fastjet::ClosestPair2D::size | ( | ) | [inline, virtual] |
Implements fastjet::ClosestPair2DBase.
Definition at line 226 of file ClosestPair2D.hh.
{ return _points.size() - _available_points.size(); }
std::stack<Point *> fastjet::ClosestPair2D::_available_points [private] |
Definition at line 135 of file ClosestPair2D.hh.
unsigned int fastjet::ClosestPair2D::_cp_search_range [private] |
Definition at line 179 of file ClosestPair2D.hh.
std::auto_ptr<MinHeap> fastjet::ClosestPair2D::_heap [private] |
Definition at line 133 of file ClosestPair2D.hh.
Coord2D fastjet::ClosestPair2D::_left_corner [private] |
pieces needed for converting coordinates to integer
Definition at line 171 of file ClosestPair2D.hh.
const unsigned int fastjet::ClosestPair2D::_nshift = 3 [static, private] |
Definition at line 103 of file ClosestPair2D.hh.
std::vector<Point> fastjet::ClosestPair2D::_points [private] |
Definition at line 134 of file ClosestPair2D.hh.
std::vector<Point *> fastjet::ClosestPair2D::_points_under_review [private] |
points that are "under review" in some way
Definition at line 138 of file ClosestPair2D.hh.
double fastjet::ClosestPair2D::_range [private] |
Definition at line 172 of file ClosestPair2D.hh.
triplet<unsigned int> fastjet::ClosestPair2D::_rel_shifts [private] |
Definition at line 177 of file ClosestPair2D.hh.
const unsigned int fastjet::ClosestPair2D::_remove_heap_entry = 1 [static, private] |
Definition at line 141 of file ClosestPair2D.hh.
const unsigned int fastjet::ClosestPair2D::_review_heap_entry = 2 [static, private] |
Definition at line 142 of file ClosestPair2D.hh.
const unsigned int fastjet::ClosestPair2D::_review_neighbour = 4 [static, private] |
Definition at line 143 of file ClosestPair2D.hh.
triplet<unsigned int> fastjet::ClosestPair2D::_shifts [private] |
Definition at line 176 of file ClosestPair2D.hh.
triplet<std::auto_ptr<Tree> > fastjet::ClosestPair2D::_trees [private] |
Definition at line 132 of file ClosestPair2D.hh.