ClusterSequence Class Reference

deals with clustering More...

#include <ClusterSequence.hh>

Inheritance diagram for ClusterSequence:

[legend]

Collaboration diagram for ClusterSequence:

[legend]

List of all members.

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...
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 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
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: decide the range in eta allocate the tiles set up the cross-referencing info between 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 clustering

Definition at line 66 of file ClusterSequence.hh.

---

Member Typedef Documentation

typedef std::pair<double,TwoVertices> ClusterSequence::DijEntry [private]

Definition at line 601 of file ClusterSequence.hh.

typedef std::multimap<double,TwoVertices> ClusterSequence::DistMap [private]

Definition at line 602 of file ClusterSequence.hh.

typedef std::pair<int,int> ClusterSequence::TwoVertices [private]

Definition at line 600 of file ClusterSequence.hh.

---

Member Enumeration Documentation

enum ClusterSequence::JetType

Enumerator:

Invalid
InexistentParent
BeamJet

Definition at line 408 of file ClusterSequence.hh.

{Invalid=-3, InexistentParent = -2, BeamJet = -1};

---

Constructor & Destructor Documentation

ClusterSequence::ClusterSequence

(

)

[inline]

default constructor

Definition at line 72 of file ClusterSequence.hh.

{}

template<class L >

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().

                                                                      {

  // transfer the initial jets (type L) into our own array
  _transfer_input_jets(pseudojets);

  // run the clustering
  _initialise_and_run(R,strategy,writeout_combinations);
}

template<class L >

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().

                                                                      {

  // transfer the initial jets (type L) into our own array
  _transfer_input_jets(pseudojets);

  // run the clustering
  _initialise_and_run(jet_def,writeout_combinations);
}

ClusterSequence::~ClusterSequence

(

)

[virtual]

Definition at line 54 of file ClusterSequence.cc.

{}

---

Member Function Documentation

void 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, and jet_scale_for_algorithm().

Referenced by _delaunay_cluster().

                                                                {
  
  double yiB = jet_scale_for_algorithm(_jets[ii]);
  if (yiB == 0.0) {
    // in this case convention is that we do not worry about distances
    // but directly state that nearest neighbour is beam
    DijMap.insert(DijEntry(yiB,  TwoVertices(ii,-1)));
  } else {
    double DeltaR2 = DNN->NearestNeighbourDistance(ii) * _invR2;
    // Logic of following bit is: only add point to map if it has
    // smaller kt2 than nearest neighbour j (if it has larger kt,
    // then: either it is j's nearest neighbour and then we will
    // include dij when we come to j; or it is not j's nearest
    // neighbour and j will recombine with someone else).
    
    // If DeltaR2 > 1.0 then in any case it will recombine with beam rather
    // than with any neighbours.
    // (put general normalisation here at some point)
    if (DeltaR2 > 1.0) {
      DijMap.insert(DijEntry(yiB,  TwoVertices(ii,-1)));
    } else {
      double kt2i = jet_scale_for_algorithm(_jets[ii]);
      int jj = DNN->NearestNeighbourIndex(ii);
      if (kt2i <= jet_scale_for_algorithm(_jets[jj])) {
        double dij = DeltaR2 * kt2i;
        DijMap.insert(DijEntry(dij, TwoVertices(ii,jj)));
      }
    }
  }
}

void 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 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 938 of file ClusterSequence.cc.

References _history, _jets, _writeout_combinations, ClusterSequence::history_element::child, ClusterSequence::history_element::dij, Invalid, ClusterSequence::history_element::jetp_index, ClusterSequence::history_element::max_dij_so_far, ClusterSequence::history_element::parent1, and ClusterSequence::history_element::parent2.

Referenced by _do_iB_recombination_step(), and _do_ij_recombination_step().

                                   {

  history_element element;
  element.parent1 = parent1;
  element.parent2 = parent2;
  element.jetp_index = jetp_index;
  element.child = Invalid;
  element.dij   = dij;
  element.max_dij_so_far = max(dij,_history[_history.size()-1].max_dij_so_far);
  _history.push_back(element);

  int local_step = _history.size()-1;
  assert(local_step == step_number);

  assert(parent1 >= 0);
  _history[parent1].child = local_step;
  if (parent2 >= 0) {_history[parent2].child = local_step;}

  // get cross-referencing right from PseudoJets
  if (jetp_index != Invalid) {
    assert(jetp_index >= 0);
    //cout << _jets.size() <<" "<<jetp_index<<"\n";
    _jets[jetp_index].set_cluster_hist_index(local_step);
  }

  if (_writeout_combinations) {
    cout << local_step << ": " 
         << parent1 << " with " << parent2
         << "; y = "<< dij<<endl;
  }

}

void 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().

template<class J >

double 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().

                                                                                   {
  double kt2 = jet->kt2;
  if (jet->NN != NULL) {if (jet->NN->kt2 < kt2) {kt2 = jet->NN->kt2;}}
  return jet->NN_dist * kt2;
}

template<>

double ClusterSequence::_bj_dist	(	const EEBriefJet *const	jeta,
		const EEBriefJet *const	jetb
	)			const [inline]

Definition at line 95 of file ClusterSequence_N2.cc.

References ClusterSequence::EEBriefJet::nx, ClusterSequence::EEBriefJet::ny, and ClusterSequence::EEBriefJet::nz.

                                                     {
  double dist = 1.0 
    - jeta->nx*jetb->nx
    - jeta->ny*jetb->ny
    - jeta->nz*jetb->nz;
  dist *= 2; // distance is _2_*min(Ei^2,Ej^2)*(1-cos theta)
  //cout << "Dist = " << dist << ": " 
  //     << jeta->nx << " "
  //     << jeta->ny << " "
  //     << jeta->nz << " "
  //     <<endl;
  return dist;
}

template<class J >

double 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 twopi.

Referenced by _bj_set_NN_crosscheck(), _bj_set_NN_nocross(), _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().

                                                                  {
  double dphi = std::abs(jetA->phi - jetB->phi);
  double deta = (jetA->eta - jetB->eta);
  if (dphi > pi) {dphi = twopi - dphi;}
  return dphi*dphi + deta*deta;
}

template<class J >

J* 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.

References _jets.

          {
    J * res;
    for(res = head; res<tail; res++) {
      if (_jets[res->_jets_index].cluster_hist_index() == hist_index) {break;}
    }
    return res;
  }

void ClusterSequence::_bj_remove_from_tiles

(

TiledJet *const

jet

)

[private]

Definition at line 50 of file ClusterSequence_TiledN2.cc.

References _tiles, ClusterSequence::Tile::head, ClusterSequence::TiledJet::next, ClusterSequence::TiledJet::previous, and ClusterSequence::TiledJet::tile_index.

                                                                {
  Tile * tile = & _tiles[jet->tile_index];

  if (jet->previous == NULL) {
    // we are at head of the tile, so reset it.
    // If this was the only jet on the tile then tile->head will now be NULL
    tile->head = jet->next;
  } else {
    // adjust link from previous jet in this tile
    jet->previous->next = jet->next;
  }
  if (jet->next != NULL) {
    // adjust backwards-link from next jet in this tile
    jet->next->previous = jet->previous;
  }
}

void 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().

template<>

void 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, ClusterSequence::EEBriefJet::_jets_index, _R2, _Rparam, ee_genkt_algorithm, ee_kt_algorithm, JetDefinition::extra_param(), jet_def(), ClusterSequence::EEBriefJet::kt2, ClusterSequence::EEBriefJet::NN, ClusterSequence::EEBriefJet::NN_dist, fastjet::norm(), ClusterSequence::EEBriefJet::nx, ClusterSequence::EEBriefJet::ny, and ClusterSequence::EEBriefJet::nz.

                                                                                 {

  double E = _jets[_jets_index].E();
  double scale = E*E; // the default energy scale for the kt alg
  double p  = jet_def().extra_param(); // in case we're ee_genkt
  switch (_jet_algorithm) {
  case ee_kt_algorithm:
    assert(_Rparam > 2.0); // force this to be true! [not best place, but works]
    // recall that _invR2 is artificially set to 1 for this alg
    // so that we automatically have dij = scale * 2(1-cos theta_ij)
    // Normally, _Rparam should be automatically set to 4 from JetDefinition
    break; 
  case ee_genkt_algorithm:
    if (p <= 0 && scale < 1e-300) scale = 1e-300; // same dodgy safety as genkt
    scale = pow(scale,p);
    break;
  default:
    throw Error("Unrecognised jet algorithm");
  }
  jetA->kt2  = scale; // "kt2" might one day be renamed as "scale" or some such

  double norm = _jets[_jets_index].modp2();
  if (norm > 0) {
    norm = 1.0/sqrt(norm);
    jetA->nx = norm * _jets[_jets_index].px();
    jetA->ny = norm * _jets[_jets_index].py();
    jetA->nz = norm * _jets[_jets_index].pz();
  } else {
    jetA->nx = 0.0;
    jetA->ny = 0.0;
    jetA->nz = 1.0;
  }
  jetA->_jets_index = _jets_index;
  // initialise NN info as well
  jetA->NN_dist = _R2;
  jetA->NN      = NULL;
}

template<class J >

void 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().

                                                                         {
    jetA->eta  = _jets[_jets_index].rap();
    jetA->phi  = _jets[_jets_index].phi_02pi();
    jetA->kt2  = jet_scale_for_algorithm(_jets[_jets_index]);
    jetA->_jets_index = _jets_index;
    // initialise NN info as well
    jetA->NN_dist = _R2;
    jetA->NN      = NULL;
}

template<class J >

void 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.

                                                                {
  double NN_dist = _R2;
  J * NN  = NULL;
  for (J * jetB = head; jetB != tail; jetB++) {
    double dist = _bj_dist(jet,jetB);
    if (dist < NN_dist) {
      NN_dist = dist;
      NN = jetB;
    }
    if (dist < jetB->NN_dist) {
      jetB->NN_dist = dist;
      jetB->NN = jet;
    }
  }
  jet->NN = NN;
  jet->NN_dist = NN_dist;
}

template<class J >

void 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.

                                                                            {
  double NN_dist = _R2;
  J * NN  = NULL;
  if (head < jet) {
    for (J * jetB = head; jetB != jet; jetB++) {
      double dist = _bj_dist(jet,jetB);
      if (dist < NN_dist) {
        NN_dist = dist;
        NN = jetB;
      }
    }
  }
  if (tail > jet) {
    for (J * jetB = jet+1; jetB != tail; jetB++) {
      double dist = _bj_dist(jet,jetB);
      if (dist < NN_dist) {
        NN_dist = dist;
        NN = jetB;
      }
    }
  }
  jet->NN = NN;
  jet->NN_dist = NN_dist;
}

void 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, cambridge_algorithm, Invalid, fastjet::d0::inline_maths::min(), and twopi.

