|
|
|
#include <EmissionVetoHook1.h>
Public Member Functions | |
| bool | canVetoFSREmission () |
| bool | canVetoISREmission () |
| bool | canVetoMPIEmission () |
| bool | canVetoMPIStep () |
| bool | doVetoFSREmission (int, const Pythia8::Event &e, int iSys, bool) |
| bool | doVetoISREmission (int, const Pythia8::Event &e, int iSys) |
| bool | doVetoMPIEmission (int, const Pythia8::Event &e) |
| bool | doVetoMPIStep (int nMPI, const Pythia8::Event &e) |
| EmissionVetoHook1 (int nFinalIn, bool vetoOnIn, int vetoCountIn, int pThardModeIn, int pTemtModeIn, int emittedModeIn, int pTdefModeIn, bool MPIvetoOnIn, int VerbosityIn) | |
| void | fatalEmissionVeto (string message) |
| int | numberVetoMPIStep () |
| double | pTcalc (const Pythia8::Event &e, int i, int j, int k, int r, int xSRin) |
| double | pTpowheg (const Pythia8::Event &e, int i, int j, bool FSR) |
| double | pTpythia (const Pythia8::Event &e, int RadAfterBranch, int EmtAfterBranch, int RecAfterBranch, bool FSR) |
| ~EmissionVetoHook1 () | |
Private Attributes | |
| bool | accepted |
| int | emittedMode |
| int | MPIvetoOn |
| int | nAcceptSeq |
| int | nFinal |
| int | nFinalExt |
| unsigned long int | nFSRveto |
| unsigned long int | nISRveto |
| int | pTdefMode |
| int | pTemtMode |
| double | pThard |
| int | pThardMode |
| double | pTMPI |
| int | Verbosity |
| int | vetoCount |
| int | vetoOn |
Definition at line 3 of file EmissionVetoHook1.h.
| EmissionVetoHook1::EmissionVetoHook1 | ( | int | nFinalIn, |
| bool | vetoOnIn, | ||
| int | vetoCountIn, | ||
| int | pThardModeIn, | ||
| int | pTemtModeIn, | ||
| int | emittedModeIn, | ||
| int | pTdefModeIn, | ||
| bool | MPIvetoOnIn, | ||
| int | VerbosityIn | ||
| ) | [inline] |
Definition at line 8 of file EmissionVetoHook1.h.
| EmissionVetoHook1::~EmissionVetoHook1 | ( | ) | [inline] |
Definition at line 17 of file EmissionVetoHook1.h.
References gather_cfg::cout, nFSRveto, and nISRveto.
| bool EmissionVetoHook1::canVetoFSREmission | ( | ) | [inline] |
| bool EmissionVetoHook1::canVetoISREmission | ( | ) | [inline] |
| bool EmissionVetoHook1::canVetoMPIEmission | ( | ) | [inline] |
| bool EmissionVetoHook1::canVetoMPIStep | ( | ) | [inline] |
Definition at line 24 of file EmissionVetoHook1.h.
{ return true; }
| bool EmissionVetoHook1::doVetoFSREmission | ( | int | , |
| const Pythia8::Event & | e, | ||
| int | iSys, | ||
| bool | |||
| ) |
Definition at line 384 of file EmissionVetoHook1.cc.
References accepted, gather_cfg::cout, emittedMode, fatalEmissionVeto(), i, j, gen::k, min, nAcceptSeq, nFSRveto, pTcalc(), pTemtMode, pThard, alignCSCRings::r, ntuplemaker::status, swap(), vetoCount, and vetoOn.
{
// Must be radiation from the hard system
if (iSys != 0) return false;
// If we already have accepted 'vetoCount' emissions in a row, do nothing
if (vetoOn && nAcceptSeq >= vetoCount) return false;
// Pythia radiator (before and after), emitted and recoiler (after)
int iRecAft = e.size() - 1;
int iEmt = e.size() - 2;
int iRadAft = e.size() - 3;
int iRadBef = e[iEmt].mother1();
if ( (e[iRecAft].status() != 52 && e[iRecAft].status() != -53) ||
e[iEmt].status() != 51 || e[iRadAft].status() != 51) {
e.list();
fatalEmissionVeto(string("Couldn't find Pythia FSR emission"));
}
// Behaviour based on pTemtMode:
// 0 - pT of emitted w.r.t. radiator before
// 1 - min(pT of emitted w.r.t. all incoming/outgoing)
// 2 - min(pT of all outgoing w.r.t. all incoming/outgoing)
int xSR = (pTemtMode == 0) ? 1 : -1;
int i = (pTemtMode == 0) ? iRadBef : -1;
int k = (pTemtMode == 0) ? iRadAft : -1;
int r = (pTemtMode == 0) ? iRecAft : -1;
// When pTemtMode is 0 or 1, iEmt has been selected
double pTemt = -1.;
if (pTemtMode == 0 || pTemtMode == 1) {
// Which parton is emitted, based on emittedMode:
// 0 - Pythia definition of emitted
// 1 - Pythia definition of radiated after emission
// 2 - Random selection of emitted or radiated after emission
// 3 - Try both emitted and radiated after emission
int j = iRadAft;
if (emittedMode == 0 || (emittedMode == 2 && rndmPtr->flat() < 0.5)) j++;
for (int jLoop = 0; jLoop < 2; jLoop++) {
if (jLoop == 0) pTemt = pTcalc(e, i, j, k, r, xSR);
else if (jLoop == 1) pTemt = min(pTemt, pTcalc(e, i, j, k, r, xSR));
// For emittedMode == 3, have tried iRadAft, now try iEmt
if (emittedMode != 3) break;
if (k != -1) swap(j, k); else j = iEmt;
}
// If pTemtMode is 2, then try all final-state partons as emitted
} else if (pTemtMode == 2) {
pTemt = pTcalc(e, i, -1, k, r, xSR);
}
#ifdef DBGOUTPUT
cout << "doVetoFSREmission: pTemt = " << pTemt << endl << endl;
#endif
// Veto if pTemt > pThard
if (pTemt > pThard) {
nAcceptSeq = 0;
nFSRveto++;
return true;
}
// Else mark that an emission has been accepted and continue
nAcceptSeq++;
accepted = true;
return false;
}
| bool EmissionVetoHook1::doVetoISREmission | ( | int | , |
| const Pythia8::Event & | e, | ||
| int | iSys | ||
| ) |
Definition at line 333 of file EmissionVetoHook1.cc.
References accepted, gather_cfg::cout, fatalEmissionVeto(), i, j, gen::k, nAcceptSeq, nISRveto, pTcalc(), pTemtMode, pThard, alignCSCRings::r, ntuplemaker::status, vetoCount, and vetoOn.
{
// Must be radiation from the hard system
if (iSys != 0) return false;
// If we already have accepted 'vetoCount' emissions in a row, do nothing
if (vetoOn && nAcceptSeq >= vetoCount) return false;
// Pythia radiator after, emitted and recoiler after.
int iRadAft = -1, iEmt = -1, iRecAft = -1;
for (int i = e.size() - 1; i > 0; i--) {
if (iRadAft == -1 && e[i].status() == -41) iRadAft = i;
else if (iEmt == -1 && e[i].status() == 43) iEmt = i;
else if (iRecAft == -1 && e[i].status() == -42) iRecAft = i;
if (iRadAft != -1 && iEmt != -1 && iRecAft != -1) break;
}
if (iRadAft == -1 || iEmt == -1 || iRecAft == -1) {
e.list();
fatalEmissionVeto(string("Couldn't find Pythia ISR emission"));
}
// pTemtMode == 0: pT of emitted w.r.t. radiator
// pTemtMode == 1: min(pT of emitted w.r.t. all incoming/outgoing)
// pTemtMode == 2: min(pT of all outgoing w.r.t. all incoming/outgoing)
int xSR = (pTemtMode == 0) ? 0 : -1;
int i = (pTemtMode == 0) ? iRadAft : -1;
int j = (pTemtMode != 2) ? iEmt : -1;
int k = -1;
int r = (pTemtMode == 0) ? iRecAft : -1;
double pTemt = pTcalc(e, i, j, k, r, xSR);
#ifdef DBGOUTPUT
cout << "doVetoISREmission: pTemt = " << pTemt << endl << endl;
#endif
// Veto if pTemt > pThard
if (pTemt > pThard) {
nAcceptSeq = 0;
nISRveto++;
return true;
}
// Else mark that an emission has been accepted and continue
nAcceptSeq++;
accepted = true;
return false;
}
| bool EmissionVetoHook1::doVetoMPIEmission | ( | int | , |
| const Pythia8::Event & | e | ||
| ) |
Definition at line 458 of file EmissionVetoHook1.cc.
References gather_cfg::cout, MPIvetoOn, and pTMPI.
| bool EmissionVetoHook1::doVetoMPIStep | ( | int | nMPI, |
| const Pythia8::Event & | e | ||
| ) |
Definition at line 248 of file EmissionVetoHook1.cc.
References abs, accepted, prof2calltree::count, gather_cfg::cout, fatalEmissionVeto(), first, i, prof2calltree::last, MPIvetoOn, nAcceptSeq, nFinal, nFinalExt, pTcalc(), pThard, pThardMode, pTMPI, and Verbosity.
{
// Extra check on nMPI
if (nMPI > 1) return false;
// Find if there is a POWHEG emission. Go backwards through the
// event record until there is a non-final particle. Also sum pT and
// find pT_1 for possible MPI vetoing
int count = 0, inonfinal = 0;
double pT1 = 0., pTsum = 0.;
for (int i = e.size() - 1; i > 0; i--) {
inonfinal = i;
if (e[i].isFinal()) {
count++;
pT1 = e[i].pT();
pTsum += e[i].pT();
} else break;
}
nFinal = nFinalExt;
if (nFinal < 0) { // nFinal is not specified from external, try to find out
int first = -1, myid;
int last = -1;
for(int ip = 2; ip < e.size(); ip++) {
myid = e[ip].id();
if(abs(myid) < 6 || abs(myid) == 21) continue;
first = ip;
break;
}
if(first < 0) fatalEmissionVeto(string("signal particles not found"));
for(int ip = first; ip < e.size(); ip++) {
myid = e[ip].id();
if(abs(myid) < 6 || abs(myid) == 21) continue;
last = ip;
}
nFinal = last - inonfinal;
}
// Extra check that we have the correct final state
if (count != nFinal && count != nFinal + 1)
fatalEmissionVeto(string("Wrong number of final state particles in event"));
// Flag if POWHEG radiation present and index
bool isEmt = (count == nFinal) ? false : true;
int iEmt = (isEmt) ? e.size() - 1 : -1;
// If there is no radiation or if pThardMode is 0 then set pThard to QRen.
if (!isEmt || pThardMode == 0) {
pThard = infoPtr->QRen();
// If pThardMode is 1 then the pT of the POWHEG emission is checked against
// all other incoming and outgoing partons, with the minimal value taken
} else if (pThardMode == 1) {
pThard = pTcalc(e, -1, iEmt, -1, -1, -1);
// If pThardMode is 2, then the pT of all final-state partons is checked
// against all other incoming and outgoing partons, with the minimal value
// taken
} else if (pThardMode == 2) {
pThard = pTcalc(e, -1, -1, -1, -1, -1);
}
// Find MPI veto pT if necessary
if (MPIvetoOn) {
pTMPI = (isEmt) ? pTsum / 2. : pT1;
}
if(Verbosity)
cout << "doVetoMPIStep: QFac = " << infoPtr->QFac()
<< ", QRen = " << infoPtr->QRen()
<< ", pThard = " << pThard << endl << endl;
// Initialise other variables
accepted = false;
nAcceptSeq = 0;
// Do not veto the event
return false;
}
| void EmissionVetoHook1::fatalEmissionVeto | ( | string | message | ) |
Definition at line 4 of file EmissionVetoHook1.cc.
References edm::errors::Configuration, and Exception.
Referenced by doVetoFSREmission(), doVetoISREmission(), and doVetoMPIStep().
{
throw edm::Exception(edm::errors::Configuration,"Pythia8Interface")
<< "EmissionVeto: " << message << endl;
}
| int EmissionVetoHook1::numberVetoMPIStep | ( | ) | [inline] |
Definition at line 25 of file EmissionVetoHook1.h.
{ return 1; }
| double EmissionVetoHook1::pTcalc | ( | const Pythia8::Event & | e, |
| int | i, | ||
| int | j, | ||
| int | k, | ||
| int | r, | ||
| int | xSRin | ||
| ) |
Definition at line 117 of file EmissionVetoHook1.cc.
References gather_cfg::cout, min, pTdefMode, pTpowheg(), and pTpythia().
Referenced by doVetoFSREmission(), doVetoISREmission(), and doVetoMPIStep().
{
// Loop over ISR and FSR if necessary
double pTemt = -1., pTnow;
int xSR1 = (xSRin == -1) ? 0 : xSRin;
int xSR2 = (xSRin == -1) ? 2 : xSRin + 1;
for (int xSR = xSR1; xSR < xSR2; xSR++) {
// FSR flag
bool FSR = (xSR == 0) ? false : true;
// If all necessary arguments have been given, then directly calculate.
// POWHEG ISR and FSR, need i and j.
if ((pTdefMode == 0 || pTdefMode == 1) && i > 0 && j > 0) {
pTemt = pTpowheg(e, i, j, (pTdefMode == 0) ? false : FSR);
// Pythia ISR, need i, j and r.
} else if (!FSR && pTdefMode == 2 && i > 0 && j > 0 && r > 0) {
pTemt = pTpythia(e, i, j, r, FSR);
// Pythia FSR, need k, j and r.
} else if (FSR && pTdefMode == 2 && j > 0 && k > 0 && r > 0) {
pTemt = pTpythia(e, k, j, r, FSR);
// Otherwise need to try all possible combintations.
} else {
// Start by finding incoming legs to the hard system after
// branching (radiator after branching, i for ISR).
// Use partonSystemsPtr to find incoming just prior to the
// branching and track mothers.
int iInA = partonSystemsPtr->getInA(0);
int iInB = partonSystemsPtr->getInB(0);
while (e[iInA].mother1() != 1) { iInA = e[iInA].mother1(); }
while (e[iInB].mother1() != 2) { iInB = e[iInB].mother1(); }
// If we do not have j, then try all final-state partons
int jNow = (j > 0) ? j : 0;
int jMax = (j > 0) ? j + 1 : e.size();
for (; jNow < jMax; jNow++) {
// Final-state and coloured jNow or photon only
if (!e[jNow].isFinal()) continue;
if (e[jNow].colType() == 0 && e[jNow].id() != 22) continue;
// POWHEG
if (pTdefMode == 0 || pTdefMode == 1) {
// ISR - only done once as just kinematical pT
if (!FSR) {
pTnow = pTpowheg(e, iInA, jNow, (pTdefMode == 0) ? false : FSR);
if (pTnow > 0.) pTemt = (pTemt < 0) ? pTnow : min(pTemt, pTnow);
// FSR - try all outgoing partons from system before branching
// as i. Note that for the hard system, there is no
// "before branching" information.
} else {
int outSize = partonSystemsPtr->sizeOut(0);
for (int iMem = 0; iMem < outSize; iMem++) {
int iNow = partonSystemsPtr->getOut(0, iMem);
// Coloured only, i != jNow and no carbon copies
if (iNow == jNow || e[iNow].colType() == 0) continue;
if (jNow == e[iNow].daughter1()
&& jNow == e[iNow].daughter2()) continue;
pTnow = pTpowheg(e, iNow, jNow, (pTdefMode == 0)
? false : FSR);
if (pTnow > 0.) pTemt = (pTemt < 0)
? pTnow : min(pTemt, pTnow);
} // for (iMem)
} // if (!FSR)
// Pythia
} else if (pTdefMode == 2) {
// ISR - other incoming as recoiler
if (!FSR) {
pTnow = pTpythia(e, iInA, jNow, iInB, FSR);
if (pTnow > 0.) pTemt = (pTemt < 0) ? pTnow : min(pTemt, pTnow);
pTnow = pTpythia(e, iInB, jNow, iInA, FSR);
if (pTnow > 0.) pTemt = (pTemt < 0) ? pTnow : min(pTemt, pTnow);
// FSR - try all final-state coloured partons as radiator
// after emission (k).
} else {
for (int kNow = 0; kNow < e.size(); kNow++) {
if (kNow == jNow || !e[kNow].isFinal() ||
e[kNow].colType() == 0) continue;
// For this kNow, need to have a recoiler.
// Try two incoming.
pTnow = pTpythia(e, kNow, jNow, iInA, FSR);
if (pTnow > 0.) pTemt = (pTemt < 0)
? pTnow : min(pTemt, pTnow);
pTnow = pTpythia(e, kNow, jNow, iInB, FSR);
if (pTnow > 0.) pTemt = (pTemt < 0)
? pTnow : min(pTemt, pTnow);
// Try all other outgoing.
for (int rNow = 0; rNow < e.size(); rNow++) {
if (rNow == kNow || rNow == jNow ||
!e[rNow].isFinal() || e[rNow].colType() == 0) continue;
pTnow = pTpythia(e, kNow, jNow, rNow, FSR);
if (pTnow > 0.) pTemt = (pTemt < 0)
? pTnow : min(pTemt, pTnow);
} // for (rNow)
} // for (kNow)
} // if (!FSR)
} // if (pTdefMode)
} // for (j)
}
} // for (xSR)
#ifdef DBGOUTPUT
cout << "pTcalc: i = " << i << ", j = " << j << ", k = " << k
<< ", r = " << r << ", xSR = " << xSRin
<< ", pTemt = " << pTemt << endl;
#endif
return pTemt;
}
| double EmissionVetoHook1::pTpowheg | ( | const Pythia8::Event & | e, |
| int | i, | ||
| int | j, | ||
| bool | FSR | ||
| ) |
Definition at line 74 of file EmissionVetoHook1.cc.
References gather_cfg::cout, alignCSCRings::e, j, AlCaHLTBitMon_ParallelJobs::p, and mathSSE::sqrt().
Referenced by pTcalc().
{
// pT value for FSR and ISR
double pTnow = 0.;
if (FSR) {
// POWHEG d_ij (in CM frame). Note that the incoming beams have not
// been updated in the parton systems pointer yet (i.e. prior to any
// potential recoil).
int iInA = partonSystemsPtr->getInA(0);
int iInB = partonSystemsPtr->getInB(0);
double betaZ = - ( e[iInA].pz() + e[iInB].pz() ) /
( e[iInA].e() + e[iInB].e() );
Pythia8::Vec4 iVecBst(e[i].p()), jVecBst(e[j].p());
iVecBst.bst(0., 0., betaZ);
jVecBst.bst(0., 0., betaZ);
pTnow = sqrt( (iVecBst + jVecBst).m2Calc() *
iVecBst.e() * jVecBst.e() /
Pythia8::pow2(iVecBst.e() + jVecBst.e()) );
} else {
// POWHEG pT_ISR is just kinematic pT
pTnow = e[j].pT();
}
// Check result
if (pTnow < 0.) {
cout << "Warning: pTpowheg was negative" << endl;
return -1.;
}
#ifdef DBGOUTPUT
cout << "pTpowheg: i = " << i << ", j = " << j
<< ", pTnow = " << pTnow << endl;
#endif
return pTnow;
}
| double EmissionVetoHook1::pTpythia | ( | const Pythia8::Event & | e, |
| int | RadAfterBranch, | ||
| int | EmtAfterBranch, | ||
| int | RecAfterBranch, | ||
| bool | FSR | ||
| ) |
Definition at line 17 of file EmissionVetoHook1.cc.
References abs, gather_cfg::cout, mathSSE::sqrt(), and z.
Referenced by pTcalc().
{
// Convenient shorthands for later
Pythia8::Vec4 radVec = e[RadAfterBranch].p();
Pythia8::Vec4 emtVec = e[EmtAfterBranch].p();
Pythia8::Vec4 recVec = e[RecAfterBranch].p();
int radID = e[RadAfterBranch].id();
// Calculate virtuality of splitting
double sign = (FSR) ? 1. : -1.;
Pythia8::Vec4 Q(radVec + sign * emtVec);
double Qsq = sign * Q.m2Calc();
// Mass term of radiator
double m2Rad = (abs(radID) >= 4 && abs(radID) < 7) ?
Pythia8::pow2(particleDataPtr->m0(radID)) : 0.;
// z values for FSR and ISR
double z, pTnow;
if (FSR) {
// Construct 2 -> 3 variables
Pythia8::Vec4 sum = radVec + recVec + emtVec;
double m2Dip = sum.m2Calc();
double x1 = 2. * (sum * radVec) / m2Dip;
double x3 = 2. * (sum * emtVec) / m2Dip;
z = x1 / (x1 + x3);
pTnow = z * (1. - z);
} else {
// Construct dipoles before/after splitting
Pythia8::Vec4 qBR(radVec - emtVec + recVec);
Pythia8::Vec4 qAR(radVec + recVec);
z = qBR.m2Calc() / qAR.m2Calc();
pTnow = (1. - z);
}
// Virtuality with correct sign
pTnow *= (Qsq - sign * m2Rad);
// Can get negative pT for massive splittings
if (pTnow < 0.) {
cout << "Warning: pTpythia was negative" << endl;
return -1.;
}
#ifdef DBGOUTPUT
cout << "pTpythia: rad = " << RadAfterBranch << ", emt = "
<< EmtAfterBranch << ", rec = " << RecAfterBranch
<< ", pTnow = " << sqrt(pTnow) << endl;
#endif
// Return pT
return sqrt(pTnow);
}
bool EmissionVetoHook1::accepted [private] |
Definition at line 53 of file EmissionVetoHook1.h.
Referenced by doVetoFSREmission(), doVetoISREmission(), and doVetoMPIStep().
int EmissionVetoHook1::emittedMode [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by doVetoFSREmission().
int EmissionVetoHook1::MPIvetoOn [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by canVetoMPIEmission(), doVetoMPIEmission(), and doVetoMPIStep().
int EmissionVetoHook1::nAcceptSeq [private] |
Definition at line 55 of file EmissionVetoHook1.h.
Referenced by doVetoFSREmission(), doVetoISREmission(), and doVetoMPIStep().
int EmissionVetoHook1::nFinal [private] |
Definition at line 51 of file EmissionVetoHook1.h.
Referenced by doVetoMPIStep().
int EmissionVetoHook1::nFinalExt [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by doVetoMPIStep().
unsigned long int EmissionVetoHook1::nFSRveto [private] |
Definition at line 57 of file EmissionVetoHook1.h.
Referenced by doVetoFSREmission(), and ~EmissionVetoHook1().
unsigned long int EmissionVetoHook1::nISRveto [private] |
Definition at line 57 of file EmissionVetoHook1.h.
Referenced by doVetoISREmission(), and ~EmissionVetoHook1().
int EmissionVetoHook1::pTdefMode [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by pTcalc().
int EmissionVetoHook1::pTemtMode [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by doVetoFSREmission(), and doVetoISREmission().
double EmissionVetoHook1::pThard [private] |
Definition at line 52 of file EmissionVetoHook1.h.
Referenced by doVetoFSREmission(), doVetoISREmission(), and doVetoMPIStep().
int EmissionVetoHook1::pThardMode [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by doVetoMPIStep().
double EmissionVetoHook1::pTMPI [private] |
Definition at line 52 of file EmissionVetoHook1.h.
Referenced by doVetoMPIEmission(), and doVetoMPIStep().
int EmissionVetoHook1::Verbosity [private] |
Definition at line 58 of file EmissionVetoHook1.h.
Referenced by doVetoMPIStep().
int EmissionVetoHook1::vetoCount [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by doVetoFSREmission(), and doVetoISREmission().
int EmissionVetoHook1::vetoOn [private] |
Definition at line 49 of file EmissionVetoHook1.h.
Referenced by canVetoFSREmission(), canVetoISREmission(), doVetoFSREmission(), and doVetoISREmission().
1.7.3