Main Page | Namespace List | Class Hierarchy | Class List | Directories | File List | Namespace Members | Class Members | File Members

fastjet::ClusterSequenceActiveAreaExplicitGhosts Class Reference

Class that behaves essentially like ClusterSequence except that it also provides access to the area of a jet (which will be a random quantity. More...

#include <ClusterSequenceActiveAreaExplicitGhosts.hh>

Inheritance diagram for fastjet::ClusterSequenceActiveAreaExplicitGhosts:

Inheritance graph
[legend]
Collaboration diagram for fastjet::ClusterSequenceActiveAreaExplicitGhosts:

Collaboration graph
[legend]
List of all members.

Public Member Functions

template<class L>
 ClusterSequenceActiveAreaExplicitGhosts (const std::vector< L > &pseudojets, const JetDefinition &jet_def, const ActiveAreaSpec &area_spec, const bool &writeout_combinations=false)
 constructor using a ActiveAreaSpec to specify how the area is to be measured
template<class L>
void _initialise (const std::vector< L > &pseudojets, const JetDefinition &jet_def, const ActiveAreaSpec &area_spec, const bool &writeout_combinations)
 does the actual work of initialisation
unsigned int n_hard_particles () const
 returns the number of hard particles (i.e. those supplied by the user).
virtual double area (const PseudoJet &jet) const
 returns the area of a jet
virtual PseudoJet area_4vector (const PseudoJet &jet) const
 returns a four vector corresponding to the sum (E-scheme) of the ghost four-vectors composing the jet area, normalised such that for a small contiguous area the p_t of the extended_area jet is equal to area of the jet.
bool is_pure_ghost (const PseudoJet &jet) const
 true if a jet is made exclusively of ghosts
bool is_pure_ghost (int history_index) const
 true if the entry in the history index corresponds to a ghost; if hist_ix does not correspond to an actual particle (i.e.
double total_area () const
 returns the total area under study

Private Member Functions

void _add_ghosts (const ActiveAreaSpec &area_spec)
void _post_process ()
 routine to be called after the processing is done so as to establish summary information on all the jets (areas, whether pure ghost, etc.

Private Attributes

int _n_ghosts
double _ghost_area
std::vector< bool > _is_pure_ghost
std::vector< double > _areas
std::vector< PseudoJet_area_4vectors
unsigned int _initial_hard_n

Detailed Description

Class that behaves essentially like ClusterSequence except that it also provides access to the area of a jet (which will be a random quantity.

.. Figure out what to do about seeds later...)

Definition at line 47 of file ClusterSequenceActiveAreaExplicitGhosts.hh.


Constructor & Destructor Documentation

template<class L>
fastjet::ClusterSequenceActiveAreaExplicitGhosts::ClusterSequenceActiveAreaExplicitGhosts const std::vector< L > &  pseudojets,
const JetDefinition jet_def,
const ActiveAreaSpec area_spec,
const bool &  writeout_combinations = false
[inline]
 

constructor using a ActiveAreaSpec to specify how the area is to be measured

Definition at line 53 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

00057            : ClusterSequenceWithArea() {
00058            _initialise(pseudojets,jet_def,area_spec,writeout_combinations); }


Member Function Documentation

void fastjet::ClustSeqActAreaEG::_add_ghosts const ActiveAreaSpec area_spec  )  [private]
 

Definition at line 44 of file ClusterSequenceActiveAreaExplicitGhosts.cc.

References _ghost_area, _initial_hard_n, _is_pure_ghost, fastjet::ClusterSequence::_jets, _n_ghosts, fastjet::ActiveAreaSpec::actual_ghost_area(), fastjet::ActiveAreaSpec::add_ghosts(), and fastjet::ActiveAreaSpec::n_ghosts().

00045                                                            {
00046 
00047   // add the ghosts to the jets
00048   area_spec.add_ghosts(_jets);
00049 
00050   // now add labelling...
00051   for (unsigned i = _initial_hard_n; i < _jets.size(); i++) {
00052     //_jets[i].set_user_index(1);
00053     _is_pure_ghost.push_back(true);
00054   }
00055 
00056   // and record some info from the area_spec
00057   _ghost_area = area_spec.actual_ghost_area();
00058   _n_ghosts   = area_spec.n_ghosts();
00059 }

template<class L>
void fastjet::ClusterSequenceActiveAreaExplicitGhosts::_initialise const std::vector< L > &  pseudojets,
const JetDefinition jet_def,
const ActiveAreaSpec area_spec,
const bool &  writeout_combinations
 

does the actual work of initialisation

Definition at line 118 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

00121                                               {
00122   // don't reserve space yet -- will be done below
00123 
00124   // insert initial jets this way so that any type L that can be
00125   // converted to a pseudojet will work fine (basically PseudoJet
00126   // and any type that has [] subscript access to the momentum
00127   // components, such as CLHEP HepLorentzVector).
00128   for (unsigned int i = 0; i < pseudojets.size(); i++) {
00129     PseudoJet mom(pseudojets[i]);
00130     //mom.set_user_index(0); // for user's particles (user index now lost...)
00131     _jets.push_back(mom);
00132     _is_pure_ghost.push_back(false);
00133   }
00134 
00135   _initial_hard_n = _jets.size();
00136 
00137   _add_ghosts(area_spec);
00138 
00139   if (writeout_combinations) {
00140     std::cout << "# Printing particles including ghosts\n";
00141     for (unsigned j = 0; j < _jets.size(); j++) {
00142       printf("%5u %20.13f %20.13f %20.13e\n",
00143                j,_jets[j].rap(),_jets[j].phi_02pi(),_jets[j].kt2());
00144     }
00145     std::cout << "# Finished printing particles including ghosts\n";
00146   }
00147 
00148   // this will ensure that we can still point to jets without
00149   // difficulties arising!
00150   _jets.reserve(_jets.size()*2);
00151 
00152   // run the clustering
00153   _initialise_and_run(jet_def,writeout_combinations);
00154 
00155   // set up all other information
00156   _post_process();
00157 }

void fastjet::ClustSeqActAreaEG::_post_process  )  [private]
 

routine to be called after the processing is done so as to establish summary information on all the jets (areas, whether pure ghost, etc.

)

Definition at line 97 of file ClusterSequenceActiveAreaExplicitGhosts.cc.

References _area_4vectors, _areas, _ghost_area, fastjet::ClusterSequence::_history, fastjet::ClusterSequence::_initial_n, _is_pure_ghost, fastjet::ClusterSequence::_jet_def, fastjet::ClusterSequence::_jets, fastjet::ClusterSequence::BeamJet, and fastjet::JetDefinition::recombiner().

00097                                       {
00098 
00099   // sort out sizes
00100   _areas.resize(_history.size());
00101   _area_4vectors.resize(_history.size());
00102   _is_pure_ghost.resize(_history.size());
00103   
00104   // First set up areas for the initial particles (ghost=_ghost_area,
00105   // real particles = 0); recall that _initial_n here is the number of
00106   // particles including ghosts
00107   for (int i = 0; i < _initial_n; i++) {
00108     if (_is_pure_ghost[i]) {
00109       _areas[i] = _ghost_area;
00110       // normalise pt to be _ghost_area (NB we make use of fact that
00111       // for initial particles, jet and clust_hist index are the same).
00112       _area_4vectors[i] = (_ghost_area/_jets[i].perp()) * _jets[i];
00113     } else {
00114       _areas[i] = 0;
00115       _area_4vectors[i] = PseudoJet(0.0,0.0,0.0,0.0);
00116     }
00117   }
00118   
00119   // next follow the branching through and set up the areas 
00120   // and ghost-nature at each step of the clustering (rather than
00121   // each jet).
00122   for (unsigned i = _initial_n; i < _history.size(); i++) {
00123     if (_history[i].parent2 == BeamJet) {
00124       _is_pure_ghost[i]  = _is_pure_ghost[_history[i].parent1];
00125       _areas[i]          = _areas[_history[i].parent1];
00126       _area_4vectors[i] = _area_4vectors[_history[i].parent1];
00127     } else {
00128       _is_pure_ghost[i]  = _is_pure_ghost[_history[i].parent1] && 
00129                            _is_pure_ghost[_history[i].parent2]   ;
00130       _areas[i]          = _areas[_history[i].parent1] + 
00131                            _areas[_history[i].parent2]  ;                          
00132       _jet_def.recombiner()->recombine(_area_4vectors[_history[i].parent1], 
00133                            _area_4vectors[_history[i].parent2],
00134                            _area_4vectors[i]);  
00135 //      _area_4vectors[i] = _area_4vectors[_history[i].parent1] + 
00136 //                          _area_4vectors[_history[i].parent2]  ;
00137     }
00138 
00139   }
00140   
00141 }

double fastjet::ClustSeqActAreaEG::area const PseudoJet jet  )  const [virtual]
 

