CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
BPHWriteSpecificDecay.cc
Go to the documentation of this file.
2 
4 
12 
20 
27 
30 
33 
36 
39 
40 #include <set>
41 #include <string>
42 #include <iostream>
43 #include <fstream>
44 
45 using namespace std;
46 
47 #define SET_PAR(TYPE, NAME, PSET) (NAME = PSET.getParameter<TYPE>(#NAME))
48 // SET_PAR(string,xyz,ps);
49 // is equivalent to
50 // ( xyz = ps.getParameter< string >( "xyx" ) )
51 
53  usePV = (!SET_PAR(string, pVertexLabel, ps).empty());
54  usePM = (!SET_PAR(string, patMuonLabel, ps).empty());
55  useCC = (!SET_PAR(string, ccCandsLabel, ps).empty());
56  usePF = (!SET_PAR(string, pfCandsLabel, ps).empty());
57  usePC = (!SET_PAR(string, pcCandsLabel, ps).empty());
58  useGP = (!SET_PAR(string, gpCandsLabel, ps).empty());
59  SET_PAR(string, oniaName, ps);
60  SET_PAR(string, sdName, ps);
61  SET_PAR(string, ssName, ps);
62  SET_PAR(string, buName, ps);
63  SET_PAR(string, bdName, ps);
64  SET_PAR(string, bsName, ps);
65 
66  SET_PAR(bool, writeMomentum, ps);
67  SET_PAR(bool, writeVertex, ps);
68 
69  rMap["Onia"] = Onia;
70  rMap["PHiMuMu"] = Pmm;
71  rMap["Psi1"] = Psi1;
72  rMap["Psi2"] = Psi2;
73  rMap["Ups"] = Ups;
74  rMap["Ups1"] = Ups1;
75  rMap["Ups2"] = Ups2;
76  rMap["Ups3"] = Ups3;
77  rMap["Kx0"] = Kx0;
78  rMap["PhiKK"] = Pkk;
79  rMap["Bu"] = Bu;
80  rMap["Bd"] = Bd;
81  rMap["Bs"] = Bs;
82 
83  pMap["ptMin"] = ptMin;
84  pMap["etaMax"] = etaMax;
85  pMap["mJPsiMin"] = mPsiMin;
86  pMap["mJPsiMax"] = mPsiMax;
87  pMap["mKx0Min"] = mKx0Min;
88  pMap["mKx0Max"] = mKx0Max;
89  pMap["mPhiMin"] = mPhiMin;
90  pMap["mPhiMax"] = mPhiMax;
91  pMap["massMin"] = massMin;
92  pMap["massMax"] = massMax;
93  pMap["probMin"] = probMin;
94  pMap["massFitMin"] = mFitMin;
95  pMap["massFitMax"] = mFitMax;
96  pMap["constrMass"] = constrMass;
97  pMap["constrSigma"] = constrSigma;
98 
99  fMap["constrMJPsi"] = constrMJPsi;
100  fMap["writeCandidate"] = writeCandidate;
101 
102  recoOnia = recoKx0 = writeKx0 = recoPkk = writePkk = recoBu = writeBu = recoBd = writeBd = recoBs = writeBs = false;
103 
104  writeOnia = true;
105  const vector<edm::ParameterSet> recoSelect = ps.getParameter<vector<edm::ParameterSet> >("recoSelect");
106  int iSel;
107  int nSel = recoSelect.size();
108  for (iSel = 0; iSel < nSel; ++iSel)
109  setRecoParameters(recoSelect[iSel]);
110  if (!recoOnia)
111  writeOnia = false;
112 
113  if (recoBu)
114  recoOnia = true;
115  if (recoBd)
116  recoOnia = recoKx0 = true;
117  if (recoBs)
118  recoOnia = recoPkk = true;
119  if (writeBu)
120  writeOnia = true;
121  if (writeBd)
122  writeOnia = writeKx0 = true;
123  if (writeBs)
124  writeOnia = writePkk = true;
125 
126  if (usePV)
127  consume<vector<reco::Vertex> >(pVertexToken, pVertexLabel);
128  if (usePM)
129  consume<pat::MuonCollection>(patMuonToken, patMuonLabel);
130  if (useCC)
131  consume<vector<pat::CompositeCandidate> >(ccCandsToken, ccCandsLabel);
132  if (usePF)
133  consume<vector<reco::PFCandidate> >(pfCandsToken, pfCandsLabel);
134  if (usePC)
135  consume<vector<BPHTrackReference::candidate> >(pcCandsToken, pcCandsLabel);
136  if (useGP)
137  consume<vector<pat::GenericParticle> >(gpCandsToken, gpCandsLabel);
138 
139  if (writeOnia)
140  produces<pat::CompositeCandidateCollection>(oniaName);
141  if (writeKx0)
142  produces<pat::CompositeCandidateCollection>(sdName);
143  if (writePkk)
144  produces<pat::CompositeCandidateCollection>(ssName);
145  if (writeBu)
146  produces<pat::CompositeCandidateCollection>(buName);
147  if (writeBd)
148  produces<pat::CompositeCandidateCollection>(bdName);
149  if (writeBs)
150  produces<pat::CompositeCandidateCollection>(bsName);
151 }
152 
154 
157  desc.add<string>("pVertexLabel", "");
158  desc.add<string>("patMuonLabel", "");
159  desc.add<string>("ccCandsLabel", "");
160  desc.add<string>("pfCandsLabel", "");
161  desc.add<string>("pcCandsLabel", "");
162  desc.add<string>("gpCandsLabel", "");
163  desc.add<string>("oniaName", "oniaCand");
164  desc.add<string>("sdName", "kx0Cand");
165  desc.add<string>("ssName", "phiCand");
166  desc.add<string>("buName", "buFitted");
167  desc.add<string>("bdName", "bdFitted");
168  desc.add<string>("bsName", "bsFitted");
169  desc.add<bool>("writeVertex", true);
170  desc.add<bool>("writeMomentum", true);
172  dpar.add<string>("name");
173  dpar.add<double>("ptMin", -2.0e35);
174  dpar.add<double>("etaMax", -2.0e35);
175  dpar.add<double>("mJPsiMin", -2.0e35);
176  dpar.add<double>("mJPsiMax", -2.0e35);
177  dpar.add<double>("mKx0Min", -2.0e35);
178  dpar.add<double>("mKx0Max", -2.0e35);
179  dpar.add<double>("mPhiMin", -2.0e35);
180  dpar.add<double>("mPhiMax", -2.0e35);
181  dpar.add<double>("massMin", -2.0e35);
182  dpar.add<double>("massMax", -2.0e35);
183  dpar.add<double>("probMin", -2.0e35);
184  dpar.add<double>("massFitMin", -2.0e35);
185  dpar.add<double>("massFitMax", -2.0e35);
186  dpar.add<double>("constrMass", -2.0e35);
187  dpar.add<double>("constrSigma", -2.0e35);
188  dpar.add<bool>("constrMJPsi", true);
189  dpar.add<bool>("writeCandidate", true);
190  vector<edm::ParameterSet> rpar;
191  desc.addVPSet("recoSelect", dpar, rpar);
192  descriptions.add("bphWriteSpecificDecay", desc);
193  return;
194 }
195 
197 
199  fill(ev, es);
200  if (writeOnia)
201  write(ev, lFull, oniaName);
202  if (writeKx0)
203  write(ev, lSd, sdName);
204  if (writePkk)
205  write(ev, lSs, ssName);
206  if (writeBu)
207  write(ev, lBu, buName);
208  if (writeBd)
209  write(ev, lBd, bdName);
210  if (writeBs)
211  write(ev, lBs, bsName);
212  return;
213 }
214 
216  lFull.clear();
217  lJPsi.clear();
218  lSd.clear();
219  lSs.clear();
220  lBu.clear();
221  lBd.clear();
222  lBs.clear();
223  jPsiOMap.clear();
224  pvRefMap.clear();
225  ccRefMap.clear();
226 
227  // get magnetic field
229  es.get<IdealMagneticFieldRecord>().get(magneticField);
230 
231  // get object collections
232  // collections are got through "BPHTokenWrapper" interface to allow
233  // uniform access in different CMSSW versions
234 
236  pVertexToken.get(ev, pVertices);
237  int npv = pVertices->size();
238 
239  int nrc = 0;
240 
241  // get reco::PFCandidate collection (in full AOD )
243  if (usePF) {
244  pfCandsToken.get(ev, pfCands);
245  nrc = pfCands->size();
246  }
247 
248  // get pat::PackedCandidate collection (in MiniAOD)
249  // pat::PackedCandidate is not defined in CMSSW_5XY, so a
250  // typedef (BPHTrackReference::candidate) is used, actually referring
251  // to pat::PackedCandidate only for CMSSW versions where it's defined
253  if (usePC) {
254  pcCandsToken.get(ev, pcCands);
255  nrc = pcCands->size();
256  }
257 
258  // get pat::GenericParticle collection (in skimmed data)
260  if (useGP) {
261  gpCandsToken.get(ev, gpCands);
262  nrc = gpCands->size();
263  }
264 
265  // get pat::Muon collection (in full AOD and MiniAOD)
267  if (usePM) {
268  patMuonToken.get(ev, patMuon);
269  }
270 
271  // get muons from pat::CompositeCandidate objects describing onia;
272  // muons from all composite objects are copied to an unique std::vector
273  vector<const reco::Candidate*> muDaugs;
274  set<const pat::Muon*> muonSet;
275  typedef multimap<const reco::Candidate*, const pat::CompositeCandidate*> mu_cc_map;
276  mu_cc_map muCCMap;
277  if (useCC) {
279  ccCandsToken.get(ev, ccCands);
280  int n = ccCands->size();
281  muDaugs.clear();
282  muDaugs.reserve(n);
283  muonSet.clear();
284  set<const pat::Muon*>::const_iterator iter;
285  set<const pat::Muon*>::const_iterator iend;
286  int i;
287  for (i = 0; i < n; ++i) {
288  const pat::CompositeCandidate& cc = ccCands->at(i);
289  int j;
290  int m = cc.numberOfDaughters();
291  for (j = 0; j < m; ++j) {
292  const reco::Candidate* dp = cc.daughter(j);
293  const pat::Muon* mp = dynamic_cast<const pat::Muon*>(dp);
294  iter = muonSet.begin();
295  iend = muonSet.end();
296  bool add = (mp != nullptr) && (muonSet.find(mp) == iend);
297  while (add && (iter != iend)) {
298  if (BPHRecoBuilder::sameTrack(mp, *iter++, 1.0e-5))
299  add = false;
300  }
301  if (add)
302  muonSet.insert(mp);
303  // associate muon to the CompositeCandidate containing it
304  muCCMap.insert(pair<const reco::Candidate*, const pat::CompositeCandidate*>(dp, &cc));
305  }
306  }
307  iter = muonSet.begin();
308  iend = muonSet.end();
309  while (iter != iend)
310  muDaugs.push_back(*iter++);
311  }
312 
313  map<recoType, map<parType, double> >::const_iterator rIter = parMap.begin();
314  map<recoType, map<parType, double> >::const_iterator rIend = parMap.end();
315 
316  // reconstruct quarkonia
317 
318  BPHOniaToMuMuBuilder* onia = nullptr;
319  if (recoOnia) {
320  if (usePM)
321  onia = new BPHOniaToMuMuBuilder(
322  es, BPHRecoBuilder::createCollection(patMuon, "cfmig"), BPHRecoBuilder::createCollection(patMuon, "cfmig"));
323  else if (useCC)
324  onia = new BPHOniaToMuMuBuilder(
325  es, BPHRecoBuilder::createCollection(muDaugs, "cfmig"), BPHRecoBuilder::createCollection(muDaugs, "cfmig"));
326  }
327 
328  if (onia != nullptr) {
329  while (rIter != rIend) {
330  const map<recoType, map<parType, double> >::value_type& rEntry = *rIter++;
331  recoType rType = rEntry.first;
332  const map<parType, double>& pMap = rEntry.second;
334  switch (rType) {
335  case Pmm:
337  break;
338  case Psi1:
340  break;
341  case Psi2:
343  break;
344  case Ups:
346  break;
347  case Ups1:
349  break;
350  case Ups2:
352  break;
353  case Ups3:
355  break;
356  default:
357  continue;
358  }
359  map<parType, double>::const_iterator pIter = pMap.begin();
360  map<parType, double>::const_iterator pIend = pMap.end();
361  while (pIter != pIend) {
362  const map<parType, double>::value_type& pEntry = *pIter++;
363  parType id = pEntry.first;
364  double pv = pEntry.second;
365  switch (id) {
366  case ptMin:
367  onia->setPtMin(type, pv);
368  break;
369  case etaMax:
370  onia->setEtaMax(type, pv);
371  break;
372  case massMin:
373  onia->setMassMin(type, pv);
374  break;
375  case massMax:
376  onia->setMassMax(type, pv);
377  break;
378  case probMin:
379  onia->setProbMin(type, pv);
380  break;
381  case constrMass:
382  onia->setConstr(type, pv, onia->getConstrSigma(type));
383  break;
384  case constrSigma:
385  onia->setConstr(type, onia->getConstrMass(type), pv);
386  break;
387  default:
388  break;
389  }
390  }
391  }
392  lFull = onia->build();
393  }
394 
395  // associate onia to primary vertex
396 
397  int iFull;
398  int nFull = lFull.size();
399  map<const BPHRecoCandidate*, const reco::Vertex*> oniaVtxMap;
400 
401  typedef mu_cc_map::const_iterator mu_cc_iter;
402  for (iFull = 0; iFull < nFull; ++iFull) {
403  const reco::Vertex* pVtx = nullptr;
404  int pvId = 0;
405  const BPHPlusMinusCandidate* ptr = lFull[iFull].get();
406  const std::vector<const reco::Candidate*>& daugs = ptr->daughters();
407 
408  // try to recover primary vertex association in skim data:
409  // get the CompositeCandidate containing both muons
410  pair<mu_cc_iter, mu_cc_iter> cc0 = muCCMap.equal_range(ptr->originalReco(daugs[0]));
411  pair<mu_cc_iter, mu_cc_iter> cc1 = muCCMap.equal_range(ptr->originalReco(daugs[1]));
412  mu_cc_iter iter0 = cc0.first;
413  mu_cc_iter iend0 = cc0.second;
414  mu_cc_iter iter1 = cc1.first;
415  mu_cc_iter iend1 = cc1.second;
416  while ((iter0 != iend0) && (pVtx == nullptr)) {
417  const pat::CompositeCandidate* ccp = iter0++->second;
418  while (iter1 != iend1) {
419  if (ccp != iter1++->second)
420  continue;
421  pVtx = ccp->userData<reco::Vertex>("PVwithmuons");
422  const reco::Vertex* sVtx = nullptr;
423  const reco::Vertex::Point& pPos = pVtx->position();
424  float dMin = 999999.;
425  int ipv;
426  for (ipv = 0; ipv < npv; ++ipv) {
427  const reco::Vertex* tVtx = &pVertices->at(ipv);
428  const reco::Vertex::Point& tPos = tVtx->position();
429  float dist = pow(pPos.x() - tPos.x(), 2) + pow(pPos.y() - tPos.y(), 2) + pow(pPos.z() - tPos.z(), 2);
430  if (dist < dMin) {
431  dMin = dist;
432  sVtx = tVtx;
433  pvId = ipv;
434  }
435  }
436  pVtx = sVtx;
437  break;
438  }
439  }
440 
441  // if not found, as ofr other type of inut data,
442  // try to get the nearest primary vertex in z direction
443  if (pVtx == nullptr) {
444  const reco::Vertex::Point& sVtp = ptr->vertex().position();
445  GlobalPoint cPos(sVtp.x(), sVtp.y(), sVtp.z());
446  const pat::CompositeCandidate& sCC = ptr->composite();
447  GlobalVector cDir(sCC.px(), sCC.py(), sCC.pz());
448  GlobalPoint bPos(0.0, 0.0, 0.0);
449  GlobalVector bDir(0.0, 0.0, 1.0);
451  bool state = ttmd.calculate(GlobalTrajectoryParameters(cPos, cDir, TrackCharge(0), &(*magneticField)),
452  GlobalTrajectoryParameters(bPos, bDir, TrackCharge(0), &(*magneticField)));
453  float minDz = 999999.;
454  float extrapZ = (state ? ttmd.points().first.z() : -9e20);
455  int ipv;
456  for (ipv = 0; ipv < npv; ++ipv) {
457  const reco::Vertex& tVtx = pVertices->at(ipv);
458  float deltaZ = fabs(extrapZ - tVtx.position().z());
459  if (deltaZ < minDz) {
460  minDz = deltaZ;
461  pVtx = &tVtx;
462  pvId = ipv;
463  }
464  }
465  }
466 
467  oniaVtxMap[ptr] = pVtx;
468  pvRefMap[ptr] = vertex_ref(pVertices, pvId);
469  }
470  pVertexToken.get(ev, pVertices);
471 
472  // get JPsi subsample and associate JPsi candidate to original
473  // generic onia candidate
474  if (nFull)
475  lJPsi = onia->getList(BPHOniaToMuMuBuilder::Psi1);
476 
477  int nJPsi = lJPsi.size();
478  delete onia;
479 
480  if (!nJPsi)
481  return;
482  if (!nrc)
483  return;
484 
485  int ij;
486  int io;
487  int nj = lJPsi.size();
488  int no = lFull.size();
489  for (ij = 0; ij < nj; ++ij) {
490  const BPHRecoCandidate* jp = lJPsi[ij].get();
491  for (io = 0; io < no; ++io) {
492  const BPHRecoCandidate* oc = lFull[io].get();
493  if ((jp->originalReco(jp->getDaug("MuPos")) == oc->originalReco(oc->getDaug("MuPos"))) &&
494  (jp->originalReco(jp->getDaug("MuNeg")) == oc->originalReco(oc->getDaug("MuNeg")))) {
495  jPsiOMap[jp] = oc;
496  break;
497  }
498  }
499  }
500 
501  // build and dump Bu
502 
503  BPHBuToJPsiKBuilder* bu = nullptr;
504  if (recoBu) {
505  if (usePF)
506  bu = new BPHBuToJPsiKBuilder(es, lJPsi, BPHRecoBuilder::createCollection(pfCands));
507  else if (usePC)
508  bu = new BPHBuToJPsiKBuilder(es, lJPsi, BPHRecoBuilder::createCollection(pcCands));
509  else if (useGP)
510  bu = new BPHBuToJPsiKBuilder(es, lJPsi, BPHRecoBuilder::createCollection(gpCands));
511  }
512 
513  if (bu != nullptr) {
514  rIter = parMap.find(Bu);
515  if (rIter != rIend) {
516  const map<parType, double>& pMap = rIter->second;
517  map<parType, double>::const_iterator pIter = pMap.begin();
518  map<parType, double>::const_iterator pIend = pMap.end();
519  while (pIter != pIend) {
520  const map<parType, double>::value_type& pEntry = *pIter++;
521  parType id = pEntry.first;
522  double pv = pEntry.second;
523  switch (id) {
524  case ptMin:
525  bu->setKPtMin(pv);
526  break;
527  case etaMax:
528  bu->setKEtaMax(pv);
529  break;
530  case mPsiMin:
531  bu->setJPsiMassMin(pv);
532  break;
533  case mPsiMax:
534  bu->setJPsiMassMax(pv);
535  break;
536  case massMin:
537  bu->setMassMin(pv);
538  break;
539  case massMax:
540  bu->setMassMax(pv);
541  break;
542  case probMin:
543  bu->setProbMin(pv);
544  break;
545  case mFitMin:
546  bu->setMassFitMin(pv);
547  break;
548  case mFitMax:
549  bu->setMassFitMax(pv);
550  break;
551  case constrMJPsi:
552  bu->setConstr(pv > 0);
553  break;
554  case writeCandidate:
555  writeBu = (pv > 0);
556  break;
557  default:
558  break;
559  }
560  }
561  }
562  lBu = bu->build();
563  delete bu;
564  }
565 
566  // build and dump Kx0
567 
568  vector<BPHPlusMinusConstCandPtr> lKx0;
569  BPHKx0ToKPiBuilder* kx0 = nullptr;
570  if (recoKx0) {
571  if (usePF)
572  kx0 = new BPHKx0ToKPiBuilder(
574  else if (usePC)
575  kx0 = new BPHKx0ToKPiBuilder(
577  else if (useGP)
578  kx0 = new BPHKx0ToKPiBuilder(
580  }
581 
582  if (kx0 != nullptr) {
583  rIter = parMap.find(Kx0);
584  if (rIter != rIend) {
585  const map<parType, double>& pMap = rIter->second;
586  map<parType, double>::const_iterator pIter = pMap.begin();
587  map<parType, double>::const_iterator pIend = pMap.end();
588  while (pIter != pIend) {
589  const map<parType, double>::value_type& pEntry = *pIter++;
590  parType id = pEntry.first;
591  double pv = pEntry.second;
592  switch (id) {
593  case ptMin:
594  kx0->setPtMin(pv);
595  break;
596  case etaMax:
597  kx0->setEtaMax(pv);
598  break;
599  case massMin:
600  kx0->setMassMin(pv);
601  break;
602  case massMax:
603  kx0->setMassMax(pv);
604  break;
605  case probMin:
606  kx0->setProbMin(pv);
607  break;
608  case constrMass:
609  kx0->setConstr(pv, kx0->getConstrSigma());
610  break;
611  case constrSigma:
612  kx0->setConstr(kx0->getConstrMass(), pv);
613  break;
614  case writeCandidate:
615  writeKx0 = (pv > 0);
616  break;
617  default:
618  break;
619  }
620  }
621  }
622  lKx0 = kx0->build();
623  delete kx0;
624  }
625 
626  int nKx0 = lKx0.size();
627 
628  // build and dump Bd
629 
630  if (recoBd && nKx0) {
631  BPHBdToJPsiKxBuilder* bd = new BPHBdToJPsiKxBuilder(es, lJPsi, lKx0);
632  rIter = parMap.find(Bd);
633  if (rIter != rIend) {
634  const map<parType, double>& pMap = rIter->second;
635  map<parType, double>::const_iterator pIter = pMap.begin();
636  map<parType, double>::const_iterator pIend = pMap.end();
637  while (pIter != pIend) {
638  const map<parType, double>::value_type& pEntry = *pIter++;
639  parType id = pEntry.first;
640  double pv = pEntry.second;
641  switch (id) {
642  case mPsiMin:
643  bd->setJPsiMassMin(pv);
644  break;
645  case mPsiMax:
646  bd->setJPsiMassMax(pv);
647  break;
648  case mKx0Min:
649  bd->setKxMassMin(pv);
650  break;
651  case mKx0Max:
652  bd->setKxMassMax(pv);
653  break;
654  case massMin:
655  bd->setMassMin(pv);
656  break;
657  case massMax:
658  bd->setMassMax(pv);
659  break;
660  case probMin:
661  bd->setProbMin(pv);
662  break;
663  case mFitMin:
664  bd->setMassFitMin(pv);
665  break;
666  case mFitMax:
667  bd->setMassFitMax(pv);
668  break;
669  case constrMJPsi:
670  bd->setConstr(pv > 0);
671  break;
672  case writeCandidate:
673  writeBd = (pv > 0);
674  break;
675  default:
676  break;
677  }
678  }
679  }
680 
681  lBd = bd->build();
682  delete bd;
683 
684  set<BPHRecoConstCandPtr> sKx0;
685  int iBd;
686  int nBd = lBd.size();
687  for (iBd = 0; iBd < nBd; ++iBd)
688  sKx0.insert(lBd[iBd]->getComp("Kx0"));
689  set<BPHRecoConstCandPtr>::const_iterator iter = sKx0.begin();
690  set<BPHRecoConstCandPtr>::const_iterator iend = sKx0.end();
691  while (iter != iend)
692  lSd.push_back(*iter++);
693  }
694 
695  // build and dump Phi
696 
697  vector<BPHPlusMinusConstCandPtr> lPhi;
698  BPHPhiToKKBuilder* phi = nullptr;
699  if (recoPkk) {
700  if (usePF)
701  phi = new BPHPhiToKKBuilder(
703  else if (usePC)
704  phi = new BPHPhiToKKBuilder(
706  else if (useGP)
707  phi = new BPHPhiToKKBuilder(
709  }
710 
711  if (phi != nullptr) {
712  rIter = parMap.find(Pkk);
713  if (rIter != rIend) {
714  const map<parType, double>& pMap = rIter->second;
715  map<parType, double>::const_iterator pIter = pMap.begin();
716  map<parType, double>::const_iterator pIend = pMap.end();
717  while (pIter != pIend) {
718  const map<parType, double>::value_type& pEntry = *pIter++;
719  parType id = pEntry.first;
720  double pv = pEntry.second;
721  switch (id) {
722  case ptMin:
723  phi->setPtMin(pv);
724  break;
725  case etaMax:
726  phi->setEtaMax(pv);
727  break;
728  case massMin:
729  phi->setMassMin(pv);
730  break;
731  case massMax:
732  phi->setMassMax(pv);
733  break;
734  case probMin:
735  phi->setProbMin(pv);
736  break;
737  case constrMass:
738  phi->setConstr(pv, phi->getConstrSigma());
739  break;
740  case constrSigma:
741  phi->setConstr(phi->getConstrMass(), pv);
742  break;
743  case writeCandidate:
744  writePkk = (pv > 0);
745  break;
746  default:
747  break;
748  }
749  }
750  }
751  lPhi = phi->build();
752  delete phi;
753  }
754 
755  int nPhi = lPhi.size();
756 
757  // build and dump Bs
758 
759  if (recoBs && nPhi) {
760  BPHBsToJPsiPhiBuilder* bs = new BPHBsToJPsiPhiBuilder(es, lJPsi, lPhi);
761  rIter = parMap.find(Bs);
762  if (rIter != rIend) {
763  const map<parType, double>& pMap = rIter->second;
764  map<parType, double>::const_iterator pIter = pMap.begin();
765  map<parType, double>::const_iterator pIend = pMap.end();
766  while (pIter != pIend) {
767  const map<parType, double>::value_type& pEntry = *pIter++;
768  parType id = pEntry.first;
769  double pv = pEntry.second;
770  switch (id) {
771  case mPsiMin:
772  bs->setJPsiMassMin(pv);
773  break;
774  case mPsiMax:
775  bs->setJPsiMassMax(pv);
776  break;
777  case mPhiMin:
778  bs->setPhiMassMin(pv);
779  break;
780  case mPhiMax:
781  bs->setPhiMassMax(pv);
782  break;
783  case massMin:
784  bs->setMassMin(pv);
785  break;
786  case massMax:
787  bs->setMassMax(pv);
788  break;
789  case probMin:
790  bs->setProbMin(pv);
791  break;
792  case mFitMin:
793  bs->setMassFitMin(pv);
794  break;
795  case mFitMax:
796  bs->setMassFitMax(pv);
797  break;
798  case constrMJPsi:
799  bs->setConstr(pv > 0);
800  break;
801  case writeCandidate:
802  writeBs = (pv > 0);
803  break;
804  default:
805  break;
806  }
807  }
808  }
809 
810  lBs = bs->build();
811  delete bs;
812 
813  set<BPHRecoConstCandPtr> sPhi;
814  int iBs;
815  int nBs = lBs.size();
816  for (iBs = 0; iBs < nBs; ++iBs)
817  sPhi.insert(lBs[iBs]->getComp("Phi"));
818  set<BPHRecoConstCandPtr>::const_iterator iter = sPhi.begin();
819  set<BPHRecoConstCandPtr>::const_iterator iend = sPhi.end();
820  while (iter != iend)
821  lSs.push_back(*iter++);
822  }
823 
824  return;
825 }
826 
828 
830  const string& name = ps.getParameter<string>("name");
831  bool writeCandidate = ps.getParameter<bool>("writeCandidate");
832  switch (rMap[name]) {
833  case Onia:
834  recoOnia = true;
835  writeOnia = writeCandidate;
836  break;
837  case Pmm:
838  case Psi1:
839  case Psi2:
840  case Ups:
841  case Ups1:
842  case Ups2:
843  case Ups3:
844  recoOnia = true;
845  break;
846  case Kx0:
847  recoKx0 = true;
848  writeKx0 = writeCandidate;
849  break;
850  case Pkk:
851  recoPkk = true;
852  writePkk = writeCandidate;
853  break;
854  case Bu:
855  recoBu = true;
856  writeBu = writeCandidate;
857  break;
858  case Bd:
859  recoBd = true;
860  writeBd = writeCandidate;
861  break;
862  case Bs:
863  recoBs = true;
864  writeBs = writeCandidate;
865  break;
866  }
867 
868  map<string, parType>::const_iterator pIter = pMap.begin();
869  map<string, parType>::const_iterator pIend = pMap.end();
870  while (pIter != pIend) {
871  const map<string, parType>::value_type& entry = *pIter++;
872  const string& pn = entry.first;
873  parType id = entry.second;
874  double pv = ps.getParameter<double>(pn);
875  if (pv > -1.0e35)
876  edm::LogVerbatim("Configuration") << "BPHWriteSpecificDecay::setRecoParameters: set " << pn << " for " << name
877  << " : " << (parMap[rMap[name]][id] = pv);
878  }
879 
880  map<string, parType>::const_iterator fIter = fMap.begin();
881  map<string, parType>::const_iterator fIend = fMap.end();
882  while (fIter != fIend) {
883  const map<string, parType>::value_type& entry = *fIter++;
884  const string& fn = entry.first;
885  parType id = entry.second;
886  edm::LogVerbatim("Configuration") << "BPHWriteSpecificDecay::setRecoParameters: set " << fn << " for " << name
887  << " : " << (parMap[rMap[name]][id] = (ps.getParameter<bool>(fn) ? 1 : -1));
888  }
889 }
890 
892 
Analysis-level particle class.
type
Definition: HCALResponse.h:21
void setEtaMax(double eta)
double getConstrSigma() const
static bool sameTrack(const reco::Candidate *lCand, const reco::Candidate *rCand, double minPDifference)
T getParameter(std::string const &) const
std::vector< BPHRecoConstCandPtr > build()
build Bu candidates
ParameterDescriptionBase * addVPSet(U const &iLabel, ParameterSetDescription const &validator, std::vector< ParameterSet > const &defaults)
std::pair< GlobalPoint, GlobalPoint > points() const override
void setJPsiMassMin(double m)
set cuts
void setMassMax(oniaType type, double m)
double getConstrMass(oniaType type) const
void setPtMin(double pt)
set cuts
void setJPsiMassMin(double m)
set cuts
void setKPtMin(double pt)
set cuts
void setMassFitMax(double m)
void setKEtaMax(double eta)
void setMassMin(double m)
bool ev
void setProbMin(double p)
void setPtMin(double pt)
set cuts
double getConstrMass() const
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
const Point & position() const
position
Definition: Vertex.h:113
#define SET_PAR(TYPE, NAME, PSET)
virtual const std::vector< const reco::Candidate * > & daughters() const
virtual void fill(edm::Event &ev, const edm::EventSetup &es)
void setProbMin(oniaType type, double p)
bool calculate(const TrajectoryStateOnSurface &sta, const TrajectoryStateOnSurface &stb) override
static BPHGenericCollection * createCollection(const edm::Handle< T > &collection, const std::string &list="cfhpmig")
std::vector< BPHPlusMinusConstCandPtr > build()
build Phi candidates
U second(std::pair< T, U > const &p)
const pat::CompositeCandidate & composite() const override
get a composite by the simple sum of simple particles
int TrackCharge
Definition: TrackCharge.h:4
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
BPHWriteSpecificDecay(const edm::ParameterSet &ps)
double getConstrSigma(oniaType type) const
void setEtaMax(double eta)
def pv(vc)
Definition: MetAnalyzer.py:7
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
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
const Candidate * daughter(size_type) const override
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
void setProbMin(double p)
void setMassFitMin(double m)
ParameterDescriptionBase * add(U const &iLabel, T const &value)
virtual const reco::Candidate * originalReco(const reco::Candidate *daug) const
get the original particle from the clone
void setMassMin(double m)
void add(std::map< std::string, TH1 * > &h, TH1 *hist)
void setMassMax(double m)
double getConstrSigma() const
void setJPsiMassMin(double m)
set cuts
void add(std::string const &label, ParameterSetDescription const &psetDescription)
std::vector< BPHPlusMinusConstCandPtr > build()
build resonance candidates
void setConstr(double mass, double sigma)
void setConstr(oniaType type, double mass, double sigma)
void setEtaMax(oniaType type, double eta)
T get() const
Definition: EventSetup.h:73
std::vector< BPHPlusMinusConstCandPtr > getList(oniaType type, BPHRecoSelect *dSel=0, BPHMomentumSelect *mSel=0, BPHVertexSelect *vSel=0, BPHFitSelect *kSel=0)
const_iterator begin() const
first daughter const_iterator
Definition: Candidate.h:144
size_type numberOfDaughters() const override
number of daughters
std::vector< BPHRecoConstCandPtr > build()
build Bs candidates
void setMassMax(double m)
virtual const reco::Candidate * getDaug(const std::string &name) const
void setMassMin(oniaType type, double m)
void produce(edm::Event &ev, const edm::EventSetup &es) override
Analysis-level muon class.
Definition: Muon.h:51
void setRecoParameters(const edm::ParameterSet &ps)
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:30
void setPtMin(oniaType type, double pt)
set cuts
std::vector< BPHRecoConstCandPtr > build()
build Bs candidates
void setConstr(double mass, double sigma)
void setJPsiMassMax(double m)
double getConstrMass() const
virtual const reco::Vertex & vertex() const
get reconstructed vertex
std::vector< BPHPlusMinusConstCandPtr > build()
build Phi candidates