Rivet analyses referenceLHCB_2012_I1097092$B_c^+\to J/\psi\pi^+$ and $B_c^+\to J/\psi\pi^+\pi^+\pi^-$Experiment: LHCB (LHC) Inspire ID: 1097092 Status: VALIDATED NOHEPDATA Authors:
Beam energies: ANY Run details:
Angular distributions in between the muons $B_c^+\to J/\psi\pi^+$ and $B_c^+\to J/\psi\pi^+\pi^+\pi^-$ with the muons from $J/\psi\to\mu^+\mu^-$ and mass distributions in $B_c^+\to J/\psi\pi^+\pi^+\pi^-$. The data were read from the plots in the paper and may not be corrected for efficiency and acceptance. Source code: LHCB_2012_I1097092.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/UnstableParticles.hh"
4#include "Rivet/Projections/DecayedParticles.hh"
5
6namespace Rivet {
7
8
9 /// @brief B_c - > jpsi pi+ and pi+pi+pi-
10 class LHCB_2012_I1097092 : public Analysis {
11 public:
12
13 /// Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(LHCB_2012_I1097092);
15
16
17 /// @name Analysis methods
18 /// @{
19
20 /// Book histograms and initialise projections before the run
21 void init() {
22 UnstableParticles ufs = UnstableParticles(Cuts::abspid==541);
23 declare(ufs, "UFS");
24 DecayedParticles BC(ufs);
25 BC.addStable( PID::PI0);
26 BC.addStable( PID::K0S);
27 BC.addStable( PID::JPSI);
28 declare(BC, "BC");
29 for(unsigned int ix=0;ix<2;++ix) {
30 book(_h_mass [ix],1+ix,1,1);
31 book(_h_ctheta[ix],3,1,1+ix);
32 }
33 }
34
35
36 /// Perform the per-event analysis
37 void analyze(const Event& event) {
38 static const map<PdgId,unsigned int> & mode1 = { { 443,1}, { 211,1} };
39 static const map<PdgId,unsigned int> & mode1CC = { { 443,1}, {-211,1} };
40 static const map<PdgId,unsigned int> & mode2 = { { 443,1}, { 211,2}, {-211,1} };
41 static const map<PdgId,unsigned int> & mode2CC = { { 443,1}, {-211,2}, { 211,1} };
42 DecayedParticles BC = apply<DecayedParticles>(event, "BC");
43 // loop over particles
44 for(unsigned int ix=0;ix<BC.decaying().size();++ix) {
45 int sign = BC.decaying()[ix].pid()/BC.decaying()[ix].abspid();
46 int imode=-1;
47 if ((sign== 1 && BC.modeMatches(ix,2,mode1 )) ||
48 (sign==-1 && BC.modeMatches(ix,2,mode1CC))) {
49 imode=0;
50 }
51 else if ((sign== 1 && BC.modeMatches(ix,4,mode2 )) ||
52 (sign==-1 && BC.modeMatches(ix,4,mode2CC))) {
53 const Particle & pim = BC.decayProducts()[ix].at(-sign*211)[0];
54 const Particles & pip = BC.decayProducts()[ix].at( sign*211);
55 _h_mass[0]->fill((pim.momentum()+pip[0].momentum()+pip[1].momentum()).mass());
56 _h_mass[1]->fill((pim.momentum()+pip[0].momentum()).mass());
57 _h_mass[1]->fill((pim.momentum()+pip[1].momentum()).mass());
58 imode=1;
59 }
60 else
61 continue;
62 // check the J/psi decay mode
63 const Particle & jpsi = BC.decayProducts()[ix].at(443)[0];
64 if(jpsi.children().size()!=2 || jpsi.children()[0].pid()!=-jpsi.children()[1].pid() ||
65 jpsi.children()[0].abspid()!=13) continue;
66 Particle muon = jpsi.children()[0].pid()==13 ? jpsi.children()[0] : jpsi.children()[1];
67 LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(BC.decaying()[ix].momentum().betaVec());
68 FourMomentum pJPsi = boost1.transform(jpsi.momentum());
69 FourMomentum pMu = boost1.transform(muon.momentum());
70 // to j/psi rest frame
71 LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pJPsi.betaVec());
72 Vector3 axis = pJPsi.p3().unit();
73 FourMomentum pp = boost2.transform(pMu);
74 // calculate angle
75 double cTheta = pp.p3().unit().dot(axis);
76 _h_ctheta[imode]->fill(cTheta);
77 }
78 }
79
80
81 /// Normalise histograms etc., after the run
82 void finalize() {
83 for(unsigned int ix=0;ix<2;++ix) {
84 normalize(_h_mass [ix]);
85 normalize(_h_ctheta[ix]);
86 }
87 }
88
89 /// @}
90
91
92 /// @name Histograms
93 /// @{
94 Histo1DPtr _h_mass[2];
95 Histo1DPtr _h_ctheta[2];
96 /// @}
97
98
99 };
100
101
102 RIVET_DECLARE_PLUGIN(LHCB_2012_I1097092);
103
104}
|