File indexing completed on 2023-03-17 11:27:17
0001
0002
0003
0004
0005
0006
0007 #include "Validation/EventGenerator/interface/BasicHepMCValidation.h"
0008
0009 #include "CLHEP/Units/defs.h"
0010 #include "CLHEP/Units/PhysicalConstants.h"
0011 #include "Validation/EventGenerator/interface/DQMHelper.h"
0012 using namespace edm;
0013
0014 BasicHepMCValidation::BasicHepMCValidation(const edm::ParameterSet &iPSet)
0015 : wmanager_(iPSet, consumesCollector()), hepmcCollection_(iPSet.getParameter<edm::InputTag>("hepmcCollection")) {
0016 hepmcCollectionToken_ = consumes<HepMCProduct>(hepmcCollection_);
0017 fPDGTableToken = esConsumes<edm::Transition::BeginRun>();
0018 }
0019
0020 BasicHepMCValidation::~BasicHepMCValidation() {}
0021
0022 void BasicHepMCValidation::dqmBeginRun(const edm::Run &r, const edm::EventSetup &c) {
0023 fPDGTable = c.getHandle(fPDGTableToken);
0024 }
0025
0026 namespace {
0027
0028 constexpr double logPdfMax = 3.0;
0029 constexpr double logPdfMin = -3.0;
0030 constexpr int logPdfNbin = 150;
0031 constexpr double logPdfBinsize = (logPdfMax - logPdfMin) / (double)logPdfNbin;
0032 constexpr double logQScaleMax = 4.0;
0033 constexpr double logQScaleMin = -1.0;
0034 constexpr int logQScaleNbin = 500;
0035 constexpr double logQScaleBinsize = (logQScaleMax - logQScaleMin) / (double)logQScaleNbin;
0036 }
0037
0038 void BasicHepMCValidation::bookHistograms(DQMStore::IBooker &i, edm::Run const &, edm::EventSetup const &) {
0039
0040 DQMHelper dqm(&i);
0041 i.setCurrentFolder("Generator/Particles");
0042
0043
0044 nEvt = dqm.book1dHisto("nEvt", "n analyzed Events", 1, 0., 1.);
0045
0046
0047
0048
0049 particles.push_back(ParticleMonitor("u", 1, i));
0050 particles.push_back(ParticleMonitor("ubar", -1, i));
0051 particles.push_back(ParticleMonitor("d", 2, i));
0052 particles.push_back(ParticleMonitor("dbar", -2, i));
0053 particles.push_back(ParticleMonitor("s", 3, i));
0054 particles.push_back(ParticleMonitor("sbar", -3, i));
0055 particles.push_back(ParticleMonitor("c", 4, i));
0056 particles.push_back(ParticleMonitor("cbar", -4, i));
0057 particles.push_back(ParticleMonitor("b", 5, i));
0058 particles.push_back(ParticleMonitor("bbar", -5, i));
0059 particles.push_back(ParticleMonitor("t", 6, i));
0060 particles.push_back(ParticleMonitor("tbar", -6, i));
0061
0062
0063 particles.push_back(ParticleMonitor("eminus", 11, i));
0064 particles.push_back(ParticleMonitor("eplus", -11, i));
0065 particles.push_back(ParticleMonitor("nue", 12, i));
0066 particles.push_back(ParticleMonitor("nuebar", -12, i));
0067 particles.push_back(ParticleMonitor("muminus", 13, i));
0068 particles.push_back(ParticleMonitor("muplus", -13, i));
0069 particles.push_back(ParticleMonitor("numu", 14, i));
0070 particles.push_back(ParticleMonitor("numubar", -14, i));
0071 particles.push_back(ParticleMonitor("tauminus", 15, i));
0072 particles.push_back(ParticleMonitor("tauplus", -15, i));
0073 particles.push_back(ParticleMonitor("nutau", 16, i));
0074 particles.push_back(ParticleMonitor("nutaubar", -16, i));
0075
0076
0077 particles.push_back(ParticleMonitor("Wplus", 24, i));
0078 particles.push_back(ParticleMonitor("Wminus", -24, i));
0079 particles.push_back(ParticleMonitor("Z", 23, i));
0080 particles.push_back(ParticleMonitor("gamma", 22, i));
0081 particles.push_back(ParticleMonitor("gluon", 21, i));
0082
0083
0084 particles.push_back(ParticleMonitor("piplus", 211, i, true));
0085 particles.push_back(ParticleMonitor("piminus", -211, i, true));
0086 particles.push_back(ParticleMonitor("pizero", 111, i, true));
0087 particles.push_back(ParticleMonitor("Kplus", 321, i));
0088 particles.push_back(ParticleMonitor("Kminus", -321, i));
0089 particles.push_back(ParticleMonitor("Klzero", 130, i));
0090 particles.push_back(ParticleMonitor("Kszero", 310, i));
0091
0092
0093 particles.push_back(ParticleMonitor("p", 2212, i, true));
0094 particles.push_back(ParticleMonitor("pbar", -2212, i, true));
0095 particles.push_back(ParticleMonitor("n", 2112, i, true));
0096 particles.push_back(ParticleMonitor("nbar", -2112, i, true));
0097 particles.push_back(ParticleMonitor("lambda0", 3122, i));
0098 particles.push_back(ParticleMonitor("lambda0bar", -3122, i));
0099
0100
0101 particles.push_back(ParticleMonitor("Dplus", 411, i));
0102 particles.push_back(ParticleMonitor("Dminus", -411, i));
0103 particles.push_back(ParticleMonitor("Dzero", 421, i));
0104 particles.push_back(ParticleMonitor("Dzerobar", -421, i));
0105
0106
0107 particles.push_back(ParticleMonitor("Bplus", 521, i));
0108 particles.push_back(ParticleMonitor("Bminus", -521, i));
0109 particles.push_back(ParticleMonitor("Bzero", 511, i));
0110 particles.push_back(ParticleMonitor("Bzerobar", -511, i));
0111 particles.push_back(ParticleMonitor("Bszero", 531, i));
0112 particles.push_back(ParticleMonitor("Bszerobar", -531, i));
0113
0114
0115 otherPtclNumber = dqm.book1dHisto(
0116 "otherPtclNumber", "Log10(No. other ptcls)", 60, -1, 5, "log_{10}(No. other ptcls)", "Number of Events");
0117 otherPtclMomentum = dqm.book1dHisto("otherPtclMomentum",
0118 "Log10(p) other ptcls",
0119 60,
0120 -2,
0121 4,
0122 "log10(P^{other ptcls}) (log_{10}(GeV))",
0123 "Number of Events");
0124
0125
0126 genPtclNumber = dqm.book1dHisto(
0127 "genPtclNumber", "Log10(No. all particles)", 60, -1, 5, "log10(No. all particles)", "Number of Events");
0128 genVrtxNumber = dqm.book1dHisto(
0129 "genVrtxNumber", "Log10(No. all vertexs)", 60, -1, 5, "log10(No. all vertexs)", "Number of Events");
0130
0131 stablePtclNumber = dqm.book1dHisto("stablePtclNumber",
0132 "Log10(No. stable particles)",
0133 50,
0134 0,
0135 5,
0136 "log10(No. stable particles)",
0137 "Number of Events");
0138 stablePtclPhi = dqm.book1dHisto(
0139 "stablePtclPhi", "stable Ptcl Phi", 360, -180, 180, "#phi^{stable Ptcl} (rad)", "Number of Events");
0140 stablePtclEta = dqm.book1dHisto(
0141 "stablePtclEta", "stable Ptcl Eta (pseudo rapidity)", 220, -11, 11, "#eta^{stable Ptcl}", "Number of Events");
0142 stablePtclCharge =
0143 dqm.book1dHisto("stablePtclCharge", "stablePtclCharge", 5, -2, 2, "Charge^{stable ptcls}", "Number of Events");
0144 stableChaNumber = dqm.book1dHisto("stableChaNumber",
0145 "Log10(No. stable charged particles)",
0146 50,
0147 0,
0148 5,
0149 "log_{10}(No. stable charged particles)",
0150 "Number of Events");
0151 stablePtclp = dqm.book1dHisto("stablePtclp",
0152 "Log10(p) stable ptcl p",
0153 80,
0154 -4,
0155 4,
0156 "log_{10}(P^{stable ptcl}) (log_{10}(GeV))",
0157 "Number of Events");
0158 stablePtclpT = dqm.book1dHisto("stablePtclpT",
0159 "Log10(pT) stable ptcl pT",
0160 80,
0161 -4,
0162 4,
0163 "log_{10}(P_{t}^{stable ptcl}) (log_{10}(GeV))",
0164 "Number of Events");
0165 partonNumber =
0166 dqm.book1dHisto("partonNumber", "number of partons", 100, 0, 100, "number of partons", "Number of Events");
0167 partonpT =
0168 dqm.book1dHisto("partonpT", "Log10(pT) parton pT", 80, -4, 4, "Log10(P_{t}^{parton})", "Number of Events");
0169 outVrtxStablePtclNumber = dqm.book1dHisto("outVrtxStablePtclNumber",
0170 "No. outgoing stable ptcls from vrtx",
0171 10,
0172 0,
0173 10,
0174 "No. outgoing stable ptcls from vrtx",
0175 "Number of Events");
0176
0177 outVrtxPtclNumber = dqm.book1dHisto("outVrtxPtclNumber",
0178 "No. outgoing ptcls from vrtx",
0179 30,
0180 0,
0181 30,
0182 "No. outgoing ptcls from vrtx",
0183 "Number of Events");
0184 vrtxZ = dqm.book1dHisto("VrtxZ", "VrtxZ", 50, -250, 250, "Z_{Vtx}", "Number of Events");
0185 vrtxRadius = dqm.book1dHisto("vrtxRadius", "vrtxRadius", 50, 0, 50, "R_{vtx}", "Number of Events");
0186
0187 unknownPDTNumber = dqm.book1dHisto("unknownPDTNumber",
0188 "Log10(No. unknown ptcls PDT)",
0189 60,
0190 -1,
0191 5,
0192 "log_{10}(No. unknown ptcls PDT)",
0193 "Number of Events");
0194 genPtclStatus = dqm.book1dHisto("genPtclStatus", "Status of genParticle", 200, 0, 200., "", "Number of Events");
0195
0196 Bjorken_x = dqm.book1dHisto("Bjorken_x", "Bjorken_x", 1000, 0.0, 1.0, "Bjorken_{x}", "Number of Events");
0197 pdf_u = dqm.book1dHisto(
0198 "pdf_u", "Log10(PDF(u,x,Q))", logPdfNbin, logPdfMin, logPdfMax, "log_{10}(x*f(x))", "Number of Events");
0199 pdf_ubar = dqm.book1dHisto("pdf_ubar",
0200 "Log10(PDF(ubar,x,Q))",
0201 logPdfNbin,
0202 logPdfMin,
0203 logPdfMax,
0204 "log_{10}(x*f(x))",
0205 "Number of Events");
0206 pdf_d = dqm.book1dHisto(
0207 "pdf_d", "Log10(PDF(d,x,Q))", logPdfNbin, logPdfMin, logPdfMax, "log_{10}(x*f(x))", "Number of Events");
0208 pdf_dbar = dqm.book1dHisto("pdf_dbar",
0209 "Log10(PDF(dbar,x,Q))",
0210 logPdfNbin,
0211 logPdfMin,
0212 logPdfMax,
0213 "log_{10}(x*f(x))",
0214 "Number of Events");
0215 pdf_ssbar = dqm.book1dHisto("pdf_ssbar",
0216 "Log10(PDF(ssbar,x,Q))",
0217 logPdfNbin,
0218 logPdfMin,
0219 logPdfMax,
0220 "log_{10}(x*f(x))",
0221 "Number of Events");
0222 pdf_ccbar = dqm.book1dHisto("pdf_ccbar",
0223 "Log10(PDF(ccbar,x,Q))",
0224 logPdfNbin,
0225 logPdfMin,
0226 logPdfMax,
0227 "log_{10}(x*f(x))",
0228 "Number of Events");
0229 pdf_bbbar = dqm.book1dHisto("pdf_bbbar",
0230 "Log10(PDF(bbbar,x,Q))",
0231 logPdfNbin,
0232 logPdfMin,
0233 logPdfMax,
0234 "log_{10}(x*f(x))",
0235 "Number of Events");
0236 pdf_g = dqm.book1dHisto(
0237 "pdf_g", "Log10(PDF(g,x,Q))", logPdfNbin, logPdfMin, logPdfMax, "log_{10}(x*f(x))", "Number of Events");
0238 scalePDF = dqm.book1dHisto("scalePDF",
0239 "Log10(Q-scale(GeV))",
0240 logQScaleNbin,
0241 logQScaleMin,
0242 logQScaleMax,
0243 "log_{10}(Q-scale(GeV))",
0244 "Number of Events");
0245 parton1Id = dqm.book1dHisto("parton1Id", "ID of parton 1", 45, -14.5, 30.5, "ID", "Number of Events");
0246 parton2Id = dqm.book1dHisto("parton2Id", "ID of parton 2", 45, -14.5, 30.5, "ID", "Number of Events");
0247
0248 status1ShortLived = dqm.book1dHisto("status1ShortLived", "Status 1 short lived", 11, 0, 11, "", "Number of Events");
0249 status1ShortLived->setBinLabel(1, "d/dbar");
0250 status1ShortLived->setBinLabel(2, "u/ubar");
0251 status1ShortLived->setBinLabel(3, "s/sbar");
0252 status1ShortLived->setBinLabel(4, "c/cbar");
0253 status1ShortLived->setBinLabel(5, "b/bbar");
0254 status1ShortLived->setBinLabel(6, "t/tbar");
0255 status1ShortLived->setBinLabel(7, "g");
0256 status1ShortLived->setBinLabel(8, "tau-/tau+");
0257 status1ShortLived->setBinLabel(9, "Z0");
0258 status1ShortLived->setBinLabel(10, "W-/W+");
0259 status1ShortLived->setBinLabel(11, "PDG = 7,8,17,25-99");
0260
0261 log10DeltaEcms = dqm.book1dHisto("DeltaEcms1log10",
0262 "log_{10} of deviation from nominal Ecms",
0263 200,
0264 -1.,
0265 5.,
0266 "log_{10}(#DeltaE) (log_{10}(GeV))",
0267 "Number of Events");
0268 DeltaEcms =
0269 dqm.book1dHisto("DeltaEcms1", "deviation from nominal Ecms", 200, -1., 1., "#DeltaE (GeV)", "Number of Events");
0270 DeltaPx =
0271 dqm.book1dHisto("DeltaPx1", "deviation from nominal Px", 200, -1., 1., "#DeltaP_{x} (GeV)", "Number of Events");
0272 DeltaPy =
0273 dqm.book1dHisto("DeltaPy1", "deviation from nominal Py", 200, -1., 1., "#DeltaP_{y} (GeV)", "Number of Events");
0274 DeltaPz =
0275 dqm.book1dHisto("DeltaPz1", "deviation from nominal Pz", 200, -1., 1., "#DeltaP_{z} (GeV)", "Number of Events");
0276 return;
0277 }
0278
0279 void BasicHepMCValidation::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0280
0281 int partonNum = 0;
0282
0283 int outVrtxStablePtclNum = 0;
0284 int stablePtclNum = 0;
0285 int otherPtclNum = 0;
0286 int unknownPDTNum = 0;
0287 int stableChaNum = 0;
0288
0289 double bjorken = 0.;
0290 double logPdf1 = 0.;
0291 double logPdf2 = 0.;
0292 double logQScale = 0.;
0293
0294 double etotal = 0.;
0295 double pxtotal = 0.;
0296 double pytotal = 0.;
0297 double pztotal = 0.;
0298
0299
0300 edm::Handle<HepMCProduct> evt;
0301 iEvent.getByToken(hepmcCollectionToken_, evt);
0302
0303
0304 HepMC::GenEvent const *myGenEvent = evt->GetEvent();
0305
0306 double weight = wmanager_.weight(iEvent);
0307
0308 nEvt->Fill(0.5, weight);
0309
0310 genPtclNumber->Fill(log10(myGenEvent->particles_size()), weight);
0311 genVrtxNumber->Fill(log10(myGenEvent->vertices_size()), weight);
0312
0313
0314 HepMC::PdfInfo const *pdf = myGenEvent->pdf_info();
0315 if (pdf) {
0316 bjorken = ((pdf->x1()) / ((pdf->x1()) + (pdf->x2())));
0317 logQScale = log10(pdf->scalePDF());
0318 if (logQScale > logQScaleMax)
0319 logQScale = logQScaleMax - 0.5 * logQScaleBinsize;
0320 if (logQScale < logQScaleMin)
0321 logQScale = logQScaleMin + 0.5 * logQScaleBinsize;
0322 logPdf1 = log10(pdf->pdf1());
0323 if (logPdf1 > logPdfMax)
0324 logPdf1 = logPdfMax - 0.5 * logPdfBinsize;
0325 if (logPdf1 < logPdfMin)
0326 logPdf1 = logPdfMin + 0.5 * logPdfBinsize;
0327 logPdf2 = log10(pdf->pdf2());
0328 if (logPdf2 > logPdfMax)
0329 logPdf2 = logPdfMax - 0.5 * logPdfBinsize;
0330 if (logPdf2 < logPdfMin)
0331 logPdf2 = logPdfMin + 0.5 * logPdfBinsize;
0332 Bjorken_x->Fill(bjorken, weight);
0333 scalePDF->Fill(logQScale, weight);
0334 parton1Id->Fill((double)pdf->id1(), weight);
0335 parton2Id->Fill((double)pdf->id2(), weight);
0336 if (pdf->id1() == 2)
0337 pdf_u->Fill(logPdf1, weight);
0338 if (pdf->id2() == 2)
0339 pdf_u->Fill(logPdf2, weight);
0340 if (pdf->id1() == -2)
0341 pdf_ubar->Fill(logPdf1, weight);
0342 if (pdf->id2() == -2)
0343 pdf_ubar->Fill(logPdf2, weight);
0344 if (pdf->id1() == 1)
0345 pdf_d->Fill(logPdf1, weight);
0346 if (pdf->id2() == 1)
0347 pdf_d->Fill(logPdf2, weight);
0348 if (pdf->id1() == -1)
0349 pdf_dbar->Fill(logPdf1, weight);
0350 if (pdf->id2() == -1)
0351 pdf_dbar->Fill(logPdf2, weight);
0352 if (std::abs(pdf->id1()) == 3)
0353 pdf_ssbar->Fill(logPdf1, weight);
0354 if (std::abs(pdf->id2()) == 3)
0355 pdf_ssbar->Fill(logPdf2, weight);
0356 if (std::abs(pdf->id1()) == 4)
0357 pdf_ccbar->Fill(logPdf1, weight);
0358 if (std::abs(pdf->id2()) == 4)
0359 pdf_ccbar->Fill(logPdf2, weight);
0360 if (std::abs(pdf->id1()) == 5)
0361 pdf_bbbar->Fill(logPdf1, weight);
0362 if (std::abs(pdf->id2()) == 5)
0363 pdf_bbbar->Fill(logPdf2, weight);
0364 if (std::abs(pdf->id1()) == 21)
0365 pdf_g->Fill(logPdf1, weight);
0366 if (std::abs(pdf->id2()) == 21)
0367 pdf_g->Fill(logPdf2, weight);
0368 }
0369
0370
0371 HepMC::GenEvent::vertex_const_iterator vrtxBegin = myGenEvent->vertices_begin();
0372 HepMC::GenEvent::vertex_const_iterator vrtxEnd = myGenEvent->vertices_end();
0373 unsigned int nvtx(0);
0374 for (HepMC::GenEvent::vertex_const_iterator vrtxIt = vrtxBegin; vrtxIt != vrtxEnd; ++vrtxIt) {
0375
0376 HepMC::GenVertex const *vrtx = *vrtxIt;
0377 outVrtxPtclNumber->Fill(vrtx->particles_out_size(), weight);
0378
0379
0380 if (nvtx == 0) {
0381 vrtxZ->Fill(vrtx->point3d().z(), weight);
0382 vrtxRadius->Fill(vrtx->point3d().perp(), weight);
0383 }
0384
0385 HepMC::GenVertex::particles_out_const_iterator vrtxPtclBegin = vrtx->particles_out_const_begin();
0386 HepMC::GenVertex::particles_out_const_iterator vrtxPtclEnd = vrtx->particles_out_const_end();
0387 outVrtxStablePtclNum = 0;
0388 for (HepMC::GenVertex::particles_out_const_iterator vrtxPtclIt = vrtxPtclBegin; vrtxPtclIt != vrtxPtclEnd;
0389 ++vrtxPtclIt) {
0390 HepMC::GenParticle const *vrtxPtcl = *vrtxPtclIt;
0391 if (vrtxPtcl->status() == 1) {
0392 ++outVrtxStablePtclNum;
0393
0394 }
0395 }
0396 outVrtxStablePtclNumber->Fill(outVrtxStablePtclNum, weight);
0397 nvtx++;
0398 }
0399
0400
0401 HepMC::GenEvent::particle_const_iterator ptclBegin = myGenEvent->particles_begin();
0402 HepMC::GenEvent::particle_const_iterator ptclEnd = myGenEvent->particles_end();
0403 for (HepMC::GenEvent::particle_const_iterator ptclIt = ptclBegin; ptclIt != ptclEnd; ++ptclIt) {
0404
0405 HepMC::GenParticle const *ptcl = *ptclIt;
0406 int Id = ptcl->pdg_id();
0407 float Log_p = log10(ptcl->momentum().rho());
0408 double charge = 999.;
0409 int status = ptcl->status();
0410 const HepPDT::ParticleData *PData = fPDGTable->particle(HepPDT::ParticleID(Id));
0411 if (PData == nullptr) {
0412
0413 ++unknownPDTNum;
0414 } else
0415 charge = PData->charge();
0416
0417
0418 genPtclStatus->Fill((float)status, weight);
0419
0420
0421 if (ptcl->status() == 1) {
0422 ++stablePtclNum;
0423 stablePtclPhi->Fill(ptcl->momentum().phi() / CLHEP::degree,
0424 weight);
0425 stablePtclEta->Fill(ptcl->momentum().pseudoRapidity(), weight);
0426 stablePtclCharge->Fill(charge, weight);
0427 stablePtclp->Fill(Log_p, weight);
0428 stablePtclpT->Fill(log10(ptcl->momentum().perp()), weight);
0429 if (charge != 0. && charge != 999.)
0430 ++stableChaNum;
0431 if (std::abs(Id) == 1)
0432 status1ShortLived->Fill(1, weight);
0433 if (std::abs(Id) == 2)
0434 status1ShortLived->Fill(2, weight);
0435 if (std::abs(Id) == 3)
0436 status1ShortLived->Fill(3, weight);
0437 if (std::abs(Id) == 4)
0438 status1ShortLived->Fill(4, weight);
0439 if (std::abs(Id) == 5)
0440 status1ShortLived->Fill(5, weight);
0441 if (std::abs(Id) == 6)
0442 status1ShortLived->Fill(6, weight);
0443 if (Id == 21)
0444 status1ShortLived->Fill(7, weight);
0445 if (std::abs(Id) == 15)
0446 status1ShortLived->Fill(8, weight);
0447 if (Id == 23)
0448 status1ShortLived->Fill(9, weight);
0449 if (std::abs(Id) == 24)
0450 status1ShortLived->Fill(10, weight);
0451 if (std::abs(Id) == 7 || std::abs(Id) == 8 || std::abs(Id) == 17 || (std::abs(Id) >= 25 && std::abs(Id) <= 99))
0452 status1ShortLived->Fill(11, weight);
0453 etotal += ptcl->momentum().e();
0454 pxtotal += ptcl->momentum().px();
0455 pytotal += ptcl->momentum().py();
0456 pztotal += ptcl->momentum().pz();
0457 }
0458
0459 if (abs(Id) < 6 || abs(Id) == 22) {
0460 ++partonNum;
0461 partonpT->Fill(Log_p, weight);
0462 }
0463
0464 bool indentified = false;
0465 for (unsigned int i = 0; i < particles.size(); i++) {
0466 if (particles.at(i).Fill(ptcl, weight)) {
0467 indentified = true;
0468 break;
0469 }
0470 }
0471 if (!indentified) {
0472 ++otherPtclNum;
0473 otherPtclMomentum->Fill(Log_p, weight);
0474 }
0475 }
0476
0477
0478 double ecms = 13000.;
0479 if (myGenEvent->valid_beam_particles()) {
0480 ecms = myGenEvent->beam_particles().first->momentum().e() + myGenEvent->beam_particles().second->momentum().e();
0481 }
0482 log10DeltaEcms->Fill(log10(fabs(etotal - ecms)), weight);
0483 DeltaEcms->Fill(etotal - ecms, weight);
0484 DeltaPx->Fill(pxtotal, weight);
0485 DeltaPy->Fill(pytotal, weight);
0486 DeltaPz->Fill(pztotal, weight);
0487
0488
0489 stablePtclNumber->Fill(log10(stablePtclNum + 0.1), weight);
0490 stableChaNumber->Fill(log10(stableChaNum + 0.1), weight);
0491 otherPtclNumber->Fill(log10(otherPtclNum + 0.1), weight);
0492 unknownPDTNumber->Fill(log10(unknownPDTNum + 0.1), weight);
0493
0494 partonNumber->Fill(partonNum, weight);
0495 for (unsigned int i = 0; i < particles.size(); i++) {
0496 particles.at(i).FillCount(weight);
0497 };
0498
0499 }