Referenced by _initialise_and_run().

                                         {

  if (_jet_algorithm != cambridge_algorithm) throw Error("_CP2DChan_cluster called for a jet-finder that is not the cambridge algorithm");

  unsigned int n = _jets.size();

  vector<MirrorInfo>   coordIDs(2*n);  // link from original to mirror indices
  vector<int>          jetIDs(2*n);     // link from mirror to original indices
  vector<Coord2D>      coords(2*n);   // our coordinates (and copies)

  // start things off...
  double minrap = numeric_limits<double>::max();
  double maxrap = -minrap;
  int coord_index = 0;
  for (unsigned i = 0; i < n; i++) {
    // separate out points with infinite rapidity
    if (_jets[i].E() == abs(_jets[i].pz()) && _jets[i].perp2() == 0.0) {
      coordIDs[i] = MirrorInfo(BeamJet,BeamJet);
    } else {
      coordIDs[i].orig   = coord_index;
      coordIDs[i].mirror = coord_index+1;
      coords[coord_index]   = Coord2D(_jets[i].rap(), _jets[i].phi_02pi());
      coords[coord_index+1] = Coord2D(_jets[i].rap(), _jets[i].phi_02pi()+twopi);
      jetIDs[coord_index]   = i;
      jetIDs[coord_index+1] = i;
      minrap = min(coords[coord_index].x,minrap);
      maxrap = max(coords[coord_index].x,maxrap);
      coord_index += 2;
    }
  }
  // label remaining "mirror info" as saying that there are no jets
  for (unsigned i = n; i < 2*n; i++) {coordIDs[i].orig = Invalid;}

  // set them to their actual size...
  coords.resize(coord_index);

  // establish limits (with some leeway on rapidity)
  Coord2D left_edge(minrap-1.0, 0.0);
  Coord2D right_edge(maxrap+1.0, 2*twopi);

  // now create the closest pair search object
  ClosestPair2D cp(coords, left_edge, right_edge);

  // and start iterating...
  vector<Coord2D> new_points(2);
  vector<unsigned int> cIDs_to_remove(4);
  vector<unsigned int> new_cIDs(2);
  do {
    // find out which pair we recombine
    unsigned int cID1, cID2;
    double distance2;
    cp.closest_pair(cID1,cID2,distance2);
    distance2 *= _invR2;

    // if the closest distance > 1, we exit and go to "inclusive" stage
    if (distance2 > 1.0) {break;}

    // do the recombination...
    int jet_i = jetIDs[cID1];
    int jet_j = jetIDs[cID2];
    int newjet_k;
    _do_ij_recombination_step(jet_i, jet_j, distance2, newjet_k);

    // now prepare operations on CP structure
    cIDs_to_remove[0] = coordIDs[jet_i].orig;
    cIDs_to_remove[1] = coordIDs[jet_i].mirror;
    cIDs_to_remove[2] = coordIDs[jet_j].orig;
    cIDs_to_remove[3] = coordIDs[jet_j].mirror;
    new_points[0] = Coord2D(_jets[newjet_k].rap(),_jets[newjet_k].phi_02pi());
    new_points[1] = Coord2D(_jets[newjet_k].rap(),_jets[newjet_k].phi_02pi()+twopi);
    // carry out the CP operation
    //cp.replace_many(cIDs_to_remove, new_points, new_cIDs);
    // remarkable the following is faster...
    new_cIDs[0] = cp.replace(cIDs_to_remove[0], cIDs_to_remove[2], new_points[0]);
    new_cIDs[1] = cp.replace(cIDs_to_remove[1], cIDs_to_remove[3], new_points[1]);
    // signal that the following jets no longer exist
    coordIDs[jet_i].orig = Invalid;
    coordIDs[jet_j].orig = Invalid;
    // and do bookkeeping for new jet
    coordIDs[newjet_k] = MirrorInfo(new_cIDs[0], new_cIDs[1]);
    jetIDs[new_cIDs[0]] = newjet_k;
    jetIDs[new_cIDs[1]] = newjet_k;

    // if we've reached one jet we should exit...
    n--;
    if (n == 1) {break;}

  } while(true);
  

  // now do "beam" recombinations 
  //for (unsigned int jet_i = 0; jet_i < coordIDs.size(); jet_i++) {
  //  // if we have a predeclared beam jet or a valid set of mirror
  //  // coordinates, recombine then recombine this jet with the beam
  //  if (coordIDs[jet_i].orig == BeamJet || coordIDs[jet_i].orig > 0) {
  //    // diB = 1 _always_ (for the cambridge alg.)
  //    _do_iB_recombination_step(jet_i, 1.0);
  //  }
  //}

  _do_Cambridge_inclusive_jets();

  //n = _history.size();
  //for (unsigned int hist_i = 0; hist_i < n; hist_i++) {
  //  if (_history[hist_i].child == Invalid) {
  //    _do_iB_recombination_step(_history[hist_i].jetp_index, 1.0);
  //  }
  //}

}

void 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 cambridge_algorithm.

Referenced by _CP2DChan_cluster_2piMultD(), and _initialise_and_run().

                                               {

  if (_jet_algorithm != cambridge_algorithm) throw Error("CP2DChan clustering method called for a jet-finder that is not the cambridge algorithm");

  // run the clustering with mirror copies kept such that only things
  // within _Rparam of a border are mirrored
  _CP2DChan_limited_cluster(_Rparam);

  //
  _do_Cambridge_inclusive_jets();
}

void 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().

                                                  {

  // do a first run of clustering up to a small distance parameter,
  if (_Rparam >= 0.39) {
    _CP2DChan_limited_cluster(min(_Rparam/2,0.3));
  }

  // and then the final round of clustering
  _CP2DChan_cluster_2pi2R ();
}

void 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, Invalid, Private::make_mirror(), fastjet::d0::inline_maths::min(), and fastjet::pi.

Referenced by _CP2DChan_cluster_2pi2R(), and _CP2DChan_cluster_2piMultD().

                                                            {
  
  unsigned int n = _initial_n;

  vector<MirrorInfo>   coordIDs(2*n); // coord IDs of a given jetID
  vector<int>          jetIDs(2*n);   // jet ID for a given coord ID
  vector<Coord2D>      coords(2*n);   // our coordinates (and copies)

  // start things off...
  double minrap = numeric_limits<double>::max();
  double maxrap = -minrap;
  int coord_index = -1;
  int n_active = 0;
  for (unsigned jet_i = 0; jet_i < _jets.size(); jet_i++) {

    // skip jets that already have children or that have infinite
    // rapidity
    if (_history[_jets[jet_i].cluster_hist_index()].child != Invalid ||
        (_jets[jet_i].E() == abs(_jets[jet_i].pz()) && 
         _jets[jet_i].perp2() == 0.0)
        ) {continue;}

    n_active++;

    coordIDs[jet_i].orig = ++coord_index;
    coords[coord_index]  = Coord2D(_jets[jet_i].rap(), _jets[jet_i].phi_02pi());
    jetIDs[coord_index]  = jet_i;
    minrap = min(coords[coord_index].x,minrap);
    maxrap = max(coords[coord_index].x,maxrap);

    Coord2D mirror_point(coords[coord_index]);
    if (make_mirror(mirror_point, Dlim)) {
      coordIDs[jet_i].mirror = ++coord_index;
      coords[coord_index] = mirror_point;
      jetIDs[coord_index] = jet_i;
    } else {
      coordIDs[jet_i].mirror = Invalid;
    }
  }

  // set them to their actual size...
  coords.resize(coord_index+1);

  // establish limits (with some leeway on rapidity)
  Coord2D left_edge(minrap-1.0, -pi);
  Coord2D right_edge(maxrap+1.0, 3*pi);

  //cerr << "minrap, maxrap = " << minrap << " " << maxrap << endl;

  // now create the closest pair search object
  ClosestPair2D cp(coords, left_edge, right_edge);

  // cross check to see what's going on...
  //cerr << "Tree depths before: ";
  //cp.print_tree_depths(cerr);

  // and start iterating...
  vector<Coord2D> new_points(2);
  vector<unsigned int> cIDs_to_remove(4);
  vector<unsigned int> new_cIDs(2);

  do {
    // find out which pair we recombine
    unsigned int cID1, cID2;
    double distance2;
    cp.closest_pair(cID1,cID2,distance2);

    // if the closest distance > Dlim, we exit and go to "inclusive" stage
    if (distance2 > Dlim*Dlim) {break;}

    // normalise distance as we like it
    distance2 *= _invR2;

    // do the recombination...
    int jet_i = jetIDs[cID1];
    int jet_j = jetIDs[cID2];
    int newjet_k;
    _do_ij_recombination_step(jet_i, jet_j, distance2, newjet_k);

    // don't bother with any further action if only one active particle
    // is left (also avoid closest-pair error [cannot remove last particle]).
    if (--n_active == 1) {break;}

    // now prepare operations on CP structure
    cIDs_to_remove.resize(0);
    cIDs_to_remove.push_back(coordIDs[jet_i].orig);
    cIDs_to_remove.push_back(coordIDs[jet_j].orig);
    if (coordIDs[jet_i].mirror != Invalid) 
      cIDs_to_remove.push_back(coordIDs[jet_i].mirror);
    if (coordIDs[jet_j].mirror != Invalid) 
      cIDs_to_remove.push_back(coordIDs[jet_j].mirror);

    Coord2D new_point(_jets[newjet_k].rap(),_jets[newjet_k].phi_02pi());
    new_points.resize(0);
    new_points.push_back(new_point);
    if (make_mirror(new_point, Dlim)) new_points.push_back(new_point);
    
    // carry out actions on search tree
    cp.replace_many(cIDs_to_remove, new_points, new_cIDs);

    // now fill in info for new points...
    coordIDs[newjet_k].orig = new_cIDs[0];
    jetIDs[new_cIDs[0]]       = newjet_k;
    if (new_cIDs.size() == 2) {
      coordIDs[newjet_k].mirror = new_cIDs[1];
      jetIDs[new_cIDs[1]]         = newjet_k;
    } else {coordIDs[newjet_k].mirror = Invalid;}
    
    //n_active--;
    //if (n_active == 1) {break;}

  } while(true);
  
  // cross check to see what's going on...
  //cerr << "Max tree depths after: ";
  //cp.print_tree_depths(cerr);

}

void 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 228 of file ClusterSequence.cc.

References _invR2, _jet_algorithm, _jet_def, _plugin_activated, _print_banner(), _R2, _Rparam, _strategy, _writeout_combinations, JetDefinition::jet_algorithm(), JetDefinition::R(), and JetDefinition::strategy().

Referenced by ClusterSequenceActiveArea::_initialise_AA(), and _initialise_and_run().

                                                                          {

  // let the user know what's going on
  _print_banner();

  // make a local copy of the jet definition (for future use?)
  _jet_def = jet_def;
  
  _writeout_combinations = writeout_combinations;
  _jet_algorithm = jet_def.jet_algorithm();
  _Rparam = jet_def.R();  _R2 = _Rparam*_Rparam; _invR2 = 1.0/_R2;
  _strategy = jet_def.strategy();

  // disallow interference from the plugin
  _plugin_activated = false;
  
}

void 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, NlnN, NlnN3pi, NlnN4pi, strategy_string(), and twopi.

Referenced by _initialise_and_run().

                                         {

  int n = _jets.size();

  vector<EtaPhi> points(n); // recall EtaPhi is just a typedef'd pair<double>
  for (int i = 0; i < n; i++) {
    points[i] = EtaPhi(_jets[i].rap(),_jets[i].phi_02pi());
    points[i].sanitize(); // make sure things are in the right range
  }

  // initialise our DNN structure with the set of points
  DynamicNearestNeighbours * DNN;
#ifndef DROP_CGAL // strategy = NlnN* are not supported if we drop CGAL...
  bool verbose = false;
  bool ignore_nearest_is_mirror = (_Rparam < twopi);
  if (_strategy == NlnN4pi) {
    DNN = new Dnn4piCylinder(points,verbose);
  } else if (_strategy == NlnN3pi) {
    DNN = new Dnn3piCylinder(points,ignore_nearest_is_mirror,verbose);
  } else if (_strategy == NlnN) {
    DNN = new Dnn2piCylinder(points,ignore_nearest_is_mirror,verbose);
  } else 
#else
  if (_strategy == NlnN4pi || _strategy == NlnN3pi || _strategy == NlnN) {
    ostringstream err;
    err << "ERROR: Requested strategy "<<strategy_string()<<" but it is not"<<endl;
    err << "       supported because FastJet was compiled without CGAL"<<endl;
    throw Error(err.str());
    //assert(false);
  }
#endif // DROP_CGAL
  {
    ostringstream err;
    err << "ERROR: Unrecognized value for strategy: "<<_strategy<<endl;
    assert(false);
    throw Error(err.str());
  }

  // We will find nearest neighbour for each vertex, and include
  // distance in map (NB DistMap is a typedef given in the .h file)
  DistMap DijMap;

  // fill the map with the minimal (as far as we know) subset of Dij
  // distances (i.e. nearest neighbour ones).
  for (int ii = 0; ii < n; ii++) {
    _add_ktdistance_to_map(ii, DijMap, DNN);
  }

  // run the clustering (go up to i=n-1, but then will stop half-way down,
  // when we reach that point -- it will be the final beam jet and there
  // will be no nearest neighbours to find).
  for (int i=0;i<n;i++) {
    // find nearest vertices
    // NB: skip cases where the point is not there anymore!
    TwoVertices SmallestDijPair;
    int jet_i, jet_j;
    double SmallestDij;
    bool Valid2;
    bool recombine_with_beam;
    do { 
      SmallestDij = DijMap.begin()->first;
      SmallestDijPair = DijMap.begin()->second;
      jet_i = SmallestDijPair.first;
      jet_j = SmallestDijPair.second;
      // distance is immediately removed regardless of whether or not
      // it is used.
      // Some temporary testing code relating to problems with the gcc-3.2 compiler
      //cout << "got here and size is "<< DijMap.size()<< " and it is "<<SmallestDij <<"\n";
      //cout <<  jet_i << " "<< jet_j<<"\n";
      DijMap.erase(DijMap.begin());
      //cout << "got beyond here\n";

      // need to "prime" the validity of jet_j in such a way that 
      // if it corresponds to the beam then it is automatically valid.
      recombine_with_beam = (jet_j == BeamJet);
      if (!recombine_with_beam) {Valid2 = DNN->Valid(jet_j);} 
      else {Valid2 = true;}

    } while ( !DNN->Valid(jet_i) || !Valid2);


    // The following part acts just on jet momenta and on the history.
    // The action on the nearest-neighbour structures takes place
    // later (only if at least 2 jets are around).
    if (! recombine_with_beam) {
      int nn; // will be index of new jet
      _do_ij_recombination_step(jet_i, jet_j, SmallestDij, nn);
      //OBS // merge the two jets, add new jet, remove old ones
      //OBS _jets.push_back(_jets[jet_i] + _jets[jet_j]);
      //OBS 
      //OBS int nn = _jets.size()-1;
      //OBS _jets[nn].set_cluster_hist_index(n+i);
      //OBS 
      //OBS // get corresponding indices in history structure
      //OBS int hist_i = _jets[jet_i].cluster_hist_index();
      //OBS int hist_j = _jets[jet_j].cluster_hist_index();
      //OBS 
      //OBS 
      //OBS _add_step_to_history(n+i,min(hist_i,hist_j), max(hist_i,hist_j),
      //OBS                   _jets.size()-1, SmallestDij);

      // add new point to points vector
      EtaPhi newpoint(_jets[nn].rap(), _jets[nn].phi_02pi());
      newpoint.sanitize(); // make sure it is in correct range
      points.push_back(newpoint);
    } else {
      // recombine the jet with the beam
      _do_iB_recombination_step(jet_i, SmallestDij);
      //OBS _add_step_to_history(n+i,_jets[jet_i].cluster_hist_index(),BeamJet,
      //OBS                        Invalid, SmallestDij);
    }

    // exit the loop because we do not want to look for nearest neighbours
    // etc. of zero partons
    if (i == n-1) {break;}

    vector<int> updated_neighbours;
    if (! recombine_with_beam) {
      // update DNN
      int point3;
      DNN->RemoveCombinedAddCombination(jet_i, jet_j, 
                                       points[points.size()-1], point3,
                                       updated_neighbours);
      // C++ beginners' comment: static_cast to unsigned int is necessary
      // to do away with warnings about type mismatch between point3 (int) 
      // and points.size (unsigned int)
      if (static_cast<unsigned int> (point3) != points.size()-1) {
        throw Error("INTERNAL ERROR: point3 != points.size()-1");}
    } else {
      // update DNN
      DNN->RemovePoint(jet_i, updated_neighbours);
    }

    // update map
    vector<int>::iterator it = updated_neighbours.begin();
    for (; it != updated_neighbours.end(); ++it) {
      int ii = *it;
      _add_ktdistance_to_map(ii, DijMap, DNN);
    }
      
  } // end clustering loop 
  
  // remember to clean up!
  delete DNN;
}

