fastjet::ClusterSequence Class Reference
deals with clustering More...
#include <ClusterSequence.hh>
Classes | |
| struct | BriefJet |
| the fundamental structure which contains the minimal info about a jet, as needed for our plain N^2 algorithm -- the idea is to put all info that will be accessed N^2 times into an array of BriefJets. More... | |
| struct | EEBriefJet |
| fundamental structure for e+e- clustering More... | |
| class | Extras |
| a class intended to serve as a base in case a plugin needs to associate extra information with a ClusterSequence (see SISConePlugin. More... | |
| class | GenBriefJet |
| class to help with a generic clustering sequence More... | |
| struct | history_element |
| a single element in the clustering history (see vector _history below). More... | |
| struct | Tile |
| The fundamental structures to be used for the tiled N^2 algorithm (see CCN27-44 for some discussion of pattern of tiling). More... | |
| class | TiledJet |
| structure analogous to BriefJet, but with the extra information needed for dealing with tiles More... | |
Public Types | |
| enum | JetType { Invalid = -3, InexistentParent = -2, BeamJet = -1 } |
Public Member Functions | |
| ClusterSequence () | |
| default constructor | |
| template<class L > | |
| ClusterSequence (const std::vector< L > &pseudojets, const double &R=1.0, const Strategy &strategy=Best, const bool &writeout_combinations=false) | |
| create a clustersequence starting from the supplied set of pseudojets and clustering them with the long-invariant kt algorithm (E-scheme recombination) with the supplied value for R. | |
| template<class L > | |
| ClusterSequence (const std::vector< L > &pseudojets, const JetDefinition &jet_def, const bool &writeout_combinations=false) | |
| create a clustersequence starting from the supplied set of pseudojets and clustering them with jet definition specified by jet_def (which also specifies the clustering strategy) | |
| virtual | ~ClusterSequence () |
| std::vector< PseudoJet > | inclusive_jets (const double &ptmin=0.0) const |
| return a vector of all jets (in the sense of the inclusive algorithm) with pt >= ptmin. | |
| int | n_exclusive_jets (const double &dcut) const |
| return the number of jets (in the sense of the exclusive algorithm) that would be obtained when running the algorithm with the given dcut. | |
| std::vector< PseudoJet > | exclusive_jets (const double &dcut) const |
| return a vector of all jets (in the sense of the exclusive algorithm) that would be obtained when running the algorithm with the given dcut. | |
| std::vector< PseudoJet > | exclusive_jets (const int &njets) const |
| return a vector of all jets when the event is clustered (in the exclusive sense) to exactly njets. | |
| double | exclusive_dmerge (const int &njets) const |
| return the dmin corresponding to the recombination that went from n+1 to n jets (sometimes known as d_{n n+1}). | |
| double | exclusive_dmerge_max (const int &njets) const |
| return the maximum of the dmin encountered during all recombinations up to the one that led to an n-jet final state; identical to exclusive_dmerge, except in cases where the dmin do not increase monotonically. | |
| double | exclusive_ymerge (int njets) const |
| return the ymin corresponding to the recombination that went from n+1 to n jets (sometimes known as y_{n n+1}). | |
| double | exclusive_ymerge_max (int njets) const |
| same as exclusive_dmerge_max, but normalised to squared total energy | |
| int | n_exclusive_jets_ycut (double ycut) const |
| the number of exclusive jets at the given ycut | |
| std::vector< PseudoJet > | exclusive_jets_ycut (double ycut) const |
| the exclusive jets obtained at the given ycut | |
| std::vector< PseudoJet > | exclusive_subjets (const PseudoJet &jet, const double &dcut) const |
| return a vector of all subjets of the current jet (in the sense of the exclusive algorithm) that would be obtained when running the algorithm with the given dcut. | |
| int | n_exclusive_subjets (const PseudoJet &jet, const double &dcut) const |
| return the size of exclusive_subjets(. | |
| std::vector< PseudoJet > | exclusive_subjets (const PseudoJet &jet, int nsub) const |
| return the list of subjets obtained by unclustering the supplied jet down to n subjets (or all constituents if there are fewer than n). | |
| double | exclusive_subdmerge (const PseudoJet &jet, int nsub) const |
| return the dij that was present in the merging nsub+1 -> nsub subjets inside this jet. | |
| double | exclusive_subdmerge_max (const PseudoJet &jet, int nsub) const |
| return the maximum dij that occurred in the whole event at the stage that the nsub+1 -> nsub merge of subjets occurred inside this jet. | |
| double | Q () const |
| returns the sum of all energies in the event (relevant mainly for e+e-) | |
| double | Q2 () const |
| return Q()^2 | |
| bool | object_in_jet (const PseudoJet &object, const PseudoJet &jet) const |
| returns true iff the object is included in the jet. | |
| bool | has_parents (const PseudoJet &jet, PseudoJet &parent1, PseudoJet &parent2) const |
| if the jet has parents in the clustering, it returns true and sets parent1 and parent2 equal to them. | |
| bool | has_child (const PseudoJet &jet, PseudoJet &child) const |
| if the jet has a child then return true and give the child jet otherwise return false and set the child to zero | |
| bool | has_child (const PseudoJet &jet, const PseudoJet *&childp) const |
| Version of has_child that sets a pointer to the child if the child exists;. | |
| bool | has_partner (const PseudoJet &jet, PseudoJet &partner) const |
| if this jet has a child (and so a partner) return true and give the partner, otherwise return false and set the partner to zero | |
| std::vector< PseudoJet > | constituents (const PseudoJet &jet) const |
| return a vector of the particles that make up jet | |
| void | print_jets_for_root (const std::vector< PseudoJet > &jets, std::ostream &ostr=std::cout) const |
| output the supplied vector of jets in a format that can be read by an appropriate root script; the format is: jet-n jet-px jet-py jet-pz jet-E particle-n particle-rap particle-phi particle-pt particle-n particle-rap particle-phi particle-pt | |
| void | print_jets_for_root (const std::vector< PseudoJet > &jets, const std::string &filename, const std::string &comment="") const |
| print jets for root to the file labelled filename, with an optional comment at the beginning | |
| void | add_constituents (const PseudoJet &jet, std::vector< PseudoJet > &subjet_vector) const |
| add on to subjet_vector the constituents of jet (for internal use mainly) | |
| Strategy | strategy_used () const |
| return the enum value of the strategy used to cluster the event | |
| std::string | strategy_string () const |
| const JetDefinition & | jet_def () const |
| return a reference to the jet definition | |
| double | jet_scale_for_algorithm (const PseudoJet &jet) const |
| returns the scale associated with a jet as required for this clustering algorithm (kt^2 for the kt-algorithm, 1 for the Cambridge algorithm). | |
| void | plugin_record_ij_recombination (int jet_i, int jet_j, double dij, int &newjet_k) |
| record the fact that there has been a recombination between jets()[jet_i] and jets()[jet_k], with the specified dij, and return the index (newjet_k) allocated to the new jet, whose momentum is assumed to be the 4-vector sum of that of jet_i and jet_j | |
| void | plugin_record_ij_recombination (int jet_i, int jet_j, double dij, const PseudoJet &newjet, int &newjet_k) |
| as for the simpler variant of plugin_record_ij_recombination, except that the new jet is attributed the momentum and user_index of newjet | |
| void | plugin_record_iB_recombination (int jet_i, double diB) |
| record the fact that there has been a recombination between jets()[jet_i] and the beam, with the specified diB; when looking for inclusive jets, any iB recombination will returned to the user as a jet. | |
| void | plugin_associate_extras (std::auto_ptr< Extras > extras_in) |
| the plugin can associated some extra information with the ClusterSequence object by calling this function | |
| bool | plugin_activated () const |
| returns true when the plugin is allowed to run the show. | |
| const Extras * | extras () const |
| returns a pointer to the extras object (may be null) | |
| template<class GBJ > | |
| void | plugin_simple_N2_cluster () |
| allows a plugin to run a templated clustering (nearest-neighbour heuristic) | |
| const std::vector< PseudoJet > & | jets () const |
| allow the user to access the internally stored _jets() array, which contains both the initial particles and the various intermediate and final stages of recombination. | |
| const std::vector < history_element > & | history () const |
| allow the user to access the raw internal history. | |
| unsigned int | n_particles () const |
| returns the number of particles that were provided to the clustering algorithm (helps the user find their way around the history and jets objects if they weren't paying attention beforehand). | |
| std::vector< int > | particle_jet_indices (const std::vector< PseudoJet > &) const |
| returns a vector of size n_particles() which indicates, for each of the initial particles (in the order in which they were supplied), which of the supplied jets it belongs to; if it does not belong to any of the supplied jets, the index is set to -1; | |
| std::vector< int > | unique_history_order () const |
| routine that returns an order in which to read the history such that clusterings that lead to identical jet compositions but different histories (because of degeneracies in the clustering order) will have matching constituents for each matching entry in the unique_history_order. | |
| std::vector< PseudoJet > | unclustered_particles () const |
| return the set of particles that have not been clustered. | |
| void | transfer_from_sequence (ClusterSequence &from_seq) |
| transfer the sequence contained in other_seq into our own; any plugin "extras" contained in the from_seq will be lost from there. | |
| template<> | |
| void | _bj_set_jetinfo (EEBriefJet *const jetA, const int _jets_index) const |
| template<> | |
| double | _bj_dist (const EEBriefJet *const jeta, const EEBriefJet *const jetb) const |
Static Public Member Functions | |
| static void | set_jet_algorithm (JetAlgorithm jet_algorithm) |
| set the default (static) jet finder across all current and future ClusterSequence objects -- deprecated and obsolescent (i.e. | |
| static void | set_jet_finder (JetAlgorithm jet_algorithm) |
| same as above for backward compatibility | |
Protected Member Functions | |
| bool | _potentially_valid (const PseudoJet &jet) const |
| returns true if the jet has a history index contained within the range of this CS | |
| template<class L > | |
| void | _transfer_input_jets (const std::vector< L > &pseudojets) |
| transfer the vector<L> of input jets into our own vector<PseudoJet> _jets (with some reserved space for future growth). | |
| void | _initialise_and_run (const JetDefinition &jet_def, const bool &writeout_combinations) |
| This is the routine that will do all the initialisation and then run the clustering (may be called by various constructors). | |
| void | _initialise_and_run (const double &R, const Strategy &strategy, const bool &writeout_combinations) |
| This is an alternative routine for initialising and running the clustering, provided for legacy purposes. | |
| void | _decant_options (const JetDefinition &jet_def, const bool &writeout_combinations) |
| fills in the various member variables with "decanted" options from the jet_definition and writeout_combinations variables | |
| void | _fill_initial_history () |
| fill out the history (and jet cross refs) related to the initial set of jets (assumed already to have been "transferred"), without any clustering | |
| void | _do_ij_recombination_step (const int &jet_i, const int &jet_j, const double &dij, int &newjet_k) |
| carry out the recombination between the jets numbered jet_i and jet_j, at distance scale dij; return the index newjet_k of the result of the recombination of i and j. | |
| void | _do_iB_recombination_step (const int &jet_i, const double &diB) |
| carry out an recombination step in which _jets[jet_i] merges with the beam, | |
| void | get_subhist_set (std::set< const history_element * > &subhist, const PseudoJet &jet, double dcut, int maxjet) const |
set subhist to be a set pointers to history entries corresponding to the subjets of this jet; one stops going working down through the subjets either when
| |
Protected Attributes | |
| JetDefinition | _jet_def |
| std::vector< PseudoJet > | _jets |
| This contains the physical PseudoJets; for each PseudoJet one can find the corresponding position in the _history by looking at _jets[i].cluster_hist_index(). | |
| std::vector< history_element > | _history |
| this vector will contain the branching history; for each stage, _history[i].jetp_index indicates where to look in the _jets vector to get the physical PseudoJet. | |
| bool | _writeout_combinations |
| int | _initial_n |
| double | _Rparam |
| double | _R2 |
| double | _invR2 |
| double | _Qtot |
| Strategy | _strategy |
| JetAlgorithm | _jet_algorithm |
Static Protected Attributes | |
| static JetAlgorithm | _default_jet_algorithm = kt_algorithm |
Private Types | |
| typedef std::pair< int, int > | TwoVertices |
| typedef std::pair< double, TwoVertices > | DijEntry |
| typedef std::multimap< double, TwoVertices > | DistMap |
Private Member Functions | |
| void | _really_dumb_cluster () |
| Run the clustering in a very slow variant of the N^3 algorithm. | |
| void | _delaunay_cluster () |
| Run the clustering using a Hierarchical Delaunay triangulation and STL maps to achieve O(N*ln N) behaviour. | |
| template<class BJ > | |
| void | _simple_N2_cluster () |
| void | _tiled_N2_cluster () |
| run a tiled clustering | |
| void | _faster_tiled_N2_cluster () |
| run a tiled clustering | |
| void | _minheap_faster_tiled_N2_cluster () |
| run a tiled clustering, with our minheap for keeping track of the smallest dij | |
| void | _CP2DChan_cluster () |
| a 4pi variant of the closest pair clustering | |
| void | _CP2DChan_cluster_2pi2R () |
| a variant of the closest pair clustering which uses a region of size 2pi+2R in phi. | |
| void | _CP2DChan_cluster_2piMultD () |
| a variant of the closest pair clustering which uses a region of size 2pi + 2*0.3 and then carries on with 2pi+2*R | |
| void | _CP2DChan_limited_cluster (double D) |
| clusters only up to a distance Dlim -- does not deal with "inclusive" jets -- these are left to some other part of the program | |
| void | _do_Cambridge_inclusive_jets () |
| void | _add_step_to_history (const int &step_number, const int &parent1, const int &parent2, const int &jetp_index, const double &dij) |
| void | _extract_tree_children (int pos, std::valarray< bool > &, const std::valarray< int > &, std::vector< int > &) const |
| internal routine associated with the construction of the unique history order (following children in the tree) | |
| void | _extract_tree_parents (int pos, std::valarray< bool > &, const std::valarray< int > &, std::vector< int > &) const |
| internal routine associated with the construction of the unique history order (following parents in the tree) | |
| void | _add_ktdistance_to_map (const int &ii, DistMap &DijMap, const DynamicNearestNeighbours *DNN) |
| currently used only in the Voronoi based code | |
| void | _print_banner () |
| for making sure the user knows what it is they're running... | |
| template<class J > | |
| void | _bj_set_jetinfo (J *const jet, const int _jets_index) const |
| set the kinematic and labelling info for jeta so that it corresponds to _jets[_jets_index] | |
| void | _bj_remove_from_tiles (TiledJet *const jet) const |
| "remove" this jet, which implies updating links of neighbours and perhaps modifying the tile structure | |
| template<class J > | |
| double | _bj_dist (const J *const jeta, const J *const jetb) const |
| return the distance between two BriefJet objects | |
| template<class J > | |
| double | _bj_diJ (const J *const jeta) const |
| template<class J > | |
| J * | _bj_of_hindex (const int hist_index, J *const head, J *const tail) const |
| for testing purposes only: if in the range head--tail-1 there is a a jet which corresponds to hist_index in the history, then return a pointer to that jet; otherwise return tail. | |
| template<class J > | |
| void | _bj_set_NN_nocross (J *const jeta, J *const head, const J *const tail) const |
| updates (only towards smaller distances) the NN for jeta without checking whether in the process jeta itself might be a new NN of one of the jets being scanned -- span the range head to tail-1 with assumption that jeta is not contained in that range | |
| template<class J > | |
| void | _bj_set_NN_crosscheck (J *const jeta, J *const head, const J *const tail) const |
| reset the NN for jeta and DO check whether in the process jeta itself might be a new NN of one of the jets being scanned -- span the range head to tail-1 with assumption that jeta is not contained in that range | |
| int | _tile_index (int ieta, int iphi) const |
| int | _tile_index (const double &eta, const double &phi) const |
| return the tile index corresponding to the given eta,phi point | |
| void | _tj_set_jetinfo (TiledJet *const jet, const int _jets_index) |
| void | _bj_remove_from_tiles (TiledJet *const jet) |
| void | _initialise_tiles () |
Set up the tiles:
| |
| void | _print_tiles (TiledJet *briefjets) const |
| output the contents of the tiles | |
| void | _add_neighbours_to_tile_union (const int tile_index, std::vector< int > &tile_union, int &n_near_tiles) const |
| void | _add_untagged_neighbours_to_tile_union (const int tile_index, std::vector< int > &tile_union, int &n_near_tiles) |
| void | _simple_N2_cluster_BriefJet () |
| to help instantiation (fj 2.4.0; did not quite work on gcc 33 and os x 10.3?) | |
| void | _simple_N2_cluster_EEBriefJet () |
| to avoid issues with template instantiation (OS X 10.3, gcc 3.3) | |
Private Attributes | |
| bool | _plugin_activated |
| std::auto_ptr< Extras > | _extras |
| std::vector< Tile > | _tiles |
| double | _tiles_eta_min |
| double | _tiles_eta_max |
| double | _tile_size_eta |
| double | _tile_size_phi |
| int | _n_tiles_phi |
| int | _tiles_ieta_min |
| int | _tiles_ieta_max |
Static Private Attributes | |
| static bool | _first_time = true |
| will be set by default to be true for the first run | |
| static int | _n_exclusive_warnings = 0 |
| record the number of warnings provided about the exclusive algorithm -- so that we don't print it out more than a few times. | |
| static const int | n_tile_neighbours = 9 |
| number of neighbours that a tile will have (rectangular geometry gives 9 neighbours). | |
Detailed Description
deals with clusteringDefinition at line 66 of file ClusterSequence.hh.
Member Typedef Documentation
typedef std::pair<double,TwoVertices> fastjet::ClusterSequence::DijEntry [private] |
Definition at line 601 of file ClusterSequence.hh.
typedef std::multimap<double,TwoVertices> fastjet::ClusterSequence::DistMap [private] |
Definition at line 602 of file ClusterSequence.hh.
typedef std::pair<int,int> fastjet::ClusterSequence::TwoVertices [private] |
Definition at line 600 of file ClusterSequence.hh.
Member Enumeration Documentation
Definition at line 408 of file ClusterSequence.hh.
00408 {Invalid=-3, InexistentParent = -2, BeamJet = -1};
Constructor & Destructor Documentation
| fastjet::ClusterSequence::ClusterSequence | ( | ) | [inline] |
| fastjet::ClusterSequence::ClusterSequence | ( | const std::vector< L > & | pseudojets, | |
| const double & | R = 1.0, |
|||
| const Strategy & | strategy = Best, |
|||
| const bool & | writeout_combinations = false | |||
| ) | [inline] |
create a clustersequence starting from the supplied set of pseudojets and clustering them with the long-invariant kt algorithm (E-scheme recombination) with the supplied value for R.
If strategy=DumbN3 a very stupid N^3 algorithm is used for the clustering; otherwise strategy = NlnN* uses cylinders algorithms with some number of pi coverage. If writeout_combinations=true a summary of the recombination sequence is written out
Definition at line 796 of file ClusterSequence.hh.
References _initialise_and_run(), and _transfer_input_jets().
00800 { 00801 00802 // transfer the initial jets (type L) into our own array 00803 _transfer_input_jets(pseudojets); 00804 00805 // run the clustering 00806 _initialise_and_run(R,strategy,writeout_combinations); 00807 }
| fastjet::ClusterSequence::ClusterSequence | ( | const std::vector< L > & | pseudojets, | |
| const JetDefinition & | jet_def, | |||
| const bool & | writeout_combinations = false | |||
| ) | [inline] |
create a clustersequence starting from the supplied set of pseudojets and clustering them with jet definition specified by jet_def (which also specifies the clustering strategy)
constructor of a jet-clustering sequence from a vector of four-momenta, with the jet definition specified by jet_def
Definition at line 813 of file ClusterSequence.hh.
References _initialise_and_run(), and _transfer_input_jets().
00816 { 00817 00818 // transfer the initial jets (type L) into our own array 00819 _transfer_input_jets(pseudojets); 00820 00821 // run the clustering 00822 _initialise_and_run(jet_def,writeout_combinations); 00823 }
| fastjet::ClusterSequence::~ClusterSequence | ( | ) | [virtual] |
Member Function Documentation
| void fastjet::ClusterSequence::_add_ktdistance_to_map | ( | const int & | ii, | |
| DistMap & | DijMap, | |||
| const DynamicNearestNeighbours * | DNN | |||
| ) | [private] |
currently used only in the Voronoi based code
Add the current kt distance for particle ii to the map (DijMap) using information from the DNN object.
Work as follows:
. if the kt is zero then it's nearest neighbour is taken to be the the beam jet and the distance is zero.
. if cylinder distance to nearest neighbour > _Rparam then it is yiB that is smallest and this is added to map.
. otherwise if the nearest neighbour jj has a larger kt then add dij to the map.
. otherwise do nothing
Definition at line 220 of file ClusterSequence_Delaunay.cc.
References _invR2, _jets, jet_scale_for_algorithm(), fastjet::DynamicNearestNeighbours::NearestNeighbourDistance(), and fastjet::DynamicNearestNeighbours::NearestNeighbourIndex().
Referenced by _delaunay_cluster().
00223 { 00224 00225 double yiB = jet_scale_for_algorithm(_jets[ii]); 00226 if (yiB == 0.0) { 00227 // in this case convention is that we do not worry about distances 00228 // but directly state that nearest neighbour is beam 00229 DijMap.insert(DijEntry(yiB, TwoVertices(ii,-1))); 00230 } else { 00231 double DeltaR2 = DNN->NearestNeighbourDistance(ii) * _invR2; 00232 // Logic of following bit is: only add point to map if it has 00233 // smaller kt2 than nearest neighbour j (if it has larger kt, 00234 // then: either it is j's nearest neighbour and then we will 00235 // include dij when we come to j; or it is not j's nearest 00236 // neighbour and j will recombine with someone else). 00237 00238 // If DeltaR2 > 1.0 then in any case it will recombine with beam rather 00239 // than with any neighbours. 00240 // (put general normalisation here at some point) 00241 if (DeltaR2 > 1.0) { 00242 DijMap.insert(DijEntry(yiB, TwoVertices(ii,-1))); 00243 } else { 00244 double kt2i = jet_scale_for_algorithm(_jets[ii]); 00245 int jj = DNN->NearestNeighbourIndex(ii); 00246 if (kt2i <= jet_scale_for_algorithm(_jets[jj])) { 00247 double dij = DeltaR2 * kt2i; 00248 DijMap.insert(DijEntry(dij, TwoVertices(ii,jj))); 00249 } 00250 } 00251 } 00252 }
| void fastjet::ClusterSequence::_add_neighbours_to_tile_union | ( | const int | tile_index, | |
| std::vector< int > & | tile_union, | |||
| int & | n_near_tiles | |||
| ) | const [private] |
Referenced by _tiled_N2_cluster().
| void fastjet::ClusterSequence::_add_step_to_history | ( | const int & | step_number, | |
| const int & | parent1, | |||
| const int & | parent2, | |||
| const int & | jetp_index, | |||
| const double & | dij | |||
| ) | [private] |
Definition at line 936 of file ClusterSequence.cc.
References _history, _jets, _writeout_combinations, fastjet::ClusterSequence::history_element::child, fastjet::ClusterSequence::history_element::dij, Invalid, fastjet::ClusterSequence::history_element::jetp_index, fastjet::ClusterSequence::history_element::max_dij_so_far, fastjet::ClusterSequence::history_element::parent1, and fastjet::ClusterSequence::history_element::parent2.
Referenced by _do_iB_recombination_step(), and _do_ij_recombination_step().
00939 { 00940 00941 history_element element; 00942 element.parent1 = parent1; 00943 element.parent2 = parent2; 00944 element.jetp_index = jetp_index; 00945 element.child = Invalid; 00946 element.dij = dij; 00947 element.max_dij_so_far = max(dij,_history[_history.size()-1].max_dij_so_far); 00948 _history.push_back(element); 00949 00950 int local_step = _history.size()-1; 00951 assert(local_step == step_number); 00952 00953 assert(parent1 >= 0); 00954 _history[parent1].child = local_step; 00955 if (parent2 >= 0) {_history[parent2].child = local_step;} 00956 00957 // get cross-referencing right from PseudoJets 00958 if (jetp_index != Invalid) { 00959 assert(jetp_index >= 0); 00960 //cout << _jets.size() <<" "<<jetp_index<<"\n"; 00961 _jets[jetp_index].set_cluster_hist_index(local_step); 00962 } 00963 00964 if (_writeout_combinations) { 00965 cout << local_step << ": " 00966 << parent1 << " with " << parent2 00967 << "; y = "<< dij<<endl; 00968 } 00969 00970 }
| void fastjet::ClusterSequence::_add_untagged_neighbours_to_tile_union | ( | const int | tile_index, | |
| std::vector< int > & | tile_union, | |||
| int & | n_near_tiles | |||
| ) | [private] |
Referenced by _faster_tiled_N2_cluster(), and _minheap_faster_tiled_N2_cluster().
| double fastjet::ClusterSequence::_bj_diJ | ( | const J *const | jeta | ) | const [inline, private] |
Definition at line 863 of file ClusterSequence.hh.
Referenced by _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().
00863 { 00864 double kt2 = jet->kt2; 00865 if (jet->NN != NULL) {if (jet->NN->kt2 < kt2) {kt2 = jet->NN->kt2;}} 00866 return jet->NN_dist * kt2; 00867 }
| double fastjet::ClusterSequence::_bj_dist | ( | const EEBriefJet *const | jeta, | |
| const EEBriefJet *const | jetb | |||
| ) | const [inline] |
Definition at line 95 of file ClusterSequence_N2.cc.
References fastjet::ClusterSequence::EEBriefJet::nx, fastjet::ClusterSequence::EEBriefJet::ny, and fastjet::ClusterSequence::EEBriefJet::nz.
00097 { 00098 double dist = 1.0 00099 - jeta->nx*jetb->nx 00100 - jeta->ny*jetb->ny 00101 - jeta->nz*jetb->nz; 00102 dist *= 2; // distance is _2_*min(Ei^2,Ej^2)*(1-cos theta) 00103 //cout << "Dist = " << dist << ": " 00104 // << jeta->nx << " " 00105 // << jeta->ny << " " 00106 // << jeta->nz << " " 00107 // <<endl; 00108 return dist; 00109 }
| double fastjet::ClusterSequence::_bj_dist | ( | const J *const | jeta, | |
| const J *const | jetb | |||
| ) | const [inline, private] |
return the distance between two BriefJet objects
Definition at line 854 of file ClusterSequence.hh.
References fastjet::pi, and fastjet::twopi.
Referenced by _bj_set_NN_crosscheck(), _bj_set_NN_nocross(), _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().
00855 { 00856 double dphi = std::abs(jetA->phi - jetB->phi); 00857 double deta = (jetA->eta - jetB->eta); 00858 if (dphi > pi) {dphi = twopi - dphi;} 00859 return dphi*dphi + deta*deta; 00860 }
| J* fastjet::ClusterSequence::_bj_of_hindex | ( | const int | hist_index, | |
| J *const | head, | |||
| J *const | tail | |||
| ) | const [inline, private] |
for testing purposes only: if in the range head--tail-1 there is a a jet which corresponds to hist_index in the history, then return a pointer to that jet; otherwise return tail.
Definition at line 670 of file ClusterSequence.hh.
00673 { 00674 J * res; 00675 for(res = head; res<tail; res++) { 00676 if (_jets[res->_jets_index].cluster_hist_index() == hist_index) {break;} 00677 } 00678 return res; 00679 }
| void fastjet::ClusterSequence::_bj_remove_from_tiles | ( | TiledJet *const | jet | ) | [private] |
Definition at line 50 of file ClusterSequence_TiledN2.cc.
References _tiles, fastjet::ClusterSequence::Tile::head, fastjet::ClusterSequence::TiledJet::next, fastjet::ClusterSequence::TiledJet::previous, and fastjet::ClusterSequence::TiledJet::tile_index.
00050 { 00051 Tile * tile = & _tiles[jet->tile_index]; 00052 00053 if (jet->previous == NULL) { 00054 // we are at head of the tile, so reset it. 00055 // If this was the only jet on the tile then tile->head will now be NULL 00056 tile->head = jet->next; 00057 } else { 00058 // adjust link from previous jet in this tile 00059 jet->previous->next = jet->next; 00060 } 00061 if (jet->next != NULL) { 00062 // adjust backwards-link from next jet in this tile 00063 jet->next->previous = jet->previous; 00064 } 00065 }
| void fastjet::ClusterSequence::_bj_remove_from_tiles | ( | TiledJet *const | jet | ) | const [private] |
"remove" this jet, which implies updating links of neighbours and perhaps modifying the tile structure
Referenced by _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().
| void fastjet::ClusterSequence::_bj_set_jetinfo | ( | EEBriefJet *const | jetA, | |
| const int | _jets_index | |||
| ) | const [inline] |
Definition at line 54 of file ClusterSequence_N2.cc.
References _jet_algorithm, _jets, fastjet::ClusterSequence::EEBriefJet::_jets_index, _R2, _Rparam, fastjet::ee_genkt_algorithm, fastjet::ee_kt_algorithm, fastjet::JetDefinition::extra_param(), jet_def(), fastjet::ClusterSequence::EEBriefJet::kt2, fastjet::ClusterSequence::EEBriefJet::NN, fastjet::ClusterSequence::EEBriefJet::NN_dist, fastjet::norm(), fastjet::ClusterSequence::EEBriefJet::nx, fastjet::ClusterSequence::EEBriefJet::ny, and fastjet::ClusterSequence::EEBriefJet::nz.
00055 { 00056 00057 double E = _jets[_jets_index].E(); 00058 double scale = E*E; // the default energy scale for the kt alg 00059 double p = jet_def().extra_param(); // in case we're ee_genkt 00060 switch (_jet_algorithm) { 00061 case ee_kt_algorithm: 00062 assert(_Rparam > 2.0); // force this to be true! [not best place, but works] 00063 // recall that _invR2 is artificially set to 1 for this alg 00064 // so that we automatically have dij = scale * 2(1-cos theta_ij) 00065 // Normally, _Rparam should be automatically set to 4 from JetDefinition 00066 break; 00067 case ee_genkt_algorithm: 00068 if (p <= 0 && scale < 1e-300) scale = 1e-300; // same dodgy safety as genkt 00069 scale = pow(scale,p); 00070 break; 00071 default: 00072 throw Error("Unrecognised jet algorithm"); 00073 } 00074 jetA->kt2 = scale; // "kt2" might one day be renamed as "scale" or some such 00075 00076 double norm = _jets[_jets_index].modp2(); 00077 if (norm > 0) { 00078 norm = 1.0/sqrt(norm); 00079 jetA->nx = norm * _jets[_jets_index].px(); 00080 jetA->ny = norm * _jets[_jets_index].py(); 00081 jetA->nz = norm * _jets[_jets_index].pz(); 00082 } else { 00083 jetA->nx = 0.0; 00084 jetA->ny = 0.0; 00085 jetA->nz = 1.0; 00086 } 00087 jetA->_jets_index = _jets_index; 00088 // initialise NN info as well 00089 jetA->NN_dist = _R2; 00090 jetA->NN = NULL; 00091 }
| void fastjet::ClusterSequence::_bj_set_jetinfo | ( | J *const | jet, | |
| const int | _jets_index | |||
| ) | const [inline, private] |
set the kinematic and labelling info for jeta so that it corresponds to _jets[_jets_index]
Definition at line 839 of file ClusterSequence.hh.
References _jets, _R2, and jet_scale_for_algorithm().
00840 { 00841 jetA->eta = _jets[_jets_index].rap(); 00842 jetA->phi = _jets[_jets_index].phi_02pi(); 00843 jetA->kt2 = jet_scale_for_algorithm(_jets[_jets_index]); 00844 jetA->_jets_index = _jets_index; 00845 // initialise NN info as well 00846 jetA->NN_dist = _R2; 00847 jetA->NN = NULL; 00848 }
| void fastjet::ClusterSequence::_bj_set_NN_crosscheck | ( | J *const | jeta, | |
| J *const | head, | |||
| const J *const | tail | |||
| ) | const [inline, private] |
reset the NN for jeta and DO check whether in the process jeta itself might be a new NN of one of the jets being scanned -- span the range head to tail-1 with assumption that jeta is not contained in that range
Definition at line 901 of file ClusterSequence.hh.
References _bj_dist(), and _R2.
00902 { 00903 double NN_dist = _R2; 00904 J * NN = NULL; 00905 for (J * jetB = head; jetB != tail; jetB++) { 00906 double dist = _bj_dist(jet,jetB); 00907 if (dist < NN_dist) { 00908 NN_dist = dist; 00909 NN = jetB; 00910 } 00911 if (dist < jetB->NN_dist) { 00912 jetB->NN_dist = dist; 00913 jetB->NN = jet; 00914 } 00915 } 00916 jet->NN = NN; 00917 jet->NN_dist = NN_dist; 00918 }
| void fastjet::ClusterSequence::_bj_set_NN_nocross | ( | J *const | jeta, | |
| J *const | head, | |||
| const J *const | tail | |||
| ) | const [inline, private] |
updates (only towards smaller distances) the NN for jeta without checking whether in the process jeta itself might be a new NN of one of the jets being scanned -- span the range head to tail-1 with assumption that jeta is not contained in that range
Definition at line 873 of file ClusterSequence.hh.
References _bj_dist(), and _R2.
00874 { 00875 double NN_dist = _R2; 00876 J * NN = NULL; 00877 if (head < jet) { 00878 for (J * jetB = head; jetB != jet; jetB++) { 00879 double dist = _bj_dist(jet,jetB); 00880 if (dist < NN_dist) { 00881 NN_dist = dist; 00882 NN = jetB; 00883 } 00884 } 00885 } 00886 if (tail > jet) { 00887 for (J * jetB = jet+1; jetB != tail; jetB++) { 00888 double dist = _bj_dist(jet,jetB); 00889 if (dist < NN_dist) { 00890 NN_dist = dist; 00891 NN = jetB; 00892 } 00893 } 00894 } 00895 jet->NN = NN; 00896 jet->NN_dist = NN_dist; 00897 }
| void fastjet::ClusterSequence::_CP2DChan_cluster | ( | ) | [private] |
a 4pi variant of the closest pair clustering
Definition at line 223 of file ClusterSequence_CP2DChan.cc.
References _do_Cambridge_inclusive_jets(), _do_ij_recombination_step(), _invR2, _jet_algorithm, _jets, BeamJet, fastjet::cambridge_algorithm, fastjet::ClosestPair2D::closest_pair(), Invalid, fastjet::d0::inline_maths::min(), fastjet::ClosestPair2D::replace(), and fastjet::twopi.
Referenced by _initialise_and_run().
00223 { 00224 00225 if (_jet_algorithm != cambridge_algorithm) throw Error("_CP2DChan_cluster called for a jet-finder that is not the cambridge algorithm"); 00226 00227 unsigned int n = _jets.size(); 00228 00229 vector<MirrorInfo> coordIDs(2*n); // link from original to mirror indices 00230 vector<int> jetIDs(2*n); // link from mirror to original indices 00231 vector<Coord2D> coords(2*n); // our coordinates (and copies) 00232 00233 // start things off... 00234 double minrap = numeric_limits<double>::max(); 00235 double maxrap = -minrap; 00236 int coord_index = 0; 00237 for (unsigned i = 0; i < n; i++) { 00238 // separate out points with infinite rapidity 00239 if (_jets[i].E() == abs(_jets[i].pz()) && _jets[i].perp2() == 0.0) { 00240 coordIDs[i] = MirrorInfo(BeamJet,BeamJet); 00241 } else { 00242 coordIDs[i].orig = coord_index; 00243 coordIDs[i].mirror = coord_index+1; 00244 coords[coord_index] = Coord2D(_jets[i].rap(), _jets[i].phi_02pi()); 00245 coords[coord_index+1] = Coord2D(_jets[i].rap(), _jets[i].phi_02pi()+twopi); 00246 jetIDs[coord_index] = i; 00247 jetIDs[coord_index+1] = i; 00248 minrap = min(coords[coord_index].x,minrap); 00249 maxrap = max(coords[coord_index].x,maxrap); 00250 coord_index += 2; 00251 } 00252 } 00253 // label remaining "mirror info" as saying that there are no jets 00254 for (unsigned i = n; i < 2*n; i++) {coordIDs[i].orig = Invalid;} 00255 00256 // set them to their actual size... 00257 coords.resize(coord_index); 00258 00259 // establish limits (with some leeway on rapidity) 00260 Coord2D left_edge(minrap-1.0, 0.0); 00261 Coord2D right_edge(maxrap+1.0, 2*twopi); 00262 00263 // now create the closest pair search object 00264 ClosestPair2D cp(coords, left_edge, right_edge); 00265 00266 // and start iterating... 00267 vector<Coord2D> new_points(2); 00268 vector<unsigned int> cIDs_to_remove(4); 00269 vector<unsigned int> new_cIDs(2); 00270 do { 00271 // find out which pair we recombine 00272 unsigned int cID1, cID2; 00273 double distance2; 00274 cp.closest_pair(cID1,cID2,distance2); 00275 distance2 *= _invR2; 00276 00277 // if the closest distance > 1, we exit and go to "inclusive" stage 00278 if (distance2 > 1.0) {break;} 00279 00280 // do the recombination... 00281 int jet_i = jetIDs[cID1]; 00282 int jet_j = jetIDs[cID2]; 00283 int newjet_k; 00284 _do_ij_recombination_step(jet_i, jet_j, distance2, newjet_k); 00285 00286 // now prepare operations on CP structure 00287 cIDs_to_remove[0] = coordIDs[jet_i].orig; 00288 cIDs_to_remove[1] = coordIDs[jet_i].mirror; 00289 cIDs_to_remove[2] = coordIDs[jet_j].orig; 00290 cIDs_to_remove[3] = coordIDs[jet_j].mirror; 00291 new_points[0] = Coord2D(_jets[newjet_k].rap(),_jets[newjet_k].phi_02pi()); 00292 new_points[1] = Coord2D(_jets[newjet_k].rap(),_jets[newjet_k].phi_02pi()+twopi); 00293 // carry out the CP operation 00294 //cp.replace_many(cIDs_to_remove, new_points, new_cIDs); 00295 // remarkable the following is faster... 00296 new_cIDs[0] = cp.replace(cIDs_to_remove[0], cIDs_to_remove[2], new_points[0]); 00297 new_cIDs[1] = cp.replace(cIDs_to_remove[1], cIDs_to_remove[3], new_points[1]); 00298 // signal that the following jets no longer exist 00299 coordIDs[jet_i].orig = Invalid; 00300 coordIDs[jet_j].orig = Invalid; 00301 // and do bookkeeping for new jet 00302 coordIDs[newjet_k] = MirrorInfo(new_cIDs[0], new_cIDs[1]); 00303 jetIDs[new_cIDs[0]] = newjet_k; 00304 jetIDs[new_cIDs[1]] = newjet_k; 00305 00306 // if we've reached one jet we should exit... 00307 n--; 00308 if (n == 1) {break;} 00309 00310 } while(true); 00311 00312 00313 // now do "beam" recombinations 00314 //for (unsigned int jet_i = 0; jet_i < coordIDs.size(); jet_i++) { 00315 // // if we have a predeclared beam jet or a valid set of mirror 00316 // // coordinates, recombine then recombine this jet with the beam 00317 // if (coordIDs[jet_i].orig == BeamJet || coordIDs[jet_i].orig > 0) { 00318 // // diB = 1 _always_ (for the cambridge alg.) 00319 // _do_iB_recombination_step(jet_i, 1.0); 00320 // } 00321 //} 00322 00323 _do_Cambridge_inclusive_jets(); 00324 00325 //n = _history.size(); 00326 //for (unsigned int hist_i = 0; hist_i < n; hist_i++) { 00327 // if (_history[hist_i].child == Invalid) { 00328 // _do_iB_recombination_step(_history[hist_i].jetp_index, 1.0); 00329 // } 00330 //} 00331 00332 }
| void fastjet::ClusterSequence::_CP2DChan_cluster_2pi2R | ( | ) | [private] |
a variant of the closest pair clustering which uses a region of size 2pi+2R in phi.
Definition at line 193 of file ClusterSequence_CP2DChan.cc.
References _CP2DChan_limited_cluster(), _do_Cambridge_inclusive_jets(), _jet_algorithm, _Rparam, and fastjet::cambridge_algorithm.
Referenced by _CP2DChan_cluster_2piMultD(), and _initialise_and_run().
00193 { 00194 00195 if (_jet_algorithm != cambridge_algorithm) throw Error("CP2DChan clustering method called for a jet-finder that is not the cambridge algorithm"); 00196 00197 // run the clustering with mirror copies kept such that only things 00198 // within _Rparam of a border are mirrored 00199 _CP2DChan_limited_cluster(_Rparam); 00200 00201 // 00202 _do_Cambridge_inclusive_jets(); 00203 }
| void fastjet::ClusterSequence::_CP2DChan_cluster_2piMultD | ( | ) | [private] |
a variant of the closest pair clustering which uses a region of size 2pi + 2*0.3 and then carries on with 2pi+2*R
Definition at line 209 of file ClusterSequence_CP2DChan.cc.
References _CP2DChan_cluster_2pi2R(), _CP2DChan_limited_cluster(), _Rparam, and fastjet::d0::inline_maths::min().
Referenced by _initialise_and_run().
00209 { 00210 00211 // do a first run of clustering up to a small distance parameter, 00212 if (_Rparam >= 0.39) { 00213 _CP2DChan_limited_cluster(min(_Rparam/2,0.3)); 00214 } 00215 00216 // and then the final round of clustering 00217 _CP2DChan_cluster_2pi2R (); 00218 }
| void fastjet::ClusterSequence::_CP2DChan_limited_cluster | ( | double | D | ) | [private] |
clusters only up to a distance Dlim -- does not deal with "inclusive" jets -- these are left to some other part of the program
Definition at line 69 of file ClusterSequence_CP2DChan.cc.
References _do_ij_recombination_step(), _history, _initial_n, _invR2, _jets, fastjet::ClosestPair2D::closest_pair(), Invalid, fastjet::Private::make_mirror(), fastjet::d0::inline_maths::min(), fastjet::pi, and fastjet::ClosestPair2D::replace_many().
Referenced by _CP2DChan_cluster_2pi2R(), and _CP2DChan_cluster_2piMultD().
00069 { 00070 00071 unsigned int n = _initial_n; 00072 00073 vector<MirrorInfo> coordIDs(2*n); // coord IDs of a given jetID 00074 vector<int> jetIDs(2*n); // jet ID for a given coord ID 00075 vector<Coord2D> coords(2*n); // our coordinates (and copies) 00076 00077 // start things off... 00078 double minrap = numeric_limits<double>::max(); 00079 double maxrap = -minrap; 00080 int coord_index = -1; 00081 int n_active = 0; 00082 for (unsigned jet_i = 0; jet_i < _jets.size(); jet_i++) { 00083 00084 // skip jets that already have children or that have infinite 00085 // rapidity 00086 if (_history[_jets[jet_i].cluster_hist_index()].child != Invalid || 00087 (_jets[jet_i].E() == abs(_jets[jet_i].pz()) && 00088 _jets[jet_i].perp2() == 0.0) 00089 ) {continue;} 00090 00091 n_active++; 00092 00093 coordIDs[jet_i].orig = ++coord_index; 00094 coords[coord_index] = Coord2D(_jets[jet_i].rap(), _jets[jet_i].phi_02pi()); 00095 jetIDs[coord_index] = jet_i; 00096 minrap = min(coords[coord_index].x,minrap); 00097 maxrap = max(coords[coord_index].x,maxrap); 00098 00099 Coord2D mirror_point(coords[coord_index]); 00100 if (make_mirror(mirror_point, Dlim)) { 00101 coordIDs[jet_i].mirror = ++coord_index; 00102 coords[coord_index] = mirror_point; 00103 jetIDs[coord_index] = jet_i; 00104 } else { 00105 coordIDs[jet_i].mirror = Invalid; 00106 } 00107 } 00108 00109 // set them to their actual size... 00110 coords.resize(coord_index+1); 00111 00112 // establish limits (with some leeway on rapidity) 00113 Coord2D left_edge(minrap-1.0, -pi); 00114 Coord2D right_edge(maxrap+1.0, 3*pi); 00115 00116 //cerr << "minrap, maxrap = " << minrap << " " << maxrap << endl; 00117 00118 // now create the closest pair search object 00119 ClosestPair2D cp(coords, left_edge, right_edge); 00120 00121 // cross check to see what's going on... 00122 //cerr << "Tree depths before: "; 00123 //cp.print_tree_depths(cerr); 00124 00125 // and start iterating... 00126 vector<Coord2D> new_points(2); 00127 vector<unsigned int> cIDs_to_remove(4); 00128 vector<unsigned int> new_cIDs(2); 00129 00130 do { 00131 // find out which pair we recombine 00132 unsigned int cID1, cID2; 00133 double distance2; 00134 cp.closest_pair(cID1,cID2,distance2); 00135 00136 // if the closest distance > Dlim, we exit and go to "inclusive" stage 00137 if (distance2 > Dlim*Dlim) {break;} 00138 00139 // normalise distance as we like it 00140 distance2 *= _invR2; 00141 00142 // do the recombination... 00143 int jet_i = jetIDs[cID1]; 00144 int jet_j = jetIDs[cID2]; 00145 int newjet_k; 00146 _do_ij_recombination_step(jet_i, jet_j, distance2, newjet_k); 00147 00148 // don't bother with any further action if only one active particle 00149 // is left (also avoid closest-pair error [cannot remove last particle]). 00150 if (--n_active == 1) {break;} 00151 00152 // now prepare operations on CP structure 00153 cIDs_to_remove.resize(0); 00154 cIDs_to_remove.push_back(coordIDs[jet_i].orig); 00155 cIDs_to_remove.push_back(coordIDs[jet_j].orig); 00156 if (coordIDs[jet_i].mirror != Invalid) 00157 cIDs_to_remove.push_back(coordIDs[jet_i].mirror); 00158 if (coordIDs[jet_j].mirror != Invalid) 00159 cIDs_to_remove.push_back(coordIDs[jet_j].mirror); 00160 00161 Coord2D new_point(_jets[newjet_k].rap(),_jets[newjet_k].phi_02pi()); 00162 new_points.resize(0); 00163 new_points.push_back(new_point); 00164 if (make_mirror(new_point, Dlim)) new_points.push_back(new_point); 00165 00166 // carry out actions on search tree 00167 cp.replace_many(cIDs_to_remove, new_points, new_cIDs); 00168 00169 // now fill in info for new points... 00170 coordIDs[newjet_k].orig = new_cIDs[0]; 00171 jetIDs[new_cIDs[0]] = newjet_k; 00172 if (new_cIDs.size() == 2) { 00173 coordIDs[newjet_k].mirror = new_cIDs[1]; 00174 jetIDs[new_cIDs[1]] = newjet_k; 00175 } else {coordIDs[newjet_k].mirror = Invalid;} 00176 00178 //n_active--; 00179 //if (n_active == 1) {break;} 00180 00181 } while(true); 00182 00183 // cross check to see what's going on... 00184 //cerr << "Max tree depths after: "; 00185 //cp.print_tree_depths(cerr); 00186 00187 }
| void fastjet::ClusterSequence::_decant_options | ( | const JetDefinition & | jet_def, | |
| const bool & | writeout_combinations | |||
| ) | [protected] |
fills in the various member variables with "decanted" options from the jet_definition and writeout_combinations variables
Definition at line 226 of file ClusterSequence.cc.
References _invR2, _jet_algorithm, _jet_def, _plugin_activated, _print_banner(), _R2, _Rparam, _strategy, _writeout_combinations, fastjet::JetDefinition::jet_algorithm(), fastjet::JetDefinition::R(), and fastjet::JetDefinition::strategy().
Referenced by fastjet::ClusterSequenceActiveArea::_initialise_AA(), and _initialise_and_run().
00227 { 00228 00229 // let the user know what's going on 00230 _print_banner(); 00231 00232 // make a local copy of the jet definition (for future use?) 00233 _jet_def = jet_def; 00234 00235 _writeout_combinations = writeout_combinations; 00236 _jet_algorithm = jet_def.jet_algorithm(); 00237 _Rparam = jet_def.R(); _R2 = _Rparam*_Rparam; _invR2 = 1.0/_R2; 00238 _strategy = jet_def.strategy(); 00239 00240 // disallow interference from the plugin 00241 _plugin_activated = false; 00242 00243 }
| void fastjet::ClusterSequence::_delaunay_cluster | ( | ) | [private] |
Run the clustering using a Hierarchical Delaunay triangulation and STL maps to achieve O(N*ln N) behaviour.
There may be internally asserted assumptions about absence of points with coincident eta-phi coordinates.
Definition at line 58 of file ClusterSequence_Delaunay.cc.
References _add_ktdistance_to_map(), _do_iB_recombination_step(), _do_ij_recombination_step(), _jets, _Rparam, _strategy, BeamJet, fastjet::NlnN, fastjet::NlnN3pi, fastjet::NlnN4pi, fastjet::DynamicNearestNeighbours::RemoveCombinedAddCombination(), fastjet::DynamicNearestNeighbours::RemovePoint(), fastjet::EtaPhi::sanitize(), strategy_string(), fastjet::twopi, and fastjet::DynamicNearestNeighbours::Valid().
Referenced by _initialise_and_run().
00058 { 00059 00060 int n = _jets.size(); 00061 00062 vector<EtaPhi> points(n); // recall EtaPhi is just a typedef'd pair<double> 00063 for (int i = 0; i < n; i++) { 00064 points[i] = EtaPhi(_jets[i].rap(),_jets[i].phi_02pi()); 00065 points[i].sanitize(); // make sure things are in the right range 00066 } 00067 00068 // initialise our DNN structure with the set of points 00069 DynamicNearestNeighbours * DNN; 00070 #ifndef DROP_CGAL // strategy = NlnN* are not supported if we drop CGAL... 00071 bool verbose = false; 00072 bool ignore_nearest_is_mirror = (_Rparam < twopi); 00073 if (_strategy == NlnN4pi) { 00074 DNN = new Dnn4piCylinder(points,verbose); 00075 } else if (_strategy == NlnN3pi) { 00076 DNN = new Dnn3piCylinder(points,ignore_nearest_is_mirror,verbose); 00077 } else if (_strategy == NlnN) { 00078 DNN = new Dnn2piCylinder(points,ignore_nearest_is_mirror,verbose); 00079 } else 00080 #else 00081 if (_strategy == NlnN4pi || _strategy == NlnN3pi || _strategy == NlnN) { 00082 ostringstream err; 00083 err << "ERROR: Requested strategy "<<strategy_string()<<" but it is not"<<endl; 00084 err << " supported because FastJet was compiled without CGAL"<<endl; 00085 throw Error(err.str()); 00086 //assert(false); 00087 } 00088 #endif // DROP_CGAL 00089 { 00090 ostringstream err; 00091 err << "ERROR: Unrecognized value for strategy: "<<_strategy<<endl; 00092 assert(false); 00093 throw Error(err.str()); 00094 } 00095 00096 // We will find nearest neighbour for each vertex, and include 00097 // distance in map (NB DistMap is a typedef given in the .h file) 00098 DistMap DijMap; 00099 00100 // fill the map with the minimal (as far as we know) subset of Dij 00101 // distances (i.e. nearest neighbour ones). 00102 for (int ii = 0; ii < n; ii++) { 00103 _add_ktdistance_to_map(ii, DijMap, DNN); 00104 } 00105 00106 // run the clustering (go up to i=n-1, but then will stop half-way down, 00107 // when we reach that point -- it will be the final beam jet and there 00108 // will be no nearest neighbours to find). 00109 for (int i=0;i<n;i++) { 00110 // find nearest vertices 00111 // NB: skip cases where the point is not there anymore! 00112 TwoVertices SmallestDijPair; 00113 int jet_i, jet_j; 00114 double SmallestDij; 00115 bool Valid2; 00116 bool recombine_with_beam; 00117 do { 00118 SmallestDij = DijMap.begin()->first; 00119 SmallestDijPair = DijMap.begin()->second; 00120 jet_i = SmallestDijPair.first; 00121 jet_j = SmallestDijPair.second; 00122 // distance is immediately removed regardless of whether or not 00123 // it is used. 00124 // Some temporary testing code relating to problems with the gcc-3.2 compiler 00125 //cout << "got here and size is "<< DijMap.size()<< " and it is "<<SmallestDij <<"\n"; 00126 //cout << jet_i << " "<< jet_j<<"\n"; 00127 DijMap.erase(DijMap.begin()); 00128 //cout << "got beyond here\n"; 00129 00130 // need to "prime" the validity of jet_j in such a way that 00131 // if it corresponds to the beam then it is automatically valid. 00132 recombine_with_beam = (jet_j == BeamJet); 00133 if (!recombine_with_beam) {Valid2 = DNN->Valid(jet_j);} 00134 else {Valid2 = true;} 00135 00136 } while ( !DNN->Valid(jet_i) || !Valid2); 00137 00138 00139 // The following part acts just on jet momenta and on the history. 00140 // The action on the nearest-neighbour structures takes place 00141 // later (only if at least 2 jets are around). 00142 if (! recombine_with_beam) { 00143 int nn; // will be index of new jet 00144 _do_ij_recombination_step(jet_i, jet_j, SmallestDij, nn); 00145 //OBS // merge the two jets, add new jet, remove old ones 00146 //OBS _jets.push_back(_jets[jet_i] + _jets[jet_j]); 00147 //OBS 00148 //OBS int nn = _jets.size()-1; 00149 //OBS _jets[nn].set_cluster_hist_index(n+i); 00150 //OBS 00151 //OBS // get corresponding indices in history structure 00152 //OBS int hist_i = _jets[jet_i].cluster_hist_index(); 00153 //OBS int hist_j = _jets[jet_j].cluster_hist_index(); 00154 //OBS 00155 //OBS 00156 //OBS _add_step_to_history(n+i,min(hist_i,hist_j), max(hist_i,hist_j), 00157 //OBS _jets.size()-1, SmallestDij); 00158 00159 // add new point to points vector 00160 EtaPhi newpoint(_jets[nn].rap(), _jets[nn].phi_02pi()); 00161 newpoint.sanitize(); // make sure it is in correct range 00162 points.push_back(newpoint); 00163 } else { 00164 // recombine the jet with the beam 00165 _do_iB_recombination_step(jet_i, SmallestDij); 00166 //OBS _add_step_to_history(n+i,_jets[jet_i].cluster_hist_index(),BeamJet, 00167 //OBS Invalid, SmallestDij); 00168 } 00169 00170 // exit the loop because we do not want to look for nearest neighbours 00171 // etc. of zero partons 00172 if (i == n-1) {break;} 00173 00174 vector<int> updated_neighbours; 00175 if (! recombine_with_beam) { 00176 // update DNN 00177 int point3; 00178 DNN->RemoveCombinedAddCombination(jet_i, jet_j, 00179 points[points.size()-1], point3, 00180 updated_neighbours); 00181 // C++ beginners' comment: static_cast to unsigned int is necessary 00182 // to do away with warnings about type mismatch between point3 (int) 00183 // and points.size (unsigned int) 00184 if (static_cast<unsigned int> (point3) != points.size()-1) { 00185 throw Error("INTERNAL ERROR: point3 != points.size()-1");} 00186 } else { 00187 // update DNN 00188 DNN->RemovePoint(jet_i, updated_neighbours); 00189 } 00190 00191 // update map 00192 vector<int>::iterator it = updated_neighbours.begin(); 00193 for (; it != updated_neighbours.end(); ++it) { 00194 int ii = *it; 00195 _add_ktdistance_to_map(ii, DijMap, DNN); 00196 } 00197 00198 } // end clustering loop 00199 00200 // remember to clean up! 00201 delete DNN; 00202 }
| void fastjet::ClusterSequence::_do_Cambridge_inclusive_jets | ( | ) | [private] |
Definition at line 336 of file ClusterSequence_CP2DChan.cc.
References _do_iB_recombination_step(), _history, and Invalid.
Referenced by _CP2DChan_cluster(), and _CP2DChan_cluster_2pi2R().
00336 { 00337 unsigned int n = _history.size(); 00338 for (unsigned int hist_i = 0; hist_i < n; hist_i++) { 00339 if (_history[hist_i].child == Invalid) { 00340 _do_iB_recombination_step(_history[hist_i].jetp_index, 1.0); 00341 } 00342 } 00343 }
| void fastjet::ClusterSequence::_do_iB_recombination_step | ( | const int & | jet_i, | |
| const double & | diB | |||
| ) | [protected] |
carry out an recombination step in which _jets[jet_i] merges with the beam,
carries out the bookkeeping associated with the step of recombining jet_i with the beam
Definition at line 1113 of file ClusterSequence.cc.
References _add_step_to_history(), _history, _jets, BeamJet, and Invalid.
Referenced by _delaunay_cluster(), _do_Cambridge_inclusive_jets(), _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), _really_dumb_cluster(), _tiled_N2_cluster(), and fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history().
01114 { 01115 // get history index 01116 int newstep_k = _history.size(); 01117 01118 // recombine the jet with the beam 01119 _add_step_to_history(newstep_k,_jets[jet_i].cluster_hist_index(),BeamJet, 01120 Invalid, diB); 01121 01122 }
| void fastjet::ClusterSequence::_do_ij_recombination_step | ( | const int & | jet_i, | |
| const int & | jet_j, | |||
| const double & | dij, | |||
| int & | newjet_k | |||
| ) | [protected] |
carry out the recombination between the jets numbered jet_i and jet_j, at distance scale dij; return the index newjet_k of the result of the recombination of i and j.
carries out the bookkeeping associated with the step of recombining jet_i and jet_j (assuming a distance dij) and returns the index of the recombined jet, newjet_k.
Definition at line 1080 of file ClusterSequence.cc.
References _add_step_to_history(), _history, _jet_def, _jets, fastjet::d0::inline_maths::min(), and fastjet::JetDefinition::recombiner().
Referenced by _CP2DChan_cluster(), _CP2DChan_limited_cluster(), _delaunay_cluster(), _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), _really_dumb_cluster(), _tiled_N2_cluster(), and fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history().
01083 { 01084 01085 // create the new jet by recombining the first two 01086 PseudoJet newjet; 01087 _jet_def.recombiner()->recombine(_jets[jet_i], _jets[jet_j], newjet); 01088 _jets.push_back(newjet); 01089 // original version... 01090 //_jets.push_back(_jets[jet_i] + _jets[jet_j]); 01091 01092 // get its index 01093 newjet_k = _jets.size()-1; 01094 01095 // get history index 01096 int newstep_k = _history.size(); 01097 // and provide jet with the info 01098 _jets[newjet_k].set_cluster_hist_index(newstep_k); 01099 01100 // finally sort out the history 01101 int hist_i = _jets[jet_i].cluster_hist_index(); 01102 int hist_j = _jets[jet_j].cluster_hist_index(); 01103 01104 _add_step_to_history(newstep_k, min(hist_i, hist_j), max(hist_i,hist_j), 01105 newjet_k, dij); 01106 01107 }
| void fastjet::ClusterSequence::_extract_tree_children | ( | int | pos, | |
| std::valarray< bool > & | , | |||
| const std::valarray< int > & | , | |||
| std::vector< int > & | ||||
| ) | const [private] |
internal routine associated with the construction of the unique history order (following children in the tree)
Reimplemented in fastjet::ClusterSequenceActiveArea.
Referenced by unique_history_order().
| void fastjet::ClusterSequence::_extract_tree_parents | ( | int | pos, | |
| std::valarray< bool > & | , | |||
| const std::valarray< int > & | , | |||
| std::vector< int > & | ||||
| ) | const [private] |
internal routine associated with the construction of the unique history order (following parents in the tree)
Reimplemented in fastjet::ClusterSequenceActiveArea.
| void fastjet::ClusterSequence::_faster_tiled_N2_cluster | ( | ) | [private] |
run a tiled clustering
Definition at line 510 of file ClusterSequence_TiledN2.cc.
References _add_untagged_neighbours_to_tile_union(), _bj_diJ(), _bj_dist(), _bj_remove_from_tiles(), _do_iB_recombination_step(), _do_ij_recombination_step(), _initialise_tiles(), _invR2, _jets, fastjet::ClusterSequence::TiledJet::_jets_index, _R2, _tiles, _tj_set_jetinfo(), fastjet::ClusterSequence::Tile::begin_tiles, fastjet::ClusterSequence::TiledJet::diJ_posn, fastjet::ClusterSequence::Tile::end_tiles, fastjet::ClusterSequence::Tile::head, n_tile_neighbours, fastjet::ClusterSequence::TiledJet::next, fastjet::ClusterSequence::TiledJet::NN, fastjet::ClusterSequence::TiledJet::NN_dist, fastjet::ClusterSequence::Tile::tagged, and fastjet::ClusterSequence::TiledJet::tile_index.
Referenced by _initialise_and_run().
00510 { 00511 00512 _initialise_tiles(); 00513 00514 int n = _jets.size(); 00515 TiledJet * briefjets = new TiledJet[n]; 00516 TiledJet * jetA = briefjets, * jetB; 00517 TiledJet oldB; 00518 00519 00520 // will be used quite deep inside loops, but declare it here so that 00521 // memory (de)allocation gets done only once 00522 vector<int> tile_union(3*n_tile_neighbours); 00523 00524 // initialise the basic jet info 00525 for (int i = 0; i< n; i++) { 00526 _tj_set_jetinfo(jetA, i); 00527 //cout << i<<": "<<jetA->tile_index<<"\n"; 00528 jetA++; // move on to next entry of briefjets 00529 } 00530 TiledJet * head = briefjets; // a nicer way of naming start 00531 00532 // set up the initial nearest neighbour information 00533 vector<Tile>::const_iterator tile; 00534 for (tile = _tiles.begin(); tile != _tiles.end(); tile++) { 00535 // first do it on this tile 00536 for (jetA = tile->head; jetA != NULL; jetA = jetA->next) { 00537 for (jetB = tile->head; jetB != jetA; jetB = jetB->next) { 00538 double dist = _bj_dist(jetA,jetB); 00539 if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;} 00540 if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;} 00541 } 00542 } 00543 // then do it for RH tiles 00544 for (Tile ** RTile = tile->RH_tiles; RTile != tile->end_tiles; RTile++) { 00545 for (jetA = tile->head; jetA != NULL; jetA = jetA->next) { 00546 for (jetB = (*RTile)->head; jetB != NULL; jetB = jetB->next) { 00547 double dist = _bj_dist(jetA,jetB); 00548 if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;} 00549 if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;} 00550 } 00551 } 00552 } 00553 // no need to do it for LH tiles, since they are implicitly done 00554 // when we set NN for both jetA and jetB on the RH tiles. 00555 } 00556 00557 00558 // now create the diJ (where J is i's NN) table -- remember that 00559 // we differ from standard normalisation here by a factor of R2 00560 // (corrected for at the end). 00561 struct diJ_plus_link { 00562 double diJ; // the distance 00563 TiledJet * jet; // the jet (i) for which we've found this distance 00564 // (whose NN will the J). 00565 }; 00566 diJ_plus_link * diJ = new diJ_plus_link[n]; 00567 jetA = head; 00568 for (int i = 0; i < n; i++) { 00569 diJ[i].diJ = _bj_diJ(jetA); // kt distance * R^2 00570 diJ[i].jet = jetA; // our compact diJ table will not be in 00571 jetA->diJ_posn = i; // one-to-one corresp. with non-compact jets, 00572 // so set up bi-directional correspondence here. 00573 jetA++; // have jetA follow i 00574 } 00575 00576 // now run the recombination loop 00577 int history_location = n-1; 00578 while (n > 0) { 00579 00580 // find the minimum of the diJ on this round 00581 diJ_plus_link * best, *stop; // pointers a bit faster than indices 00582 // could use best to keep track of diJ min, but it turns out to be 00583 // marginally faster to have a separate variable (avoids n 00584 // dereferences at the expense of n/2 assignments). 00585 double diJ_min = diJ[0].diJ; // initialise the best one here. 00586 best = diJ; // and here 00587 stop = diJ+n; 00588 for (diJ_plus_link * here = diJ+1; here != stop; here++) { 00589 if (here->diJ < diJ_min) {best = here; diJ_min = here->diJ;} 00590 } 00591 00592 // do the recombination between A and B 00593 history_location++; 00594 jetA = best->jet; 00595 jetB = jetA->NN; 00596 // put the normalisation back in 00597 diJ_min *= _invR2; 00598 00599 if (jetB != NULL) { 00600 // jet-jet recombination 00601 // If necessary relabel A & B to ensure jetB < jetA, that way if 00602 // the larger of them == newtail then that ends up being jetA and 00603 // the new jet that is added as jetB is inserted in a position that 00604 // has a future! 00605 if (jetA < jetB) {swap(jetA,jetB);} 00606 00607 int nn; // new jet index 00608 _do_ij_recombination_step(jetA->_jets_index, jetB->_jets_index, diJ_min, nn); 00609 00610 //OBS// get the two history indices 00611 //OBSint ihstry_a = _jets[jetA->_jets_index].cluster_hist_index(); 00612 //OBSint ihstry_b = _jets[jetB->_jets_index].cluster_hist_index(); 00613 //OBS// create the recombined jet 00614 //OBS_jets.push_back(_jets[jetA->_jets_index] + _jets[jetB->_jets_index]); 00615 //OBSint nn = _jets.size() - 1; 00616 //OBS_jets[nn].set_cluster_hist_index(history_location); 00617 //OBS// update history 00618 //OBS//cout <<n-1<<" "<<jetA-head<<" "<<jetB-head<<"; "; 00619 //OBS_add_step_to_history(history_location, 00620 //OBS min(ihstry_a,ihstry_b),max(ihstry_a,ihstry_b), 00621 //OBS nn, diJ_min); 00622 // what was jetB will now become the new jet 00623 _bj_remove_from_tiles(jetA); 00624 oldB = * jetB; // take a copy because we will need it... 00625 _bj_remove_from_tiles(jetB); 00626 _tj_set_jetinfo(jetB, nn); // cause jetB to become _jets[nn] 00627 // (also registers the jet in the tiling) 00628 } else { 00629 // jet-beam recombination 00630 // get the hist_index 00631 _do_iB_recombination_step(jetA->_jets_index, diJ_min); 00632 //OBSint ihstry_a = _jets[jetA->_jets_index].cluster_hist_index(); 00633 //OBS//cout <<n-1<<" "<<jetA-head<<" "<<-1<<"; "; 00634 //OBS_add_step_to_history(history_location,ihstry_a,BeamJet,Invalid,diJ_min); 00635 _bj_remove_from_tiles(jetA); 00636 } 00637 00638 // first establish the set of tiles over which we are going to 00639 // have to run searches for updated and new nearest-neighbours -- 00640 // basically a combination of vicinity of the tiles of the two old 00641 // and one new jet. 00642 int n_near_tiles = 0; 00643 _add_untagged_neighbours_to_tile_union(jetA->tile_index, 00644 tile_union, n_near_tiles); 00645 if (jetB != NULL) { 00646 if (jetB->tile_index != jetA->tile_index) { 00647 _add_untagged_neighbours_to_tile_union(jetB->tile_index, 00648 tile_union,n_near_tiles); 00649 } 00650 if (oldB.tile_index != jetA->tile_index && 00651 oldB.tile_index != jetB->tile_index) { 00652 _add_untagged_neighbours_to_tile_union(oldB.tile_index, 00653 tile_union,n_near_tiles); 00654 } 00655 } 00656 00657 // now update our nearest neighbour info and diJ table 00658 // first reduce size of table 00659 n--; 00660 // then compactify the diJ by taking the last of the diJ and copying 00661 // it to the position occupied by the diJ for jetA 00662 diJ[n].jet->diJ_posn = jetA->diJ_posn; 00663 diJ[jetA->diJ_posn] = diJ[n]; 00664 00665 // Initialise jetB's NN distance as well as updating it for 00666 // other particles. 00667 // Run over all tiles in our union 00668 for (int itile = 0; itile < n_near_tiles; itile++) { 00669 Tile * tile = &_tiles[tile_union[itile]]; 00670 tile->tagged = false; // reset tag, since we're done with unions 00671 // run over all jets in the current tile 00672 for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) { 00673 // see if jetI had jetA or jetB as a NN -- if so recalculate the NN 00674 if (jetI->NN == jetA || (jetI->NN == jetB && jetB != NULL)) { 00675 jetI->NN_dist = _R2; 00676 jetI->NN = NULL; 00677 // now go over tiles that are neighbours of I (include own tile) 00678 for (Tile ** near_tile = tile->begin_tiles; 00679 near_tile != tile->end_tiles; near_tile++) { 00680 // and then over the contents of that tile 00681 for (TiledJet * jetJ = (*near_tile)->head; 00682 jetJ != NULL; jetJ = jetJ->next) { 00683 double dist = _bj_dist(jetI,jetJ); 00684 if (dist < jetI->NN_dist && jetJ != jetI) { 00685 jetI->NN_dist = dist; jetI->NN = jetJ; 00686 } 00687 } 00688 } 00689 diJ[jetI->diJ_posn].diJ = _bj_diJ(jetI); // update diJ kt-dist 00690 } 00691 // check whether new jetB is closer than jetI's current NN and 00692 // if jetI is closer than jetB's current (evolving) nearest 00693 // neighbour. Where relevant update things 00694 if (jetB != NULL) { 00695 double dist = _bj_dist(jetI,jetB); 00696 if (dist < jetI->NN_dist) { 00697 if (jetI != jetB) { 00698 jetI->NN_dist = dist; 00699 jetI->NN = jetB; 00700 diJ[jetI->diJ_posn].diJ = _bj_diJ(jetI); // update diJ... 00701 } 00702 } 00703 if (dist < jetB->NN_dist) { 00704 if (jetI != jetB) { 00705 jetB->NN_dist = dist; 00706 jetB->NN = jetI;} 00707 } 00708 } 00709 } 00710 } 00711 00712 // finally, register the updated kt distance for B 00713 if (jetB != NULL) {diJ[jetB->diJ_posn].diJ = _bj_diJ(jetB);} 00714 00715 } 00716 00717 // final cleaning up; 00718 delete[] diJ; 00719 delete[] briefjets; 00720 }
| void fastjet::ClusterSequence::_fill_initial_history | ( | ) | [protected] |
fill out the history (and jet cross refs) related to the initial set of jets (assumed already to have been "transferred"), without any clustering
Definition at line 248 of file ClusterSequence.cc.
References _history, _initial_n, _jet_def, _jets, _Qtot, fastjet::ClusterSequence::history_element::child, fastjet::ClusterSequence::history_element::dij, InexistentParent, Invalid, fastjet::ClusterSequence::history_element::jetp_index, fastjet::ClusterSequence::history_element::max_dij_so_far, fastjet::ClusterSequence::history_element::parent1, fastjet::ClusterSequence::history_element::parent2, and fastjet::JetDefinition::recombiner().
Referenced by fastjet::ClusterSequenceActiveArea::_initialise_AA(), and _initialise_and_run().
00248 { 00249 00250 //if (_jets.size() == 0) {throw Error("Cannot run jet-finder on empty event");} 00251 00252 // reserve sufficient space for everything 00253 _jets.reserve(_jets.size()*2); 00254 _history.reserve(_jets.size()*2); 00255 00256 _Qtot = 0; 00257 00258 for (int i = 0; i < static_cast<int>(_jets.size()) ; i++) { 00259 history_element element; 00260 element.parent1 = InexistentParent; 00261 element.parent2 = InexistentParent; 00262 element.child = Invalid; 00263 element.jetp_index = i; 00264 element.dij = 0.0; 00265 element.max_dij_so_far = 0.0; 00266 00267 _history.push_back(element); 00268 00269 // do any momentum preprocessing needed by the recombination scheme 00270 _jet_def.recombiner()->preprocess(_jets[i]); 00271 00272 // get cross-referencing right from PseudoJets 00273 _jets[i].set_cluster_hist_index(i); 00274 00275 // determine the total energy in the event 00276 _Qtot += _jets[i].E(); 00277 } 00278 _initial_n = _jets.size(); 00279 }
| void fastjet::ClusterSequence::_initialise_and_run | ( | const double & | R, | |
| const Strategy & | strategy, | |||
| const bool & | writeout_combinations | |||
| ) | [protected] |
This is an alternative routine for initialising and running the clustering, provided for legacy purposes.
The jet finder is that specified in the static member _default_jet_algorithm.
Definition at line 57 of file ClusterSequence.cc.
References _default_jet_algorithm, _initialise_and_run(), and jet_def().
00060 { 00061 00062 JetDefinition jet_def(_default_jet_algorithm, R, strategy); 00063 _initialise_and_run(jet_def, writeout_combinations); 00064 }
| void fastjet::ClusterSequence::_initialise_and_run | ( | const JetDefinition & | jet_def, | |
| const bool & | writeout_combinations | |||
| ) | [protected] |
This is the routine that will do all the initialisation and then run the clustering (may be called by various constructors).
It assumes _jets contains the momenta to be clustered.
Definition at line 68 of file ClusterSequence.cc.
References _CP2DChan_cluster(), _CP2DChan_cluster_2pi2R(), _CP2DChan_cluster_2piMultD(), _decant_options(), _delaunay_cluster(), _faster_tiled_N2_cluster(), _fill_initial_history(), _invR2, _jet_algorithm, _jet_def, _jets, _minheap_faster_tiled_N2_cluster(), _plugin_activated, _R2, _really_dumb_cluster(), _Rparam, _simple_N2_cluster_BriefJet(), _simple_N2_cluster_EEBriefJet(), _strategy, _tiled_N2_cluster(), fastjet::antikt_algorithm, fastjet::Best, fastjet::cambridge_algorithm, fastjet::ee_genkt_algorithm, fastjet::ee_kt_algorithm, fastjet::JetDefinition::jet_algorithm(), fastjet::d0::inline_maths::min(), fastjet::N2MinHeapTiled, fastjet::N2Plain, fastjet::N2PoorTiled, fastjet::N2Tiled, fastjet::N3Dumb, n_particles(), fastjet::NlnN, fastjet::NlnN3pi, fastjet::NlnN4pi, fastjet::NlnNCam, fastjet::NlnNCam2pi2R, fastjet::NlnNCam4pi, fastjet::pi, fastjet::JetDefinition::plugin(), and fastjet::plugin_algorithm.
Referenced by fastjet::ClusterSequenceActiveArea::_initialise_AA(), _initialise_and_run(), and ClusterSequence().
00070 { 00071 00072 // transfer all relevant info into internal variables 00073 _decant_options(jet_def, writeout_combinations); 00074 00075 // set up the history entries for the initial particles (those 00076 // currently in _jets) 00077 _fill_initial_history(); 00078 00079 // don't run anything if the event is empty 00080 if (n_particles() == 0) return; 00081 00082 // ----- deal with special cases: plugins & e+e- ------ 00083 if (_jet_algorithm == plugin_algorithm) { 00084 // allows plugin_xyz() functions to modify cluster sequence 00085 _plugin_activated = true; 00086 // let the plugin do its work here 00087 _jet_def.plugin()->run_clustering( (*this) ); 00088 _plugin_activated = false; 00089 return; 00090 } else if (_jet_algorithm == ee_kt_algorithm || 00091 _jet_algorithm == ee_genkt_algorithm) { 00092 // ignore requested strategy 00093 _strategy = N2Plain; 00094 if (_jet_algorithm == ee_kt_algorithm) { 00095 // make sure that R is large enough so that "beam" recomb only 00096 // occurs when a single particle is left 00097 // Normally, this should be automatically set to 4 from JetDefinition 00098 assert(_Rparam > 2.0); 00099 // this is used to renormalise the dij to get a "standard" form 00100 // and our convention in e+e- will be different from that 00101 // in long.inv case; NB: _invR2 name should be changed -> _renorm_dij? 00102 _invR2 = 1.0; 00103 } else { 00104 // as of 2009-01-09, choose R to be an angular distance, in 00105 // radians. Since the algorithm uses 2(1-cos(theta)) as its 00106 // squared angular measure, make sure that the _R2 is defined 00107 // in a similar way. 00108 if (_Rparam > pi) { 00109 // choose a value that ensures that back-to-back particles will 00110 // always recombine 00111 //_R2 = 4.0000000000001; 00112 _R2 = 2 * ( 3.0 + cos(_Rparam) ); 00113 } else { 00114 _R2 = 2 * ( 1.0 - cos(_Rparam) ); 00115 } 00116 _invR2 = 1.0/_R2; 00117 } 00118 //_simple_N2_cluster<EEBriefJet>(); 00119 _simple_N2_cluster_EEBriefJet(); 00120 return; 00121 } 00122 00123 00124 // automatically redefine the strategy according to N if that is 00125 // what the user requested -- transition points (and especially 00126 // their R-dependence) are based on empirical observations for a 00127 // R=0.4, 0.7 and 1.0, running on toth (3.4GHz, Pentium IV D [dual 00128 // core] with 2MB of cache). 00129 if (_strategy == Best) { 00130 int N = _jets.size(); 00131 if (N > 6200/pow(_Rparam,2.0) 00132 && jet_def.jet_algorithm() == cambridge_algorithm) { 00133 _strategy = NlnNCam;} 00134 else 00135 #ifndef DROP_CGAL 00136 if ((N > 16000/pow(_Rparam,1.15) && jet_def.jet_algorithm() != antikt_algorithm) 00137 || N > 35000/pow(_Rparam,1.15)) { 00138 _strategy = NlnN; } 00139 else 00140 #endif // DROP_CGAL 00141 if (N > 450) { 00142 _strategy = N2MinHeapTiled; 00143 } 00144 else if (N > 55*max(0.5,min(1.0,_Rparam))) {// empirical scaling with R 00145 _strategy = N2Tiled; 00146 } else { 00147 _strategy = N2Plain; 00148 } 00149 } 00150 00151 00152 // run the code containing the selected strategy 00153 if (_strategy == NlnN || _strategy == NlnN3pi 00154 || _strategy == NlnN4pi ) { 00155 this->_delaunay_cluster(); 00156 } else if (_strategy == N3Dumb ) { 00157 this->_really_dumb_cluster(); 00158 } else if (_strategy == N2Tiled) { 00159 this->_faster_tiled_N2_cluster(); 00160 } else if (_strategy == N2PoorTiled) { 00161 this->_tiled_N2_cluster(); 00162 } else if (_strategy == N2Plain) { 00163 // BriefJet provides standard long.invariant kt alg. 00164 //this->_simple_N2_cluster<BriefJet>(); 00165 this->_simple_N2_cluster_BriefJet(); 00166 } else if (_strategy == N2MinHeapTiled) { 00167 this->_minheap_faster_tiled_N2_cluster(); 00168 } else if (_strategy == NlnNCam4pi) { 00169 this->_CP2DChan_cluster(); 00170 } else if (_strategy == NlnNCam2pi2R) { 00171 this->_CP2DChan_cluster_2pi2R(); 00172 } else if (_strategy == NlnNCam) { 00173 this->_CP2DChan_cluster_2piMultD(); 00174 } else { 00175 ostringstream err; 00176 err << "Unrecognised value for strategy: "<<_strategy; 00177 throw Error(err.str()); 00178 //assert(false); 00179 } 00180 }
| void fastjet::ClusterSequence::_initialise_tiles | ( | ) | [private] |
Set up the tiles:
- decide the range in eta
- allocate the tiles
- set up the cross-referencing info between tiles.
The neighbourhood of a tile is set up as follows
LRR LXR LLR
such that tiles is an array containing XLLLLRRRR with pointers | \ RH_tiles \ surrounding_tiles
with appropriate precautions when close to the edge of the tiled region.
Definition at line 86 of file ClusterSequence_TiledN2.cc.
References _jets, _n_tiles_phi, _Rparam, _tile_index(), _tile_size_eta, _tile_size_phi, _tiles, _tiles_eta_max, _tiles_eta_min, _tiles_ieta_max, _tiles_ieta_min, fastjet::ClusterSequence::Tile::begin_tiles, fastjet::ClusterSequence::Tile::end_tiles, fastjet::ClusterSequence::Tile::head, fastjet::ClusterSequence::Tile::RH_tiles, fastjet::ClusterSequence::Tile::surrounding_tiles, fastjet::ClusterSequence::Tile::tagged, and fastjet::twopi.
Referenced by _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().
00086 { 00087 00088 // first decide tile sizes (with a lower bound to avoid huge memory use with 00089 // very small R) 00090 double default_size = max(0.1,_Rparam); 00091 _tile_size_eta = default_size; 00092 _n_tiles_phi = int(floor(twopi/default_size)); 00093 _tile_size_phi = twopi / _n_tiles_phi; // >= _Rparam and fits in 2pi 00094 00095 // always include zero rapidity in the tiling region 00096 _tiles_eta_min = 0.0; 00097 _tiles_eta_max = 0.0; 00098 // but go no further than following 00099 const double maxrap = 7.0; 00100 00101 // and find out how much further one should go 00102 for(unsigned int i = 0; i < _jets.size(); i++) { 00103 double eta = _jets[i].rap(); 00104 // first check if eta is in range -- to avoid taking into account 00105 // very spurious rapidities due to particles with near-zero kt. 00106 if (abs(eta) < maxrap) { 00107 if (eta < _tiles_eta_min) {_tiles_eta_min = eta;} 00108 if (eta > _tiles_eta_max) {_tiles_eta_max = eta;} 00109 } 00110 } 00111 00112 // now adjust the values 00113 _tiles_ieta_min = int(floor(_tiles_eta_min/_tile_size_eta)); 00114 _tiles_ieta_max = int(floor( _tiles_eta_max/_tile_size_eta)); 00115 _tiles_eta_min = _tiles_ieta_min * _tile_size_eta; 00116 _tiles_eta_max = _tiles_ieta_max * _tile_size_eta; 00117 00118 // allocate the tiles 00119 _tiles.resize((_tiles_ieta_max-_tiles_ieta_min+1)*_n_tiles_phi); 00120 00121 // now set up the cross-referencing between tiles 00122 for (int ieta = _tiles_ieta_min; ieta <= _tiles_ieta_max; ieta++) { 00123 for (int iphi = 0; iphi < _n_tiles_phi; iphi++) { 00124 Tile * tile = & _tiles[_tile_index(ieta,iphi)]; 00125 // no jets in this tile yet 00126 tile->head = NULL; // first element of tiles points to itself 00127 tile->begin_tiles[0] = tile; 00128 Tile ** pptile = & (tile->begin_tiles[0]); 00129 pptile++; 00130 // 00131 // set up L's in column to the left of X 00132 tile->surrounding_tiles = pptile; 00133 if (ieta > _tiles_ieta_min) { 00134 // with the itile subroutine, we can safely run tiles from 00135 // idphi=-1 to idphi=+1, because it takes care of 00136 // negative and positive boundaries 00137 for (int idphi = -1; idphi <=+1; idphi++) { 00138 *pptile = & _tiles[_tile_index(ieta-1,iphi+idphi)]; 00139 pptile++; 00140 } 00141 } 00142 // now set up last L (below X) 00143 *pptile = & _tiles[_tile_index(ieta,iphi-1)]; 00144 pptile++; 00145 // set up first R (above X) 00146 tile->RH_tiles = pptile; 00147 *pptile = & _tiles[_tile_index(ieta,iphi+1)]; 00148 pptile++; 00149 // set up remaining R's, to the right of X 00150 if (ieta < _tiles_ieta_max) { 00151 for (int idphi = -1; idphi <= +1; idphi++) { 00152 *pptile = & _tiles[_tile_index(ieta+1,iphi+idphi)]; 00153 pptile++; 00154 } 00155 } 00156 // now put semaphore for end tile 00157 tile->end_tiles = pptile; 00158 // finally make sure tiles are untagged 00159 tile->tagged = false; 00160 } 00161 } 00162 00163 }
| void fastjet::ClusterSequence::_minheap_faster_tiled_N2_cluster | ( | ) | [private] |
run a tiled clustering, with our minheap for keeping track of the smallest dij
Definition at line 727 of file ClusterSequence_TiledN2.cc.
References _add_untagged_neighbours_to_tile_union(), _bj_diJ(), _bj_dist(), _bj_remove_from_tiles(), _do_iB_recombination_step(), _do_ij_recombination_step(), _initialise_tiles(), _invR2, _jets, fastjet::ClusterSequence::TiledJet::_jets_index, _R2, _tiles, _tj_set_jetinfo(), fastjet::ClusterSequence::Tile::begin_tiles, fastjet::ClusterSequence::Tile::end_tiles, fastjet::ClusterSequence::Tile::head, fastjet::ClusterSequence::TiledJet::label_minheap_update_done(), fastjet::MinHeap::minloc(), fastjet::MinHeap::minval(), n_tile_neighbours, fastjet::ClusterSequence::TiledJet::next, fastjet::ClusterSequence::TiledJet::NN, fastjet::ClusterSequence::TiledJet::NN_dist, fastjet::MinHeap::remove(), fastjet::ClusterSequence::Tile::tagged, fastjet::ClusterSequence::TiledJet::tile_index, and fastjet::MinHeap::update().
Referenced by _initialise_and_run().
00727 { 00728 00729 _initialise_tiles(); 00730 00731 int n = _jets.size(); 00732 TiledJet * briefjets = new TiledJet[n]; 00733 TiledJet * jetA = briefjets, * jetB; 00734 TiledJet oldB; 00735 00736 00737 // will be used quite deep inside loops, but declare it here so that 00738 // memory (de)allocation gets done only once 00739 vector<int> tile_union(3*n_tile_neighbours); 00740 00741 // initialise the basic jet info 00742 for (int i = 0; i< n; i++) { 00743 _tj_set_jetinfo(jetA, i); 00744 //cout << i<<": "<<jetA->tile_index<<"\n"; 00745 jetA++; // move on to next entry of briefjets 00746 } 00747 TiledJet * head = briefjets; // a nicer way of naming start 00748 00749 // set up the initial nearest neighbour information 00750 vector<Tile>::const_iterator tile; 00751 for (tile = _tiles.begin(); tile != _tiles.end(); tile++) { 00752 // first do it on this tile 00753 for (jetA = tile->head; jetA != NULL; jetA = jetA->next) { 00754 for (jetB = tile->head; jetB != jetA; jetB = jetB->next) { 00755 double dist = _bj_dist(jetA,jetB); 00756 if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;} 00757 if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;} 00758 } 00759 } 00760 // then do it for RH tiles 00761 for (Tile ** RTile = tile->RH_tiles; RTile != tile->end_tiles; RTile++) { 00762 for (jetA = tile->head; jetA != NULL; jetA = jetA->next) { 00763 for (jetB = (*RTile)->head; jetB != NULL; jetB = jetB->next) { 00764 double dist = _bj_dist(jetA,jetB); 00765 if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;} 00766 if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;} 00767 } 00768 } 00769 } 00770 // no need to do it for LH tiles, since they are implicitly done 00771 // when we set NN for both jetA and jetB on the RH tiles. 00772 } 00773 00774 00778 //struct diJ_plus_link { 00779 // double diJ; // the distance 00780 // TiledJet * jet; // the jet (i) for which we've found this distance 00781 // // (whose NN will the J). 00782 //}; 00783 //diJ_plus_link * diJ = new diJ_plus_link[n]; 00784 //jetA = head; 00785 //for (int i = 0; i < n; i++) { 00786 // diJ[i].diJ = _bj_diJ(jetA); // kt distance * R^2 00787 // diJ[i].jet = jetA; // our compact diJ table will not be in 00788 // jetA->diJ_posn = i; // one-to-one corresp. with non-compact jets, 00789 // // so set up bi-directional correspondence here. 00790 // jetA++; // have jetA follow i 00791 //} 00792 00793 vector<double> diJs(n); 00794 for (int i = 0; i < n; i++) { 00795 diJs[i] = _bj_diJ(&briefjets[i]); 00796 briefjets[i].label_minheap_update_done(); 00797 } 00798 MinHeap minheap(diJs); 00799 // have a stack telling us which jets we'll have to update on the heap 00800 vector<TiledJet *> jets_for_minheap; 00801 jets_for_minheap.reserve(n); 00802 00803 // now run the recombination loop 00804 int history_location = n-1; 00805 while (n > 0) { 00806 00807 double diJ_min = minheap.minval() *_invR2; 00808 jetA = head + minheap.minloc(); 00809 00810 // do the recombination between A and B 00811 history_location++; 00812 jetB = jetA->NN; 00813 00814 if (jetB != NULL) { 00815 // jet-jet recombination 00816 // If necessary relabel A & B to ensure jetB < jetA, that way if 00817 // the larger of them == newtail then that ends up being jetA and 00818 // the new jet that is added as jetB is inserted in a position that 00819 // has a future! 00820 if (jetA < jetB) {swap(jetA,jetB);} 00821 00822 int nn; // new jet index 00823 _do_ij_recombination_step(jetA->_jets_index, jetB->_jets_index, diJ_min, nn); 00824 00825 // what was jetB will now become the new jet 00826 _bj_remove_from_tiles(jetA); 00827 oldB = * jetB; // take a copy because we will need it... 00828 _bj_remove_from_tiles(jetB); 00829 _tj_set_jetinfo(jetB, nn); // cause jetB to become _jets[nn] 00830 // (also registers the jet in the tiling) 00831 } else { 00832 // jet-beam recombination 00833 // get the hist_index 00834 _do_iB_recombination_step(jetA->_jets_index, diJ_min); 00835 _bj_remove_from_tiles(jetA); 00836 } 00837 00838 // remove the minheap entry for jetA 00839 minheap.remove(jetA-head); 00840 00841 // first establish the set of tiles over which we are going to 00842 // have to run searches for updated and new nearest-neighbours -- 00843 // basically a combination of vicinity of the tiles of the two old 00844 // and one new jet. 00845 int n_near_tiles = 0; 00846 _add_untagged_neighbours_to_tile_union(jetA->tile_index, 00847 tile_union, n_near_tiles); 00848 if (jetB != NULL) { 00849 if (jetB->tile_index != jetA->tile_index) { 00850 _add_untagged_neighbours_to_tile_union(jetB->tile_index, 00851 tile_union,n_near_tiles); 00852 } 00853 if (oldB.tile_index != jetA->tile_index && 00854 oldB.tile_index != jetB->tile_index) { 00855 _add_untagged_neighbours_to_tile_union(oldB.tile_index, 00856 tile_union,n_near_tiles); 00857 } 00858 // indicate that we'll have to update jetB in the minheap 00859 jetB->label_minheap_update_needed(); 00860 jets_for_minheap.push_back(jetB); 00861 } 00862 00863 00864 // Initialise jetB's NN distance as well as updating it for 00865 // other particles. 00866 // Run over all tiles in our union 00867 for (int itile = 0; itile < n_near_tiles; itile++) { 00868 Tile * tile = &_tiles[tile_union[itile]]; 00869 tile->tagged = false; // reset tag, since we're done with unions 00870 // run over all jets in the current tile 00871 for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) { 00872 // see if jetI had jetA or jetB as a NN -- if so recalculate the NN 00873 if (jetI->NN == jetA || (jetI->NN == jetB && jetB != NULL)) { 00874 jetI->NN_dist = _R2; 00875 jetI->NN = NULL; 00876 // label jetI as needing heap action... 00877 if (!jetI->minheap_update_needed()) { 00878 jetI->label_minheap_update_needed(); 00879 jets_for_minheap.push_back(jetI);} 00880 // now go over tiles that are neighbours of I (include own tile) 00881 for (Tile ** near_tile = tile->begin_tiles; 00882 near_tile != tile->end_tiles; near_tile++) { 00883 // and then over the contents of that tile 00884 for (TiledJet * jetJ = (*near_tile)->head; 00885 jetJ != NULL; jetJ = jetJ->next) { 00886 double dist = _bj_dist(jetI,jetJ); 00887 if (dist < jetI->NN_dist && jetJ != jetI) { 00888 jetI->NN_dist = dist; jetI->NN = jetJ; 00889 } 00890 } 00891 } 00892 } 00893 // check whether new jetB is closer than jetI's current NN and 00894 // if jetI is closer than jetB's current (evolving) nearest 00895 // neighbour. Where relevant update things 00896 if (jetB != NULL) { 00897 double dist = _bj_dist(jetI,jetB); 00898 if (dist < jetI->NN_dist) { 00899 if (jetI != jetB) { 00900 jetI->NN_dist = dist; 00901 jetI->NN = jetB; 00902 // label jetI as needing heap action... 00903 if (!jetI->minheap_update_needed()) { 00904 jetI->label_minheap_update_needed(); 00905 jets_for_minheap.push_back(jetI);} 00906 } 00907 } 00908 if (dist < jetB->NN_dist) { 00909 if (jetI != jetB) { 00910 jetB->NN_dist = dist; 00911 jetB->NN = jetI;} 00912 } 00913 } 00914 } 00915 } 00916 00917 // deal with jets whose minheap entry needs updating 00918 while (jets_for_minheap.size() > 0) { 00919 TiledJet * jetI = jets_for_minheap.back(); 00920 jets_for_minheap.pop_back(); 00921 minheap.update(jetI-head, _bj_diJ(jetI)); 00922 jetI->label_minheap_update_done(); 00923 } 00924 n--; 00925 } 00926 00927 // final cleaning up; 00928 delete[] briefjets; 00929 }
| bool fastjet::ClusterSequence::_potentially_valid | ( | const PseudoJet & | jet | ) | const [inline, protected] |
returns true if the jet has a history index contained within the range of this CS
Definition at line 487 of file ClusterSequence.hh.
References fastjet::PseudoJet::cluster_hist_index().
Referenced by object_in_jet().
00487 { 00488 return jet.cluster_hist_index() >= 0 00489 && jet.cluster_hist_index() < int(_history.size()); 00490 }
| void fastjet::ClusterSequence::_print_banner | ( | ) | [private] |
for making sure the user knows what it is they're running...
Definition at line 197 of file ClusterSequence.cc.
References _first_time, and fastjet_version.
Referenced by _decant_options().
00197 { 00198 00199 if (!_first_time) {return;} 00200 _first_time = false; 00201 00202 00203 //Symp. Discr. Alg, p.472 (2002) and CGAL (http://www.cgal.org); 00204 00205 cout << "#--------------------------------------------------------------------------\n"; 00206 cout << "# FastJet release " << fastjet_version << endl; 00207 cout << "# Written by M. Cacciari, G.P. Salam and G. Soyez \n"; 00208 cout << "# http://www.fastjet.fr \n"; 00209 cout << "# \n"; 00210 cout << "# Longitudinally invariant Kt, anti-Kt, and inclusive Cambridge/Aachen \n"; 00211 cout << "# clustering using fast geometric algorithms, with area measures and optional\n"; 00212 cout << "# external jet-finder plugins. \n"; 00213 cout << "# Please cite Phys. Lett. B641 (2006) [hep-ph/0512210] if you use this code.\n"; 00214 cout << "# \n"; 00215 cout << "# This package uses T.Chan's closest pair algorithm, Proc.13th ACM-SIAM \n"; 00216 cout << "# Symp. Discr. Alg, p.472 (2002), S.Fortune's Voronoi algorithm and code " ; 00217 #ifndef DROP_CGAL 00218 cout << endl << "# and CGAL: http://www.cgal.org/"; 00219 #endif // DROP_CGAL 00220 cout << ".\n"; 00221 cout << "#-------------------------------------------------------------------------\n"; 00222 }
| void fastjet::ClusterSequence::_print_tiles | ( | TiledJet * | briefjets | ) | const [private] |
output the contents of the tiles
Definition at line 212 of file ClusterSequence_TiledN2.cc.
References _tiles.
00212 { 00213 for (vector<Tile>::const_iterator tile = _tiles.begin(); 00214 tile < _tiles.end(); tile++) { 00215 cout << "Tile " << tile - _tiles.begin()<<" = "; 00216 vector<int> list; 00217 for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) { 00218 list.push_back(jetI-briefjets); 00219 //cout <<" "<<jetI-briefjets; 00220 } 00221 sort(list.begin(),list.end()); 00222 for (unsigned int i = 0; i < list.size(); i++) {cout <<" "<<list[i];} 00223 cout <<"\n"; 00224 } 00225 }
| void fastjet::ClusterSequence::_really_dumb_cluster | ( | ) | [private] |
Run the clustering in a very slow variant of the N^3 algorithm.
The only thing this routine has going for it is that memory usage is O(N)!
Definition at line 50 of file ClusterSequence_DumbN3.cc.
References _do_iB_recombination_step(), _do_ij_recombination_step(), _invR2, _jets, BeamJet, jet_scale_for_algorithm(), fastjet::d0::inline_maths::min(), and fastjet::d0::inline_maths::y().
Referenced by _initialise_and_run().
00050 { 00051 00052 // the array that will be overwritten here will be one 00053 // of pointers to jets. 00054 vector<PseudoJet *> jetsp(_jets.size()); 00055 vector<int> indices(_jets.size()); 00056 00057 for (size_t i = 0; i<_jets.size(); i++) { 00058 jetsp[i] = & _jets[i]; 00059 indices[i] = i; 00060 } 00061 00062 for (int n = jetsp.size(); n > 0; n--) { 00063 int ii, jj; 00064 // find smallest beam distance [remember jet_scale_for_algorithm 00065 // will return kt^2 for the kt-algorithm and 1 for the Cambridge/Aachen] 00066 double ymin = jet_scale_for_algorithm(*(jetsp[0])); 00067 ii = 0; jj = -2; 00068 for (int i = 0; i < n; i++) { 00069 double yiB = jet_scale_for_algorithm(*(jetsp[i])); 00070 if (yiB < ymin) { 00071 ymin = yiB; ii = i; jj = -2;} 00072 } 00073 00074 // find smallest distance between pair of jetsp 00075 for (int i = 0; i < n-1; i++) { 00076 for (int j = i+1; j < n; j++) { 00077 //double y = jetsp[i]->kt_distance(*jetsp[j])*_invR2; 00078 double y = min(jet_scale_for_algorithm(*(jetsp[i])), 00079 jet_scale_for_algorithm(*(jetsp[j]))) 00080 * jetsp[i]->plain_distance(*jetsp[j])*_invR2; 00081 if (y < ymin) {ymin = y; ii = i; jj = j;} 00082 } 00083 } 00084 00085 // output recombination sequence 00086 // old "ktclus" way of labelling 00087 //cout <<n<< " "<< ii+1 << " with " << jj+1 << "; y = "<< ymin<<endl; 00088 // new delaunay way of labelling 00089 int jjindex_or_beam, iiindex; 00090 if (jj < 0) {jjindex_or_beam = BeamJet; iiindex = indices[ii];} 00091 else { 00092 jjindex_or_beam = max(indices[ii],indices[jj]); 00093 iiindex = min(indices[ii],indices[jj]); 00094 } 00095 00096 // now recombine 00097 int newn = 2*jetsp.size() - n; 00098 if (jj >= 0) { 00099 // combine pair 00100 int nn; // new jet index 00101 _do_ij_recombination_step(jetsp[ii]-&_jets[0], 00102 jetsp[jj]-&_jets[0], ymin, nn); 00103 00104 // sort out internal bookkeeping 00105 jetsp[ii] = &_jets[nn]; 00106 // have jj point to jet that was pointed at by n-1 00107 // (since original jj is no longer current, so put n-1 into jj) 00108 jetsp[jj] = jetsp[n-1]; 00109 indices[ii] = newn; 00110 indices[jj] = indices[n-1]; 00111 00112 //OBS_jets.push_back(*jetsp[ii] + *jetsp[jj]); 00113 //OBSjetsp[ii] = &_jets[_jets.size()-1]; 00114 //OBS// have jj point to jet that was pointed at by n-1 00115 //OBS// (since original jj is no longer current, so put n-1 into jj) 00116 //OBSjetsp[jj] = jetsp[n-1]; 00117 //OBS 00118 //OBSindices[ii] = newn; 00119 //OBSindices[jj] = indices[n-1]; 00120 //OBS_add_step_to_history(newn,iiindex, 00121 //OBS jjindex_or_beam,_jets.size()-1,ymin); 00122 } else { 00123 // combine ii with beam 00124 _do_iB_recombination_step(jetsp[ii]-&_jets[0], ymin); 00125 // put last jet (pointer) in place of ii (which has disappeared) 00126 jetsp[ii] = jetsp[n-1]; 00127 indices[ii] = indices[n-1]; 00128 //OBS_add_step_to_history(newn,iiindex,jjindex_or_beam,Invalid, ymin); 00129 } 00130 } 00131 00132 }
| void fastjet::ClusterSequence::_simple_N2_cluster | ( | ) | [inline, private] |
| void fastjet::ClusterSequence::_simple_N2_cluster_BriefJet | ( | ) | [private] |
to help instantiation (fj 2.4.0; did not quite work on gcc 33 and os x 10.3?)
to avoid issues with template instantiation (OS X 10.3, gcc 3.3)
Definition at line 115 of file ClusterSequence_N2.cc.
Referenced by _initialise_and_run().
| void fastjet::ClusterSequence::_simple_N2_cluster_EEBriefJet | ( | ) | [private] |
to avoid issues with template instantiation (OS X 10.3, gcc 3.3)
Definition at line 121 of file ClusterSequence_N2.cc.
Referenced by _initialise_and_run().
| int fastjet::ClusterSequence::_tile_index | ( | const double & | eta, | |
| const double & | phi | |||
| ) | const [private] |
return the tile index corresponding to the given eta,phi point
Definition at line 168 of file ClusterSequence_TiledN2.cc.
References _n_tiles_phi, _tile_size_eta, _tile_size_phi, _tiles_eta_max, _tiles_eta_min, _tiles_ieta_max, _tiles_ieta_min, and fastjet::twopi.
00168 { 00169 int ieta, iphi; 00170 if (eta <= _tiles_eta_min) {ieta = 0;} 00171 else if (eta >= _tiles_eta_max) {ieta = _tiles_ieta_max-_tiles_ieta_min;} 00172 else { 00173 //ieta = int(floor((eta - _tiles_eta_min) / _tile_size_eta)); 00174 ieta = int(((eta - _tiles_eta_min) / _tile_size_eta)); 00175 // following needed in case of rare but nasty rounding errors 00176 if (ieta > _tiles_ieta_max-_tiles_ieta_min) { 00177 ieta = _tiles_ieta_max-_tiles_ieta_min;} 00178 } 00179 // allow for some extent of being beyond range in calculation of phi 00180 // as well 00181 //iphi = (int(floor(phi/_tile_size_phi)) + _n_tiles_phi) % _n_tiles_phi; 00182 // with just int and no floor, things run faster but beware 00183 iphi = int((phi+twopi)/_tile_size_phi) % _n_tiles_phi; 00184 return (iphi + ieta * _n_tiles_phi); 00185 }
| int fastjet::ClusterSequence::_tile_index | ( | int | ieta, | |
| int | iphi | |||
| ) | const [inline, private] |
Definition at line 728 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tj_set_jetinfo().
00728 { 00729 // note that (-1)%n = -1 so that we have to add _n_tiles_phi 00730 // before performing modulo operation 00731 return (ieta-_tiles_ieta_min)*_n_tiles_phi 00732 + (iphi+_n_tiles_phi) % _n_tiles_phi; 00733 }
| void fastjet::ClusterSequence::_tiled_N2_cluster | ( | ) | [private] |
run a tiled clustering
Definition at line 267 of file ClusterSequence_TiledN2.cc.
References _add_neighbours_to_tile_union(), _bj_diJ(), _bj_dist(), _bj_remove_from_tiles(), _do_iB_recombination_step(), _do_ij_recombination_step(), _initialise_tiles(), _invR2, _jets, fastjet::ClusterSequence::TiledJet::_jets_index, _R2, _tiles, _tj_set_jetinfo(), fastjet::ClusterSequence::Tile::begin_tiles, fastjet::ClusterSequence::Tile::end_tiles, fastjet::ClusterSequence::Tile::head, n_tile_neighbours, fastjet::ClusterSequence::TiledJet::next, fastjet::ClusterSequence::TiledJet::NN, fastjet::ClusterSequence::TiledJet::NN_dist, fastjet::ClusterSequence::TiledJet::previous, and fastjet::ClusterSequence::TiledJet::tile_index.
Referenced by _initialise_and_run().
00267 { 00268 00269 _initialise_tiles(); 00270 00271 int n = _jets.size(); 00272 TiledJet * briefjets = new TiledJet[n]; 00273 TiledJet * jetA = briefjets, * jetB; 00274 TiledJet oldB; 00275 00276 00277 // will be used quite deep inside loops, but declare it here so that 00278 // memory (de)allocation gets done only once 00279 vector<int> tile_union(3*n_tile_neighbours); 00280 00281 // initialise the basic jet info 00282 for (int i = 0; i< n; i++) { 00283 _tj_set_jetinfo(jetA, i); 00284 //cout << i<<": "<<jetA->tile_index<<"\n"; 00285 jetA++; // move on to next entry of briefjets 00286 } 00287 TiledJet * tail = jetA; // a semaphore for the end of briefjets 00288 TiledJet * head = briefjets; // a nicer way of naming start 00289 00290 // set up the initial nearest neighbour information 00291 vector<Tile>::const_iterator tile; 00292 for (tile = _tiles.begin(); tile != _tiles.end(); tile++) { 00293 // first do it on this tile 00294 for (jetA = tile->head; jetA != NULL; jetA = jetA->next) { 00295 for (jetB = tile->head; jetB != jetA; jetB = jetB->next) { 00296 double dist = _bj_dist(jetA,jetB); 00297 if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;} 00298 if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;} 00299 } 00300 } 00301 // then do it for RH tiles 00302 for (Tile ** RTile = tile->RH_tiles; RTile != tile->end_tiles; RTile++) { 00303 for (jetA = tile->head; jetA != NULL; jetA = jetA->next) { 00304 for (jetB = (*RTile)->head; jetB != NULL; jetB = jetB->next) { 00305 double dist = _bj_dist(jetA,jetB); 00306 if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;} 00307 if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;} 00308 } 00309 } 00310 } 00311 } 00312 00313 // now create the diJ (where J is i's NN) table -- remember that 00314 // we differ from standard normalisation here by a factor of R2 00315 double * diJ = new double[n]; 00316 jetA = head; 00317 for (int i = 0; i < n; i++) { 00318 diJ[i] = _bj_diJ(jetA); 00319 jetA++; // have jetA follow i 00320 } 00321 00322 // now run the recombination loop 00323 int history_location = n-1; 00324 while (tail != head) { 00325 00326 // find the minimum of the diJ on this round 00327 double diJ_min = diJ[0]; 00328 int diJ_min_jet = 0; 00329 for (int i = 1; i < n; i++) { 00330 if (diJ[i] < diJ_min) {diJ_min_jet = i; diJ_min = diJ[i];} 00331 } 00332 00333 // do the recombination between A and B 00334 history_location++; 00335 jetA = & briefjets[diJ_min_jet]; 00336 jetB = jetA->NN; 00337 // put the normalisation back in 00338 diJ_min *= _invR2; 00339 00340 //if (n == 19) {cout << "Hello "<<jetA-head<<" "<<jetB-head<<"\n";} 00341 00342 //cout <<" WILL RECOMBINE "<< jetA-briefjets<<" "<<jetB-briefjets<<"\n"; 00343 00344 if (jetB != NULL) { 00345 // jet-jet recombination 00346 // If necessary relabel A & B to ensure jetB < jetA, that way if 00347 // the larger of them == newtail then that ends up being jetA and 00348 // the new jet that is added as jetB is inserted in a position that 00349 // has a future! 00350 if (jetA < jetB) {swap(jetA,jetB);} 00351 00352 int nn; // new jet index 00353 _do_ij_recombination_step(jetA->_jets_index, jetB->_jets_index, diJ_min, nn); 00354 00355 //OBS// get the two history indices 00356 //OBSint hist_a = _jets[jetA->_jets_index].cluster_hist_index(); 00357 //OBSint hist_b = _jets[jetB->_jets_index].cluster_hist_index(); 00358 //OBS// create the recombined jet 00359 //OBS_jets.push_back(_jets[jetA->_jets_index] + _jets[jetB->_jets_index]); 00360 //OBSint nn = _jets.size() - 1; 00361 //OBS_jets[nn].set_cluster_hist_index(history_location); 00362 //OBS// update history 00363 //OBS//cout <<n-1<<" "<<jetA-head<<" "<<jetB-head<<"; "; 00364 //OBS_add_step_to_history(history_location, 00365 //OBS min(hist_a,hist_b),max(hist_a,hist_b), 00366 //OBS nn, diJ_min); 00367 00368 // what was jetB will now become the new jet 00369 _bj_remove_from_tiles(jetA); 00370 oldB = * jetB; // take a copy because we will need it... 00371 _bj_remove_from_tiles(jetB); 00372 _tj_set_jetinfo(jetB, nn); // also registers the jet in the tiling 00373 } else { 00374 // jet-beam recombination 00375 _do_iB_recombination_step(jetA->_jets_index, diJ_min); 00376 00377 //OBS// get the hist_index 00378 //OBSint hist_a = _jets[jetA->_jets_index].cluster_hist_index(); 00379 //OBS//cout <<n-1<<" "<<jetA-head<<" "<<-1<<"; "; 00380 //OBS_add_step_to_history(history_location,hist_a,BeamJet,Invalid,diJ_min); 00381 _bj_remove_from_tiles(jetA); 00382 } 00383 00384 // first establish the set of tiles over which we are going to 00385 // have to run searches for updated and new nearest-neighbours 00386 int n_near_tiles = 0; 00387 _add_neighbours_to_tile_union(jetA->tile_index, tile_union, n_near_tiles); 00388 if (jetB != NULL) { 00389 bool sort_it = false; 00390 if (jetB->tile_index != jetA->tile_index) { 00391 sort_it = true; 00392 _add_neighbours_to_tile_union(jetB->tile_index,tile_union,n_near_tiles); 00393 } 00394 if (oldB.tile_index != jetA->tile_index && 00395 oldB.tile_index != jetB->tile_index) { 00396 sort_it = true; 00397 _add_neighbours_to_tile_union(oldB.tile_index,tile_union,n_near_tiles); 00398 } 00399 00400 if (sort_it) { 00401 // sort the tiles before then compressing the list 00402 sort(tile_union.begin(), tile_union.begin()+n_near_tiles); 00403 // and now condense the list 00404 int nnn = 1; 00405 for (int i = 1; i < n_near_tiles; i++) { 00406 if (tile_union[i] != tile_union[nnn-1]) { 00407 tile_union[nnn] = tile_union[i]; 00408 nnn++; 00409 } 00410 } 00411 n_near_tiles = nnn; 00412 } 00413 } 00414 00415 // now update our nearest neighbour info and diJ table 00416 // first reduce size of table 00417 tail--; n--; 00418 if (jetA == tail) { 00419 // there is nothing to be done 00420 } else { 00421 // Copy last jet contents and diJ info into position of jetA 00422 *jetA = *tail; 00423 diJ[jetA - head] = diJ[tail-head]; 00424 // IN the tiling fix pointers to tail and turn them into 00425 // pointers to jetA (from predecessors, successors and the tile 00426 // head if need be) 00427 if (jetA->previous == NULL) { 00428 _tiles[jetA->tile_index].head = jetA; 00429 } else { 00430 jetA->previous->next = jetA; 00431 } 00432 if (jetA->next != NULL) {jetA->next->previous = jetA;} 00433 } 00434 00435 // Initialise jetB's NN distance as well as updating it for 00436 // other particles. 00437 for (int itile = 0; itile < n_near_tiles; itile++) { 00438 Tile * tile = &_tiles[tile_union[itile]]; 00439 for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) { 00440 // see if jetI had jetA or jetB as a NN -- if so recalculate the NN 00441 if (jetI->NN == jetA || (jetI->NN == jetB && jetB != NULL)) { 00442 jetI->NN_dist = _R2; 00443 jetI->NN = NULL; 00444 // now go over tiles that are neighbours of I (include own tile) 00445 for (Tile ** near_tile = tile->begin_tiles; 00446 near_tile != tile->end_tiles; near_tile++) { 00447 // and then over the contents of that tile 00448 for (TiledJet * jetJ = (*near_tile)->head; 00449 jetJ != NULL; jetJ = jetJ->next) { 00450 double dist = _bj_dist(jetI,jetJ); 00451 if (dist < jetI->NN_dist && jetJ != jetI) { 00452 jetI->NN_dist = dist; jetI->NN = jetJ; 00453 } 00454 } 00455 } 00456 diJ[jetI-head] = _bj_diJ(jetI); // update diJ 00457 } 00458 // check whether new jetB is closer than jetI's current NN and 00459 // if need to update things 00460 if (jetB != NULL) { 00461 double dist = _bj_dist(jetI,jetB); 00462 if (dist < jetI->NN_dist) { 00463 if (jetI != jetB) { 00464 jetI->NN_dist = dist; 00465 jetI->NN = jetB; 00466 diJ[jetI-head] = _bj_diJ(jetI); // update diJ... 00467 } 00468 } 00469 if (dist < jetB->NN_dist) { 00470 if (jetI != jetB) { 00471 jetB->NN_dist = dist; 00472 jetB->NN = jetI;} 00473 } 00474 } 00475 } 00476 } 00477 00478 00479 if (jetB != NULL) {diJ[jetB-head] = _bj_diJ(jetB);} 00480 //cout << n<<" "<<briefjets[95].NN-briefjets<<" "<<briefjets[95].NN_dist <<"\n"; 00481 00482 // remember to update pointers to tail 00483 for (Tile ** near_tile = _tiles[tail->tile_index].begin_tiles; 00484 near_tile!= _tiles[tail->tile_index].end_tiles; near_tile++){ 00485 // and then the contents of that tile 00486 for (TiledJet * jetJ = (*near_tile)->head; 00487 jetJ != NULL; jetJ = jetJ->next) { 00488 if (jetJ->NN == tail) {jetJ->NN = jetA;} 00489 } 00490 } 00491 00492 //for (int i = 0; i < n; i++) { 00493 // if (briefjets[i].NN-briefjets >= n && briefjets[i].NN != NULL) {cout <<"YOU MUST BE CRAZY for n ="<<n<<", i = "<<i<<", NN = "<<briefjets[i].NN-briefjets<<"\n";} 00494 //} 00495 00496 00497 if (jetB != NULL) {diJ[jetB-head] = _bj_diJ(jetB);} 00498 //cout << briefjets[95].NN-briefjets<<" "<<briefjets[95].NN_dist <<"\n"; 00499 00500 } 00501 00502 // final cleaning up; 00503 delete[] diJ; 00504 delete[] briefjets; 00505 }
| void fastjet::ClusterSequence::_tj_set_jetinfo | ( | TiledJet *const | jet, | |
| const int | _jets_index | |||
| ) | [inline, private] |
Definition at line 191 of file ClusterSequence_TiledN2.cc.
References _tile_index(), _tiles, fastjet::ClusterSequence::Tile::head, fastjet::ClusterSequence::TiledJet::next, and fastjet::ClusterSequence::TiledJet::previous.
Referenced by _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().
00192 { 00193 // first call the generic setup 00194 _bj_set_jetinfo<>(jet, _jets_index); 00195 00196 // Then do the setup specific to the tiled case. 00197 00198 // Find out which tile it belonds to 00199 jet->tile_index = _tile_index(jet->eta, jet->phi); 00200 00201 // Insert it into the tile's linked list of jets 00202 Tile * tile = &_tiles[jet->tile_index]; 00203 jet->previous = NULL; 00204 jet->next = tile->head; 00205 if (jet->next != NULL) {jet->next->previous = jet;} 00206 tile->head = jet; 00207 }
| void fastjet::ClusterSequence::_transfer_input_jets | ( | const std::vector< L > & | pseudojets | ) | [inline, protected] |
transfer the vector<L> of input jets into our own vector<PseudoJet> _jets (with some reserved space for future growth).
Definition at line 777 of file ClusterSequence.hh.
References _jets.
Referenced by ClusterSequence().
00778 { 00779 00780 // this will ensure that we can point to jets without difficulties 00781 // arising. 00782 _jets.reserve(pseudojets.size()*2); 00783 00784 // insert initial jets this way so that any type L that can be 00785 // converted to a pseudojet will work fine (basically PseudoJet 00786 // and any type that has [] subscript access to the momentum 00787 // components, such as CLHEP HepLorentzVector). 00788 for (unsigned int i = 0; i < pseudojets.size(); i++) { 00789 _jets.push_back(pseudojets[i]);} 00790 00791 }
| void fastjet::ClusterSequence::add_constituents | ( | const PseudoJet & | jet, | |
| std::vector< PseudoJet > & | subjet_vector | |||
| ) | const |
add on to subjet_vector the constituents of jet (for internal use mainly)
Referenced by constituents().
return a vector of the particles that make up jet
Definition at line 815 of file ClusterSequence.cc.
References add_constituents().
Referenced by main(), print_jet(), print_jets(), and print_jets_and_sub().
00815 { 00816 vector<PseudoJet> subjets; 00817 add_constituents(jet, subjets); 00818 return subjets; 00819 }
| double fastjet::ClusterSequence::exclusive_dmerge | ( | const int & | njets | ) | const |
return the dmin corresponding to the recombination that went from n+1 to n jets (sometimes known as d_{n n+1}).
return the dmin corresponding to the recombination that went from n+1 to n jets
Definition at line 535 of file ClusterSequence.cc.
References _history, and _initial_n.
Referenced by main().
00535 { 00536 assert(njets >= 0); 00537 if (njets >= _initial_n) {return 0.0;} 00538 return _history[2*_initial_n-njets-1].dij; 00539 }
| double fastjet::ClusterSequence::exclusive_dmerge_max | ( | const int & | njets | ) | const |
return the maximum of the dmin encountered during all recombinations up to the one that led to an n-jet final state; identical to exclusive_dmerge, except in cases where the dmin do not increase monotonically.
return the maximum of the dmin encountered during all recombinations up to the one that led to an n-jet final state; identical to exclusive_dmerge, except in cases where the dmin do not increase monotonically.
Definition at line 547 of file ClusterSequence.cc.
References _history, and _initial_n.
00547 { 00548 assert(njets >= 0); 00549 if (njets >= _initial_n) {return 0.0;} 00550 return _history[2*_initial_n-njets-1].max_dij_so_far; 00551 }
| vector< PseudoJet > fastjet::ClusterSequence::exclusive_jets | ( | const int & | njets | ) | const |
return a vector of all jets when the event is clustered (in the exclusive sense) to exactly njets.
Definition at line 463 of file ClusterSequence.cc.
References _history, _initial_n, _jet_def, _jets, _n_exclusive_warnings, fastjet::cambridge_algorithm, fastjet::ee_genkt_algorithm, fastjet::ee_kt_algorithm, fastjet::JetDefinition::extra_param(), fastjet::genkt_algorithm, fastjet::JetDefinition::jet_algorithm(), jets(), fastjet::kt_algorithm, fastjet::JetDefinition::plugin(), and fastjet::plugin_algorithm.
00463 { 00464 00465 // make sure the user does not ask for more than jets than there 00466 // were particles in the first place. 00467 assert (njets <= _initial_n); 00468 00469 // provide a warning when extracting exclusive jets for algorithms 00470 // that does not support it explicitly. 00471 // Native algorithm that support it are: kt, ee_kt, cambridge, 00472 // genkt and ee_genkt (both with p>=0) 00473 // For plugins, we check Plugin::exclusive_sequence_meaningful() 00474 if (( _jet_def.jet_algorithm() != kt_algorithm) && 00475 ( _jet_def.jet_algorithm() != cambridge_algorithm) && 00476 ( _jet_def.jet_algorithm() != ee_kt_algorithm) && 00477 (((_jet_def.jet_algorithm() != genkt_algorithm) && 00478 (_jet_def.jet_algorithm() != ee_genkt_algorithm)) || 00479 (_jet_def.extra_param() <0)) && 00480 ((_jet_def.jet_algorithm() != plugin_algorithm) || 00481 (!_jet_def.plugin()->exclusive_sequence_meaningful())) && 00482 (_n_exclusive_warnings < 5)) { 00483 _n_exclusive_warnings++; 00484 cerr << "FastJet WARNING: dcut and exclusive jets for jet-finders other than kt should be interpreted with care." << endl; 00485 } 00486 00487 00488 // calculate the point where we have to stop the clustering. 00489 // relation between stop_point, njets assumes one extra jet disappears 00490 // at each clustering. 00491 int stop_point = 2*_initial_n - njets; 00492 00493 // some sanity checking to make sure that e+e- does not give us 00494 // surprises (should we ever implement e+e-)... 00495 if (2*_initial_n != static_cast<int>(_history.size())) { 00496 ostringstream err; 00497 err << "2*_initial_n != _history.size() -- this endangers internal assumptions!\n"; 00498 throw Error(err.str()); 00499 //assert(false); 00500 } 00501 00502 // now go forwards and reconstitute the jets that we have -- 00503 // basically for any history element, see if the parent jets to 00504 // which it refers were created before the stopping point -- if they 00505 // were then add them to the list, otherwise they are subsequent 00506 // recombinations of the jets that we are looking for. 00507 vector<PseudoJet> jets; 00508 for (unsigned int i = stop_point; i < _history.size(); i++) { 00509 int parent1 = _history[i].parent1; 00510 if (parent1 < stop_point) { 00511 jets.push_back(_jets[_history[parent1].jetp_index]); 00512 } 00513 int parent2 = _history[i].parent2; 00514 if (parent2 < stop_point && parent2 > 0) { 00515 jets.push_back(_jets[_history[parent2].jetp_index]); 00516 } 00517 00518 } 00519 00520 // sanity check... 00521 if (static_cast<int>(jets.size()) != njets) { 00522 ostringstream err; 00523 err << "ClusterSequence::exclusive_jets: size of returned vector (" 00524 <<jets.size()<<") does not coincide with requested number of jets (" 00525 <<njets<<")"; 00526 throw Error(err.str()); 00527 } 00528 00529 return jets; 00530 }
| vector< PseudoJet > fastjet::ClusterSequence::exclusive_jets | ( | const double & | dcut | ) | const |
return a vector of all jets (in the sense of the exclusive algorithm) that would be obtained when running the algorithm with the given dcut.
Definition at line 455 of file ClusterSequence.cc.
References n_exclusive_jets().
Referenced by main(), and print_jets_and_sub().
00455 { 00456 int njets = n_exclusive_jets(dcut); 00457 return exclusive_jets(njets); 00458 }
| std::vector<PseudoJet> fastjet::ClusterSequence::exclusive_jets_ycut | ( | double | ycut | ) | const [inline] |
the exclusive jets obtained at the given ycut
Definition at line 145 of file ClusterSequence.hh.
Referenced by main().
00145 { 00146 int njets = n_exclusive_jets_ycut(ycut); 00147 return exclusive_jets(njets); 00148 }
| double fastjet::ClusterSequence::exclusive_subdmerge | ( | const PseudoJet & | jet, | |
| int | nsub | |||
| ) | const |
return the dij that was present in the merging nsub+1 -> nsub subjets inside this jet.
return the dij that was present in the merging nsub+1 -> nsub subjets inside this jet.
If the jet has nsub or fewer constituents, it will return 0.
will be zero if nconst <= nsub, since highest will be an original particle have zero dij
Definition at line 619 of file ClusterSequence.cc.
References get_subhist_set().
00619 { 00620 set<const history_element*> subhist; 00621 00622 // get the set of history elements that correspond to subjets at 00623 // scale dcut 00624 get_subhist_set(subhist, jet, -1.0, nsub); 00625 00626 set<const history_element*>::iterator highest = subhist.end(); 00627 highest--; 00630 return (*highest)->dij; 00631 }
| double fastjet::ClusterSequence::exclusive_subdmerge_max | ( | const PseudoJet & | jet, | |
| int | nsub | |||
| ) | const |
return the maximum dij that occurred in the whole event at the stage that the nsub+1 -> nsub merge of subjets occurred inside this jet.
return the maximum dij that occurred in the whole event at the stage that the nsub+1 -> nsub merge of subjets occurred inside this jet.
If the jet has nsub or fewer constituents, it will return 0.
will be zero if nconst <= nsub, since highest will be an original particle have zero dij
Definition at line 640 of file ClusterSequence.cc.
References get_subhist_set().
Referenced by print_jets_and_sub().
00640 { 00641 00642 set<const history_element*> subhist; 00643 00644 // get the set of history elements that correspond to subjets at 00645 // scale dcut 00646 get_subhist_set(subhist, jet, -1.0, nsub); 00647 00648 set<const history_element*>::iterator highest = subhist.end(); 00649 highest--; 00652 return (*highest)->max_dij_so_far; 00653 }
| std::vector< PseudoJet > fastjet::ClusterSequence::exclusive_subjets | ( | const PseudoJet & | jet, | |
| int | n | |||
| ) | const |
return the list of subjets obtained by unclustering the supplied jet down to n subjets (or all constituents if there are fewer than n).
return the list of subjets obtained by unclustering the supplied jet down to n subjets (or all constituents if there are fewer than n).
requires n ln n time
Definition at line 595 of file ClusterSequence.cc.
00595 { 00596 00597 set<const history_element*> subhist; 00598 00599 // get the set of history elements that correspond to subjets at 00600 // scale dcut 00601 get_subhist_set(subhist, jet, -1.0, n); 00602 00603 // now transfer this into a sequence of jets 00604 vector<PseudoJet> subjets; 00605 subjets.reserve(subhist.size()); 00606 for (set<const history_element*>::iterator elem = subhist.begin(); 00607 elem != subhist.end(); elem++) { 00608 subjets.push_back(_jets[(*elem)->jetp_index]); 00609 } 00610 return subjets; 00611 }
| std::vector< PseudoJet > fastjet::ClusterSequence::exclusive_subjets | ( | const PseudoJet & | jet, | |
| const double & | dcut | |||
| ) | const |
return a vector of all subjets of the current jet (in the sense of the exclusive algorithm) that would be obtained when running the algorithm with the given dcut.
return a vector of all subjets of the current jet (in the sense of the exclusive algorithm) that would be obtained when running the algorithm with the given dcut.
Time taken is O(m ln m), where m is the number of subjets that are found. If m gets to be of order of the total number of constituents in the jet, this could be substantially slower than just getting that list of constituents.
Definition at line 559 of file ClusterSequence.cc.
Referenced by main(), and print_jets_and_sub().
00559 { 00560 00561 set<const history_element*> subhist; 00562 00563 // get the set of history elements that correspond to subjets at 00564 // scale dcut 00565 get_subhist_set(subhist, jet, dcut, 0); 00566 00567 // now transfer this into a sequence of jets 00568 vector<PseudoJet> subjets; 00569 subjets.reserve(subhist.size()); 00570 for (set<const history_element*>::iterator elem = subhist.begin(); 00571 elem != subhist.end(); elem++) { 00572 subjets.push_back(_jets[(*elem)->jetp_index]); 00573 } 00574 return subjets; 00575 }
| double fastjet::ClusterSequence::exclusive_ymerge | ( | int | njets | ) | const [inline] |
return the ymin corresponding to the recombination that went from n+1 to n jets (sometimes known as y_{n n+1}).
Definition at line 136 of file ClusterSequence.hh.
Referenced by main().
00136 {return exclusive_dmerge(njets) / Q2();}
| double fastjet::ClusterSequence::exclusive_ymerge_max | ( | int | njets | ) | const [inline] |
same as exclusive_dmerge_max, but normalised to squared total energy
Definition at line 139 of file ClusterSequence.hh.
References exclusive_ymerge_max().
Referenced by exclusive_ymerge_max().
00139 {return exclusive_ymerge_max(njets)/Q2();}
| const Extras* fastjet::ClusterSequence::extras | ( | ) | const [inline] |
returns a pointer to the extras object (may be null)
Definition at line 320 of file ClusterSequence.hh.
Referenced by main().
00320 {return _extras.get();}
| void fastjet::ClusterSequence::get_subhist_set | ( | std::set< const history_element * > & | subhist, | |
| const PseudoJet & | jet, | |||
| double | dcut, | |||
| int | maxjet | |||
| ) | const [protected] |
set subhist to be a set pointers to history entries corresponding to the subjets of this jet; one stops going working down through the subjets either when
- there is no further to go
- one has found maxjet entries
- max_dij_so_far <= dcut By setting maxjet=0 one can use just dcut; by setting dcut<0 one can use jet maxjet
Referenced by exclusive_subdmerge(), exclusive_subdmerge_max(), and n_exclusive_subjets().
Version of has_child that sets a pointer to the child if the child exists;.
Definition at line 767 of file ClusterSequence.cc.
References _history, _jets, fastjet::ClusterSequence::history_element::child, and fastjet::PseudoJet::cluster_hist_index().
00767 { 00768 00769 const history_element & hist = _history[jet.cluster_hist_index()]; 00770 00771 // check that this jet has a child and that the child corresponds to 00772 // a true jet [RETHINK-IF-CHANGE-NUMBERING: what is the right 00773 // behaviour if the child is the same jet but made inclusive...?] 00774 if (hist.child >= 0 && _history[hist.child].jetp_index >= 0) { 00775 childp = &(_jets[_history[hist.child].jetp_index]); 00776 return true; 00777 } else { 00778 childp = NULL; 00779 return false; 00780 } 00781 }
if the jet has a child then return true and give the child jet otherwise return false and set the child to zero
if the jet has a child then return true and give the child jet otherwise return false and set the child to zero
Definition at line 745 of file ClusterSequence.cc.
Referenced by object_in_jet().
00745 { 00746 00747 //const history_element & hist = _history[jet.cluster_hist_index()]; 00748 // 00749 //if (hist.child >= 0) { 00750 // child = _jets[_history[hist.child].jetp_index]; 00751 // return true; 00752 //} else { 00753 // child = PseudoJet(0.0,0.0,0.0,0.0); 00754 // return false; 00755 //} 00756 const PseudoJet * childp; 00757 bool res = has_child(jet, childp); 00758 if (res) { 00759 child = *childp; 00760 return true; 00761 } else { 00762 child = PseudoJet(0.0,0.0,0.0,0.0); 00763 return false; 00764 } 00765 }
| bool fastjet::ClusterSequence::has_parents | ( | const PseudoJet & | jet, | |
| PseudoJet & | parent1, | |||
| PseudoJet & | parent2 | |||
| ) | const |
if the jet has parents in the clustering, it returns true and sets parent1 and parent2 equal to them.
if the jet has parents in the clustering, it returns true and sets parent1 and parent2 equal to them.
if it has no parents it returns false and sets parent1 and parent2 to zero
if it has no parents it returns false and sets parent1 and parent2 to zero
Definition at line 719 of file ClusterSequence.cc.
References _history, _jets, fastjet::PseudoJet::cluster_hist_index(), fastjet::ClusterSequence::history_element::parent1, fastjet::ClusterSequence::history_element::parent2, and fastjet::PseudoJet::perp2().
Referenced by main().
00720 { 00721 00722 const history_element & hist = _history[jet.cluster_hist_index()]; 00723 00724 // make sure we do not run into any unexpected situations -- 00725 // i.e. both parents valid, or neither 00726 assert ((hist.parent1 >= 0 && hist.parent2 >= 0) || 00727 (hist.parent1 < 0 && hist.parent2 < 0)); 00728 00729 if (hist.parent1 < 0) { 00730 parent1 = PseudoJet(0.0,0.0,0.0,0.0); 00731 parent2 = parent1; 00732 return false; 00733 } else { 00734 parent1 = _jets[_history[hist.parent1].jetp_index]; 00735 parent2 = _jets[_history[hist.parent2].jetp_index]; 00736 // order the parents in decreasing pt 00737 if (parent1.perp2() < parent2.perp2()) swap(parent1,parent2); 00738 return true; 00739 } 00740 }
if this jet has a child (and so a partner) return true and give the partner, otherwise return false and set the partner to zero
if this jet has a child (and so a partner) return true and give the partner, otherwise return false and set the partner to zero
Definition at line 788 of file ClusterSequence.cc.
References _history, _jets, fastjet::ClusterSequence::history_element::child, fastjet::PseudoJet::cluster_hist_index(), fastjet::ClusterSequence::history_element::parent1, and fastjet::ClusterSequence::history_element::parent2.
00789 { 00790 00791 const history_element & hist = _history[jet.cluster_hist_index()]; 00792 00793 // make sure we have a child and that the child does not correspond 00794 // to a clustering with the beam (or some other invalid quantity) 00795 if (hist.child >= 0 && _history[hist.child].parent2 >= 0) { 00796 const history_element & child_hist = _history[hist.child]; 00797 if (child_hist.parent1 == jet.cluster_hist_index()) { 00798 // partner will be child's parent2 -- for iB clustering 00799 // parent2 will not be valid 00800 partner = _jets[_history[child_hist.parent2].jetp_index]; 00801 } else { 00802 // partner will be child's parent1 00803 partner = _jets[_history[child_hist.parent1].jetp_index]; 00804 } 00805 return true; 00806 } else { 00807 partner = PseudoJet(0.0,0.0,0.0,0.0); 00808 return false; 00809 } 00810 }
| const std::vector< ClusterSequence::history_element > & fastjet::ClusterSequence::history | ( | ) | const [inline] |
allow the user to access the raw internal history.
This is present (as for jets()) in part so that protected derived classes can access this information about other ClusterSequences.
A user who wishes to follow the details of the ClusterSequence can also make use of this information (and should consult the history_element documentation for more information), but should be aware that these internal structures may evolve in future FastJet versions.
Definition at line 830 of file ClusterSequence.hh.
References _history.
Referenced by fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history(), main(), and fastjet::NestedDefsPlugin::run_clustering().
00830 { 00831 return _history; 00832 }
| vector< PseudoJet > fastjet::ClusterSequence::inclusive_jets | ( | const double & | ptmin = 0.0 |
) | const |
return a vector of all jets (in the sense of the inclusive algorithm) with pt >= ptmin.
Time taken should be of the order of the number of jets returned.
Definition at line 384 of file ClusterSequence.cc.
References _history, _jet_algorithm, _jets, fastjet::antikt_algorithm, BeamJet, fastjet::cambridge_algorithm, fastjet::cambridge_for_passive_algorithm, fastjet::ee_genkt_algorithm, fastjet::ee_kt_algorithm, fastjet::genkt_algorithm, jets(), fastjet::kt_algorithm, fastjet::PseudoJet::perp2(), and fastjet::plugin_algorithm.
Referenced by fastjet::ClusterSequenceAreaBase::empty_area(), fastjet::ClusterSequenceAreaBase::get_median_rho_and_sigma(), main(), fastjet::ClusterSequenceAreaBase::parabolic_pt_per_unit_area(), print_jets_and_sub(), fastjet::ClusterSequenceActiveArea::pt_per_unit_area(), run_jet_finder(), and fastjet::ClusterSequenceAreaBase::subtracted_jets().
00384 { 00385 double dcut = ptmin*ptmin; 00386 int i = _history.size() - 1; // last jet 00387 vector<PseudoJet> jets; 00388 if (_jet_algorithm == kt_algorithm) { 00389 while (i >= 0) { 00390 // with our specific definition of dij and diB (i.e. R appears only in 00391 // dij), then dij==diB is the same as the jet.perp2() and we can exploit 00392 // this in selecting the jets... 00393 if (_history[i].max_dij_so_far < dcut) {break;} 00394 if (_history[i].parent2 == BeamJet && _history[i].dij >= dcut) { 00395 // for beam jets 00396 int parent1 = _history[i].parent1; 00397 jets.push_back(_jets[_history[parent1].jetp_index]);} 00398 i--; 00399 } 00400 } else if (_jet_algorithm == cambridge_algorithm) { 00401 while (i >= 0) { 00402 // inclusive jets are all at end of clustering sequence in the 00403 // Cambridge algorithm -- so if we find a non-exclusive jet, then 00404 // we can exit 00405 if (_history[i].parent2 != BeamJet) {break;} 00406 int parent1 = _history[i].parent1; 00407 const PseudoJet & jet = _jets[_history[parent1].jetp_index]; 00408 if (jet.perp2() >= dcut) {jets.push_back(jet);} 00409 i--; 00410 } 00411 } else if (_jet_algorithm == plugin_algorithm 00412 || _jet_algorithm == ee_kt_algorithm 00413 || _jet_algorithm == antikt_algorithm 00414 || _jet_algorithm == genkt_algorithm 00415 || _jet_algorithm == ee_genkt_algorithm 00416 || _jet_algorithm == cambridge_for_passive_algorithm) { 00417 // for inclusive jets with a plugin algorithm, we make no 00418 // assumptions about anything (relation of dij to momenta, 00419 // ordering of the dij, etc.) 00420 while (i >= 0) { 00421 if (_history[i].parent2 == BeamJet) { 00422 int parent1 = _history[i].parent1; 00423 const PseudoJet & jet = _jets[_history[parent1].jetp_index]; 00424 if (jet.perp2() >= dcut) {jets.push_back(jet);} 00425 } 00426 i--; 00427 } 00428 } else {throw Error("cs::inclusive_jets(...): Unrecognized jet algorithm");} 00429 return jets; 00430 }
| const JetDefinition& fastjet::ClusterSequence::jet_def | ( | ) | const [inline] |
return a reference to the jet definition
Definition at line 264 of file ClusterSequence.hh.
Referenced by _bj_set_jetinfo(), fastjet::ClusterSequenceAreaBase::_check_jet_alg_good_for_median(), _initialise_and_run(), fastjet::ClusterSequence1GhostPassiveArea::_run_1GPA(), fastjet::ClusterSequenceActiveArea::_run_AA(), fastjet::ClusterSequenceArea::_warn_if_range_unsuitable(), fastjet::ClusterSequencePassiveArea::empty_area(), jet_scale_for_algorithm(), and print_jets_and_sub().
00264 {return _jet_def;}
| double fastjet::ClusterSequence::jet_scale_for_algorithm | ( | const PseudoJet & | jet | ) | const |
returns the scale associated with a jet as required for this clustering algorithm (kt^2 for the kt-algorithm, 1 for the Cambridge algorithm).
[May become virtual at some point]
Definition at line 319 of file ClusterSequence.cc.
References _jet_algorithm, _jet_def, fastjet::antikt_algorithm, fastjet::cambridge_algorithm, fastjet::cambridge_for_passive_algorithm, fastjet::JetDefinition::extra_param(), fastjet::genkt_algorithm, jet_def(), fastjet::PseudoJet::kt2(), and fastjet::kt_algorithm.
Referenced by _add_ktdistance_to_map(), _bj_set_jetinfo(), and _really_dumb_cluster().
00320 { 00321 if (_jet_algorithm == kt_algorithm) {return jet.kt2();} 00322 else if (_jet_algorithm == cambridge_algorithm) {return 1.0;} 00323 else if (_jet_algorithm == antikt_algorithm) { 00324 double kt2=jet.kt2(); 00325 return kt2 > 1e-300 ? 1.0/kt2 : 1e300; 00326 } else if (_jet_algorithm == genkt_algorithm) { 00327 double kt2 = jet.kt2(); 00328 double p = jet_def().extra_param(); 00329 if (p <= 0 && kt2 < 1e-300) kt2 = 1e-300; // dodgy safety check 00330 return pow(kt2, p); 00331 } else if (_jet_algorithm == cambridge_for_passive_algorithm) { 00332 double kt2 = jet.kt2(); 00333 double lim = _jet_def.extra_param(); 00334 if (kt2 < lim*lim && kt2 != 0.0) { 00335 return 1.0/kt2; 00336 } else {return 1.0;} 00337 } else {throw Error("Unrecognised jet algorithm");} 00338 }
| const std::vector< PseudoJet > & fastjet::ClusterSequence::jets | ( | ) | const [inline] |
allow the user to access the internally stored _jets() array, which contains both the initial particles and the various intermediate and final stages of recombination.
The first n_particles() entries are the original particles, in the order in which they were supplied to the ClusterSequence constructor. It can be useful to access them for example when examining whether a given input object is part of a specific jet, via the objects_in_jet(...) member function (which only takes PseudoJets that are registered in the ClusterSequence).
One of the other (internal uses) is related to the fact because we don't seem to be able to access protected elements of the class for an object that is not "this" (at least in case where "this" is of a slightly different kind from the object, both derived from ClusterSequence).
Definition at line 826 of file ClusterSequence.hh.
References _jets.
Referenced by exclusive_jets(), inclusive_jets(), fastjet::TrackJetPlugin::run_clustering(), fastjet::SISConeSphericalPlugin::run_clustering(), fastjet::SISConePlugin::run_clustering(), fastjet::PxConePlugin::run_clustering(), fastjet::NestedDefsPlugin::run_clustering(), fastjet::JadePlugin::run_clustering(), fastjet::EECambridgePlugin::run_clustering(), fastjet::D0RunIIConePlugin::run_clustering(), fastjet::CMSIterativeConePlugin::run_clustering(), fastjet::CDFMidPointPlugin::run_clustering(), fastjet::CDFJetCluPlugin::run_clustering(), fastjet::ATLASConePlugin::run_clustering(), and fastjet::ClusterSequenceAreaBase::subtracted_jets().
00826 { 00827 return _jets; 00828 }
| int fastjet::ClusterSequence::n_exclusive_jets | ( | const double & | dcut | ) | const |
return the number of jets (in the sense of the exclusive algorithm) that would be obtained when running the algorithm with the given dcut.
Definition at line 436 of file ClusterSequence.cc.
References _history, and _initial_n.
Referenced by exclusive_jets().
00436 { 00437 00438 // first locate the point where clustering would have stopped (i.e. the 00439 // first time max_dij_so_far > dcut) 00440 int i = _history.size() - 1; // last jet 00441 while (i >= 0) { 00442 if (_history[i].max_dij_so_far <= dcut) {break;} 00443 i--; 00444 } 00445 int stop_point = i + 1; 00446 // relation between stop_point, njets assumes one extra jet disappears 00447 // at each clustering. 00448 int njets = 2*_initial_n - stop_point; 00449 return njets; 00450 }
| int fastjet::ClusterSequence::n_exclusive_jets_ycut | ( | double | ycut | ) | const [inline] |
the number of exclusive jets at the given ycut
Definition at line 142 of file ClusterSequence.hh.
00142 {return n_exclusive_jets(ycut*Q2());}
| int fastjet::ClusterSequence::n_exclusive_subjets | ( | const PseudoJet & | jet, | |
| const double & | dcut | |||
| ) | const |
return the size of exclusive_subjets(.
..); still n ln n with same coefficient, but marginally more efficient than manually taking exclusive_subjets.size()
..); still n ln n with same coefficient, but marginally more efficient than manually taking exclusive_subjets.size()
Definition at line 581 of file ClusterSequence.cc.
References get_subhist_set().
00582 { 00583 set<const history_element*> subhist; 00584 // get the set of history elements that correspond to subjets at 00585 // scale dcut 00586 get_subhist_set(subhist, jet, dcut, 0); 00587 return subhist.size(); 00588 }
| unsigned int fastjet::ClusterSequence::n_particles | ( | ) | const [inline] |
returns the number of particles that were provided to the clustering algorithm (helps the user find their way around the history and jets objects if they weren't paying attention beforehand).
Definition at line 834 of file ClusterSequence.hh.
References _initial_n.
Referenced by _initialise_and_run(), fastjet::ClusterSequenceVoronoiArea::_initializeVA(), unclustered_particles(), and unique_history_order().
00834 {return _initial_n;}
| bool fastjet::ClusterSequence::object_in_jet | ( | const PseudoJet & | object, | |
| const PseudoJet & | jet | |||
| ) | const |
returns true iff the object is included in the jet.
NB: this is only sensible if the object is already registered within the cluster sequence, so you cannot use it with an input particle to the CS (since the particle won't have the history index set properly).
For nice clustering structures it should run in O(ln(N)) time but in worst cases (certain cone plugins) it can take O(n) time, where n is the number of particles in the jet.
Definition at line 696 of file ClusterSequence.cc.
References _potentially_valid(), fastjet::PseudoJet::cluster_hist_index(), and has_child().
00697 { 00698 00699 // make sure the object conceivably belongs to this clustering 00700 // sequence 00701 assert(_potentially_valid(object) && _potentially_valid(jet)); 00702 00703 const PseudoJet * this_object = &object; 00704 const PseudoJet * childp; 00705 while(true) { 00706 if (this_object->cluster_hist_index() == jet.cluster_hist_index()) { 00707 return true; 00708 } else if (has_child(*this_object, childp)) {this_object = childp;} 00709 else {return false;} 00710 } 00711 }
| std::vector<int> fastjet::ClusterSequence::particle_jet_indices | ( | const std::vector< PseudoJet > & | ) | const |
returns a vector of size n_particles() which indicates, for each of the initial particles (in the order in which they were supplied), which of the supplied jets it belongs to; if it does not belong to any of the supplied jets, the index is set to -1;
| bool fastjet::ClusterSequence::plugin_activated | ( | ) | const [inline] |
returns true when the plugin is allowed to run the show.
Definition at line 317 of file ClusterSequence.hh.
00317 {return _plugin_activated;}
| void fastjet::ClusterSequence::plugin_associate_extras | ( | std::auto_ptr< Extras > | extras_in | ) | [inline] |
the plugin can associated some extra information with the ClusterSequence object by calling this function
Definition at line 312 of file ClusterSequence.hh.
Referenced by fastjet::SISConeSphericalPlugin::run_clustering(), and fastjet::SISConePlugin::run_clustering().
00312 { 00313 _extras = extras_in; 00314 }
| void fastjet::ClusterSequence::plugin_record_iB_recombination | ( | int | jet_i, | |
| double | diB | |||
| ) | [inline] |
record the fact that there has been a recombination between jets()[jet_i] and the beam, with the specified diB; when looking for inclusive jets, any iB recombination will returned to the user as a jet.
Definition at line 296 of file ClusterSequence.hh.
Referenced by fastjet::TrackJetPlugin::run_clustering(), fastjet::SISConeSphericalPlugin::run_clustering(), fastjet::SISConePlugin::run_clustering(), fastjet::PxConePlugin::run_clustering(), fastjet::NestedDefsPlugin::run_clustering(), fastjet::JadePlugin::run_clustering(), fastjet::EECambridgePlugin::run_clustering(), fastjet::D0RunIIConePlugin::run_clustering(), fastjet::CMSIterativeConePlugin::run_clustering(), fastjet::CDFMidPointPlugin::run_clustering(), fastjet::CDFJetCluPlugin::run_clustering(), and fastjet::ATLASConePlugin::run_clustering().
00296 { 00297 assert(plugin_activated()); 00298 _do_iB_recombination_step(jet_i, diB); 00299 }
| void fastjet::ClusterSequence::plugin_record_ij_recombination | ( | int | jet_i, | |
| int | jet_j, | |||
| double | dij, | |||
| const PseudoJet & | newjet, | |||
| int & | newjet_k | |||
| ) |
as for the simpler variant of plugin_record_ij_recombination, except that the new jet is attributed the momentum and user_index of newjet
Definition at line 369 of file ClusterSequence.cc.
References _jets, and plugin_record_ij_recombination().
00371 { 00372 00373 plugin_record_ij_recombination(jet_i, jet_j, dij, newjet_k); 00374 00375 // now transfer newjet into place 00376 int tmp_index = _jets[newjet_k].cluster_hist_index(); 00377 _jets[newjet_k] = newjet; 00378 _jets[newjet_k].set_cluster_hist_index(tmp_index); 00379 }
| void fastjet::ClusterSequence::plugin_record_ij_recombination | ( | int | jet_i, | |
| int | jet_j, | |||
| double | dij, | |||
| int & | newjet_k | |||
| ) | [inline] |
record the fact that there has been a recombination between jets()[jet_i] and jets()[jet_k], with the specified dij, and return the index (newjet_k) allocated to the new jet, whose momentum is assumed to be the 4-vector sum of that of jet_i and jet_j
Definition at line 279 of file ClusterSequence.hh.
Referenced by plugin_record_ij_recombination(), fastjet::TrackJetPlugin::run_clustering(), fastjet::SISConeSphericalPlugin::run_clustering(), fastjet::SISConePlugin::run_clustering(), fastjet::PxConePlugin::run_clustering(), fastjet::NestedDefsPlugin::run_clustering(), fastjet::JadePlugin::run_clustering(), fastjet::EECambridgePlugin::run_clustering(), fastjet::D0RunIIConePlugin::run_clustering(), fastjet::CMSIterativeConePlugin::run_clustering(), fastjet::CDFMidPointPlugin::run_clustering(), fastjet::CDFJetCluPlugin::run_clustering(), and fastjet::ATLASConePlugin::run_clustering().
00280 { 00281 assert(plugin_activated()); 00282 _do_ij_recombination_step(jet_i, jet_j, dij, newjet_k); 00283 }
| void fastjet::ClusterSequence::plugin_simple_N2_cluster | ( | ) | [inline] |
allows a plugin to run a templated clustering (nearest-neighbour heuristic)
This has N^2 behaviour on "good" distance, but a worst case behaviour of N^3 (and many algs trigger the worst case behaviour)
For more details on how this works, see GenBriefJet below
Definition at line 329 of file ClusterSequence.hh.
00329 { 00330 assert(plugin_activated()); 00331 _simple_N2_cluster<GBJ>(); 00332 }
| void fastjet::ClusterSequence::print_jets_for_root | ( | const std::vector< PseudoJet > & | jets, | |
| const std::string & | filename, | |||
| const std::string & | comment = "" | |||
| ) | const |
print jets for root to the file labelled filename, with an optional comment at the beginning
Definition at line 849 of file ClusterSequence.cc.
References print_jets_for_root().
00851 { 00852 std::ofstream ostr(filename.c_str()); 00853 if (comment != "") ostr << "# " << comment << endl; 00854 print_jets_for_root(jets, ostr); 00855 }
| void fastjet::ClusterSequence::print_jets_for_root | ( | const std::vector< PseudoJet > & | jets, | |
| std::ostream & | ostr = std::cout | |||
| ) | const |
output the supplied vector of jets in a format that can be read by an appropriate root script; the format is: jet-n jet-px jet-py jet-pz jet-E particle-n particle-rap particle-phi particle-pt particle-n particle-rap particle-phi particle-pt
.. END ... [i.e. above repeated]
Referenced by print_jets(), and print_jets_for_root().
| double fastjet::ClusterSequence::Q | ( | ) | const [inline] |
returns the sum of all energies in the event (relevant mainly for e+e-)
Definition at line 192 of file ClusterSequence.hh.
00192 {return _Qtot;}
| double fastjet::ClusterSequence::Q2 | ( | ) | const [inline] |
return Q()^2
Definition at line 194 of file ClusterSequence.hh.
Referenced by fastjet::EECambridgePlugin::run_clustering().
| static void fastjet::ClusterSequence::set_jet_algorithm | ( | JetAlgorithm | jet_algorithm | ) | [inline, static] |
set the default (static) jet finder across all current and future ClusterSequence objects -- deprecated and obsolescent (i.e.
may be suppressed in a future release).
Definition at line 371 of file ClusterSequence.hh.
00371 {_default_jet_algorithm = jet_algorithm;}
| static void fastjet::ClusterSequence::set_jet_finder | ( | JetAlgorithm | jet_algorithm | ) | [inline, static] |
same as above for backward compatibility
Definition at line 373 of file ClusterSequence.hh.
00373 {_default_jet_algorithm = jet_algorithm;}
| string fastjet::ClusterSequence::strategy_string | ( | ) | const |
Definition at line 285 of file ClusterSequence.cc.
References _strategy, fastjet::N2MinHeapTiled, fastjet::N2Plain, fastjet::N2PoorTiled, fastjet::N2Tiled, fastjet::N3Dumb, fastjet::NlnN, fastjet::NlnN3pi, fastjet::NlnN4pi, fastjet::NlnNCam, fastjet::NlnNCam2pi2R, fastjet::NlnNCam4pi, and fastjet::plugin_strategy.
Referenced by _delaunay_cluster(), and main().
00285 { 00286 string strategy; 00287 switch(_strategy) { 00288 case NlnN: 00289 strategy = "NlnN"; break; 00290 case NlnN3pi: 00291 strategy = "NlnN3pi"; break; 00292 case NlnN4pi: 00293 strategy = "NlnN4pi"; break; 00294 case N2Plain: 00295 strategy = "N2Plain"; break; 00296 case N2Tiled: 00297 strategy = "N2Tiled"; break; 00298 case N2MinHeapTiled: 00299 strategy = "N2MinHeapTiled"; break; 00300 case N2PoorTiled: 00301 strategy = "N2PoorTiled"; break; 00302 case N3Dumb: 00303 strategy = "N3Dumb"; break; 00304 case NlnNCam4pi: 00305 strategy = "NlnNCam4pi"; break; 00306 case NlnNCam2pi2R: 00307 strategy = "NlnNCam2pi2R"; break; 00308 case NlnNCam: 00309 strategy = "NlnNCam"; break; // 2piMultD 00310 case plugin_strategy: 00311 strategy = "plugin strategy"; break; 00312 default: 00313 strategy = "Unrecognized"; 00314 } 00315 return strategy; 00316 }
| Strategy fastjet::ClusterSequence::strategy_used | ( | ) | const [inline] |
return the enum value of the strategy used to cluster the event
Definition at line 260 of file ClusterSequence.hh.
Referenced by fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history().
00260 {return _strategy;}
| void fastjet::ClusterSequence::transfer_from_sequence | ( | ClusterSequence & | from_seq | ) |
transfer the sequence contained in other_seq into our own; any plugin "extras" contained in the from_seq will be lost from there.
transfer the sequence contained in other_seq into our own; any plugin "extras" contained in the from_seq will be lost from there.
Definition at line 345 of file ClusterSequence.cc.
References _extras, _history, _initial_n, _invR2, _jet_algorithm, _jet_def, _jets, _plugin_activated, _R2, _Rparam, _strategy, and _writeout_combinations.
Referenced by fastjet::ClusterSequencePassiveArea::_initialise_and_run_PA(), and fastjet::ClusterSequenceArea::initialize_and_run_cswa().
00345 { 00346 00347 // the metadata 00348 _jet_def = from_seq._jet_def ; 00349 _writeout_combinations = from_seq._writeout_combinations ; 00350 _initial_n = from_seq._initial_n ; 00351 _Rparam = from_seq._Rparam ; 00352 _R2 = from_seq._R2 ; 00353 _invR2 = from_seq._invR2 ; 00354 _strategy = from_seq._strategy ; 00355 _jet_algorithm = from_seq._jet_algorithm ; 00356 _plugin_activated = from_seq._plugin_activated ; 00357 00358 // the data 00359 _jets = from_seq._jets; 00360 _history = from_seq._history; 00361 // the following transferse ownership of the extras from the from_seq 00362 _extras = from_seq._extras; 00363 00364 }
| vector< PseudoJet > fastjet::ClusterSequence::unclustered_particles | ( | ) | const |
return the set of particles that have not been clustered.
For kt and cam/aachen algorithms this should always be null, but for cone type algorithms it can be non-null;
Definition at line 1035 of file ClusterSequence.cc.
References _history, _jets, Invalid, and n_particles().
Referenced by fastjet::ClusterSequenceActiveAreaExplicitGhosts::empty_area(), main(), and run_jet_finder().
01035 { 01036 vector<PseudoJet> unclustered; 01037 for (unsigned i = 0; i < n_particles() ; i++) { 01038 if (_history[i].child == Invalid) 01039 unclustered.push_back(_jets[_history[i].jetp_index]); 01040 } 01041 return unclustered; 01042 }
| vector< int > fastjet::ClusterSequence::unique_history_order | ( | ) | const |
routine that returns an order in which to read the history such that clusterings that lead to identical jet compositions but different histories (because of degeneracies in the clustering order) will have matching constituents for each matching entry in the unique_history_order.
The order has the property that an entry's parents will always appear prior to that entry itself.
Roughly speaking the order is such that we first provide all steps that lead to the final jet containing particle 1; then we have the steps that lead to reconstruction of the jet containing the next-lowest-numbered unclustered particle, etc... [see GPS CCN28-12 for more info -- of course a full explanation here would be better...]
Definition at line 980 of file ClusterSequence.cc.
References _extract_tree_children(), _history, fastjet::d0::inline_maths::min(), and n_particles().
Referenced by fastjet::ClusterSequence1GhostPassiveArea::_run_1GPA(), fastjet::ClusterSequenceActiveArea::_run_AA(), and main().
00980 { 00981 00982 // first construct an array that will tell us the lowest constituent 00983 // of a given jet -- this will always be one of the original 00984 // particles, whose order is well defined and so will help us to 00985 // follow the tree in a unique manner. 00986 valarray<int> lowest_constituent(_history.size()); 00987 int hist_n = _history.size(); 00988 lowest_constituent = hist_n; // give it a large number 00989 for (int i = 0; i < hist_n; i++) { 00990 // sets things up for the initial partons 00991 lowest_constituent[i] = min(lowest_constituent[i],i); 00992 // propagates them through to the children of this parton 00993 if (_history[i].child > 0) lowest_constituent[_history[i].child] 00994 = min(lowest_constituent[_history[i].child],lowest_constituent[i]); 00995 } 00996 00997 // establish an array for what we have and have not extracted so far 00998 valarray<bool> extracted(_history.size()); extracted = false; 00999 vector<int> unique_tree; 01000 unique_tree.reserve(_history.size()); 01001 01002 // now work our way through the tree 01003 for (unsigned i = 0; i < n_particles(); i++) { 01004 if (!extracted[i]) { 01005 unique_tree.push_back(i); 01006 extracted[i] = true; 01007 _extract_tree_children(i, extracted, lowest_constituent, unique_tree); 01008 } 01009 } 01010 01011 return unique_tree; 01012 }
Member Data Documentation
JetAlgorithm fastjet::ClusterSequence::_default_jet_algorithm = kt_algorithm [static, protected] |
std::auto_ptr<Extras> fastjet::ClusterSequence::_extras [private] |
bool fastjet::ClusterSequence::_first_time = true [static, private] |
will be set by default to be true for the first run
Definition at line 612 of file ClusterSequence.hh.
Referenced by _print_banner().
std::vector<history_element> fastjet::ClusterSequence::_history [protected] |
this vector will contain the branching history; for each stage, _history[i].jetp_index indicates where to look in the _jets vector to get the physical PseudoJet.
Definition at line 540 of file ClusterSequence.hh.
Referenced by _add_step_to_history(), _CP2DChan_limited_cluster(), _do_Cambridge_inclusive_jets(), _do_iB_recombination_step(), _do_ij_recombination_step(), _fill_initial_history(), fastjet::ClusterSequencePassiveArea::_initialise_and_run_PA(), fastjet::ClusterSequenceVoronoiArea::_initializeVA(), fastjet::ClusterSequenceActiveAreaExplicitGhosts::_post_process(), fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history(), exclusive_dmerge(), exclusive_dmerge_max(), exclusive_jets(), has_child(), has_parents(), has_partner(), history(), inclusive_jets(), n_exclusive_jets(), transfer_from_sequence(), unclustered_particles(), and unique_history_order().
int fastjet::ClusterSequence::_initial_n [protected] |
Definition at line 554 of file ClusterSequence.hh.
Referenced by _CP2DChan_limited_cluster(), _fill_initial_history(), fastjet::ClusterSequenceActiveAreaExplicitGhosts::_post_process(), exclusive_dmerge(), exclusive_dmerge_max(), exclusive_jets(), n_exclusive_jets(), n_particles(), and transfer_from_sequence().
double fastjet::ClusterSequence::_invR2 [protected] |
Definition at line 555 of file ClusterSequence.hh.
Referenced by _add_ktdistance_to_map(), _CP2DChan_cluster(), _CP2DChan_limited_cluster(), _decant_options(), _faster_tiled_N2_cluster(), _initialise_and_run(), _minheap_faster_tiled_N2_cluster(), _really_dumb_cluster(), _tiled_N2_cluster(), and transfer_from_sequence().
JetAlgorithm fastjet::ClusterSequence::_jet_algorithm [protected] |
Definition at line 558 of file ClusterSequence.hh.
Referenced by _bj_set_jetinfo(), _CP2DChan_cluster(), _CP2DChan_cluster_2pi2R(), _decant_options(), _initialise_and_run(), inclusive_jets(), jet_scale_for_algorithm(), and transfer_from_sequence().
JetDefinition fastjet::ClusterSequence::_jet_def [protected] |
Definition at line 483 of file ClusterSequence.hh.
Referenced by _decant_options(), _do_ij_recombination_step(), _fill_initial_history(), _initialise_and_run(), fastjet::ClusterSequencePassiveArea::_initialise_and_run_PA(), fastjet::ClusterSequenceVoronoiArea::_initializeVA(), fastjet::ClusterSequenceActiveAreaExplicitGhosts::_post_process(), exclusive_jets(), jet_scale_for_algorithm(), and transfer_from_sequence().
std::vector<PseudoJet> fastjet::ClusterSequence::_jets [protected] |
This contains the physical PseudoJets; for each PseudoJet one can find the corresponding position in the _history by looking at _jets[i].cluster_hist_index().
Definition at line 534 of file ClusterSequence.hh.
Referenced by fastjet::ClusterSequenceActiveAreaExplicitGhosts::_add_ghosts(), _add_ktdistance_to_map(), _add_step_to_history(), _bj_set_jetinfo(), _CP2DChan_cluster(), _CP2DChan_limited_cluster(), _delaunay_cluster(), _do_iB_recombination_step(), _do_ij_recombination_step(), _faster_tiled_N2_cluster(), _fill_initial_history(), _initialise_and_run(), fastjet::ClusterSequencePassiveArea::_initialise_and_run_PA(), _initialise_tiles(), fastjet::ClusterSequenceVoronoiArea::_initializeVA(), _minheap_faster_tiled_N2_cluster(), fastjet::ClusterSequenceActiveAreaExplicitGhosts::_post_process(), _really_dumb_cluster(), fastjet::ClusterSequenceActiveArea::_resize_and_zero_AA(), fastjet::ClusterSequence1GhostPassiveArea::_run_1GPA(), fastjet::ClusterSequenceActiveArea::_run_AA(), _tiled_N2_cluster(), _transfer_input_jets(), exclusive_jets(), has_child(), has_parents(), has_partner(), inclusive_jets(), jets(), plugin_record_ij_recombination(), transfer_from_sequence(), and unclustered_particles().
int fastjet::ClusterSequence::_n_exclusive_warnings = 0 [static, private] |
record the number of warnings provided about the exclusive algorithm -- so that we don't print it out more than a few times.
Definition at line 617 of file ClusterSequence.hh.
Referenced by exclusive_jets().
int fastjet::ClusterSequence::_n_tiles_phi [private] |
Definition at line 724 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tile_index().
bool fastjet::ClusterSequence::_plugin_activated [private] |
Definition at line 563 of file ClusterSequence.hh.
Referenced by _decant_options(), _initialise_and_run(), and transfer_from_sequence().
double fastjet::ClusterSequence::_Qtot [protected] |
double fastjet::ClusterSequence::_R2 [protected] |
Definition at line 555 of file ClusterSequence.hh.
Referenced by _bj_set_jetinfo(), _bj_set_NN_crosscheck(), _bj_set_NN_nocross(), _decant_options(), _faster_tiled_N2_cluster(), _initialise_and_run(), _minheap_faster_tiled_N2_cluster(), _tiled_N2_cluster(), and transfer_from_sequence().
double fastjet::ClusterSequence::_Rparam [protected] |
Definition at line 555 of file ClusterSequence.hh.
Referenced by _bj_set_jetinfo(), _CP2DChan_cluster_2pi2R(), _CP2DChan_cluster_2piMultD(), _decant_options(), _delaunay_cluster(), _initialise_and_run(), _initialise_tiles(), and transfer_from_sequence().
Strategy fastjet::ClusterSequence::_strategy [protected] |
Definition at line 557 of file ClusterSequence.hh.
Referenced by _decant_options(), _delaunay_cluster(), _initialise_and_run(), fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history(), strategy_string(), and transfer_from_sequence().
double fastjet::ClusterSequence::_tile_size_eta [private] |
Definition at line 723 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tile_index().
double fastjet::ClusterSequence::_tile_size_phi [private] |
Definition at line 723 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tile_index().
std::vector<Tile> fastjet::ClusterSequence::_tiles [private] |
Definition at line 721 of file ClusterSequence.hh.
Referenced by _bj_remove_from_tiles(), _faster_tiled_N2_cluster(), _initialise_tiles(), _minheap_faster_tiled_N2_cluster(), _print_tiles(), _tiled_N2_cluster(), and _tj_set_jetinfo().
double fastjet::ClusterSequence::_tiles_eta_max [private] |
Definition at line 722 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tile_index().
double fastjet::ClusterSequence::_tiles_eta_min [private] |
Definition at line 722 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tile_index().
int fastjet::ClusterSequence::_tiles_ieta_max [private] |
Definition at line 724 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tile_index().
int fastjet::ClusterSequence::_tiles_ieta_min [private] |
Definition at line 724 of file ClusterSequence.hh.
Referenced by _initialise_tiles(), and _tile_index().
bool fastjet::ClusterSequence::_writeout_combinations [protected] |
Definition at line 553 of file ClusterSequence.hh.
Referenced by _add_step_to_history(), _decant_options(), and transfer_from_sequence().
const int fastjet::ClusterSequence::n_tile_neighbours = 9 [static, private] |
number of neighbours that a tile will have (rectangular geometry gives 9 neighbours).
Definition at line 703 of file ClusterSequence.hh.
Referenced by _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().
The documentation for this class was generated from the following files:
