FastJet 3.4.1
PseudoJet.hh
1//FJSTARTHEADER
2// $Id$
3//
4// Copyright (c) 2005-2023, 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_PSEUDOJET_HH__
33#define __FASTJET_PSEUDOJET_HH__
34
35#include<valarray>
36#include<vector>
37#include<cassert>
38#include<cmath>
39#include<iostream>
40#include "fastjet/config.h"
41#include "fastjet/internal/numconsts.hh"
42#include "fastjet/internal/IsBase.hh"
43#include "fastjet/SharedPtr.hh"
44#include "fastjet/Error.hh"
45#include "fastjet/PseudoJetStructureBase.hh"
46
47FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
48
49//using namespace std;
50
51/// Used to protect against parton-level events where pt can be zero
52/// for some partons, giving rapidity=infinity. KtJet fails in those cases.
53const double MaxRap = 1e5;
54
55/// default value for phi, meaning it (and rapidity) have yet to be calculated)
56const double pseudojet_invalid_phi = -100.0;
57const double pseudojet_invalid_rap = -1e200;
58
59#ifndef __FJCORE__
60// forward definition
62#endif // __FJCORE__
63
64/// @ingroup basic_classes
65/// \class PseudoJet
66/// Class to contain pseudojets, including minimal information of use to
67/// jet-clustering routines.
68class PseudoJet {
69
70 public:
71 //----------------------------------------------------------------------
72 /// @name Constructors and destructor
73 //\{
74 /// default constructor, which as of FJ3.0 provides an object for
75 /// which all operations are now valid and which has zero momentum
76 ///
77 // (cf. this is actually OK from a timing point of view and in some
78 // cases better than just having the default constructor for the
79 // internal shared pointer: see PJtiming.cc and the notes therein)
80 //
81 // note: no reset of shared pointers needed
82 PseudoJet() : _px(0), _py(0), _pz(0), _E(0) {_finish_init(); _reset_indices();}
83 //PseudoJet() : _px(0), _py(0), _pz(0), _E(0), _phi(pseudojet_invalid_phi), _rap(pseudojet_invalid_rap), _kt2(0) {_reset_indices();}
84 /// construct a pseudojet from explicit components
85 PseudoJet(const double px, const double py, const double pz, const double E);
86
87 /// constructor from any object that has px,py,pz,E = some_four_vector[0--3],
88 #ifndef SWIG
89 template <class L> PseudoJet(const L & some_four_vector);
90 #endif
91
92 // Constructor that performs minimal initialisation (only that of
93 // the shared pointers), of use in certain speed-critical contexts
94 //
95 // NB: "dummy" is commented to avoid unused-variable compiler warnings
96 PseudoJet(bool /* dummy */) {}
97
98#ifdef FASTJET_HAVE_THREAD_SAFETY
99 PseudoJet(const PseudoJet &other){ (*this)=other; }
100 PseudoJet& operator=(const PseudoJet& other);
101#endif
102
103 /// default (virtual) destructor
104 virtual ~PseudoJet(){}
105//std::shared_ptr: #ifdef FASTJET_HAVE_THREAD_SAFETY
106//std::shared_ptr: ;
107//std::shared_ptr: #else
108//std::shared_ptr: {}
109//std::shared_ptr: #endif
110 //\} ---- end of constructors and destructors --------------------------
111
112 //----------------------------------------------------------------------
113 /// @name Kinematic access functions
114 //\{
115 //----------------------------------------------------------------------
116 inline double E() const {return _E;}
117 inline double e() const {return _E;} // like CLHEP
118 inline double px() const {return _px;}
119 inline double py() const {return _py;}
120 inline double pz() const {return _pz;}
121
122 /// returns phi (in the range 0..2pi)
123 inline double phi() const {return phi_02pi();}
124
125 /// returns phi in the range -pi..pi
126 inline double phi_std() const {
127 _ensure_valid_rap_phi();
128 return _phi > pi ? _phi-twopi : _phi;}
129
130 /// returns phi in the range 0..2pi
131 inline double phi_02pi() const {
132 _ensure_valid_rap_phi();
133 return _phi;
134 }
135
136 /// returns the rapidity or some large value when the rapidity
137 /// is infinite
138 inline double rap() const {
139 _ensure_valid_rap_phi();
140 return _rap;
141 }
142
143 /// the same as rap()
144 inline double rapidity() const {return rap();} // like CLHEP
145
146 /// returns the pseudo-rapidity or some large value when the
147 /// rapidity is infinite
148 double pseudorapidity() const;
149 double eta() const {return pseudorapidity();}
150
151 /// returns the squared transverse momentum
152 inline double pt2() const {return _kt2;}
153 /// returns the scalar transverse momentum
154 inline double pt() const {return sqrt(_kt2);}
155 /// returns the squared transverse momentum
156 inline double perp2() const {return _kt2;} // like CLHEP
157 /// returns the scalar transverse momentum
158 inline double perp() const {return sqrt(_kt2);} // like CLHEP
159 /// returns the squared transverse momentum
160 inline double kt2() const {return _kt2;} // for bkwds compatibility
161
162 /// returns the squared invariant mass // like CLHEP
163 inline double m2() const {return (_E+_pz)*(_E-_pz)-_kt2;}
164 /// returns the invariant mass
165 /// (If m2() is negative then -sqrt(-m2()) is returned, as in CLHEP)
166 inline double m() const;
167
168 /// returns the squared transverse mass = kt^2+m^2
169 inline double mperp2() const {return (_E+_pz)*(_E-_pz);}
170 /// returns the transverse mass = sqrt(kt^2+m^2)
171 inline double mperp() const {return sqrt(std::abs(mperp2()));}
172 /// returns the squared transverse mass = kt^2+m^2
173 inline double mt2() const {return (_E+_pz)*(_E-_pz);}
174 /// returns the transverse mass = sqrt(kt^2+m^2)
175 inline double mt() const {return sqrt(std::abs(mperp2()));}
176
177 /// return the squared 3-vector modulus = px^2+py^2+pz^2
178 inline double modp2() const {return _kt2+_pz*_pz;}
179 /// return the 3-vector modulus = sqrt(px^2+py^2+pz^2)
180 inline double modp() const {return sqrt(_kt2+_pz*_pz);}
181
182 /// return the transverse energy
183 inline double Et() const {return (_kt2==0) ? 0.0 : _E/sqrt(1.0+_pz*_pz/_kt2);}
184 /// return the transverse energy squared
185 inline double Et2() const {return (_kt2==0) ? 0.0 : _E*_E/(1.0+_pz*_pz/_kt2);}
186
187 /// cos of the polar angle
188 /// should we have: min(1.0,max(-1.0,_pz/sqrt(modp2())));
189 inline double cos_theta() const {
190 return std::min(1.0, std::max(-1.0, _pz/sqrt(modp2())));
191 }
192 /// polar angle
193 inline double theta() const { return acos(cos_theta()); }
194
195 /// returns component i, where X==0, Y==1, Z==2, E==3
196 double operator () (int i) const ;
197 /// returns component i, where X==0, Y==1, Z==2, E==3
198 inline double operator [] (int i) const { return (*this)(i); }; // this too
199
200
201
202 /// returns kt distance (R=1) between this jet and another
203 double kt_distance(const PseudoJet & other) const;
204
205 /// returns squared cylinder (rap-phi) distance between this jet and another
206 double plain_distance(const PseudoJet & other) const;
207 /// returns squared cylinder (rap-phi) distance between this jet and
208 /// another
209 inline double squared_distance(const PseudoJet & other) const {
210 return plain_distance(other);}
211
212 /// return the cylinder (rap-phi) distance between this jet and another,
213 /// \f$\Delta_R = \sqrt{\Delta y^2 + \Delta \phi^2}\f$.
214 inline double delta_R(const PseudoJet & other) const {
215 return sqrt(squared_distance(other));
216 }
217
218 /// returns other.phi() - this.phi(), constrained to be in
219 /// range -pi .. pi
220 double delta_phi_to(const PseudoJet & other) const;
221
222 //// this seemed to compile except if it was used
223 //friend inline double
224 // kt_distance(const PseudoJet & jet1, const PseudoJet & jet2) {
225 // return jet1.kt_distance(jet2);}
226
227 /// returns distance between this jet and the beam
228 inline double beam_distance() const {return _kt2;}
229
230 /// return a valarray containing the four-momentum (components 0-2
231 /// are 3-mom, component 3 is energy).
232 std::valarray<double> four_mom() const;
233
234 //\} ------- end of kinematic access functions
235
236 // taken from CLHEP
237 enum { X=0, Y=1, Z=2, T=3, NUM_COORDINATES=4, SIZE=NUM_COORDINATES };
238
239
240 //----------------------------------------------------------------------
241 /// @name Kinematic modification functions
242 //\{
243 //----------------------------------------------------------------------
244 /// transform this jet (given in the rest frame of prest) into a jet
245 /// in the lab frame
246 PseudoJet & boost(const PseudoJet & prest);
247 /// transform this jet (given in lab) into a jet in the rest
248 /// frame of prest
249 PseudoJet & unboost(const PseudoJet & prest);
250
251 PseudoJet & operator*=(double);
252 PseudoJet & operator/=(double);
253 PseudoJet & operator+=(const PseudoJet &);
254 PseudoJet & operator-=(const PseudoJet &);
255
256//std::shared_ptr: #ifdef FASTJET_HAVE_THREAD_SAFETY
257//std::shared_ptr: /// overload the assignment through the = operator
258//std::shared_ptr: ///
259//std::shared_ptr: /// this is needed to make sure that release_from_cs is called
260//std::shared_ptr: /// before copying the structure shared pointer!
261//std::shared_ptr: void operator=(const PseudoJet &other){
262//std::shared_ptr: _release_jet_from_cs();
263//std::shared_ptr: _structure = other._structure;
264//std::shared_ptr: _user_info = other._user_info;
265//std::shared_ptr:
266//std::shared_ptr: _px = other._px;
267//std::shared_ptr: _py = other._py;
268//std::shared_ptr: _pz = other._pz;
269//std::shared_ptr: _E = other._E;
270//std::shared_ptr:
271//std::shared_ptr: _phi = other._phi;
272//std::shared_ptr: _rap = other._rap;
273//std::shared_ptr: _kt2 = other._kt2;
274//std::shared_ptr:
275//std::shared_ptr: _cluster_hist_index = other._cluster_hist_index;
276//std::shared_ptr: _user_index = other._user_index;
277//std::shared_ptr: }
278//std::shared_ptr:
279//std::shared_ptr: /// force resetting the structure pointer to an empty structure
280//std::shared_ptr: ///
281//std::shared_ptr: /// THIS IS STRICTLY MEANT FOR INTERNAL USAGE IN SELF-DELETING
282//std::shared_ptr: /// ClusterSequences. IT SHOULD NOT BE USED BY END-USERS.
283//std::shared_ptr: void force_reset_structure(){
284//std::shared_ptr: _structure.reset();
285//std::shared_ptr: }
286//std::shared_ptr: #endif
287
288 /// reset the 4-momentum according to the supplied components and
289 /// put the user and history indices back to their default values
290 inline void reset(double px, double py, double pz, double E);
291
292 /// reset the PseudoJet to be equal to psjet (including its
293 /// indices); NB if the argument is derived from a PseudoJet then
294 /// the "reset" used will be the templated version
295 ///
296 /// Note: this is included on top of the templated version because
297 /// PseudoJet is not "derived" from PseudoJet, so the templated
298 /// reset would not handle this case properly.
299 inline void reset(const PseudoJet & psjet) {
300 (*this) = psjet;
301 }
302
303 /// reset the 4-momentum according to the supplied generic 4-vector
304 /// (accessible via indexing, [0]==px,...[3]==E) and put the user
305 /// and history indices back to their default values.
306#ifndef SWIG
307 template <class L> inline void reset(const L & some_four_vector) {
308 // check if some_four_vector can be cast to a PseudoJet
309 //
310 // Note that a regular dynamic_cast would not work here because
311 // there is no guarantee that L is polymorphic. We use a more
312 // complex construct here that works also in such a case. As for
313 // dynamic_cast, NULL is returned if L is not derived from
314 // PseudoJet
315 //
316 // Note the explicit request for fastjet::cast_if_derived; when
317 // combining fastjet and fjcore, this avoids ambiguity in which of
318 // the two cast_if_derived calls to use.
319 const PseudoJet * pj = fastjet::cast_if_derived<const PseudoJet>(&some_four_vector);
320
321 if (pj){
322 (*this) = *pj;
323 } else {
324 reset(some_four_vector[0], some_four_vector[1],
325 some_four_vector[2], some_four_vector[3]);
326 }
327 }
328#endif // SWIG
329
330 /// reset the PseudoJet according to the specified pt, rapidity,
331 /// azimuth and mass (also resetting indices, etc.)
332 /// (phi should satisfy -2pi<phi<4pi)
333 inline void reset_PtYPhiM(double pt_in, double y_in, double phi_in, double m_in=0.0) {
334 reset_momentum_PtYPhiM(pt_in, y_in, phi_in, m_in);
335 _reset_indices();
336 //std::shared_ptr: _reset_shared_pointers();
337 }
338
339 /// reset the 4-momentum according to the supplied components
340 /// but leave all other information (indices, user info, etc.)
341 /// untouched
342 inline void reset_momentum(double px, double py, double pz, double E);
343
344 /// reset the 4-momentum according to the components of the supplied
345 /// PseudoJet, including cached components; note that the template
346 /// version (below) will be called for classes derived from PJ.
347 inline void reset_momentum(const PseudoJet & pj);
348
349 /// reset the 4-momentum according to the specified pt, rapidity,
350 /// azimuth and mass (phi should satisfy -2pi<phi<4pi)
351 void reset_momentum_PtYPhiM(double pt, double y, double phi, double m=0.0);
352
353 /// reset the 4-momentum according to the supplied generic 4-vector
354 /// (accessible via indexing, [0]==px,...[3]==E), but leave all
355 /// other information (indices, user info, etc.) untouched
356 template <class L> inline void reset_momentum(const L & some_four_vector) {
357 reset_momentum(some_four_vector[0], some_four_vector[1],
358 some_four_vector[2], some_four_vector[3]);
359 }
360
361 /// in some cases when setting a 4-momentum, the user/program knows
362 /// what rapidity and azimuth are associated with that 4-momentum;
363 /// by calling this routine the user can provide the information
364 /// directly to the PseudoJet and avoid expensive rap-phi
365 /// recalculations.
366 ///
367 /// - \param rap rapidity
368 /// - \param phi (in range -twopi...4*pi)
369 ///
370 /// USE WITH CAUTION: there are no checks that the rapidity and
371 /// azimuth supplied are sensible, nor does this reset the
372 /// 4-momentum components if things don't match.
373 void set_cached_rap_phi(double rap, double phi);
374
375
376 //\} --- end of kin mod functions ------------------------------------
377
378 //----------------------------------------------------------------------
379 /// @name User index functions
380 ///
381 /// To allow the user to set and access an integer index which can
382 /// be exploited by the user to associate extra information with a
383 /// particle/jet (for example pdg id, or an indication of a
384 /// particle's origin within the user's analysis)
385 //
386 //\{
387
388 /// return the user_index,
389 inline int user_index() const {return _user_index;}
390 /// set the user_index, intended to allow the user to add simple
391 /// identifying information to a particle/jet
392 inline void set_user_index(const int index) {_user_index = index;}
393
394 //\} ----- end of use index functions ---------------------------------
395
396 //----------------------------------------------------------------------
397 /// @name User information types and functions
398 ///
399 /// Allows PseudoJet to carry extra user info (as an object derived from
400 /// UserInfoBase).
401 //\{
402
403 /// @ingroup user_info
404 /// \class UserInfoBase
405 /// a base class to hold extra user information in a PseudoJet
406 ///
407 /// This is a base class to help associate extra user information
408 /// with a jet. The user should store their information in a class
409 /// derived from this. This allows information of arbitrary
410 /// complexity to be easily associated with a PseudoJet (in contrast
411 /// to the user index). For example, in a Monte Carlo simulation,
412 /// the user information might include the PDG ID, and the position
413 /// of the production vertex for the particle.
414 ///
415 /// The PseudoJet is able to store a shared pointer to any object
416 /// derived from UserInfo. The use of a shared pointer frees the
417 /// user of the need to handle the memory management associated with
418 /// the information.
419 ///
420 /// Having the user information derive from a common base class also
421 /// facilitates dynamic casting, etc.
422 ///
424 public:
425 // dummy ctor
426 UserInfoBase(){};
427
428 // dummy virtual dtor
429 // makes it polymorphic to allow for dynamic_cast
430 virtual ~UserInfoBase(){};
431 };
432
433 /// error class to be thrown if accessing user info when it doesn't
434 /// exist
435 class InexistentUserInfo : public Error {
436 public:
438 };
439
440 /// sets the internal shared pointer to the user information.
441 ///
442 /// Note that the PseudoJet will now _own_ the pointer, and delete
443 /// the corresponding object when it (the jet, and any copies of the jet)
444 /// goes out of scope.
445 void set_user_info(UserInfoBase * user_info_in) {
446 _user_info.reset(user_info_in);
447 }
448
449 /// returns a reference to the dynamic cast conversion of user_info
450 /// to type L.
451 ///
452 /// Usage: suppose you have previously set the user info with a pointer
453 /// to an object of type MyInfo,
454 ///
455 /// class MyInfo: public PseudoJet::UserInfoBase {
456 /// MyInfo(int id) : _pdg_id(id);
457 /// int pdg_id() const {return _pdg_id;}
458 /// int _pdg_id;
459 /// };
460 ///
461 /// PseudoJet particle(...);
462 /// particle.set_user_info(new MyInfo(its_pdg_id));
463 ///
464 /// Then you would access that pdg_id() as
465 ///
466 /// particle.user_info<MyInfo>().pdg_id();
467 ///
468 /// It's overkill for just a single integer, but scales easily to
469 /// more extensive information.
470 ///
471 /// Note that user_info() throws an InexistentUserInfo() error if
472 /// there is no user info; throws a std::bad_cast if the conversion
473 /// doesn't work
474 ///
475 /// If this behaviour does not fit your needs, use instead the the
476 /// user_info_ptr() or user_info_shared_ptr() member functions.
477 template<class L>
478 const L & user_info() const{
479 if (_user_info.get() == 0) throw InexistentUserInfo();
480 return dynamic_cast<const L &>(* _user_info.get());
481 }
482
483 /// returns true if the PseudoJet has user information
484 bool has_user_info() const{
485 return _user_info.get();
486 }
487
488 /// returns true if the PseudoJet has user information than can be
489 /// cast to the template argument type.
490 template<class L>
491 bool has_user_info() const{
492 return _user_info.get() && dynamic_cast<const L *>(_user_info.get());
493 }
494
495 /// retrieve a pointer to the (const) user information
497 // the line below is not needed since the next line would anyway
498 // return NULL in that case
499 //if (!_user_info) return NULL;
500 return _user_info.get();
501 }
502
503
504 /// retrieve a (const) shared pointer to the user information
506 return _user_info;
507 }
508
509 /// retrieve a (non-const) shared pointer to the user information;
510 /// you can use this, for example, to set the shared pointer, eg
511 ///
512 /// \code
513 /// p2.user_info_shared_ptr() = p1.user_info_shared_ptr();
514 /// \endcode
515 ///
516 /// or
517 ///
518 /// \code
519 /// SharedPtr<PseudoJet::UserInfoBase> info_shared(new MyInfo(...));
520 /// p2.user_info_shared_ptr() = info_shared;
521 /// \endcode
523 return _user_info;
524 }
525
526 void set_user_info_shared_ptr(const SharedPtr<UserInfoBase> & user_info_in) {
527 _user_info = user_info_in;
528 }
529
530 // \} --- end of extra info functions ---------------------------------
531
532 //----------------------------------------------------------------------
533 /// @name Description
534 ///
535 /// Since a PseudoJet can have a structure that contains a variety
536 /// of information, we provide a description that allows one to check
537 /// exactly what kind of PseudoJet we are dealing with
538 //
539 //\{
540
541 /// return a string describing what kind of PseudoJet we are dealing with
542 std::string description() const;
543
544 //\} ----- end of description functions ---------------------------------
545
546 //-------------------------------------------------------------
547 /// @name Access to the associated ClusterSequence object.
548 ///
549 /// In addition to having kinematic information, jets may contain a
550 /// reference to an associated ClusterSequence (this is the case,
551 /// for example, if the jet has been returned by a ClusterSequence
552 /// member function).
553 //\{
554 //-------------------------------------------------------------
555 /// returns true if this PseudoJet has an associated ClusterSequence.
556 bool has_associated_cluster_sequence() const;
557 /// shorthand for has_associated_cluster_sequence()
558 bool has_associated_cs() const {return has_associated_cluster_sequence();}
559
560 /// returns true if this PseudoJet has an associated and still
561 /// valid(ated) ClusterSequence.
562 bool has_valid_cluster_sequence() const;
563 /// shorthand for has_valid_cluster_sequence()
564 bool has_valid_cs() const {return has_valid_cluster_sequence();}
565
566 /// get a (const) pointer to the parent ClusterSequence (NULL if
567 /// inexistent)
568 const ClusterSequence* associated_cluster_sequence() const;
569 // shorthand for associated_cluster_sequence()
570 const ClusterSequence* associated_cs() const {return associated_cluster_sequence();}
571
572 /// if the jet has a valid associated cluster sequence then return a
573 /// pointer to it; otherwise throw an error
575 return validated_cs();
576 }
577 /// shorthand for validated_cluster_sequence()
578 const ClusterSequence * validated_cs() const;
579
580#ifndef __FJCORE__
581 /// if the jet has valid area information then return a pointer to
582 /// the associated ClusterSequenceAreaBase object; otherwise throw an error
584 return validated_csab();
585 }
586
587 /// shorthand for validated_cluster_sequence_area_base()
588 const ClusterSequenceAreaBase * validated_csab() const;
589#endif // __FJCORE__
590
591 //\}
592
593 //-------------------------------------------------------------
594 /// @name Access to the associated PseudoJetStructureBase object.
595 ///
596 /// In addition to having kinematic information, jets may contain a
597 /// reference to an associated ClusterSequence (this is the case,
598 /// for example, if the jet has been returned by a ClusterSequence
599 /// member function).
600 //\{
601 //-------------------------------------------------------------
602
603 /// set the associated structure
604 void set_structure_shared_ptr(const SharedPtr<PseudoJetStructureBase> &structure_in);
605
606 /// return true if there is some structure associated with this PseudoJet
607 bool has_structure() const;
608
609 /// return a pointer to the structure (of type
610 /// PseudoJetStructureBase*) associated with this PseudoJet.
611 ///
612 /// return NULL if there is no associated structure
613 const PseudoJetStructureBase* structure_ptr() const;
614
615 /// return a non-const pointer to the structure (of type
616 /// PseudoJetStructureBase*) associated with this PseudoJet.
617 ///
618 /// return NULL if there is no associated structure
619 ///
620 /// Only use this if you know what you are doing. In any case,
621 /// prefer the 'structure_ptr()' (the const version) to this method,
622 /// unless you really need a write access to the PseudoJet's
623 /// underlying structure.
624 PseudoJetStructureBase* structure_non_const_ptr();
625
626 /// return a pointer to the structure (of type
627 /// PseudoJetStructureBase*) associated with this PseudoJet.
628 ///
629 /// throw an error if there is no associated structure
630 const PseudoJetStructureBase* validated_structure_ptr() const;
631
632 /// return a reference to the shared pointer to the
633 /// PseudoJetStructureBase associated with this PseudoJet
634 const SharedPtr<PseudoJetStructureBase> & structure_shared_ptr() const;
635
636 /// returns a reference to the structure casted to the requested
637 /// structure type
638 ///
639 /// If there is no structure associated, an Error is thrown.
640 /// If the type is not met, a std::bad_cast error is thrown.
641 template<typename StructureType>
642 const StructureType & structure() const;
643
644 /// check if the PseudoJet has the structure resulting from a Transformer
645 /// (that is, its structure is compatible with a Transformer::StructureType).
646 /// If there is no structure, false is returned.
647 template<typename TransformerType>
648 bool has_structure_of() const;
649
650 /// this is a helper to access any structure created by a Transformer
651 /// (that is, of type Transformer::StructureType).
652 ///
653 /// If there is no structure, or if the structure is not compatible
654 /// with TransformerType, an error is thrown.
655 template<typename TransformerType>
656 const typename TransformerType::StructureType & structure_of() const;
657
658 //\}
659
660 //-------------------------------------------------------------
661 /// @name Methods for access to information about jet structure
662 ///
663 /// These allow access to jet constituents, and other jet
664 /// subtructure information. They only work if the jet is associated
665 /// with a ClusterSequence.
666 //-------------------------------------------------------------
667 //\{
668
669 /// check if it has been recombined with another PseudoJet in which
670 /// case, return its partner through the argument. Otherwise,
671 /// 'partner' is set to 0.
672 ///
673 /// an Error is thrown if this PseudoJet has no currently valid
674 /// associated ClusterSequence
675 virtual bool has_partner(PseudoJet &partner) const;
676
677 /// check if it has been recombined with another PseudoJet in which
678 /// case, return its child through the argument. Otherwise, 'child'
679 /// is set to 0.
680 ///
681 /// an Error is thrown if this PseudoJet has no currently valid
682 /// associated ClusterSequence
683 virtual bool has_child(PseudoJet &child) const;
684
685 /// check if it is the product of a recombination, in which case
686 /// return the 2 parents through the 'parent1' and 'parent2'
687 /// arguments. Otherwise, set these to 0.
688 ///
689 /// an Error is thrown if this PseudoJet has no currently valid
690 /// associated ClusterSequence
691 virtual bool has_parents(PseudoJet &parent1, PseudoJet &parent2) const;
692
693 /// check if the current PseudoJet contains the one passed as
694 /// argument.
695 ///
696 /// an Error is thrown if this PseudoJet has no currently valid
697 /// associated ClusterSequence
698 virtual bool contains(const PseudoJet &constituent) const;
699
700 /// check if the current PseudoJet is contained in the one passed as
701 /// argument.
702 ///
703 /// an Error is thrown if this PseudoJet has no currently valid
704 /// associated ClusterSequence
705 virtual bool is_inside(const PseudoJet &jet) const;
706
707
708 /// returns true if the PseudoJet has constituents
709 virtual bool has_constituents() const;
710
711 /// retrieve the constituents.
712 ///
713 /// an Error is thrown if this PseudoJet has no currently valid
714 /// associated ClusterSequence or other substructure information
715 virtual std::vector<PseudoJet> constituents() const;
716
717
718 /// returns true if the PseudoJet has support for exclusive subjets
719 virtual bool has_exclusive_subjets() const;
720
721 /// return a vector of all subjets of the current jet (in the sense
722 /// of the exclusive algorithm) that would be obtained when running
723 /// the algorithm with the given dcut.
724 ///
725 /// Time taken is O(m ln m), where m is the number of subjets that
726 /// are found. If m gets to be of order of the total number of
727 /// constituents in the jet, this could be substantially slower than
728 /// just getting that list of constituents.
729 ///
730 /// an Error is thrown if this PseudoJet has no currently valid
731 /// associated ClusterSequence
732 std::vector<PseudoJet> exclusive_subjets (const double dcut) const;
733
734 /// return the size of exclusive_subjets(...); still n ln n with same
735 /// coefficient, but marginally more efficient than manually taking
736 /// exclusive_subjets.size()
737 ///
738 /// an Error is thrown if this PseudoJet has no currently valid
739 /// associated ClusterSequence
740 int n_exclusive_subjets(const double dcut) const;
741
742 /// return the list of subjets obtained by unclustering the supplied
743 /// jet down to nsub subjets. Throws an error if there are fewer than
744 /// nsub particles in the jet.
745 ///
746 /// For ClusterSequence type jets, requires nsub ln nsub time
747 ///
748 /// An Error is thrown if this PseudoJet has no currently valid
749 /// associated ClusterSequence
750 std::vector<PseudoJet> exclusive_subjets (int nsub) const;
751
752 /// return the list of subjets obtained by unclustering the supplied
753 /// jet down to nsub subjets (or all constituents if there are fewer
754 /// than nsub).
755 ///
756 /// For ClusterSequence type jets, requires nsub ln nsub time
757 ///
758 /// An Error is thrown if this PseudoJet has no currently valid
759 /// associated ClusterSequence
760 std::vector<PseudoJet> exclusive_subjets_up_to (int nsub) const;
761
762 /// Returns the dij that was present in the merging nsub+1 -> nsub
763 /// subjets inside this jet.
764 ///
765 /// Returns 0 if there were nsub or fewer constituents in the jet.
766 ///
767 /// an Error is thrown if this PseudoJet has no currently valid
768 /// associated ClusterSequence
769 double exclusive_subdmerge(int nsub) const;
770
771 /// Returns the maximum dij that occurred in the whole event at the
772 /// stage that the nsub+1 -> nsub merge of subjets occurred inside
773 /// this jet.
774 ///
775 /// Returns 0 if there were nsub or fewer constituents in the jet.
776 ///
777 /// an Error is thrown if this PseudoJet has no currently valid
778 /// associated ClusterSequence
779 double exclusive_subdmerge_max(int nsub) const;
780
781
782 /// returns true if a jet has pieces
783 ///
784 /// By default a single particle or a jet coming from a
785 /// ClusterSequence have no pieces and this methos will return false.
786 ///
787 /// In practice, this is equivalent to have an structure of type
788 /// CompositeJetStructure.
789 virtual bool has_pieces() const;
790
791
792 /// retrieve the pieces that make up the jet.
793 ///
794 /// If the jet does not support pieces, an error is throw
795 virtual std::vector<PseudoJet> pieces() const;
796
797
798 // the following ones require a computation of the area in the
799 // parent ClusterSequence (See ClusterSequenceAreaBase for details)
800 //------------------------------------------------------------------
801#ifndef __FJCORE__
802
803 /// check if it has a defined area
804 virtual bool has_area() const;
805
806 /// return the jet (scalar) area.
807 /// throws an Error if there is no support for area in the parent CS
808 virtual double area() const;
809
810 /// return the error (uncertainty) associated with the determination
811 /// of the area of this jet.
812 /// throws an Error if there is no support for area in the parent CS
813 virtual double area_error() const;
814
815 /// return the jet 4-vector area.
816 /// throws an Error if there is no support for area in the parent CS
817 virtual PseudoJet area_4vector() const;
818
819 /// true if this jet is made exclusively of ghosts.
820 /// throws an Error if there is no support for area in the parent CS
821 virtual bool is_pure_ghost() const;
822
823#endif // __FJCORE__
824 //\} --- end of jet structure -------------------------------------
825
826
827
828 //----------------------------------------------------------------------
829 /// @name Members mainly intended for internal use
830 //----------------------------------------------------------------------
831 //\{
832 /// return the cluster_hist_index, intended to be used by clustering
833 /// routines.
834 inline int cluster_hist_index() const {return _cluster_hist_index;}
835 /// set the cluster_hist_index, intended to be used by clustering routines.
836 inline void set_cluster_hist_index(const int index) {_cluster_hist_index = index;}
837
838 /// alternative name for cluster_hist_index() [perhaps more meaningful]
840 return cluster_hist_index();}
841 /// alternative name for set_cluster_hist_index(...) [perhaps more
842 /// meaningful]
843 inline void set_cluster_sequence_history_index(const int index) {
844 set_cluster_hist_index(index);}
845
846 //\} ---- end of internal use functions ---------------------------
847
848 protected:
849
851 SharedPtr<UserInfoBase> _user_info;
852
853
854 private:
855 // NB: following order must be kept for things to behave sensibly...
856 double _px,_py,_pz,_E;
857 mutable double _phi, _rap;
858 double _kt2;
859 int _cluster_hist_index, _user_index;
860
861#ifdef FASTJET_HAVE_THREAD_SAFETY
862 enum {
863 Init_Done=1,
864 Init_NotDone=0,
865 Init_InProgress=-1
866 };
867
868 mutable std::atomic<int> _init_status;
869#endif
870
871 /// calculate phi, rap, kt2 based on the 4-momentum components
872 void _finish_init();
873 /// set the indices to default values
874 void _reset_indices();
875
876//std::shared_ptr: /// reset the shared pointers to empty ones
877//std::shared_ptr: void _reset_shared_pointers();
878//std::shared_ptr:
879//std::shared_ptr: #ifdef FASTJET_HAVE_THREAD_SAFETY
880//std::shared_ptr: /// For jets associated with a ClusterSequence, this will "free" the
881//std::shared_ptr: /// jet from the ClusterSequence. This means that, for self-deleting
882//std::shared_ptr: /// cluster sequences, it will reset the structure pointer and check
883//std::shared_ptr: /// if the CS needs to be deleted.
884//std::shared_ptr: ///
885//std::shared_ptr: /// This is the replacement mechanism for self-deleting cs (with
886//std::shared_ptr: /// set_count not working with std shared_ptr) and...
887//std::shared_ptr: ///
888//std::shared_ptr: /// IT HAS TO BE CALLED BEFORE ANY CHANGE OF THE JET STRUCTURE
889//std::shared_ptr: /// POINTER
890//std::shared_ptr: void _release_jet_from_cs();
891//std::shared_ptr: #endif //FASTJET_HAVE_THREAD_SAFETY
892
893 /// ensure that the internal values for rapidity and phi
894 /// correspond to 4-momentum structure
895#ifdef FASTJET_HAVE_THREAD_SAFETY
896 void _ensure_valid_rap_phi() const;
897#else
898 inline void _ensure_valid_rap_phi() const{
899 if (_phi == pseudojet_invalid_phi) _set_rap_phi();
900 }
901#endif
902
903 /// set cached rapidity and phi values
904 void _set_rap_phi() const;
905
906 // needed for operator* to have access to _ensure_valid_rap_phi()
907 friend PseudoJet operator*(double, const PseudoJet &);
908};
909
910
911//----------------------------------------------------------------------
912// routines for basic binary operations
913
914PseudoJet operator+(const PseudoJet &, const PseudoJet &);
915PseudoJet operator-(const PseudoJet &, const PseudoJet &);
916PseudoJet operator*(double, const PseudoJet &);
917PseudoJet operator*(const PseudoJet &, double);
918PseudoJet operator/(const PseudoJet &, double);
919
920/// returns true if the 4 momentum components of the two PseudoJets
921/// are identical and all the internal indices (user, cluster_history)
922/// + structure and user-info shared pointers are too
923bool operator==(const PseudoJet &, const PseudoJet &);
924
925/// inequality test which is exact opposite of operator==
926inline bool operator!=(const PseudoJet & a, const PseudoJet & b) {return !(a==b);}
927
928/// Can only be used with val=0 and tests whether all four
929/// momentum components are equal to val (=0.0)
930bool operator==(const PseudoJet & jet, const double val);
931inline bool operator==(const double val, const PseudoJet & jet) {return jet == val;}
932
933/// Can only be used with val=0 and tests whether at least one of the
934/// four momentum components is different from val (=0.0)
935inline bool operator!=(const PseudoJet & a, const double val) {return !(a==val);}
936inline bool operator!=( const double val, const PseudoJet & a) {return !(a==val);}
937
938/// returns the 4-vector dot product of a and b
939inline double dot_product(const PseudoJet & a, const PseudoJet & b) {
940 return a.E()*b.E() - a.px()*b.px() - a.py()*b.py() - a.pz()*b.pz();
941}
942
943/// returns the cosine of the angle between a and b
944inline double cos_theta(const PseudoJet & a, const PseudoJet & b) {
945 double dot_3d = a.px()*b.px() + a.py()*b.py() + a.pz()*b.pz();
946 return std::min(1.0, std::max(-1.0, dot_3d/sqrt(a.modp2()*b.modp2())));
947}
948
949/// returns the angle between a and b
950inline double theta(const PseudoJet & a, const PseudoJet & b) {
951 return acos(cos_theta(a,b));
952}
953
954/// returns true if the momenta of the two input jets are identical
955bool have_same_momentum(const PseudoJet &, const PseudoJet &);
956
957/// return a pseudojet with the given pt, y, phi and mass
958/// (phi should satisfy -2pi<phi<4pi)
959PseudoJet PtYPhiM(double pt, double y, double phi, double m = 0.0);
960
961//----------------------------------------------------------------------
962// Routines to do with providing sorted arrays of vectors.
963
964/// return a vector of jets sorted into decreasing transverse momentum
965std::vector<PseudoJet> sorted_by_pt(const std::vector<PseudoJet> & jets);
966
967/// return a vector of jets sorted into increasing rapidity
968std::vector<PseudoJet> sorted_by_rapidity(const std::vector<PseudoJet> & jets);
969
970/// return a vector of jets sorted into decreasing energy
971std::vector<PseudoJet> sorted_by_E(const std::vector<PseudoJet> & jets);
972
973/// return a vector of jets sorted into increasing pz
974std::vector<PseudoJet> sorted_by_pz(const std::vector<PseudoJet> & jets);
975
976//----------------------------------------------------------------------
977// some code to help sorting
978
979/// sort the indices so that values[indices[0->n-1]] is sorted
980/// into increasing order
981void sort_indices(std::vector<int> & indices,
982 const std::vector<double> & values);
983
984/// given a vector of values with a one-to-one correspondence with the
985/// vector of objects, sort objects into an order such that the
986/// associated values would be in increasing order (but don't actually
987/// touch the values vector in the process).
988template<class T> std::vector<T> objects_sorted_by_values(const std::vector<T> & objects,
989 const std::vector<double> & values) {
990 //assert(objects.size() == values.size());
991 if (objects.size() != values.size()){
992 throw Error("fastjet::objects_sorted_by_values(...): the size of the 'objects' vector must match the size of the 'values' vector");
993 }
994
995 // get a vector of indices
996 std::vector<int> indices(values.size());
997 for (size_t i = 0; i < indices.size(); i++) {indices[i] = i;}
998
999 // sort the indices
1000 sort_indices(indices, values);
1001
1002 // copy the objects
1003 std::vector<T> objects_sorted(objects.size());
1004
1005 // place the objects in the correct order
1006 for (size_t i = 0; i < indices.size(); i++) {
1007 objects_sorted[i] = objects[indices[i]];
1008 }
1009
1010 return objects_sorted;
1011}
1012
1013/// \if internal_doc
1014/// @ingroup internal
1015/// \class IndexedSortHelper
1016/// a class that helps us carry out indexed sorting.
1017/// \endif
1019public:
1020 inline IndexedSortHelper (const std::vector<double> * reference_values) {
1021 _ref_values = reference_values;
1022 };
1023 inline int operator() (const int i1, const int i2) const {
1024 return (*_ref_values)[i1] < (*_ref_values)[i2];
1025 };
1026private:
1027 const std::vector<double> * _ref_values;
1028};
1029
1030
1031//----------------------------------------------------------------------
1032/// constructor from any object that has px,py,pz,E = some_four_vector[0--3],
1033// NB: do not know if it really needs to be inline, but when it wasn't
1034// linking failed with g++ (who knows what was wrong...)
1035#ifndef SWIG
1036template <class L> inline PseudoJet::PseudoJet(const L & some_four_vector) {
1037 reset(some_four_vector);
1038}
1039#endif
1040
1041//----------------------------------------------------------------------
1042inline void PseudoJet::_reset_indices() {
1043 set_cluster_hist_index(-1);
1044 set_user_index(-1);
1045
1046 _structure.reset();
1047 _user_info.reset();
1048}
1049
1050//std::shared_ptr: inline void PseudoJet::_reset_shared_pointers() {
1051//std::shared_ptr: #ifdef FASTJET_HAVE_THREAD_SAFETY
1052//std::shared_ptr: // if the jet currently belongs to a cs, we need to release it before any chenge
1053//std::shared_ptr: _release_jet_from_cs();
1054//std::shared_ptr: #endif // FASTJET_HAVE_THREAD_SAFETY
1055//std::shared_ptr: _structure.reset();
1056//std::shared_ptr: _user_info.reset();
1057//std::shared_ptr: // and do not forget to remove it in reset_indices above
1058//std::shared_ptr: }
1059
1060
1061
1062
1063// taken literally from CLHEP
1064inline double PseudoJet::m() const {
1065 double mm = m2();
1066 return mm < 0.0 ? -std::sqrt(-mm) : std::sqrt(mm);
1067}
1068
1069
1070inline void PseudoJet::reset(double px_in, double py_in, double pz_in, double E_in) {
1071 _px = px_in;
1072 _py = py_in;
1073 _pz = pz_in;
1074 _E = E_in;
1075 _finish_init();
1076 _reset_indices();
1077 //std::shared_ptr: _reset_shared_pointers();
1078}
1079
1080inline void PseudoJet::reset_momentum(double px_in, double py_in, double pz_in, double E_in) {
1081 _px = px_in;
1082 _py = py_in;
1083 _pz = pz_in;
1084 _E = E_in;
1085 _finish_init();
1086}
1087
1088inline void PseudoJet::reset_momentum(const PseudoJet & pj) {
1089 _px = pj._px ;
1090 _py = pj._py ;
1091 _pz = pj._pz ;
1092 _E = pj._E ;
1093#ifdef FASTJET_HAVE_THREAD_SAFETY
1094 // the following lines are here to address the situation where
1095 // the current form of the jet has Init_Done, but the other jet does
1096 // not (or is in the process of being initialized, so that copying
1097 // is dangerous).
1098 // We take the approach of verifying if the pj has been initialised
1099 // and if it has we take its cached phi & rapidity, otherwise
1100 // we call _finish_init(), which sets Init_NotDone and puts
1101 // rapidity & phi to invalid values.
1102 int expected = Init_Done;
1103 if (pj._init_status.compare_exchange_weak(expected, Init_Done)) {
1104 _init_status = Init_Done;
1105 _phi = pj._phi;
1106 _rap = pj._rap;
1107 _kt2 = pj._kt2;
1108 } else {
1109 _finish_init();
1110 }
1111#else
1112 _phi = pj._phi;
1113 _rap = pj._rap;
1114 _kt2 = pj._kt2;
1115#endif
1116}
1117
1118//-------------------------------------------------------------------------------
1119// implementation of the templated accesses to the underlying structyre
1120//-------------------------------------------------------------------------------
1121
1122// returns a reference to the structure casted to the requested
1123// structure type
1124//
1125// If there is no sructure associated, an Error is thrown.
1126// If the type is not met, a std::bad_cast error is thrown.
1127template<typename StructureType>
1128const StructureType & PseudoJet::structure() const{
1129 return dynamic_cast<const StructureType &>(* validated_structure_ptr());
1130
1131}
1132
1133// check if the PseudoJet has the structure resulting from a Transformer
1134// (that is, its structure is compatible with a Transformer::StructureType)
1135template<typename TransformerType>
1136bool PseudoJet::has_structure_of() const{
1137 if (!_structure) return false;
1138
1139 return dynamic_cast<const typename TransformerType::StructureType *>(_structure.get()) != 0;
1140}
1141
1142// this is a helper to access a structure created by a Transformer
1143// (that is, of type Transformer::StructureType)
1144// NULL is returned if the corresponding type is not met
1145template<typename TransformerType>
1146const typename TransformerType::StructureType & PseudoJet::structure_of() const{
1147 if (!_structure)
1148 throw Error("Trying to access the structure of a PseudoJet without an associated structure");
1149
1150 return dynamic_cast<const typename TransformerType::StructureType &>(*_structure);
1151}
1152
1153
1154
1155//-------------------------------------------------------------------------------
1156// helper functions to build a jet made of pieces
1157//
1158// Note that there are more complete versions of these functions, with
1159// an additional argument for a recombination scheme, in
1160// JetDefinition.hh
1161// -------------------------------------------------------------------------------
1162
1163/// build a "CompositeJet" from the vector of its pieces
1164///
1165/// In this case, E-scheme recombination is assumed to compute the
1166/// total momentum
1167PseudoJet join(const std::vector<PseudoJet> & pieces);
1168
1169/// build a MergedJet from a single PseudoJet
1170PseudoJet join(const PseudoJet & j1);
1171
1172/// build a MergedJet from 2 PseudoJet
1173PseudoJet join(const PseudoJet & j1, const PseudoJet & j2);
1174
1175/// build a MergedJet from 3 PseudoJet
1176PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3);
1177
1178/// build a MergedJet from 4 PseudoJet
1179PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3, const PseudoJet & j4);
1180
1181
1182
1183FASTJET_END_NAMESPACE
1184
1185#endif // __FASTJET_PSEUDOJET_HH__
base class that sets interface for extensions of ClusterSequence that provide information about the a...
deals with clustering
base class corresponding to errors that can be thrown by FastJet
Definition: Error.hh:52
Contains any information related to the clustering that should be directly accessible to PseudoJet.
error class to be thrown if accessing user info when it doesn't exist
Definition: PseudoJet.hh:435
a base class to hold extra user information in a PseudoJet
Definition: PseudoJet.hh:423
Class to contain pseudojets, including minimal information of use to jet-clustering routines.
Definition: PseudoJet.hh:68
PseudoJet()
default constructor, which as of FJ3.0 provides an object for which all operations are now valid and ...
Definition: PseudoJet.hh:82
double cos_theta() const
cos of the polar angle should we have: min(1.0,max(-1.0,_pz/sqrt(modp2())));
Definition: PseudoJet.hh:189
const ClusterSequenceAreaBase * validated_cluster_sequence_area_base() const
if the jet has valid area information then return a pointer to the associated ClusterSequenceAreaBase...
Definition: PseudoJet.hh:583
double mt() const
returns the transverse mass = sqrt(kt^2+m^2)
Definition: PseudoJet.hh:175
int cluster_sequence_history_index() const
alternative name for cluster_hist_index() [perhaps more meaningful]
Definition: PseudoJet.hh:839
double rap() const
returns the rapidity or some large value when the rapidity is infinite
Definition: PseudoJet.hh:138
double Et() const
return the transverse energy
Definition: PseudoJet.hh:183
double delta_R(const PseudoJet &other) const
return the cylinder (rap-phi) distance between this jet and another, .
Definition: PseudoJet.hh:214
const L & user_info() const
returns a reference to the dynamic cast conversion of user_info to type L.
Definition: PseudoJet.hh:478
double mperp() const
returns the transverse mass = sqrt(kt^2+m^2)
Definition: PseudoJet.hh:171
double modp() const
return the 3-vector modulus = sqrt(px^2+py^2+pz^2)
Definition: PseudoJet.hh:180
void set_user_index(const int index)
set the user_index, intended to allow the user to add simple identifying information to a particle/je...
Definition: PseudoJet.hh:392
const UserInfoBase * user_info_ptr() const
retrieve a pointer to the (const) user information
Definition: PseudoJet.hh:496
void reset_momentum(const L &some_four_vector)
reset the 4-momentum according to the supplied generic 4-vector (accessible via indexing,...
Definition: PseudoJet.hh:356
bool has_user_info() const
returns true if the PseudoJet has user information
Definition: PseudoJet.hh:484
double squared_distance(const PseudoJet &other) const
returns squared cylinder (rap-phi) distance between this jet and another
Definition: PseudoJet.hh:209
bool has_valid_cs() const
shorthand for has_valid_cluster_sequence()
Definition: PseudoJet.hh:564
double phi() const
returns phi (in the range 0..2pi)
Definition: PseudoJet.hh:123
virtual ~PseudoJet()
default (virtual) destructor
Definition: PseudoJet.hh:104
double perp() const
returns the scalar transverse momentum
Definition: PseudoJet.hh:158
double perp2() const
returns the squared transverse momentum
Definition: PseudoJet.hh:156
double m2() const
returns the squared invariant mass // like CLHEP
Definition: PseudoJet.hh:163
void set_cluster_hist_index(const int index)
set the cluster_hist_index, intended to be used by clustering routines.
Definition: PseudoJet.hh:836
double kt2() const
returns the squared transverse momentum
Definition: PseudoJet.hh:160
bool has_associated_cs() const
shorthand for has_associated_cluster_sequence()
Definition: PseudoJet.hh:558
bool has_user_info() const
returns true if the PseudoJet has user information than can be cast to the template argument type.
Definition: PseudoJet.hh:491
double theta() const
polar angle
Definition: PseudoJet.hh:193
double phi_std() const
returns phi in the range -pi..pi
Definition: PseudoJet.hh:126
double phi_02pi() const
returns phi in the range 0..2pi
Definition: PseudoJet.hh:131
int cluster_hist_index() const
return the cluster_hist_index, intended to be used by clustering routines.
Definition: PseudoJet.hh:834
void reset_PtYPhiM(double pt_in, double y_in, double phi_in, double m_in=0.0)
reset the PseudoJet according to the specified pt, rapidity, azimuth and mass (also resetting indices...
Definition: PseudoJet.hh:333
SharedPtr< UserInfoBase > & user_info_shared_ptr()
retrieve a (non-const) shared pointer to the user information; you can use this, for example,...
Definition: PseudoJet.hh:522
double pt() const
returns the scalar transverse momentum
Definition: PseudoJet.hh:154
void set_cluster_sequence_history_index(const int index)
alternative name for set_cluster_hist_index(...) [perhaps more meaningful]
Definition: PseudoJet.hh:843
double mperp2() const
returns the squared transverse mass = kt^2+m^2
Definition: PseudoJet.hh:169
void reset(const L &some_four_vector)
reset the 4-momentum according to the supplied generic 4-vector (accessible via indexing,...
Definition: PseudoJet.hh:307
int user_index() const
return the user_index,
Definition: PseudoJet.hh:389
const SharedPtr< UserInfoBase > & user_info_shared_ptr() const
retrieve a (const) shared pointer to the user information
Definition: PseudoJet.hh:505
const ClusterSequence * validated_cluster_sequence() const
if the jet has a valid associated cluster sequence then return a pointer to it; otherwise throw an er...
Definition: PseudoJet.hh:574
double modp2() const
return the squared 3-vector modulus = px^2+py^2+pz^2
Definition: PseudoJet.hh:178
double Et2() const
return the transverse energy squared
Definition: PseudoJet.hh:185
double rapidity() const
the same as rap()
Definition: PseudoJet.hh:144
double beam_distance() const
returns distance between this jet and the beam
Definition: PseudoJet.hh:228
double pt2() const
returns the squared transverse momentum
Definition: PseudoJet.hh:152
void set_user_info(UserInfoBase *user_info_in)
sets the internal shared pointer to the user information.
Definition: PseudoJet.hh:445
void reset(const PseudoJet &psjet)
reset the PseudoJet to be equal to psjet (including its indices); NB if the argument is derived from ...
Definition: PseudoJet.hh:299
double mt2() const
returns the squared transverse mass = kt^2+m^2
Definition: PseudoJet.hh:173
An implementation of shared pointers that is broadly similar to C++11 shared_ptr (https://en....
Definition: SharedPtr.hh:341
T * get() const
get the stored pointer
Definition: SharedPtr.hh:473
void reset()
reset the pointer to default value (NULL)
Definition: SharedPtr.hh:381
Selector operator*(const Selector &s1, const Selector &s2)
successive application of 2 selectors
Definition: Selector.cc:559
bool operator!=(const PseudoJet &a, const PseudoJet &b)
inequality test which is exact opposite of operator==
Definition: PseudoJet.hh:926
vector< PseudoJet > sorted_by_rapidity(const vector< PseudoJet > &jets)
return a vector of jets sorted into increasing rapidity
Definition: PseudoJet.cc:879
const double MaxRap
Used to protect against parton-level events where pt can be zero for some partons,...
Definition: PseudoJet.hh:53
vector< PseudoJet > sorted_by_E(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing energy
Definition: PseudoJet.cc:887
bool operator==(const PseudoJet &a, const PseudoJet &b)
returns true if the 4 momentum components of the two PseudoJets are identical and all the internal in...
Definition: PseudoJet.cc:345
std::vector< T > objects_sorted_by_values(const std::vector< T > &objects, const std::vector< double > &values)
given a vector of values with a one-to-one correspondence with the vector of objects,...
Definition: PseudoJet.hh:988
double dot_product(const PseudoJet &a, const PseudoJet &b)
returns the 4-vector dot product of a and b
Definition: PseudoJet.hh:939
const double pseudojet_invalid_phi
default value for phi, meaning it (and rapidity) have yet to be calculated)
Definition: PseudoJet.hh:56
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
Definition: PseudoJet.cc:871
double theta(const PseudoJet &a, const PseudoJet &b)
returns the angle between a and b
Definition: PseudoJet.hh:950
double cos_theta(const PseudoJet &a, const PseudoJet &b)
returns the cosine of the angle between a and b
Definition: PseudoJet.hh:944
vector< PseudoJet > sorted_by_pz(const vector< PseudoJet > &jets)
return a vector of jets sorted into increasing pz
Definition: PseudoJet.cc:895