FastJet 3.4.1
MassDropTagger.cc
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#include "fastjet/tools/MassDropTagger.hh"
32#include "fastjet/ClusterSequence.hh"
33#include <sstream>
34
35FASTJET_BEGIN_NAMESPACE
36
37LimitedWarning MassDropTagger::_warnings_nonca;
38LimitedWarning MassDropTagger::_negative_mass_warning;
39
40using namespace std;
41
42//----------------------------------------------------------------------
43// MassDropTagger class implementation
44//----------------------------------------------------------------------
45
46//------------------------------------------------------------------------
47// description of the tagger
48string MassDropTagger::description() const{
49 ostringstream oss;
50 oss << "MassDropTagger with mu=" << _mu << " and ycut=" << _ycut;
51 return oss.str();
52}
53
54//------------------------------------------------------------------------
55// returns the tagged PseudoJet if successful, 0 otherwise
56// - jet the PseudoJet to tag
57PseudoJet MassDropTagger::result(const PseudoJet & jet) const{
58 PseudoJet j = jet;
59
60 // issue a warning if the jet is not obtained through a C/A
61 // clustering
64 _warnings_nonca.warn("MassDropTagger should only be applied on jets from a Cambridge/Aachen clustering; use it with other algorithms at your own risk.");
65
66
67 PseudoJet j1, j2;
68 bool had_parents;
69
70 // we just ask that we can "walk" in the cluster sequence.
71 // appropriate errors will be thrown automatically if this is not
72 // the case
73 while ((had_parents = j.has_parents(j1,j2))) {
74 if (j.m2() <= 0) {
75 _negative_mass_warning.warn(
76 "MassDropTagger: parent (sub)jet has mass^2<=0; returning null jet");
77 return PseudoJet();
78 }
79 // make parent1 the more massive jet
80 if (j1.m2() < j2.m2()) std::swap(j1,j2);
81
82 // if we pass the conditions on the mass drop and its degree of
83 // asymmetry (kt_dist/m^2 > rtycut [where kt_dist/m^2 \sim
84 // z/(1-z)), then we've found something interesting, so exit the
85 // loop
86 if ( (j1.m2() < _mu*_mu*j.m2()) && (j1.kt_distance(j2) > _ycut*j.m2()) )
87 break;
88 else
89 j = j1;
90 }
91
92 if (!had_parents)
93 // no Higgs found, return an empty PseudoJet
94 return PseudoJet();
95
96 // create the result and its structure
97 PseudoJet result_local = j;
99 s->_mu = j1.m() / j.m();
100 s->_y = j1.kt_distance(j2)/j.m2();
101
103
104 return result_local;
105}
106
107FASTJET_END_NAMESPACE
108
const JetDefinition & jet_def() const
return a reference to the jet definition
JetAlgorithm jet_algorithm() const
return information about the definition...
the structure returned by the MassDropTagger transformer.
double _y
the value of the asymmetry parameter
double _mu
the value of the mass-drop parameter
Class to contain pseudojets, including minimal information of use to jet-clustering routines.
Definition: PseudoJet.hh:68
double m2() const
returns the squared invariant mass // like CLHEP
Definition: PseudoJet.hh:163
double kt_distance(const PseudoJet &other) const
returns kt distance (R=1) between this jet and another
Definition: PseudoJet.cc:475
virtual bool has_parents(PseudoJet &parent1, PseudoJet &parent2) const
check if it is the product of a recombination, in which case return the 2 parents through the 'parent...
Definition: PseudoJet.cc:648
void set_structure_shared_ptr(const SharedPtr< PseudoJetStructureBase > &structure_in)
set the associated structure
Definition: PseudoJet.cc:561
double m() const
returns the invariant mass (If m2() is negative then -sqrt(-m2()) is returned, as in CLHEP)
Definition: PseudoJet.hh:1064
const ClusterSequence * validated_cs() const
shorthand for validated_cluster_sequence()
Definition: PseudoJet.cc:554
bool has_associated_cluster_sequence() const
returns true if this PseudoJet has an associated ClusterSequence.
Definition: PseudoJet.cc:527
An implementation of shared pointers that is broadly similar to C++11 shared_ptr (https://en....
Definition: SharedPtr.hh:341
@ cambridge_algorithm
the longitudinally invariant variant of the cambridge algorithm (aka Aachen algoithm).