File indexing completed on 2022-04-07 05:50:33
0001 #include "DataFormats/PatCandidates/interface/CompositeCandidate.h"
0002 #include "DataFormats/PatCandidates/interface/GenericParticle.h"
0003 #include "DataFormats/PatCandidates/interface/Muon.h"
0004 #include "DataFormats/TrackReco/interface/Track.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHMomentumSelect.h"
0008 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHPlusMinusCandidate.h"
0009 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoBuilder.h"
0010 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoCandidate.h"
0011 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoSelect.h"
0012 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHTrackReference.h"
0013 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHVertexSelect.h"
0014 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHBcToJPsiPiBuilder.h"
0015 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHBdToJPsiKsBuilder.h"
0016 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHBdToJPsiKxBuilder.h"
0017 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHBsToJPsiPhiBuilder.h"
0018 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHBuToJPsiKBuilder.h"
0019 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHChi2Select.h"
0020 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHK0sToPiPiBuilder.h"
0021 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHKx0ToKPiBuilder.h"
0022 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHLambda0ToPPiBuilder.h"
0023 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHLbToJPsiL0Builder.h"
0024 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHMassSelect.h"
0025 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHMuonEtaSelect.h"
0026 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHMuonPtSelect.h"
0027 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHOniaToMuMuBuilder.h"
0028 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticleMasses.h"
0029 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticleNeutralVeto.h"
0030 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticlePtSelect.h"
0031 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHPhiToKKBuilder.h"
0032 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHX3872ToJPsiPiPiBuilder.h"
0033 #include "HeavyFlavorAnalysis/SpecificDecay/plugins/BPHWriteSpecificDecay.h"
0034 #include "TrackingTools/PatternTools/interface/TwoTrackMinimumDistance.h"
0035
0036 #include <set>
0037 #include <string>
0038
0039 using namespace std;
0040
0041 #define SET_PAR(TYPE, NAME, PSET) (NAME = PSET.getParameter<TYPE>(#NAME))
0042
0043
0044
0045
0046 BPHWriteSpecificDecay::BPHWriteSpecificDecay(const edm::ParameterSet& ps) : bFieldToken(esConsumes()) {
0047 usePV = (!SET_PAR(string, pVertexLabel, ps).empty());
0048 usePM = (!SET_PAR(string, patMuonLabel, ps).empty());
0049 useCC = (!SET_PAR(string, ccCandsLabel, ps).empty());
0050 usePF = (!SET_PAR(string, pfCandsLabel, ps).empty());
0051 usePC = (!SET_PAR(string, pcCandsLabel, ps).empty());
0052 useGP = (!SET_PAR(string, gpCandsLabel, ps).empty());
0053 useK0 = (!SET_PAR(string, k0CandsLabel, ps).empty());
0054 useL0 = (!SET_PAR(string, l0CandsLabel, ps).empty());
0055 useKS = (!SET_PAR(string, kSCandsLabel, ps).empty());
0056 useLS = (!SET_PAR(string, lSCandsLabel, ps).empty());
0057 SET_PAR(string, oniaName, ps);
0058 SET_PAR(string, sdName, ps);
0059 SET_PAR(string, ssName, ps);
0060 SET_PAR(string, buName, ps);
0061 SET_PAR(string, bdName, ps);
0062 SET_PAR(string, bsName, ps);
0063 SET_PAR(string, k0Name, ps);
0064 SET_PAR(string, l0Name, ps);
0065 SET_PAR(string, b0Name, ps);
0066 SET_PAR(string, lbName, ps);
0067 SET_PAR(string, bcName, ps);
0068 SET_PAR(string, x3872Name, ps);
0069
0070 SET_PAR(bool, writeMomentum, ps);
0071 SET_PAR(bool, writeVertex, ps);
0072
0073 rMap["Onia"] = Onia;
0074 rMap["PHiMuMu"] = Pmm;
0075 rMap["Psi1"] = Psi1;
0076 rMap["Psi2"] = Psi2;
0077 rMap["Ups"] = Ups;
0078 rMap["Ups1"] = Ups1;
0079 rMap["Ups2"] = Ups2;
0080 rMap["Ups3"] = Ups3;
0081 rMap["Kx0"] = Kx0;
0082 rMap["PhiKK"] = Pkk;
0083 rMap["Bu"] = Bu;
0084 rMap["Bd"] = Bd;
0085 rMap["Bs"] = Bs;
0086 rMap["K0s"] = K0s;
0087 rMap["Lambda0"] = Lambda0;
0088 rMap["B0"] = B0;
0089 rMap["Lambdab"] = Lambdab;
0090 rMap["Bc"] = Bc;
0091 rMap["X3872"] = X3872;
0092
0093 pMap["ptMin"] = ptMin;
0094 pMap["etaMax"] = etaMax;
0095 pMap["mJPsiMin"] = mPsiMin;
0096 pMap["mJPsiMax"] = mPsiMax;
0097 pMap["mKx0Min"] = mKx0Min;
0098 pMap["mKx0Max"] = mKx0Max;
0099 pMap["mPhiMin"] = mPhiMin;
0100 pMap["mPhiMax"] = mPhiMax;
0101 pMap["mK0sMin"] = mK0sMin;
0102 pMap["mK0sMax"] = mK0sMax;
0103 pMap["mLambda0Min"] = mLambda0Min;
0104 pMap["mLambda0Max"] = mLambda0Max;
0105 pMap["massMin"] = massMin;
0106 pMap["massMax"] = massMax;
0107 pMap["probMin"] = probMin;
0108 pMap["massFitMin"] = mFitMin;
0109 pMap["massFitMax"] = mFitMax;
0110 pMap["constrMass"] = constrMass;
0111 pMap["constrSigma"] = constrSigma;
0112
0113 fMap["constrMJPsi"] = constrMJPsi;
0114 fMap["writeCandidate"] = writeCandidate;
0115
0116 recoOnia = recoKx0 = writeKx0 = recoPkk = writePkk = recoBu = writeBu = recoBd = writeBd = recoBs = writeBs =
0117 recoK0s = writeK0s = recoLambda0 = writeLambda0 = recoB0 = writeB0 = recoLambdab = writeLambdab = recoBc =
0118 writeBc = recoX3872 = writeX3872 = false;
0119
0120 writeOnia = true;
0121 const vector<edm::ParameterSet> recoSelect = ps.getParameter<vector<edm::ParameterSet>>("recoSelect");
0122 int iSel;
0123 int nSel = recoSelect.size();
0124 for (iSel = 0; iSel < nSel; ++iSel)
0125 setRecoParameters(recoSelect[iSel]);
0126 if (!recoOnia)
0127 writeOnia = false;
0128
0129 if (recoBu)
0130 recoOnia = true;
0131 if (recoBd)
0132 recoOnia = recoKx0 = true;
0133 if (recoBs)
0134 recoOnia = recoPkk = true;
0135 if (recoB0)
0136 recoOnia = recoK0s = true;
0137 if (recoLambdab)
0138 recoOnia = recoLambda0 = true;
0139 if (recoBc)
0140 recoOnia = true;
0141 if (recoX3872)
0142 recoOnia = true;
0143 if (writeBu)
0144 writeOnia = true;
0145 if (writeBd)
0146 writeOnia = writeKx0 = true;
0147 if (writeBs)
0148 writeOnia = writePkk = true;
0149 if (writeB0)
0150 writeOnia = writeK0s = true;
0151 if (writeLambdab)
0152 writeOnia = writeLambda0 = true;
0153 if (writeBc)
0154 writeOnia = true;
0155 if (writeX3872)
0156 writeOnia = true;
0157
0158 if (usePV)
0159 consume<vector<reco::Vertex>>(pVertexToken, pVertexLabel);
0160 if (usePM)
0161 consume<pat::MuonCollection>(patMuonToken, patMuonLabel);
0162 if (useCC)
0163 consume<vector<pat::CompositeCandidate>>(ccCandsToken, ccCandsLabel);
0164 if (usePF)
0165 consume<vector<reco::PFCandidate>>(pfCandsToken, pfCandsLabel);
0166 if (usePC)
0167 consume<vector<BPHTrackReference::candidate>>(pcCandsToken, pcCandsLabel);
0168 if (useGP)
0169 consume<vector<pat::GenericParticle>>(gpCandsToken, gpCandsLabel);
0170 if (useK0)
0171 consume<vector<reco::VertexCompositeCandidate>>(k0CandsToken, k0CandsLabel);
0172 if (useL0)
0173 consume<vector<reco::VertexCompositeCandidate>>(l0CandsToken, l0CandsLabel);
0174 if (useKS)
0175 consume<vector<reco::VertexCompositePtrCandidate>>(kSCandsToken, kSCandsLabel);
0176 if (useLS)
0177 consume<vector<reco::VertexCompositePtrCandidate>>(lSCandsToken, lSCandsLabel);
0178
0179 if (writeOnia)
0180 produces<pat::CompositeCandidateCollection>(oniaName);
0181 if (writeKx0)
0182 produces<pat::CompositeCandidateCollection>(sdName);
0183 if (writePkk)
0184 produces<pat::CompositeCandidateCollection>(ssName);
0185 if (writeBu)
0186 produces<pat::CompositeCandidateCollection>(buName);
0187 if (writeBd)
0188 produces<pat::CompositeCandidateCollection>(bdName);
0189 if (writeBs)
0190 produces<pat::CompositeCandidateCollection>(bsName);
0191 if (writeK0s)
0192 produces<pat::CompositeCandidateCollection>(k0Name);
0193 if (writeLambda0)
0194 produces<pat::CompositeCandidateCollection>(l0Name);
0195 if (writeB0)
0196 produces<pat::CompositeCandidateCollection>(b0Name);
0197 if (writeLambdab)
0198 produces<pat::CompositeCandidateCollection>(lbName);
0199 if (writeBc)
0200 produces<pat::CompositeCandidateCollection>(bcName);
0201 if (writeX3872)
0202 produces<pat::CompositeCandidateCollection>(x3872Name);
0203 }
0204
0205 void BPHWriteSpecificDecay::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0206 edm::ParameterSetDescription desc;
0207 desc.add<string>("pVertexLabel", "");
0208 desc.add<string>("patMuonLabel", "");
0209 desc.add<string>("ccCandsLabel", "");
0210 desc.add<string>("pfCandsLabel", "");
0211 desc.add<string>("pcCandsLabel", "");
0212 desc.add<string>("gpCandsLabel", "");
0213 desc.add<string>("k0CandsLabel", "");
0214 desc.add<string>("l0CandsLabel", "");
0215 desc.add<string>("kSCandsLabel", "");
0216 desc.add<string>("lSCandsLabel", "");
0217 desc.add<string>("oniaName", "oniaCand");
0218 desc.add<string>("sdName", "kx0Cand");
0219 desc.add<string>("ssName", "phiCand");
0220 desc.add<string>("buName", "buFitted");
0221 desc.add<string>("bdName", "bdFitted");
0222 desc.add<string>("bsName", "bsFitted");
0223 desc.add<string>("k0Name", "k0Fitted");
0224 desc.add<string>("l0Name", "l0Fitted");
0225 desc.add<string>("b0Name", "b0Fitted");
0226 desc.add<string>("lbName", "lbFitted");
0227 desc.add<string>("bcName", "bcFitted");
0228 desc.add<string>("x3872Name", "x3872Fitted");
0229 desc.add<bool>("writeVertex", true);
0230 desc.add<bool>("writeMomentum", true);
0231 edm::ParameterSetDescription dpar;
0232 dpar.add<string>("name");
0233 dpar.add<double>("ptMin", -2.0e35);
0234 dpar.add<double>("etaMax", -2.0e35);
0235 dpar.add<double>("mJPsiMin", -2.0e35);
0236 dpar.add<double>("mJPsiMax", -2.0e35);
0237 dpar.add<double>("mKx0Min", -2.0e35);
0238 dpar.add<double>("mKx0Max", -2.0e35);
0239 dpar.add<double>("mPhiMin", -2.0e35);
0240 dpar.add<double>("mPhiMax", -2.0e35);
0241 dpar.add<double>("mK0sMin", -2.0e35);
0242 dpar.add<double>("mK0sMax", -2.0e35);
0243 dpar.add<double>("mLambda0Min", -2.0e35);
0244 dpar.add<double>("mLambda0Max", -2.0e35);
0245 dpar.add<double>("massMin", -2.0e35);
0246 dpar.add<double>("massMax", -2.0e35);
0247 dpar.add<double>("probMin", -2.0e35);
0248 dpar.add<double>("massFitMin", -2.0e35);
0249 dpar.add<double>("massFitMax", -2.0e35);
0250 dpar.add<double>("constrMass", -2.0e35);
0251 dpar.add<double>("constrSigma", -2.0e35);
0252 dpar.add<bool>("constrMJPsi", true);
0253 dpar.add<bool>("writeCandidate", true);
0254 vector<edm::ParameterSet> rpar;
0255 desc.addVPSet("recoSelect", dpar, rpar);
0256 descriptions.add("bphWriteSpecificDecay", desc);
0257 return;
0258 }
0259
0260 void BPHWriteSpecificDecay::beginJob() { return; }
0261
0262 void BPHWriteSpecificDecay::produce(edm::Event& ev, const edm::EventSetup& es) {
0263 fill(ev, es);
0264 if (writeOnia)
0265 write(ev, lFull, oniaName);
0266 if (writeKx0)
0267 write(ev, lSd, sdName);
0268 if (writePkk)
0269 write(ev, lSs, ssName);
0270 if (writeBu)
0271 write(ev, lBu, buName);
0272 if (writeBd)
0273 write(ev, lBd, bdName);
0274 if (writeBs)
0275 write(ev, lBs, bsName);
0276 if (writeK0s)
0277 write(ev, lK0, k0Name);
0278 if (writeLambda0)
0279 write(ev, lL0, l0Name);
0280 if (writeB0)
0281 write(ev, lB0, b0Name);
0282 if (writeLambdab)
0283 write(ev, lLb, lbName);
0284 if (writeBc)
0285 write(ev, lBc, bcName);
0286 if (writeX3872)
0287 write(ev, lX3872, x3872Name);
0288 return;
0289 }
0290
0291 void BPHWriteSpecificDecay::fill(edm::Event& ev, const edm::EventSetup& es) {
0292 lFull.clear();
0293 lJPsi.clear();
0294 lSd.clear();
0295 lSs.clear();
0296 lBu.clear();
0297 lBd.clear();
0298 lBs.clear();
0299 lK0.clear();
0300 lL0.clear();
0301 lB0.clear();
0302 lLb.clear();
0303 lBc.clear();
0304 lX3872.clear();
0305 jPsiOMap.clear();
0306 daughMap.clear();
0307 pvRefMap.clear();
0308 ccRefMap.clear();
0309
0310
0311 const MagneticField* magneticField = &es.getData(bFieldToken);
0312
0313
0314
0315
0316
0317 edm::Handle<std::vector<reco::Vertex>> pVertices;
0318 pVertexToken.get(ev, pVertices);
0319 int npv = pVertices->size();
0320
0321 int nrc = 0;
0322
0323
0324 edm::Handle<vector<reco::PFCandidate>> pfCands;
0325 if (usePF) {
0326 pfCandsToken.get(ev, pfCands);
0327 nrc = pfCands->size();
0328 }
0329
0330
0331
0332
0333
0334 edm::Handle<vector<BPHTrackReference::candidate>> pcCands;
0335 if (usePC) {
0336 pcCandsToken.get(ev, pcCands);
0337 nrc = pcCands->size();
0338 }
0339
0340
0341 edm::Handle<vector<pat::GenericParticle>> gpCands;
0342 if (useGP) {
0343 gpCandsToken.get(ev, gpCands);
0344 nrc = gpCands->size();
0345 }
0346
0347
0348 edm::Handle<pat::MuonCollection> patMuon;
0349 if (usePM) {
0350 patMuonToken.get(ev, patMuon);
0351 }
0352
0353
0354 edm::Handle<std::vector<reco::VertexCompositeCandidate>> k0Cand;
0355 if (useK0) {
0356 k0CandsToken.get(ev, k0Cand);
0357 }
0358
0359
0360 edm::Handle<std::vector<reco::VertexCompositeCandidate>> l0Cand;
0361 if (useL0) {
0362 l0CandsToken.get(ev, l0Cand);
0363 }
0364
0365
0366 edm::Handle<std::vector<reco::VertexCompositePtrCandidate>> kSCand;
0367 if (useKS) {
0368 kSCandsToken.get(ev, kSCand);
0369 }
0370
0371
0372 edm::Handle<std::vector<reco::VertexCompositePtrCandidate>> lSCand;
0373 if (useLS) {
0374 lSCandsToken.get(ev, lSCand);
0375 }
0376
0377
0378
0379 vector<const reco::Candidate*> muDaugs;
0380 set<const pat::Muon*> muonSet;
0381 typedef multimap<const reco::Candidate*, const pat::CompositeCandidate*> mu_cc_map;
0382 mu_cc_map muCCMap;
0383 if (useCC) {
0384 edm::Handle<vector<pat::CompositeCandidate>> ccCands;
0385 ccCandsToken.get(ev, ccCands);
0386 int n = ccCands->size();
0387 muDaugs.clear();
0388 muDaugs.reserve(n);
0389 muonSet.clear();
0390 set<const pat::Muon*>::const_iterator iter;
0391 set<const pat::Muon*>::const_iterator iend;
0392 int i;
0393 for (i = 0; i < n; ++i) {
0394 const pat::CompositeCandidate& cc = ccCands->at(i);
0395 int j;
0396 int m = cc.numberOfDaughters();
0397 for (j = 0; j < m; ++j) {
0398 const reco::Candidate* dp = cc.daughter(j);
0399 const pat::Muon* mp = dynamic_cast<const pat::Muon*>(dp);
0400 iter = muonSet.begin();
0401 iend = muonSet.end();
0402 bool add = (mp != nullptr) && (muonSet.find(mp) == iend);
0403 while (add && (iter != iend)) {
0404 if (BPHRecoBuilder::sameTrack(mp, *iter++, 1.0e-5))
0405 add = false;
0406 }
0407 if (add)
0408 muonSet.insert(mp);
0409
0410 muCCMap.insert(pair<const reco::Candidate*, const pat::CompositeCandidate*>(dp, &cc));
0411 }
0412 }
0413 iter = muonSet.begin();
0414 iend = muonSet.end();
0415 while (iter != iend)
0416 muDaugs.push_back(*iter++);
0417 }
0418
0419 map<recoType, map<parType, double>>::const_iterator rIter = parMap.begin();
0420 map<recoType, map<parType, double>>::const_iterator rIend = parMap.end();
0421
0422
0423
0424 BPHOniaToMuMuBuilder* onia = nullptr;
0425 if (recoOnia) {
0426 if (usePM)
0427 onia = new BPHOniaToMuMuBuilder(
0428 es, BPHRecoBuilder::createCollection(patMuon, "ingmcf"), BPHRecoBuilder::createCollection(patMuon, "ingmcf"));
0429 else if (useCC)
0430 onia = new BPHOniaToMuMuBuilder(
0431 es, BPHRecoBuilder::createCollection(muDaugs, "ingmcf"), BPHRecoBuilder::createCollection(muDaugs, "ingmcf"));
0432 }
0433
0434 if (onia != nullptr) {
0435 while (rIter != rIend) {
0436 const map<recoType, map<parType, double>>::value_type& rEntry = *rIter++;
0437 recoType rType = rEntry.first;
0438 const map<parType, double>& pMap = rEntry.second;
0439 BPHOniaToMuMuBuilder::oniaType type;
0440 switch (rType) {
0441 case Pmm:
0442 type = BPHOniaToMuMuBuilder::Phi;
0443 break;
0444 case Psi1:
0445 type = BPHOniaToMuMuBuilder::Psi1;
0446 break;
0447 case Psi2:
0448 type = BPHOniaToMuMuBuilder::Psi2;
0449 break;
0450 case Ups:
0451 type = BPHOniaToMuMuBuilder::Ups;
0452 break;
0453 case Ups1:
0454 type = BPHOniaToMuMuBuilder::Ups1;
0455 break;
0456 case Ups2:
0457 type = BPHOniaToMuMuBuilder::Ups2;
0458 break;
0459 case Ups3:
0460 type = BPHOniaToMuMuBuilder::Ups3;
0461 break;
0462 default:
0463 continue;
0464 }
0465 map<parType, double>::const_iterator pIter = pMap.begin();
0466 map<parType, double>::const_iterator pIend = pMap.end();
0467 while (pIter != pIend) {
0468 const map<parType, double>::value_type& pEntry = *pIter++;
0469 parType id = pEntry.first;
0470 double pv = pEntry.second;
0471 switch (id) {
0472 case ptMin:
0473 onia->setPtMin(type, pv);
0474 break;
0475 case etaMax:
0476 onia->setEtaMax(type, pv);
0477 break;
0478 case massMin:
0479 onia->setMassMin(type, pv);
0480 break;
0481 case massMax:
0482 onia->setMassMax(type, pv);
0483 break;
0484 case probMin:
0485 onia->setProbMin(type, pv);
0486 break;
0487 case constrMass:
0488 onia->setConstr(type, pv, onia->getConstrSigma(type));
0489 break;
0490 case constrSigma:
0491 onia->setConstr(type, onia->getConstrMass(type), pv);
0492 break;
0493 default:
0494 break;
0495 }
0496 }
0497 }
0498 lFull = onia->build();
0499 }
0500
0501
0502
0503 int iFull;
0504 int nFull = lFull.size();
0505 map<const BPHRecoCandidate*, const reco::Vertex*> oniaVtxMap;
0506
0507 typedef mu_cc_map::const_iterator mu_cc_iter;
0508 for (iFull = 0; iFull < nFull; ++iFull) {
0509 const reco::Vertex* pVtx = nullptr;
0510 int pvId = 0;
0511 const BPHPlusMinusCandidate* ptr = lFull[iFull].get();
0512 const std::vector<const reco::Candidate*>& daugs = ptr->daughters();
0513
0514
0515
0516 pair<mu_cc_iter, mu_cc_iter> cc0 = muCCMap.equal_range(ptr->originalReco(daugs[0]));
0517 pair<mu_cc_iter, mu_cc_iter> cc1 = muCCMap.equal_range(ptr->originalReco(daugs[1]));
0518 mu_cc_iter iter0 = cc0.first;
0519 mu_cc_iter iend0 = cc0.second;
0520 mu_cc_iter iter1 = cc1.first;
0521 mu_cc_iter iend1 = cc1.second;
0522 while ((iter0 != iend0) && (pVtx == nullptr)) {
0523 const pat::CompositeCandidate* ccp = iter0++->second;
0524 while (iter1 != iend1) {
0525 if (ccp != iter1++->second)
0526 continue;
0527 pVtx = ccp->userData<reco::Vertex>("PVwithmuons");
0528 const reco::Vertex* sVtx = nullptr;
0529 const reco::Vertex::Point& pPos = pVtx->position();
0530 float dMin = 999999.;
0531 int ipv;
0532 for (ipv = 0; ipv < npv; ++ipv) {
0533 const reco::Vertex* tVtx = &pVertices->at(ipv);
0534 const reco::Vertex::Point& tPos = tVtx->position();
0535 float dist = pow(pPos.x() - tPos.x(), 2) + pow(pPos.y() - tPos.y(), 2) + pow(pPos.z() - tPos.z(), 2);
0536 if (dist < dMin) {
0537 dMin = dist;
0538 sVtx = tVtx;
0539 pvId = ipv;
0540 }
0541 }
0542 pVtx = sVtx;
0543 break;
0544 }
0545 }
0546
0547
0548
0549 if (pVtx == nullptr) {
0550 const reco::Vertex::Point& sVtp = ptr->vertex().position();
0551 GlobalPoint cPos(sVtp.x(), sVtp.y(), sVtp.z());
0552 const pat::CompositeCandidate& sCC = ptr->composite();
0553 GlobalVector cDir(sCC.px(), sCC.py(), sCC.pz());
0554 GlobalPoint bPos(0.0, 0.0, 0.0);
0555 GlobalVector bDir(0.0, 0.0, 1.0);
0556 TwoTrackMinimumDistance ttmd;
0557 bool state = ttmd.calculate(GlobalTrajectoryParameters(cPos, cDir, TrackCharge(0), &(*magneticField)),
0558 GlobalTrajectoryParameters(bPos, bDir, TrackCharge(0), &(*magneticField)));
0559 float minDz = 999999.;
0560 float extrapZ = (state ? ttmd.points().first.z() : -9e20);
0561 int ipv;
0562 for (ipv = 0; ipv < npv; ++ipv) {
0563 const reco::Vertex& tVtx = pVertices->at(ipv);
0564 float deltaZ = fabs(extrapZ - tVtx.position().z());
0565 if (deltaZ < minDz) {
0566 minDz = deltaZ;
0567 pVtx = &tVtx;
0568 pvId = ipv;
0569 }
0570 }
0571 }
0572
0573 oniaVtxMap[ptr] = pVtx;
0574 pvRefMap[ptr] = vertex_ref(pVertices, pvId);
0575 }
0576 pVertexToken.get(ev, pVertices);
0577
0578
0579
0580 if (nFull)
0581 lJPsi = onia->getList(BPHOniaToMuMuBuilder::Psi1);
0582
0583 int nJPsi = lJPsi.size();
0584 delete onia;
0585
0586 if (!nJPsi)
0587 return;
0588 if (!nrc)
0589 return;
0590
0591 int ij;
0592 int io;
0593 int nj = lJPsi.size();
0594 int no = lFull.size();
0595 for (ij = 0; ij < nj; ++ij) {
0596 const BPHRecoCandidate* jp = lJPsi[ij].get();
0597 for (io = 0; io < no; ++io) {
0598 const BPHRecoCandidate* oc = lFull[io].get();
0599 if ((jp->originalReco(jp->getDaug("MuPos")) == oc->originalReco(oc->getDaug("MuPos"))) &&
0600 (jp->originalReco(jp->getDaug("MuNeg")) == oc->originalReco(oc->getDaug("MuNeg")))) {
0601 jPsiOMap[jp] = oc;
0602 break;
0603 }
0604 }
0605 }
0606
0607
0608
0609 BPHBuToJPsiKBuilder* bu = nullptr;
0610 if (recoBu) {
0611 if (usePF)
0612 bu = new BPHBuToJPsiKBuilder(es, lJPsi, BPHRecoBuilder::createCollection(pfCands));
0613 else if (usePC)
0614 bu = new BPHBuToJPsiKBuilder(es, lJPsi, BPHRecoBuilder::createCollection(pcCands));
0615 else if (useGP)
0616 bu = new BPHBuToJPsiKBuilder(es, lJPsi, BPHRecoBuilder::createCollection(gpCands));
0617 }
0618
0619 if (bu != nullptr) {
0620 rIter = parMap.find(Bu);
0621 if (rIter != rIend) {
0622 const map<parType, double>& pMap = rIter->second;
0623 map<parType, double>::const_iterator pIter = pMap.begin();
0624 map<parType, double>::const_iterator pIend = pMap.end();
0625 while (pIter != pIend) {
0626 const map<parType, double>::value_type& pEntry = *pIter++;
0627 parType id = pEntry.first;
0628 double pv = pEntry.second;
0629 switch (id) {
0630 case ptMin:
0631 bu->setKPtMin(pv);
0632 break;
0633 case etaMax:
0634 bu->setKEtaMax(pv);
0635 break;
0636 case mPsiMin:
0637 bu->setJPsiMassMin(pv);
0638 break;
0639 case mPsiMax:
0640 bu->setJPsiMassMax(pv);
0641 break;
0642 case massMin:
0643 bu->setMassMin(pv);
0644 break;
0645 case massMax:
0646 bu->setMassMax(pv);
0647 break;
0648 case probMin:
0649 bu->setProbMin(pv);
0650 break;
0651 case mFitMin:
0652 bu->setMassFitMin(pv);
0653 break;
0654 case mFitMax:
0655 bu->setMassFitMax(pv);
0656 break;
0657 case constrMJPsi:
0658 bu->setConstr(pv > 0);
0659 break;
0660 case writeCandidate:
0661 writeBu = (pv > 0);
0662 break;
0663 default:
0664 break;
0665 }
0666 }
0667 }
0668 lBu = bu->build();
0669 delete bu;
0670 }
0671
0672
0673
0674 vector<BPHPlusMinusConstCandPtr> lKx0;
0675 BPHKx0ToKPiBuilder* kx0 = nullptr;
0676 if (recoKx0) {
0677 if (usePF)
0678 kx0 = new BPHKx0ToKPiBuilder(
0679 es, BPHRecoBuilder::createCollection(pfCands), BPHRecoBuilder::createCollection(pfCands));
0680 else if (usePC)
0681 kx0 = new BPHKx0ToKPiBuilder(
0682 es, BPHRecoBuilder::createCollection(pcCands), BPHRecoBuilder::createCollection(pcCands));
0683 else if (useGP)
0684 kx0 = new BPHKx0ToKPiBuilder(
0685 es, BPHRecoBuilder::createCollection(gpCands), BPHRecoBuilder::createCollection(gpCands));
0686 }
0687
0688 if (kx0 != nullptr) {
0689 rIter = parMap.find(Kx0);
0690 if (rIter != rIend) {
0691 const map<parType, double>& pMap = rIter->second;
0692 map<parType, double>::const_iterator pIter = pMap.begin();
0693 map<parType, double>::const_iterator pIend = pMap.end();
0694 while (pIter != pIend) {
0695 const map<parType, double>::value_type& pEntry = *pIter++;
0696 parType id = pEntry.first;
0697 double pv = pEntry.second;
0698 switch (id) {
0699 case ptMin:
0700 kx0->setPtMin(pv);
0701 break;
0702 case etaMax:
0703 kx0->setEtaMax(pv);
0704 break;
0705 case massMin:
0706 kx0->setMassMin(pv);
0707 break;
0708 case massMax:
0709 kx0->setMassMax(pv);
0710 break;
0711 case probMin:
0712 kx0->setProbMin(pv);
0713 break;
0714 case writeCandidate:
0715 writeKx0 = (pv > 0);
0716 break;
0717 default:
0718 break;
0719 }
0720 }
0721 }
0722 lKx0 = kx0->build();
0723 delete kx0;
0724 }
0725
0726 int nKx0 = lKx0.size();
0727
0728
0729
0730 if (recoBd && nKx0) {
0731 BPHBdToJPsiKxBuilder* bd = new BPHBdToJPsiKxBuilder(es, lJPsi, lKx0);
0732 rIter = parMap.find(Bd);
0733 if (rIter != rIend) {
0734 const map<parType, double>& pMap = rIter->second;
0735 map<parType, double>::const_iterator pIter = pMap.begin();
0736 map<parType, double>::const_iterator pIend = pMap.end();
0737 while (pIter != pIend) {
0738 const map<parType, double>::value_type& pEntry = *pIter++;
0739 parType id = pEntry.first;
0740 double pv = pEntry.second;
0741 switch (id) {
0742 case mPsiMin:
0743 bd->setJPsiMassMin(pv);
0744 break;
0745 case mPsiMax:
0746 bd->setJPsiMassMax(pv);
0747 break;
0748 case mKx0Min:
0749 bd->setKxMassMin(pv);
0750 break;
0751 case mKx0Max:
0752 bd->setKxMassMax(pv);
0753 break;
0754 case massMin:
0755 bd->setMassMin(pv);
0756 break;
0757 case massMax:
0758 bd->setMassMax(pv);
0759 break;
0760 case probMin:
0761 bd->setProbMin(pv);
0762 break;
0763 case mFitMin:
0764 bd->setMassFitMin(pv);
0765 break;
0766 case mFitMax:
0767 bd->setMassFitMax(pv);
0768 break;
0769 case constrMJPsi:
0770 bd->setConstr(pv > 0);
0771 break;
0772 case writeCandidate:
0773 writeBd = (pv > 0);
0774 break;
0775 default:
0776 break;
0777 }
0778 }
0779 }
0780
0781 lBd = bd->build();
0782 delete bd;
0783
0784 set<BPHRecoConstCandPtr> sKx0;
0785 int iBd;
0786 int nBd = lBd.size();
0787 for (iBd = 0; iBd < nBd; ++iBd)
0788 sKx0.insert(lBd[iBd]->getComp("Kx0"));
0789 set<BPHRecoConstCandPtr>::const_iterator iter = sKx0.begin();
0790 set<BPHRecoConstCandPtr>::const_iterator iend = sKx0.end();
0791 while (iter != iend)
0792 lSd.push_back(*iter++);
0793 }
0794
0795
0796
0797 vector<BPHPlusMinusConstCandPtr> lPhi;
0798 BPHPhiToKKBuilder* phi = nullptr;
0799 if (recoPkk) {
0800 if (usePF)
0801 phi = new BPHPhiToKKBuilder(
0802 es, BPHRecoBuilder::createCollection(pfCands), BPHRecoBuilder::createCollection(pfCands));
0803 else if (usePC)
0804 phi = new BPHPhiToKKBuilder(
0805 es, BPHRecoBuilder::createCollection(pcCands), BPHRecoBuilder::createCollection(pcCands));
0806 else if (useGP)
0807 phi = new BPHPhiToKKBuilder(
0808 es, BPHRecoBuilder::createCollection(gpCands), BPHRecoBuilder::createCollection(gpCands));
0809 }
0810
0811 if (phi != nullptr) {
0812 rIter = parMap.find(Pkk);
0813 if (rIter != rIend) {
0814 const map<parType, double>& pMap = rIter->second;
0815 map<parType, double>::const_iterator pIter = pMap.begin();
0816 map<parType, double>::const_iterator pIend = pMap.end();
0817 while (pIter != pIend) {
0818 const map<parType, double>::value_type& pEntry = *pIter++;
0819 parType id = pEntry.first;
0820 double pv = pEntry.second;
0821 switch (id) {
0822 case ptMin:
0823 phi->setPtMin(pv);
0824 break;
0825 case etaMax:
0826 phi->setEtaMax(pv);
0827 break;
0828 case massMin:
0829 phi->setMassMin(pv);
0830 break;
0831 case massMax:
0832 phi->setMassMax(pv);
0833 break;
0834 case probMin:
0835 phi->setProbMin(pv);
0836 break;
0837 case writeCandidate:
0838 writePkk = (pv > 0);
0839 break;
0840 default:
0841 break;
0842 }
0843 }
0844 }
0845 lPhi = phi->build();
0846 delete phi;
0847 }
0848
0849 int nPhi = lPhi.size();
0850
0851
0852
0853 if (recoBs && nPhi) {
0854 BPHBsToJPsiPhiBuilder* bs = new BPHBsToJPsiPhiBuilder(es, lJPsi, lPhi);
0855 rIter = parMap.find(Bs);
0856 if (rIter != rIend) {
0857 const map<parType, double>& pMap = rIter->second;
0858 map<parType, double>::const_iterator pIter = pMap.begin();
0859 map<parType, double>::const_iterator pIend = pMap.end();
0860 while (pIter != pIend) {
0861 const map<parType, double>::value_type& pEntry = *pIter++;
0862 parType id = pEntry.first;
0863 double pv = pEntry.second;
0864 switch (id) {
0865 case mPsiMin:
0866 bs->setJPsiMassMin(pv);
0867 break;
0868 case mPsiMax:
0869 bs->setJPsiMassMax(pv);
0870 break;
0871 case mPhiMin:
0872 bs->setPhiMassMin(pv);
0873 break;
0874 case mPhiMax:
0875 bs->setPhiMassMax(pv);
0876 break;
0877 case massMin:
0878 bs->setMassMin(pv);
0879 break;
0880 case massMax:
0881 bs->setMassMax(pv);
0882 break;
0883 case probMin:
0884 bs->setProbMin(pv);
0885 break;
0886 case mFitMin:
0887 bs->setMassFitMin(pv);
0888 break;
0889 case mFitMax:
0890 bs->setMassFitMax(pv);
0891 break;
0892 case constrMJPsi:
0893 bs->setConstr(pv > 0);
0894 break;
0895 case writeCandidate:
0896 writeBs = (pv > 0);
0897 break;
0898 default:
0899 break;
0900 }
0901 }
0902 }
0903
0904 lBs = bs->build();
0905 delete bs;
0906
0907 set<BPHRecoConstCandPtr> sPhi;
0908 int iBs;
0909 int nBs = lBs.size();
0910 for (iBs = 0; iBs < nBs; ++iBs)
0911 sPhi.insert(lBs[iBs]->getComp("Phi"));
0912 set<BPHRecoConstCandPtr>::const_iterator iter = sPhi.begin();
0913 set<BPHRecoConstCandPtr>::const_iterator iend = sPhi.end();
0914 while (iter != iend)
0915 lSs.push_back(*iter++);
0916 }
0917
0918
0919
0920 BPHK0sToPiPiBuilder* k0s = nullptr;
0921 if (recoK0s) {
0922 if (useK0)
0923 k0s = new BPHK0sToPiPiBuilder(es, k0Cand.product(), "cfp");
0924 else if (useKS)
0925 k0s = new BPHK0sToPiPiBuilder(es, kSCand.product(), "cfp");
0926 }
0927 if (k0s != nullptr) {
0928 rIter = parMap.find(K0s);
0929 if (rIter != rIend) {
0930 const map<parType, double>& pMap = rIter->second;
0931 map<parType, double>::const_iterator pIter = pMap.begin();
0932 map<parType, double>::const_iterator pIend = pMap.end();
0933 while (pIter != pIend) {
0934 const map<parType, double>::value_type& pEntry = *pIter++;
0935 parType id = pEntry.first;
0936 double pv = pEntry.second;
0937 switch (id) {
0938 case ptMin:
0939 k0s->setPtMin(pv);
0940 break;
0941 case etaMax:
0942 k0s->setEtaMax(pv);
0943 break;
0944 case massMin:
0945 k0s->setMassMin(pv);
0946 break;
0947 case massMax:
0948 k0s->setMassMax(pv);
0949 break;
0950 case probMin:
0951 k0s->setProbMin(pv);
0952 break;
0953 case writeCandidate:
0954 writeK0s = (pv > 0);
0955 break;
0956 default:
0957 break;
0958 }
0959 }
0960 }
0961 lK0 = k0s->build();
0962 delete k0s;
0963 }
0964
0965 int nK0 = lK0.size();
0966
0967
0968
0969 BPHLambda0ToPPiBuilder* l0s = nullptr;
0970 if (recoLambda0) {
0971 if (useL0)
0972 l0s = new BPHLambda0ToPPiBuilder(es, l0Cand.product(), "cfp");
0973 else if (useLS)
0974 l0s = new BPHLambda0ToPPiBuilder(es, lSCand.product(), "cfp");
0975 }
0976 if (l0s != nullptr) {
0977 rIter = parMap.find(Lambda0);
0978 if (rIter != rIend) {
0979 const map<parType, double>& pMap = rIter->second;
0980 map<parType, double>::const_iterator pIter = pMap.begin();
0981 map<parType, double>::const_iterator pIend = pMap.end();
0982 while (pIter != pIend) {
0983 const map<parType, double>::value_type& pEntry = *pIter++;
0984 parType id = pEntry.first;
0985 double pv = pEntry.second;
0986 switch (id) {
0987 case ptMin:
0988 l0s->setPtMin(pv);
0989 break;
0990 case etaMax:
0991 l0s->setEtaMax(pv);
0992 break;
0993 case massMin:
0994 l0s->setMassMin(pv);
0995 break;
0996 case massMax:
0997 l0s->setMassMax(pv);
0998 break;
0999 case probMin:
1000 l0s->setProbMin(pv);
1001 break;
1002 case writeCandidate:
1003 writeLambda0 = (pv > 0);
1004 break;
1005 default:
1006 break;
1007 }
1008 }
1009 }
1010 lL0 = l0s->build();
1011 delete l0s;
1012 }
1013
1014 int nL0 = lL0.size();
1015
1016
1017
1018 if (recoB0 && nK0) {
1019 BPHBdToJPsiKsBuilder* b0 = new BPHBdToJPsiKsBuilder(es, lJPsi, lK0);
1020 rIter = parMap.find(B0);
1021 if (rIter != rIend) {
1022 const map<parType, double>& pMap = rIter->second;
1023 map<parType, double>::const_iterator pIter = pMap.begin();
1024 map<parType, double>::const_iterator pIend = pMap.end();
1025 while (pIter != pIend) {
1026 const map<parType, double>::value_type& pEntry = *pIter++;
1027 parType id = pEntry.first;
1028 double pv = pEntry.second;
1029 switch (id) {
1030 case mPsiMin:
1031 b0->setJPsiMassMin(pv);
1032 break;
1033 case mPsiMax:
1034 b0->setJPsiMassMax(pv);
1035 break;
1036 case mK0sMin:
1037 b0->setK0MassMin(pv);
1038 break;
1039 case mK0sMax:
1040 b0->setK0MassMax(pv);
1041 break;
1042 case massMin:
1043 b0->setMassMin(pv);
1044 break;
1045 case massMax:
1046 b0->setMassMax(pv);
1047 break;
1048 case probMin:
1049 b0->setProbMin(pv);
1050 break;
1051 case mFitMin:
1052 b0->setMassFitMin(pv);
1053 break;
1054 case mFitMax:
1055 b0->setMassFitMax(pv);
1056 break;
1057 case constrMJPsi:
1058 b0->setConstr(pv > 0);
1059 break;
1060 case writeCandidate:
1061 writeB0 = (pv > 0);
1062 break;
1063 default:
1064 break;
1065 }
1066 }
1067 }
1068
1069 lB0 = b0->build();
1070 const map<const BPHRecoCandidate*, const BPHRecoCandidate*>& b0Map = b0->daughMap();
1071 daughMap.insert(b0Map.begin(), b0Map.end());
1072 delete b0;
1073 }
1074
1075
1076
1077 if (recoLambdab && nL0) {
1078 BPHLbToJPsiL0Builder* lb = new BPHLbToJPsiL0Builder(es, lJPsi, lL0);
1079 rIter = parMap.find(Lambdab);
1080 if (rIter != rIend) {
1081 const map<parType, double>& pMap = rIter->second;
1082 map<parType, double>::const_iterator pIter = pMap.begin();
1083 map<parType, double>::const_iterator pIend = pMap.end();
1084 while (pIter != pIend) {
1085 const map<parType, double>::value_type& pEntry = *pIter++;
1086 parType id = pEntry.first;
1087 double pv = pEntry.second;
1088 switch (id) {
1089 case mPsiMin:
1090 lb->setJPsiMassMin(pv);
1091 break;
1092 case mPsiMax:
1093 lb->setJPsiMassMax(pv);
1094 break;
1095 case mLambda0Min:
1096 lb->setLambda0MassMin(pv);
1097 break;
1098 case mLambda0Max:
1099 lb->setLambda0MassMax(pv);
1100 break;
1101 case massMin:
1102 lb->setMassMin(pv);
1103 break;
1104 case massMax:
1105 lb->setMassMax(pv);
1106 break;
1107 case probMin:
1108 lb->setProbMin(pv);
1109 break;
1110 case mFitMin:
1111 lb->setMassFitMin(pv);
1112 break;
1113 case mFitMax:
1114 lb->setMassFitMax(pv);
1115 break;
1116 case constrMJPsi:
1117 lb->setConstr(pv > 0);
1118 break;
1119 case writeCandidate:
1120 writeLambdab = (pv > 0);
1121 break;
1122 default:
1123 break;
1124 }
1125 }
1126 }
1127
1128 lLb = lb->build();
1129 const map<const BPHRecoCandidate*, const BPHRecoCandidate*>& ldMap = lb->daughMap();
1130 daughMap.insert(ldMap.begin(), ldMap.end());
1131 delete lb;
1132 }
1133
1134
1135
1136 BPHBcToJPsiPiBuilder* bc = nullptr;
1137 if (recoBc) {
1138 if (usePF)
1139 bc = new BPHBcToJPsiPiBuilder(es, lJPsi, BPHRecoBuilder::createCollection(pfCands));
1140 else if (usePC)
1141 bc = new BPHBcToJPsiPiBuilder(es, lJPsi, BPHRecoBuilder::createCollection(pcCands));
1142 else if (useGP)
1143 bc = new BPHBcToJPsiPiBuilder(es, lJPsi, BPHRecoBuilder::createCollection(gpCands));
1144 }
1145
1146 if (bc != nullptr) {
1147 rIter = parMap.find(Bc);
1148 if (rIter != rIend) {
1149 const map<parType, double>& pMap = rIter->second;
1150 map<parType, double>::const_iterator pIter = pMap.begin();
1151 map<parType, double>::const_iterator pIend = pMap.end();
1152 while (pIter != pIend) {
1153 const map<parType, double>::value_type& pEntry = *pIter++;
1154 parType id = pEntry.first;
1155 double pv = pEntry.second;
1156 switch (id) {
1157 case ptMin:
1158 bc->setPiPtMin(pv);
1159 break;
1160 case etaMax:
1161 bc->setPiEtaMax(pv);
1162 break;
1163 case mPsiMin:
1164 bc->setJPsiMassMin(pv);
1165 break;
1166 case mPsiMax:
1167 bc->setJPsiMassMax(pv);
1168 break;
1169 case massMin:
1170 bc->setMassMin(pv);
1171 break;
1172 case massMax:
1173 bc->setMassMax(pv);
1174 break;
1175 case probMin:
1176 bc->setProbMin(pv);
1177 break;
1178 case mFitMin:
1179 bc->setMassFitMin(pv);
1180 break;
1181 case mFitMax:
1182 bc->setMassFitMax(pv);
1183 break;
1184 case constrMJPsi:
1185 bc->setConstr(pv > 0);
1186 break;
1187 case writeCandidate:
1188 writeBc = (pv > 0);
1189 break;
1190 default:
1191 break;
1192 }
1193 }
1194 }
1195 lBc = bc->build();
1196 delete bc;
1197 }
1198
1199
1200
1201 BPHX3872ToJPsiPiPiBuilder* x3872 = nullptr;
1202 if (recoX3872) {
1203 if (usePF)
1204 x3872 = new BPHX3872ToJPsiPiPiBuilder(
1205 es, lJPsi, BPHRecoBuilder::createCollection(pfCands), BPHRecoBuilder::createCollection(pfCands));
1206 else if (usePC)
1207 x3872 = new BPHX3872ToJPsiPiPiBuilder(
1208 es, lJPsi, BPHRecoBuilder::createCollection(pcCands), BPHRecoBuilder::createCollection(pcCands));
1209 else if (useGP)
1210 x3872 = new BPHX3872ToJPsiPiPiBuilder(
1211 es, lJPsi, BPHRecoBuilder::createCollection(gpCands), BPHRecoBuilder::createCollection(gpCands));
1212 }
1213
1214 if (x3872 != nullptr) {
1215 rIter = parMap.find(X3872);
1216 if (rIter != rIend) {
1217 const map<parType, double>& pMap = rIter->second;
1218 map<parType, double>::const_iterator pIter = pMap.begin();
1219 map<parType, double>::const_iterator pIend = pMap.end();
1220 while (pIter != pIend) {
1221 const map<parType, double>::value_type& pEntry = *pIter++;
1222 parType id = pEntry.first;
1223 double pv = pEntry.second;
1224 switch (id) {
1225 case ptMin:
1226 x3872->setPiPtMin(pv);
1227 break;
1228 case etaMax:
1229 x3872->setPiEtaMax(pv);
1230 break;
1231 case mPsiMin:
1232 x3872->setJPsiMassMin(pv);
1233 break;
1234 case mPsiMax:
1235 x3872->setJPsiMassMax(pv);
1236 break;
1237 case massMin:
1238 x3872->setMassMin(pv);
1239 break;
1240 case massMax:
1241 x3872->setMassMax(pv);
1242 break;
1243 case probMin:
1244 x3872->setProbMin(pv);
1245 break;
1246 case mFitMin:
1247 x3872->setMassFitMin(pv);
1248 break;
1249 case mFitMax:
1250 x3872->setMassFitMax(pv);
1251 break;
1252 case constrMJPsi:
1253 x3872->setConstr(pv > 0);
1254 break;
1255 case writeCandidate:
1256 writeX3872 = (pv > 0);
1257 break;
1258 default:
1259 break;
1260 }
1261 }
1262 }
1263 lX3872 = x3872->build();
1264 delete x3872;
1265 }
1266
1267 return;
1268 }
1269
1270 void BPHWriteSpecificDecay::endJob() { return; }
1271
1272 void BPHWriteSpecificDecay::setRecoParameters(const edm::ParameterSet& ps) {
1273 const string& name = ps.getParameter<string>("name");
1274 bool writeCandidate = ps.getParameter<bool>("writeCandidate");
1275 switch (rMap[name]) {
1276 case Onia:
1277 recoOnia = true;
1278 writeOnia = writeCandidate;
1279 break;
1280 case Pmm:
1281 case Psi1:
1282 case Psi2:
1283 case Ups:
1284 case Ups1:
1285 case Ups2:
1286 case Ups3:
1287 recoOnia = true;
1288 break;
1289 case Kx0:
1290 recoKx0 = true;
1291 writeKx0 = writeCandidate;
1292 break;
1293 case Pkk:
1294 recoPkk = true;
1295 writePkk = writeCandidate;
1296 break;
1297 case Bu:
1298 recoBu = true;
1299 writeBu = writeCandidate;
1300 break;
1301 case Bd:
1302 recoBd = true;
1303 writeBd = writeCandidate;
1304 break;
1305 case Bs:
1306 recoBs = true;
1307 writeBs = writeCandidate;
1308 break;
1309 case K0s:
1310 recoK0s = true;
1311 writeK0s = writeCandidate;
1312 break;
1313 case Lambda0:
1314 recoLambda0 = true;
1315 writeLambda0 = writeCandidate;
1316 break;
1317 case B0:
1318 recoB0 = true;
1319 writeB0 = writeCandidate;
1320 break;
1321 case Lambdab:
1322 recoLambdab = true;
1323 writeLambdab = writeCandidate;
1324 break;
1325 case Bc:
1326 recoBc = true;
1327 writeBc = writeCandidate;
1328 break;
1329 case X3872:
1330 recoX3872 = true;
1331 writeX3872 = writeCandidate;
1332 break;
1333 }
1334
1335 map<string, parType>::const_iterator pIter = pMap.begin();
1336 map<string, parType>::const_iterator pIend = pMap.end();
1337 while (pIter != pIend) {
1338 const map<string, parType>::value_type& entry = *pIter++;
1339 const string& pn = entry.first;
1340 parType id = entry.second;
1341 double pv = ps.getParameter<double>(pn);
1342 if (pv > -1.0e35)
1343 edm::LogVerbatim("Configuration") << "BPHWriteSpecificDecay::setRecoParameters: set " << pn << " for " << name
1344 << " : " << (parMap[rMap[name]][id] = pv);
1345 }
1346
1347 map<string, parType>::const_iterator fIter = fMap.begin();
1348 map<string, parType>::const_iterator fIend = fMap.end();
1349 while (fIter != fIend) {
1350 const map<string, parType>::value_type& entry = *fIter++;
1351 const string& fn = entry.first;
1352 parType id = entry.second;
1353 edm::LogVerbatim("Configuration") << "BPHWriteSpecificDecay::setRecoParameters: set " << fn << " for " << name
1354 << " : " << (parMap[rMap[name]][id] = (ps.getParameter<bool>(fn) ? 1 : -1));
1355 }
1356 }
1357
1358 #include "FWCore/Framework/interface/MakerMacros.h"
1359
1360 DEFINE_FWK_MODULE(BPHWriteSpecificDecay);