void 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().

                                                    {
  unsigned int n = _history.size();
  for (unsigned int hist_i = 0; hist_i < n; hist_i++) {
    if (_history[hist_i].child == Invalid) {
      _do_iB_recombination_step(_history[hist_i].jetp_index, 1.0);
    }
  }
}

void 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 1115 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(), ClusterSequenceActiveArea::_transfer_ghost_free_history(), and plugin_record_iB_recombination().

                                                                         {
  // get history index
  int newstep_k = _history.size();

  // recombine the jet with the beam
  _add_step_to_history(newstep_k,_jets[jet_i].cluster_hist_index(),BeamJet,
                       Invalid, diB);

}

void 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 1082 of file ClusterSequence.cc.

References _add_step_to_history(), _history, _jet_def, _jets, fastjet::d0::inline_maths::min(), and 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(), ClusterSequenceActiveArea::_transfer_ghost_free_history(), and plugin_record_ij_recombination().

                                               {

  // create the new jet by recombining the first two
  PseudoJet newjet;
  _jet_def.recombiner()->recombine(_jets[jet_i], _jets[jet_j], newjet);
  _jets.push_back(newjet);
  // original version...
  //_jets.push_back(_jets[jet_i] + _jets[jet_j]);

  // get its index
  newjet_k = _jets.size()-1;

  // get history index
  int newstep_k = _history.size();
  // and provide jet with the info
  _jets[newjet_k].set_cluster_hist_index(newstep_k);

  // finally sort out the history 
  int hist_i = _jets[jet_i].cluster_hist_index();
  int hist_j = _jets[jet_j].cluster_hist_index();

  _add_step_to_history(newstep_k, min(hist_i, hist_j), max(hist_i,hist_j),
                       newjet_k, dij);

}

void 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 ClusterSequenceActiveArea.

Referenced by unique_history_order().

void 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 ClusterSequenceActiveArea.

void 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, ClusterSequence::TiledJet::_jets_index, _R2, _tiles, _tj_set_jetinfo(), ClusterSequence::Tile::begin_tiles, ClusterSequence::TiledJet::diJ_posn, ClusterSequence::Tile::end_tiles, ClusterSequence::Tile::head, n_tile_neighbours, ClusterSequence::TiledJet::next, ClusterSequence::TiledJet::NN, ClusterSequence::TiledJet::NN_dist, ClusterSequence::Tile::tagged, and ClusterSequence::TiledJet::tile_index.

Referenced by _initialise_and_run().

                                               {

  _initialise_tiles();

  int n = _jets.size();
  TiledJet * briefjets = new TiledJet[n];
  TiledJet * jetA = briefjets, * jetB;
  TiledJet oldB;
  

  // will be used quite deep inside loops, but declare it here so that
  // memory (de)allocation gets done only once
  vector<int> tile_union(3*n_tile_neighbours);
  
  // initialise the basic jet info 
  for (int i = 0; i< n; i++) {
    _tj_set_jetinfo(jetA, i);
    //cout << i<<": "<<jetA->tile_index<<"\n";
    jetA++; // move on to next entry of briefjets
  }
  TiledJet * head = briefjets; // a nicer way of naming start

  // set up the initial nearest neighbour information
  vector<Tile>::const_iterator tile;
  for (tile = _tiles.begin(); tile != _tiles.end(); tile++) {
    // first do it on this tile
    for (jetA = tile->head; jetA != NULL; jetA = jetA->next) {
      for (jetB = tile->head; jetB != jetA; jetB = jetB->next) {
        double dist = _bj_dist(jetA,jetB);
        if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;}
        if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;}
      }
    }
    // then do it for RH tiles
    for (Tile ** RTile = tile->RH_tiles; RTile != tile->end_tiles; RTile++) {
      for (jetA = tile->head; jetA != NULL; jetA = jetA->next) {
        for (jetB = (*RTile)->head; jetB != NULL; jetB = jetB->next) {
          double dist = _bj_dist(jetA,jetB);
          if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;}
          if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;}
        }
      }
    }
    // no need to do it for LH tiles, since they are implicitly done
    // when we set NN for both jetA and jetB on the RH tiles.
  }

  
  // now create the diJ (where J is i's NN) table -- remember that 
  // we differ from standard normalisation here by a factor of R2
  // (corrected for at the end). 
  struct diJ_plus_link {
    double     diJ; // the distance
    TiledJet * jet; // the jet (i) for which we've found this distance
                    // (whose NN will the J).
  };
  diJ_plus_link * diJ = new diJ_plus_link[n];
  jetA = head;
  for (int i = 0; i < n; i++) {
    diJ[i].diJ = _bj_diJ(jetA); // kt distance * R^2
    diJ[i].jet = jetA;  // our compact diJ table will not be in      
    jetA->diJ_posn = i; // one-to-one corresp. with non-compact jets,
                        // so set up bi-directional correspondence here.
    jetA++; // have jetA follow i 
  }

  // now run the recombination loop
  int history_location = n-1;
  while (n > 0) {

    // find the minimum of the diJ on this round
    diJ_plus_link * best, *stop; // pointers a bit faster than indices
    // could use best to keep track of diJ min, but it turns out to be
    // marginally faster to have a separate variable (avoids n
    // dereferences at the expense of n/2 assignments).
    double diJ_min = diJ[0].diJ; // initialise the best one here.
    best = diJ;                  // and here
    stop = diJ+n;
    for (diJ_plus_link * here = diJ+1; here != stop; here++) {
      if (here->diJ < diJ_min) {best = here; diJ_min  = here->diJ;}
    }

    // do the recombination between A and B
    history_location++;
    jetA = best->jet;
    jetB = jetA->NN;
    // put the normalisation back in
    diJ_min *= _invR2; 

    if (jetB != NULL) {
      // jet-jet recombination
      // If necessary relabel A & B to ensure jetB < jetA, that way if
      // the larger of them == newtail then that ends up being jetA and 
      // the new jet that is added as jetB is inserted in a position that
      // has a future!
      if (jetA < jetB) {swap(jetA,jetB);}

      int nn; // new jet index
      _do_ij_recombination_step(jetA->_jets_index, jetB->_jets_index, diJ_min, nn);
      
      //OBS// get the two history indices
      //OBSint ihstry_a = _jets[jetA->_jets_index].cluster_hist_index();
      //OBSint ihstry_b = _jets[jetB->_jets_index].cluster_hist_index();
      //OBS// create the recombined jet
      //OBS_jets.push_back(_jets[jetA->_jets_index] + _jets[jetB->_jets_index]);
      //OBSint nn = _jets.size() - 1;
      //OBS_jets[nn].set_cluster_hist_index(history_location);
      //OBS// update history
      //OBS//cout <<n-1<<" "<<jetA-head<<" "<<jetB-head<<"; ";
      //OBS_add_step_to_history(history_location, 
      //OBS                min(ihstry_a,ihstry_b),max(ihstry_a,ihstry_b),
      //OBS                        nn, diJ_min);
      // what was jetB will now become the new jet
      _bj_remove_from_tiles(jetA);
      oldB = * jetB;  // take a copy because we will need it...
      _bj_remove_from_tiles(jetB);
      _tj_set_jetinfo(jetB, nn); // cause jetB to become _jets[nn]
                                 // (also registers the jet in the tiling)
    } else {
      // jet-beam recombination
      // get the hist_index
      _do_iB_recombination_step(jetA->_jets_index, diJ_min);
      //OBSint ihstry_a = _jets[jetA->_jets_index].cluster_hist_index();
      //OBS//cout <<n-1<<" "<<jetA-head<<" "<<-1<<"; ";
      //OBS_add_step_to_history(history_location,ihstry_a,BeamJet,Invalid,diJ_min); 
      _bj_remove_from_tiles(jetA);
    }

    // first establish the set of tiles over which we are going to
    // have to run searches for updated and new nearest-neighbours --
    // basically a combination of vicinity of the tiles of the two old
    // and one new jet.
    int n_near_tiles = 0;
    _add_untagged_neighbours_to_tile_union(jetA->tile_index, 
                                           tile_union, n_near_tiles);
    if (jetB != NULL) {
      if (jetB->tile_index != jetA->tile_index) {
        _add_untagged_neighbours_to_tile_union(jetB->tile_index,
                                               tile_union,n_near_tiles);
      }
      if (oldB.tile_index != jetA->tile_index && 
          oldB.tile_index != jetB->tile_index) {
        _add_untagged_neighbours_to_tile_union(oldB.tile_index,
                                               tile_union,n_near_tiles);
      }
    }

    // now update our nearest neighbour info and diJ table
    // first reduce size of table
    n--;
    // then compactify the diJ by taking the last of the diJ and copying
    // it to the position occupied by the diJ for jetA
    diJ[n].jet->diJ_posn = jetA->diJ_posn;
    diJ[jetA->diJ_posn] = diJ[n];

    // Initialise jetB's NN distance as well as updating it for 
    // other particles.
    // Run over all tiles in our union 
    for (int itile = 0; itile < n_near_tiles; itile++) {
      Tile * tile = &_tiles[tile_union[itile]];
      tile->tagged = false; // reset tag, since we're done with unions
      // run over all jets in the current tile
      for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) {
        // see if jetI had jetA or jetB as a NN -- if so recalculate the NN
        if (jetI->NN == jetA || (jetI->NN == jetB && jetB != NULL)) {
          jetI->NN_dist = _R2;
          jetI->NN      = NULL;
          // now go over tiles that are neighbours of I (include own tile)
          for (Tile ** near_tile  = tile->begin_tiles; 
                       near_tile != tile->end_tiles; near_tile++) {
            // and then over the contents of that tile
            for (TiledJet * jetJ  = (*near_tile)->head; 
                            jetJ != NULL; jetJ = jetJ->next) {
              double dist = _bj_dist(jetI,jetJ);
              if (dist < jetI->NN_dist && jetJ != jetI) {
                jetI->NN_dist = dist; jetI->NN = jetJ;
              }
            }
          }
          diJ[jetI->diJ_posn].diJ = _bj_diJ(jetI); // update diJ kt-dist
        }
        // check whether new jetB is closer than jetI's current NN and
        // if jetI is closer than jetB's current (evolving) nearest
        // neighbour. Where relevant update things
        if (jetB != NULL) {
          double dist = _bj_dist(jetI,jetB);
          if (dist < jetI->NN_dist) {
            if (jetI != jetB) {
              jetI->NN_dist = dist;
              jetI->NN = jetB;
              diJ[jetI->diJ_posn].diJ = _bj_diJ(jetI); // update diJ...
            }
          }
          if (dist < jetB->NN_dist) {
            if (jetI != jetB) {
              jetB->NN_dist = dist;
              jetB->NN      = jetI;}
          }
        }
      }
    }

    // finally, register the updated kt distance for B
    if (jetB != NULL) {diJ[jetB->diJ_posn].diJ = _bj_diJ(jetB);}

  }

  // final cleaning up;
  delete[] diJ;
  delete[] briefjets;
}

void 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 250 of file ClusterSequence.cc.

