1 | #ifndef __FASTJET_TOOLS_FILTER_HH__
|
---|
2 | #define __FASTJET_TOOLS_FILTER_HH__
|
---|
3 |
|
---|
4 | //FJSTARTHEADER
|
---|
5 | // $Id: Filter.hh 4442 2020-05-05 07:50:11Z soyez $
|
---|
6 | //
|
---|
7 | // Copyright (c) 2005-2020, 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
|
---|
161 | typedef FilterStructure StructureType;
|
---|
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.
|
---|
203 | class FilterStructure : public CompositeJetStructure {
|
---|
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__
|
---|