29 #include <fastjet/tools/JHTopTagger.hh>
30 #include <fastjet/Error.hh>
31 #include <fastjet/JetDefinition.hh>
32 #include <fastjet/ClusterSequence.hh>
36 FASTJET_BEGIN_NAMESPACE
44 LimitedWarning JHTopTagger::_warnings_nonca;
48 string JHTopTagger::description()
const{
50 oss <<
"JHTopTagger with delta_p=" << _delta_p <<
", delta_r=" << _delta_r
51 <<
", cos_theta_W_max=" << _cos_theta_W_max
52 <<
" and mW = " << _mW;
53 oss << description_of_selectors();
65 throw Error(
"JHTopTagger can only be applied on jets having an associated (and valid) ClusterSequence");
71 _warnings_nonca.warn(
"JHTopTagger should only be applied on jets from a Cambridge/Aachen clustering; use it with other algorithms at your own risk.");
75 vector<PseudoJet> split0 = _split_once(jet, jet);
76 if (split0.size() == 0)
return PseudoJet();
79 vector<PseudoJet> subjets;
80 for (
unsigned i = 0; i < 2; i++) {
81 vector<PseudoJet> split1 = _split_once(split0[i], jet);
82 if (split1.size() > 0) {
83 subjets.push_back(split1[0]);
84 subjets.push_back(split1[1]);
86 subjets.push_back(split0[i]);
91 if (subjets.size() < 3)
return PseudoJet();
94 double dmW_min = numeric_limits<double>::max();
96 for (
unsigned i = 0 ; i < subjets.size()-1; i++) {
97 for (
unsigned j = i+1 ; j < subjets.size(); j++) {
98 double dmW = abs(_mW - (subjets[i]+subjets[j]).m());
100 dmW_min = dmW; ii = i; jj = j;
110 if (ii>0) std::swap(subjets[ii], subjets[0]);
111 if (jj>1) std::swap(subjets[jj], subjets[1]);
112 if (subjets[0].perp2() < subjets[1].perp2()) std::swap(subjets[0], subjets[1]);
113 if ((subjets.size()>3) && (subjets[2].perp2() < subjets[3].perp2()))
114 std::swap(subjets[2], subjets[3]);
120 PseudoJet W = join(subjets[0], subjets[1], *rec);
122 if (subjets.size()>3) {
123 non_W = join(subjets[2], subjets[3], *rec);
125 non_W = join(subjets[2], *rec);
127 PseudoJet result_local = join<JHTopTaggerStructure>(W, non_W, *rec);
129 s->_cos_theta_w = _cos_theta_W(result_local);
137 if (s->_cos_theta_w >= _cos_theta_W_max ||
138 ! _top_selector.pass(result_local) || ! _W_selector.pass(W)
169 vector<PseudoJet> JHTopTagger::_split_once(
const PseudoJet & jet_to_split,
173 vector<PseudoJet> result_local;
175 if (p2.
perp2() > p1.
perp2()) std::swap(p1,p2);
176 if (p1.
perp() < _delta_p * reference_jet.
perp())
break;
178 if (p2.
perp() < _delta_p * reference_jet.
perp()) {
183 result_local.push_back(p1);
184 result_local.push_back(p2);
193 FASTJET_END_NAMESPACE
double rap() const
returns the rapidity or some large value when the rapidity is infinite
JetAlgorithm jet_algorithm() const
return information about the definition...
PseudoJetStructureBase * structure_non_const_ptr()
return a non-const pointer to the structure (of type PseudoJetStructureBase*) associated with this Ps...
const ClusterSequence * validated_cs() const
shorthand for validated_cluster_sequence()
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...
const JetDefinition & jet_def() const
return a reference to the jet definition
bool has_valid_cluster_sequence() const
returns true if this PseudoJet has an associated and still valid(ated) ClusterSequence.
base class corresponding to errors that can be thrown by FastJet
const Recombiner * recombiner() const
return a pointer to the currently defined recombiner.
const ClusterSequence * associated_cluster_sequence() const
get a (const) pointer to the parent ClusterSequence (NULL if inexistent)
double perp() const
returns the scalar transverse momentum
the longitudinally invariant variant of the cambridge algorithm (aka Aachen algoithm).
An abstract base class that will provide the recombination scheme facilities and/or allow a user to e...
double delta_phi_to(const PseudoJet &other) const
returns other.phi() - this.phi(), constrained to be in range -pi .
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
double perp2() const
returns the squared transverse momentum
the structure returned by the JHTopTagger transformer.