CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_8_patch3/src/GeneratorInterface/Pythia8Interface/src/EmissionVetoHook.cc

Go to the documentation of this file.
00001 #include "GeneratorInterface/Pythia8Interface/interface/EmissionVetoHook.h"
00002 
00003 #include "FWCore/ServiceRegistry/interface/Service.h"
00004 void EmissionVetoHook::fatalEmissionVeto(string message) {
00005   throw edm::Exception(edm::errors::Configuration,"Pythia8Interface")
00006     << "EmissionVeto: " << message << endl;
00007 }
00008 
00009 // Use VetoMIStep to analyse the incoming LHEF event and
00010 // extract the veto scale
00011 bool EmissionVetoHook::doVetoMIStep(int, const Pythia8::Event &e) {
00012   int first=-1, myid;
00013   last = -1;
00014   for(int ip = 2; ip < e.size(); ip++) {
00015     myid = e[ip].id();
00016     if(abs(myid) < 6 || abs(myid) == 21) continue;
00017     first = ip;
00018     break;
00019   }
00020   if(first < 0) fatalEmissionVeto(string("signal particles not found"));
00021   for(int ip = first; ip < e.size(); ip++) {
00022     myid = e[ip].id();
00023     if(abs(myid) < 6 || abs(myid) == 21) continue;
00024     last = ip;
00025   }
00026   if(Verbosity)
00027     cout << "last before powheg emission = " << last << " , id = "
00028          << e[last].id() << endl;
00029 
00030   // Some events may not have radiation from POWHEG
00031   switch (e.size() - 1 - last) {
00032   case 0:
00033     // No radiation - veto scale is given by the factorisation scale
00034     pTpowheg = -1.;
00035     pTveto   = infoPtr->QFac();
00036     noRad    = true;
00037 
00038     // If this is the first no radiation event, then print scale
00039     if (firstNoRad) {
00040       if(Verbosity)
00041         cout << "Info: no POWHEG radiation, Q = " << pTveto
00042              << " GeV" << endl;
00043       firstNoRad = false;
00044     }
00045     break;
00046 
00047   case 1:
00048     // Radiation is parton last+1 - first check that it is as expected
00049     if (e[last+1].id() != 21 && e[last+1].idAbs() > 5) {
00050       cout << endl << "Emergency dump of the intermediate event: " << endl;
00051       e.list();
00052       fatalEmissionVeto(string("Error: jet is not quark/gluon"));
00053     }
00054     // Veto scale is given by jet pT
00055     pTveto = pTpowheg = e[last+1].pT();
00056     noRad  = false;
00057     noRad  = false;
00058     break;
00059   }
00060 
00061   if(Verbosity) cout << "veto pT = " << pTveto << endl;
00062 
00063   // Initialise other variables
00064   pTshower = -1.;
00065 
00066   // Do not veto the event
00067   return false;
00068 }
00069 
00070 // For subsequent ISR/FSR emissions, find the pT of the shower
00071 // emission and veto as necessary
00072 bool EmissionVetoHook::doVetoISREmission(int, const Pythia8::Event &e) {
00073   // ISR - next shower emission is given status 43
00074   int i;
00075   for (i = e.size() - 1; i > last; i--)
00076     if (e[i].isFinal() && e[i].status() == 43) break;
00077   if (i == last) {
00078     cout << endl << "Emergency dump of the intermediate event: " << endl;
00079     e.list();
00080     fatalEmissionVeto(string("Error: couldn't find ISR emission"));
00081   }
00082      
00083   // Veto if above the POWHEG scale
00084   if (e[i].pT() > pTveto) {
00085     nISRveto++;
00086     if(Verbosity) cout << "ISR vetoed" << endl;
00087     return true; 
00088   }
00089   // Store the first shower emission pT
00090   if (pTshower < 0.) pTshower = e[i].pT();
00091    
00092   return false;
00093 }
00094 
00095 bool EmissionVetoHook::doVetoFSREmission(int, const Pythia8::Event &e) {
00096   // FSR - shower emission will have status 51 and not be t/tbar
00097   int i;
00098   for (i = e.size() - 1; i > last; i--)
00099     if (e[i].isFinal() && e[i].status() == 51 &&
00100         e[i].idAbs() != 6) break;
00101   if (i == last) {
00102     cout << endl << "Emergency dump of the intermediate event: " << endl;
00103     e.list();
00104     fatalEmissionVeto(string("Error: couldn't find FSR emission"));
00105   }
00106 
00107   // Veto if above the POWHEG scale
00108   if (e[i].pT() > pTveto) {
00109     nFSRveto++;  
00110     if(Verbosity) cout << "FSR vetoed" << endl;
00111     return true;
00112   }
00113   // Store the first shower emission pT   
00114   if (pTshower < 0.) pTshower = e[i].pT();
00115 
00116   return false;
00117 }