References _history, _initial_n, _jet_def, _jets, _Qtot, ClusterSequence::history_element::child, ClusterSequence::history_element::dij, InexistentParent, Invalid, ClusterSequence::history_element::jetp_index, ClusterSequence::history_element::max_dij_so_far, ClusterSequence::history_element::parent1, ClusterSequence::history_element::parent2, and JetDefinition::recombiner().

Referenced by ClusterSequenceActiveArea::_initialise_AA(), and _initialise_and_run().

                                             {

  //if (_jets.size() == 0) {throw Error("Cannot run jet-finder on empty event");}

  // reserve sufficient space for everything
  _jets.reserve(_jets.size()*2);
  _history.reserve(_jets.size()*2);

  _Qtot = 0;

  for (int i = 0; i < static_cast<int>(_jets.size()) ; i++) {
    history_element element;
    element.parent1 = InexistentParent;
    element.parent2 = InexistentParent;
    element.child   = Invalid;
    element.jetp_index = i;
    element.dij     = 0.0;
    element.max_dij_so_far = 0.0;

    _history.push_back(element);
    
    // do any momentum preprocessing needed by the recombination scheme
    _jet_def.recombiner()->preprocess(_jets[i]);

    // get cross-referencing right from PseudoJets
    _jets[i].set_cluster_hist_index(i);

    // determine the total energy in the event
    _Qtot += _jets[i].E();
  }
  _initial_n = _jets.size();
}

void 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().

                                                                      {

  JetDefinition jet_def(_default_jet_algorithm, R, strategy);
  _initialise_and_run(jet_def, writeout_combinations);
}

void 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(), antikt_algorithm, Best, cambridge_algorithm, ee_genkt_algorithm, ee_kt_algorithm, JetDefinition::jet_algorithm(), fastjet::d0::inline_maths::min(), N2MinHeapTiled, N2Plain, N2PoorTiled, N2Tiled, N3Dumb, n_particles(), NlnN, NlnN3pi, NlnN4pi, NlnNCam, NlnNCam2pi2R, NlnNCam4pi, fastjet::pi, JetDefinition::plugin(), and plugin_algorithm.

Referenced by ClusterSequenceActiveArea::_initialise_AA(), _initialise_and_run(), and ClusterSequence().

                                                                      {

  // transfer all relevant info into internal variables
  _decant_options(jet_def, writeout_combinations);

  // set up the history entries for the initial particles (those
  // currently in _jets)
  _fill_initial_history();

  // don't run anything if the event is empty
  if (n_particles() == 0) return;

  // ----- deal with special cases: plugins & e+e- ------
  if (_jet_algorithm == plugin_algorithm) {
    // allows plugin_xyz() functions to modify cluster sequence
    _plugin_activated = true;
    // let the plugin do its work here
    _jet_def.plugin()->run_clustering( (*this) );
    _plugin_activated = false;
    return;
  } else if (_jet_algorithm == ee_kt_algorithm ||
             _jet_algorithm == ee_genkt_algorithm) {
    // ignore requested strategy
    _strategy = N2Plain;
    if (_jet_algorithm == ee_kt_algorithm) {
      // make sure that R is large enough so that "beam" recomb only
      // occurs when a single particle is left
      // Normally, this should be automatically set to 4 from JetDefinition
      assert(_Rparam > 2.0); 
      // this is used to renormalise the dij to get a "standard" form
      // and our convention in e+e- will be different from that
      // in long.inv case; NB: _invR2 name should be changed -> _renorm_dij?
      _invR2 = 1.0;
    } else {
      // as of 2009-01-09, choose R to be an angular distance, in
      // radians.  Since the algorithm uses 2(1-cos(theta)) as its
      // squared angular measure, make sure that the _R2 is defined
      // in a similar way.
      if (_Rparam > pi) {
        // choose a value that ensures that back-to-back particles will
        // always recombine 
        //_R2 = 4.0000000000001;
        _R2 = 2 * ( 3.0 + cos(_Rparam) );
      } else {
        _R2    = 2 * ( 1.0 - cos(_Rparam) );
      }
      _invR2 = 1.0/_R2;
    }
    //_simple_N2_cluster<EEBriefJet>();
    _simple_N2_cluster_EEBriefJet();
    return;
  }


  // automatically redefine the strategy according to N if that is
  // what the user requested -- transition points (and especially
  // their R-dependence) are based on empirical observations for a
  // R=0.4, 0.7 and 1.0, running on toth (3.4GHz, Pentium IV D [dual
  // core] with 2MB of cache).
  if (_strategy == Best) {
    int N = _jets.size();
    if (N > 6200/pow(_Rparam,2.0) 
        && jet_def.jet_algorithm() == cambridge_algorithm) {
      _strategy = NlnNCam;}
    else
#ifndef DROP_CGAL
      if ((N > 16000/pow(_Rparam,1.15) && jet_def.jet_algorithm() != antikt_algorithm)
        || N > 35000/pow(_Rparam,1.15)) {
      _strategy = NlnN; }   
    else                    
#endif  // DROP_CGAL
      if (N > 450) {
      _strategy = N2MinHeapTiled;
    }
    else if (N > 55*max(0.5,min(1.0,_Rparam))) {// empirical scaling with R
      _strategy = N2Tiled;
    } else {
      _strategy = N2Plain;
    }
  }


  // run the code containing the selected strategy
  if (_strategy == NlnN || _strategy == NlnN3pi 
      || _strategy == NlnN4pi ) {
    this->_delaunay_cluster();
  } else if (_strategy ==  N3Dumb ) {
    this->_really_dumb_cluster();
  } else if (_strategy == N2Tiled) {
    this->_faster_tiled_N2_cluster();
  } else if (_strategy == N2PoorTiled) {
    this->_tiled_N2_cluster();
  } else if (_strategy == N2Plain) {
    // BriefJet provides standard long.invariant kt alg.
    //this->_simple_N2_cluster<BriefJet>();
    this->_simple_N2_cluster_BriefJet();
  } else if (_strategy == N2MinHeapTiled) {
    this->_minheap_faster_tiled_N2_cluster();
  } else if (_strategy == NlnNCam4pi) {
    this->_CP2DChan_cluster();
  } else if (_strategy == NlnNCam2pi2R) {
    this->_CP2DChan_cluster_2pi2R();
  } else if (_strategy == NlnNCam) {
    this->_CP2DChan_cluster_2piMultD();
  } else {
    ostringstream err;
    err << "Unrecognised value for strategy: "<<_strategy;
    throw Error(err.str());
    //assert(false);
  }
}

void 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, ClusterSequence::Tile::begin_tiles, ClusterSequence::Tile::end_tiles, ClusterSequence::Tile::head, ClusterSequence::Tile::RH_tiles, ClusterSequence::Tile::surrounding_tiles, ClusterSequence::Tile::tagged, and twopi.

Referenced by _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().

                                        {

  // first decide tile sizes (with a lower bound to avoid huge memory use with
  // very small R)
  double default_size = max(0.1,_Rparam);
  _tile_size_eta = default_size;
  _n_tiles_phi   = int(floor(twopi/default_size));
  _tile_size_phi = twopi / _n_tiles_phi; // >= _Rparam and fits in 2pi

  // always include zero rapidity in the tiling region
  _tiles_eta_min = 0.0;
  _tiles_eta_max = 0.0;
  // but go no further than following
  const double maxrap = 7.0;

  // and find out how much further one should go
  for(unsigned int i = 0; i < _jets.size(); i++) {
    double eta = _jets[i].rap();
    // first check if eta is in range -- to avoid taking into account
    // very spurious rapidities due to particles with near-zero kt.
    if (abs(eta) < maxrap) {
      if (eta < _tiles_eta_min) {_tiles_eta_min = eta;}
      if (eta > _tiles_eta_max) {_tiles_eta_max = eta;}
    }
  }

  // now adjust the values
  _tiles_ieta_min = int(floor(_tiles_eta_min/_tile_size_eta));
  _tiles_ieta_max = int(floor( _tiles_eta_max/_tile_size_eta));
  _tiles_eta_min = _tiles_ieta_min * _tile_size_eta;
  _tiles_eta_max = _tiles_ieta_max * _tile_size_eta;

  // allocate the tiles
  _tiles.resize((_tiles_ieta_max-_tiles_ieta_min+1)*_n_tiles_phi);

  // now set up the cross-referencing between tiles
  for (int ieta = _tiles_ieta_min; ieta <= _tiles_ieta_max; ieta++) {
    for (int iphi = 0; iphi < _n_tiles_phi; iphi++) {
      Tile * tile = & _tiles[_tile_index(ieta,iphi)];
      // no jets in this tile yet
      tile->head = NULL; // first element of tiles points to itself
      tile->begin_tiles[0] =  tile;
      Tile ** pptile = & (tile->begin_tiles[0]);
      pptile++;
      //
      // set up L's in column to the left of X
      tile->surrounding_tiles = pptile;
      if (ieta > _tiles_ieta_min) {
        // with the itile subroutine, we can safely run tiles from
        // idphi=-1 to idphi=+1, because it takes care of
        // negative and positive boundaries
        for (int idphi = -1; idphi <=+1; idphi++) {
          *pptile = & _tiles[_tile_index(ieta-1,iphi+idphi)];
          pptile++;
        }       
      }
      // now set up last L (below X)
      *pptile = & _tiles[_tile_index(ieta,iphi-1)];
      pptile++;
      // set up first R (above X)
      tile->RH_tiles = pptile;
      *pptile = & _tiles[_tile_index(ieta,iphi+1)];
      pptile++;
      // set up remaining R's, to the right of X
      if (ieta < _tiles_ieta_max) {
        for (int idphi = -1; idphi <= +1; idphi++) {
          *pptile = & _tiles[_tile_index(ieta+1,iphi+idphi)];
          pptile++;
        }       
      }
      // now put semaphore for end tile
      tile->end_tiles = pptile;
      // finally make sure tiles are untagged
      tile->tagged = false;
    }
  }

}

void 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, ClusterSequence::TiledJet::_jets_index, _R2, _tiles, _tj_set_jetinfo(), ClusterSequence::Tile::begin_tiles, ClusterSequence::Tile::end_tiles, ClusterSequence::Tile::head, ClusterSequence::TiledJet::label_minheap_update_done(), n_tile_neighbours, ClusterSequence::TiledJet::next, ClusterSequence::TiledJet::NN, ClusterSequence::TiledJet::NN_dist, ClusterSequence::Tile::tagged, and ClusterSequence::TiledJet::tile_index.

