19 #include "Math/VectorUtil.h" 33 #define SET_LABEL(NAME, PSET) (NAME = PSET.getParameter<string>(#NAME)) 42 if (cand.hasUserData(
name))
43 return cand.userData<
T>(
name);
65 int n = cand.numberOfDaughters();
66 vector<const reco::Candidate*>
v;
68 for (i = 0; i <
n; ++
i) {
71 if ((mass > massMin) && (mass <
massMax))
81 int cutPixelLayers = 0,
93 if (p->
innerTrack()->hitPattern().trackerLayersWithMeasurement() <= cutTL)
95 if (p->
innerTrack()->hitPattern().pixelLayersWithMeasurement() <= cutPL)
126 : pLMin(ptMinLoose),
pTMin(ptMinTight), eLMax(etaMaxLoose), eTMax(etaMaxTight), sms(softMuonselector) {}
130 if (dptr0 ==
nullptr)
132 if (dptr1 ==
nullptr)
134 float pt0 = dptr0->
pt();
135 float pt1 = dptr1->
pt();
136 if ((pt0 < pLMin) || (pt1 < pLMin))
140 float eta0 = fabs(dptr0->
eta());
141 float eta1 = fabs(dptr1->
eta());
142 if ((eLMax > 0) && ((eta0 > eLMax) || (eta1 > eLMax)))
144 if ((eTMax > 0) && ((eta0 > eTMax) && (eta1 > eTMax)))
146 if (sms !=
nullptr) {
147 const reco::Vertex* pvtx = BPHUserData::getByRef<reco::Vertex>(cand,
"primaryVertex");
150 if (!sms->accept(*dptr0, pvtx))
152 if (!sms->accept(*dptr1, pvtx))
172 if (((mMin > 0) && (mMax < 0)) || ((mMin < 0) && (mMax > 0)) || ((mMin > 0) && (mMax > 0) && (mMin < mMax))) {
176 if ((mMax > 0) && (mass > mMax))
204 if (((mMin > 0) && (mMax < 0)) || ((mMin < 0) && (mMax > 0)) || ((mMin > 0) && (mMax > 0) && (mMin < mMax))) {
207 if ((mMax > 0) && (mass > mMax))
225 float e =
sqrt((x * x) + (y * y) + (z * z) + (mass * mass));
226 float r =
log((e + z) / (e - z)) / 2;
244 :
pMin(probMin), cMin(cosMin), sMin(sigMin) {}
246 const reco::Vertex* svtx = BPHUserData::get<reco::Vertex>(cand,
"vertex");
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());
263 float mass = cand.mass();
267 double ctauPV = distXY.
value() * cosAlpha * mass / cmom.Perp();
271 double ctauErrPV =
sqrt(ROOT::Math::Similarity(vmom, vXYe)) * mass / cmom.Perp2();
272 if ((ctauPV / ctauErrPV) < sMin)
287 :
pMin(probMin), cMin(cosMin), sMin(sigMin) {}
289 const reco::Vertex* svtx = BPHUserData::get<reco::Vertex>(cand,
"fitVertex");
296 if ((cMin > 0) || (sMin > 0)) {
297 TVector3 disp(svtx->
x() - pvtx->
x(), svtx->
y() - pvtx->
y(), 0);
299 BPHUserData::get<Vector3DBase<float, GlobalTag> >(cand,
"fitMomentum");
302 TVector3 cmom(fmom->x(), fmom->y(), 0);
303 float cosAlpha = disp.Dot(cmom) / (disp.Perp() * cmom.Perp());
308 if (!cand.hasUserFloat(
"fitMass"))
310 float mass = cand.userFloat(
"fitMass");
314 double ctauPV = distXY.
value() * cosAlpha * mass / cmom.Perp();
318 double ctauErrPV =
sqrt(ROOT::Math::Similarity(vmom, vXYe)) * mass / cmom.Perp2();
319 if ((ctauPV / ctauErrPV) < sMin)
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());
339 consume<vector<pat::CompositeCandidate> >(oniaCandsToken, oniaCandsLabel);
341 consume<vector<pat::CompositeCandidate> >(sdCandsToken, sdCandsLabel);
343 consume<vector<pat::CompositeCandidate> >(ssCandsToken, ssCandsLabel);
345 consume<vector<pat::CompositeCandidate> >(buCandsToken, buCandsLabel);
347 consume<vector<pat::CompositeCandidate> >(bdCandsToken, bdCandsLabel);
349 consume<vector<pat::CompositeCandidate> >(bsCandsToken, bsCandsLabel);
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;
374 double oniaProbMin = 0.005;
375 double oniaCosMin = -1.0;
376 double oniaSigMin = -1.0;
378 double oniaMuPtMinLoose = -1.0;
379 double oniaMuPtMinTight = -1.0;
380 double oniaMuEtaMaxLoose = -1.0;
381 double oniaMuEtaMaxTight = -1.0;
389 new BPHDaughterSelect(oniaMuPtMinLoose, oniaMuPtMinTight, oniaMuEtaMaxLoose, oniaMuEtaMaxTight, &sms);
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;
406 buJPsiBasicSelect =
new BPHCompositeBasicSelect(buJPsiMassMin, buJPsiMassMax, buJPsiPtMin, buJPsiEtaMax, buJPsiRMax);
408 buJPsiDaughterSelect =
new BPHDaughterSelect(buMuPtMinLoose, buMuPtMinTight, buMuEtaMaxLoose, buMuEtaMaxTight, &sms);
412 double bdJPsiPtMin = 8.0;
413 double bdJPsiEtaMax = -1.0;
414 double bdJPsiRMax = -1.0;
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;
428 bdJPsiBasicSelect =
new BPHCompositeBasicSelect(bdJPsiMassMin, bdJPsiMassMax, bdJPsiPtMin, bdJPsiEtaMax, bdJPsiRMax);
429 bdKx0BasicSelect =
new BPHCompositeBasicSelect(bdKx0MassMin, bdKx0MassMax, bdKx0PtMin, bdKx0EtaMax, bdKx0RMax);
431 bdJPsiDaughterSelect =
new BPHDaughterSelect(bdMuPtMinLoose, bdMuPtMinTight, bdMuEtaMaxLoose, bdMuEtaMaxTight, &sms);
435 double bsJPsiPtMin = 8.0;
436 double bsJPsiEtaMax = -1.0;
437 double bsJPsiRMax = -1.0;
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);
453 bsJPsiDaughterSelect =
new BPHDaughterSelect(bsMuPtMinLoose, bsMuPtMinTight, bsMuEtaMaxLoose, bsMuEtaMaxTight, &sms);
457 delete phiBasicSelect;
458 delete jPsiBasicSelect;
459 delete psi2BasicSelect;
460 delete upsBasicSelect;
461 delete oniaVertexSelect;
462 delete oniaDaughterSelect;
464 delete buJPsiBasicSelect;
465 delete buVertexSelect;
466 delete buJPsiDaughterSelect;
468 delete bdJPsiBasicSelect;
469 delete bdKx0BasicSelect;
470 delete bdVertexSelect;
471 delete bdJPsiDaughterSelect;
473 delete bsJPsiBasicSelect;
474 delete bsPhiBasicSelect;
475 delete bsVertexSelect;
476 delete bsJPsiDaughterSelect;
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);
492 createHisto(
"massPhi", 35, 0.85, 1.20);
493 createHisto(
"massJPsi", 35, 2.95, 3.30);
494 createHisto(
"massPsi2", 60, 3.40, 4.00);
495 createHisto(
"massUps123", 125, 8.50, 11.0);
496 createHisto(
"massBu", 50, 5.00, 6.00);
497 createHisto(
"massBd", 50, 5.00, 6.00);
498 createHisto(
"massBs", 50, 5.00, 6.00);
499 createHisto(
"mfitBu", 50, 5.00, 6.00);
500 createHisto(
"mfitBd", 50, 5.00, 6.00);
501 createHisto(
"mfitBs", 50, 5.00, 6.00);
502 createHisto(
"massBuJPsi", 35, 2.95, 3.30);
503 createHisto(
"massBdJPsi", 35, 2.95, 3.30);
504 createHisto(
"massBsJPsi", 35, 2.95, 3.30);
505 createHisto(
"massBsPhi", 50, 1.01, 1.03);
506 createHisto(
"massBdKx0", 50, 0.80, 1.05);
508 createHisto(
"massFull", 200, 2.00, 12.0);
528 oniaCandsToken.get(ev, oniaCands);
529 nqo = oniaCands->size();
532 for (iqo = 0; iqo < nqo; ++iqo) {
533 LogTrace(
"DataDump") <<
"*********** quarkonium " << iqo <<
"/" << nqo <<
" ***********";
535 if (!oniaVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(cand,
"primaryVertex")))
537 if (!oniaDaughterSelect->accept(cand))
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);
556 buCandsToken.get(ev, buCands);
557 nbu = buCands->size();
560 for (ibu = 0; ibu < nbu; ++ibu) {
561 LogTrace(
"DataDump") <<
"*********** Bu " << ibu <<
"/" << nbu <<
" ***********";
564 LogTrace(
"DataDump") <<
"JPsi: " << jPsi;
567 if (!buJPsiBasicSelect->accept(*jPsi))
569 if (!buJPsiDaughterSelect->accept(*jPsi))
571 if (!buVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi,
"primaryVertex")))
576 if (kptr->
pt() < buKPtMin)
578 fillHisto(
"Bu", cand);
579 fillHisto(
"BuJPsi", *jPsi);
588 bdCandsToken.get(ev, bdCands);
589 nbd = bdCands->size();
592 for (ibd = 0; ibd < nbd; ++ibd) {
593 LogTrace(
"DataDump") <<
"*********** Bd " << ibd <<
"/" << nbd <<
" ***********";
596 LogTrace(
"DataDump") <<
"JPsi: " << jPsi;
600 LogTrace(
"DataDump") <<
"Kx0: " << kx0;
603 if (!bdJPsiBasicSelect->accept(*jPsi))
605 if (!bdKx0BasicSelect->accept(*kx0))
607 if (!bdJPsiDaughterSelect->accept(*jPsi))
609 if (!bdVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi,
"primaryVertex")))
611 fillHisto(
"Bd", cand);
612 fillHisto(
"BdJPsi", *jPsi);
613 fillHisto(
"BdKx0", *kx0);
622 bsCandsToken.get(ev, bsCands);
623 nbs = bsCands->size();
626 for (ibs = 0; ibs < nbs; ++ibs) {
627 LogTrace(
"DataDump") <<
"*********** Bs " << ibs <<
"/" << nbs <<
" ***********";
630 LogTrace(
"DataDump") <<
"JPsi: " << jPsi;
634 LogTrace(
"DataDump") <<
"Phi: " << phi;
637 if (!bsJPsiBasicSelect->accept(*jPsi))
639 if (!bsPhiBasicSelect->accept(*phi))
641 if (!bsJPsiDaughterSelect->accept(*jPsi))
643 if (!bsVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi,
"primaryVertex")))
645 fillHisto(
"Bs", cand);
646 fillHisto(
"BsJPsi", *jPsi);
647 fillHisto(
"BsPhi", *phi);
657 fillHisto(
"mass" +
name, cand.
mass());
658 fillHisto(
"mfit" +
name, mass);
663 map<string, TH1F*>::iterator iter = histoMap.find(
name);
664 map<string, TH1F*>::iterator iend = histoMap.end();
667 iter->second->Fill(x);
672 histoMap[
name] = fs->make<TH1F>(name.c_str(), name.c_str(),
nbin, hmin, hmax);
Analysis-level particle class.
BPHDaughterSelect(float ptMinLoose, float ptMinTight, float etaMaxLoose, float etaMaxTight, const BPHSoftMuonSelect *softMuonselector=0)
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
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.
double pt() const final
transverse momentum
bool accept(const pat::CompositeCandidate &cand, const reco::Vertex *pvtx=0) const override
static const double jPsiMass
const Point & position() const
position
float userFloat(const std::string &key) const
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)
static const double kx0Mass
reco::TrackRef innerTrack() const override
reference to Track reconstructed in the tracker only (reimplemented from reco::Muon) ...
void createHisto(const std::string &name, int nbin, float hmin, float hmax)
double chi2() const
chi-squares
math::XYZPoint Point
point in the space
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...
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.
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
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)
Error error() const
return SMatrix
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
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
const BPHSoftMuonSelect * sms
Analysis-level muon class.
#define SET_LABEL(NAME, PSET)
double mass() const final
mass
void fillHisto(const std::string &name, const pat::CompositeCandidate &cand)
~BPHHistoSpecificDecay() override