00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #ifndef __FASTJET_CLOSESTPAIR2D__HH__
00032 #define __FASTJET_CLOSESTPAIR2D__HH__
00033
00034 #include<vector>
00035 #include<stack>
00036 #include<iostream>
00037 #include "fastjet/internal/ClosestPair2DBase.hh"
00038 #include "fastjet/internal/SearchTree.hh"
00039 #include "fastjet/internal/MinHeap.hh"
00040
00041 FASTJET_BEGIN_NAMESPACE
00042
00043
00046 class ClosestPair2D : public ClosestPair2DBase {
00047 public:
00052 ClosestPair2D(const std::vector<Coord2D> & positions,
00053 const Coord2D & left_corner, const Coord2D & right_corner) {
00054 _initialize(positions, left_corner, right_corner, positions.size());
00055 };
00056
00059 ClosestPair2D(const std::vector<Coord2D> & positions,
00060 const Coord2D & left_corner, const Coord2D & right_corner,
00061 const unsigned int max_size) {
00062 _initialize(positions, left_corner, right_corner, max_size);
00063 };
00064
00067 void closest_pair(unsigned int & ID1, unsigned int & ID2,
00068 double & distance2) const;
00069
00071 void remove(unsigned int ID);
00072
00075 unsigned int insert(const Coord2D &);
00076
00079 virtual unsigned int replace(unsigned int ID1, unsigned int ID2,
00080 const Coord2D & position);
00081
00084 virtual void replace_many(const std::vector<unsigned int> & IDs_to_remove,
00085 const std::vector<Coord2D> & new_positions,
00086 std::vector<unsigned int> & new_IDs);
00087
00088
00089 inline void print_tree_depths(std::ostream & outdev) const {
00090 outdev << _trees[0]->max_depth() << " "
00091 << _trees[1]->max_depth() << " "
00092 << _trees[2]->max_depth() << "\n";
00093 };
00094
00095 unsigned int size();
00096
00097 private:
00098
00099 void _initialize(const std::vector<Coord2D> & positions,
00100 const Coord2D & left_corner, const Coord2D & right_corner,
00101 const unsigned int max_size);
00102
00103 static const unsigned int _nshift = 3;
00104
00105 class Point;
00106
00109 template<class T> class triplet {
00110 public:
00111 inline const T & operator[](unsigned int i) const {return _contents[i];};
00112 inline T & operator[](unsigned int i) {return _contents[i];};
00113 private:
00114 T _contents[_nshift];
00115 };
00116
00117
00119 class Shuffle {
00120 public:
00121 unsigned int x, y;
00122 Point * point;
00123 bool operator<(const Shuffle &) const;
00124 void operator+=(unsigned int shift) {x += shift; y+= shift;};
00125 };
00126
00127 typedef SearchTree<Shuffle> Tree;
00128 typedef Tree::circulator circulator;
00129 typedef Tree::const_circulator const_circulator;
00130
00131
00132 triplet<std::auto_ptr<Tree> > _trees;
00133 std::auto_ptr<MinHeap> _heap;
00134 std::vector<Point> _points;
00135 std::stack<Point *> _available_points;
00136
00138 std::vector<Point *> _points_under_review;
00139
00140
00141 static const unsigned int _remove_heap_entry = 1;
00142 static const unsigned int _review_heap_entry = 2;
00143 static const unsigned int _review_neighbour = 4;
00144
00148 void _add_label(Point * point, unsigned int review_flag);
00149
00153 void _set_label(Point * point, unsigned int review_flag);
00154
00159 void _deal_with_points_to_review();
00160
00162 void _remove_from_search_tree(Point * point_to_remove);
00163
00165 void _insert_into_search_tree(Point * new_point);
00166
00168 void _point2shuffle(Point & , Shuffle & , unsigned int shift);
00169
00171 Coord2D _left_corner;
00172 double _range;
00173
00174 int _ID(const Point *) const;
00175
00176 triplet<unsigned int> _shifts;
00177 triplet<unsigned int> _rel_shifts;
00178
00179 unsigned int _cp_search_range;
00180 };
00181
00182
00183
00185 class ClosestPair2D::Point {
00186 public:
00188 Coord2D coord;
00190 Point * neighbour;
00192 double neighbour_dist2;
00194 triplet<circulator> circ;
00195
00197 unsigned int review_flag;
00198
00200 double distance2(const Point & other) const {
00201 return coord.distance2(other.coord);
00202 };
00203
00205
00206 };
00207
00208
00209
00212 inline bool floor_ln2_less(unsigned x, unsigned y) {
00213 if (x>y) return false;
00214 return (x < (x^y));
00215 }
00216
00217
00218
00220 inline int ClosestPair2D::_ID(const Point * point) const {
00221 return point - &(_points[0]);
00222 }
00223
00224
00225
00226 inline unsigned int ClosestPair2D::size() {
00227 return _points.size() - _available_points.size();
00228 }
00229
00230
00231
00232 FASTJET_END_NAMESPACE
00233
00234 #endif // __FASTJET_CLOSESTPAIR2D__HH__