1 | //STARTHEADER |
---|
2 | // $Id: RestFrameNSubjettinessTagger.cc 999 2013-03-04 11:48:06Z pavel $ |
---|
3 | // |
---|
4 | // Copyright (c) 2005-2011, 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 and are described in hep-ph/0512210. If you use |
---|
16 | // FastJet as part of work towards a scientific publication, please |
---|
17 | // include a citation to the FastJet paper. |
---|
18 | // |
---|
19 | // FastJet is distributed in the hope that it will be useful, |
---|
20 | // but WITHOUT ANY WARRANTY; without even the implied warranty of |
---|
21 | // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
---|
22 | // GNU General Public License for more details. |
---|
23 | // |
---|
24 | // You should have received a copy of the GNU General Public License |
---|
25 | // along with FastJet. If not, see <http://www.gnu.org/licenses/>. |
---|
26 | //---------------------------------------------------------------------- |
---|
27 | //ENDHEADER |
---|
28 | |
---|
29 | #include <fastjet/tools/RestFrameNSubjettinessTagger.hh> |
---|
30 | #include <fastjet/tools/Boost.hh> |
---|
31 | #include <fastjet/ClusterSequence.hh> |
---|
32 | #include <sstream> |
---|
33 | |
---|
34 | using namespace std; |
---|
35 | |
---|
36 | FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh |
---|
37 | |
---|
38 | //------------------------------------------------------------------------ |
---|
39 | // RestFrameNSubjettinessTagger class implementation |
---|
40 | //------------------------------------------------------------------------ |
---|
41 | |
---|
42 | //------------------------------------------------------------------------ |
---|
43 | // tagger description |
---|
44 | string RestFrameNSubjettinessTagger::description() const{ |
---|
45 | ostringstream oss; |
---|
46 | oss << "RestFrameNSubjettiness tagger that performs clustering in the jet rest frame with " |
---|
47 | << _subjet_def.description() |
---|
48 | << ", supplemented with cuts tau_2 < " << _t2cut |
---|
49 | << " and cos(theta_s) < " << _costscut; |
---|
50 | return oss.str(); |
---|
51 | } |
---|
52 | |
---|
53 | |
---|
54 | //------------------------------------------------------------------------ |
---|
55 | // action on a single jet |
---|
56 | // returns the tagged PseudoJet if successful, 0 otherwise |
---|
57 | PseudoJet RestFrameNSubjettinessTagger::result(const PseudoJet & jet) const{ |
---|
58 | // make sure that the jet has constituents |
---|
59 | if (!jet.has_constituents()) |
---|
60 | throw("The jet you try to tag needs to have accessible constituents"); |
---|
61 | |
---|
62 | // get the constituents and boost them into the rest frame of the jet |
---|
63 | vector<PseudoJet> rest_input = jet.constituents(); |
---|
64 | for (unsigned int i=0; i<rest_input.size(); i++) |
---|
65 | rest_input[i].unboost(jet); |
---|
66 | |
---|
67 | ClusterSequence cs_rest(rest_input, _subjet_def); |
---|
68 | vector<PseudoJet> subjets = (_use_exclusive) |
---|
69 | ? cs_rest.exclusive_jets(2) |
---|
70 | : sorted_by_E(cs_rest.inclusive_jets()); |
---|
71 | |
---|
72 | // impose the cuts in the rest-frame |
---|
73 | if (subjets.size()<2) return PseudoJet(); |
---|
74 | |
---|
75 | const PseudoJet &j0 = subjets[0]; |
---|
76 | const PseudoJet &j1 = subjets[1]; |
---|
77 | |
---|
78 | /// impose the cut on cos(theta_s) |
---|
79 | double ct0 = (j0.px()*jet.px() + j0.py()*jet.py() + j0.pz()*jet.pz()) |
---|
80 | /sqrt(j0.modp2()*jet.modp2()); |
---|
81 | double ct1 = (j1.px()*jet.px() + j1.py()*jet.py() + j1.pz()*jet.pz()) |
---|
82 | /sqrt(j1.modp2()*jet.modp2()); |
---|
83 | if ((ct0 > _costscut) || (ct1 > _costscut)) return PseudoJet(); |
---|
84 | |
---|
85 | // ccompute the 2-subjettiness and impose the coresponding cut |
---|
86 | double tau2 = 0.0; |
---|
87 | for (unsigned int i=0; i<rest_input.size(); i++) |
---|
88 | tau2 += min(dot_product(rest_input[i], j0), |
---|
89 | dot_product(rest_input[i], j1)); |
---|
90 | |
---|
91 | tau2 *= (2.0/jet.m2()); |
---|
92 | |
---|
93 | if (tau2 > _t2cut) return PseudoJet(); |
---|
94 | |
---|
95 | // We have a positive tag, |
---|
96 | // - boost everything back into the lab frame |
---|
97 | // - record the info in the interface |
---|
98 | // Note that in order to point to the correct Clustersequence, the |
---|
99 | // subjets must be taken from the boosted one. We extract that |
---|
100 | // through the history index of the rest-frame subjets |
---|
101 | ClusterSequence * cs_structure = new ClusterSequence(); |
---|
102 | Boost boost(jet); |
---|
103 | cs_structure->transfer_from_sequence(cs_rest, &boost); |
---|
104 | PseudoJet subjet_lab1 = cs_structure->jets()[cs_rest.history()[subjets[0].cluster_hist_index()].jetp_index]; |
---|
105 | PseudoJet subjet_lab2 = cs_structure->jets()[cs_rest.history()[subjets[0].cluster_hist_index()].jetp_index]; |
---|
106 | |
---|
107 | PseudoJet result_local = join<StructureType>(subjet_lab1,subjet_lab2); |
---|
108 | StructureType * s = (StructureType *) result_local.structure_non_const_ptr(); |
---|
109 | // s->_original_jet = jet; |
---|
110 | s->_tau2 = tau2; |
---|
111 | s->_costhetas = max(ct0, ct1); |
---|
112 | |
---|
113 | // keep the rest-frame CS alive |
---|
114 | cs_structure->delete_self_when_unused(); |
---|
115 | |
---|
116 | return result_local; |
---|
117 | } |
---|
118 | |
---|
119 | |
---|
120 | FASTJET_END_NAMESPACE |
---|