FastJet  3.1.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
MassDropTagger.cc
1 //FJSTARTHEADER
2 // $Id: MassDropTagger.cc 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 #include <fastjet/tools/MassDropTagger.hh>
32 #include <fastjet/ClusterSequence.hh>
33 #include <sstream>
34 
35 FASTJET_BEGIN_NAMESPACE
36 
37 LimitedWarning MassDropTagger::_warnings_nonca;
38 
39 using namespace std;
40 
41 //----------------------------------------------------------------------
42 // MassDropTagger class implementation
43 //----------------------------------------------------------------------
44 
45 //------------------------------------------------------------------------
46 // description of the tagger
47 string MassDropTagger::description() const{
48  ostringstream oss;
49  oss << "MassDropTagger with mu=" << _mu << " and ycut=" << _ycut;
50  return oss.str();
51 }
52 
53 //------------------------------------------------------------------------
54 // returns the tagged PseudoJet if successful, 0 otherwise
55 // - jet the PseudoJet to tag
56 PseudoJet MassDropTagger::result(const PseudoJet & jet) const{
57  PseudoJet j = jet;
58 
59  // issue a warning if the jet is not obtained through a C/A
60  // clustering
61  if ((! j.has_associated_cluster_sequence()) ||
63  _warnings_nonca.warn("MassDropTagger should only be applied on jets from a Cambridge/Aachen clustering; use it with other algorithms at your own risk.");
64 
65 
66  PseudoJet j1, j2;
67  bool had_parents;
68 
69  // we just ask that we can "walk" in the cluster sequence.
70  // appropriate errors will be thrown automatically if this is not
71  // the case
72  while ((had_parents = j.has_parents(j1,j2))) {
73  // make parent1 the more massive jet
74  if (j1.m2() < j2.m2()) std::swap(j1,j2);
75 
76  // if we pass the conditions on the mass drop and its degree of
77  // asymmetry (kt_dist/m^2 > rtycut [where kt_dist/m^2 \sim
78  // z/(1-z)), then we've found something interesting, so exit the
79  // loop
80  if ( (j1.m2() < _mu*_mu*j.m2()) && (j1.kt_distance(j2) > _ycut*j.m2()) )
81  break;
82  else
83  j = j1;
84  }
85 
86  if (!had_parents)
87  // no Higgs found, return an empty PseudoJet
88  return PseudoJet();
89 
90  // create the result and its structure
91  PseudoJet result_local = j;
92  MassDropTaggerStructure * s = new MassDropTaggerStructure(result_local);
93 // s->_original_jet = jet;
94  s->_mu = (j.m2()!=0.0) ? sqrt(j1.m2()/j.m2()) : 0.0;
95  s->_y = (j.m2()!=0.0) ? j1.kt_distance(j2)/j.m2() : 0.0;
96 
98 
99  return result_local;
100 }
101 
102 FASTJET_END_NAMESPACE
103 
JetAlgorithm jet_algorithm() const
return information about the definition...
const ClusterSequence * validated_cs() const
shorthand for validated_cluster_sequence()
Definition: PseudoJet.cc:457
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:547
double m2() const
returns the squared invariant mass // like CLHEP
Definition: PseudoJet.hh:148
const JetDefinition & jet_def() const
return a reference to the jet definition
the structure returned by the MassDropTagger transformer.
double _y
the value of the asymmetry parameter
double kt_distance(const PseudoJet &other) const
returns kt distance (R=1) between this jet and another
Definition: PseudoJet.cc:378
an implementation of C++0x shared pointers (or boost's)
Definition: SharedPtr.hh:116
void set_structure_shared_ptr(const SharedPtr< PseudoJetStructureBase > &structure)
set the associated structure
Definition: PseudoJet.cc:464
the longitudinally invariant variant of the cambridge algorithm (aka Aachen algoithm).
bool has_associated_cluster_sequence() const
returns true if this PseudoJet has an associated ClusterSequence.
Definition: PseudoJet.cc:430
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
Definition: PseudoJet.hh:67
double _mu
the value of the mass-drop parameter