Referenced by _initialise_and_run().

                                                       {

  _initialise_tiles();

  int n = _jets.size();
  TiledJet * briefjets = new TiledJet[n];
  TiledJet * jetA = briefjets, * jetB;
  TiledJet oldB;
  

  // will be used quite deep inside loops, but declare it here so that
  // memory (de)allocation gets done only once
  vector<int> tile_union(3*n_tile_neighbours);
  
  // initialise the basic jet info 
  for (int i = 0; i< n; i++) {
    _tj_set_jetinfo(jetA, i);
    //cout << i<<": "<<jetA->tile_index<<"\n";
    jetA++; // move on to next entry of briefjets
  }
  TiledJet * head = briefjets; // a nicer way of naming start

  // set up the initial nearest neighbour information
  vector<Tile>::const_iterator tile;
  for (tile = _tiles.begin(); tile != _tiles.end(); tile++) {
    // first do it on this tile
    for (jetA = tile->head; jetA != NULL; jetA = jetA->next) {
      for (jetB = tile->head; jetB != jetA; jetB = jetB->next) {
        double dist = _bj_dist(jetA,jetB);
        if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;}
        if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;}
      }
    }
    // then do it for RH tiles
    for (Tile ** RTile = tile->RH_tiles; RTile != tile->end_tiles; RTile++) {
      for (jetA = tile->head; jetA != NULL; jetA = jetA->next) {
        for (jetB = (*RTile)->head; jetB != NULL; jetB = jetB->next) {
          double dist = _bj_dist(jetA,jetB);
          if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;}
          if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;}
        }
      }
    }
    // no need to do it for LH tiles, since they are implicitly done
    // when we set NN for both jetA and jetB on the RH tiles.
  }

  
  //struct diJ_plus_link {
  //  double     diJ; // the distance
  //  TiledJet * jet; // the jet (i) for which we've found this distance
  //                  // (whose NN will the J).
  //};
  //diJ_plus_link * diJ = new diJ_plus_link[n];
  //jetA = head;
  //for (int i = 0; i < n; i++) {
  //  diJ[i].diJ = _bj_diJ(jetA); // kt distance * R^2
  //  diJ[i].jet = jetA;  // our compact diJ table will not be in            
  //  jetA->diJ_posn = i; // one-to-one corresp. with non-compact jets,
  //                      // so set up bi-directional correspondence here.
  //  jetA++; // have jetA follow i 
  //}

  vector<double> diJs(n);
  for (int i = 0; i < n; i++) {
    diJs[i] = _bj_diJ(&briefjets[i]);
    briefjets[i].label_minheap_update_done();
  }
  MinHeap minheap(diJs);
  // have a stack telling us which jets we'll have to update on the heap
  vector<TiledJet *> jets_for_minheap;
  jets_for_minheap.reserve(n); 

  // now run the recombination loop
  int history_location = n-1;
  while (n > 0) {

    double diJ_min = minheap.minval() *_invR2;
    jetA = head + minheap.minloc();

    // do the recombination between A and B
    history_location++;
    jetB = jetA->NN;

    if (jetB != NULL) {
      // jet-jet recombination
      // If necessary relabel A & B to ensure jetB < jetA, that way if
      // the larger of them == newtail then that ends up being jetA and 
      // the new jet that is added as jetB is inserted in a position that
      // has a future!
      if (jetA < jetB) {swap(jetA,jetB);}

      int nn; // new jet index
      _do_ij_recombination_step(jetA->_jets_index, jetB->_jets_index, diJ_min, nn);
      
      // what was jetB will now become the new jet
      _bj_remove_from_tiles(jetA);
      oldB = * jetB;  // take a copy because we will need it...
      _bj_remove_from_tiles(jetB);
      _tj_set_jetinfo(jetB, nn); // cause jetB to become _jets[nn]
                                 // (also registers the jet in the tiling)
    } else {
      // jet-beam recombination
      // get the hist_index
      _do_iB_recombination_step(jetA->_jets_index, diJ_min);
      _bj_remove_from_tiles(jetA);
    }

    // remove the minheap entry for jetA
    minheap.remove(jetA-head);

    // first establish the set of tiles over which we are going to
    // have to run searches for updated and new nearest-neighbours --
    // basically a combination of vicinity of the tiles of the two old
    // and one new jet.
    int n_near_tiles = 0;
    _add_untagged_neighbours_to_tile_union(jetA->tile_index, 
                                           tile_union, n_near_tiles);
    if (jetB != NULL) {
      if (jetB->tile_index != jetA->tile_index) {
        _add_untagged_neighbours_to_tile_union(jetB->tile_index,
                                               tile_union,n_near_tiles);
      }
      if (oldB.tile_index != jetA->tile_index && 
          oldB.tile_index != jetB->tile_index) {
        _add_untagged_neighbours_to_tile_union(oldB.tile_index,
                                               tile_union,n_near_tiles);
      }
      // indicate that we'll have to update jetB in the minheap
      jetB->label_minheap_update_needed();
      jets_for_minheap.push_back(jetB);
    }


    // Initialise jetB's NN distance as well as updating it for 
    // other particles.
    // Run over all tiles in our union 
    for (int itile = 0; itile < n_near_tiles; itile++) {
      Tile * tile = &_tiles[tile_union[itile]];
      tile->tagged = false; // reset tag, since we're done with unions
      // run over all jets in the current tile
      for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) {
        // see if jetI had jetA or jetB as a NN -- if so recalculate the NN
        if (jetI->NN == jetA || (jetI->NN == jetB && jetB != NULL)) {
          jetI->NN_dist = _R2;
          jetI->NN      = NULL;
          // label jetI as needing heap action...
          if (!jetI->minheap_update_needed()) {
            jetI->label_minheap_update_needed();
            jets_for_minheap.push_back(jetI);}
          // now go over tiles that are neighbours of I (include own tile)
          for (Tile ** near_tile  = tile->begin_tiles; 
                       near_tile != tile->end_tiles; near_tile++) {
            // and then over the contents of that tile
            for (TiledJet * jetJ  = (*near_tile)->head; 
                            jetJ != NULL; jetJ = jetJ->next) {
              double dist = _bj_dist(jetI,jetJ);
              if (dist < jetI->NN_dist && jetJ != jetI) {
                jetI->NN_dist = dist; jetI->NN = jetJ;
              }
            }
          }
        }
        // check whether new jetB is closer than jetI's current NN and
        // if jetI is closer than jetB's current (evolving) nearest
        // neighbour. Where relevant update things
        if (jetB != NULL) {
          double dist = _bj_dist(jetI,jetB);
          if (dist < jetI->NN_dist) {
            if (jetI != jetB) {
              jetI->NN_dist = dist;
              jetI->NN = jetB;
              // label jetI as needing heap action...
              if (!jetI->minheap_update_needed()) {
                jetI->label_minheap_update_needed();
                jets_for_minheap.push_back(jetI);}
            }
          }
          if (dist < jetB->NN_dist) {
            if (jetI != jetB) {
              jetB->NN_dist = dist;
              jetB->NN      = jetI;}
          }
        }
      }
    }

    // deal with jets whose minheap entry needs updating
    while (jets_for_minheap.size() > 0) {
      TiledJet * jetI = jets_for_minheap.back(); 
      jets_for_minheap.pop_back();
      minheap.update(jetI-head, _bj_diJ(jetI));
      jetI->label_minheap_update_done();
    }
    n--;
  }

  // final cleaning up;
  delete[] briefjets;
}

bool 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 _history, and PseudoJet::cluster_hist_index().

Referenced by object_in_jet().

                                                       {
    return jet.cluster_hist_index() >= 0 
      && jet.cluster_hist_index() < int(_history.size());
  }

void 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().

                                    {

  if (!_first_time) {return;}
  _first_time = false;
  
  
  //Symp. Discr. Alg, p.472 (2002) and  CGAL (http://www.cgal.org);

  cout << "#--------------------------------------------------------------------------\n";
  cout << "#                      FastJet release " << fastjet_version << endl;
  cout << "#            Written by M. Cacciari, G.P. Salam and G. Soyez            \n"; 
  cout << "#                         http://www.fastjet.fr                         \n"; 
  cout << "#                                                                       \n";
  cout << "# Longitudinally invariant Kt, anti-Kt, and inclusive Cambridge/Aachen  \n";
  cout << "# clustering using fast geometric algorithms, with area measures and optional\n";
  cout << "# external jet-finder plugins.                                          \n";
  cout << "# Please cite Phys. Lett. B641 (2006) [hep-ph/0512210] if you use this code.\n";
  cout << "#                                                                       \n";
  cout << "# This package uses T.Chan's closest pair algorithm, Proc.13th ACM-SIAM \n";
  cout << "# Symp. Discr. Alg, p.472 (2002), S.Fortune's Voronoi algorithm and code " ;
#ifndef DROP_CGAL
  cout << endl << "# and CGAL: http://www.cgal.org/";
#endif  // DROP_CGAL
  cout << ".\n";
  cout << "#-------------------------------------------------------------------------\n";
  // make sure we really have the output done.
  cout.flush();
}

void ClusterSequence::_print_tiles

(

TiledJet *

briefjets

)

const [private]

output the contents of the tiles

Definition at line 212 of file ClusterSequence_TiledN2.cc.

References _tiles.

                                                              {
  for (vector<Tile>::const_iterator tile = _tiles.begin(); 
       tile < _tiles.end(); tile++) {
    cout << "Tile " << tile - _tiles.begin()<<" = ";
    vector<int> list;
    for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) {
      list.push_back(jetI-briefjets);
      //cout <<" "<<jetI-briefjets;
    }
    sort(list.begin(),list.end());
    for (unsigned int i = 0; i < list.size(); i++) {cout <<" "<<list[i];}
    cout <<"\n";
  }
}

void 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().

                                            {

  // the array that will be overwritten here will be one
  // of pointers to jets.
  vector<PseudoJet *> jetsp(_jets.size());
  vector<int>         indices(_jets.size());

  for (size_t i = 0; i<_jets.size(); i++) {
    jetsp[i] = & _jets[i];
    indices[i] = i;
  }

  for (int n = jetsp.size(); n > 0; n--) {
    int ii, jj;
    // find smallest beam distance [remember jet_scale_for_algorithm 
    // will return kt^2 for the kt-algorithm and 1 for the Cambridge/Aachen]
    double ymin = jet_scale_for_algorithm(*(jetsp[0]));
    ii = 0; jj = -2;
    for (int i = 0; i < n; i++) {
      double yiB = jet_scale_for_algorithm(*(jetsp[i]));
      if (yiB < ymin) {
        ymin = yiB; ii = i; jj = -2;}
    }

    // find smallest distance between pair of jetsp
    for (int i = 0; i < n-1; i++) {
      for (int j = i+1; j < n; j++) {
        //double y = jetsp[i]->kt_distance(*jetsp[j])*_invR2;
        double y = min(jet_scale_for_algorithm(*(jetsp[i])), 
                       jet_scale_for_algorithm(*(jetsp[j])))
                    * jetsp[i]->plain_distance(*jetsp[j])*_invR2;
        if (y < ymin) {ymin = y; ii = i; jj = j;}
      }
    }

    // output recombination sequence
    // old "ktclus" way of labelling
    //cout <<n<< " "<< ii+1 << " with " << jj+1 << "; y = "<< ymin<<endl;
    // new delaunay way of labelling
    int jjindex_or_beam, iiindex;
    if (jj < 0) {jjindex_or_beam = BeamJet; iiindex = indices[ii];} 
    else {
      jjindex_or_beam = max(indices[ii],indices[jj]);
      iiindex =         min(indices[ii],indices[jj]);
    }

    // now recombine
    int newn = 2*jetsp.size() - n;
    if (jj >= 0) {
      // combine pair
      int nn; // new jet index
      _do_ij_recombination_step(jetsp[ii]-&_jets[0], 
                                jetsp[jj]-&_jets[0], ymin, nn);
      
      // sort out internal bookkeeping
      jetsp[ii] = &_jets[nn];
      // have jj point to jet that was pointed at by n-1 
      // (since original jj is no longer current, so put n-1 into jj)
      jetsp[jj] = jetsp[n-1];
      indices[ii] = newn;
      indices[jj] = indices[n-1];

      //OBS_jets.push_back(*jetsp[ii] + *jetsp[jj]);
      //OBSjetsp[ii] = &_jets[_jets.size()-1];
      //OBS// have jj point to jet that was pointed at by n-1 
      //OBS// (since original jj is no longer current, so put n-1 into jj)
      //OBSjetsp[jj] = jetsp[n-1];
      //OBS
      //OBSindices[ii] = newn;
      //OBSindices[jj] = indices[n-1];
      //OBS_add_step_to_history(newn,iiindex,
      //OBS                           jjindex_or_beam,_jets.size()-1,ymin);
    } else {
      // combine ii with beam
      _do_iB_recombination_step(jetsp[ii]-&_jets[0], ymin);
      // put last jet (pointer) in place of ii (which has disappeared)
      jetsp[ii] = jetsp[n-1];
      indices[ii] = indices[n-1];
      //OBS_add_step_to_history(newn,iiindex,jjindex_or_beam,Invalid, ymin);
    }
  }

}

template<class BJ >

void ClusterSequence::_simple_N2_cluster

(

)

[inline, private]

void 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().

                                                  {  
  _simple_N2_cluster<BriefJet>();
}

void 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().

                                                    {  
  _simple_N2_cluster<EEBriefJet>();
}

int 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 twopi.

                                                                             {
  int ieta, iphi;
  if      (eta <= _tiles_eta_min) {ieta = 0;}
  else if (eta >= _tiles_eta_max) {ieta = _tiles_ieta_max-_tiles_ieta_min;}
  else {
    //ieta = int(floor((eta - _tiles_eta_min) / _tile_size_eta));
    ieta = int(((eta - _tiles_eta_min) / _tile_size_eta));
    // following needed in case of rare but nasty rounding errors
    if (ieta > _tiles_ieta_max-_tiles_ieta_min) {
      ieta = _tiles_ieta_max-_tiles_ieta_min;} 
  }
  // allow for some extent of being beyond range in calculation of phi
  // as well
  //iphi = (int(floor(phi/_tile_size_phi)) + _n_tiles_phi) % _n_tiles_phi;
  // with just int and no floor, things run faster but beware
  iphi = int((phi+twopi)/_tile_size_phi) % _n_tiles_phi;
  return (iphi + ieta * _n_tiles_phi);
}

int ClusterSequence::_tile_index	(	int	ieta,
		int	iphi
	)			const [inline, private]

Definition at line 728 of file ClusterSequence.hh.

References _n_tiles_phi, and _tiles_ieta_min.

Referenced by _initialise_tiles(), and _tj_set_jetinfo().

                                                    {
    // note that (-1)%n = -1 so that we have to add _n_tiles_phi
    // before performing modulo operation
    return (ieta-_tiles_ieta_min)*_n_tiles_phi
                  + (iphi+_n_tiles_phi) % _n_tiles_phi;
  }

void 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, ClusterSequence::TiledJet::_jets_index, _R2, _tiles, _tj_set_jetinfo(), ClusterSequence::Tile::begin_tiles, ClusterSequence::Tile::end_tiles, ClusterSequence::Tile::head, n_tile_neighbours, ClusterSequence::TiledJet::next, ClusterSequence::TiledJet::NN, ClusterSequence::TiledJet::NN_dist, ClusterSequence::TiledJet::previous, and ClusterSequence::TiledJet::tile_index.

