CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
BPHHistoSpecificDecay.cc
Go to the documentation of this file.
2 
4 
6 
8 
11 
14 
16 
18 #include "TMath.h"
19 #include "Math/VectorUtil.h"
20 #include "TVector3.h"
21 
22 #include <TH1.h>
23 #include <TFile.h>
24 
25 #include <set>
26 #include <string>
27 #include <iostream>
28 #include <fstream>
29 #include <cmath>
30 
31 using namespace std;
32 
33 #define SET_LABEL(NAME, PSET) (NAME = PSET.getParameter<string>(#NAME))
34 // SET_LABEL(xyz,ps);
35 // is equivalent to
36 // xyz = ps.getParameter<string>( "xyx" )
37 
38 class BPHUserData {
39 public:
40  template <class T>
41  static const T* get(const pat::CompositeCandidate& cand, const string& name) {
42  if (cand.hasUserData(name))
43  return cand.userData<T>(name);
44  return nullptr;
45  }
46  template <class T>
47  static const T* getByRef(const pat::CompositeCandidate& cand, const string& name) {
48  if (cand.hasUserData(name)) {
49  typedef edm::Ref<std::vector<T> > objRef;
50  const objRef* ref = cand.userData<objRef>(name);
51  if (ref == nullptr)
52  return nullptr;
53  if (ref->isNull())
54  return nullptr;
55  return ref->get();
56  }
57  return nullptr;
58  }
59 };
60 
61 class BPHDaughters {
62 public:
63  static vector<const reco::Candidate*> get(const pat::CompositeCandidate& cand, float massMin, float massMax) {
64  int i;
65  int n = cand.numberOfDaughters();
66  vector<const reco::Candidate*> v;
67  v.reserve(n);
68  for (i = 0; i < n; ++i) {
69  const reco::Candidate* dptr = cand.daughter(i);
70  float mass = dptr->mass();
71  if ((mass > massMin) && (mass < massMax))
72  v.push_back(dptr);
73  }
74  return v;
75  }
76 };
77 
79 public:
80  BPHSoftMuonSelect(int cutTrackerLayers = 5,
81  int cutPixelLayers = 0,
82  float maxDxy = 0.3,
83  float maxDz = 20.0,
84  bool goodMuon = true,
85  bool highPurity = true)
86  : cutTL(cutTrackerLayers), cutPL(cutPixelLayers), maxXY(maxDxy), maxZ(maxDz), gM(goodMuon), hP(highPurity) {}
87  bool accept(const reco::Candidate& cand, const reco::Vertex* pv) const {
88  const pat::Muon* p = reinterpret_cast<const pat::Muon*>(&cand);
89  if (p == nullptr)
90  return false;
92  return false;
93  if (p->innerTrack()->hitPattern().trackerLayersWithMeasurement() <= cutTL)
94  return false;
95  if (p->innerTrack()->hitPattern().pixelLayersWithMeasurement() <= cutPL)
96  return false;
97  if (hP && !p->innerTrack()->quality(reco::TrackBase::highPurity))
98  return false;
99  if (pv == nullptr)
100  return true;
101  const reco::Vertex::Point& pos = pv->position();
102  if (fabs(p->innerTrack()->dxy(pos)) >= maxXY)
103  return false;
104  if (fabs(p->innerTrack()->dz(pos)) >= maxZ)
105  return false;
106  return true;
107  }
108 
109 private:
110  const reco::Vertex* pv;
111  int cutTL;
112  int cutPL;
113  float maxXY;
114  float maxZ;
115  bool gM;
116  bool hP;
117 };
118 
120 public:
121  BPHDaughterSelect(float ptMinLoose,
122  float ptMinTight,
123  float etaMaxLoose,
124  float etaMaxTight,
125  const BPHSoftMuonSelect* softMuonselector = nullptr)
126  : pLMin(ptMinLoose), pTMin(ptMinTight), eLMax(etaMaxLoose), eTMax(etaMaxTight), sms(softMuonselector) {}
127  bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pv = nullptr) const override {
128  const reco::Candidate* dptr0 = cand.daughter(0);
129  const reco::Candidate* dptr1 = cand.daughter(1);
130  if (dptr0 == nullptr)
131  return false;
132  if (dptr1 == nullptr)
133  return false;
134  float pt0 = dptr0->pt();
135  float pt1 = dptr1->pt();
136  if ((pt0 < pLMin) || (pt1 < pLMin))
137  return false;
138  if ((pt0 < pTMin) && (pt1 < pTMin))
139  return false;
140  float eta0 = fabs(dptr0->eta());
141  float eta1 = fabs(dptr1->eta());
142  if ((eLMax > 0) && ((eta0 > eLMax) || (eta1 > eLMax)))
143  return false;
144  if ((eTMax > 0) && ((eta0 > eTMax) && (eta1 > eTMax)))
145  return false;
146  if (sms != nullptr) {
147  const reco::Vertex* pvtx = BPHUserData::getByRef<reco::Vertex>(cand, "primaryVertex");
148  if (pvtx == nullptr)
149  return false;
150  if (!sms->accept(*dptr0, pvtx))
151  return false;
152  if (!sms->accept(*dptr1, pvtx))
153  return false;
154  }
155  return true;
156  }
157 
158 private:
159  float pLMin;
160  float pTMin;
161  float eLMax;
162  float eTMax;
164 };
165 
167 public:
169  float massMin, float massMax, float ptMin = -1.0, float etaMax = -1.0, float rapidityMax = -1.0)
170  : mMin(massMin), mMax(massMax), pMin(ptMin), eMax(etaMax), yMax(rapidityMax) {}
171  bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pv = nullptr) const override {
172  if (((mMin > 0) && (mMax < 0)) || ((mMin < 0) && (mMax > 0)) || ((mMin > 0) && (mMax > 0) && (mMin < mMax))) {
173  float mass = cand.mass();
174  if (mass < mMin)
175  return false;
176  if ((mMax > 0) && (mass > mMax))
177  return false;
178  }
179  if (cand.pt() < pMin)
180  return false;
181  if ((eMax > 0) && (fabs(cand.eta()) > eMax))
182  return false;
183  if ((yMax > 0) && (fabs(cand.rapidity()) > yMax))
184  return false;
185  return true;
186  }
187 
188 private:
189  float mMin;
190  float mMax;
191  float pMin;
192  float eMax;
193  float yMax;
194 };
195 
197 public:
198  BPHFittedBasicSelect(float massMin, float massMax, float ptMin = -1.0, float etaMax = -1.0, float rapidityMax = -1.0)
199  : mMin(massMin), mMax(massMax), pMin(ptMin), eMax(etaMax), rMax(rapidityMax) {}
200  bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pv = nullptr) const override {
201  if (!cand.hasUserFloat("fitMass"))
202  return false;
203  float mass = cand.userFloat("fitMass");
204  if (((mMin > 0) && (mMax < 0)) || ((mMin < 0) && (mMax > 0)) || ((mMin > 0) && (mMax > 0) && (mMin < mMax))) {
205  if (mass < mMin)
206  return false;
207  if ((mMax > 0) && (mass > mMax))
208  return false;
209  }
210  const Vector3DBase<float, GlobalTag>* fmom = BPHUserData::get<Vector3DBase<float, GlobalTag> >(cand, "fitMomentum");
211  if (fmom == nullptr)
212  return false;
213  if (pMin > 0) {
214  if (fmom->transverse() < pMin)
215  return false;
216  }
217  if (eMax > 0) {
218  if (fabs(fmom->eta()) > eMax)
219  return false;
220  }
221  if (rMax > 0) {
222  float x = fmom->x();
223  float y = fmom->y();
224  float z = fmom->z();
225  float e = sqrt((x * x) + (y * y) + (z * z) + (mass * mass));
226  float r = log((e + z) / (e - z)) / 2;
227  if (fabs(r) > rMax)
228  return false;
229  }
230  return true;
231  }
232 
233 private:
234  float mMin;
235  float mMax;
236  float pMin;
237  float eMax;
238  float rMax;
239 };
240 
242 public:
243  BPHCompositeVertexSelect(float probMin, float cosMin = -1.0, float sigMin = -1.0)
244  : pMin(probMin), cMin(cosMin), sMin(sigMin) {}
245  bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pvtx = nullptr) const override {
246  const reco::Vertex* svtx = BPHUserData::get<reco::Vertex>(cand, "vertex");
247  if (svtx == nullptr)
248  return false;
249  if (pvtx == nullptr)
250  return false;
251  if (pMin > 0) {
252  if (ChiSquaredProbability(svtx->chi2(), svtx->ndof()) < pMin)
253  return false;
254  }
255  if ((cMin > 0) || (sMin > 0)) {
256  TVector3 disp(svtx->x() - pvtx->x(), svtx->y() - pvtx->y(), 0);
257  TVector3 cmom(cand.px(), cand.py(), 0);
258  float cosAlpha = disp.Dot(cmom) / (disp.Perp() * cmom.Perp());
259  if (cosAlpha < cMin)
260  return false;
261  if (sMin < 0)
262  return true;
263  float mass = cand.mass();
264  AlgebraicVector3 vmom(cand.px(), cand.py(), 0);
265  VertexDistanceXY vdistXY;
266  Measurement1D distXY = vdistXY.distance(*svtx, *pvtx);
267  double ctauPV = distXY.value() * cosAlpha * mass / cmom.Perp();
268  GlobalError sve = svtx->error();
269  GlobalError pve = pvtx->error();
270  AlgebraicSymMatrix33 vXYe = sve.matrix() + pve.matrix();
271  double ctauErrPV = sqrt(ROOT::Math::Similarity(vmom, vXYe)) * mass / cmom.Perp2();
272  if ((ctauPV / ctauErrPV) < sMin)
273  return false;
274  }
275  return true;
276  }
277 
278 private:
279  float pMin;
280  float cMin;
281  float sMin;
282 };
283 
285 public:
286  BPHFittedVertexSelect(float probMin, float cosMin = -1.0, float sigMin = -1.0)
287  : pMin(probMin), cMin(cosMin), sMin(sigMin) {}
288  bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pvtx) const override {
289  const reco::Vertex* svtx = BPHUserData::get<reco::Vertex>(cand, "fitVertex");
290  if (svtx == nullptr)
291  return false;
292  if (pMin > 0) {
293  if (ChiSquaredProbability(svtx->chi2(), svtx->ndof()) < pMin)
294  return false;
295  }
296  if ((cMin > 0) || (sMin > 0)) {
297  TVector3 disp(svtx->x() - pvtx->x(), svtx->y() - pvtx->y(), 0);
298  const Vector3DBase<float, GlobalTag>* fmom =
299  BPHUserData::get<Vector3DBase<float, GlobalTag> >(cand, "fitMomentum");
300  if (fmom == nullptr)
301  return false;
302  TVector3 cmom(fmom->x(), fmom->y(), 0);
303  float cosAlpha = disp.Dot(cmom) / (disp.Perp() * cmom.Perp());
304  if (cosAlpha < cMin)
305  return false;
306  if (sMin < 0)
307  return true;
308  if (!cand.hasUserFloat("fitMass"))
309  return false;
310  float mass = cand.userFloat("fitMass");
311  AlgebraicVector3 vmom(fmom->x(), fmom->y(), 0);
312  VertexDistanceXY vdistXY;
313  Measurement1D distXY = vdistXY.distance(*svtx, *pvtx);
314  double ctauPV = distXY.value() * cosAlpha * mass / cmom.Perp();
315  GlobalError sve = svtx->error();
316  GlobalError pve = pvtx->error();
317  AlgebraicSymMatrix33 vXYe = sve.matrix() + pve.matrix();
318  double ctauErrPV = sqrt(ROOT::Math::Similarity(vmom, vXYe)) * mass / cmom.Perp2();
319  if ((ctauPV / ctauErrPV) < sMin)
320  return false;
321  }
322  return true;
323  }
324 
325 private:
326  float pMin;
327  float cMin;
328  float sMin;
329 };
330 
332  useOnia = (!SET_LABEL(oniaCandsLabel, ps).empty());
333  useSd = (!SET_LABEL(sdCandsLabel, ps).empty());
334  useSs = (!SET_LABEL(ssCandsLabel, ps).empty());
335  useBu = (!SET_LABEL(buCandsLabel, ps).empty());
336  useBd = (!SET_LABEL(bdCandsLabel, ps).empty());
337  useBs = (!SET_LABEL(bsCandsLabel, ps).empty());
338  if (useOnia)
339  consume<vector<pat::CompositeCandidate> >(oniaCandsToken, oniaCandsLabel);
340  if (useSd)
341  consume<vector<pat::CompositeCandidate> >(sdCandsToken, sdCandsLabel);
342  if (useSs)
343  consume<vector<pat::CompositeCandidate> >(ssCandsToken, ssCandsLabel);
344  if (useBu)
345  consume<vector<pat::CompositeCandidate> >(buCandsToken, buCandsLabel);
346  if (useBd)
347  consume<vector<pat::CompositeCandidate> >(bdCandsToken, bdCandsLabel);
348  if (useBs)
349  consume<vector<pat::CompositeCandidate> >(bsCandsToken, bsCandsLabel);
350 
351  static const BPHSoftMuonSelect sms;
352 
353  double phiMassMin = 0.85;
354  double phiMassMax = 3.30;
355  double phiPtMin = 16.0;
356  double phiEtaMax = -1.0;
357  double phiRMax = -1.0;
358  double jPsiMassMin = 2.95;
359  double jPsiMassMax = 3.30;
360  double jPsiPtMin = 16.0;
361  double jPsiEtaMax = -1.0;
362  double jPsiRMax = -1.0;
363  double psi2MassMin = 3.40;
364  double psi2MassMax = 4.00;
365  double psi2PtMin = 13.0;
366  double psi2EtaMax = -1.0;
367  double psi2RMax = -1.0;
368  double upsMassMin = 8.50;
369  double upsMassMax = 11.0;
370  double upsPtMin = 13.0;
371  double upsEtaMax = -1.0;
372  double upsRMax = -1.0;
373 
374  double oniaProbMin = 0.005;
375  double oniaCosMin = -1.0;
376  double oniaSigMin = -1.0;
377 
378  double oniaMuPtMinLoose = -1.0;
379  double oniaMuPtMinTight = -1.0;
380  double oniaMuEtaMaxLoose = -1.0;
381  double oniaMuEtaMaxTight = -1.0;
382 
383  phiBasicSelect = new BPHCompositeBasicSelect(phiMassMin, phiMassMax, phiPtMin, phiEtaMax, phiRMax);
384  jPsiBasicSelect = new BPHCompositeBasicSelect(jPsiMassMin, jPsiMassMax, jPsiPtMin, jPsiEtaMax, jPsiRMax);
385  psi2BasicSelect = new BPHCompositeBasicSelect(psi2MassMin, psi2MassMax, psi2PtMin, psi2EtaMax, psi2RMax);
386  upsBasicSelect = new BPHCompositeBasicSelect(upsMassMin, upsMassMax, upsPtMin, upsEtaMax, upsRMax);
387  oniaVertexSelect = new BPHCompositeVertexSelect(oniaProbMin, oniaCosMin, oniaSigMin);
388  oniaDaughterSelect =
389  new BPHDaughterSelect(oniaMuPtMinLoose, oniaMuPtMinTight, oniaMuEtaMaxLoose, oniaMuEtaMaxTight, &sms);
390 
391  double buJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
392  double buJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
393  double buJPsiPtMin = 8.0;
394  double buJPsiEtaMax = -1.0;
395  double buJPsiRMax = -1.0;
396  double buProbMin = 0.10;
397  double buCosMin = 0.99;
398  double buSigMin = 3.0;
399  double buMuPtMinLoose = 4.0;
400  double buMuPtMinTight = 4.0;
401  double buMuEtaMaxLoose = 2.2;
402  double buMuEtaMaxTight = 2.2;
403 
404  buKPtMin = 1.6;
405 
406  buJPsiBasicSelect = new BPHCompositeBasicSelect(buJPsiMassMin, buJPsiMassMax, buJPsiPtMin, buJPsiEtaMax, buJPsiRMax);
407  buVertexSelect = new BPHFittedVertexSelect(buProbMin, buCosMin, buSigMin);
408  buJPsiDaughterSelect = new BPHDaughterSelect(buMuPtMinLoose, buMuPtMinTight, buMuEtaMaxLoose, buMuEtaMaxTight, &sms);
409 
410  double bdJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
411  double bdJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
412  double bdJPsiPtMin = 8.0;
413  double bdJPsiEtaMax = -1.0;
414  double bdJPsiRMax = -1.0;
415  double bdKx0MassMin = BPHParticleMasses::kx0Mass - 0.075;
416  double bdKx0MassMax = BPHParticleMasses::kx0Mass + 0.075;
417  double bdKx0PtMin = -1.0;
418  double bdKx0EtaMax = -1.0;
419  double bdKx0RMax = -1.0;
420  double bdProbMin = 0.10;
421  double bdCosMin = 0.99;
422  double bdSigMin = 3.0;
423  double bdMuPtMinLoose = 4.0;
424  double bdMuPtMinTight = 4.0;
425  double bdMuEtaMaxLoose = 2.2;
426  double bdMuEtaMaxTight = 2.2;
427 
428  bdJPsiBasicSelect = new BPHCompositeBasicSelect(bdJPsiMassMin, bdJPsiMassMax, bdJPsiPtMin, bdJPsiEtaMax, bdJPsiRMax);
429  bdKx0BasicSelect = new BPHCompositeBasicSelect(bdKx0MassMin, bdKx0MassMax, bdKx0PtMin, bdKx0EtaMax, bdKx0RMax);
430  bdVertexSelect = new BPHFittedVertexSelect(bdProbMin, bdCosMin, bdSigMin);
431  bdJPsiDaughterSelect = new BPHDaughterSelect(bdMuPtMinLoose, bdMuPtMinTight, bdMuEtaMaxLoose, bdMuEtaMaxTight, &sms);
432 
433  double bsJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
434  double bsJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
435  double bsJPsiPtMin = 8.0;
436  double bsJPsiEtaMax = -1.0;
437  double bsJPsiRMax = -1.0;
438  double bsPhiMassMin = BPHParticleMasses::phiMass - 0.010;
439  double bsPhiMassMax = BPHParticleMasses::phiMass + 0.010;
440  double bsPhiPtMin = -1.0;
441  double bsPhiEtaMax = -1.0;
442  double bsPhiRMax = -1.0;
443  double bsProbMin = 0.10;
444  double bsCosMin = 0.99;
445  double bsSigMin = 3.0;
446  double bsMuPtMinLoose = 4.0;
447  double bsMuPtMinTight = 4.0;
448  double bsMuEtaMaxLoose = 2.2;
449  double bsMuEtaMaxTight = 2.2;
450  bsJPsiBasicSelect = new BPHCompositeBasicSelect(bsJPsiMassMin, bsJPsiMassMax, bsJPsiPtMin, bsJPsiEtaMax, bsJPsiRMax);
451  bsPhiBasicSelect = new BPHCompositeBasicSelect(bsPhiMassMin, bsPhiMassMax, bsPhiPtMin, bsPhiEtaMax, bsPhiRMax);
452  bsVertexSelect = new BPHFittedVertexSelect(bsProbMin, bsCosMin, bsSigMin);
453  bsJPsiDaughterSelect = new BPHDaughterSelect(bsMuPtMinLoose, bsMuPtMinTight, bsMuEtaMaxLoose, bsMuEtaMaxTight, &sms);
454 }
455 
457  delete phiBasicSelect;
458  delete jPsiBasicSelect;
459  delete psi2BasicSelect;
460  delete upsBasicSelect;
461  delete oniaVertexSelect;
462  delete oniaDaughterSelect;
463 
464  delete buJPsiBasicSelect;
465  delete buVertexSelect;
466  delete buJPsiDaughterSelect;
467 
468  delete bdJPsiBasicSelect;
469  delete bdKx0BasicSelect;
470  delete bdVertexSelect;
471  delete bdJPsiDaughterSelect;
472 
473  delete bsJPsiBasicSelect;
474  delete bsPhiBasicSelect;
475  delete bsVertexSelect;
476  delete bsJPsiDaughterSelect;
477 }
478 
481  desc.add<string>("oniaCandsLabel", "");
482  desc.add<string>("sdCandsLabel", "");
483  desc.add<string>("ssCandsLabel", "");
484  desc.add<string>("buCandsLabel", "");
485  desc.add<string>("bdCandsLabel", "");
486  desc.add<string>("bsCandsLabel", "");
487  descriptions.add("bphHistoSpecificDecay", desc);
488  return;
489 }
490 
492  createHisto("massPhi", 35, 0.85, 1.20); // Phi mass
493  createHisto("massJPsi", 35, 2.95, 3.30); // JPsi mass
494  createHisto("massPsi2", 60, 3.40, 4.00); // Psi2 mass
495  createHisto("massUps123", 125, 8.50, 11.0); // Ups mass
496  createHisto("massBu", 50, 5.00, 6.00); // Bu mass
497  createHisto("massBd", 50, 5.00, 6.00); // Bd mass
498  createHisto("massBs", 50, 5.00, 6.00); // Bs mass
499  createHisto("mfitBu", 50, 5.00, 6.00); // Bu mass, with constraint
500  createHisto("mfitBd", 50, 5.00, 6.00); // Bd mass, with constraint
501  createHisto("mfitBs", 50, 5.00, 6.00); // Bs mass, with constraint
502  createHisto("massBuJPsi", 35, 2.95, 3.30); // JPsi mass in Bu decay
503  createHisto("massBdJPsi", 35, 2.95, 3.30); // JPsi mass in Bd decay
504  createHisto("massBsJPsi", 35, 2.95, 3.30); // JPsi mass in Bs decay
505  createHisto("massBsPhi", 50, 1.01, 1.03); // Phi mass in Bs decay
506  createHisto("massBdKx0", 50, 0.80, 1.05); // Kx0 mass in Bd decay
507 
508  createHisto("massFull", 200, 2.00, 12.0); // Full onia mass
509 
510  return;
511 }
512 
514  // get magnetic field
516  es.get<IdealMagneticFieldRecord>().get(magneticField);
517 
518  // get object collections
519  // collections are got through "BPHTokenWrapper" interface to allow
520  // uniform access in different CMSSW versions
521 
523 
525  int iqo;
526  int nqo = 0;
527  if (useOnia) {
528  oniaCandsToken.get(ev, oniaCands);
529  nqo = oniaCands->size();
530  }
531 
532  for (iqo = 0; iqo < nqo; ++iqo) {
533  LogTrace("DataDump") << "*********** quarkonium " << iqo << "/" << nqo << " ***********";
534  const pat::CompositeCandidate& cand = oniaCands->at(iqo);
535  if (!oniaVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(cand, "primaryVertex")))
536  continue;
537  if (!oniaDaughterSelect->accept(cand))
538  continue;
539  fillHisto("Full", cand);
540  if (phiBasicSelect->accept(cand))
541  fillHisto("Phi", cand);
542  if (jPsiBasicSelect->accept(cand))
543  fillHisto("JPsi", cand);
544  if (psi2BasicSelect->accept(cand))
545  fillHisto("Psi2", cand);
546  if (upsBasicSelect->accept(cand))
547  fillHisto("Ups123", cand);
548  }
549 
551 
553  int ibu;
554  int nbu = 0;
555  if (useBu) {
556  buCandsToken.get(ev, buCands);
557  nbu = buCands->size();
558  }
559 
560  for (ibu = 0; ibu < nbu; ++ibu) {
561  LogTrace("DataDump") << "*********** Bu " << ibu << "/" << nbu << " ***********";
562  const pat::CompositeCandidate& cand = buCands->at(ibu);
563  const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
564  LogTrace("DataDump") << "JPsi: " << jPsi;
565  if (jPsi == nullptr)
566  continue;
567  if (!buJPsiBasicSelect->accept(*jPsi))
568  continue;
569  if (!buJPsiDaughterSelect->accept(*jPsi))
570  continue;
571  if (!buVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")))
572  continue;
573  const reco::Candidate* kptr = BPHDaughters::get(cand, 0.49, 0.50).front();
574  if (kptr == nullptr)
575  continue;
576  if (kptr->pt() < buKPtMin)
577  continue;
578  fillHisto("Bu", cand);
579  fillHisto("BuJPsi", *jPsi);
580  }
581 
583 
585  int ibd;
586  int nbd = 0;
587  if (useBd) {
588  bdCandsToken.get(ev, bdCands);
589  nbd = bdCands->size();
590  }
591 
592  for (ibd = 0; ibd < nbd; ++ibd) {
593  LogTrace("DataDump") << "*********** Bd " << ibd << "/" << nbd << " ***********";
594  const pat::CompositeCandidate& cand = bdCands->at(ibd);
595  const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
596  LogTrace("DataDump") << "JPsi: " << jPsi;
597  if (jPsi == nullptr)
598  continue;
599  const pat::CompositeCandidate* kx0 = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToKx0");
600  LogTrace("DataDump") << "Kx0: " << kx0;
601  if (kx0 == nullptr)
602  continue;
603  if (!bdJPsiBasicSelect->accept(*jPsi))
604  continue;
605  if (!bdKx0BasicSelect->accept(*kx0))
606  continue;
607  if (!bdJPsiDaughterSelect->accept(*jPsi))
608  continue;
609  if (!bdVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")))
610  continue;
611  fillHisto("Bd", cand);
612  fillHisto("BdJPsi", *jPsi);
613  fillHisto("BdKx0", *kx0);
614  }
615 
617 
619  int ibs;
620  int nbs = 0;
621  if (useBs) {
622  bsCandsToken.get(ev, bsCands);
623  nbs = bsCands->size();
624  }
625 
626  for (ibs = 0; ibs < nbs; ++ibs) {
627  LogTrace("DataDump") << "*********** Bs " << ibs << "/" << nbs << " ***********";
628  const pat::CompositeCandidate& cand = bsCands->at(ibs);
629  const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
630  LogTrace("DataDump") << "JPsi: " << jPsi;
631  if (jPsi == nullptr)
632  continue;
633  const pat::CompositeCandidate* phi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToPhi");
634  LogTrace("DataDump") << "Phi: " << phi;
635  if (phi == nullptr)
636  continue;
637  if (!bsJPsiBasicSelect->accept(*jPsi))
638  continue;
639  if (!bsPhiBasicSelect->accept(*phi))
640  continue;
641  if (!bsJPsiDaughterSelect->accept(*jPsi))
642  continue;
643  if (!bsVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")))
644  continue;
645  fillHisto("Bs", cand);
646  fillHisto("BsJPsi", *jPsi);
647  fillHisto("BsPhi", *phi);
648  }
649 
650  return;
651 }
652 
654 
656  float mass = (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1);
657  fillHisto("mass" + name, cand.mass());
658  fillHisto("mfit" + name, mass);
659  return;
660 }
661 
662 void BPHHistoSpecificDecay::fillHisto(const string& name, float x) {
663  map<string, TH1F*>::iterator iter = histoMap.find(name);
664  map<string, TH1F*>::iterator iend = histoMap.end();
665  if (iter == iend)
666  return;
667  iter->second->Fill(x);
668  return;
669 }
670 
671 void BPHHistoSpecificDecay::createHisto(const string& name, int nbin, float hmin, float hmax) {
672  histoMap[name] = fs->make<TH1F>(name.c_str(), name.c_str(), nbin, hmin, hmax);
673  return;
674 }
675 
677 
Analysis-level particle class.
BPHDaughterSelect(float ptMinLoose, float ptMinTight, float etaMaxLoose, float etaMaxTight, const BPHSoftMuonSelect *softMuonselector=0)
const reco::Vertex * pv
double eta() const final
momentum pseudorapidity
BPHHistoSpecificDecay(const edm::ParameterSet &ps)
BPHCompositeBasicSelect(float massMin, float massMax, float ptMin=-1.0, float etaMax=-1.0, float rapidityMax=-1.0)
BPHFittedVertexSelect(float probMin, float cosMin=-1.0, float sigMin=-1.0)
static const T * getByRef(const pat::CompositeCandidate &cand, const string &name)
bool accept(const pat::CompositeCandidate &cand, const reco::Vertex *pv=0) const override
const AlgebraicSymMatrix33 matrix() const
double y() const
y coordinate
Definition: Vertex.h:117
virtual const Candidate * daughter(size_type i) const =0
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
double rapidity() const final
rapidity
bool hasUserFloat(const std::string &key) const
Return true if there is a user-defined float with a given name.
Definition: PATObject.h:378
double pt() const final
transverse momentum
bool ev
bool accept(const pat::CompositeCandidate &cand, const reco::Vertex *pvtx=0) const override
static const double jPsiMass
const Point & position() const
position
Definition: Vertex.h:113
float userFloat(const std::string &key) const
Definition: PATObject.h:876
bool accept(const pat::CompositeCandidate &cand, const reco::Vertex *pv=0) const override
void analyze(const edm::Event &ev, const edm::EventSetup &es) override
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
static const double kx0Mass
reco::TrackRef innerTrack() const override
reference to Track reconstructed in the tracker only (reimplemented from reco::Muon) ...
Definition: Muon.h:72
T sqrt(T t)
Definition: SSEVec.h:19
void createHisto(const std::string &name, int nbin, float hmin, float hmax)
def pv(vc)
Definition: MetAnalyzer.py:7
double chi2() const
chi-squares
Definition: Vertex.h:102
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
float ChiSquaredProbability(double chiSquared, double nrDOF)
BPHFittedBasicSelect(float massMin, float massMax, float ptMin=-1.0, float etaMax=-1.0, float rapidityMax=-1.0)
const Candidate * daughter(size_type) const override
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
const T * userData(const std::string &key) const
Returns user-defined data. Returns NULL if the data is not present, or not of type T...
Definition: PATObject.h:327
ParameterDescriptionBase * add(U const &iLabel, T const &value)
bool hasUserData(const std::string &key) const
Check if user data with a specific type is present.
Definition: PATObject.h:332
#define LogTrace(id)
double ndof() const
Definition: Vertex.h:109
bool isGoodMuon(const reco::Muon &muon, SelectionType type, reco::Muon::ArbitrationType arbitrationType=reco::Muon::SegmentAndTrackArbitration)
main GoodMuon wrapper call
bool accept(const reco::Candidate &cand, const reco::Vertex *pv) const
static const double phiMass
bool accept(const pat::CompositeCandidate &cand, const reco::Vertex *pvtx) const override
double x() const
x coordinate
Definition: Vertex.h:115
bool accept(const pat::CompositeCandidate &cand, const reco::Vertex *pv=0) const override
BPHCompositeVertexSelect(float probMin, float cosMin=-1.0, float sigMin=-1.0)
virtual double eta() const =0
momentum pseudorapidity
virtual double pt() const =0
transverse momentum
virtual double mass() const =0
mass
void add(std::string const &label, ParameterSetDescription const &psetDescription)
double value() const
Definition: Measurement1D.h:25
Error error() const
return SMatrix
Definition: Vertex.h:149
T get() const
Definition: EventSetup.h:73
ROOT::Math::SVector< double, 3 > AlgebraicVector3
BPHSoftMuonSelect(int cutTrackerLayers=5, int cutPixelLayers=0, float maxDxy=0.3, float maxDz=20.0, bool goodMuon=true, bool highPurity=true)
static vector< const reco::Candidate * > get(const pat::CompositeCandidate &cand, float massMin, float massMax)
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
long double T
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
const BPHSoftMuonSelect * sms
Analysis-level muon class.
Definition: Muon.h:51
#define SET_LABEL(NAME, PSET)
double mass() const final
mass
void fillHisto(const std::string &name, const pat::CompositeCandidate &cand)