returns the area of a jet

Reimplemented from fastjet::ClusterSequenceWithArea.

Definition at line 64 of file ClusterSequenceActiveAreaExplicitGhosts.cc.

References _areas.

Referenced by fastjet::ClusterSequenceActiveArea::_transfer_areas().

00064                                                            {
00065   return _areas[jet.cluster_hist_index()];
00066 }

PseudoJet fastjet::ClustSeqActAreaEG::area_4vector const PseudoJet jet  )  const [virtual]
 

returns a four vector corresponding to the sum (E-scheme) of the ghost four-vectors composing the jet area, normalised such that for a small contiguous area the p_t of the extended_area jet is equal to area of the jet.

Reimplemented from fastjet::ClusterSequenceWithArea.

Definition at line 78 of file ClusterSequenceActiveAreaExplicitGhosts.cc.

References _area_4vectors, and fastjet::PseudoJet::cluster_hist_index().

Referenced by fastjet::ClusterSequenceActiveArea::_transfer_areas().

00078                                                                       {
00079   return _area_4vectors[jet.cluster_hist_index()];
00080 }

bool fastjet::ClustSeqActAreaEG::is_pure_ghost int  history_index  )  const
 

true if the entry in the history index corresponds to a ghost; if hist_ix does not correspond to an actual particle (i.e.

hist_ix < 0), then the result is false.