Referenced by _initialise_and_run().

                                        {

  _initialise_tiles();

  int n = _jets.size();
  TiledJet * briefjets = new TiledJet[n];
  TiledJet * jetA = briefjets, * jetB;
  TiledJet oldB;
  

  // will be used quite deep inside loops, but declare it here so that
  // memory (de)allocation gets done only once
  vector<int> tile_union(3*n_tile_neighbours);
  
  // initialise the basic jet info 
  for (int i = 0; i< n; i++) {
    _tj_set_jetinfo(jetA, i);
    //cout << i<<": "<<jetA->tile_index<<"\n";
    jetA++; // move on to next entry of briefjets
  }
  TiledJet * tail = jetA; // a semaphore for the end of briefjets
  TiledJet * head = briefjets; // a nicer way of naming start

  // set up the initial nearest neighbour information
  vector<Tile>::const_iterator tile;
  for (tile = _tiles.begin(); tile != _tiles.end(); tile++) {
    // first do it on this tile
    for (jetA = tile->head; jetA != NULL; jetA = jetA->next) {
      for (jetB = tile->head; jetB != jetA; jetB = jetB->next) {
        double dist = _bj_dist(jetA,jetB);
        if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;}
        if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;}
      }
    }
    // then do it for RH tiles
    for (Tile ** RTile = tile->RH_tiles; RTile != tile->end_tiles; RTile++) {
      for (jetA = tile->head; jetA != NULL; jetA = jetA->next) {
        for (jetB = (*RTile)->head; jetB != NULL; jetB = jetB->next) {
          double dist = _bj_dist(jetA,jetB);
          if (dist < jetA->NN_dist) {jetA->NN_dist = dist; jetA->NN = jetB;}
          if (dist < jetB->NN_dist) {jetB->NN_dist = dist; jetB->NN = jetA;}
        }
      }
    }
  }
  
  // now create the diJ (where J is i's NN) table -- remember that 
  // we differ from standard normalisation here by a factor of R2
  double * diJ = new double[n];
  jetA = head;
  for (int i = 0; i < n; i++) {
    diJ[i] = _bj_diJ(jetA);
    jetA++; // have jetA follow i
  }

  // now run the recombination loop
  int history_location = n-1;
  while (tail != head) {

    // find the minimum of the diJ on this round
    double diJ_min = diJ[0];
    int diJ_min_jet = 0;
    for (int i = 1; i < n; i++) {
      if (diJ[i] < diJ_min) {diJ_min_jet = i; diJ_min  = diJ[i];}
    }

    // do the recombination between A and B
    history_location++;
    jetA = & briefjets[diJ_min_jet];
    jetB = jetA->NN;
    // put the normalisation back in
    diJ_min *= _invR2; 

    //if (n == 19) {cout << "Hello "<<jetA-head<<" "<<jetB-head<<"\n";}

    //cout <<" WILL RECOMBINE "<< jetA-briefjets<<" "<<jetB-briefjets<<"\n";

    if (jetB != NULL) {
      // jet-jet recombination
      // If necessary relabel A & B to ensure jetB < jetA, that way if
      // the larger of them == newtail then that ends up being jetA and 
      // the new jet that is added as jetB is inserted in a position that
      // has a future!
      if (jetA < jetB) {swap(jetA,jetB);}

      int nn; // new jet index
      _do_ij_recombination_step(jetA->_jets_index, jetB->_jets_index, diJ_min, nn);

      //OBS// get the two history indices
      //OBSint hist_a = _jets[jetA->_jets_index].cluster_hist_index();
      //OBSint hist_b = _jets[jetB->_jets_index].cluster_hist_index();
      //OBS// create the recombined jet
      //OBS_jets.push_back(_jets[jetA->_jets_index] + _jets[jetB->_jets_index]);
      //OBSint nn = _jets.size() - 1;
      //OBS_jets[nn].set_cluster_hist_index(history_location);
      //OBS// update history
      //OBS//cout <<n-1<<" "<<jetA-head<<" "<<jetB-head<<"; ";
      //OBS_add_step_to_history(history_location, 
      //OBS                min(hist_a,hist_b),max(hist_a,hist_b),
      //OBS                        nn, diJ_min);

      // what was jetB will now become the new jet
      _bj_remove_from_tiles(jetA);
      oldB = * jetB;  // take a copy because we will need it...
      _bj_remove_from_tiles(jetB);
      _tj_set_jetinfo(jetB, nn); // also registers the jet in the tiling
    } else {
      // jet-beam recombination
      _do_iB_recombination_step(jetA->_jets_index, diJ_min);
            
      //OBS// get the hist_index
      //OBSint hist_a = _jets[jetA->_jets_index].cluster_hist_index();
      //OBS//cout <<n-1<<" "<<jetA-head<<" "<<-1<<"; ";
      //OBS_add_step_to_history(history_location,hist_a,BeamJet,Invalid,diJ_min); 
      _bj_remove_from_tiles(jetA);
    }

    // first establish the set of tiles over which we are going to 
    // have to run searches for updated and new nearest-neighbours
    int n_near_tiles = 0;
    _add_neighbours_to_tile_union(jetA->tile_index, tile_union, n_near_tiles);
    if (jetB != NULL) {
      bool sort_it = false;
      if (jetB->tile_index != jetA->tile_index) {
        sort_it = true;
        _add_neighbours_to_tile_union(jetB->tile_index,tile_union,n_near_tiles);
      }
      if (oldB.tile_index != jetA->tile_index && 
          oldB.tile_index != jetB->tile_index) {
        sort_it = true;
        _add_neighbours_to_tile_union(oldB.tile_index,tile_union,n_near_tiles);
      }

      if (sort_it) {
        // sort the tiles before then compressing the list
        sort(tile_union.begin(), tile_union.begin()+n_near_tiles);
        // and now condense the list
        int nnn = 1;
        for (int i = 1; i < n_near_tiles; i++) {
          if (tile_union[i] != tile_union[nnn-1]) {
            tile_union[nnn] = tile_union[i]; 
            nnn++;
          }
        }
        n_near_tiles = nnn;
      }
    }

    // now update our nearest neighbour info and diJ table
    // first reduce size of table
    tail--; n--;
    if (jetA == tail) {
      // there is nothing to be done
    } else {
      // Copy last jet contents and diJ info into position of jetA
      *jetA = *tail;
      diJ[jetA - head] = diJ[tail-head];
      // IN the tiling fix pointers to tail and turn them into
      // pointers to jetA (from predecessors, successors and the tile
      // head if need be)
      if (jetA->previous == NULL) {
        _tiles[jetA->tile_index].head = jetA;
      } else {
        jetA->previous->next = jetA;
      }
      if (jetA->next != NULL) {jetA->next->previous = jetA;}
    }

    // Initialise jetB's NN distance as well as updating it for 
    // other particles.
    for (int itile = 0; itile < n_near_tiles; itile++) {
      Tile * tile = &_tiles[tile_union[itile]];
      for (TiledJet * jetI = tile->head; jetI != NULL; jetI = jetI->next) {
        // see if jetI had jetA or jetB as a NN -- if so recalculate the NN
        if (jetI->NN == jetA || (jetI->NN == jetB && jetB != NULL)) {
          jetI->NN_dist = _R2;
          jetI->NN      = NULL;
          // now go over tiles that are neighbours of I (include own tile)
          for (Tile ** near_tile  = tile->begin_tiles; 
                       near_tile != tile->end_tiles; near_tile++) {
            // and then over the contents of that tile
            for (TiledJet * jetJ  = (*near_tile)->head; 
                            jetJ != NULL; jetJ = jetJ->next) {
              double dist = _bj_dist(jetI,jetJ);
              if (dist < jetI->NN_dist && jetJ != jetI) {
                jetI->NN_dist = dist; jetI->NN = jetJ;
              }
            }
          }
          diJ[jetI-head] = _bj_diJ(jetI); // update diJ 
        }
        // check whether new jetB is closer than jetI's current NN and
        // if need to update things
        if (jetB != NULL) {
          double dist = _bj_dist(jetI,jetB);
          if (dist < jetI->NN_dist) {
            if (jetI != jetB) {
              jetI->NN_dist = dist;
              jetI->NN = jetB;
              diJ[jetI-head] = _bj_diJ(jetI); // update diJ...
            }
          }
          if (dist < jetB->NN_dist) {
            if (jetI != jetB) {
              jetB->NN_dist = dist;
              jetB->NN      = jetI;}
          }
        }
      }
    }


    if (jetB != NULL) {diJ[jetB-head] = _bj_diJ(jetB);}
    //cout << n<<" "<<briefjets[95].NN-briefjets<<" "<<briefjets[95].NN_dist <<"\n";

    // remember to update pointers to tail
    for (Tile ** near_tile = _tiles[tail->tile_index].begin_tiles; 
                 near_tile!= _tiles[tail->tile_index].end_tiles; near_tile++){
      // and then the contents of that tile
      for (TiledJet * jetJ = (*near_tile)->head; 
                     jetJ != NULL; jetJ = jetJ->next) {
        if (jetJ->NN == tail) {jetJ->NN = jetA;}
      }
    }

    //for (int i = 0; i < n; i++) {
    //  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";}
    //}


    if (jetB != NULL) {diJ[jetB-head] = _bj_diJ(jetB);}
    //cout << briefjets[95].NN-briefjets<<" "<<briefjets[95].NN_dist <<"\n";

  }

  // final cleaning up;
  delete[] diJ;
  delete[] briefjets;
}

void 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, ClusterSequence::Tile::head, ClusterSequence::TiledJet::next, and ClusterSequence::TiledJet::previous.

Referenced by _faster_tiled_N2_cluster(), _minheap_faster_tiled_N2_cluster(), and _tiled_N2_cluster().

                                                                     {
  // first call the generic setup
  _bj_set_jetinfo<>(jet, _jets_index);

  // Then do the setup specific to the tiled case.

  // Find out which tile it belonds to
  jet->tile_index = _tile_index(jet->eta, jet->phi);

  // Insert it into the tile's linked list of jets
  Tile * tile = &_tiles[jet->tile_index];
  jet->previous   = NULL;
  jet->next       = tile->head;
  if (jet->next != NULL) {jet->next->previous = jet;}
  tile->head      = jet;
}

template<class L >

void 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().

                                                                        {

  // this will ensure that we can point to jets without difficulties
  // arising.
  _jets.reserve(pseudojets.size()*2);

  // insert initial jets this way so that any type L that can be
  // converted to a pseudojet will work fine (basically PseudoJet
  // and any type that has [] subscript access to the momentum
  // components, such as CLHEP HepLorentzVector).
  for (unsigned int i = 0; i < pseudojets.size(); i++) {
    _jets.push_back(pseudojets[i]);}
  
}

void 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().

vector< PseudoJet > ClusterSequence::constituents

(

const PseudoJet &

jet

)

const

return a vector of the particles that make up jet

Definition at line 817 of file ClusterSequence.cc.

References add_constituents().

                                                                            {
  vector<PseudoJet> subjets;
  add_constituents(jet, subjets);
  return subjets;
}

double 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 537 of file ClusterSequence.cc.

References _history, and _initial_n.

Referenced by exclusive_ymerge().

                                                                 {
  assert(njets >= 0);
  if (njets >= _initial_n) {return 0.0;}
  return _history[2*_initial_n-njets-1].dij;
}

double 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.

Definition at line 549 of file ClusterSequence.cc.

References _history, and _initial_n.

Referenced by exclusive_ymerge_max().

                                                                     {
  assert(njets >= 0);
  if (njets >= _initial_n) {return 0.0;}
  return _history[2*_initial_n-njets-1].max_dij_so_far;
}

vector< PseudoJet > 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 465 of file ClusterSequence.cc.

