File indexing completed on 2022-04-07 05:50:32
0001 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0002 #include "DataFormats/PatCandidates/interface/Muon.h"
0003 #include "FWCore/Common/interface/TriggerNames.h"
0004 #include "FWCore/Framework/interface/ESHandle.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticleMasses.h"
0008 #include "HeavyFlavorAnalysis/SpecificDecay/plugins/BPHHistoSpecificDecay.h"
0009 #include "RecoVertex/VertexTools/interface/VertexDistanceXY.h"
0010
0011 #include "TFile.h"
0012 #include "TH1.h"
0013 #include "TMath.h"
0014 #include "TTree.h"
0015 #include "TVector3.h"
0016
0017 #include "Math/VectorUtil.h"
0018 #include <set>
0019 #include <string>
0020 #include <iostream>
0021 #include <fstream>
0022 #include <cmath>
0023
0024 using namespace std;
0025
0026 #define SET_LABEL(NAME, PSET) (NAME = PSET.getParameter<string>(#NAME))
0027
0028
0029
0030
0031 #define CAT3(A, B, C) A##B##C
0032 #define STRING_NX(A) #A
0033 #define STRING(A) STRING_NX(A)
0034 #define CHK_TRIG(RESULTS, NAMES, INDEX, PATH) \
0035 if (NAMES[INDEX].find(STRING(CAT3(HLT_, PATH, _v))) == 0) { \
0036 flag_##PATH = RESULTS->accept(INDEX); \
0037 continue; \
0038 }
0039
0040
0041
0042
0043 class VertexAnalysis {
0044 public:
0045 static double cAlpha(const reco::Vertex* pvtx, const reco::Vertex* svtx, float px, float py) {
0046 TVector3 disp(svtx->x() - pvtx->x(), svtx->y() - pvtx->y(), 0);
0047 TVector3 cmom(px, py, 0);
0048 return disp.Dot(cmom) / (disp.Perp() * cmom.Perp());
0049 }
0050 static double cAlpha(const reco::Vertex* pvtx, const reco::Vertex* svtx, const TVector3& cmom) {
0051 TVector3 disp(svtx->x() - pvtx->x(), svtx->y() - pvtx->y(), 0);
0052 return disp.Dot(cmom) / (disp.Perp() * cmom.Perp());
0053 }
0054 static void dist2D(const reco::Vertex* pvtx,
0055 const reco::Vertex* svtx,
0056 float px,
0057 float py,
0058 float mass,
0059 double& ctauPV,
0060 double& ctauErrPV) {
0061 dist2D(pvtx, svtx, px, py, cAlpha(pvtx, svtx, px, py), mass, ctauPV, ctauErrPV);
0062 return;
0063 }
0064 static void dist2D(const reco::Vertex* pvtx,
0065 const reco::Vertex* svtx,
0066 float px,
0067 float py,
0068 double cosAlpha,
0069 float mass,
0070 double& ctauPV,
0071 double& ctauErrPV) {
0072 TVector3 cmom(px, py, 0);
0073 AlgebraicVector3 vmom(px, py, 0);
0074 VertexDistanceXY vdistXY;
0075 Measurement1D distXY = vdistXY.distance(*svtx, *pvtx);
0076 ctauPV = distXY.value() * cosAlpha * mass / cmom.Perp();
0077 GlobalError sve = svtx->error();
0078 GlobalError pve = pvtx->error();
0079 AlgebraicSymMatrix33 vXYe = sve.matrix() + pve.matrix();
0080 ctauErrPV = sqrt(ROOT::Math::Similarity(vmom, vXYe)) * mass / cmom.Perp2();
0081 return;
0082 }
0083 };
0084
0085 class BPHUserData {
0086 public:
0087 template <class T>
0088 static const T* get(const pat::CompositeCandidate& cand, const string& name) {
0089 if (cand.hasUserData(name))
0090 return cand.userData<T>(name);
0091 return nullptr;
0092 }
0093 static float get(const pat::CompositeCandidate& cand, const string& name, float d = 0.0) {
0094 if (cand.hasUserFloat(name))
0095 return cand.userFloat(name);
0096 return d;
0097 }
0098 template <class T>
0099 static const T* getByRef(const pat::CompositeCandidate& cand, const string& name) {
0100 if (cand.hasUserData(name)) {
0101 typedef edm::Ref<std::vector<T>> objRef;
0102 const objRef* ref = cand.userData<objRef>(name);
0103 if (ref == nullptr)
0104 return nullptr;
0105 if (ref->isNull())
0106 return nullptr;
0107 return ref->get();
0108 }
0109 return nullptr;
0110 }
0111 };
0112
0113 class BPHDaughters {
0114 public:
0115 static vector<const reco::Candidate*> get(const pat::CompositeCandidate& cand, float massMin, float massMax) {
0116 int i;
0117 int n = cand.numberOfDaughters();
0118 vector<const reco::Candidate*> v;
0119 v.reserve(n);
0120 for (i = 0; i < n; ++i) {
0121 const reco::Candidate* dptr = cand.daughter(i);
0122 float mass = dptr->mass();
0123 if ((mass > massMin) && (mass < massMax))
0124 v.push_back(dptr);
0125 }
0126 return v;
0127 }
0128 };
0129
0130 class BPHSoftMuonSelect {
0131 public:
0132 BPHSoftMuonSelect(int cutTrackerLayers = 5,
0133 int cutPixelLayers = 0,
0134 float maxDxy = 0.3,
0135 float maxDz = 20.0,
0136 bool goodMuon = true,
0137 bool highPurity = true)
0138 : cutTL(cutTrackerLayers), cutPL(cutPixelLayers), maxXY(maxDxy), maxZ(maxDz), gM(goodMuon), hP(highPurity) {}
0139 ~BPHSoftMuonSelect() {}
0140 bool accept(const reco::Candidate& cand, const reco::Vertex* pv) const {
0141 const pat::Muon* p = dynamic_cast<const pat::Muon*>(&cand);
0142 if (p == nullptr)
0143 return false;
0144 if (gM && !muon::isGoodMuon(*p, muon::TMOneStationTight))
0145 return false;
0146 if (p->innerTrack()->hitPattern().trackerLayersWithMeasurement() <= cutTL)
0147 return false;
0148 if (p->innerTrack()->hitPattern().pixelLayersWithMeasurement() <= cutPL)
0149 return false;
0150 if (hP && !p->innerTrack()->quality(reco::TrackBase::highPurity))
0151 return false;
0152 if (pv == nullptr)
0153 return true;
0154 const reco::Vertex::Point& pos = pv->position();
0155 if (fabs(p->innerTrack()->dxy(pos)) >= maxXY)
0156 return false;
0157 if (fabs(p->innerTrack()->dz(pos)) >= maxZ)
0158 return false;
0159 return true;
0160 }
0161
0162 private:
0163 const reco::Vertex* pv;
0164 int cutTL;
0165 int cutPL;
0166 float maxXY;
0167 float maxZ;
0168 bool gM;
0169 bool hP;
0170 };
0171
0172 class BPHDaughterSelect : public BPHHistoSpecificDecay::CandidateSelect {
0173 public:
0174 BPHDaughterSelect(float ptMinLoose,
0175 float ptMinTight,
0176 float etaMaxLoose,
0177 float etaMaxTight,
0178 const BPHSoftMuonSelect* softMuonselector = nullptr)
0179 : pLMin(ptMinLoose), pTMin(ptMinTight), eLMax(etaMaxLoose), eTMax(etaMaxTight), sms(softMuonselector) {}
0180 ~BPHDaughterSelect() override {}
0181 bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pv = nullptr) const override {
0182 return accept(cand, pLMin, pTMin, eLMax, eTMax, pv, sms);
0183 }
0184 static bool accept(const pat::CompositeCandidate& cand,
0185 float ptMinLoose,
0186 float ptMinTight,
0187 float etaMaxLoose,
0188 float etaMaxTight,
0189 const reco::Vertex* pv = nullptr,
0190 const BPHSoftMuonSelect* softMuonselector = nullptr) {
0191 const reco::Candidate* dptr0 = cand.daughter(0);
0192 const reco::Candidate* dptr1 = cand.daughter(1);
0193 if (dptr0 == nullptr)
0194 return false;
0195 if (dptr1 == nullptr)
0196 return false;
0197 float pt0 = dptr0->pt();
0198 float pt1 = dptr1->pt();
0199 if ((pt0 < ptMinLoose) || (pt1 < ptMinLoose))
0200 return false;
0201 if ((pt0 < ptMinTight) && (pt1 < ptMinTight))
0202 return false;
0203 float eta0 = fabs(dptr0->eta());
0204 float eta1 = fabs(dptr1->eta());
0205 if ((etaMaxLoose > 0) && ((eta0 > etaMaxLoose) || (eta1 > etaMaxLoose)))
0206 return false;
0207 if ((etaMaxTight > 0) && ((eta0 > etaMaxTight) && (eta1 > etaMaxTight)))
0208 return false;
0209 if (softMuonselector != nullptr) {
0210 const reco::Vertex* pvtx = BPHUserData::getByRef<reco::Vertex>(cand, "primaryVertex");
0211 if (pvtx == nullptr)
0212 return false;
0213 if (!softMuonselector->accept(*dptr0, pvtx))
0214 return false;
0215 if (!softMuonselector->accept(*dptr1, pvtx))
0216 return false;
0217 }
0218 return true;
0219 }
0220
0221 private:
0222 float pLMin;
0223 float pTMin;
0224 float eLMax;
0225 float eTMax;
0226 const BPHSoftMuonSelect* sms;
0227 };
0228
0229 class BPHCompositeBasicSelect : public BPHHistoSpecificDecay::CandidateSelect {
0230 public:
0231 BPHCompositeBasicSelect(
0232 float massMin, float massMax, float ptMin = -1.0, float etaMax = -1.0, float rapidityMax = -1.0)
0233 : mMin(massMin), mMax(massMax), pMin(ptMin), eMax(etaMax), yMax(rapidityMax) {}
0234 ~BPHCompositeBasicSelect() override {}
0235 bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pv = nullptr) const override {
0236 if (((mMin > 0) && (mMax < 0)) || ((mMin < 0) && (mMax > 0)) || ((mMin > 0) && (mMax > 0) && (mMin < mMax))) {
0237 float mass = cand.mass();
0238 if (mass < mMin)
0239 return false;
0240 if ((mMax > 0) && (mass > mMax))
0241 return false;
0242 }
0243 if (cand.pt() < pMin)
0244 return false;
0245 if ((eMax > 0) && (fabs(cand.eta()) > eMax))
0246 return false;
0247 if ((yMax > 0) && (fabs(cand.rapidity()) > yMax))
0248 return false;
0249 return true;
0250 }
0251
0252 private:
0253 float mMin;
0254 float mMax;
0255 float pMin;
0256 float eMax;
0257 float yMax;
0258 };
0259
0260 class BPHFittedBasicSelect : public BPHHistoSpecificDecay::CandidateSelect {
0261 public:
0262 BPHFittedBasicSelect(float massMin, float massMax, float ptMin = -1.0, float etaMax = -1.0, float rapidityMax = -1.0)
0263 : mMin(massMin), mMax(massMax), pMin(ptMin), eMax(etaMax), yMax(rapidityMax) {}
0264 ~BPHFittedBasicSelect() override {}
0265 bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pv = nullptr) const override {
0266 if (!cand.hasUserFloat("fitMass"))
0267 return false;
0268 float mass = cand.userFloat("fitMass");
0269 if (((mMin > 0) && (mMax < 0)) || ((mMin < 0) && (mMax > 0)) || ((mMin > 0) && (mMax > 0) && (mMin < mMax))) {
0270 if (mass < mMin)
0271 return false;
0272 if ((mMax > 0) && (mass > mMax))
0273 return false;
0274 }
0275 const Vector3DBase<float, GlobalTag>* fmom = BPHUserData::get<Vector3DBase<float, GlobalTag>>(cand, "fitMomentum");
0276 if (fmom == nullptr)
0277 return false;
0278 if (pMin > 0) {
0279 if (fmom->transverse() < pMin)
0280 return false;
0281 }
0282 if (eMax > 0) {
0283 if (fabs(fmom->eta()) > eMax)
0284 return false;
0285 }
0286 if (yMax > 0) {
0287 float x = fmom->x();
0288 float y = fmom->y();
0289 float z = fmom->z();
0290 float e = sqrt((x * x) + (y * y) + (z * z) + (mass * mass));
0291 float r = log((e + z) / (e - z)) / 2;
0292 if (fabs(r) > yMax)
0293 return false;
0294 }
0295 return true;
0296 }
0297
0298 private:
0299 float mMin;
0300 float mMax;
0301 float pMin;
0302 float eMax;
0303 float yMax;
0304 };
0305
0306 class BPHGenericVertexSelect : public BPHHistoSpecificDecay::CandidateSelect {
0307 public:
0308 BPHGenericVertexSelect(char vType, float probMin, float cosMin = -2.0, float sigMin = -1.0, char dMode = 'r')
0309 : type(vType), pMin(probMin), cMin(cosMin), sMin(sigMin), mode(dMode) {}
0310 ~BPHGenericVertexSelect() override {}
0311 bool accept(const pat::CompositeCandidate& cand, const reco::Vertex* pvtx) const override {
0312 if (pvtx == nullptr)
0313 return false;
0314 const reco::Vertex* svtx = nullptr;
0315 float px;
0316 float py;
0317 float mass;
0318 switch (type) {
0319 case 'c':
0320 svtx = BPHUserData::get<reco::Vertex>(cand, "vertex");
0321 px = cand.px();
0322 py = cand.py();
0323 mass = cand.mass();
0324 break;
0325 case 'f':
0326 svtx = BPHUserData::get<reco::Vertex>(cand, "fitVertex");
0327 {
0328 const Vector3DBase<float, GlobalTag>* fmom =
0329 BPHUserData::get<Vector3DBase<float, GlobalTag>>(cand, "fitMomentum");
0330 if (fmom == nullptr)
0331 return false;
0332 px = fmom->x();
0333 py = fmom->y();
0334 }
0335 if (!cand.hasUserFloat("fitMass"))
0336 return false;
0337 mass = cand.userFloat("fitMass");
0338 break;
0339 default:
0340 return false;
0341 }
0342 if (svtx == nullptr)
0343 return false;
0344 if (pMin > 0) {
0345 if (ChiSquaredProbability(svtx->chi2(), svtx->ndof()) < pMin)
0346 return false;
0347 }
0348 if ((cMin > -1.0) || (sMin > 0)) {
0349 float cosAlpha = VertexAnalysis::cAlpha(pvtx, svtx, px, py);
0350 if (cosAlpha < cMin)
0351 return false;
0352 if (sMin < 0)
0353 return true;
0354 double ctauPV;
0355 double ctauErrPV;
0356 VertexAnalysis::dist2D(pvtx, svtx, px, py, cosAlpha, mass, ctauPV, ctauErrPV);
0357 float dTest;
0358 switch (mode) {
0359 case 'a':
0360 case 'd':
0361 dTest = ctauPV;
0362 break;
0363 case 'r':
0364 case 's':
0365 default:
0366 dTest = ctauPV / ctauErrPV;
0367 break;
0368 }
0369 if (dTest < sMin)
0370 return false;
0371 }
0372 return true;
0373 }
0374
0375 private:
0376 char type;
0377 float pMin;
0378 float cMin;
0379 float sMin;
0380 char mode;
0381 };
0382
0383 BPHHistoSpecificDecay::BPHHistoSpecificDecay(const edm::ParameterSet& ps) {
0384 useTrig = (!SET_LABEL(trigResultsLabel, ps).empty());
0385 useOnia = (!SET_LABEL(oniaCandsLabel, ps).empty());
0386 useSd = (!SET_LABEL(sdCandsLabel, ps).empty());
0387 useSs = (!SET_LABEL(ssCandsLabel, ps).empty());
0388 useBu = (!SET_LABEL(buCandsLabel, ps).empty());
0389 useBd = (!SET_LABEL(bdCandsLabel, ps).empty());
0390 useBs = (!SET_LABEL(bsCandsLabel, ps).empty());
0391 useK0 = (!SET_LABEL(k0CandsLabel, ps).empty());
0392 useL0 = (!SET_LABEL(l0CandsLabel, ps).empty());
0393 useB0 = (!SET_LABEL(b0CandsLabel, ps).empty());
0394 useLb = (!SET_LABEL(lbCandsLabel, ps).empty());
0395 useBc = (!SET_LABEL(bcCandsLabel, ps).empty());
0396 useX3872 = (!SET_LABEL(x3872CandsLabel, ps).empty());
0397 if (useTrig)
0398 consume<edm::TriggerResults>(trigResultsToken, trigResultsLabel);
0399 if (useOnia)
0400 consume<vector<pat::CompositeCandidate>>(oniaCandsToken, oniaCandsLabel);
0401 if (useSd)
0402 consume<vector<pat::CompositeCandidate>>(sdCandsToken, sdCandsLabel);
0403 if (useSs)
0404 consume<vector<pat::CompositeCandidate>>(ssCandsToken, ssCandsLabel);
0405 if (useBu)
0406 consume<vector<pat::CompositeCandidate>>(buCandsToken, buCandsLabel);
0407 if (useBd)
0408 consume<vector<pat::CompositeCandidate>>(bdCandsToken, bdCandsLabel);
0409 if (useBs)
0410 consume<vector<pat::CompositeCandidate>>(bsCandsToken, bsCandsLabel);
0411 if (useK0)
0412 consume<vector<pat::CompositeCandidate>>(k0CandsToken, k0CandsLabel);
0413 if (useL0)
0414 consume<vector<pat::CompositeCandidate>>(l0CandsToken, l0CandsLabel);
0415 if (useB0)
0416 consume<vector<pat::CompositeCandidate>>(b0CandsToken, b0CandsLabel);
0417 if (useLb)
0418 consume<vector<pat::CompositeCandidate>>(lbCandsToken, lbCandsLabel);
0419 if (useBc)
0420 consume<vector<pat::CompositeCandidate>>(bcCandsToken, bcCandsLabel);
0421 if (useX3872)
0422 consume<vector<pat::CompositeCandidate>>(x3872CandsToken, x3872CandsLabel);
0423
0424 static const BPHSoftMuonSelect* sms = new BPHSoftMuonSelect;
0425
0426
0427
0428 double phiIMassMin = 0.85;
0429 double phiIMassMax = 1.20;
0430 double phiIPtMin = 18.0;
0431 double phiIEtaMax = -1.0;
0432 double phiIYMax = -1.0;
0433 double phiBMassMin = 0.85;
0434 double phiBMassMax = 1.20;
0435 double phiBPtMin = 14.0;
0436 double phiBEtaMax = -1.0;
0437 double phiBYMax = 1.25;
0438 double jPsiIMassMin = 2.80;
0439 double jPsiIMassMax = 3.40;
0440 double jPsiIPtMin = 25.0;
0441 double jPsiIEtaMax = -1.0;
0442 double jPsiIYMax = -1.0;
0443 double jPsiBMassMin = 2.80;
0444 double jPsiBMassMax = 3.40;
0445 double jPsiBPtMin = 20.0;
0446 double jPsiBEtaMax = -1.0;
0447 double jPsiBYMax = 1.25;
0448 double psi2IMassMin = 3.40;
0449 double psi2IMassMax = 4.00;
0450 double psi2IPtMin = 18.0;
0451 double psi2IEtaMax = -1.0;
0452 double psi2IYMax = -1.0;
0453 double psi2BMassMin = 3.40;
0454 double psi2BMassMax = 4.00;
0455 double psi2BPtMin = 10.0;
0456 double psi2BEtaMax = -1.0;
0457 double psi2BYMax = 1.25;
0458 double upsIMassMin = 8.50;
0459 double upsIMassMax = 11.0;
0460 double upsIPtMin = 15.0;
0461 double upsIEtaMax = -1.0;
0462 double upsIYMax = -1.0;
0463 double upsBMassMin = 8.50;
0464 double upsBMassMax = 11.0;
0465 double upsBPtMin = 12.0;
0466
0467
0468 double upsBEtaMax = -1.0;
0469 double upsBYMax = 1.4;
0470
0471 double oniaProbMin = 0.005;
0472 double oniaCosMin = -2.0;
0473 double oniaSigMin = -1.0;
0474
0475 double oniaMuPtMinLoose = 2.0;
0476 double oniaMuPtMinTight = -1.0;
0477 double oniaMuEtaMaxLoose = -1.0;
0478 double oniaMuEtaMaxTight = -1.0;
0479
0480 phiIBasicSelect = new BPHCompositeBasicSelect(phiIMassMin, phiIMassMax, phiIPtMin, phiIEtaMax, phiIYMax);
0481 phiBBasicSelect = new BPHCompositeBasicSelect(phiBMassMin, phiBMassMax, phiBPtMin, phiBEtaMax, phiBYMax);
0482 jPsiIBasicSelect = new BPHCompositeBasicSelect(jPsiIMassMin, jPsiIMassMax, jPsiIPtMin, jPsiIEtaMax, jPsiIYMax);
0483 jPsiBBasicSelect = new BPHCompositeBasicSelect(jPsiBMassMin, jPsiBMassMax, jPsiBPtMin, jPsiBEtaMax, jPsiBYMax);
0484 psi2IBasicSelect = new BPHCompositeBasicSelect(psi2IMassMin, psi2IMassMax, psi2IPtMin, psi2IEtaMax, psi2IYMax);
0485 psi2BBasicSelect = new BPHCompositeBasicSelect(psi2BMassMin, psi2BMassMax, psi2BPtMin, psi2BEtaMax, psi2BYMax);
0486 upsIBasicSelect = new BPHCompositeBasicSelect(upsIMassMin, upsIMassMax, upsIPtMin, upsIEtaMax, upsIYMax);
0487 upsBBasicSelect = new BPHCompositeBasicSelect(upsBMassMin, upsBMassMax, upsBPtMin, upsBEtaMax, upsBYMax);
0488
0489 oniaVertexSelect = new BPHGenericVertexSelect('c', oniaProbMin, oniaCosMin, oniaSigMin);
0490 oniaDaughterSelect =
0491 new BPHDaughterSelect(oniaMuPtMinLoose, oniaMuPtMinTight, oniaMuEtaMaxLoose, oniaMuEtaMaxTight, sms);
0492
0493 double npJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0494 double npJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0495 double npJPsiPtMin = 8.0;
0496 double npJPsiEtaMax = -1.0;
0497 double npJPsiYMax = -1.0;
0498 double npMuPtMinLoose = 4.0;
0499 double npMuPtMinTight = -1.0;
0500 double npMuEtaMaxLoose = 2.2;
0501 double npMuEtaMaxTight = -1.0;
0502
0503 npJPsiBasicSelect = new BPHCompositeBasicSelect(npJPsiMassMin, npJPsiMassMax, npJPsiPtMin, npJPsiEtaMax, npJPsiYMax);
0504 npJPsiDaughterSelect = new BPHDaughterSelect(npMuPtMinLoose, npMuPtMinTight, npMuEtaMaxLoose, npMuEtaMaxTight, sms);
0505
0506
0507
0508 double buIMassMin = 0.0;
0509 double buIMassMax = 999999.0;
0510 double buIPtMin = 27.0;
0511 double buIEtaMax = -1.0;
0512 double buIYMax = -1.0;
0513 double buIJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0514 double buIJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0515 double buIJPsiPtMin = 25.0;
0516 double buIJPsiEtaMax = -1.0;
0517 double buIJPsiYMax = -1.0;
0518 double buIProbMin = 0.15;
0519 double buICosMin = -2.0;
0520 double buISigMin = -1.0;
0521
0522
0523
0524
0525
0526 buIKPtMin = 2.0;
0527
0528 buIBasicSelect = new BPHFittedBasicSelect(buIMassMin, buIMassMax, buIPtMin, buIEtaMax, buIYMax);
0529 buIJPsiBasicSelect =
0530 new BPHCompositeBasicSelect(buIJPsiMassMin, buIJPsiMassMax, buIJPsiPtMin, buIJPsiEtaMax, buIJPsiYMax);
0531 buIVertexSelect = new BPHGenericVertexSelect('f', buIProbMin, buICosMin, buISigMin);
0532 buIJPsiDaughterSelect = nullptr;
0533
0534
0535
0536
0537 double buDMassMin = 0.0;
0538 double buDMassMax = 999999.0;
0539 double buDPtMin = 10.0;
0540 double buDEtaMax = -1.0;
0541 double buDYMax = -1.0;
0542 double buDJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0543 double buDJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0544 double buDJPsiPtMin = 8.0;
0545 double buDJPsiEtaMax = -1.0;
0546 double buDJPsiYMax = -1.0;
0547 double buDProbMin = 0.10;
0548 double buDCosMin = 0.99;
0549 double buDSigMin = 3.0;
0550
0551
0552
0553
0554
0555 buDKPtMin = 1.6;
0556
0557 buDBasicSelect = new BPHFittedBasicSelect(buDMassMin, buDMassMax, buDPtMin, buDEtaMax, buDYMax);
0558 buDJPsiBasicSelect =
0559 new BPHCompositeBasicSelect(buDJPsiMassMin, buDJPsiMassMax, buDJPsiPtMin, buDJPsiEtaMax, buDJPsiYMax);
0560 buDVertexSelect = new BPHGenericVertexSelect('f', buDProbMin, buDCosMin, buDSigMin);
0561 buDJPsiDaughterSelect = nullptr;
0562
0563
0564
0565
0566
0567
0568 double bdIMassMin = 0.0;
0569 double bdIMassMax = 999999.0;
0570 double bdIPtMin = 27.0;
0571 double bdIEtaMax = -1.0;
0572 double bdIYMax = -1.0;
0573 double bdIJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0574 double bdIJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0575 double bdIJPsiPtMin = 25.0;
0576 double bdIJPsiEtaMax = -1.0;
0577 double bdIJPsiYMax = -1.0;
0578 double bdIKx0MassMin = BPHParticleMasses::kx0Mass - 0.100;
0579 double bdIKx0MassMax = BPHParticleMasses::kx0Mass + 0.100;
0580 double bdIKx0PtMin = -1.0;
0581 double bdIKx0EtaMax = -1.0;
0582 double bdIKx0YMax = -1.0;
0583 double bdIProbMin = 0.15;
0584 double bdICosMin = -2.0;
0585 double bdISigMin = -1.0;
0586
0587
0588
0589
0590
0591 bdIBasicSelect = new BPHFittedBasicSelect(bdIMassMin, bdIMassMax, bdIPtMin, bdIEtaMax, bdIYMax);
0592 bdIJPsiBasicSelect =
0593 new BPHCompositeBasicSelect(bdIJPsiMassMin, bdIJPsiMassMax, bdIJPsiPtMin, bdIJPsiEtaMax, bdIJPsiYMax);
0594 bdIKx0BasicSelect = new BPHCompositeBasicSelect(bdIKx0MassMin, bdIKx0MassMax, bdIKx0PtMin, bdIKx0EtaMax, bdIKx0YMax);
0595 bdIVertexSelect = new BPHGenericVertexSelect('f', bdIProbMin, bdICosMin, bdISigMin);
0596 bdIJPsiDaughterSelect = nullptr;
0597
0598
0599
0600
0601 double bdDMassMin = 0.0;
0602 double bdDMassMax = 999999.0;
0603 double bdDPtMin = 10.0;
0604 double bdDEtaMax = -1.0;
0605 double bdDYMax = -1.0;
0606 double bdDJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0607 double bdDJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0608 double bdDJPsiPtMin = 8.0;
0609 double bdDJPsiEtaMax = -1.0;
0610 double bdDJPsiYMax = -1.0;
0611 double bdDKx0MassMin = BPHParticleMasses::kx0Mass - 0.100;
0612 double bdDKx0MassMax = BPHParticleMasses::kx0Mass + 0.100;
0613 double bdDKx0PtMin = -1.0;
0614 double bdDKx0EtaMax = -1.0;
0615 double bdDKx0YMax = -1.0;
0616 double bdDProbMin = 0.10;
0617 double bdDCosMin = 0.99;
0618 double bdDSigMin = 3.0;
0619
0620
0621
0622
0623
0624 bdDBasicSelect = new BPHFittedBasicSelect(bdDMassMin, bdDMassMax, bdDPtMin, bdDEtaMax, bdDYMax);
0625 bdDJPsiBasicSelect =
0626 new BPHCompositeBasicSelect(bdDJPsiMassMin, bdDJPsiMassMax, bdDJPsiPtMin, bdDJPsiEtaMax, bdDJPsiYMax);
0627 bdDKx0BasicSelect = new BPHCompositeBasicSelect(bdDKx0MassMin, bdDKx0MassMax, bdDKx0PtMin, bdDKx0EtaMax, bdDKx0YMax);
0628 bdDVertexSelect = new BPHGenericVertexSelect('f', bdDProbMin, bdDCosMin, bdDSigMin);
0629 bdDJPsiDaughterSelect = nullptr;
0630
0631
0632
0633
0634
0635
0636 double bsIMassMin = 0.0;
0637 double bsIMassMax = 999999.0;
0638 double bsIPtMin = 27.0;
0639 double bsIEtaMax = -1.0;
0640 double bsIYMax = -1.0;
0641 double bsIJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0642 double bsIJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0643 double bsIJPsiPtMin = 25.0;
0644 double bsIJPsiEtaMax = -1.0;
0645 double bsIJPsiYMax = -1.0;
0646 double bsIPhiMassMin = BPHParticleMasses::phiMass - 0.010;
0647 double bsIPhiMassMax = BPHParticleMasses::phiMass + 0.010;
0648 double bsIPhiPtMin = -1.0;
0649 double bsIPhiEtaMax = -1.0;
0650 double bsIPhiYMax = -1.0;
0651 double bsIProbMin = 0.15;
0652 double bsICosMin = -2.0;
0653 double bsISigMin = -1.0;
0654
0655
0656
0657
0658
0659 bsIBasicSelect = new BPHFittedBasicSelect(bsIMassMin, bsIMassMax, bsIPtMin, bsIEtaMax, bsIYMax);
0660 bsIJPsiBasicSelect =
0661 new BPHCompositeBasicSelect(bsIJPsiMassMin, bsIJPsiMassMax, bsIJPsiPtMin, bsIJPsiEtaMax, bsIJPsiYMax);
0662 bsIPhiBasicSelect = new BPHCompositeBasicSelect(bsIPhiMassMin, bsIPhiMassMax, bsIPhiPtMin, bsIPhiEtaMax, bsIPhiYMax);
0663 bsIVertexSelect = new BPHGenericVertexSelect('f', bsIProbMin, bsICosMin, bsISigMin);
0664 bsIJPsiDaughterSelect = nullptr;
0665
0666
0667
0668
0669 double bsDMassMin = 0.0;
0670 double bsDMassMax = 999999.0;
0671 double bsDPtMin = 10.0;
0672 double bsDEtaMax = -1.0;
0673 double bsDYMax = -1.0;
0674 double bsDJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0675 double bsDJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0676 double bsDJPsiPtMin = 8.0;
0677 double bsDJPsiEtaMax = -1.0;
0678 double bsDJPsiYMax = -1.0;
0679 double bsDPhiMassMin = BPHParticleMasses::phiMass - 0.010;
0680 double bsDPhiMassMax = BPHParticleMasses::phiMass + 0.010;
0681 double bsDPhiPtMin = -1.0;
0682 double bsDPhiEtaMax = -1.0;
0683 double bsDPhiYMax = -1.0;
0684 double bsDProbMin = 0.10;
0685 double bsDCosMin = 0.99;
0686 double bsDSigMin = 3.0;
0687
0688
0689
0690
0691
0692 bsDBasicSelect = new BPHFittedBasicSelect(bsDMassMin, bsDMassMax, bsDPtMin, bsDEtaMax, bsDYMax);
0693 bsDJPsiBasicSelect =
0694 new BPHCompositeBasicSelect(bsDJPsiMassMin, bsDJPsiMassMax, bsDJPsiPtMin, bsDJPsiEtaMax, bsDJPsiYMax);
0695 bsDPhiBasicSelect = new BPHCompositeBasicSelect(bsDPhiMassMin, bsDPhiMassMax, bsDPhiPtMin, bsDPhiEtaMax, bsDPhiYMax);
0696 bsDVertexSelect = new BPHGenericVertexSelect('f', bsDProbMin, bsDCosMin, bsDSigMin);
0697 bsDJPsiDaughterSelect = nullptr;
0698
0699
0700
0701
0702
0703
0704 double b0IMassMin = 0.0;
0705 double b0IMassMax = 999999.0;
0706 double b0IPtMin = 27.0;
0707 double b0IEtaMax = -1.0;
0708 double b0IYMax = -1.0;
0709 double b0IJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0710 double b0IJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0711 double b0IJPsiPtMin = 25.0;
0712 double b0IJPsiEtaMax = -1.0;
0713 double b0IJPsiYMax = -1.0;
0714 double b0IK0sMassMin = BPHParticleMasses::k0sMass - 0.010;
0715 double b0IK0sMassMax = BPHParticleMasses::k0sMass + 0.010;
0716 double b0IK0sPtMin = -1.0;
0717 double b0IK0sEtaMax = -1.0;
0718 double b0IK0sYMax = -1.0;
0719 double b0IProbMin = 0.15;
0720 double b0ICosMin = -2.0;
0721 double b0ISigMin = -1.0;
0722
0723
0724
0725
0726
0727 b0IBasicSelect = new BPHFittedBasicSelect(b0IMassMin, b0IMassMax, b0IPtMin, b0IEtaMax, b0IYMax);
0728 b0IJPsiBasicSelect =
0729 new BPHCompositeBasicSelect(b0IJPsiMassMin, b0IJPsiMassMax, b0IJPsiPtMin, b0IJPsiEtaMax, b0IJPsiYMax);
0730 b0IK0sBasicSelect = new BPHFittedBasicSelect(b0IK0sMassMin, b0IK0sMassMax, b0IK0sPtMin, b0IK0sEtaMax, b0IK0sYMax);
0731 b0IVertexSelect = new BPHGenericVertexSelect('f', b0IProbMin, b0ICosMin, b0ISigMin);
0732 b0IJPsiDaughterSelect = nullptr;
0733
0734
0735
0736
0737 double b0DMassMin = 0.0;
0738 double b0DMassMax = 999999.0;
0739 double b0DPtMin = 10.0;
0740 double b0DEtaMax = -1.0;
0741 double b0DYMax = -1.0;
0742 double b0DJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0743 double b0DJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0744 double b0DJPsiPtMin = 8.0;
0745 double b0DJPsiEtaMax = -1.0;
0746 double b0DJPsiYMax = -1.0;
0747 double b0DK0sMassMin = BPHParticleMasses::k0sMass - 0.010;
0748 double b0DK0sMassMax = BPHParticleMasses::k0sMass + 0.010;
0749 double b0DK0sPtMin = -1.0;
0750 double b0DK0sEtaMax = -1.0;
0751 double b0DK0sYMax = -1.0;
0752 double b0DProbMin = 0.10;
0753 double b0DCosMin = 0.99;
0754 double b0DSigMin = 3.0;
0755
0756
0757
0758
0759
0760 b0DBasicSelect = new BPHFittedBasicSelect(b0DMassMin, b0DMassMax, b0DPtMin, b0DEtaMax, b0DYMax);
0761 b0DJPsiBasicSelect =
0762 new BPHCompositeBasicSelect(b0DJPsiMassMin, b0DJPsiMassMax, b0DJPsiPtMin, b0DJPsiEtaMax, b0DJPsiYMax);
0763 b0DK0sBasicSelect = new BPHFittedBasicSelect(b0DK0sMassMin, b0DK0sMassMax, b0DK0sPtMin, b0DK0sEtaMax, b0DK0sYMax);
0764 b0DVertexSelect = new BPHGenericVertexSelect('f', b0DProbMin, b0DCosMin, b0DSigMin);
0765 b0DJPsiDaughterSelect = nullptr;
0766
0767
0768
0769
0770
0771
0772 double lbIMassMin = 0.0;
0773 double lbIMassMax = 999999.0;
0774 double lbIPtMin = 27.0;
0775 double lbIEtaMax = -1.0;
0776 double lbIYMax = -1.0;
0777 double lbIJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0778 double lbIJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0779 double lbIJPsiPtMin = 25.0;
0780 double lbIJPsiEtaMax = -1.0;
0781 double lbIJPsiYMax = -1.0;
0782 double lbILambda0MassMin = BPHParticleMasses::lambda0Mass - 0.006;
0783 double lbILambda0MassMax = BPHParticleMasses::lambda0Mass + 0.006;
0784 double lbILambda0PtMin = -1.0;
0785 double lbILambda0EtaMax = -1.0;
0786 double lbILambda0YMax = -1.0;
0787 double lbIProbMin = 0.10;
0788 double lbICosMin = -2.0;
0789 double lbISigMin = -1.0;
0790
0791
0792
0793
0794
0795 lbIBasicSelect = new BPHFittedBasicSelect(lbIMassMin, lbIMassMax, lbIPtMin, lbIEtaMax, lbIYMax);
0796 lbIJPsiBasicSelect =
0797 new BPHCompositeBasicSelect(lbIJPsiMassMin, lbIJPsiMassMax, lbIJPsiPtMin, lbIJPsiEtaMax, lbIJPsiYMax);
0798 lbILambda0BasicSelect =
0799 new BPHFittedBasicSelect(lbILambda0MassMin, lbILambda0MassMax, lbILambda0PtMin, lbILambda0EtaMax, lbILambda0YMax);
0800 lbIVertexSelect = new BPHGenericVertexSelect('f', lbIProbMin, lbICosMin, lbISigMin);
0801 lbIJPsiDaughterSelect = nullptr;
0802
0803
0804
0805
0806 double lbDMassMin = 0.0;
0807 double lbDMassMax = 999999.0;
0808 double lbDPtMin = 10.0;
0809 double lbDEtaMax = -1.0;
0810 double lbDYMax = -1.0;
0811 double lbDJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0812 double lbDJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0813 double lbDJPsiPtMin = 8.0;
0814 double lbDJPsiEtaMax = -1.0;
0815 double lbDJPsiYMax = -1.0;
0816 double lbDLambda0MassMin = BPHParticleMasses::lambda0Mass - 0.006;
0817 double lbDLambda0MassMax = BPHParticleMasses::lambda0Mass + 0.006;
0818 double lbDLambda0PtMin = -1.0;
0819 double lbDLambda0EtaMax = -1.0;
0820 double lbDLambda0YMax = -1.0;
0821 double lbDProbMin = 0.10;
0822 double lbDCosMin = 0.99;
0823 double lbDSigMin = 3.0;
0824
0825
0826
0827
0828
0829 lbDBasicSelect = new BPHFittedBasicSelect(lbDMassMin, lbDMassMax, lbDPtMin, lbDEtaMax, lbDYMax);
0830 lbDJPsiBasicSelect =
0831 new BPHCompositeBasicSelect(lbDJPsiMassMin, lbDJPsiMassMax, lbDJPsiPtMin, lbDJPsiEtaMax, lbDJPsiYMax);
0832 lbDLambda0BasicSelect =
0833 new BPHFittedBasicSelect(lbDLambda0MassMin, lbDLambda0MassMax, lbDLambda0PtMin, lbDLambda0EtaMax, lbDLambda0YMax);
0834 lbDVertexSelect = new BPHGenericVertexSelect('f', lbDProbMin, lbDCosMin, lbDSigMin);
0835 lbDJPsiDaughterSelect = nullptr;
0836
0837
0838
0839
0840
0841
0842 double bcIMassMin = 0.0;
0843 double bcIMassMax = 999999.0;
0844 double bcIPtMin = 27.0;
0845 double bcIEtaMax = -1.0;
0846 double bcIYMax = -1.0;
0847 double bcIJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0848 double bcIJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0849 double bcIJPsiPtMin = 25.0;
0850 double bcIJPsiEtaMax = -1.0;
0851 double bcIJPsiYMax = -1.0;
0852 double bcIJPsiProbMin = 0.005;
0853 double bcIProbMin = 0.10;
0854 double bcICosMin = -2.0;
0855 double bcISigMin = -1.0;
0856 double bcIDistMin = 0.01;
0857
0858
0859
0860
0861
0862 bcIPiPtMin = 3.5;
0863
0864 bcIBasicSelect = new BPHFittedBasicSelect(bcIMassMin, bcIMassMax, bcIPtMin, bcIEtaMax, bcIYMax);
0865 bcIJPsiBasicSelect =
0866 new BPHCompositeBasicSelect(bcIJPsiMassMin, bcIJPsiMassMax, bcIJPsiPtMin, bcIJPsiEtaMax, bcIJPsiYMax);
0867 bcIJPsiVertexSelect = new BPHGenericVertexSelect('c', bcIJPsiProbMin);
0868 bcIVertexSelect = new BPHGenericVertexSelect('f', bcIProbMin, bcICosMin, bcISigMin, bcIDistMin);
0869 bcIJPsiDaughterSelect = nullptr;
0870
0871
0872
0873
0874 double bcDMassMin = 0.0;
0875 double bcDMassMax = 999999.0;
0876 double bcDPtMin = 8.0;
0877 double bcDEtaMax = -1.0;
0878 double bcDYMax = -1.0;
0879 double bcDJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0880 double bcDJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0881 double bcDJPsiPtMin = 7.0;
0882 double bcDJPsiEtaMax = -1.0;
0883 double bcDJPsiYMax = -1.0;
0884 double bcDJPsiProbMin = 0.005;
0885 double bcDProbMin = 0.10;
0886 double bcDCosMin = 0.99;
0887 double bcDSigMin = 3.0;
0888
0889
0890
0891
0892
0893 bcJPsiDcaMax = 0.5;
0894 bcDPiPtMin = 3.5;
0895
0896 bcDBasicSelect = new BPHFittedBasicSelect(bcDMassMin, bcDMassMax, bcDPtMin, bcDEtaMax, bcDYMax);
0897 bcDJPsiBasicSelect =
0898 new BPHCompositeBasicSelect(bcDJPsiMassMin, bcDJPsiMassMax, bcDJPsiPtMin, bcDJPsiEtaMax, bcDJPsiYMax);
0899 bcDJPsiVertexSelect = new BPHGenericVertexSelect('c', bcDJPsiProbMin);
0900 bcDVertexSelect = new BPHGenericVertexSelect('f', bcDProbMin, bcDCosMin, bcDSigMin);
0901 bcDJPsiDaughterSelect = nullptr;
0902
0903
0904
0905
0906
0907
0908 double x3872IMassMin = 0.0;
0909 double x3872IMassMax = 999999.0;
0910 double x3872IPtMin = 27.0;
0911 double x3872IEtaMax = -1.0;
0912 double x3872IYMax = -1.0;
0913 double x3872IJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0914 double x3872IJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0915 double x3872IJPsiPtMin = 25.0;
0916 double x3872IJPsiEtaMax = -1.0;
0917 double x3872IJPsiYMax = -1.0;
0918 double x3872IJPsiProbMin = 0.10;
0919 double x3872IProbMin = 0.10;
0920 double x3872ICosMin = -2.0;
0921 double x3872ISigMin = -1.0;
0922 double x3872IDistMin = 0.01;
0923
0924
0925
0926
0927
0928 x3872JPsiDcaMax = 0.5;
0929 x3872IPiPtMin = 1.2;
0930
0931 x3872IBasicSelect = new BPHFittedBasicSelect(x3872IMassMin, x3872IMassMax, x3872IPtMin, x3872IEtaMax, x3872IYMax);
0932 x3872IJPsiBasicSelect = new BPHCompositeBasicSelect(
0933 x3872IJPsiMassMin, x3872IJPsiMassMax, x3872IJPsiPtMin, x3872IJPsiEtaMax, x3872IJPsiYMax);
0934 x3872IJPsiVertexSelect = new BPHGenericVertexSelect('c', x3872IJPsiProbMin);
0935 x3872IVertexSelect = new BPHGenericVertexSelect('f', x3872IProbMin, x3872ICosMin, x3872ISigMin, x3872IDistMin);
0936 x3872IJPsiDaughterSelect = nullptr;
0937
0938
0939
0940
0941
0942 double x3872DMassMin = 0.0;
0943 double x3872DMassMax = 999999.0;
0944 double x3872DPtMin = 8.0;
0945 double x3872DEtaMax = -1.0;
0946 double x3872DYMax = -1.0;
0947 double x3872DJPsiMassMin = BPHParticleMasses::jPsiMass - 0.150;
0948 double x3872DJPsiMassMax = BPHParticleMasses::jPsiMass + 0.150;
0949 double x3872DJPsiPtMin = 7.0;
0950 double x3872DJPsiEtaMax = -1.0;
0951 double x3872DJPsiYMax = -1.0;
0952 double x3872DJPsiProbMin = 0.10;
0953 double x3872DProbMin = 0.10;
0954 double x3872DCosMin = 0.99;
0955 double x3872DSigMin = 3.0;
0956
0957
0958
0959
0960
0961 x3872DPiPtMin = 1.2;
0962
0963 x3872DBasicSelect = new BPHFittedBasicSelect(x3872DMassMin, x3872DMassMax, x3872DPtMin, x3872DEtaMax, x3872DYMax);
0964 x3872DJPsiBasicSelect = new BPHCompositeBasicSelect(
0965 x3872DJPsiMassMin, x3872DJPsiMassMax, x3872DJPsiPtMin, x3872DJPsiEtaMax, x3872DJPsiYMax);
0966 x3872DJPsiVertexSelect = new BPHGenericVertexSelect('c', x3872DJPsiProbMin);
0967 x3872DVertexSelect = new BPHGenericVertexSelect('f', x3872DProbMin, x3872DCosMin, x3872DSigMin);
0968 x3872DJPsiDaughterSelect = nullptr;
0969
0970
0971
0972
0973 }
0974
0975 BPHHistoSpecificDecay::~BPHHistoSpecificDecay() {
0976 delete phiIBasicSelect;
0977 delete phiBBasicSelect;
0978 delete jPsiIBasicSelect;
0979 delete jPsiBBasicSelect;
0980 delete psi2IBasicSelect;
0981 delete psi2BBasicSelect;
0982 delete upsIBasicSelect;
0983 delete upsBBasicSelect;
0984 delete oniaVertexSelect;
0985 delete oniaDaughterSelect;
0986
0987 delete npJPsiBasicSelect;
0988 delete npJPsiDaughterSelect;
0989
0990 delete buIBasicSelect;
0991 delete buIJPsiBasicSelect;
0992 delete buIVertexSelect;
0993 delete buIJPsiDaughterSelect;
0994 delete buDBasicSelect;
0995 delete buDJPsiBasicSelect;
0996 delete buDVertexSelect;
0997 delete buDJPsiDaughterSelect;
0998
0999 delete bdIBasicSelect;
1000 delete bdIJPsiBasicSelect;
1001 delete bdIKx0BasicSelect;
1002 delete bdIVertexSelect;
1003 delete bdIJPsiDaughterSelect;
1004 delete bdDBasicSelect;
1005 delete bdDJPsiBasicSelect;
1006 delete bdDKx0BasicSelect;
1007 delete bdDVertexSelect;
1008 delete bdDJPsiDaughterSelect;
1009
1010 delete bsIBasicSelect;
1011 delete bsIJPsiBasicSelect;
1012 delete bsIPhiBasicSelect;
1013 delete bsIVertexSelect;
1014 delete bsIJPsiDaughterSelect;
1015 delete bsDBasicSelect;
1016 delete bsDJPsiBasicSelect;
1017 delete bsDPhiBasicSelect;
1018 delete bsDVertexSelect;
1019 delete bsDJPsiDaughterSelect;
1020
1021 delete b0IBasicSelect;
1022 delete b0IJPsiBasicSelect;
1023 delete b0IK0sBasicSelect;
1024 delete b0IVertexSelect;
1025 delete b0IJPsiDaughterSelect;
1026 delete b0DBasicSelect;
1027 delete b0DJPsiBasicSelect;
1028 delete b0DK0sBasicSelect;
1029 delete b0DVertexSelect;
1030 delete b0DJPsiDaughterSelect;
1031
1032 delete lbIBasicSelect;
1033 delete lbIJPsiBasicSelect;
1034 delete lbILambda0BasicSelect;
1035 delete lbIVertexSelect;
1036 delete lbIJPsiDaughterSelect;
1037 delete lbDBasicSelect;
1038 delete lbDJPsiBasicSelect;
1039 delete lbDLambda0BasicSelect;
1040 delete lbDVertexSelect;
1041 delete lbDJPsiDaughterSelect;
1042
1043 delete bcIBasicSelect;
1044 delete bcIJPsiBasicSelect;
1045 delete bcIJPsiVertexSelect;
1046 delete bcIVertexSelect;
1047 delete bcIJPsiDaughterSelect;
1048 delete bcDBasicSelect;
1049 delete bcDJPsiBasicSelect;
1050 delete bcDJPsiVertexSelect;
1051 delete bcDVertexSelect;
1052 delete bcDJPsiDaughterSelect;
1053
1054 delete x3872IBasicSelect;
1055 delete x3872IJPsiBasicSelect;
1056 delete x3872IJPsiVertexSelect;
1057 delete x3872IVertexSelect;
1058 delete x3872IJPsiDaughterSelect;
1059 delete x3872DBasicSelect;
1060 delete x3872DJPsiBasicSelect;
1061 delete x3872DJPsiVertexSelect;
1062 delete x3872DVertexSelect;
1063 delete x3872DJPsiDaughterSelect;
1064 }
1065
1066 void BPHHistoSpecificDecay::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
1067 edm::ParameterSetDescription desc;
1068 desc.add<string>("trigResultsLabel", "");
1069 desc.add<string>("oniaCandsLabel", "");
1070 desc.add<string>("sdCandsLabel", "");
1071 desc.add<string>("ssCandsLabel", "");
1072 desc.add<string>("buCandsLabel", "");
1073 desc.add<string>("bdCandsLabel", "");
1074 desc.add<string>("bsCandsLabel", "");
1075 desc.add<string>("k0CandsLabel", "");
1076 desc.add<string>("l0CandsLabel", "");
1077 desc.add<string>("b0CandsLabel", "");
1078 desc.add<string>("lbCandsLabel", "");
1079 desc.add<string>("bcCandsLabel", "");
1080 desc.add<string>("x3872CandsLabel", "");
1081 descriptions.add("process.bphHistoSpecificDecay", desc);
1082 return;
1083 }
1084
1085 void BPHHistoSpecificDecay::beginJob() {
1086 createHisto("massDIPhi", 50, 0.90, 1.15);
1087 createHisto("massTIPhi", 50, 0.90, 1.15);
1088 createHisto("massDBPhi", 50, 0.90, 1.15);
1089 createHisto("massTBPhi", 50, 0.90, 1.15);
1090 createHisto("massDIJPsi", 35, 2.95, 3.30);
1091 createHisto("massTIJPsi", 35, 2.95, 3.30);
1092 createHisto("massDBJPsi", 35, 2.95, 3.30);
1093 createHisto("massTBJPsi", 35, 2.95, 3.30);
1094 createHisto("massDIPsi2", 60, 3.40, 4.00);
1095 createHisto("massTIPsi2", 60, 3.40, 4.00);
1096 createHisto("massDBPsi2", 60, 3.40, 4.00);
1097 createHisto("massTBPsi2", 60, 3.40, 4.00);
1098 createHisto("massDIUps123", 115, 8.70, 11.0);
1099 createHisto("massTIUps123", 115, 8.70, 11.0);
1100 createHisto("massDBUps123", 115, 8.70, 11.0);
1101 createHisto("massTBUps123", 115, 8.70, 11.0);
1102 createHisto("massDIBu", 100, 5.00, 6.00);
1103 createHisto("massTIBu", 100, 5.00, 6.00);
1104 createHisto("massDDBu", 100, 5.00, 6.00);
1105 createHisto("massTDBu", 100, 5.00, 6.00);
1106 createHisto("massDIBd", 100, 5.00, 6.00);
1107 createHisto("massTIBd", 100, 5.00, 6.00);
1108 createHisto("massDDBd", 100, 5.00, 6.00);
1109 createHisto("massTDBd", 100, 5.00, 6.00);
1110 createHisto("massDIBs", 100, 5.00, 6.00);
1111 createHisto("massTIBs", 100, 5.00, 6.00);
1112 createHisto("massDDBs", 100, 5.00, 6.00);
1113 createHisto("massTDBs", 100, 5.00, 6.00);
1114 createHisto("massDIBc", 100, 6.00, 7.00);
1115 createHisto("massTIBc", 100, 6.00, 7.00);
1116 createHisto("massDDBc", 100, 6.00, 7.00);
1117 createHisto("massTDBc", 100, 6.00, 7.00);
1118 createHisto("massDIX3872", 40, 3.60, 4.00);
1119 createHisto("massTIX3872", 40, 3.60, 4.00);
1120 createHisto("massDDX3872", 40, 3.60, 4.00);
1121 createHisto("massTDX3872", 40, 3.60, 4.00);
1122 createHisto("mfitDIBu", 100, 5.00, 6.00);
1123 createHisto("mfitTIBu", 100, 5.00, 6.00);
1124 createHisto("mfitDDBu", 100, 5.00, 6.00);
1125 createHisto("mfitTDBu", 100, 5.00, 6.00);
1126 createHisto("mfitDIBd", 100, 5.00, 6.00);
1127 createHisto("mfitTIBd", 100, 5.00, 6.00);
1128 createHisto("mfitDDBd", 100, 5.00, 6.00);
1129 createHisto("mfitTDBd", 100, 5.00, 6.00);
1130 createHisto("mfitDIBs", 100, 5.00, 6.00);
1131 createHisto("mfitTIBs", 100, 5.00, 6.00);
1132 createHisto("mfitDDBs", 100, 5.00, 6.00);
1133 createHisto("mfitTDBs", 100, 5.00, 6.00);
1134 createHisto("mfitDIBc", 100, 6.00, 7.00);
1135 createHisto("mfitTIBc", 100, 6.00, 7.00);
1136 createHisto("mfitDDBc", 100, 6.00, 7.00);
1137 createHisto("mfitTDBc", 100, 6.00, 7.00);
1138 createHisto("mfitDIX3872", 40, 3.60, 4.00);
1139 createHisto("mfitTIX3872", 40, 3.60, 4.00);
1140 createHisto("mfitDDX3872", 40, 3.60, 4.00);
1141 createHisto("mfitTDX3872", 40, 3.60, 4.00);
1142 createHisto("massDIBuJPsi", 35, 2.95, 3.30);
1143 createHisto("massTIBuJPsi", 35, 2.95, 3.30);
1144 createHisto("massDDBuJPsi", 35, 2.95, 3.30);
1145 createHisto("massTDBuJPsi", 35, 2.95, 3.30);
1146 createHisto("massDIBdJPsi", 35, 2.95, 3.30);
1147 createHisto("massTIBdJPsi", 35, 2.95, 3.30);
1148 createHisto("massDDBdJPsi", 35, 2.95, 3.30);
1149 createHisto("massTDBdJPsi", 35, 2.95, 3.30);
1150 createHisto("massDIBdKx0", 50, 0.80, 1.05);
1151 createHisto("massTIBdKx0", 50, 0.80, 1.05);
1152 createHisto("massDDBdKx0", 50, 0.80, 1.05);
1153 createHisto("massTDBdKx0", 50, 0.80, 1.05);
1154 createHisto("massDIBsJPsi", 35, 2.95, 3.30);
1155 createHisto("massTIBsJPsi", 35, 2.95, 3.30);
1156 createHisto("massDDBsJPsi", 35, 2.95, 3.30);
1157 createHisto("massTDBsJPsi", 35, 2.95, 3.30);
1158 createHisto("massDIBsPhi", 50, 1.01, 1.03);
1159 createHisto("massTIBsPhi", 50, 1.01, 1.03);
1160 createHisto("massDDBsPhi", 50, 1.01, 1.03);
1161 createHisto("massTDBsPhi", 50, 1.01, 1.03);
1162 createHisto("massDIBcJPsi", 35, 2.95, 3.30);
1163 createHisto("massTIBcJPsi", 35, 2.95, 3.30);
1164 createHisto("massDDBcJPsi", 35, 2.95, 3.30);
1165 createHisto("massTDBcJPsi", 35, 2.95, 3.30);
1166 createHisto("massDIX3JPsi", 35, 2.95, 3.30);
1167 createHisto("massTIX3JPsi", 35, 2.95, 3.30);
1168 createHisto("massDDX3JPsi", 35, 2.95, 3.30);
1169 createHisto("massTDX3JPsi", 35, 2.95, 3.30);
1170 createHisto("massDK0s", 50, 0.40, 0.60);
1171 createHisto("mfitDK0s", 50, 0.40, 0.60);
1172 createHisto("massDLambda0", 60, 1.00, 1.30);
1173 createHisto("mfitDLambda0", 60, 1.00, 1.30);
1174 createHisto("massDIB0", 50, 5.00, 6.00);
1175 createHisto("massTIB0", 50, 5.00, 6.00);
1176 createHisto("massDDB0", 50, 5.00, 6.00);
1177 createHisto("massTDB0", 50, 5.00, 6.00);
1178 createHisto("mfitDIB0", 50, 5.00, 6.00);
1179 createHisto("mfitTIB0", 50, 5.00, 6.00);
1180 createHisto("mfitDDB0", 50, 5.00, 6.00);
1181 createHisto("mfitTDB0", 50, 5.00, 6.00);
1182 createHisto("massDIB0JPsi", 35, 2.95, 3.30);
1183 createHisto("massTIB0JPsi", 35, 2.95, 3.30);
1184 createHisto("massDDB0JPsi", 35, 2.95, 3.30);
1185 createHisto("massTDB0JPsi", 35, 2.95, 3.30);
1186 createHisto("massDIB0K0s", 50, 0.40, 0.60);
1187 createHisto("massTIB0K0s", 50, 0.40, 0.60);
1188 createHisto("massDDB0K0s", 50, 0.40, 0.60);
1189 createHisto("massTDB0K0s", 50, 0.40, 0.60);
1190 createHisto("mfitDIB0K0s", 50, 0.40, 0.60);
1191 createHisto("mfitTIB0K0s", 50, 0.40, 0.60);
1192 createHisto("mfitDDB0K0s", 50, 0.40, 0.60);
1193 createHisto("mfitTDB0K0s", 50, 0.40, 0.60);
1194 createHisto("massDILambdab", 25, 5.00, 6.00);
1195 createHisto("massTILambdab", 25, 5.00, 6.00);
1196 createHisto("massDDLambdab", 25, 5.00, 6.00);
1197 createHisto("massTDLambdab", 25, 5.00, 6.00);
1198 createHisto("mfitDILambdab", 25, 5.00, 6.00);
1199 createHisto("mfitTILambdab", 25, 5.00, 6.00);
1200 createHisto("mfitDDLambdab", 25, 5.00, 6.00);
1201 createHisto("mfitTDLambdab", 25, 5.00, 6.00);
1202 createHisto("massDILbJPsi", 35, 2.95, 3.30);
1203 createHisto("massTILbJPsi", 35, 2.95, 3.30);
1204 createHisto("massDDLbJPsi", 35, 2.95, 3.30);
1205 createHisto("massTDLbJPsi", 35, 2.95, 3.30);
1206 createHisto("massDILbL0", 60, 1.00, 1.30);
1207 createHisto("massTILbL0", 60, 1.00, 1.30);
1208 createHisto("massDDLbL0", 60, 1.00, 1.30);
1209 createHisto("massTDLbL0", 60, 1.00, 1.30);
1210 createHisto("mfitDILbL0", 60, 1.00, 1.30);
1211 createHisto("mfitTILbL0", 60, 1.00, 1.30);
1212 createHisto("mfitDDLbL0", 60, 1.00, 1.30);
1213 createHisto("mfitTDLbL0", 60, 1.00, 1.30);
1214
1215 createHisto("massFull", 200, 2.00, 12.0);
1216
1217 createHisto("ctauDIJPsi", 60, -0.05, 0.25);
1218 createHisto("ctauTIJPsi", 60, -0.05, 0.25);
1219 createHisto("ctauDBJPsi", 60, -0.05, 0.25);
1220 createHisto("ctauTBJPsi", 60, -0.05, 0.25);
1221 createHisto("ctauDIBu", 60, -0.05, 0.25);
1222 createHisto("ctauTIBu", 60, -0.05, 0.25);
1223 createHisto("ctauDDBu", 60, -0.05, 0.25);
1224 createHisto("ctauTDBu", 60, -0.05, 0.25);
1225 createHisto("ctauDIBd", 60, -0.05, 0.25);
1226 createHisto("ctauTIBd", 60, -0.05, 0.25);
1227 createHisto("ctauDDBd", 60, -0.05, 0.25);
1228 createHisto("ctauTDBd", 60, -0.05, 0.25);
1229 createHisto("ctauDIBs", 60, -0.05, 0.25);
1230 createHisto("ctauTIBs", 60, -0.05, 0.25);
1231 createHisto("ctauDDBs", 60, -0.05, 0.25);
1232 createHisto("ctauTDBs", 60, -0.05, 0.25);
1233 createHisto("ctauDIB0", 60, -0.05, 0.25);
1234 createHisto("ctauTIB0", 60, -0.05, 0.25);
1235 createHisto("ctauDDB0", 60, -0.05, 0.25);
1236 createHisto("ctauTDB0", 60, -0.05, 0.25);
1237 createHisto("ctauDILambdab", 60, -0.05, 0.25);
1238 createHisto("ctauTILambdab", 60, -0.05, 0.25);
1239 createHisto("ctauDDLambdab", 60, -0.05, 0.25);
1240 createHisto("ctauTDLambdab", 60, -0.05, 0.25);
1241
1242 recoName = new string;
1243 tree = fs->make<TTree>("BPHReco", "BPHReco");
1244 b_runNumber = tree->Branch("runNumber", &runNumber, "runNumber/i");
1245 b_lumiSection = tree->Branch("lumiSection", &lumiSection, "lumiSection/i");
1246 b_eventNumber = tree->Branch("eventNumber", &eventNumber, "eventNumber/i");
1247 b_recoName = tree->Branch("recoName", &recoName, 8192, 99);
1248 b_recoMass = tree->Branch("recoMass", &recoMass, "recoMass/F");
1249 b_recoTime = tree->Branch("recoTime", &recoTime, "recoTime/F");
1250 b_recoErrT = tree->Branch("recoErrT", &recoErrT, "recoErrT/F");
1251
1252 return;
1253 }
1254
1255 void BPHHistoSpecificDecay::analyze(const edm::Event& ev, const edm::EventSetup& es) {
1256
1257
1258 static map<string, ofstream*> ofMap;
1259 if (ofMap.empty()) {
1260 ofMap["BarPhi"] = nullptr;
1261 ofMap["IncJPsi"] = nullptr;
1262 ofMap["BarJPsi"] = nullptr;
1263 ofMap["IncPsi2"] = nullptr;
1264 ofMap["BarPsi2"] = nullptr;
1265 ofMap["BarUpsilon123"] = nullptr;
1266 ofMap["InclusiveBu"] = nullptr;
1267 ofMap["DisplacedBu"] = nullptr;
1268 ofMap["InclusiveBd"] = nullptr;
1269 ofMap["DisplacedBd"] = nullptr;
1270 ofMap["InclusiveBs"] = nullptr;
1271 ofMap["DisplacedBs"] = nullptr;
1272 ofMap["K0s"] = nullptr;
1273 ofMap["Lambda0"] = nullptr;
1274 ofMap["InclusiveB0"] = nullptr;
1275 ofMap["DisplacedB0"] = nullptr;
1276 ofMap["InclusiveLambdab"] = nullptr;
1277 ofMap["DisplacedLambdab"] = nullptr;
1278 ofMap["InclusiveBc"] = nullptr;
1279 ofMap["DisplacedBc"] = nullptr;
1280 ofMap["InclusiveX3872"] = nullptr;
1281 ofMap["DisplacedX3872"] = nullptr;
1282 map<string, ofstream*>::iterator iter = ofMap.begin();
1283 map<string, ofstream*>::iterator iend = ofMap.end();
1284 string name = "list";
1285 while (iter != iend) {
1286 iter->second = new ofstream(name + iter->first);
1287 ++iter;
1288 }
1289 }
1290
1291
1292 runNumber = ev.id().run();
1293 lumiSection = ev.id().luminosityBlock();
1294 eventNumber = ev.id().event();
1295
1296
1297
1298
1299
1300
1301
1302 edm::Handle<edm::TriggerResults> trigResults;
1303 const edm::TriggerNames* trigNames = nullptr;
1304 if (useTrig) {
1305 trigResultsToken.get(ev, trigResults);
1306 if (trigResults.isValid())
1307 trigNames = &ev.triggerNames(*trigResults);
1308 }
1309
1310 bool flag_Dimuon25_Jpsi = false;
1311 bool flag_Dimuon20_Jpsi_Barrel_Seagulls = false;
1312 bool flag_Dimuon14_Phi_Barrel_Seagulls = false;
1313 bool flag_Dimuon18_PsiPrime = false;
1314 bool flag_Dimuon10_PsiPrime_Barrel_Seagulls = false;
1315 bool flag_Dimuon12_Upsilon_eta1p5 = false;
1316 bool flag_Dimuon12_Upsilon_y1p4 = false;
1317 bool flag_DoubleMu4_JpsiTrk_Displaced = false;
1318 if (trigNames != nullptr) {
1319 const edm::TriggerNames::Strings& names = trigNames->triggerNames();
1320 int iObj;
1321 int nObj = names.size();
1322 for (iObj = 0; iObj < nObj; ++iObj) {
1323
1324 CHK_TRIG(trigResults, names, iObj, Dimuon25_Jpsi)
1325 CHK_TRIG(trigResults, names, iObj, Dimuon20_Jpsi_Barrel_Seagulls)
1326 CHK_TRIG(trigResults, names, iObj, Dimuon14_Phi_Barrel_Seagulls)
1327 CHK_TRIG(trigResults, names, iObj, Dimuon18_PsiPrime)
1328 CHK_TRIG(trigResults, names, iObj, Dimuon10_PsiPrime_Barrel_Seagulls)
1329 CHK_TRIG(trigResults, names, iObj, Dimuon12_Upsilon_eta1p5)
1330 CHK_TRIG(trigResults, names, iObj, Dimuon12_Upsilon_y1p4)
1331 CHK_TRIG(trigResults, names, iObj, DoubleMu4_JpsiTrk_Displaced)
1332 }
1333 }
1334
1335
1336
1337
1338
1339
1340
1341
1342
1343
1344
1345
1346
1347
1348
1349
1350
1351
1352 edm::Handle<vector<pat::CompositeCandidate>> oniaCands;
1353 int iqo;
1354 int nqo = 0;
1355 if (useOnia) {
1356 oniaCandsToken.get(ev, oniaCands);
1357 nqo = oniaCands->size();
1358 }
1359
1360 for (iqo = 0; iqo < nqo; ++iqo) {
1361 LogTrace("DataDump") << "*********** quarkonium " << iqo << "/" << nqo << " ***********";
1362 const pat::CompositeCandidate& cand = oniaCands->at(iqo);
1363 if (!oniaVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(cand, "primaryVertex")))
1364 continue;
1365 if (!oniaDaughterSelect->accept(cand))
1366 continue;
1367 fillHisto("Full", cand, 'c');
1368 if (phiBBasicSelect->accept(cand)) {
1369 fillHisto("DBPhi", cand, 'c');
1370 if (flag_Dimuon14_Phi_Barrel_Seagulls)
1371 fillHisto("TBPhi", cand, 'c');
1372 if (flag_Dimuon14_Phi_Barrel_Seagulls)
1373 *ofMap["BarPhi"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1374 << cand.mass() << endl;
1375 }
1376 if (jPsiIBasicSelect->accept(cand)) {
1377 fillHisto("DIJPsi", cand, 'c');
1378 if (flag_Dimuon25_Jpsi)
1379 fillHisto("TIJPsi", cand, 'c');
1380 if (flag_Dimuon25_Jpsi)
1381 *ofMap["IncJPsi"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1382 << cand.mass() << endl;
1383 }
1384 if (jPsiBBasicSelect->accept(cand)) {
1385 fillHisto("DBJPsi", cand, 'c');
1386 if (flag_Dimuon20_Jpsi_Barrel_Seagulls)
1387 fillHisto("TBJPsi", cand, 'c');
1388 if (flag_Dimuon20_Jpsi_Barrel_Seagulls)
1389 *ofMap["BarJPsi"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1390 << cand.mass() << endl;
1391 }
1392 if (psi2IBasicSelect->accept(cand)) {
1393 fillHisto("DIPsi2", cand, 'c');
1394 if (flag_Dimuon18_PsiPrime)
1395 fillHisto("TIPsi2", cand, 'c');
1396 if (flag_Dimuon18_PsiPrime)
1397 *ofMap["IncPsi2"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1398 << cand.mass() << endl;
1399 }
1400 if (psi2BBasicSelect->accept(cand)) {
1401 fillHisto("DBPsi2", cand, 'c');
1402 if (flag_Dimuon10_PsiPrime_Barrel_Seagulls)
1403 fillHisto("TBPsi2", cand, 'c');
1404 if (flag_Dimuon10_PsiPrime_Barrel_Seagulls)
1405 *ofMap["BarPsi2"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1406 << cand.mass() << endl;
1407 }
1408 if (upsBBasicSelect->accept(cand)) {
1409 fillHisto("DBUps123", cand, 'c');
1410 if (flag_Dimuon12_Upsilon_eta1p5 || flag_Dimuon12_Upsilon_y1p4)
1411 fillHisto("TBUps123", cand, 'c');
1412 if (flag_Dimuon12_Upsilon_eta1p5 || flag_Dimuon12_Upsilon_y1p4)
1413 *ofMap["BarUpsilon123"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1414 << cand.mass() << endl;
1415 }
1416 }
1417
1418
1419
1420 edm::Handle<vector<pat::CompositeCandidate>> buCands;
1421 int ibu;
1422 int nbu = 0;
1423 if (useBu) {
1424 buCandsToken.get(ev, buCands);
1425 nbu = buCands->size();
1426 }
1427
1428 for (ibu = 0; ibu < nbu; ++ibu) {
1429 LogTrace("DataDump") << "*********** Bu " << ibu << "/" << nbu << " ***********";
1430 const pat::CompositeCandidate& cand = buCands->at(ibu);
1431 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1432 LogTrace("DataDump") << "JPsi: " << jPsi;
1433 if (jPsi == nullptr)
1434 continue;
1435 if (!npJPsiBasicSelect->accept(*jPsi))
1436 continue;
1437 if (!npJPsiDaughterSelect->accept(*jPsi))
1438 continue;
1439 const reco::Candidate* kptr = BPHDaughters::get(cand, 0.49, 0.50).front();
1440 if (kptr == nullptr)
1441 continue;
1442 if (buIBasicSelect->accept(cand) && buIJPsiBasicSelect->accept(*jPsi) &&
1443
1444 buIVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1445 (kptr->pt() > buIKPtMin)) {
1446 fillHisto("DIBu", cand, 'f');
1447 fillHisto("DIBuJPsi", *jPsi, 'c');
1448 if (flag_Dimuon25_Jpsi) {
1449 fillHisto("TIBu", cand, 'f');
1450 fillHisto("TIBuJPsi", *jPsi, 'c');
1451 *ofMap["InclusiveBu"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1452 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1453 }
1454 }
1455 if (buDBasicSelect->accept(cand) && buDJPsiBasicSelect->accept(*jPsi) &&
1456
1457 buDVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1458 (kptr->pt() > buDKPtMin)) {
1459 fillHisto("DDBu", cand, 'f');
1460 fillHisto("DDBuJPsi", *jPsi, 'c');
1461 if (flag_DoubleMu4_JpsiTrk_Displaced) {
1462 fillHisto("TDBu", cand, 'f');
1463 fillHisto("TDBuJPsi", *jPsi, 'c');
1464 *ofMap["DisplacedBu"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1465 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1466 }
1467 }
1468 }
1469
1470
1471
1472 edm::Handle<vector<pat::CompositeCandidate>> bdCands;
1473 int ibd;
1474 int nbd = 0;
1475 if (useBd) {
1476 bdCandsToken.get(ev, bdCands);
1477 nbd = bdCands->size();
1478 }
1479
1480 for (ibd = 0; ibd < nbd; ++ibd) {
1481 LogTrace("DataDump") << "*********** Bd " << ibd << "/" << nbd << " ***********";
1482 const pat::CompositeCandidate& cand = bdCands->at(ibd);
1483 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1484 LogTrace("DataDump") << "JPsi: " << jPsi;
1485 if (jPsi == nullptr)
1486 continue;
1487 if (!npJPsiBasicSelect->accept(*jPsi))
1488 continue;
1489 if (!npJPsiDaughterSelect->accept(*jPsi))
1490 continue;
1491 const pat::CompositeCandidate* kx0 = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToKx0");
1492 LogTrace("DataDump") << "Kx0: " << kx0;
1493 if (kx0 == nullptr)
1494 continue;
1495 if (bdIBasicSelect->accept(cand) && bdIJPsiBasicSelect->accept(*jPsi) && bdIKx0BasicSelect->accept(*kx0) &&
1496
1497 bdIVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1498 fillHisto("DIBd", cand, 'f');
1499 fillHisto("DIBdJPsi", *jPsi, 'c');
1500 fillHisto("DIBdKx0", *kx0, 'c');
1501 if (flag_Dimuon25_Jpsi) {
1502 fillHisto("TIBd", cand, 'f');
1503 fillHisto("TIBdJPsi", *jPsi, 'c');
1504 fillHisto("TIBdKx0", *kx0, 'c');
1505 *ofMap["InclusiveBd"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1506 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1507 }
1508 }
1509 if (bdDBasicSelect->accept(cand) && bdDJPsiBasicSelect->accept(*jPsi) && bdDKx0BasicSelect->accept(*kx0) &&
1510
1511 bdDVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1512 fillHisto("DDBd", cand, 'f');
1513 fillHisto("DDBdJPsi", *jPsi, 'c');
1514 fillHisto("DDBdKx0", *kx0, 'c');
1515 if (flag_DoubleMu4_JpsiTrk_Displaced) {
1516 fillHisto("TDBd", cand, 'f');
1517 fillHisto("TDBdJPsi", *jPsi, 'c');
1518 fillHisto("TDBdKx0", *kx0, 'c');
1519 *ofMap["DisplacedBd"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1520 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1521 }
1522 }
1523 }
1524
1525
1526
1527 edm::Handle<vector<pat::CompositeCandidate>> bsCands;
1528 int ibs;
1529 int nbs = 0;
1530 if (useBs) {
1531 bsCandsToken.get(ev, bsCands);
1532 nbs = bsCands->size();
1533 }
1534
1535 for (ibs = 0; ibs < nbs; ++ibs) {
1536 LogTrace("DataDump") << "*********** Bs " << ibs << "/" << nbs << " ***********";
1537 const pat::CompositeCandidate& cand = bsCands->at(ibs);
1538 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1539 LogTrace("DataDump") << "JPsi: " << jPsi;
1540 if (jPsi == nullptr)
1541 continue;
1542 if (!npJPsiBasicSelect->accept(*jPsi))
1543 continue;
1544 if (!npJPsiDaughterSelect->accept(*jPsi))
1545 continue;
1546 const pat::CompositeCandidate* phi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToPhi");
1547 LogTrace("DataDump") << "Phi: " << phi;
1548 if (phi == nullptr)
1549 continue;
1550 if (bsIBasicSelect->accept(cand) && bsIJPsiBasicSelect->accept(*jPsi) && bsIPhiBasicSelect->accept(*phi) &&
1551
1552 bsIVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1553 fillHisto("DIBs", cand, 'f');
1554 fillHisto("DIBsJPsi", *jPsi, 'c');
1555 fillHisto("DIBsPhi", *phi, 'c');
1556 if (flag_Dimuon25_Jpsi) {
1557 fillHisto("TIBs", cand, 'f');
1558 fillHisto("TIBsJPsi", *jPsi, 'c');
1559 fillHisto("TIBsPhi", *phi, 'c');
1560 *ofMap["InclusiveBs"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1561 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1562 }
1563 }
1564 if (bsDBasicSelect->accept(cand) && bsDJPsiBasicSelect->accept(*jPsi) && bsDPhiBasicSelect->accept(*phi) &&
1565
1566 bsDVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1567 fillHisto("DDBs", cand, 'f');
1568 fillHisto("DDBsJPsi", *jPsi, 'c');
1569 fillHisto("DDBsPhi", *phi, 'c');
1570 if (flag_DoubleMu4_JpsiTrk_Displaced) {
1571 fillHisto("TDBs", cand, 'f');
1572 fillHisto("TDBsJPsi", *jPsi, 'c');
1573 fillHisto("TDBsPhi", *phi, 'c');
1574 *ofMap["DisplacedBs"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1575 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1576 }
1577 }
1578 }
1579
1580
1581
1582 edm::Handle<vector<pat::CompositeCandidate>> k0Cands;
1583 int ik0;
1584 int nk0 = 0;
1585 if (useK0) {
1586 k0CandsToken.get(ev, k0Cands);
1587 nk0 = k0Cands->size();
1588 }
1589
1590 for (ik0 = 0; ik0 < nk0; ++ik0) {
1591 LogTrace("DataDump") << "*********** K0 " << ik0 << "/" << nk0 << " ***********";
1592 const pat::CompositeCandidate& cand = k0Cands->at(ik0);
1593 fillHisto("DK0s", cand, 'f');
1594 *ofMap["K0s"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1595 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1596 }
1597
1598
1599
1600 edm::Handle<vector<pat::CompositeCandidate>> l0Cands;
1601 int il0;
1602 int nl0 = 0;
1603 if (useL0) {
1604 l0CandsToken.get(ev, l0Cands);
1605 nl0 = l0Cands->size();
1606 }
1607
1608 for (il0 = 0; il0 < nl0; ++il0) {
1609 LogTrace("DataDump") << "*********** Lambda0 " << il0 << "/" << nl0 << " ***********";
1610 const pat::CompositeCandidate& cand = l0Cands->at(il0);
1611 fillHisto("DLambda0", cand, 'f');
1612 *ofMap["Lambda0"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1613 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1614 }
1615
1616
1617
1618 edm::Handle<vector<pat::CompositeCandidate>> b0Cands;
1619 int ib0;
1620 int nb0 = 0;
1621 if (useB0) {
1622 b0CandsToken.get(ev, b0Cands);
1623 nb0 = b0Cands->size();
1624 }
1625
1626
1627
1628
1629
1630 for (ib0 = 0; ib0 < nb0; ++ib0) {
1631 LogTrace("DataDump") << "*********** B0 " << ib0 << "/" << nb0 << " ***********";
1632 const pat::CompositeCandidate& cand = b0Cands->at(ib0);
1633 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1634 LogTrace("DataDump") << "JPsi: " << jPsi;
1635 if (jPsi == nullptr)
1636 continue;
1637 if (!npJPsiBasicSelect->accept(*jPsi))
1638 continue;
1639 if (!npJPsiDaughterSelect->accept(*jPsi))
1640 continue;
1641 const pat::CompositeCandidate* k0s = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToK0s");
1642 LogTrace("DataDump") << "K0s: " << k0s;
1643 if (k0s == nullptr)
1644 continue;
1645 if (b0IBasicSelect->accept(cand) && b0IJPsiBasicSelect->accept(*jPsi) && b0IK0sBasicSelect->accept(*k0s) &&
1646
1647 b0IVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1648 fillHisto("DIB0", cand, 'f');
1649 fillHisto("DIB0JPsi", *jPsi, 'c');
1650 fillHisto("DIB0K0s", *k0s, 'c');
1651 if (flag_Dimuon25_Jpsi) {
1652 fillHisto("TIB0", cand, 'f');
1653 fillHisto("TIB0JPsi", *jPsi, 'c');
1654 fillHisto("TIB0K0s", *k0s, 'c');
1655 *ofMap["InclusiveB0"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1656 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1657 }
1658 }
1659 if (b0DBasicSelect->accept(cand) && b0DJPsiBasicSelect->accept(*jPsi) && b0DK0sBasicSelect->accept(*k0s) &&
1660
1661 b0DVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1662 fillHisto("DDB0", cand, 'f');
1663 fillHisto("DDB0JPsi", *jPsi, 'c');
1664 fillHisto("DDB0K0s", *k0s, 'c');
1665 if (flag_DoubleMu4_JpsiTrk_Displaced) {
1666 fillHisto("TDB0", cand, 'f');
1667 fillHisto("TDB0JPsi", *jPsi, 'c');
1668 fillHisto("TDB0K0s", *k0s, 'c');
1669 *ofMap["DisplacedB0"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1670 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1671 }
1672 }
1673 }
1674
1675
1676
1677 edm::Handle<vector<pat::CompositeCandidate>> lbCands;
1678 int ilb;
1679 int nlb = 0;
1680 if (useLb) {
1681 lbCandsToken.get(ev, lbCands);
1682 nlb = lbCands->size();
1683 }
1684
1685 for (ilb = 0; ilb < nlb; ++ilb) {
1686 LogTrace("DataDump") << "*********** Lambdab " << ilb << "/" << nlb << " ***********";
1687 const pat::CompositeCandidate& cand = lbCands->at(ilb);
1688 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1689 LogTrace("DataDump") << "JPsi: " << jPsi;
1690 if (jPsi == nullptr)
1691 continue;
1692 if (!npJPsiBasicSelect->accept(*jPsi))
1693 continue;
1694 if (!npJPsiDaughterSelect->accept(*jPsi))
1695 continue;
1696 const pat::CompositeCandidate* l0 = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToLambda0");
1697 LogTrace("DataDump") << "Lambda0: " << l0;
1698 if (l0 == nullptr)
1699 continue;
1700 if (lbIBasicSelect->accept(cand) && lbIJPsiBasicSelect->accept(*jPsi) && lbILambda0BasicSelect->accept(*l0) &&
1701
1702 lbIVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1703 fillHisto("DILambdab", cand, 'f');
1704 fillHisto("DILbJPsi", *jPsi, 'c');
1705 fillHisto("DILbL0", *l0, 'c');
1706 if (flag_Dimuon25_Jpsi) {
1707 fillHisto("TILambdab", cand, 'f');
1708 fillHisto("TILbJPsi", *jPsi, 'c');
1709 fillHisto("TILbL0", *l0, 'c');
1710 *ofMap["InclusiveLambdab"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1711 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1712 }
1713 }
1714 if (lbDBasicSelect->accept(cand) && lbDJPsiBasicSelect->accept(*jPsi) && lbDLambda0BasicSelect->accept(*l0) &&
1715
1716 lbDVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex"))) {
1717 fillHisto("DDLambdab", cand, 'f');
1718 fillHisto("DDLbJPsi", *jPsi, 'c');
1719 fillHisto("DDLbL0", *l0, 'c');
1720 if (flag_DoubleMu4_JpsiTrk_Displaced) {
1721 fillHisto("TDLambdab", cand, 'f');
1722 fillHisto("TDLbJPsi", *jPsi, 'c');
1723 fillHisto("TDLbL0", *l0, 'c');
1724 *ofMap["DisplacedLambdab"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1725 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1726 }
1727 }
1728 }
1729
1730
1731
1732 edm::Handle<vector<pat::CompositeCandidate>> bcCands;
1733 int ibc;
1734 int nbc = 0;
1735 if (useBc) {
1736 bcCandsToken.get(ev, bcCands);
1737 nbc = bcCands->size();
1738 }
1739
1740 for (ibc = 0; ibc < nbc; ++ibc) {
1741 LogTrace("DataDump") << "*********** Bc " << ibc << "/" << nbc << " ***********";
1742 const pat::CompositeCandidate& cand = bcCands->at(ibc);
1743 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1744 LogTrace("DataDump") << "JPsi: " << jPsi;
1745 if (jPsi == nullptr)
1746 continue;
1747
1748 if (!npJPsiBasicSelect->accept(*jPsi))
1749 continue;
1750 if (!npJPsiDaughterSelect->accept(*jPsi))
1751 continue;
1752 const reco::Candidate* pptr = BPHDaughters::get(cand, 0.13, 0.14).front();
1753 if (pptr == nullptr)
1754 continue;
1755
1756 if (bcIBasicSelect->accept(cand) && bcIJPsiBasicSelect->accept(*jPsi) &&
1757
1758 bcIJPsiVertexSelect->accept(*jPsi, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1759 bcIVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1760 (pptr->pt() > bcIPiPtMin)) {
1761 fillHisto("DIBc", cand, 'f');
1762 fillHisto("DIBcJPsi", *jPsi, 'c');
1763 if (flag_Dimuon25_Jpsi) {
1764 fillHisto("TIBc", cand, 'f');
1765 fillHisto("TIBcJPsi", *jPsi, 'c');
1766 *ofMap["InclusiveBc"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1767 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1768 }
1769 }
1770 if (bcDBasicSelect->accept(cand) && bcDJPsiBasicSelect->accept(*jPsi) &&
1771
1772 bcDVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1773 bcDVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1774 (pptr->pt() > bcDPiPtMin)) {
1775 fillHisto("DDBc", cand, 'f');
1776 fillHisto("DDBcJPsi", *jPsi, 'c');
1777 if (flag_DoubleMu4_JpsiTrk_Displaced) {
1778 fillHisto("TDBc", cand, 'f');
1779 fillHisto("TDBcJPsi", *jPsi, 'c');
1780 *ofMap["DisplacedBc"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1781 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1782 }
1783 }
1784 }
1785
1786
1787
1788 edm::Handle<vector<pat::CompositeCandidate>> x3872Cands;
1789 int ix3872;
1790 int nx3872 = 0;
1791 if (useX3872) {
1792 x3872CandsToken.get(ev, x3872Cands);
1793 nx3872 = x3872Cands->size();
1794 }
1795
1796 for (ix3872 = 0; ix3872 < nx3872; ++ix3872) {
1797 LogTrace("DataDump") << "*********** X3872 " << ix3872 << "/" << nx3872 << " ***********";
1798 const pat::CompositeCandidate& cand = x3872Cands->at(ix3872);
1799 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1800 LogTrace("DataDump") << "JPsi: " << jPsi;
1801 if (jPsi == nullptr)
1802 continue;
1803
1804 if (!npJPsiBasicSelect->accept(*jPsi))
1805 continue;
1806 if (!npJPsiDaughterSelect->accept(*jPsi))
1807 continue;
1808 const reco::Candidate* ppt1 = BPHDaughters::get(cand, 0.13, 0.14)[0];
1809 const reco::Candidate* ppt2 = BPHDaughters::get(cand, 0.13, 0.14)[1];
1810 if (ppt1 == nullptr)
1811 continue;
1812 if (ppt2 == nullptr)
1813 continue;
1814 if (x3872IBasicSelect->accept(cand) && x3872IJPsiBasicSelect->accept(*jPsi) &&
1815
1816 x3872IJPsiVertexSelect->accept(*jPsi, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1817 x3872IVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1818 (ppt1->pt() > x3872IPiPtMin) && (ppt2->pt() > x3872IPiPtMin)) {
1819 fillHisto("DIX3872", cand, 'f');
1820 fillHisto("DIX3872JPsi", *jPsi, 'c');
1821 if (flag_Dimuon25_Jpsi) {
1822 fillHisto("TIX3872", cand, 'f');
1823 fillHisto("TIX3872JPsi", *jPsi, 'c');
1824 *ofMap["InclusiveX3872"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1825 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1826 }
1827 }
1828 if (x3872DBasicSelect->accept(cand) && x3872DJPsiBasicSelect->accept(*jPsi) &&
1829
1830 x3872DVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1831 x3872DVertexSelect->accept(cand, BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex")) &&
1832 (ppt1->pt() > x3872DPiPtMin) && (ppt2->pt() > x3872DPiPtMin)) {
1833 fillHisto("DDX3872", cand, 'f');
1834 fillHisto("DDX3872JPsi", *jPsi, 'c');
1835 if (flag_DoubleMu4_JpsiTrk_Displaced) {
1836 fillHisto("TDX3872", cand, 'f');
1837 fillHisto("TDX3872JPsi", *jPsi, 'c');
1838 *ofMap["DisplacedX3872"] << ev.id().run() << ' ' << ev.id().luminosityBlock() << ' ' << ev.id().event() << ' '
1839 << (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1) << endl;
1840 }
1841 }
1842 }
1843
1844 return;
1845 }
1846
1847 void BPHHistoSpecificDecay::endJob() {
1848
1849 return;
1850 }
1851
1852 void BPHHistoSpecificDecay::fillHisto(const string& name, const pat::CompositeCandidate& cand, char svType) {
1853 *recoName = name;
1854 float mass = (cand.hasUserFloat("fitMass") ? cand.userFloat("fitMass") : -1);
1855 fillHisto("mass" + name, cand.mass());
1856 fillHisto("mfit" + name, mass);
1857 recoMass = mass;
1858 recoTime = -999999.0;
1859 recoErrT = -999999.0;
1860
1861 const reco::Vertex* pvtx = BPHUserData::getByRef<reco::Vertex>(cand, "primaryVertex");
1862 if (pvtx == nullptr) {
1863 const pat::CompositeCandidate* jPsi = BPHUserData::getByRef<pat::CompositeCandidate>(cand, "refToJPsi");
1864 if (jPsi == nullptr)
1865 return;
1866 pvtx = BPHUserData::getByRef<reco::Vertex>(*jPsi, "primaryVertex");
1867 }
1868
1869
1870 if (pvtx != nullptr) {
1871 const reco::Vertex* svtx = nullptr;
1872 if (svType == 'f')
1873 svtx = BPHUserData::get<reco::Vertex>(cand, "fitVertex");
1874 if (svtx == nullptr)
1875 svtx = BPHUserData::get<reco::Vertex>(cand, "vertex");
1876
1877 if (svtx != nullptr) {
1878 float px;
1879 float py;
1880 const Vector3DBase<float, GlobalTag>* fmom =
1881 BPHUserData::get<Vector3DBase<float, GlobalTag>>(cand, "fitMomentum");
1882 if (fmom != nullptr) {
1883 px = fmom->x();
1884 py = fmom->y();
1885 } else {
1886 px = cand.px();
1887 py = cand.py();
1888 }
1889 double ctauPV2;
1890 double ctauErrPV2;
1891 VertexAnalysis::dist2D(pvtx, svtx, px, py, mass, ctauPV2, ctauErrPV2);
1892 recoTime = ctauPV2;
1893 recoErrT = ctauErrPV2;
1894 fillHisto("ctau" + name, ctauPV2);
1895 }
1896 }
1897 tree->Fill();
1898 return;
1899 }
1900
1901 void BPHHistoSpecificDecay::fillHisto(const string& name, float x) {
1902 map<string, TH1F*>::iterator iter = histoMap.find(name);
1903 map<string, TH1F*>::iterator iend = histoMap.end();
1904 if (iter == iend)
1905 return;
1906 iter->second->Fill(x);
1907 return;
1908 }
1909
1910 void BPHHistoSpecificDecay::createHisto(const string& name, int nbin, float hmin, float hmax) {
1911 histoMap[name] = fs->make<TH1F>(name.c_str(), name.c_str(), nbin, hmin, hmax);
1912 return;
1913 }
1914
1915 #include "FWCore/Framework/interface/MakerMacros.h"
1916
1917 DEFINE_FWK_MODULE(BPHHistoSpecificDecay);