00001 //STARTHEADER 00002 // $Id: ClusterSequenceActiveAreaExplicitGhosts.cc 482 2007-02-28 17:07:24Z cacciari $ 00003 // 00004 // Copyright (c) 2005-2006, Matteo Cacciari and Gavin Salam 00005 // 00006 //---------------------------------------------------------------------- 00007 // This file is part of FastJet. 00008 // 00009 // FastJet is free software; you can redistribute it and/or modify 00010 // it under the terms of the GNU General Public License as published by 00011 // the Free Software Foundation; either version 2 of the License, or 00012 // (at your option) any later version. 00013 // 00014 // The algorithms that underlie FastJet have required considerable 00015 // development and are described in hep-ph/0512210. If you use 00016 // FastJet as part of work towards a scientific publication, please 00017 // include a citation to the FastJet paper. 00018 // 00019 // FastJet is distributed in the hope that it will be useful, 00020 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00021 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00022 // GNU General Public License for more details. 00023 // 00024 // You should have received a copy of the GNU General Public License 00025 // along with FastJet; if not, write to the Free Software 00026 // Foundation, Inc.: 00027 // 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00028 //---------------------------------------------------------------------- 00029 //ENDHEADER 00030 00031 #include "fastjet/ClusterSequenceActiveAreaExplicitGhosts.hh" 00032 00033 using namespace std; 00034 00035 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh 00036 00037 // save some typing 00038 typedef ClusterSequenceActiveAreaExplicitGhosts ClustSeqActAreaEG; 00039 00040 00041 00042 //---------------------------------------------------------------------- 00044 void ClustSeqActAreaEG::_add_ghosts ( 00045 const ActiveAreaSpec & area_spec) { 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 } 00060 00061 00062 //---------------------------------------------------------------------- 00063 // return the area of a jet 00064 double ClustSeqActAreaEG::area (const PseudoJet & jet) const { 00065 return _areas[jet.cluster_hist_index()]; 00066 } 00067 00068 00069 //---------------------------------------------------------------------- 00070 // return the total area 00071 double ClustSeqActAreaEG::total_area () const { 00072 return _n_ghosts * _ghost_area; 00073 } 00074 00075 00076 //---------------------------------------------------------------------- 00077 // return the extended area of a jet 00078 PseudoJet ClustSeqActAreaEG::area_4vector (const PseudoJet & jet) const { 00079 return _area_4vectors[jet.cluster_hist_index()]; 00080 } 00081 00082 //---------------------------------------------------------------------- 00083 bool ClustSeqActAreaEG::is_pure_ghost(const PseudoJet & jet) const 00084 { 00085 return _is_pure_ghost[jet.cluster_hist_index()]; 00086 } 00087 00088 //---------------------------------------------------------------------- 00089 bool ClustSeqActAreaEG::is_pure_ghost(int hist_ix) const 00090 { 00091 return hist_ix >= 0 ? _is_pure_ghost[hist_ix] : false; 00092 } 00093 00094 00095 //====================================================================== 00096 // sort out the areas 00097 void ClustSeqActAreaEG::_post_process() { 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 } 00142 00143 FASTJET_END_NAMESPACE 00144