Reference Manual

Go to the documentation of this file.
00001 /*class BasicHepMCValidation
00002  *  
00003  *  Class to fill dqm monitor elements from existing EDM file
00004  *
00005  *  $Date: 2011/10/30 09:16:22 $
00006  *  $Revision: 1.7 $
00007  */
00008  
00009 #include "Validation/EventGenerator/interface/BasicHepMCValidation.h"
00010 
00011 #include "CLHEP/Units/defs.h"
00012 #include "CLHEP/Units/PhysicalConstants.h"
00013 
00014 using namespace edm;
00015 
00016 BasicHepMCValidation::BasicHepMCValidation(const edm::ParameterSet& iPSet): 
00017   _wmanager(iPSet),
00018   hepmcCollection_(iPSet.getParameter<edm::InputTag>("hepmcCollection"))
00019 {    
00020   dbe = 0;
00021   dbe = edm::Service<DQMStore>().operator->();
00022 }
00023 
00024 BasicHepMCValidation::~BasicHepMCValidation() {}
00025 
00026 void BasicHepMCValidation::beginJob()
00027 {
00028   if(dbe){
00030         dbe->setCurrentFolder("Generator/Particles");
00031     
00032     // Number of analyzed events
00033     nEvt = dbe->book1D("nEvt", "n analyzed Events", 1, 0., 1.);
00034         
00037         uNumber = dbe->book1D("uNumber", "No. u", 20, 0, 20);
00038         dNumber = dbe->book1D("dNumber", "No. d", 20, 0, 20);
00039         sNumber = dbe->book1D("sNumber", "No. s", 20, 0, 20);
00040     cNumber = dbe->book1D("cNumber", "No. c", 20, 0, 20);
00041         bNumber = dbe->book1D("bNumber", "No. b", 20, 0, 20);
00042         tNumber = dbe->book1D("tNumber", "No. t", 20, 0, 20);
00043         //
00044         ubarNumber = dbe->book1D("ubarNumber", "No. ubar", 20, 0, 20);
00045         dbarNumber = dbe->book1D("dbarNumber", "No. dbar", 20, 0, 20);
00046         sbarNumber = dbe->book1D("sbarNumber", "No. sbar", 20, 0, 20);
00047     cbarNumber = dbe->book1D("cbarNumber", "No. cbar", 20, 0, 20);
00048         bbarNumber = dbe->book1D("bbarNumber", "No. bbar", 20, 0, 20);
00049         tbarNumber = dbe->book1D("tbarNumber", "No. tbar", 20, 0, 20);
00050         //
00051         eminusNumber = dbe->book1D("eminusNumber", "No. e-", 20, 0, 20);
00052         nueNumber = dbe->book1D("nueNumber", "No. nu_e", 20, 0, 20);
00053         muminusNumber = dbe->book1D("muminusNumber", "No. mu-", 20, 0, 20);
00054         numuNumber = dbe->book1D("numuNumber", "No. nu_mu", 20, 0, 20);
00055         tauminusNumber = dbe->book1D("tauminusNumber", "No. tau-", 20, 0, 20);
00056         nutauNumber = dbe->book1D("nutauNumber", "No. nu_tau", 20, 0, 20);
00057         //
00058         eplusNumber = dbe->book1D("eplusNumber", "No. e+", 20, 0, 20);
00059         nuebarNumber = dbe->book1D("nuebarNumber", "No. nu_e_bar", 20, 0, 20);
00060         muplusNumber = dbe->book1D("muplusNumber", "No. mu+", 20, 0, 20);
00061         numubarNumber = dbe->book1D("numuNumber", "No. nu_mu_bar", 20, 0, 20);
00062         tauplusNumber = dbe->book1D("tauplusNumber", "No. tau+", 20, 0, 20);
00063         nutaubarNumber = dbe->book1D("nutauNumber", "No. nu_tau_bar", 20, 0, 20);
00064         //
00065         WplusNumber = dbe->book1D("WplusNumber", "No. W+", 20, 0, 20);
00066         WminusNumber = dbe->book1D("WminusNumber", "No. W-", 20, 0, 20);
00067         ZNumber = dbe->book1D("ZNumber", "No. Z", 20, 0, 20);
00068         gammaNumber = dbe->book1D("gammaNumber", "Log10(No. gamma)", 60, -1, 5); //Log
00069         gluNumber = dbe->book1D("gluonNumber", "Log10(No. gluons)", 60, -1, 5); //Log
00070         //
00071         piplusNumber = dbe->book1D("piplusNumber", "Log10(No. pi+)", 60, -1, 5); //Log
00072         piminusNumber = dbe->book1D("piminusNumber", "Log10(No. pi-)", 60, -1, 5); //Log
00073         pizeroNumber = dbe->book1D("pizeroNumber", "Log10(No. pi_0)", 60, -1, 5); //Log
00074         KplusNumber = dbe->book1D("KplusNumber", "No. K+", 100, 0, 100);
00075         KminusNumber = dbe->book1D("KminusNumber", "No. K-", 100, 0, 100);
00076         KlzeroNumber = dbe->book1D("KlzeroNumber", "No. K_l^0", 100, 0, 100);
00077         KszeroNumber = dbe->book1D("KszeroNumber", "No. K_s^0", 100, 0, 100);
00078         //
00079         pNumber = dbe->book1D("pNumber", "No. p", 100, 0, 100);
00080         pbarNumber = dbe->book1D("pbarNumber", "No. pbar", 100, 0, 100);
00081         nNumber = dbe->book1D("nNumber", "No. n", 100, 0, 100);
00082         nbarNumber = dbe->book1D("nbarNumber", "No. nbar", 100, 0, 100);
00083         l0Number = dbe->book1D("l0Number", "No. Lambda0", 100, 0, 100);
00084         l0barNumber = dbe->book1D("l0barNumber", "No. Lambda0bar", 100, 0, 100);
00085         //
00086         DplusNumber = dbe->book1D("DplusNumber", "No. D+", 20, 0, 20);
00087         DminusNumber = dbe->book1D("DminusNumber", "No. D-", 20, 0, 20);
00088         DzeroNumber = dbe->book1D("DzeroNumber", "No. D^0", 20, 0, 20);
00089         //
00090         BplusNumber = dbe->book1D("BplusNumber", "No. B+", 20, 0, 20);
00091         BminusNumber = dbe->book1D("BminusNumber", "No. B-", 20, 0, 20);
00092         BzeroNumber = dbe->book1D("BzeroNumber", "No. B^0", 20, 0, 20);
00093         BszeroNumber = dbe->book1D("BszeroNumber", "No. B^0_s", 20, 0, 20);
00094         //
00095         otherPtclNumber = dbe->book1D("otherPtclNumber", "Log10(No. other ptcls)", 60, -1, 5); //Log
00096 
00097         //Momentum 
00098         uMomentum = dbe->book1D("uMomentum", "Log10(p) u", 60, -2, 4);
00099         dMomentum = dbe->book1D("dMomentum", "Log10(p) d", 60, -2, 4);
00100         sMomentum = dbe->book1D("sMomentum", "Log10(p) s", 60, -2, 4);
00101     cMomentum = dbe->book1D("cMomentum", "Log10(p) c", 60, -2, 4);
00102         bMomentum = dbe->book1D("bMomentum", "Log10(p) b", 60, -2, 4);
00103         tMomentum = dbe->book1D("tMomentum", "Log10(p) t", 60, -2, 4);
00104         //
00105         ubarMomentum = dbe->book1D("ubarMomentum", "Log10(p) ubar", 60, -2, 4);
00106         dbarMomentum = dbe->book1D("dbarMomentum", "Log10(p) dbar", 60, -2, 4);
00107         sbarMomentum = dbe->book1D("sbarMomentum", "Log10(p) sbar", 60, -2, 4);
00108     cbarMomentum = dbe->book1D("cbarMomentum", "Log10(p) cbar", 60, -2, 4);
00109         bbarMomentum = dbe->book1D("bbarMomentum", "Log10(p) bbar", 60, -2, 4);
00110         tbarMomentum = dbe->book1D("tbarMomentum", "Log10(p) tbar", 60, -2, 4);
00111         //
00112         eminusMomentum = dbe->book1D("eminusMomentum", "Log10(p) e-", 60, -2, 4);
00113         nueMomentum = dbe->book1D("nueMomentum", "Log10(Momentum) nue", 60, -2, 4);
00114         muminusMomentum = dbe->book1D("muminusMomentum", "Log10(p) mu-", 60, -2, 4);
00115         numuMomentum = dbe->book1D("numuMomentum", "Log10(p) numu", 60, -2, 4);
00116         tauminusMomentum = dbe->book1D("tauminusMomentum", "Log10(p) tau-", 60, -2, 4);
00117         nutauMomentum = dbe->book1D("nutauMomentum", "Log10(p) nutau", 60, -2, 4);
00118         //
00119         eplusMomentum = dbe->book1D("eplusMomentum", "Log10(p) e+", 60, -2, 4);
00120         nuebarMomentum = dbe->book1D("nuebarMomentum", "Log10(p) nuebar", 60, -2, 4);
00121         muplusMomentum = dbe->book1D("muplusMomentum", "Log10(p) mu+", 60, -2, 4);
00122         numubarMomentum = dbe->book1D("numuMomentum", "Log10(p) numubar", 60, -2, 4);
00123         tauplusMomentum = dbe->book1D("tauplusMomentum", "Log10(p) tau+", 60, -2, 4);
00124         nutaubarMomentum = dbe->book1D("nutauMomentum", "Log10(p) nutaubar", 60, -2, 4);
00125         //
00126         gluMomentum = dbe->book1D("gluonMomentum", "Log10(p) gluons", 70, -3, 4);
00127         WplusMomentum = dbe->book1D("WplusMomentum", "Log10(p) W+", 60, -2, 4);
00128         WminusMomentum = dbe->book1D("WminusMomentum", "Log10(p) W-", 60, -2, 4);
00129         ZMomentum = dbe->book1D("ZMomentum", "Log10(p) Z", 60, -2, 4);
00130         gammaMomentum = dbe->book1D("gammaMomentum", "Log10(p) gamma", 70, -3, 4);
00131         //
00132         piplusMomentum = dbe->book1D("piplusMomentum", "Log10(p) pi+", 60, -2, 4);
00133         piminusMomentum = dbe->book1D("piminusMomentum", "Log10(p) pi-", 60, -2, 4);
00134         pizeroMomentum = dbe->book1D("pizeroMomentum", "Log10(p) pi_0", 60, -2, 4);
00135         KplusMomentum = dbe->book1D("KplusMomentum", "Log10(p) K+", 60, -2, 4);
00136         KminusMomentum = dbe->book1D("KminusMomentum", "Log10(p) K-", 60, -2, 4);
00137         KlzeroMomentum = dbe->book1D("KlzeroMomentum", "Log10(p) K_l^0", 60, -2, 4);
00138         KszeroMomentum = dbe->book1D("KszeroMomentum", "Log10(p) K_s^0", 60, -2, 4);
00139         //
00140         pMomentum = dbe->book1D("pMomentum", "Log10(p) p", 60, -2, 4);
00141         pbarMomentum = dbe->book1D("pbarMomentum", "Log10(p) pbar", 60, -2, 4);
00142         nMomentum = dbe->book1D("nMomentum", "Log10(p) n", 60, -2, 4);
00143         nbarMomentum = dbe->book1D("nbarMomentum", "Log10(p) nbar", 60, -2, 4);
00144         l0Momentum = dbe->book1D("l0Momentum", "Log10(p) Lambda0", 60, -2, 4);
00145         l0barMomentum = dbe->book1D("l0barMomentum", "Log10(p) Lambda0bar", 60, -2, 4);
00146         //
00147         DplusMomentum = dbe->book1D("DplusMomentum", "Log10(p) D+", 60, -2, 4);
00148         DminusMomentum = dbe->book1D("DminusMomentum", "Log10(p) D-", 60, -2, 4);
00149         DzeroMomentum = dbe->book1D("DzeroMomentum", "Log10(p) D^0", 60, -2, 4);
00150         //
00151         BplusMomentum = dbe->book1D("BplusMomentum", "Log10(p) B+", 60, -2, 4);
00152         BminusMomentum = dbe->book1D("BminusMomentum", "Log10(p) B-", 60, -2, 4);
00153         BzeroMomentum = dbe->book1D("BzeroMomentum", "Log10(p) B^0", 60, -2, 4);
00154         BszeroMomentum = dbe->book1D("BszeroMomentum", "Log10(p) B^0_s", 60, -2, 4);
00155         //
00156         otherPtclMomentum = dbe->book1D("otherPtclMomentum", "Log10(p) other ptcls", 60, -2, 4);
00157 
00159         genPtclNumber = dbe->book1D("genPtclNumber", "Log10(No. all particles)", 60, -1, 5); //Log
00160         genVrtxNumber = dbe->book1D("genVrtxNumber", "Log10(No. all vertexs)", 60, -1, 5); //Log
00161         //
00162         stablePtclNumber= dbe->book1D("stablePtclNumber", "Log10(No. stable particles)", 50, 0, 5); //Log
00163         stablePtclPhi = dbe->book1D("stablePtclPhi", "stable Ptcl Phi", 360, -180, 180);
00164         stablePtclEta = dbe->book1D("stablePtclEta", "stable Ptcl Eta (pseudo rapidity)", 220, -11, 11);
00165         stablePtclCharge = dbe->book1D("stablePtclCharge", "stablePtclCharge", 5, -2, 2);
00166         stableChaNumber= dbe->book1D("stableChaNumber", "Log10(No. stable charged particles)", 50, 0, 5); //Log
00167         stablePtclp = dbe->book1D("stablePtclp", "Log10(p) stable ptcl p", 80, -4, 4); //Log
00168         stablePtclpT = dbe->book1D("stablePtclpT", "Log10(pT) stable ptcl pT", 80, -4, 4); //Log
00169         partonNumber = dbe->book1D("partonNumber", "number of partons", 100, 0, 100);
00170         partonpT = dbe->book1D("partonpT", "Log10(pT) parton pT", 80, -4, 4); //Log
00171         outVrtxStablePtclNumber = dbe->book1D("outVrtxStablePtclNumber", "No. outgoing stable ptcls from vrtx", 10, 0, 10); 
00172         //
00173         outVrtxPtclNumber = dbe->book1D("outVrtxPtclNumber", "No. outgoing ptcls from vrtx", 30, 0, 30);
00174         vrtxZ = dbe->book1D("VrtxZ", "VrtxZ", 50 , -250, 250);
00175         vrtxRadius = dbe->book1D("vrtxRadius", "vrtxRadius", 50, 0, 50);
00176         //
00177         unknownPDTNumber = dbe->book1D("unknownPDTNumber", "Log10(No. unknown ptcls PDT)", 60, -1, 5); //Log
00178     genPtclStatus = dbe->book1D("genPtclStatus", "Status of genParticle", 200,0,200.);
00179         //
00180     Bjorken_x = dbe->book1D("Bjorken_x", "Bjorken_x", 1000, 0.0, 1.0);
00181     //
00182     status1ShortLived = dbe->book1D("status1ShortLived","Status 1 short lived", 11, 0, 11);
00183     status1ShortLived->setBinLabel(1,"d/dbar");
00184     status1ShortLived->setBinLabel(2,"u/ubar");
00185     status1ShortLived->setBinLabel(3,"s/sbar");
00186     status1ShortLived->setBinLabel(4,"c/cbar");
00187     status1ShortLived->setBinLabel(5,"b/bbar");
00188     status1ShortLived->setBinLabel(6,"t/tbar");
00189     status1ShortLived->setBinLabel(7,"g");
00190     status1ShortLived->setBinLabel(8,"tau-/tau+");
00191     status1ShortLived->setBinLabel(9,"Z0");
00192     status1ShortLived->setBinLabel(10,"W-/W+");
00193     status1ShortLived->setBinLabel(11,"PDG = 7,8,17,25-99");
00194 
00195     DeltaEcms = dbe->book1D("DeltaEcms1","deviation from nominal Ecms", 200,-1., 1.);
00196     DeltaPx = dbe->book1D("DeltaPx1","deviation from nominal Px", 200,-1., 1.);
00197     DeltaPy = dbe->book1D("DeltaPy1","deviation from nominal Py", 200,-1., 1.);
00198     DeltaPz = dbe->book1D("DeltaPz1","deviation from nominal Pz", 200,-1., 1.);
00199 
00200   }
00201   return;
00202 }
00203 
00204 void BasicHepMCValidation::endJob(){return;}
00205 void BasicHepMCValidation::beginRun(const edm::Run& iRun,const edm::EventSetup& iSetup)
00206 {
00208   iSetup.getData( fPDGTable );
00209   return;
00210 }
00211 void BasicHepMCValidation::endRun(const edm::Run& iRun,const edm::EventSetup& iSetup){return;}
00212 void BasicHepMCValidation::analyze(const edm::Event& iEvent,const edm::EventSetup& iSetup)
00213 { 
00215   int uNum = 0; int dNum = 0; int sNum = 0; int cNum = 0; int bNum = 0; int tNum = 0;
00216   int ubarNum = 0; int dbarNum = 0; int sbarNum = 0; int cbarNum = 0; int bbarNum = 0; int tbarNum = 0;
00217   int partonNum = 0;
00218   //
00219   int eminusNum = 0; int nueNum = 0; int muminusNum = 0; int numuNum = 0; int tauminusNum = 0; int nutauNum = 0;
00220   int eplusNum = 0; int nuebarNum = 0; int muplusNum = 0; int numubarNum = 0; int tauplusNum = 0; int nutaubarNum = 0;
00221   //
00222   int gluNum = 0; int WplusNum = 0; int WminusNum = 0; int ZNum = 0; int gammaNum = 0;
00223   //
00224   int piplusNum = 0; int piminusNum = 0; int pizeroNum = 0; int KplusNum = 0; int KminusNum = 0; int KlzeroNum = 0; int KszeroNum = 0;  
00225   //
00226   int pNum = 0; int pbarNum = 0; int nNum = 0; int nbarNum = 0; int l0Num = 0; int l0barNum = 0;
00227   //
00228   int DplusNum = 0; int DminusNum = 0; int DzeroNum = 0; int BplusNum = 0; int BminusNum = 0; int BzeroNum = 0; int BszeroNum = 0;
00229   //
00230   int outVrtxStablePtclNum = 0; int stablePtclNum = 0; int otherPtclNum = 0; int unknownPDTNum = 0; int stableChaNum = 0;
00231   //
00232   double bjorken = 0.;
00233   //
00234   double etotal = 0. ; double pxtotal = 0.; double pytotal = 0.; double pztotal = 0.;
00235 
00237   edm::Handle<HepMCProduct> evt;
00238   iEvent.getByLabel(hepmcCollection_, evt);
00239 
00240   //Get EVENT
00241   HepMC::GenEvent *myGenEvent = new HepMC::GenEvent(*(evt->GetEvent()));
00242 
00243   double weight = _wmanager.weight(iEvent);
00244 
00245   nEvt->Fill(0.5,weight);
00246 
00247   genPtclNumber->Fill(log10(myGenEvent->particles_size()),weight);     
00248   genVrtxNumber->Fill(log10(myGenEvent->vertices_size()),weight);
00249 
00251   HepMC::PdfInfo *pdf = myGenEvent->pdf_info();    
00252   if(pdf){
00253     bjorken = ((pdf->x1())/((pdf->x1())+(pdf->x2())));
00254   }
00255   Bjorken_x->Fill(bjorken,weight);
00256 
00257   //Looping through the VERTICES in the event
00258   HepMC::GenEvent::vertex_const_iterator vrtxBegin = myGenEvent->vertices_begin();
00259   HepMC::GenEvent::vertex_const_iterator vrtxEnd = myGenEvent->vertices_end();
00260   for(HepMC::GenEvent::vertex_const_iterator vrtxIt = vrtxBegin; vrtxIt!=vrtxEnd; ++vrtxIt)
00261     {
00263       HepMC::GenVertex *vrtx = *vrtxIt;
00264       outVrtxPtclNumber->Fill(vrtx->particles_out_size(),weight); //std::cout << "all " << vrtx->particles_out_size() << '\n';
00265       vrtxZ->Fill(vrtx->point3d().z(),weight);
00266       vrtxRadius->Fill(vrtx->point3d().perp(),weight);
00267         
00269       HepMC::GenVertex::particles_out_const_iterator vrtxPtclBegin = vrtx->particles_out_const_begin();
00270       HepMC::GenVertex::particles_out_const_iterator vrtxPtclEnd = vrtx->particles_out_const_end();
00271       outVrtxStablePtclNum = 0;
00272       for(HepMC::GenVertex::particles_out_const_iterator vrtxPtclIt = vrtxPtclBegin; vrtxPtclIt != vrtxPtclEnd; ++vrtxPtclIt)
00273         {
00274           HepMC::GenParticle *vrtxPtcl = *vrtxPtclIt;
00275           if (vrtxPtcl->status() == 1){
00276             ++outVrtxStablePtclNum; //std::cout << "stable " << outVrtxStablePtclNum << '\n';
00277           }
00278         }
00279       outVrtxStablePtclNumber->Fill(outVrtxStablePtclNum,weight);
00280     }//vertices
00281 
00282     
00284   HepMC::GenEvent::particle_const_iterator ptclBegin = myGenEvent->particles_begin();
00285   HepMC::GenEvent::particle_const_iterator ptclEnd = myGenEvent->particles_end();
00286   for(HepMC::GenEvent::particle_const_iterator ptclIt = ptclBegin; ptclIt!=ptclEnd; ++ptclIt)
00287     {
00288     
00290       HepMC::GenParticle *ptcl = *ptclIt;
00291       int Id = ptcl->pdg_id(); // std::cout << Id << '\n'; 
00292       float Log_p = log10( ptcl->momentum().rho() );
00293       double charge = 999.;     // for the charge it's needed a HepPDT method
00294       int status = ptcl->status();
00295       const HepPDT::ParticleData* PData = fPDGTable->particle(HepPDT::ParticleID(Id));
00296       if(PData==0) {
00297         //          std::cout << "Unknown id = " << Id << '\n';
00298             ++unknownPDTNum;
00299       }
00300       else
00301             charge = PData->charge();
00302 
00304       genPtclStatus->Fill((float)status,weight);
00305 
00307       if(ptcl->status() == 1){
00308             ++stablePtclNum;
00309             stablePtclPhi->Fill(ptcl->momentum().phi()/CLHEP::degree,weight); //std::cout << ptcl->polarization().phi() << '\n';
00310             stablePtclEta->Fill(ptcl->momentum().pseudoRapidity(),weight);
00311             stablePtclCharge->Fill(charge,weight); // std::cout << ptclData.charge() << '\n';
00312             stablePtclp->Fill(Log_p,weight);
00313             stablePtclpT->Fill(log10(ptcl->momentum().perp()),weight);
00314         if (charge != 0. && charge != 999.) ++stableChaNum;
00315         if ( std::abs(Id) == 1 ) status1ShortLived->Fill(1,weight);
00316         if ( std::abs(Id) == 2 ) status1ShortLived->Fill(2,weight);
00317         if ( std::abs(Id) == 3 ) status1ShortLived->Fill(3,weight);
00318         if ( std::abs(Id) == 4 ) status1ShortLived->Fill(4,weight);
00319         if ( std::abs(Id) == 5 ) status1ShortLived->Fill(5,weight);
00320         if ( std::abs(Id) == 6 ) status1ShortLived->Fill(6,weight);
00321         if ( Id == 21 ) status1ShortLived->Fill(7,weight);
00322         if ( std::abs(Id) == 15 ) status1ShortLived->Fill(8,weight);
00323         if ( Id == 23 ) status1ShortLived->Fill(9,weight);
00324         if ( std::abs(Id) == 24 ) status1ShortLived->Fill(10,weight);
00325         if ( std::abs(Id) == 7 || std::abs(Id) == 8 || std::abs(Id) == 17 || (std::abs(Id) >= 25 && std::abs(Id) <= 99) ) status1ShortLived->Fill(11,weight);
00326         etotal += ptcl->momentum().e(); 
00327         pxtotal += ptcl->momentum().px(); 
00328         pytotal += ptcl->momentum().py(); 
00329         pztotal += ptcl->momentum().pz(); 
00330       }
00331 
00332       if (abs(Id) < 6 || abs(Id) == 22){
00333         ++partonNum; partonpT->Fill(Log_p,weight);
00334       }
00335 
00337       switch(abs(Id)){
00338 
00339       case 1 : {
00340                 if(Id > 0) {
00341           ++dNum; dMomentum->Fill(Log_p,weight);}
00342                 else{
00343           ++dbarNum; dbarMomentum->Fill(Log_p,weight);}
00344       }
00345                 break;
00346                 //
00347       case 2 : {
00348                 if(Id > 0) {
00349           ++uNum; uMomentum->Fill(Log_p,weight);}
00350                 else{
00351           ++ubarNum; ubarMomentum->Fill(Log_p,weight);}
00352       }
00353                 break;
00354                 //
00355       case 3 :  {
00356                 if(Id > 0) {
00357           ++sNum; sMomentum->Fill(Log_p,weight);}
00358                 else{
00359           ++sbarNum; sbarMomentum->Fill(Log_p,weight);}
00360       }
00361                 break;
00362                 //
00363       case 4 : {
00364                 if(Id > 0) {
00365           ++cNum; cMomentum->Fill(Log_p,weight);}
00366                 else{
00367           ++cbarNum; cbarMomentum->Fill(Log_p,weight);}
00368       }
00369                 break;
00370                 //
00371       case 5 : {
00372                 if(Id > 0) {
00373           ++bNum; bMomentum->Fill(Log_p,weight);}
00374                 else{
00375           ++bbarNum; bbarMomentum->Fill(Log_p,weight);}
00376       }
00377                 break;
00378                 //
00379       case 6 : {
00380                 if(Id > 0) {
00381           ++tNum; tMomentum->Fill(Log_p,weight);}
00382                 else{
00383           ++tbarNum; tbarMomentum->Fill(Log_p,weight);}
00384       }
00385                 break;
00386                 //
00387       case 11 : {
00388                 if(Id > 0) {
00389           ++eminusNum; eminusMomentum->Fill(Log_p,weight);}
00390                 else{
00391           ++eplusNum; eplusMomentum->Fill(Log_p,weight);}
00392       }
00393                 break;
00394                 //
00395       case 12 : {
00396                 if(Id > 0) {
00397           ++nueNum; nueMomentum->Fill(Log_p, weight);}
00398                 else{
00399           ++nuebarNum; nuebarMomentum->Fill(Log_p,weight);}
00400       }
00401                 break;
00402                 //
00403       case 13 : {
00404                 if(Id > 0) {
00405           ++muminusNum; muminusMomentum->Fill(Log_p,weight);}
00406                 else{
00407           ++muplusNum; muplusMomentum->Fill(Log_p,weight);}
00408       }
00409                 break;
00410                 //
00411       case 14 : {
00412                 if(Id > 0) {
00413           ++numuNum; numuMomentum->Fill(Log_p,weight);}
00414                 else{
00415           ++numubarNum; numubarMomentum->Fill(Log_p,weight);}
00416       }
00417                 break;
00418                 //
00419       case 15 : {
00420                 if(Id > 0) {
00421           ++tauminusNum; tauminusMomentum->Fill(Log_p,weight);}
00422                 else{
00423           ++tauplusNum; tauplusMomentum->Fill(Log_p,weight);}
00424       }
00425                 break;
00426                 //
00427       case 16 : {
00428                 if(Id > 0) {
00429           ++nutauNum; nutauMomentum->Fill(Log_p,weight);}
00430                 else{
00431           ++nutaubarNum; nutaubarMomentum->Fill(Log_p,weight);}
00432       }
00433                 break;
00434                 //
00435                 //
00436       case 21 : {
00437                 ++gluNum; gluMomentum->Fill(Log_p,weight); 
00438       }
00439                 break;
00440                 //
00441       case 22 : {
00442                 ++gammaNum; gammaMomentum->Fill(Log_p,weight);
00443       }
00444                 break;
00445                 //
00446       case 23 : {
00447                 ++ZNum; ZMomentum->Fill(Log_p,weight);
00448       }
00449                 break;
00450       case 24 : {
00451                 if(Id > 0) {
00452           ++WplusNum; WplusMomentum->Fill(Log_p,weight);}
00453                 else{
00454           ++WminusNum; WminusMomentum->Fill(Log_p,weight);}
00455       }
00456                 break;
00457                 //
00458                 //
00459       case 211 : {
00460                 if(Id > 0) {
00461           ++piplusNum; piplusMomentum->Fill(Log_p,weight);}
00462                 else{
00463           ++piminusNum; piminusMomentum->Fill(Log_p,weight);}
00464       }
00465                 break;
00466                 //
00467       case 111 : {
00468                 ++pizeroNum; pizeroMomentum->Fill(Log_p,weight);
00469       }
00470                 break;
00471                 //
00472       case 321 : {
00473                 if(Id > 0) {
00474           ++KplusNum; KplusMomentum->Fill(Log_p,weight);}
00475                 else{
00476           ++KminusNum; KminusMomentum->Fill(Log_p,weight);}
00477       }
00478                 break;
00479                 //
00480       case 130 : {
00481         ++KlzeroNum; KlzeroMomentum->Fill(Log_p,weight);
00482       }
00483                 break;
00484                 //
00485       case 310 : {
00486                 ++KszeroNum; KszeroMomentum->Fill(Log_p,weight);
00487       }
00488                 break;
00489                 //
00490                 //
00491       case 2212 : {
00492                 if(Id > 0) {
00493           ++pNum; pMomentum->Fill(Log_p,weight);}
00494                 else{
00495           ++pbarNum; pbarMomentum->Fill(Log_p,weight);}
00496       }
00497                 break;
00498                 //
00499       case 2112 : {
00500                 if(Id > 0) {
00501           ++nNum; nMomentum->Fill(Log_p,weight);}
00502                 else{
00503           ++nbarNum; nbarMomentum->Fill(Log_p,weight);}
00504       }
00505                 break;
00506                 //
00507                 //
00508       case 3122 : {
00509                 if(Id > 0) {
00510           ++l0Num; l0Momentum->Fill(Log_p,weight);}
00511                 else{
00512           ++l0barNum; l0barMomentum->Fill(Log_p,weight);}
00513       }
00514         break;
00515         //
00516         //
00517       case 411 : {
00518                 if(Id > 0) {
00519           ++DplusNum; DplusMomentum->Fill(Log_p,weight);}
00520                 else{
00521           ++DminusNum; DminusMomentum->Fill(Log_p,weight);}
00522       }
00523                 break;
00524                 //
00525       case 421 : {
00526                 ++DzeroNum; DzeroMomentum->Fill(Log_p,weight);
00527       }
00528                 break;
00529                 //
00530       case 521 : {
00531                 if(Id > 0) {
00532           ++BplusNum; BplusMomentum->Fill(Log_p,weight);}
00533                 else{
00534           ++BminusNum; BminusMomentum->Fill(Log_p,weight);}
00535       }
00536                 break;
00537                 //
00538       case 511 : {
00539         ++BzeroNum; BzeroMomentum->Fill(Log_p,weight);
00540       }
00541                 break;
00542                 //
00543       case 531 : {
00544                 ++BszeroNum; BszeroMomentum->Fill(Log_p,weight);
00545       }
00546                 break;
00547                 //
00548       default : {
00549                 ++otherPtclNum; otherPtclMomentum->Fill(Log_p,weight);
00550       }
00551       }//switch
00552       //        if( 0 < Id && 100 > Id) ++part_counter[Id];
00553     }//event particles
00554 
00555 
00556   // set a default sqrt(s) and then check in the event
00557   double ecms = 7000.;
00558   if ( myGenEvent->valid_beam_particles() ) {
00559     ecms = myGenEvent->beam_particles().first->momentum().e()+myGenEvent->beam_particles().second->momentum().e();
00560   }
00561   DeltaEcms->Fill(etotal-ecms,weight);
00562   DeltaPx->Fill(pxtotal,weight);
00563   DeltaPy->Fill(pytotal,weight);
00564   DeltaPz->Fill(pztotal,weight);
00565 
00566  
00568   stablePtclNumber->Fill(log10(stablePtclNum+0.1),weight); 
00569   stableChaNumber->Fill(log10(stableChaNum+0.1),weight); 
00570   otherPtclNumber->Fill(log10(otherPtclNum+0.1),weight);
00571   unknownPDTNumber->Fill(log10(unknownPDTNum+0.1),weight);
00572   //
00573   dNumber->Fill(dNum,weight); uNumber->Fill(uNum,weight); sNumber->Fill(sNum,weight); cNumber->Fill(cNum,weight); bNumber->Fill(bNum,weight); tNumber->Fill(tNum,weight);  
00574   dbarNumber->Fill(dbarNum,weight); ubarNumber->Fill(ubarNum,weight); sbarNumber->Fill(sbarNum,weight); cbarNumber->Fill(cbarNum,weight); bbarNumber->Fill(bbarNum,weight); tbarNumber->Fill(tbarNum,weight); 
00575   partonNumber->Fill(partonNum,weight);
00576   //
00577   eminusNumber->Fill(eminusNum,weight); nueNumber->Fill(nueNum,weight); muminusNumber->Fill(muminusNum,weight); numuNumber->Fill(numuNum,weight); tauminusNumber->Fill(tauminusNum,weight); nutauNumber->Fill(nutauNum,weight);  
00578   eplusNumber->Fill(eplusNum,weight); nuebarNumber->Fill(nuebarNum,weight); muplusNumber->Fill(muplusNum,weight); numubarNumber->Fill(numubarNum,weight); tauplusNumber->Fill(tauplusNum,weight); nutaubarNumber->Fill(nutaubarNum,weight);  
00579   //
00580   ZNumber->Fill(ZNum,weight); WminusNumber->Fill(WminusNum,weight); WplusNumber->Fill(WplusNum,weight); 
00581   gammaNumber->Fill(log10(gammaNum+0.1),weight);
00582   gluNumber->Fill(log10(gluNum+0.1),weight);
00583   //
00584   piplusNumber->Fill(log10(piplusNum+0.1),weight);
00585   piminusNumber->Fill(log10(piminusNum+0.1),weight);
00586   pizeroNumber->Fill(log10(pizeroNum+0.1),weight);
00587   KplusNumber->Fill(KplusNum,weight); KminusNumber->Fill(KminusNum,weight); KlzeroNumber->Fill(KlzeroNum,weight); KszeroNumber->Fill(KszeroNum,weight); 
00588   //
00589   pNumber->Fill(pNum,weight); pbarNumber->Fill(pbarNum,weight); nNumber->Fill(nNum,weight); nbarNumber->Fill(nbarNum,weight); l0Number->Fill(l0Num); l0barNumber->Fill(l0barNum,weight);    
00590   //
00591   DplusNumber->Fill(DplusNum,weight); DminusNumber->Fill(DminusNum,weight); DzeroNumber->Fill(DzeroNum,weight); BplusNumber->Fill(BplusNum,weight); BminusNumber->Fill(BminusNum,weight); BzeroNumber->Fill(BzeroNum,weight); BszeroNumber->Fill(BszeroNum,weight); 
00592 
00593   delete myGenEvent;
00594 }//analyze
/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/Validation/EventGenerator/plugins/BasicHepMCValidation.cc
