1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
| // -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/FastJets.hh"
#define I_KNOW_THE_INITIAL_QUARKS_PROJECTION_IS_DODGY_BUT_NEED_TO_USE_IT
#include "Rivet/Projections/InitialQuarks.hh"
namespace Rivet {
/// Studies of 2 b-jet + 2 jet production in proton-proton collisions at 7 TeV
class CMS_2016_I1486238 : public Analysis {
public:
/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2016_I1486238);
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
FastJets akt(FinalState(), FastJets::ANTIKT, 0.5);
declare(akt, "antikT");
book(_h_Deltaphi_newway, 1,1,1);
book(_h_deltaphiafterlight, 9,1,1);
book(_h_SumPLight, 5,1,1);
book(_h_LeadingBJetpt, 11,1,1);
book(_h_SubleadingBJetpt, 15,1,1);
book(_h_LeadingLightJetpt, 13,1,1);
book(_h_SubleadingLightJetpt, 17,1,1);
book(_h_LeadingBJeteta, 10,1,1);
book(_h_SubleadingBJeteta, 14,1,1);
book(_h_LeadingLightJeteta, 12,1,1);
book(_h_SubleadingLightJeteta, 16,1,1);
}
/// Perform the per-event analysis
void analyze(const Event& event) {
const Jets& jets = apply<JetAlg>(event, "antikT").jetsByPt(Cuts::absrap < 4.7 && Cuts::pT > 20*GeV);
if (jets.size() < 4) vetoEvent;
// Initial quarks
/// @note Quark-level tagging...
Particles bquarks;
for (ConstGenParticlePtr p : HepMCUtils::particles(event.genEvent())) {
if (abs(p->pdg_id()) == PID::BQUARK) bquarks += Particle(p);
}
Jets bjets, ljets;
for (const Jet& j : jets) {
const bool btag = any(bquarks, deltaRLess(j, 0.3));
// for (const Particle& b : bquarks) if (deltaR(j, b) < 0.3) btag = true;
(btag && j.abseta() < 2.4 ? bjets : ljets).push_back(j);
}
// Fill histograms
const double weight = 1.0;
if (bjets.size() >= 2 && ljets.size() >= 2) {
_h_LeadingBJetpt->fill(bjets[0].pT()/GeV, weight);
_h_SubleadingBJetpt->fill(bjets[1].pT()/GeV, weight);
_h_LeadingLightJetpt->fill(ljets[0].pT()/GeV, weight);
_h_SubleadingLightJetpt->fill(ljets[1].pT()/GeV, weight);
//
_h_LeadingBJeteta->fill(bjets[0].eta(), weight);
_h_SubleadingBJeteta->fill(bjets[1].eta(), weight);
_h_LeadingLightJeteta->fill(ljets[0].eta(), weight);
_h_SubleadingLightJeteta->fill(ljets[1].eta(), weight);
const double lightdphi = deltaPhi(ljets[0], ljets[1]);
_h_deltaphiafterlight->fill(lightdphi, weight);
const double vecsumlightjets = sqrt(sqr(ljets[0].px()+ljets[1].px()) + sqr(ljets[0].py()+ljets[1].py())); //< @todo Just (lj0+lj1).pT()? Or use add_quad
const double term2 = vecsumlightjets/(sqrt(sqr(ljets[0].px()) + sqr(ljets[0].py())) + sqrt(sqr(ljets[1].px()) + sqr(ljets[1].py()))); //< @todo lj0.pT() + lj1.pT()? Or add_quad
_h_SumPLight->fill(term2, weight);
const double pxBsyst2 = bjets[0].px()+bjets[1].px(); // @todo (bj0+bj1).px()
const double pyBsyst2 = bjets[0].py()+bjets[1].py(); // @todo (bj0+bj1).py()
const double pxJetssyst2 = ljets[0].px()+ljets[1].px(); // @todo (lj0+lj1).px()
const double pyJetssyst2 = ljets[0].py()+ljets[1].py(); // @todo (lj0+lj1).py()
const double modulusB2 = sqrt(sqr(pxBsyst2)+sqr(pyBsyst2)); //< @todo add_quad
const double modulusJets2 = sqrt(sqr(pxJetssyst2)+sqr(pyJetssyst2)); //< @todo add_quad
const double cosphiBsyst2 = pxBsyst2/modulusB2;
const double cosphiJetssyst2 = pxJetssyst2/modulusJets2;
const double phiBsyst2 = ((pyBsyst2 > 0) ? 1 : -1) * acos(cosphiBsyst2); //< @todo sign(pyBsyst2)
const double phiJetssyst2 = sign(pyJetssyst2) * acos(cosphiJetssyst2);
const double Dphi2 = deltaPhi(phiBsyst2, phiJetssyst2);
_h_Deltaphi_newway->fill(Dphi2,weight);
}
}
/// Normalise histograms etc., after the run
void finalize() {
const double invlumi = crossSection()/picobarn/sumOfWeights();
normalize(_h_SumPLight); normalize(_h_deltaphiafterlight); normalize(_h_Deltaphi_newway);
scale(_h_LeadingLightJetpt, invlumi); scale(_h_SubleadingLightJetpt, invlumi);
scale(_h_LeadingBJetpt, invlumi); scale(_h_SubleadingBJetpt, invlumi);
scale(_h_LeadingLightJeteta, invlumi); scale(_h_SubleadingLightJeteta, invlumi);
scale(_h_LeadingBJeteta, invlumi); scale(_h_SubleadingBJeteta, invlumi);
}
//@}
private:
/// @name Histograms
//@{
Histo1DPtr _h_deltaphiafterlight, _h_Deltaphi_newway, _h_SumPLight;
Histo1DPtr _h_LeadingBJetpt, _h_SubleadingBJetpt, _h_LeadingLightJetpt, _h_SubleadingLightJetpt;
Histo1DPtr _h_LeadingBJeteta, _h_SubleadingBJeteta, _h_LeadingLightJeteta, _h_SubleadingLightJeteta;
};
// Hook for the plugin system
RIVET_DECLARE_PLUGIN(CMS_2016_I1486238);
}
|