29 #include "fastjet/ClusterSequenceActiveAreaExplicitGhosts.hh" 34 FASTJET_BEGIN_NAMESPACE
37 typedef ClusterSequenceActiveAreaExplicitGhosts ClustSeqActAreaEG;
40 LimitedWarning ClustSeqActAreaEG::_warnings;
44 void ClustSeqActAreaEG::_add_ghosts (
45 const GhostedAreaSpec & ghost_spec) {
48 ghost_spec.add_ghosts(_jets);
51 for (
unsigned i = _initial_hard_n; i < _jets.size(); i++) {
53 _is_pure_ghost.push_back(
true);
57 _ghost_area = ghost_spec.actual_ghost_area();
58 _n_ghosts = ghost_spec.n_ghosts();
64 double ClustSeqActAreaEG::area (
const PseudoJet & jet)
const {
71 double ClustSeqActAreaEG::total_area ()
const {
72 return _n_ghosts * _ghost_area;
83 bool ClustSeqActAreaEG::is_pure_ghost(
const PseudoJet & jet)
const 89 bool ClustSeqActAreaEG::is_pure_ghost(
int hist_ix)
const 91 return hist_ix >= 0 ? _is_pure_ghost[hist_ix] :
false;
95 double ClustSeqActAreaEG::empty_area(
const Selector & selector)
const {
98 throw Error(
"ClusterSequenceActiveAreaExplicitGhosts: empty area can only be computed from selectors applying jet by jet");
101 vector<PseudoJet> unclust = unclustered_particles();
102 double area_local = 0.0;
103 for (
unsigned iu = 0; iu < unclust.size(); iu++) {
104 if (is_pure_ghost(unclust[iu]) && selector.
pass(unclust[iu])) {
105 area_local += _ghost_area;
113 void ClustSeqActAreaEG::_post_process() {
117 _max_ghost_perp2 = 0.0;
118 for (
int i = 0; i < _initial_n; i++) {
119 if (_is_pure_ghost[i] && _jets[i].perp2() > _max_ghost_perp2)
120 _max_ghost_perp2 = _jets[i].perp2();
124 double danger_ratio = numeric_limits<double>::epsilon();
125 danger_ratio = danger_ratio * danger_ratio;
126 _has_dangerous_particles =
false;
127 for (
int i = 0; i < _initial_n; i++) {
128 if (!_is_pure_ghost[i] &&
129 danger_ratio * _jets[i].perp2() <= _max_ghost_perp2) {
130 _has_dangerous_particles =
true;
135 if (_has_dangerous_particles) _warnings.warn(
"ClusterSequenceActiveAreaExplicitGhosts: \n ghosts not sufficiently soft wrt some of the input particles\n a common cause is (unphysical?) input particles with pt=0 but finite rapidity");
138 _areas.resize(_history.size());
139 _area_4vectors.resize(_history.size());
140 _is_pure_ghost.resize(_history.size());
158 for (
int i = 0; i < _initial_n; i++) {
159 if (_is_pure_ghost[i]) {
160 _areas[i] = _ghost_area;
169 _area_4vectors[i].reset_momentum(_jets[i]);
170 _area_4vectors[i] *= (_ghost_area/_jets[i].perp());
173 _area_4vectors[i] =
PseudoJet(0.0,0.0,0.0,0.0);
180 for (
unsigned i = _initial_n; i < _history.size(); i++) {
181 if (_history[i].parent2 == BeamJet) {
182 _is_pure_ghost[i] = _is_pure_ghost[_history[i].parent1];
183 _areas[i] = _areas[_history[i].parent1];
184 _area_4vectors[i] = _area_4vectors[_history[i].parent1];
186 _is_pure_ghost[i] = _is_pure_ghost[_history[i].parent1] &&
187 _is_pure_ghost[_history[i].parent2] ;
188 _areas[i] = _areas[_history[i].parent1] +
189 _areas[_history[i].parent2] ;
190 _jet_def.recombiner()->recombine(_area_4vectors[_history[i].parent1],
191 _area_4vectors[_history[i].parent2],
201 FASTJET_END_NAMESPACE
bool applies_jet_by_jet() const
returns true if this can be applied jet by jet
bool pass(const PseudoJet &jet) const
return true if the jet passes the selection
base class corresponding to errors that can be thrown by FastJet
int cluster_hist_index() const
return the cluster_hist_index, intended to be used by clustering routines.
Class that encodes information about cuts and other selection criteria that can be applied to PseudoJ...
Class to contain pseudojets, including minimal information of use to jet-clustering routines...