References _history, _initial_n, _jet_def, _jets, _n_exclusive_warnings, cambridge_algorithm, ee_genkt_algorithm, ee_kt_algorithm, JetDefinition::extra_param(), genkt_algorithm, JetDefinition::jet_algorithm(), jets(), kt_algorithm, JetDefinition::plugin(), and plugin_algorithm.

                                                                          {

  // make sure the user does not ask for more than jets than there
  // were particles in the first place.
  assert (njets <= _initial_n);

  // provide a warning when extracting exclusive jets for algorithms 
  // that does not support it explicitly.
  // Native algorithm that support it are: kt, ee_kt, cambridge, 
  //   genkt and ee_genkt (both with p>=0)
  // For plugins, we check Plugin::exclusive_sequence_meaningful()
  if (( _jet_def.jet_algorithm() != kt_algorithm) &&
      ( _jet_def.jet_algorithm() != cambridge_algorithm) &&
      ( _jet_def.jet_algorithm() != ee_kt_algorithm) &&
      (((_jet_def.jet_algorithm() != genkt_algorithm) && 
        (_jet_def.jet_algorithm() != ee_genkt_algorithm)) || 
       (_jet_def.extra_param() <0)) &&
      ((_jet_def.jet_algorithm() != plugin_algorithm) ||
       (!_jet_def.plugin()->exclusive_sequence_meaningful())) &&
      (_n_exclusive_warnings < 5)) {
    _n_exclusive_warnings++;
    cerr << "FastJet WARNING: dcut and exclusive jets for jet-finders other than kt should be interpreted with care." << endl;
  }


  // calculate the point where we have to stop the clustering.
  // relation between stop_point, njets assumes one extra jet disappears
  // at each clustering.
  int stop_point = 2*_initial_n - njets;

  // some sanity checking to make sure that e+e- does not give us
  // surprises (should we ever implement e+e-)...
  if (2*_initial_n != static_cast<int>(_history.size())) {
    ostringstream err;
    err << "2*_initial_n != _history.size() -- this endangers internal assumptions!\n";
    throw Error(err.str());
    //assert(false);
  }

  // now go forwards and reconstitute the jets that we have --
  // basically for any history element, see if the parent jets to
  // which it refers were created before the stopping point -- if they
  // were then add them to the list, otherwise they are subsequent
  // recombinations of the jets that we are looking for.
  vector<PseudoJet> jets;
  for (unsigned int i = stop_point; i < _history.size(); i++) {
    int parent1 = _history[i].parent1;
    if (parent1 < stop_point) {
      jets.push_back(_jets[_history[parent1].jetp_index]);
    }
    int parent2 = _history[i].parent2;
    if (parent2 < stop_point && parent2 > 0) {
      jets.push_back(_jets[_history[parent2].jetp_index]);
    }
    
  }

  // sanity check...
  if (static_cast<int>(jets.size()) != njets) {
    ostringstream err;
    err << "ClusterSequence::exclusive_jets: size of returned vector ("
         <<jets.size()<<") does not coincide with requested number of jets ("
         <<njets<<")";
    throw Error(err.str());
  }

  return jets;
}

vector< PseudoJet > 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 457 of file ClusterSequence.cc.

References n_exclusive_jets().

Referenced by exclusive_jets_ycut().

                                                                            {
  int njets = n_exclusive_jets(dcut);
  return exclusive_jets(njets);
}

std::vector<PseudoJet> ClusterSequence::exclusive_jets_ycut

(

double

ycut

)

const [inline]

the exclusive jets obtained at the given ycut

Definition at line 145 of file ClusterSequence.hh.

References exclusive_jets(), and n_exclusive_jets_ycut().

                                                               {
    int njets = n_exclusive_jets_ycut(ycut);
    return exclusive_jets(njets);
  }

double 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.

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 621 of file ClusterSequence.cc.

References get_subhist_set().

                                                                                 {
  set<const history_element*> subhist;

  // get the set of history elements that correspond to subjets at
  // scale dcut
  get_subhist_set(subhist, jet, -1.0, nsub);
  
  set<const history_element*>::iterator highest = subhist.end();
  highest--;
  return (*highest)->dij;
}

double 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.

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 642 of file ClusterSequence.cc.

References get_subhist_set().

                                                                                     {

  set<const history_element*> subhist;

  // get the set of history elements that correspond to subjets at
  // scale dcut
  get_subhist_set(subhist, jet, -1.0, nsub);
  
  set<const history_element*>::iterator highest = subhist.end();
  highest--;
  return (*highest)->max_dij_so_far;
}

std::vector< PseudoJet > 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).

requires n ln n time

Definition at line 597 of file ClusterSequence.cc.

                                        {

  set<const history_element*> subhist;

  // get the set of history elements that correspond to subjets at
  // scale dcut
  get_subhist_set(subhist, jet, -1.0, n);

  // now transfer this into a sequence of jets
  vector<PseudoJet> subjets;
  subjets.reserve(subhist.size());
  for (set<const history_element*>::iterator elem = subhist.begin(); 
       elem != subhist.end(); elem++) {
    subjets.push_back(_jets[(*elem)->jetp_index]);
  }
  return subjets;
}

std::vector< PseudoJet > 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.

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 561 of file ClusterSequence.cc.

                                                      {

  set<const history_element*> subhist;

  // get the set of history elements that correspond to subjets at
  // scale dcut
  get_subhist_set(subhist, jet, dcut, 0);

  // now transfer this into a sequence of jets
  vector<PseudoJet> subjets;
  subjets.reserve(subhist.size());
  for (set<const history_element*>::iterator elem = subhist.begin(); 
       elem != subhist.end(); elem++) {
    subjets.push_back(_jets[(*elem)->jetp_index]);
  }
  return subjets;
}

double 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.

References exclusive_dmerge(), and Q2().

{return exclusive_dmerge(njets) / Q2();}

double 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_dmerge_max(), and Q2().

{return exclusive_dmerge_max(njets)/Q2();}

const Extras* ClusterSequence::extras

(

)

const [inline]

returns a pointer to the extras object (may be null)

Definition at line 320 of file ClusterSequence.hh.

References _extras.

{return _extras.get();}

void 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().

bool ClusterSequence::has_child	(	const PseudoJet &	jet,
		const PseudoJet *&	childp
	)			const

Version of has_child that sets a pointer to the child if the child exists;.

Definition at line 769 of file ClusterSequence.cc.

References _history, _jets, ClusterSequence::history_element::child, and PseudoJet::cluster_hist_index().

                                                                                       {

  const history_element & hist = _history[jet.cluster_hist_index()];

  // check that this jet has a child and that the child corresponds to
  // a true jet [RETHINK-IF-CHANGE-NUMBERING: what is the right
  // behaviour if the child is the same jet but made inclusive...?]
  if (hist.child >= 0 && _history[hist.child].jetp_index >= 0) {
    childp = &(_jets[_history[hist.child].jetp_index]);
    return true;
  } else {
    childp = NULL;
    return false;
  }
}

bool ClusterSequence::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

Definition at line 747 of file ClusterSequence.cc.

Referenced by object_in_jet().

                                                                              {

  //const history_element & hist = _history[jet.cluster_hist_index()];
  //
  //if (hist.child >= 0) {
  //  child = _jets[_history[hist.child].jetp_index];
  //  return true;
  //} else {
  //  child = PseudoJet(0.0,0.0,0.0,0.0);
  //  return false;
  //}
  const PseudoJet * childp;
  bool res = has_child(jet, childp);
  if (res) {
    child = *childp;
    return true;
  } else {
    child = PseudoJet(0.0,0.0,0.0,0.0);
    return false;
  }
}

bool 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 it has no parents it returns false and sets parent1 and parent2 to zero

Definition at line 721 of file ClusterSequence.cc.

References _history, _jets, PseudoJet::cluster_hist_index(), ClusterSequence::history_element::parent1, ClusterSequence::history_element::parent2, and PseudoJet::perp2().

                                                         {

  const history_element & hist = _history[jet.cluster_hist_index()];

  // make sure we do not run into any unexpected situations --
  // i.e. both parents valid, or neither
  assert ((hist.parent1 >= 0 && hist.parent2 >= 0) || 
          (hist.parent1 < 0 && hist.parent2 < 0));

  if (hist.parent1 < 0) {
    parent1 = PseudoJet(0.0,0.0,0.0,0.0);
    parent2 = parent1;
    return false;
  } else {
    parent1 = _jets[_history[hist.parent1].jetp_index];
    parent2 = _jets[_history[hist.parent2].jetp_index];
    // order the parents in decreasing pt
    if (parent1.perp2() < parent2.perp2()) swap(parent1,parent2);
    return true;
  }
}

bool ClusterSequence::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

Definition at line 790 of file ClusterSequence.cc.

References _history, _jets, ClusterSequence::history_element::child, PseudoJet::cluster_hist_index(), ClusterSequence::history_element::parent1, and ClusterSequence::history_element::parent2.

                                                         {

  const history_element & hist = _history[jet.cluster_hist_index()];

  // make sure we have a child and that the child does not correspond
  // to a clustering with the beam (or some other invalid quantity)
  if (hist.child >= 0 && _history[hist.child].parent2 >= 0) {
    const history_element & child_hist = _history[hist.child];
    if (child_hist.parent1 == jet.cluster_hist_index()) {
      // partner will be child's parent2 -- for iB clustering
      // parent2 will not be valid
      partner = _jets[_history[child_hist.parent2].jetp_index];
    } else {
      // partner will be child's parent1
      partner = _jets[_history[child_hist.parent1].jetp_index];
    }
    return true;
  } else {
    partner = PseudoJet(0.0,0.0,0.0,0.0);
    return false;
  }
}

const std::vector< ClusterSequence::history_element > & 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 ClusterSequenceActiveArea::_transfer_ghost_free_history(), and fastjet::NestedDefsPlugin::run_clustering().

                                                                                         {
  return _history;
}

vector< PseudoJet > 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 386 of file ClusterSequence.cc.

References _history, _jet_algorithm, _jets, antikt_algorithm, BeamJet, cambridge_algorithm, cambridge_for_passive_algorithm, ee_genkt_algorithm, ee_kt_algorithm, genkt_algorithm, jets(), kt_algorithm, PseudoJet::perp2(), and plugin_algorithm.

Referenced by ClusterSequenceAreaBase::empty_area(), ClusterSequenceAreaBase::get_median_rho_and_sigma(), ClusterSequenceAreaBase::parabolic_pt_per_unit_area(), ClusterSequenceActiveArea::pt_per_unit_area(), and ClusterSequenceAreaBase::subtracted_jets().

                                                                            {
  double dcut = ptmin*ptmin;
  int i = _history.size() - 1; // last jet
  vector<PseudoJet> jets;
  if (_jet_algorithm == kt_algorithm) {
    while (i >= 0) {
      // with our specific definition of dij and diB (i.e. R appears only in 
      // dij), then dij==diB is the same as the jet.perp2() and we can exploit
      // this in selecting the jets...
      if (_history[i].max_dij_so_far < dcut) {break;}
      if (_history[i].parent2 == BeamJet && _history[i].dij >= dcut) {
        // for beam jets
        int parent1 = _history[i].parent1;
        jets.push_back(_jets[_history[parent1].jetp_index]);}
      i--;
    }
  } else if (_jet_algorithm == cambridge_algorithm) {
    while (i >= 0) {
      // inclusive jets are all at end of clustering sequence in the
      // Cambridge algorithm -- so if we find a non-exclusive jet, then
      // we can exit
      if (_history[i].parent2 != BeamJet) {break;}
      int parent1 = _history[i].parent1;
      const PseudoJet & jet = _jets[_history[parent1].jetp_index];
      if (jet.perp2() >= dcut) {jets.push_back(jet);}
      i--;
    }
  } else if (_jet_algorithm == plugin_algorithm 
             || _jet_algorithm == ee_kt_algorithm
             || _jet_algorithm == antikt_algorithm
             || _jet_algorithm == genkt_algorithm
             || _jet_algorithm == ee_genkt_algorithm
             || _jet_algorithm == cambridge_for_passive_algorithm) {
    // for inclusive jets with a plugin algorithm, we make no
    // assumptions about anything (relation of dij to momenta,
    // ordering of the dij, etc.)
    while (i >= 0) {
      if (_history[i].parent2 == BeamJet) {
        int parent1 = _history[i].parent1;
        const PseudoJet & jet = _jets[_history[parent1].jetp_index];
        if (jet.perp2() >= dcut) {jets.push_back(jet);}
      }
      i--;
    }
  } else {throw Error("cs::inclusive_jets(...): Unrecognized jet algorithm");}
  return jets;
}

const JetDefinition& ClusterSequence::jet_def

(

)

const [inline]

return a reference to the jet definition

Definition at line 264 of file ClusterSequence.hh.

References _jet_def.

Referenced by _bj_set_jetinfo(), ClusterSequenceAreaBase::_check_jet_alg_good_for_median(), _initialise_and_run(), ClusterSequence1GhostPassiveArea::_run_1GPA(), ClusterSequenceActiveArea::_run_AA(), ClusterSequenceArea::_warn_if_range_unsuitable(), ClusterSequencePassiveArea::empty_area(), jet_scale_for_algorithm(), and ClusterSequenceAreaBase::n_empty_jets().

{return _jet_def;}

double 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 321 of file ClusterSequence.cc.

References _jet_algorithm, _jet_def, antikt_algorithm, cambridge_algorithm, cambridge_for_passive_algorithm, JetDefinition::extra_param(), genkt_algorithm, jet_def(), PseudoJet::kt2(), and kt_algorithm.

Referenced by _add_ktdistance_to_map(), _bj_set_jetinfo(), and _really_dumb_cluster().

                                                               {
  if (_jet_algorithm == kt_algorithm)             {return jet.kt2();}
  else if (_jet_algorithm == cambridge_algorithm) {return 1.0;}
  else if (_jet_algorithm == antikt_algorithm) {
    double kt2=jet.kt2();
    return kt2 > 1e-300 ? 1.0/kt2 : 1e300;
  } else if (_jet_algorithm == genkt_algorithm) {
    double kt2 = jet.kt2();
    double p   = jet_def().extra_param();
    if (p <= 0 && kt2 < 1e-300) kt2 = 1e-300; // dodgy safety check
    return pow(kt2, p);
  } else if (_jet_algorithm == cambridge_for_passive_algorithm) {
    double kt2 = jet.kt2();
    double lim = _jet_def.extra_param();
    if (kt2 < lim*lim && kt2 != 0.0) {
      return 1.0/kt2;
    } else {return 1.0;}
  } else {throw Error("Unrecognised jet algorithm");}
}

