Fork me on GitHub

source: git/external/fastjet/ClusterSequenceActiveArea.hh@ ff37d75

ImprovedOutputFile Timing dual_readout llp
Last change on this file since ff37d75 was 35cdc46, checked in by Pavel Demin <demin@…>, 10 years ago

upgrade FastJet to version 3.1.0-beta.1, upgrade Nsubjettiness to version 2.1.0, add SoftKiller version 1.0.0

  • Property mode set to 100644
File size: 8.7 KB
Line 
1//FJSTARTHEADER
2// $Id: ClusterSequenceActiveArea.hh 3619 2014-08-13 14:17:19Z salam $
3//
4// Copyright (c) 2005-2014, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
5//
6//----------------------------------------------------------------------
7// This file is part of FastJet.
8//
9// FastJet is free software; you can redistribute it and/or modify
10// it under the terms of the GNU General Public License as published by
11// the Free Software Foundation; either version 2 of the License, or
12// (at your option) any later version.
13//
14// The algorithms that underlie FastJet have required considerable
15// development. They are described in the original FastJet paper,
16// hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
17// FastJet as part of work towards a scientific publication, please
18// quote the version you use and include a citation to the manual and
19// optionally also to hep-ph/0512210.
20//
21// FastJet is distributed in the hope that it will be useful,
22// but WITHOUT ANY WARRANTY; without even the implied warranty of
23// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24// GNU General Public License for more details.
25//
26// You should have received a copy of the GNU General Public License
27// along with FastJet. If not, see <http://www.gnu.org/licenses/>.
28//----------------------------------------------------------------------
29//FJENDHEADER
30
31#ifndef __FASTJET_CLUSTERSEQUENCEACTIVEAREA_HH__
32#define __FASTJET_CLUSTERSEQUENCEACTIVEAREA_HH__
33
34
35#include "fastjet/PseudoJet.hh"
36#include "fastjet/ClusterSequenceAreaBase.hh"
37#include "fastjet/ClusterSequenceActiveAreaExplicitGhosts.hh"
38#include<iostream>
39#include<vector>
40
41//------------ backwards compatibility with version 2.1 -------------
42// for backwards compatibility make ActiveAreaSpec name available
43//#include "fastjet/ActiveAreaSpec.hh"
44//#include "fastjet/ClusterSequenceWithArea.hh"
45//--------------------------------------------------------------------
46
47
48FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
49
50//using namespace std;
51
52/// @ingroup sec_area_classes
53/// \class ClusterSequenceActiveArea
54/// Like ClusterSequence with computation of the active jet area
55///
56/// Class that behaves essentially like ClusterSequence except
57/// that it also provides access to the area of a jet (which
58/// will be a random quantity... Figure out what to do about seeds
59/// later...)
60///
61/// This class should not be used directly. Rather use
62/// ClusterSequenceArea with the appropriate AreaDefinition
63class ClusterSequenceActiveArea : public ClusterSequenceAreaBase {
64public:
65
66 /// default constructor
67 ClusterSequenceActiveArea() {}
68
69 /// constructor based on JetDefinition and GhostedAreaSpec
70 template<class L> ClusterSequenceActiveArea
71 (const std::vector<L> & pseudojets,
72 const JetDefinition & jet_def_in,
73 const GhostedAreaSpec & ghost_spec,
74 const bool & writeout_combinations = false) ;
75
76 virtual double area (const PseudoJet & jet) const {
77 return _average_area[jet.cluster_hist_index()];};
78 virtual double area_error (const PseudoJet & jet) const {
79 return _average_area2[jet.cluster_hist_index()];};
80
81 virtual PseudoJet area_4vector (const PseudoJet & jet) const {
82 return _average_area_4vector[jet.cluster_hist_index()];};
83
84 /// enum providing a variety of tentative strategies for estimating
85 /// the background (e.g. non-jet) activity in a highly populated event; the
86 /// one that has been most extensively tested is median.
87 ///
88 /// These strategies are OBSOLETE and deprecated (see comment
89 /// for pt_per_unit_area).
90 enum mean_pt_strategies{median=0, non_ghost_median, pttot_over_areatot,
91 pttot_over_areatot_cut, mean_ratio_cut, play,
92 median_4vector};
93
94 /// return the transverse momentum per unit area according to one
95 /// of the above strategies; for some strategies (those with "cut"
96 /// in their name) the parameter "range" allows one to exclude a
97 /// subset of the jets for the background estimation, those that
98 /// have pt/area > median(pt/area)*range.
99 ///
100 /// NB: This call is OBSOLETE and deprecated; use a
101 /// JetMedianBackgroundEstimator or GridMedianBackgroundEstimator
102 /// instead.
103 double pt_per_unit_area(mean_pt_strategies strat=median,
104 double range=2.0 ) const;
105
106 /// rewrite the empty area from the parent class, so as to use
107 /// all info at our disposal
108 /// return the total area, corresponding to a given Selector, that
109 /// consists of ghost jets or unclustered ghosts
110 ///
111 /// The selector passed as an argument needs to apply jet by jet.
112 virtual double empty_area(const Selector & selector) const;
113
114 /// return the true number of empty jets (replaces
115 /// ClusterSequenceAreaBase::n_empty_jets(...))
116 virtual double n_empty_jets(const Selector & selector) const;
117
118protected:
119 void _resize_and_zero_AA ();
120 void _initialise_AA(const JetDefinition & jet_def,
121 const GhostedAreaSpec & ghost_spec,
122 const bool & writeout_combinations,
123 bool & continue_running);
124
125 void _run_AA(const GhostedAreaSpec & ghost_spec);
126
127 void _postprocess_AA(const GhostedAreaSpec & ghost_spec);
128
129 /// does the initialisation and running specific to the active
130 /// areas class
131 void _initialise_and_run_AA (const JetDefinition & jet_def,
132 const GhostedAreaSpec & ghost_spec,
133 const bool & writeout_combinations = false);
134
135 /// transfer the history (and jet-momenta) from clust_seq to our
136 /// own internal structure while removing ghosts
137 void _transfer_ghost_free_history(
138 const ClusterSequenceActiveAreaExplicitGhosts & clust_seq);
139
140
141 /// transfer areas from the ClusterSequenceActiveAreaExplicitGhosts
142 /// object into our internal area bookkeeping...
143 void _transfer_areas(const std::vector<int> & unique_hist_order,
144 const ClusterSequenceActiveAreaExplicitGhosts & );
145
146 /// child classes benefit from having these at their disposal
147 std::valarray<double> _average_area, _average_area2;
148 std::valarray<PseudoJet> _average_area_4vector;
149
150 /// returns true if there are any particles whose transverse momentum
151 /// if so low that there's a risk of the ghosts having modified the
152 /// clustering sequence
153 bool has_dangerous_particles() const {return _has_dangerous_particles;}
154
155private:
156
157
158 double _non_jet_area, _non_jet_area2, _non_jet_number;
159
160 double _maxrap_for_area; // max rap where we put ghosts
161 double _safe_rap_for_area; // max rap where we trust jet areas
162
163 bool _has_dangerous_particles;
164
165
166 /// routine for extracting the tree in an order that will be independent
167 /// of any degeneracies in the recombination sequence that don't
168 /// affect the composition of the final jets
169 void _extract_tree(std::vector<int> &) const;
170 /// do the part of the extraction associated with pos, working
171 /// through its children and their parents
172 void _extract_tree_children(int pos, std::valarray<bool> &, const std::valarray<int> &, std::vector<int> &) const;
173 /// do the part of the extraction associated with the parents of pos.
174 void _extract_tree_parents (int pos, std::valarray<bool> &, const std::valarray<int> &, std::vector<int> &) const;
175
176 /// check if two jets have the same momentum to within the
177 /// tolerance (and if pt's are not the same we're forgiving and
178 /// look to see if the energy is the same)
179 void _throw_unless_jets_have_same_perp_or_E(const PseudoJet & jet,
180 const PseudoJet & refjet,
181 double tolerance,
182 const ClusterSequenceActiveAreaExplicitGhosts & jets_ghosted_seq
183 ) const;
184
185 /// since we are playing nasty games with seeds, we should warn
186 /// the user a few times
187 //static int _n_seed_warnings;
188 //const static int _max_seed_warnings = 10;
189
190 // record the number of repeats
191 int _ghost_spec_repeat;
192
193 /// a class for our internal storage of ghost jets
194 class GhostJet : public PseudoJet {
195 public:
196 GhostJet(const PseudoJet & j, double a) : PseudoJet(j), area(a){}
197 double area;
198 };
199
200 std::vector<GhostJet> _ghost_jets;
201 std::vector<GhostJet> _unclustered_ghosts;
202};
203
204
205
206
207template<class L> ClusterSequenceActiveArea::ClusterSequenceActiveArea
208(const std::vector<L> & pseudojets,
209 const JetDefinition & jet_def_in,
210 const GhostedAreaSpec & ghost_spec,
211 const bool & writeout_combinations) {
212
213 // transfer the initial jets (type L) into our own array
214 _transfer_input_jets(pseudojets);
215
216 // run the clustering for active areas
217 _initialise_and_run_AA(jet_def_in, ghost_spec, writeout_combinations);
218
219}
220
221
222
223FASTJET_END_NAMESPACE
224
225#endif // __FASTJET_CLUSTERSEQUENCEACTIVEAREA_HH__
Note: See TracBrowser for help on using the repository browser.