FastJet  3.3.0
Filter.hh
1 #ifndef __FASTJET_TOOLS_FILTER_HH__
2 #define __FASTJET_TOOLS_FILTER_HH__
3 
4 //FJSTARTHEADER
5 // $Id: Filter.hh 3846 2015-03-08 08:35:43Z soyez $
6 //
7 // Copyright (c) 2005-2014, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
8 //
9 //----------------------------------------------------------------------
10 // This file is part of FastJet.
11 //
12 // FastJet is free software; you can redistribute it and/or modify
13 // it under the terms of the GNU General Public License as published by
14 // the Free Software Foundation; either version 2 of the License, or
15 // (at your option) any later version.
16 //
17 // The algorithms that underlie FastJet have required considerable
18 // development. They are described in the original FastJet paper,
19 // hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
20 // FastJet as part of work towards a scientific publication, please
21 // quote the version you use and include a citation to the manual and
22 // optionally also to hep-ph/0512210.
23 //
24 // FastJet is distributed in the hope that it will be useful,
25 // but WITHOUT ANY WARRANTY; without even the implied warranty of
26 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
27 // GNU General Public License for more details.
28 //
29 // You should have received a copy of the GNU General Public License
30 // along with FastJet. If not, see <http://www.gnu.org/licenses/>.
31 //----------------------------------------------------------------------
32 //FJENDHEADER
33 
34 #include <fastjet/ClusterSequence.hh>
35 #include <fastjet/Selector.hh>
36 #include <fastjet/CompositeJetStructure.hh> // to derive the FilterStructure from CompositeJetStructure
37 #include <fastjet/tools/Transformer.hh> // to derive Filter from Transformer
38 #include <iostream>
39 #include <string>
40 
41 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
42 
43 // fwd declarations
44 class Filter;
45 class FilterStructure;
46 
47 //----------------------------------------------------------------------
48 /// @ingroup tools_generic
49 /// \class Filter
50 /// Class that helps perform filtering (Butterworth, Davison, Rubin
51 /// and Salam, arXiv:0802.2470) and trimming (Krohn, Thaler and Wang,
52 /// arXiv:0912.1342) on jets, optionally in conjunction with
53 /// subtraction (Cacciari and Salam, arXiv:0707.1378).
54 ///
55 /// For example, to apply filtering that reclusters a jet's
56 /// constituents with the Cambridge/Aachen jet algorithm with R=0.3
57 /// and then selects the 3 hardest subjets, one can use the following
58 /// code:
59 /// \code
60 /// Filter filter(JetDefinition(cambridge_algorithm, 0.3), SelectorNHardest(3));
61 /// PseudoJet filtered_jet = filter(original_jet);
62 /// \endcode
63 ///
64 /// To obtain trimming, involving for example the selection of all
65 /// subjets carrying at least 3% of the original jet's pt, the
66 /// selector would be replaced by SelectorPtFractionMin(0.03).
67 ///
68 /// To additionally perform subtraction on the subjets prior to
69 /// selection, either include a 3rd argument specifying the background
70 /// density rho, or call the set_subtractor(...) member function. If
71 /// subtraction is requested, the original jet must be the result of a
72 /// clustering with active area with explicit ghosts support or a
73 /// merging of such pieces.
74 ///
75 /// The information on the subjets that were kept and rejected can be
76 /// obtained using:
77 /// \code
78 /// vector<PseudoJet> kept_subjets = filtered_jet.pieces();
79 /// vector<PseudoJet> rejected_subjets = filtered_jet.structure_of<Filter>().rejected();
80 /// \endcode
81 ///
82 /// \section impl Implementation Note
83 ///
84 /// If the original jet was defined with the Cambridge/Aachen
85 /// algorithm (or is made of pieces each of which comes from the C/A
86 /// alg) and the filtering definition is C/A, then the filter does not
87 /// rerun the C/A algorithm on the constituents, but instead makes use
88 /// of the existent C/A cluster sequence in the original jet. This
89 /// increases the speed of the filter.
90 ///
91 /// See also \subpage Example11 for a further usage example.
92 ///
93 /// Support for areas, reuse of C/A cluster sequences, etc.,
94 /// considerably complicates the implementation of Filter. For an
95 /// explanation of how a simpler filter might be coded, see the
96 /// "User-defined transformers" appendix of the manual.
97 class Filter : public Transformer{
98 public:
99  /// trivial ctor
100  /// Note: this is just for derived classes
101  /// a Filter initialised through this constructor will not work!
102  Filter() : _Rfiltfunc(0), _initialised(false){};
103 
104  /// define a filter that decomposes a jet into subjets using a
105  /// generic JetDefinition and then keeps only a subset of these
106  /// subjets according to a Selector. Optionally, each subjet may be
107  /// internally bakground-subtracted prior to selection.
108  ///
109  /// \param subjet_def the jet definition applied to obtain the subjets
110  /// \param selector the Selector applied to compute the kept subjets
111  /// \param rho if non-zero, backgruond-subtract each subjet befor selection
112  ///
113  /// Note: internal subtraction only applies on jets that are
114  /// obtained with a cluster sequence with area support and explicit
115  /// ghosts
116  Filter(JetDefinition subjet_def, Selector selector, double rho = 0.0) :
117  _subjet_def(subjet_def), _Rfiltfunc(0), _Rfilt(-1), _selector(selector), _rho(rho), _subtractor(0), _initialised(true) {}
118 
119  /// Same as the full constructor (see above) but just specifying the radius
120  /// By default, Cambridge-Aachen is used
121  /// If the jet (or all its pieces) is obtained with a non-default
122  /// recombiner, that one will be used
123  /// \param Rfilt the filtering radius
124  Filter(double Rfilt, Selector selector, double rho = 0.0) :
125  _Rfiltfunc(0), _Rfilt(Rfilt), _selector(selector), _rho(rho), _subtractor(0), _initialised(true) {
126  if (_Rfilt<0)
127  throw Error("Attempt to create a Filter with a negative filtering radius");
128  }
129 
130  /// Same as the full constructor (see above) but just specifying a
131  /// filtering radius that will depend on the jet being filtered
132  /// As for the previous case, Cambridge-Aachen is used
133  /// If the jet (or all its pieces) is obtained with a non-default
134  /// recombiner, that one will be used
135  /// \param Rfilt_func the filtering radius function of a PseudoJet
136  Filter(FunctionOfPseudoJet<double> *Rfilt_func, Selector selector, double rho = 0.0) :
137  _Rfiltfunc(Rfilt_func), _Rfilt(-1), _selector(selector), _rho(rho), _subtractor(0), _initialised(true) {}
138 
139  /// default dtor
140  virtual ~Filter(){};
141 
142  /// Set a subtractor that is applied to all individual subjets before
143  /// deciding which ones to keep. It takes precedence over a non-zero rho.
144  void set_subtractor(const FunctionOfPseudoJet<PseudoJet> * subtractor_in) {_subtractor = subtractor_in;}
145 
146  /// Set a subtractor that is applied to all individual subjets before
147  /// deciding which ones to keep. It takes precedence over a non-zero rho.
148  const FunctionOfPseudoJet<PseudoJet> * subtractor() const{ return _subtractor;}
149 
150  /// runs the filtering and sets kept and rejected to be the jets of interest
151  /// (with non-zero rho, they will have been subtracted).
152  ///
153  /// \param jet the jet that gets filtered
154  /// \return the filtered jet
155  virtual PseudoJet result(const PseudoJet & jet) const;
156 
157  /// class description
158  virtual std::string description() const;
159 
160  // the type of the associated structure
162 
163 private:
164  /// Sets filtered_elements to be all the subjets on which filtering will work.
165  /// It also sets the subjet_def to be used in joining things (the bit of
166  /// subjet def that is of interest for later is the recombiner).
167  ///
168  /// this returns true if teh optimisation trick for C/A reclustering has been used
169  bool _set_filtered_elements(const PseudoJet & jet,
170  std::vector<PseudoJet> & filtered_elements) const;
171 
172  /// gather the information about what is kept and rejected under the
173  /// form of a PseudoJet with a special ClusterSequenceInfo
174  ///
175  /// The last argument (ca_optimisation_used) should be true if the
176  /// optimisation trick for C/A reclustering has been used (in which
177  /// case some extra tests have to be run for non-explicit-ghost
178  /// areas)
179  PseudoJet _finalise(const PseudoJet & jet,
180  std::vector<PseudoJet> & kept,
181  std::vector<PseudoJet> & rejected,
182  bool ca_optimisation_used) const;
183 
184  bool _uses_subtraction() const {return (_subtractor || _rho != 0);}
185 
186  JetDefinition _subjet_def; ///< the jet definition to use to extract the subjets
187  FunctionOfPseudoJet<double> *_Rfiltfunc;
188  ///< a dynamic filtering radius function of the jet being filtered
189  double _Rfilt; ///< a constant specifying the subjet radius (with C/A)
190  Selector _selector; ///< the subjet selection criterium
191  double _rho; ///< the background density (used for subtraction when possible)
192  const FunctionOfPseudoJet<PseudoJet> * _subtractor; ///< for subtracting bkgd density from subjets
193 
194  bool _initialised; ///< true when the Filter has been properly intialised
195 };
196 
197 
198 
199 //----------------------------------------------------------------------
200 /// @ingroup tools_generic
201 /// \class FilterStructure
202 /// Class to contain structure information for a filtered jet.
204 public:
205  /// constructor from an original ClusterSequenceInfo
206  /// We just share the original ClusterSequenceWrapper and initialise
207  /// the rest
208  FilterStructure(const std::vector<PseudoJet> & pieces_in,
209  const JetDefinition::Recombiner *rec = 0)
210  : CompositeJetStructure(pieces_in, rec){}
211 
212  /// virtual dtor to allow further overloading
213  virtual ~FilterStructure(){}
214 
215  /// description
216  virtual std::string description() const { return "Filtered PseudoJet"; }
217 
218  //------------------------------------------------------------------
219  /// @name The filter-specific information
220  //------------------------------------------------------------------
221 
222 // /// returns the original jet (the first of the original jets
223 // /// if you filtered a collection of jets)
224 // const PseudoJet & original() const {return _original_jet;}
225 
226  /// returns the subjets that were not kept during the filtering procedure
227  /// (subtracted if the filter requests it, and valid in the original cs)
228  const std::vector<PseudoJet> & rejected() const {return _rejected;}
229 
230  friend class Filter; // allow the filter to change the protected/private members
231 
232 protected:
233 // PseudoJet _original_jet; ///< the original jet
234  std::vector<PseudoJet> _rejected; ///< the subjets rejected by the filter
235 };
236 
237 
238 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh
239 
240 #endif // __FASTJET_TOOLS_FILTER_HH__
Filter(JetDefinition subjet_def, Selector selector, double rho=0.0)
define a filter that decomposes a jet into subjets using a generic JetDefinition and then keeps only ...
Definition: Filter.hh:116
const FunctionOfPseudoJet< PseudoJet > * subtractor() const
Set a subtractor that is applied to all individual subjets before deciding which ones to keep...
Definition: Filter.hh:148
void set_subtractor(const FunctionOfPseudoJet< PseudoJet > *subtractor_in)
Set a subtractor that is applied to all individual subjets before deciding which ones to keep...
Definition: Filter.hh:144
Class that helps perform filtering (Butterworth, Davison, Rubin and Salam, arXiv:0802.2470) and trimming (Krohn, Thaler and Wang, arXiv:0912.1342) on jets, optionally in conjunction with subtraction (Cacciari and Salam, arXiv:0707.1378).
Definition: Filter.hh:97
Base (abstract) class for a jet transformer.
Definition: Transformer.hh:71
The structure for a jet made of pieces.
FilterStructure(const std::vector< PseudoJet > &pieces_in, const JetDefinition::Recombiner *rec=0)
constructor from an original ClusterSequenceInfo We just share the original ClusterSequenceWrapper an...
Definition: Filter.hh:208
Filter(double Rfilt, Selector selector, double rho=0.0)
Same as the full constructor (see above) but just specifying the radius By default, Cambridge-Aachen is used If the jet (or all its pieces) is obtained with a non-default recombiner, that one will be used.
Definition: Filter.hh:124
Filter()
trivial ctor Note: this is just for derived classes a Filter initialised through this constructor wil...
Definition: Filter.hh:102
base class corresponding to errors that can be thrown by FastJet
Definition: Error.hh:47
Class to contain structure information for a filtered jet.
Definition: Filter.hh:203
virtual std::string description() const
description
Definition: Filter.hh:216
const std::vector< PseudoJet > & rejected() const
returns the subjets that were not kept during the filtering procedure (subtracted if the filter reque...
Definition: Filter.hh:228
virtual ~FilterStructure()
virtual dtor to allow further overloading
Definition: Filter.hh:213
Class that encodes information about cuts and other selection criteria that can be applied to PseudoJ...
Definition: Selector.hh:149
Filter(FunctionOfPseudoJet< double > *Rfilt_func, Selector selector, double rho=0.0)
Same as the full constructor (see above) but just specifying a filtering radius that will depend on t...
Definition: Filter.hh:136
virtual ~Filter()
default dtor
Definition: Filter.hh:140
std::vector< PseudoJet > _rejected
the subjets rejected by the filter
Definition: Filter.hh:234
An abstract base class that will provide the recombination scheme facilities and/or allow a user to e...
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