Definition at line 89 of file ClusterSequenceActiveAreaExplicitGhosts.cc.

References _is_pure_ghost.

00090 {
00091   return hist_ix >= 0 ? _is_pure_ghost[hist_ix] : false;
00092 }

bool fastjet::ClustSeqActAreaEG::is_pure_ghost const PseudoJet jet  )  const
 

true if a jet is made exclusively of ghosts

Definition at line 83 of file ClusterSequenceActiveAreaExplicitGhosts.cc.

References _is_pure_ghost, and fastjet::PseudoJet::cluster_hist_index().

Referenced by fastjet::ClusterSequenceActiveArea::_transfer_areas(), and fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history().

00084 {
00085   return _is_pure_ghost[jet.cluster_hist_index()];
00086 }

unsigned int fastjet::ClusterSequenceActiveAreaExplicitGhosts::n_hard_particles  )  const [inline]
 

returns the number of hard particles (i.e. those supplied by the user).

Definition at line 160 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

References _initial_hard_n.

00160 {return _initial_hard_n;}

double fastjet::ClustSeqActAreaEG::total_area  )  const
 

returns the total area under study

Definition at line 71 of file ClusterSequenceActiveAreaExplicitGhosts.cc.

References _ghost_area, and _n_ghosts.

00071                                             {
00072   return _n_ghosts * _ghost_area;
00073 }


Member Data Documentation

std::vector<PseudoJet> fastjet::ClusterSequenceActiveAreaExplicitGhosts::_area_4vectors [private]
 

Definition at line 98 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

Referenced by _post_process(), and area_4vector().

std::vector<double> fastjet::ClusterSequenceActiveAreaExplicitGhosts::_areas [private]
 

Definition at line 97 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

Referenced by _post_process(), and area().

double fastjet::ClusterSequenceActiveAreaExplicitGhosts::_ghost_area [private]
 

Definition at line 95 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

Referenced by _add_ghosts(), _post_process(), and total_area().

unsigned int fastjet::ClusterSequenceActiveAreaExplicitGhosts::_initial_hard_n [private]
 

Definition at line 100 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

Referenced by _add_ghosts(), and n_hard_particles().

std::vector<bool> fastjet::ClusterSequenceActiveAreaExplicitGhosts::_is_pure_ghost [private]
 

Definition at line 96 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

Referenced by _add_ghosts(), _post_process(), and is_pure_ghost().

int fastjet::ClusterSequenceActiveAreaExplicitGhosts::_n_ghosts [private]
 

Definition at line 94 of file ClusterSequenceActiveAreaExplicitGhosts.hh.

Referenced by _add_ghosts(), and total_area().


The documentation for this class was generated from the following files:
Generated on Mon Apr 2 20:58:18 2007 for fastjet by  doxygen 1.4.2