Rivet analyses referenceBELLE_2019_I1693396Decay kinematics of semileptonic B0→D∗− decays.Experiment: BELLE (KEKB) Inspire ID: 1693396 Status: VALIDATED NOHEPDATA Authors:
Beam energies: ANY Run details:
Measurement of recoil w, helicity and decay plane angles of semileptonc ˉB0 to D∗+ decays. N.B. the data is not unfolded but the efficiencies and response matrices used in the paper are applied to the results of the simulation. Source code: BELLE_2019_I1693396.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/UnstableParticles.hh"
5
6namespace Rivet {
7
8
9 /// @brief B0 -> D*- semileptonic
10 class BELLE_2019_I1693396 : public Analysis {
11 public:
12
13 /// Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(BELLE_2019_I1693396);
15
16
17 /// @name Analysis methods
18 /// @{
19
20 /// Book histograms and initialise projections before the run
21 void init() {
22
23 // Initialise and register projections
24 declare(UnstableParticles(Cuts::pid==511), "UFS");
25
26 // Book histograms
27 for(unsigned int ix=0;ix<2;++ix) {
28 book(_h[0][ix], "TMP/h_w_" +toString(ix+1), refData(1, 1, ix+1));
29 book(_h[1][ix], "TMP/h_costhl_"+toString(ix+1), refData(3, 1, ix+1));
30 book(_h[2][ix], "TMP/h_costhv_"+toString(ix+1), refData(2, 1, ix+1));
31 book(_h[3][ix], "TMP/h_chi_" +toString(ix+1), refData(4, 1, ix+1));
32 }
33 }
34
35 /// Perform the per-event analysis
36 bool analyzeDecay(Particle mother, vector<int> ids) {
37 // There is no point in looking for decays with less particles than to be analysed
38 if (mother.children().size() == ids.size()) {
39 bool decayfound = true;
40 for (int id : ids) {
41 if (!contains(mother, id)) decayfound = false;
42 }
43 return decayfound;
44 }
45 return false;
46 }
47
48 bool contains(Particle& mother, int id) {
49 return any(mother.children(), HasPID(id));
50 }
51
52 double recoilW(const Particle& mother) {
53 FourMomentum lepton, neutrino, meson, q;
54 for (const Particle& c : mother.children()) {
55 if (c.isNeutrino()) neutrino = c.mom();
56 else if (c.isChargedLepton()) lepton = c.mom();
57 else if (c.isHadron()) meson = c.mom();
58 }
59 q = lepton + neutrino; //no hadron before
60 double mb2= mother.mom()*mother.mom();
61 double mD2 = meson*meson;
62 return (mb2 + mD2 - q*q )/ (2. * sqrt(mb2) * sqrt(mD2) );
63 }
64
65 /// Perform the per-event analysis
66 void analyze(const Event& event) {
67 FourMomentum pl, pnu, pB, pD, pDs, ppi;
68 // Iterate of B0bar mesons
69 for(const Particle& p : apply<UnstableParticles>(event, "UFS").particles()) {
70 pB = p.momentum();
71 // Find semileptonic decays
72 int iloc=-1;
73 if (analyzeDecay(p, {PID::DSTARMINUS,12,-11})) iloc = 0;
74 else if(analyzeDecay(p, {PID::DSTARMINUS,14,-13}) ) iloc = 1;
75 else continue;
76 _h[0][iloc]->fill(recoilW(p));
77 // Get the necessary momenta for the angles
78 bool foundDdecay=false;
79 for (const Particle & c : p.children()) {
80 if ( (c.pid() == PID::DSTARMINUS) && (analyzeDecay(c, {PID::PIMINUS, PID::D0BAR}) || analyzeDecay(c, {PID::PI0, PID::DMINUS})) ) {
81 foundDdecay=true;
82 pDs = c.momentum();
83 for (const Particle & dc : c.children()) {
84 if (dc.hasCharm()) pD = dc.momentum();
85 else ppi = dc.momentum();
86 }
87 }
88 if (c.pid() == -11 || c.pid() == -13) pl = c.momentum();
89 if (c.pid() == 12 || c.pid() == 14) pnu = c.momentum();
90 }
91 // This is the angle analysis
92 if (!foundDdecay) continue;
93 // First boost all relevant momenta into the B-rest frame
94 const LorentzTransform B_boost = LorentzTransform::mkFrameTransformFromBeta(pB.betaVec());
95 // Momenta in B rest frame:
96 FourMomentum lv_brest_Dstar = B_boost.transform(pDs);
97 FourMomentum lv_brest_w = B_boost.transform(pB - pDs);
98 FourMomentum lv_brest_D = B_boost.transform(pD);
99 FourMomentum lv_brest_lep = B_boost.transform(pl);
100
101 const LorentzTransform Ds_boost = LorentzTransform::mkFrameTransformFromBeta(lv_brest_Dstar.betaVec());
102 FourMomentum lv_Dstarrest_D = Ds_boost.transform(lv_brest_D);
103 const LorentzTransform W_boost = LorentzTransform::mkFrameTransformFromBeta(lv_brest_w.betaVec());
104 FourMomentum lv_wrest_lep = W_boost.transform(lv_brest_lep);
105
106 double cos_thetaV = cos(lv_brest_Dstar.p3().angle(lv_Dstarrest_D.p3()));
107 _h[2][iloc]->fill(cos_thetaV);
108
109 double cos_thetaL = cos(lv_brest_w.p3().angle(lv_wrest_lep.p3()));
110 _h[1][iloc]->fill(cos_thetaL);
111
112 Vector3 LTrans = lv_wrest_lep.p3() - cos_thetaL*lv_wrest_lep.p3().perp()*lv_brest_w.p3().unit();
113 Vector3 VTrans = lv_Dstarrest_D.p3() - cos_thetaV*lv_Dstarrest_D.p3().perp()*lv_brest_Dstar.p3().unit();
114 float chi = atan2(LTrans.cross(VTrans).dot(lv_brest_w.p3().unit()), LTrans.dot(VTrans));
115 _h[3][iloc]->fill(chi);
116 }
117 }
118
119
120 /// Normalise histograms etc., after the run
121 void finalize() {
122 // efficiencies
123 vector<double> eff[4][2]={{{2.72,5.72,7.7,9.1,10.03,10.61,10.74,10.67,10.23,9.1},
124 {2.68,5.66,7.66,9.05,9.91,10.43,10.6,10.52,10.04,9.14}},
125 {{3.12,3.97,5.73,7.96,9.31,9.85,10.23,10.59,11.06,11.21},
126 {3.16,3.52,5.19,7.59,9.1,9.78,10.27,10.43,11.0,11.36}},
127 {{11.72,11.52,11.35,10.88,10.2,9.34,8.29,7.16,6.05,4.82},
128 {11.54,11.43,11.14,10.74,10.09,9.29,8.25,7.1,5.97,4.72}},
129 {{8.6,8.74,8.96,9.3,9.81,9.82,9.33,9.0,8.77,8.59},
130 {8.51,8.67,8.82,9.15,9.7,9.73,9.2,8.83,8.62,8.54}}};
131 vector<double> efe[4][2]={{{0.02,0.02,0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03},
132 {0.02,0.02,0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03}},
133 {{0.03,0.02,0.02,0.03,0.03,0.03,0.03,0.03,0.03,0.03},
134 {0.03,0.02,0.02,0.03,0.03,0.03,0.03,0.03,0.03,0.03}},
135 {{0.03,0.03,0.03,0.04,0.04,0.04,0.03,0.03,0.02,0.02},
136 {0.03,0.03,0.03,0.04,0.04,0.04,0.03,0.03,0.02,0.02}},
137 {{0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03},
138 {0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03,0.03}}};
139 // response matricesdouble
140 double response[4][2][10][10] = {{{{0.803,0.053,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
141 {0.197,0.778,0.098,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
142 {0.0,0.168,0.717,0.126,0.002,0.0,0.0,0.0,0.0,0.0},
143 {0.0,0.001,0.182,0.667,0.149,0.006,0.0,0.0,0.0,0.0},
144 {0.0,0.0,0.004,0.199,0.626,0.167,0.011,0.0,0.0,0.0},
145 {0.0,0.0,0.0,0.009,0.207,0.592,0.177,0.015,0.0,0.0},
146 {0.0,0.0,0.0,0.0,0.016,0.215,0.575,0.183,0.018,0.0},
147 {0.0,0.0,0.0,0.0,0.0,0.021,0.213,0.567,0.186,0.017},
148 {0.0,0.0,0.0,0.0,0.0,0.0,0.024,0.214,0.598,0.186},
149 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.022,0.198,0.797}},
150 {{0.961,0.024,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
151 {0.038,0.952,0.027,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
152 {0.0,0.021,0.948,0.041,0.001,0.0,0.0,0.0,0.0,0.0},
153 {0.0,0.001,0.023,0.918,0.067,0.003,0.001,0.001,0.001,0.0},
154 {0.0,0.001,0.001,0.04,0.871,0.097,0.005,0.001,0.001,0.0},
155 {0.0,0.0,0.0,0.001,0.06,0.817,0.129,0.006,0.001,0.0},
156 {0.0,0.0,0.0,0.0,0.001,0.082,0.758,0.164,0.007,0.001},
157 {0.0,0.0,0.0,0.0,0.0,0.001,0.106,0.698,0.196,0.008},
158 {0.0,0.0,0.0,0.0,0.0,0.0,0.001,0.128,0.657,0.212},
159 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.002,0.137,0.777}}},
160 {{{0.918,0.077,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
161 {0.082,0.806,0.095,0.001,0.0,0.0,0.0,0.0,0.0,0.0},
162 {0.0,0.115,0.761,0.101,0.002,0.0,0.0,0.0,0.0,0.0},
163 {0.0,0.001,0.141,0.735,0.105,0.002,0.0,0.0,0.0,0.0},
164 {0.0,0.0,0.002,0.16,0.719,0.1,0.001,0.0,0.0,0.0},
165 {0.0,0.0,0.0,0.003,0.17,0.722,0.093,0.001,0.0,0.0},
166 {0.0,0.0,0.0,0.0,0.003,0.173,0.738,0.08,0.001,0.0},
167 {0.0,0.0,0.0,0.0,0.0,0.002,0.166,0.771,0.072,0.0},
168 {0.0,0.0,0.0,0.0,0.0,0.0,0.001,0.147,0.819,0.064},
169 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.001,0.108,0.936}},
170 {{0.659,0.129,0.011,0.003,0.002,0.002,0.002,0.004,0.013,0.144},
171 {0.151,0.691,0.132,0.012,0.004,0.002,0.002,0.002,0.004,0.016},
172 {0.015,0.141,0.697,0.147,0.016,0.005,0.002,0.002,0.002,0.005},
173 {0.005,0.012,0.134,0.671,0.162,0.018,0.005,0.002,0.002,0.002},
174 {0.002,0.004,0.013,0.14,0.634,0.155,0.016,0.004,0.002,0.002},
175 {0.002,0.002,0.004,0.015,0.155,0.633,0.141,0.013,0.004,0.003},
176 {0.002,0.002,0.002,0.004,0.018,0.163,0.67,0.136,0.012,0.004},
177 {0.005,0.002,0.002,0.002,0.005,0.015,0.147,0.695,0.14,0.015},
178 {0.016,0.004,0.002,0.002,0.002,0.004,0.013,0.132,0.691,0.15},
179 {0.142,0.013,0.003,0.002,0.002,0.002,0.003,0.012,0.13,0.659}}},
180 {{{0.812,0.051,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
181 {0.188,0.784,0.096,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
182 {0.0,0.164,0.728,0.126,0.002,0.0,0.0,0.0,0.0,0.0},
183 {0.0,0.001,0.172,0.676,0.149,0.006,0.0,0.0,0.0,0.0},
184 {0.0,0.0,0.004,0.19,0.631,0.165,0.01,0.0,0.0,0.0},
185 {0.0,0.0,0.0,0.008,0.203,0.6,0.181,0.016,0.0,0.0},
186 {0.0,0.0,0.0,0.0,0.014,0.209,0.578,0.187,0.019,0.0},
187 {0.0,0.0,0.0,0.0,0.0,0.02,0.209,0.573,0.195,0.017},
188 {0.0,0.0,0.0,0.0,0.0,0.0,0.022,0.205,0.6,0.195},
189 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.019,0.186,0.788}},
190 {{0.959,0.022,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
191 {0.039,0.955,0.012,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
192 {0.0,0.021,0.96,0.022,0.001,0.0,0.0,0.0,0.0,0.0},
193 {0.001,0.001,0.026,0.931,0.043,0.001,0.0,0.0,0.0,0.0},
194 {0.0,0.0,0.001,0.047,0.889,0.07,0.002,0.0,0.0,0.0},
195 {0.0,0.001,0.0,0.0,0.067,0.837,0.103,0.002,0.001,0.0},
196 {0.0,0.0,0.0,0.0,0.0,0.091,0.778,0.138,0.003,0.0},
197 {0.0,0.0,0.0,0.0,0.0,0.0,0.117,0.715,0.174,0.004},
198 {0.0,0.0,0.0,0.0,0.0,0.0,0.001,0.142,0.672,0.193},
199 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.002,0.151,0.803}}},
200 {{{0.918,0.077,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0},
201 {0.082,0.805,0.091,0.001,0.0,0.0,0.0,0.0,0.0,0.0},
202 {0.0,0.117,0.763,0.101,0.002,0.0,0.0,0.0,0.0,0.0},
203 {0.0,0.001,0.142,0.735,0.103,0.002,0.0,0.0,0.0,0.0},
204 {0.0,0.0,0.003,0.159,0.723,0.098,0.001,0.0,0.0,0.0},
205 {0.0,0.0,0.0,0.004,0.169,0.726,0.091,0.001,0.0,0.0},
206 {0.0,0.0,0.0,0.0,0.004,0.172,0.745,0.082,0.001,0.0},
207 {0.0,0.0,0.0,0.0,0.0,0.002,0.161,0.771,0.074,0.0},
208 {0.0,0.0,0.0,0.0,0.0,0.0,0.001,0.145,0.817,0.066},
209 {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.107,0.934}},
210 {{0.653,0.129,0.012,0.004,0.003,0.002,0.002,0.004,0.014,0.144},
211 {0.152,0.686,0.13,0.013,0.004,0.003,0.002,0.002,0.005,0.017},
212 {0.016,0.143,0.693,0.147,0.016,0.006,0.003,0.002,0.003,0.005},
213 {0.005,0.013,0.138,0.667,0.16,0.018,0.005,0.002,0.002,0.003},
214 {0.003,0.004,0.013,0.142,0.63,0.156,0.015,0.004,0.002,0.002},
215 {0.002,0.002,0.004,0.015,0.158,0.629,0.142,0.013,0.004,0.003},
216 {0.003,0.002,0.002,0.005,0.018,0.164,0.667,0.138,0.013,0.005},
217 {0.005,0.003,0.002,0.003,0.006,0.016,0.148,0.692,0.141,0.016},
218 {0.017,0.004,0.002,0.002,0.003,0.005,0.013,0.131,0.686,0.152},
219 {0.144,0.014,0.004,0.002,0.002,0.002,0.004,0.012,0.129,0.654}}}};
220 // correct the values
221 for(unsigned int ix=0;ix<4;++ix) {
222 for(unsigned int iy=0;iy<2;++iy) {
223 Estimate1DPtr corrected;
224 book(corrected,ix+1,1,iy+1);
225 // first extract values and errors applying efficiency
226 Vector<10> val,err;
227 for(unsigned int ibin=0;ibin<_h[ix][iy]->bins().size();++ibin) {
228 val[ibin] = eff[ix][iy][ibin]/100. * _h[ix][iy]->bins()[ibin+1].sumW();
229 err[ibin] = sqr(eff[ix][iy][ibin]/100. * _h[ix][iy]->bins()[ibin+1].errW());
230 sqr(efe[ix][iy][ibin]/100. * _h[ix][iy]->bins()[ibin+1].sumW());
231 }
232 // put response into a matrix
233 Matrix<10> R,R2;
234 for(unsigned int i1=0;i1<10;++i1) {
235 for(unsigned int i2=0;i2<10;++i2) {
236 R .set(i1,i2, response[ix][iy][i1][i2]);
237 R2.set(i1,i2,sqr(response[ix][iy][i1][i2]));
238 }
239 }
240 // multiply to get value and error^2
241 val = multiply(R ,val);
242 err = multiply(R2,err);
243 // total for normalization
244 double total=0.;
245 for(unsigned int i1=0;i1<10;++i1) total+=val[i1];
246 // finally the output scatter
247 for(unsigned int ibin=0;ibin<_h[ix][iy]->bins().size();++ibin) {
248 double dx = 0.5*_h[ix][iy]->bins()[ibin+1].xWidth();
249 double dy = sqrt(err[ibin])/total/2./dx;
250 corrected->bin(ibin+1).set(val[ibin]/total/2./dx, dy);
251 }
252 }
253 }
254 }
255
256 /// @}
257
258
259 /// @name Histograms
260 /// @{
261 Histo1DPtr _h[4][2];
262 /// @}
263
264
265 };
266
267
268 RIVET_DECLARE_PLUGIN(BELLE_2019_I1693396);
269
270}
|