FastJet  3.3.1
ClusterSequenceActiveArea.hh
1 //FJSTARTHEADER
2 // $Id: ClusterSequenceActiveArea.hh 4354 2018-04-22 07:12:37Z salam $
3 //
4 // Copyright (c) 2005-2018, 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 
48 FASTJET_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
64 public:
65 
66  /// default constructor
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 FASTJET_OVERRIDE {
77  return _average_area[jet.cluster_hist_index()];};
78  virtual double area_error (const PseudoJet & jet) const FASTJET_OVERRIDE {
79  return _average_area2[jet.cluster_hist_index()];};
80 
81  virtual PseudoJet area_4vector (const PseudoJet & jet) const FASTJET_OVERRIDE {
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 FASTJET_OVERRIDE;
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 FASTJET_OVERRIDE;
117 
118 protected:
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,
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 
155 private:
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 
207 template<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 
223 FASTJET_END_NAMESPACE
224 
225 #endif // __FASTJET_CLUSTERSEQUENCEACTIVEAREA_HH__
mean_pt_strategies
enum providing a variety of tentative strategies for estimating the background (e.g.
Like ClusterSequence with computation of the active jet area.
virtual double area(const PseudoJet &jet) const
return the area associated with the given jet; this base class returns 0.
bool has_dangerous_particles() const
returns true if there are any particles whose transverse momentum if so low that there&#39;s a risk of th...
base class that sets interface for extensions of ClusterSequence that provide information about the a...
Like ClusterSequence with computation of the active jet area with the addition of explicit ghosts...
virtual PseudoJet area_4vector(const PseudoJet &jet) const
return a PseudoJet whose 4-vector is defined by the following integral
std::valarray< double > _average_area
child classes benefit from having these at their disposal
Class that encodes information about cuts and other selection criteria that can be applied to PseudoJ...
Definition: Selector.hh:149
virtual double area_error(const PseudoJet &jet) const
return the error (uncertainty) associated with the determination of the area of this jet; this base c...
Parameters to configure the computation of jet areas using ghosts.
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
Definition: PseudoJet.hh:67
class that is intended to hold a full definition of the jet clusterer