const std::vector< PseudoJet > & 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 ClusterSequenceAreaBase::subtracted_jets().

                                                                 {
  return _jets;
}

int 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 438 of file ClusterSequence.cc.

References _history, and _initial_n.

Referenced by exclusive_jets(), and n_exclusive_jets_ycut().

                                                                {

  // first locate the point where clustering would have stopped (i.e. the
  // first time max_dij_so_far > dcut)
  int i = _history.size() - 1; // last jet
  while (i >= 0) {
    if (_history[i].max_dij_so_far <= dcut) {break;}
    i--;
  }
  int stop_point = i + 1;
  // relation between stop_point, njets assumes one extra jet disappears
  // at each clustering.
  int njets = 2*_initial_n - stop_point;
  return njets;
}

int 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.

References n_exclusive_jets(), and Q2().

Referenced by exclusive_jets_ycut().

{return n_exclusive_jets(ycut*Q2());}

int 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()

Definition at line 583 of file ClusterSequence.cc.

References get_subhist_set().

                                                   {
  set<const history_element*> subhist;
  // get the set of history elements that correspond to subjets at
  // scale dcut
  get_subhist_set(subhist, jet, dcut, 0);
  return subhist.size();
}

unsigned int 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(), ClusterSequenceVoronoiArea::_initializeVA(), unclustered_particles(), and unique_history_order().

{return _initial_n;}

bool 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 698 of file ClusterSequence.cc.

References _potentially_valid(), PseudoJet::cluster_hist_index(), and has_child().

                                                                 {

  // make sure the object conceivably belongs to this clustering
  // sequence
  assert(_potentially_valid(object) && _potentially_valid(jet));

  const PseudoJet * this_object = &object;
  const PseudoJet * childp;
  while(true) {
    if (this_object->cluster_hist_index() == jet.cluster_hist_index()) {
      return true;
    } else if (has_child(*this_object, childp)) {this_object = childp;}
    else {return false;}
  }
}

std::vector<int> 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 ClusterSequence::plugin_activated

(

)

const [inline]

returns true when the plugin is allowed to run the show.

Definition at line 317 of file ClusterSequence.hh.

References _plugin_activated.

Referenced by plugin_record_iB_recombination(), plugin_record_ij_recombination(), and plugin_simple_N2_cluster().

{return _plugin_activated;}

void 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.

References _extras.

Referenced by fastjet::SISConeSphericalPlugin::run_clustering(), and fastjet::SISConePlugin::run_clustering().

                                                                     {
    _extras = extras_in;
  }

void 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.

References _do_iB_recombination_step(), and plugin_activated().

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().

                                                             {
    assert(plugin_activated());
    _do_iB_recombination_step(jet_i, diB);
  }

void 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 371 of file ClusterSequence.cc.

References _jets, and plugin_record_ij_recombination().

                                                     {

  plugin_record_ij_recombination(jet_i, jet_j, dij, newjet_k);

  // now transfer newjet into place
  int tmp_index = _jets[newjet_k].cluster_hist_index();
  _jets[newjet_k] = newjet;
  _jets[newjet_k].set_cluster_hist_index(tmp_index);
}

void 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.

References _do_ij_recombination_step(), and plugin_activated().

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().

                                                      {
    assert(plugin_activated());
    _do_ij_recombination_step(jet_i, jet_j, dij, newjet_k);
  }

template<class GBJ >

void 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.

References plugin_activated().

                                                       {
    assert(plugin_activated());
    _simple_N2_cluster<GBJ>();
  }

void 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 851 of file ClusterSequence.cc.

References print_jets_for_root().

                                                                            {
  std::ofstream ostr(filename.c_str());
  if (comment != "") ostr << "# " << comment << endl;
  print_jets_for_root(jets, ostr);
}

void 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_for_root().

double 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.

References _Qtot.

{return _Qtot;}

double ClusterSequence::Q2

(

)

const [inline]

return Q()^2

Definition at line 194 of file ClusterSequence.hh.

References _Qtot.

Referenced by exclusive_ymerge(), exclusive_ymerge_max(), n_exclusive_jets_ycut(), and fastjet::EECambridgePlugin::run_clustering().

{return _Qtot*_Qtot;}

static void 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.

References _default_jet_algorithm.

{_default_jet_algorithm = jet_algorithm;}

static void ClusterSequence::set_jet_finder

(

JetAlgorithm

jet_algorithm

)

[inline, static]

same as above for backward compatibility

Definition at line 373 of file ClusterSequence.hh.

References _default_jet_algorithm.

{_default_jet_algorithm = jet_algorithm;}

string ClusterSequence::strategy_string

(

)

const

Definition at line 287 of file ClusterSequence.cc.

References _strategy, N2MinHeapTiled, N2Plain, N2PoorTiled, N2Tiled, N3Dumb, NlnN, NlnN3pi, NlnN4pi, NlnNCam, NlnNCam2pi2R, NlnNCam4pi, and plugin_strategy.

Referenced by _delaunay_cluster().

                                                {
  string strategy;
  switch(_strategy) {
  case NlnN:
    strategy = "NlnN"; break;
  case NlnN3pi:
    strategy = "NlnN3pi"; break;
  case NlnN4pi:
    strategy = "NlnN4pi"; break;
  case N2Plain:
    strategy = "N2Plain"; break;
  case N2Tiled:
    strategy = "N2Tiled"; break;
  case N2MinHeapTiled:
    strategy = "N2MinHeapTiled"; break;
  case N2PoorTiled:
    strategy = "N2PoorTiled"; break;
  case N3Dumb:
    strategy = "N3Dumb"; break;
  case NlnNCam4pi:
    strategy = "NlnNCam4pi"; break;
  case NlnNCam2pi2R:
    strategy = "NlnNCam2pi2R"; break;
  case NlnNCam:
    strategy = "NlnNCam"; break; // 2piMultD
  case plugin_strategy:
    strategy = "plugin strategy"; break;
  default:
    strategy = "Unrecognized";
  }
  return strategy;
}  

Strategy 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.

References _strategy.

Referenced by ClusterSequenceActiveArea::_transfer_ghost_free_history().

{return _strategy;}

void 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.

Definition at line 347 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 ClusterSequencePassiveArea::_initialise_and_run_PA(), and ClusterSequenceArea::initialize_and_run_cswa().

                                                                       {

  // the metadata
  _jet_def                 = from_seq._jet_def                ;
  _writeout_combinations   = from_seq._writeout_combinations  ;
  _initial_n               = from_seq._initial_n              ;
  _Rparam                  = from_seq._Rparam                 ;
  _R2                      = from_seq._R2                     ;
  _invR2                   = from_seq._invR2                  ;
  _strategy                = from_seq._strategy               ;
  _jet_algorithm           = from_seq._jet_algorithm          ;
  _plugin_activated        = from_seq._plugin_activated       ;

  // the data
  _jets     = from_seq._jets;
  _history  = from_seq._history;
  // the following transferse ownership of the extras from the from_seq
  _extras   = from_seq._extras;

}

vector< PseudoJet > 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 1037 of file ClusterSequence.cc.

References _history, _jets, Invalid, and n_particles().

                                                               {
  vector<PseudoJet> unclustered;
  for (unsigned i = 0; i < n_particles() ; i++) {
    if (_history[i].child == Invalid) 
      unclustered.push_back(_jets[_history[i].jetp_index]);
  }
  return unclustered;
}

vector< int > 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 982 of file ClusterSequence.cc.

References _extract_tree_children(), _history, fastjet::d0::inline_maths::min(), and n_particles().

Referenced by ClusterSequence1GhostPassiveArea::_run_1GPA(), and ClusterSequenceActiveArea::_run_AA().

                                                        {

  // first construct an array that will tell us the lowest constituent
  // of a given jet -- this will always be one of the original
  // particles, whose order is well defined and so will help us to
  // follow the tree in a unique manner.
  valarray<int> lowest_constituent(_history.size());
  int hist_n = _history.size();
  lowest_constituent = hist_n; // give it a large number
  for (int i = 0; i < hist_n; i++) {
    // sets things up for the initial partons
    lowest_constituent[i] = min(lowest_constituent[i],i); 
    // propagates them through to the children of this parton
    if (_history[i].child > 0) lowest_constituent[_history[i].child] 
      = min(lowest_constituent[_history[i].child],lowest_constituent[i]);
  }

  // establish an array for what we have and have not extracted so far
  valarray<bool> extracted(_history.size()); extracted = false;
  vector<int> unique_tree;
  unique_tree.reserve(_history.size());

  // now work our way through the tree
  for (unsigned i = 0; i < n_particles(); i++) {
    if (!extracted[i]) {
      unique_tree.push_back(i);
      extracted[i] = true;
      _extract_tree_children(i, extracted, lowest_constituent, unique_tree);
    }
  }

  return unique_tree;
}

---

Member Data Documentation

JetAlgorithm ClusterSequence::_default_jet_algorithm = kt_algorithm [static, protected]

Definition at line 482 of file ClusterSequence.hh.

Referenced by _initialise_and_run(), set_jet_algorithm(), and set_jet_finder().

std::auto_ptr<Extras> ClusterSequence::_extras [private]

Definition at line 564 of file ClusterSequence.hh.

Referenced by extras(), plugin_associate_extras(), and transfer_from_sequence().

bool 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> 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(), ClusterSequencePassiveArea::_initialise_and_run_PA(), ClusterSequenceVoronoiArea::_initializeVA(), _potentially_valid(), 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 ClusterSequence::_initial_n [protected]

Definition at line 554 of file ClusterSequence.hh.

Referenced by _CP2DChan_limited_cluster(), _fill_initial_history(), exclusive_dmerge(), exclusive_dmerge_max(), exclusive_jets(), n_exclusive_jets(), n_particles(), and transfer_from_sequence().

double 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 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 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(), ClusterSequencePassiveArea::_initialise_and_run_PA(), ClusterSequenceVoronoiArea::_initializeVA(), exclusive_jets(), jet_def(), jet_scale_for_algorithm(), and transfer_from_sequence().

std::vector<PseudoJet> 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 ClusterSequenceActiveAreaExplicitGhosts::_add_ghosts(), _add_ktdistance_to_map(), _add_step_to_history(), _bj_of_hindex(), _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(), ClusterSequencePassiveArea::_initialise_and_run_PA(), _initialise_tiles(), ClusterSequenceVoronoiArea::_initializeVA(), _minheap_faster_tiled_N2_cluster(), _really_dumb_cluster(), ClusterSequenceActiveArea::_resize_and_zero_AA(), ClusterSequence1GhostPassiveArea::_run_1GPA(), 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 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 ClusterSequence::_n_tiles_phi [private]

Definition at line 724 of file ClusterSequence.hh.

Referenced by _initialise_tiles(), and _tile_index().

bool ClusterSequence::_plugin_activated [private]

Definition at line 563 of file ClusterSequence.hh.

Referenced by _decant_options(), _initialise_and_run(), plugin_activated(), and transfer_from_sequence().

double ClusterSequence::_Qtot [protected]

Definition at line 556 of file ClusterSequence.hh.

Referenced by _fill_initial_history(), Q(), and Q2().

double 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 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 ClusterSequence::_strategy [protected]

Definition at line 557 of file ClusterSequence.hh.

Referenced by _decant_options(), _delaunay_cluster(), _initialise_and_run(), ClusterSequenceActiveArea::_transfer_ghost_free_history(), strategy_string(), strategy_used(), and transfer_from_sequence().

double ClusterSequence::_tile_size_eta [private]

Definition at line 723 of file ClusterSequence.hh.

Referenced by _initialise_tiles(), and _tile_index().

double ClusterSequence::_tile_size_phi [private]

Definition at line 723 of file ClusterSequence.hh.

Referenced by _initialise_tiles(), and _tile_index().

std::vector<Tile> 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 ClusterSequence::_tiles_eta_max [private]

Definition at line 722 of file ClusterSequence.hh.

Referenced by _initialise_tiles(), and _tile_index().

double ClusterSequence::_tiles_eta_min [private]

Definition at line 722 of file ClusterSequence.hh.

Referenced by _initialise_tiles(), and _tile_index().

int ClusterSequence::_tiles_ieta_max [private]

Definition at line 724 of file ClusterSequence.hh.

Referenced by _initialise_tiles(), and _tile_index().

int ClusterSequence::_tiles_ieta_min [private]

Definition at line 724 of file ClusterSequence.hh.

Referenced by _initialise_tiles(), and _tile_index().

bool 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 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:

• include/fastjet/ClusterSequence.hh
• src/ClusterSequence.cc
• src/ClusterSequence_CP2DChan.cc
• src/ClusterSequence_Delaunay.cc
• src/ClusterSequence_DumbN3.cc
• src/ClusterSequence_N2.cc
• src/ClusterSequence_TiledN2.cc