Rivet analyses referenceCMS_2011_I890166$K_s$, $\Lambda$, and Cascade$-$ transverse momentum and rapidity spectra at 900 and 7000 GeV.Experiment: CMS (LHC) Inspire ID: 890166 Status: VALIDATED Authors:
Beam energies: (450.0, 450.0); (3500.0, 3500.0) GeV Run details:
The spectra of $K_S$, $\Lambda$, and Cascade- particles were measured versus transverse-momentum ($p_\perp$) and rapidity ($y$) in proton-proton collisions at center-of-mass energies 900 and 7000 GeV. The production is normalized to all non-single-diffractive (NSD) events using corrections for trigger and selection efficiency, acceptance, and branching ratios. The results cover a rapidity range of $|y| < 2$ and a $p_\perp$ range from 0 to 10 GeV ($K_S$ and $\Lambda$) and 0 to 6 GeV (Cascade-). Antiparticles are included in all measurements so only the sums of $\Lambda$ and $\bar{\Lambda}$, and Cascade$-$ and anti-Cascade$-$ are given. The rapidity distributions are shown versus $|y|$ but normalized to a unit of $y$. Ratios of $\Lambda$/$K_S$ and Cascade$-$/$\Lambda$ production versus $p_\perp$ and $|y|$ are also given, with somewhat smaller systematic uncertainties than obtained from taking the ratio of the individual distributions.' The data were corrected according to the SD/DD/ND content of the CMS trigger, as predicted by PYTHIA6. The uncertainties connected with correct or uncorrect modelling of diffraction were included in the systematic errors. Source code: CMS_2011_I890166.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/UnstableParticles.hh"
4
5namespace Rivet {
6
7
8 /// @brief CMS strange particle spectra (Ks, Lambda, Cascade) in pp at 900 and 7000 GeV
9 ///
10 /// @author Kevin Stenson
11 class CMS_2011_I890166 : public Analysis {
12 public:
13
14 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2011_I890166);
15
16
17 void init() {
18
19 UnstableParticles ufs(Cuts::absrap < 2);
20 declare(ufs, "UFS");
21
22 // Particle distributions versus rapidity and transverse momentum
23 for (double eVal : allowedEnergies()) {
24 const string en = toString(int(eVal));
25 if (isCompatibleWithSqrtS(eVal)) _sqs = en;
26 bool offset(en == "7000"s);
27
28 book(_h[en+"dNKshort_dy"], 1, 1, 1+offset);
29 book(_h[en+"dNKshort_dpT"], 2, 1, 1+offset);
30 book(_h[en+"dNLambda_dy"], 3, 1, 1+offset);
31 book(_h[en+"dNLambda_dpT"], 4, 1, 1+offset);
32 book(_h[en+"dNXi_dy"], 5, 1, 1+offset);
33 book(_h[en+"dNXi_dpT"], 6, 1, 1+offset);
34 //
35 book(_e[en+"LampT_KpT"], 7, 1, 1+offset);
36 book(_e[en+"XipT_LampT"], 8, 1, 1+offset);
37 book(_e[en+"Lamy_Ky"], 9, 1, 1+offset);
38 book(_e[en+"Xiy_Lamy"], 10, 1, 1+offset);
39
40 }
41 if (_sqs == "" && !merging()) {
42 throw BeamError("Invalid beam energy for " + name() + "\n");
43 }
44 }
45
46
47 void analyze(const Event& event) {
48
49 const UnstableParticles& parts = apply<UnstableParticles>(event, "UFS");
50 for (const Particle& p : parts.particles()) {
51 switch (p.abspid()) {
52 case PID::K0S:
53 _h[_sqs+"dNKshort_dy"]->fill(p.absrap());
54 _h[_sqs+"dNKshort_dpT"]->fill(p.pT()/GeV);
55 break;
56
57 case PID::LAMBDA:
58 // Lambda should not have Cascade or Omega ancestors since they should not decay. But just in case...
59 if ( !( p.hasAncestorWith(Cuts::pid == 3322) ||
60 p.hasAncestorWith(Cuts::pid == -3322) ||
61 p.hasAncestorWith(Cuts::pid == 3312) ||
62 p.hasAncestorWith(Cuts::pid == -3312) ||
63 p.hasAncestorWith(Cuts::pid == 3334) ||
64 p.hasAncestorWith(Cuts::pid == -3334) ) ) {
65 _h[_sqs+"dNLambda_dy"]->fill(p.absrap());
66 _h[_sqs+"dNLambda_dpT"]->fill(p.pT()/GeV);
67 }
68 break;
69
70 case PID::XIMINUS:
71 // Cascade should not have Omega ancestors since it should not decay. But just in case...
72 if ( !( p.hasAncestorWith(Cuts::pid == 3334) ||
73 p.hasAncestorWith(Cuts::pid == -3334) ) ) {
74 _h[_sqs+"dNXi_dy"]->fill(p.absrap());
75 _h[_sqs+"dNXi_dpT"]->fill(p.pT()/GeV);
76 }
77 break;
78 }
79 }
80 }
81
82
83 void finalize() {
84
85 for (double eVal : allowedEnergies()) {
86
87 const string en = toString(int(eVal));
88
89 divide(_h[en+"dNLambda_dpT"], _h[en+"dNKshort_dpT"], _e[en+"LampT_KpT"]);
90
91 YODA::Histo1D denom = _h[en+"dNLambda_dpT"]->clone();
92 denom.rebinXTo(_h[en+"dNXi_dpT"]->xEdges());
93 divide(*_h[en+"dNXi_dpT"], denom, _e[en+"XipT_LampT"]);
94
95 divide(_h[en+"dNLambda_dy"], _h[en+"dNKshort_dy"], _e[en+"Lamy_Ky"]);
96 divide(_h[en+"dNXi_dy"], _h[en+"dNLambda_dy"], _e[en+"Xiy_Lamy"]);
97 }
98
99 const double normpT = 1.0/sumOfWeights();
100 const double normy = 0.5*normpT; // Accounts for using |y| instead of y
101 for (auto& item : _h) {
102 if (item.first.find("_dy") != string::npos) {
103 scale(item.second, normy);
104 }
105 else {
106 scale(item.second, normpT);
107 }
108 }
109 }
110
111
112 private:
113
114 /// @name Particle distributions versus rapidity and transverse momentum
115 /// @{
116 map<string,Histo1DPtr> _h;
117 map<string,Estimate1DPtr> _e;
118
119 string _sqs = "";
120 /// @}
121
122 };
123
124
125
126 RIVET_DECLARE_ALIASED_PLUGIN(CMS_2011_I890166, CMS_2011_S8978280);
127
128}
|