Fork me on GitHub

source: git/external/fastjet/ClusterSequenceStructure.hh@ 8497ac6

ImprovedOutputFile Timing dual_readout llp
Last change on this file since 8497ac6 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: 10.7 KB
Line 
1//FJSTARTHEADER
2// $Id: ClusterSequenceStructure.hh 3433 2014-07-23 08:17:03Z 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
32#ifndef __FASTJET_CLUSTER_SEQUENCE_STRUCTURE_HH__
33#define __FASTJET_CLUSTER_SEQUENCE_STRUCTURE_HH__
34
35#include "fastjet/internal/base.hh"
36#include "fastjet/SharedPtr.hh"
37#include "fastjet/PseudoJetStructureBase.hh"
38
39#include <vector>
40
41FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
42
43/// @ingroup extra_info
44/// \class ClusterSequenceStructure
45///
46/// Contains any information related to the clustering that should be
47/// directly accessible to PseudoJet.
48///
49/// By default, this class implements basic access to the
50/// ClusterSequence related to a PseudoJet (like its constituents or
51/// its area). But it can be overloaded in order e.g. to give access
52/// to the jet substructure.
53///
54// Design question: Do we only put the methods that can be overloaded
55// or do we put everything a PJ can have access to? I think both cost
56// the same number of indirections. The first option limits the amount
57// of coding and maybe has a clearer structure. The second is more
58// consistent (everything related to the same thing is at the same
59// place) and gives better access for derived classes. We'll go for
60// the second option.
61class ClusterSequenceStructure : public PseudoJetStructureBase{
62public:
63 /// default ctor
64 ClusterSequenceStructure() : _associated_cs(NULL){}
65
66 /// ctor with initialisation to a given ClusterSequence
67 ///
68 /// In principle, this is reserved for initialisation by the parent
69 /// ClusterSequence
70 ClusterSequenceStructure(const ClusterSequence *cs){
71 set_associated_cs(cs);
72 };
73
74 /// default (virtual) dtor
75 virtual ~ClusterSequenceStructure();
76
77 /// description
78 virtual std::string description() const{ return "PseudoJet with an associated ClusterSequence"; }
79
80 //-------------------------------------------------------------
81 /// @name Direct access to the associated ClusterSequence object.
82 ///
83 /// Get access to the associated ClusterSequence (if any)
84 //\{
85 //-------------------------------------------------------------
86 /// returns true if there is an associated ClusterSequence
87 virtual bool has_associated_cluster_sequence() const{ return true;}
88
89 /// get a (const) pointer to the parent ClusterSequence (NULL if
90 /// inexistent)
91 virtual const ClusterSequence* associated_cluster_sequence() const;
92
93 /// returns true if there is a valid associated ClusterSequence
94 virtual bool has_valid_cluster_sequence() const;
95
96 /// if the jet has a valid associated cluster sequence then return a
97 /// pointer to it; otherwise throw an error
98 virtual const ClusterSequence * validated_cs() const;
99
100#ifndef __FJCORE__
101 /// if the jet has valid area information then return a pointer to
102 /// the associated ClusterSequenceAreaBase object; otherwise throw an error
103 virtual const ClusterSequenceAreaBase * validated_csab() const;
104#endif // __FJCORE__
105
106 /// set the associated csw
107 virtual void set_associated_cs(const ClusterSequence * new_cs){
108 _associated_cs = new_cs;
109 }
110 //\}
111
112 //-------------------------------------------------------------
113 /// @name Methods for access to information about jet structure
114 ///
115 /// These allow access to jet constituents, and other jet
116 /// subtructure information. They only work if the jet is associated
117 /// with a ClusterSequence.
118 //-------------------------------------------------------------
119 //\{
120
121 /// check if it has been recombined with another PseudoJet in which
122 /// case, return its partner through the argument. Otherwise,
123 /// 'partner' is set to 0.
124 ///
125 /// an Error is thrown if this PseudoJet has no currently valid
126 /// associated ClusterSequence
127 virtual bool has_partner(const PseudoJet &reference, PseudoJet &partner) const;
128
129 /// check if it has been recombined with another PseudoJet in which
130 /// case, return its child through the argument. Otherwise, 'child'
131 /// is set to 0.
132 ///
133 /// an Error is thrown if this PseudoJet has no currently valid
134 /// associated ClusterSequence
135 virtual bool has_child(const PseudoJet &reference, PseudoJet &child) const;
136
137 /// check if it is the product of a recombination, in which case
138 /// return the 2 parents through the 'parent1' and 'parent2'
139 /// arguments. Otherwise, set these to 0.
140 ///
141 /// an Error is thrown if this PseudoJet has no currently valid
142 /// associated ClusterSequence
143 virtual bool has_parents(const PseudoJet &reference, PseudoJet &parent1, PseudoJet &parent2) const;
144
145 /// check if the reference PseudoJet is contained in the second one
146 /// passed as argument.
147 ///
148 /// an Error is thrown if this PseudoJet has no currently valid
149 /// associated ClusterSequence
150 ///
151 /// false is returned if the 2 PseudoJet do not belong the same
152 /// ClusterSequence
153 virtual bool object_in_jet(const PseudoJet &reference, const PseudoJet &jet) const;
154
155 /// return true if the structure supports constituents.
156 ///
157 /// an Error is thrown if this PseudoJet has no currently valid
158 /// associated ClusterSequence
159 virtual bool has_constituents() const;
160
161 /// retrieve the constituents.
162 ///
163 /// an Error is thrown if this PseudoJet has no currently valid
164 /// associated ClusterSequence
165 virtual std::vector<PseudoJet> constituents(const PseudoJet &reference) const;
166
167
168 /// return true if the structure supports exclusive_subjets.
169 ///
170 /// an Error is thrown if this PseudoJet has no currently valid
171 /// associated ClusterSequence
172 virtual bool has_exclusive_subjets() const;
173
174 /// return a vector of all subjets of the current jet (in the sense
175 /// of the exclusive algorithm) that would be obtained when running
176 /// the algorithm with the given dcut.
177 ///
178 /// Time taken is O(m ln m), where m is the number of subjets that
179 /// are found. If m gets to be of order of the total number of
180 /// constituents in the jet, this could be substantially slower than
181 /// just getting that list of constituents.
182 ///
183 /// an Error is thrown if this PseudoJet has no currently valid
184 /// associated ClusterSequence
185 virtual std::vector<PseudoJet> exclusive_subjets(const PseudoJet &reference, const double & dcut) const;
186
187 /// return the size of exclusive_subjets(...); still n ln n with same
188 /// coefficient, but marginally more efficient than manually taking
189 /// exclusive_subjets.size()
190 ///
191 /// an Error is thrown if this PseudoJet has no currently valid
192 /// associated ClusterSequence
193 virtual int n_exclusive_subjets(const PseudoJet &reference, const double & dcut) const;
194
195 /// return the list of subjets obtained by unclustering the supplied
196 /// jet down to nsub subjets (or all constituents if there are fewer
197 /// than nsub).
198 ///
199 /// requires nsub ln nsub time
200 ///
201 /// an Error is thrown if this PseudoJet has no currently valid
202 /// associated ClusterSequence
203 virtual std::vector<PseudoJet> exclusive_subjets_up_to (const PseudoJet &reference, int nsub) const;
204
205 /// return the dij that was present in the merging nsub+1 -> nsub
206 /// subjets inside this jet.
207 ///
208 /// an Error is thrown if this PseudoJet has no currently valid
209 /// associated ClusterSequence
210 virtual double exclusive_subdmerge(const PseudoJet &reference, int nsub) const;
211
212 /// return the maximum dij that occurred in the whole event at the
213 /// stage that the nsub+1 -> nsub merge of subjets occurred inside
214 /// this jet.
215 ///
216 /// an Error is thrown if this PseudoJet has no currently valid
217 /// associated ClusterSequence
218 virtual double exclusive_subdmerge_max(const PseudoJet &reference, int nsub) const;
219
220
221 //-------------------------------------------------------------------
222 // information related to the pieces of the jet
223 //-------------------------------------------------------------------
224 /// by convention, a jet associated with a ClusterSequence will have
225 /// its parents as pieces
226 virtual bool has_pieces(const PseudoJet &reference) const;
227
228 /// by convention, a jet associated with a ClusterSequence will have
229 /// its parents as pieces
230 ///
231 /// if it has no parents, then there will only be a single piece:
232 /// itself
233 ///
234 /// Note that to answer that question, we need to access the cluster
235 /// sequence. If the cluster sequence has gone out of scope, an
236 /// error will be thrown
237 virtual std::vector<PseudoJet> pieces(const PseudoJet &reference) const;
238
239
240 // the following ones require a computation of the area in the
241 // parent ClusterSequence (See ClusterSequenceAreaBase for details)
242 //------------------------------------------------------------------
243#ifndef __FJCORE__
244
245 /// check if it has a defined area
246 virtual bool has_area() const;
247
248 /// return the jet (scalar) area.
249 /// throws an Error if there is no support for area in the parent CS
250 virtual double area(const PseudoJet &reference) const;
251
252 /// return the error (uncertainty) associated with the determination
253 /// of the area of this jet.
254 /// throws an Error if there is no support for area in the parent CS
255 virtual double area_error(const PseudoJet &reference) const;
256
257 /// return the jet 4-vector area.
258 /// throws an Error if there is no support for area in the parent CS
259 virtual PseudoJet area_4vector(const PseudoJet &reference) const;
260
261 /// true if this jet is made exclusively of ghosts.
262 /// throws an Error if there is no support for area in the parent CS
263 virtual bool is_pure_ghost(const PseudoJet &reference) const;
264
265#endif // __FJCORE__
266 //\} --- end of jet structure -------------------------------------
267
268protected:
269 const ClusterSequence *_associated_cs;
270};
271
272FASTJET_END_NAMESPACE
273
274#endif // __FASTJET_CLUSTER_SEQUENCE_STRUCTURE_HH__
Note: See TracBrowser for help on using the repository browser.