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
00010
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
00031 switch (e.size() - 1 - last) {
00032 case 0:
00033
00034 pTpowheg = -1.;
00035 pTveto = infoPtr->QFac();
00036 noRad = true;
00037
00038
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
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
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
00064 pTshower = -1.;
00065
00066
00067 return false;
00068 }
00069
00070
00071
00072 bool EmissionVetoHook::doVetoISREmission(int, const Pythia8::Event &e) {
00073
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
00084 if (e[i].pT() > pTveto) {
00085 nISRveto++;
00086 if(Verbosity) cout << "ISR vetoed" << endl;
00087 return true;
00088 }
00089
00090 if (pTshower < 0.) pTshower = e[i].pT();
00091
00092 return false;
00093 }
00094
00095 bool EmissionVetoHook::doVetoFSREmission(int, const Pythia8::Event &e) {
00096
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
00108 if (e[i].pT() > pTveto) {
00109 nFSRveto++;
00110 if(Verbosity) cout << "FSR vetoed" << endl;
00111 return true;
00112 }
00113
00114 if (pTshower < 0.) pTshower = e[i].pT();
00115
00116 return false;
00117 }