CMS 3D CMS Logo

List of all members | Public Member Functions | Private Attributes
EmissionVetoHook1 Class Reference

#include <EmissionVetoHook1.h>

Inheritance diagram for EmissionVetoHook1:

Public Member Functions

bool canVetoFSREmission () override
 
bool canVetoISREmission () override
 
bool canVetoMPIEmission () override
 
bool canVetoMPIStep () override
 
bool doVetoFSREmission (int, const Pythia8::Event &e, int iSys, bool) override
 
bool doVetoISREmission (int, const Pythia8::Event &e, int iSys) override
 
bool doVetoMPIEmission (int, const Pythia8::Event &e) override
 
bool doVetoMPIStep (int nMPI, const Pythia8::Event &e) override
 
 EmissionVetoHook1 (int nFinalIn, bool vetoOnIn, int vetoCountIn, int pThardModeIn, int pTemtModeIn, int emittedModeIn, int pTdefModeIn, bool MPIvetoOnIn, int QEDvetoModeIn, int nFinalModeIn, int VerbosityIn)
 
void fatalEmissionVeto (std::string message)
 
int numberVetoMPIStep () override
 
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 () override
 

Private Attributes

bool accepted
 
int emittedMode
 
bool isEmt
 
int MPIvetoOn
 
int nAcceptSeq
 
int nFinal
 
int nFinalExt
 
int nFinalMode
 
unsigned long int nFSRveto
 
unsigned long int nISRveto
 
int pTdefMode
 
int pTemtMode
 
double pThard
 
int pThardMode
 
double pTMPI
 
int QEDvetoMode
 
int Verbosity
 
int vetoCount
 
int vetoOn
 

Detailed Description

Definition at line 3 of file EmissionVetoHook1.h.

Constructor & Destructor Documentation

◆ EmissionVetoHook1()

EmissionVetoHook1::EmissionVetoHook1 ( int  nFinalIn,
bool  vetoOnIn,
int  vetoCountIn,
int  pThardModeIn,
int  pTemtModeIn,
int  emittedModeIn,
int  pTdefModeIn,
bool  MPIvetoOnIn,
int  QEDvetoModeIn,
int  nFinalModeIn,
int  VerbosityIn 
)
inline

Definition at line 6 of file EmissionVetoHook1.h.

17  : nFinalExt(nFinalIn),
18  vetoOn(vetoOnIn),
19  vetoCount(vetoCountIn),
20  pThardMode(pThardModeIn),
21  pTemtMode(pTemtModeIn),
22  emittedMode(emittedModeIn),
23  pTdefMode(pTdefModeIn),
24  MPIvetoOn(MPIvetoOnIn),
25  QEDvetoMode(QEDvetoModeIn),
26  nFinalMode(nFinalModeIn),
27  nISRveto(0),
28  nFSRveto(0),
29  Verbosity(VerbosityIn) {}
unsigned long int nISRveto
unsigned long int nFSRveto

◆ ~EmissionVetoHook1()

EmissionVetoHook1::~EmissionVetoHook1 ( )
inlineoverride

Definition at line 30 of file EmissionVetoHook1.h.

References gather_cfg::cout, nFSRveto, and nISRveto.

30  {
31  std::cout << "Number of ISR vetoed = " << nISRveto << std::endl;
32  std::cout << "Number of FSR vetoed = " << nFSRveto << std::endl;
33  }
unsigned long int nISRveto
unsigned long int nFSRveto

Member Function Documentation

◆ canVetoFSREmission()

bool EmissionVetoHook1::canVetoFSREmission ( )
inlineoverride

Definition at line 44 of file EmissionVetoHook1.h.

References vetoOn.

44 { return vetoOn; }

◆ canVetoISREmission()

bool EmissionVetoHook1::canVetoISREmission ( )
inlineoverride

Definition at line 41 of file EmissionVetoHook1.h.

References vetoOn.

41 { return vetoOn; }

◆ canVetoMPIEmission()

bool EmissionVetoHook1::canVetoMPIEmission ( )
inlineoverride

Definition at line 47 of file EmissionVetoHook1.h.

References MPIvetoOn.

47 { return MPIvetoOn; }

◆ canVetoMPIStep()

bool EmissionVetoHook1::canVetoMPIStep ( )
inlineoverride

Definition at line 37 of file EmissionVetoHook1.h.

37 { return true; }

◆ doVetoFSREmission()

bool EmissionVetoHook1::doVetoFSREmission ( int  ,
const Pythia8::Event &  e,
int  iSys,
bool  inResonance 
)
override

Definition at line 439 of file EmissionVetoHook1.cc.

References gather_cfg::cout, MillePedeFileConverter_cfg::e, mps_fire::i, dqmiolumiharvest::j, dqmdumpme::k, SiStripPI::min, mps_update::status, AlCaHLTBitMon_QueryRunRegistry::string, and std::swap().

439  {
440  // Must be radiation from the hard system
441  if (iSys != 0)
442  return false;
443 
444  // only use for outside resonance vetos in combination with bb4l:FSREmission:veto
445  if (!inResonance && settingsPtr->flag("POWHEG:bb4l:FSREmission:veto") == 1)
446  return false;
447 
448  // If we already have accepted 'vetoCount' emissions in a row, do nothing
449  if (vetoOn && nAcceptSeq >= vetoCount)
450  return false;
451 
452  // Pythia radiator (before and after), emitted and recoiler (after)
453  int iRecAft = e.size() - 1;
454  int iEmt = e.size() - 2;
455  int iRadAft = e.size() - 3;
456  int iRadBef = e[iEmt].mother1();
457  if ((e[iRecAft].status() != 52 && e[iRecAft].status() != -53) || e[iEmt].status() != 51 ||
458  e[iRadAft].status() != 51) {
459  e.list();
460  fatalEmissionVeto(std::string("Couldn't find Pythia FSR emission"));
461  }
462 
463  // Behaviour based on pTemtMode:
464  // 0 - pT of emitted w.r.t. radiator before
465  // 1 - std::min(pT of emitted w.r.t. all incoming/outgoing)
466  // 2 - std::min(pT of all outgoing w.r.t. all incoming/outgoing)
467  int xSR = (pTemtMode == 0) ? 1 : -1;
468  int i = (pTemtMode == 0) ? iRadBef : -1;
469  int k = (pTemtMode == 0) ? iRadAft : -1;
470  int r = (pTemtMode == 0) ? iRecAft : -1;
471 
472  // When pTemtMode is 0 or 1, iEmt has been selected
473  double pTemt = -1.;
474  if (pTemtMode == 0 || pTemtMode == 1) {
475  // Which parton is emitted, based on emittedMode:
476  // 0 - Pythia definition of emitted
477  // 1 - Pythia definition of radiated after emission
478  // 2 - Random selection of emitted or radiated after emission
479  // 3 - Try both emitted and radiated after emission
480  int j = iRadAft;
481  if (emittedMode == 0 || (emittedMode == 2 && rndmPtr->flat() < 0.5))
482  j++;
483 
484  for (int jLoop = 0; jLoop < 2; jLoop++) {
485  if (jLoop == 0)
486  pTemt = pTcalc(e, i, j, k, r, xSR);
487  else if (jLoop == 1)
488  pTemt = std::min(pTemt, pTcalc(e, i, j, k, r, xSR));
489 
490  // For emittedMode == 3, have tried iRadAft, now try iEmt
491  if (emittedMode != 3)
492  break;
493  if (k != -1)
494  std::swap(j, k);
495  else
496  j = iEmt;
497  }
498 
499  // If pTemtMode is 2, then try all final-state partons as emitted
500  } else if (pTemtMode == 2) {
501  pTemt = pTcalc(e, i, -1, k, r, xSR);
502  }
503 
504 #ifdef DBGOUTPUT
505  std::cout << "doVetoFSREmission: pTemt = " << pTemt << std::endl << std::endl;
506 #endif
507 
508  // If a Born configuration, and a colorless FS, and QEDvetoMode=2,
509  // then don't veto photons, W, or Z harder than pThard
510  bool vetoParton = (!isEmt && e[iEmt].colType() == 0 && QEDvetoMode == 2) ? false : true;
511 
512  // Veto if pTemt > pThard
513  if (pTemt > pThard) {
514  if (!vetoParton) {
515  // Don't veto ANY emissions afterwards
516  nAcceptSeq = vetoCount - 1;
517  } else {
518  nAcceptSeq = 0;
519  nFSRveto++;
520  return true;
521  }
522  }
523 
524  // Else mark that an emission has been accepted and continue
525  nAcceptSeq++;
526  accepted = true;
527  return false;
528 }
void swap(edm::DataFrameContainer &lhs, edm::DataFrameContainer &rhs)
unsigned long int nFSRveto
double pTcalc(const Pythia8::Event &e, int i, int j, int k, int r, int xSRin)
void fatalEmissionVeto(std::string message)

◆ doVetoISREmission()

bool EmissionVetoHook1::doVetoISREmission ( int  ,
const Pythia8::Event &  e,
int  iSys 
)
override

Definition at line 373 of file EmissionVetoHook1.cc.

References gather_cfg::cout, MillePedeFileConverter_cfg::e, mps_fire::i, dqmiolumiharvest::j, dqmdumpme::k, mps_update::status, and AlCaHLTBitMon_QueryRunRegistry::string.

373  {
374  // Must be radiation from the hard system
375  if (iSys != 0)
376  return false;
377 
378  // If we already have accepted 'vetoCount' emissions in a row, do nothing
379  if (vetoOn && nAcceptSeq >= vetoCount)
380  return false;
381 
382  // Pythia radiator after, emitted and recoiler after.
383  int iRadAft = -1, iEmt = -1, iRecAft = -1;
384  for (int i = e.size() - 1; i > 0; i--) {
385  if (iRadAft == -1 && e[i].status() == -41)
386  iRadAft = i;
387  else if (iEmt == -1 && e[i].status() == 43)
388  iEmt = i;
389  else if (iRecAft == -1 && e[i].status() == -42)
390  iRecAft = i;
391  if (iRadAft != -1 && iEmt != -1 && iRecAft != -1)
392  break;
393  }
394  if (iRadAft == -1 || iEmt == -1 || iRecAft == -1) {
395  e.list();
396  fatalEmissionVeto(std::string("Couldn't find Pythia ISR emission"));
397  }
398 
399  // pTemtMode == 0: pT of emitted w.r.t. radiator
400  // pTemtMode == 1: std::min(pT of emitted w.r.t. all incoming/outgoing)
401  // pTemtMode == 2: std::min(pT of all outgoing w.r.t. all incoming/outgoing)
402  int xSR = (pTemtMode == 0) ? 0 : -1;
403  int i = (pTemtMode == 0) ? iRadAft : -1;
404  int j = (pTemtMode != 2) ? iEmt : -1;
405  int k = -1;
406  int r = (pTemtMode == 0) ? iRecAft : -1;
407  double pTemt = pTcalc(e, i, j, k, r, xSR);
408 
409 #ifdef DBGOUTPUT
410  std::cout << "doVetoISREmission: pTemt = " << pTemt << std::endl << std::endl;
411 #endif
412 
413  // If a Born configuration, and a colorless FS, and QEDvetoMode=2,
414  // then don't veto photons, W, or Z harder than pThard
415  bool vetoParton = (!isEmt && e[iEmt].colType() == 0 && QEDvetoMode == 2) ? false : true;
416 
417  // Veto if pTemt > pThard
418  if (pTemt > pThard) {
419  if (!vetoParton) {
420  // Don't veto ANY emissions afterwards
421  nAcceptSeq = vetoCount - 1;
422  } else {
423  nAcceptSeq = 0;
424  nISRveto++;
425  return true;
426  }
427  }
428 
429  // Else mark that an emission has been accepted and continue
430  nAcceptSeq++;
431  accepted = true;
432  return false;
433 }
unsigned long int nISRveto
double pTcalc(const Pythia8::Event &e, int i, int j, int k, int r, int xSRin)
void fatalEmissionVeto(std::string message)

◆ doVetoMPIEmission()

bool EmissionVetoHook1::doVetoMPIEmission ( int  ,
const Pythia8::Event &  e 
)
override

Definition at line 534 of file EmissionVetoHook1.cc.

References gather_cfg::cout, and MillePedeFileConverter_cfg::e.

534  {
535  if (MPIvetoOn) {
536  if (e[e.size() - 1].pT() > pTMPI) {
537 #ifdef DBGOUTPUT
538  std::cout << "doVetoMPIEmission: pTnow = " << e[e.size() - 1].pT() << ", pTMPI = " << pTMPI << std::endl
539  << std::endl;
540 #endif
541  return true;
542  }
543  }
544  return false;
545 }

◆ doVetoMPIStep()

bool EmissionVetoHook1::doVetoMPIStep ( int  nMPI,
const Pythia8::Event &  e 
)
override

Definition at line 253 of file EmissionVetoHook1.cc.

References funct::abs(), submitPVResolutionJobs::count, gather_cfg::cout, MillePedeFileConverter_cfg::e, first, mps_fire::i, dqmdumpme::last, HLT_2022v15_cff::pT1, AlCaHLTBitMon_QueryRunRegistry::string, and remoteMonitoring_LASER_era2018_cfg::Verbosity.

253  {
254  if (nFinalMode == 3 && pThardMode != 0)
256  "When nFinalMode is set to 3, ptHardMode should be set to 0, since the emission variables in doVetoMPIStep are "
257  "not set correctly case when there are three possible particle Born particle counts."));
258 
259  // Extra check on nMPI
260  if (nMPI > 1)
261  return false;
262 
263  // Find if there is a POWHEG emission. Go backwards through the
264  // event record until there is a non-final particle. Also sum pT and
265  // find pT_1 for possible MPI vetoing
266  // Flag if POWHEG radiation is present and index at the same time
267  int count = 0, inonfinal = 0;
268  double pT1 = 0., pTsum = 0.;
269  isEmt = false;
270  int iEmt = -1;
271 
272  for (int i = e.size() - 1; i > 0; i--) {
273  inonfinal = i;
274  if (e[i].isFinal()) {
275  count++;
276  pT1 = e[i].pT();
277  pTsum += e[i].pT();
278  // the following was added for bbbar4l and will be triggered by specifying nfinalmode == 2
279  // the solution provided by Tomas may not be process independent but should work for hvq and bb4l
280  // if the particle is a light quark or a gluon and
281  // comes for a light quark or a gluon (not a resonance, not a top)
282  // then it is the POWHEG emission (hvq) or the POWHEG emission in production (bb4l)
283  if (nFinalMode == 2) {
284  if ((abs(e[i].id()) < 6 || e[i].id() == 21) &&
285  (abs(e[e[i].mother1()].id()) < 6 || e[e[i].mother1()].id() == 21)) {
286  isEmt = true;
287  iEmt = i;
288  }
289  }
290  } else
291  break;
292  }
293 
294  nFinal = nFinalExt;
295  if (nFinal < 0 || nFinalMode == 1) { // nFinal is not specified from external, try to find out
296  int first = -1, myid;
297  int last = -1;
298  for (int ip = 2; ip < e.size(); ip++) {
299  myid = e[ip].id();
300  if (abs(myid) < 6 || abs(myid) == 21)
301  continue;
302  first = ip;
303  break;
304  }
305  if (first < 0)
306  fatalEmissionVeto(std::string("signal particles not found"));
307  for (int ip = first; ip < e.size(); ip++) {
308  myid = e[ip].id();
309  if (abs(myid) < 6 || abs(myid) == 21)
310  continue;
311  last = ip;
312  }
313  nFinal = last - inonfinal;
314  }
315 
316  // don't perform a cross check in case of nfinalmode == 2
317  if (nFinalMode != 2) {
318  // Extra check that we have the correct final state
319  // In POWHEG WGamma, both w+0/1 jets and w+gamma+0/1 jets are generated at the same time, which leads to three different possible numbers of particles
320  // Normally, there would be only two possible numbers of particles for X and X+1 jet
321  if (nFinalMode == 3) {
322  if (count != nFinal && count != nFinal + 1 && count != nFinal - 1)
323  fatalEmissionVeto(std::string("Wrong number of final state particles in event"));
324  } else {
325  if (count != nFinal && count != nFinal + 1)
326  fatalEmissionVeto(std::string("Wrong number of final state particles in event"));
327  }
328  // Flag if POWHEG radiation present and index
329  if (count == nFinal + 1)
330  isEmt = true;
331  if (isEmt)
332  iEmt = e.size() - 1;
333  }
334 
335  // If there is no radiation or if pThardMode is 0 then set pThard to QRen.
336  if (!isEmt || pThardMode == 0) {
337  pThard = infoPtr->QRen();
338 
339  // If pThardMode is 1 then the pT of the POWHEG emission is checked against
340  // all other incoming and outgoing partons, with the minimal value taken
341  } else if (pThardMode == 1) {
342  pThard = pTcalc(e, -1, iEmt, -1, -1, -1);
343 
344  // If pThardMode is 2, then the pT of all final-state partons is checked
345  // against all other incoming and outgoing partons, with the minimal value
346  // taken
347  } else if (pThardMode == 2) {
348  pThard = pTcalc(e, -1, -1, -1, -1, -1);
349  }
350 
351  // Find MPI veto pT if necessary
352  if (MPIvetoOn) {
353  pTMPI = (isEmt) ? pTsum / 2. : pT1;
354  }
355 
356  if (Verbosity)
357  std::cout << "doVetoMPIStep: QFac = " << infoPtr->QFac() << ", QRen = " << infoPtr->QRen()
358  << ", pThard = " << pThard << std::endl
359  << std::endl;
360 
361  // Initialise other variables
362  accepted = false;
363  nAcceptSeq = 0; // should not reset nISRveto = nFSRveto = nFSRvetoBB4l here
364 
365  // Do not veto the event
366  return false;
367 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double pTcalc(const Pythia8::Event &e, int i, int j, int k, int r, int xSRin)
void fatalEmissionVeto(std::string message)

◆ fatalEmissionVeto()

void EmissionVetoHook1::fatalEmissionVeto ( std::string  message)

Definition at line 8 of file EmissionVetoHook1.cc.

References edm::errors::Configuration, and Exception.

8  {
9  throw edm::Exception(edm::errors::Configuration, "Pythia8Interface") << "EmissionVeto: " << message << std::endl;
10 }

◆ numberVetoMPIStep()

int EmissionVetoHook1::numberVetoMPIStep ( )
inlineoverride

Definition at line 38 of file EmissionVetoHook1.h.

38 { return 1; }

◆ pTcalc()

double EmissionVetoHook1::pTcalc ( const Pythia8::Event &  e,
int  i,
int  j,
int  k,
int  r,
int  xSRin 
)

Definition at line 112 of file EmissionVetoHook1.cc.

References gather_cfg::cout, MillePedeFileConverter_cfg::e, mps_fire::i, dqmiolumiharvest::j, dqmdumpme::k, and SiStripPI::min.

112  {
113  // Loop over ISR and FSR if necessary
114  double pTemt = -1., pTnow;
115  int xSR1 = (xSRin == -1) ? 0 : xSRin;
116  int xSR2 = (xSRin == -1) ? 2 : xSRin + 1;
117  for (int xSR = xSR1; xSR < xSR2; xSR++) {
118  // FSR flag
119  bool FSR = (xSR == 0) ? false : true;
120 
121  // If all necessary arguments have been given, then directly calculate.
122  // POWHEG ISR and FSR, need i and j.
123  if ((pTdefMode == 0 || pTdefMode == 1) && i > 0 && j > 0) {
124  pTemt = pTpowheg(e, i, j, (pTdefMode == 0) ? false : FSR);
125 
126  // Pythia ISR, need i, j and r.
127  } else if (!FSR && pTdefMode == 2 && i > 0 && j > 0 && r > 0) {
128  pTemt = pTpythia(e, i, j, r, FSR);
129 
130  // Pythia FSR, need k, j and r.
131  } else if (FSR && pTdefMode == 2 && j > 0 && k > 0 && r > 0) {
132  pTemt = pTpythia(e, k, j, r, FSR);
133 
134  // Otherwise need to try all possible combintations.
135  } else {
136  // Start by finding incoming legs to the hard system after
137  // branching (radiator after branching, i for ISR).
138  // Use partonSystemsPtr to find incoming just prior to the
139  // branching and track mothers.
140  int iInA = partonSystemsPtr->getInA(0);
141  int iInB = partonSystemsPtr->getInB(0);
142  while (e[iInA].mother1() != 1) {
143  iInA = e[iInA].mother1();
144  }
145  while (e[iInB].mother1() != 2) {
146  iInB = e[iInB].mother1();
147  }
148 
149  // If we do not have j, then try all final-state partons
150  int jNow = (j > 0) ? j : 0;
151  int jMax = (j > 0) ? j + 1 : e.size();
152  for (; jNow < jMax; jNow++) {
153  // Final-state only
154  if (!e[jNow].isFinal())
155  continue;
156  // Exclude photons (and W/Z!)
157  if (QEDvetoMode == 0 && e[jNow].colType() == 0)
158  continue;
159 
160  // POWHEG
161  if (pTdefMode == 0 || pTdefMode == 1) {
162  // ISR - only done once as just kinematical pT
163  if (!FSR) {
164  pTnow = pTpowheg(e, iInA, jNow, (pTdefMode == 0) ? false : FSR);
165  if (pTnow > 0.)
166  pTemt = (pTemt < 0) ? pTnow : std::min(pTemt, pTnow);
167 
168  // FSR - try all outgoing partons from system before branching
169  // as i. Note that for the hard system, there is no
170  // "before branching" information.
171  } else {
172  int outSize = partonSystemsPtr->sizeOut(0);
173  for (int iMem = 0; iMem < outSize; iMem++) {
174  int iNow = partonSystemsPtr->getOut(0, iMem);
175 
176  // if i != jNow and no carbon copies
177  if (iNow == jNow)
178  continue;
179  // Exlude photons (and W/Z!)
180  if (QEDvetoMode == 0 && e[iNow].colType() == 0)
181  continue;
182  if (jNow == e[iNow].daughter1() && jNow == e[iNow].daughter2())
183  continue;
184 
185  pTnow = pTpowheg(e, iNow, jNow, (pTdefMode == 0) ? false : FSR);
186  if (pTnow > 0.)
187  pTemt = (pTemt < 0) ? pTnow : std::min(pTemt, pTnow);
188  } // for (iMem)
189 
190  } // if (!FSR)
191 
192  // Pythia
193  } else if (pTdefMode == 2) {
194  // ISR - other incoming as recoiler
195  if (!FSR) {
196  pTnow = pTpythia(e, iInA, jNow, iInB, FSR);
197  if (pTnow > 0.)
198  pTemt = (pTemt < 0) ? pTnow : std::min(pTemt, pTnow);
199  pTnow = pTpythia(e, iInB, jNow, iInA, FSR);
200  if (pTnow > 0.)
201  pTemt = (pTemt < 0) ? pTnow : std::min(pTemt, pTnow);
202 
203  // FSR - try all final-state coloured partons as radiator
204  // after emission (k).
205  } else {
206  for (int kNow = 0; kNow < e.size(); kNow++) {
207  if (kNow == jNow || !e[kNow].isFinal())
208  continue;
209  if (QEDvetoMode == 0 && e[kNow].colType() == 0)
210  continue;
211 
212  // For this kNow, need to have a recoiler.
213  // Try two incoming.
214  pTnow = pTpythia(e, kNow, jNow, iInA, FSR);
215  if (pTnow > 0.)
216  pTemt = (pTemt < 0) ? pTnow : std::min(pTemt, pTnow);
217  pTnow = pTpythia(e, kNow, jNow, iInB, FSR);
218  if (pTnow > 0.)
219  pTemt = (pTemt < 0) ? pTnow : std::min(pTemt, pTnow);
220 
221  // Try all other outgoing.
222  for (int rNow = 0; rNow < e.size(); rNow++) {
223  if (rNow == kNow || rNow == jNow || !e[rNow].isFinal())
224  continue;
225  if (QEDvetoMode == 0 && e[rNow].colType() == 0)
226  continue;
227  pTnow = pTpythia(e, kNow, jNow, rNow, FSR);
228  if (pTnow > 0.)
229  pTemt = (pTemt < 0) ? pTnow : std::min(pTemt, pTnow);
230  } // for (rNow)
231 
232  } // for (kNow)
233  } // if (!FSR)
234  } // if (pTdefMode)
235  } // for (j)
236  }
237  } // for (xSR)
238 
239 #ifdef DBGOUTPUT
240  std::cout << "pTcalc: i = " << i << ", j = " << j << ", k = " << k << ", r = " << r << ", xSR = " << xSRin
241  << ", pTemt = " << pTemt << std::endl;
242 #endif
243 
244  return pTemt;
245 }
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)

◆ pTpowheg()

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, MillePedeFileConverter_cfg::e, mps_fire::i, dqmiolumiharvest::j, AlCaHLTBitMon_ParallelJobs::p, L1TauEmu::pow2(), and mathSSE::sqrt().

74  {
75  // pT value for FSR and ISR
76  double pTnow = 0.;
77  if (FSR) {
78  // POWHEG d_ij (in CM frame). Note that the incoming beams have not
79  // been updated in the parton systems pointer yet (i.e. prior to any
80  // potential recoil).
81  int iInA = partonSystemsPtr->getInA(0);
82  int iInB = partonSystemsPtr->getInB(0);
83  double betaZ = -(e[iInA].pz() + e[iInB].pz()) / (e[iInA].e() + e[iInB].e());
84  Pythia8::Vec4 iVecBst(e[i].p()), jVecBst(e[j].p());
85  iVecBst.bst(0., 0., betaZ);
86  jVecBst.bst(0., 0., betaZ);
87  pTnow = sqrt((iVecBst + jVecBst).m2Calc() * iVecBst.e() * jVecBst.e() / Pythia8::pow2(iVecBst.e() + jVecBst.e()));
88 
89  } else {
90  // POWHEG pT_ISR is just kinematic pT
91  pTnow = e[j].pT();
92  }
93 
94  // Check result
95  if (pTnow < 0.) {
96  std::cout << "Warning: pTpowheg was negative" << std::endl;
97  return -1.;
98  }
99 
100 #ifdef DBGOUTPUT
101  std::cout << "pTpowheg: i = " << i << ", j = " << j << ", pTnow = " << pTnow << std::endl;
102 #endif
103 
104  return pTnow;
105 }
constexpr int pow2(int x)
Definition: TauNNIdHW.h:59
T sqrt(T t)
Definition: SSEVec.h:19
ExtVec< T, 4 > Vec4
Definition: ExtVec.h:64

◆ pTpythia()

double EmissionVetoHook1::pTpythia ( const Pythia8::Event &  e,
int  RadAfterBranch,
int  EmtAfterBranch,
int  RecAfterBranch,
bool  FSR 
)

Definition at line 20 of file EmissionVetoHook1.cc.

References funct::abs(), gather_cfg::cout, MillePedeFileConverter_cfg::e, L1TauEmu::pow2(), Validation_hcalonly_cfi::sign, and mathSSE::sqrt().

21  {
22  // Convenient shorthands for later
23  Pythia8::Vec4 radVec = e[RadAfterBranch].p();
24  Pythia8::Vec4 emtVec = e[EmtAfterBranch].p();
25  Pythia8::Vec4 recVec = e[RecAfterBranch].p();
26  int radID = e[RadAfterBranch].id();
27 
28  // Calculate virtuality of splitting
29  double sign = (FSR) ? 1. : -1.;
30  Pythia8::Vec4 Q(radVec + sign * emtVec);
31  double Qsq = sign * Q.m2Calc();
32 
33  // Mass term of radiator
34  double m2Rad = (abs(radID) >= 4 && abs(radID) < 7) ? Pythia8::pow2(particleDataPtr->m0(radID)) : 0.;
35 
36  // z values for FSR and ISR
37  double z, pTnow;
38  if (FSR) {
39  // Construct 2 -> 3 variables
40  Pythia8::Vec4 sum = radVec + recVec + emtVec;
41  double m2Dip = sum.m2Calc();
42  double x1 = 2. * (sum * radVec) / m2Dip;
43  double x3 = 2. * (sum * emtVec) / m2Dip;
44  z = x1 / (x1 + x3);
45  pTnow = z * (1. - z);
46 
47  } else {
48  // Construct dipoles before/after splitting
49  Pythia8::Vec4 qBR(radVec - emtVec + recVec);
50  Pythia8::Vec4 qAR(radVec + recVec);
51  z = qBR.m2Calc() / qAR.m2Calc();
52  pTnow = (1. - z);
53  }
54 
55  // Virtuality with correct sign
56  pTnow *= (Qsq - sign * m2Rad);
57 
58  // Can get negative pT for massive splittings
59  if (pTnow < 0.) {
60  std::cout << "Warning: pTpythia was negative" << std::endl;
61  return -1.;
62  }
63 
64 #ifdef DBGOUTPUT
65  std::cout << "pTpythia: rad = " << RadAfterBranch << ", emt = " << EmtAfterBranch << ", rec = " << RecAfterBranch
66  << ", pTnow = " << sqrt(pTnow) << std::endl;
67 #endif
68 
69  // Return pT
70  return sqrt(pTnow);
71 }
constexpr int pow2(int x)
Definition: TauNNIdHW.h:59
T sqrt(T t)
Definition: SSEVec.h:19
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
ExtVec< T, 4 > Vec4
Definition: ExtVec.h:64

Member Data Documentation

◆ accepted

bool EmissionVetoHook1::accepted
private

Definition at line 64 of file EmissionVetoHook1.h.

◆ emittedMode

int EmissionVetoHook1::emittedMode
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ isEmt

bool EmissionVetoHook1::isEmt
private

Definition at line 64 of file EmissionVetoHook1.h.

◆ MPIvetoOn

int EmissionVetoHook1::MPIvetoOn
private

Definition at line 61 of file EmissionVetoHook1.h.

Referenced by canVetoMPIEmission().

◆ nAcceptSeq

int EmissionVetoHook1::nAcceptSeq
private

Definition at line 66 of file EmissionVetoHook1.h.

◆ nFinal

int EmissionVetoHook1::nFinal
private

Definition at line 62 of file EmissionVetoHook1.h.

◆ nFinalExt

int EmissionVetoHook1::nFinalExt
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ nFinalMode

int EmissionVetoHook1::nFinalMode
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ nFSRveto

unsigned long int EmissionVetoHook1::nFSRveto
private

Definition at line 68 of file EmissionVetoHook1.h.

Referenced by ~EmissionVetoHook1().

◆ nISRveto

unsigned long int EmissionVetoHook1::nISRveto
private

Definition at line 68 of file EmissionVetoHook1.h.

Referenced by ~EmissionVetoHook1().

◆ pTdefMode

int EmissionVetoHook1::pTdefMode
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ pTemtMode

int EmissionVetoHook1::pTemtMode
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ pThard

double EmissionVetoHook1::pThard
private

Definition at line 63 of file EmissionVetoHook1.h.

◆ pThardMode

int EmissionVetoHook1::pThardMode
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ pTMPI

double EmissionVetoHook1::pTMPI
private

Definition at line 63 of file EmissionVetoHook1.h.

◆ QEDvetoMode

int EmissionVetoHook1::QEDvetoMode
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ Verbosity

int EmissionVetoHook1::Verbosity
private

Definition at line 69 of file EmissionVetoHook1.h.

◆ vetoCount

int EmissionVetoHook1::vetoCount
private

Definition at line 61 of file EmissionVetoHook1.h.

◆ vetoOn

int EmissionVetoHook1::vetoOn
private

Definition at line 61 of file EmissionVetoHook1.h.

Referenced by canVetoFSREmission(), and canVetoISREmission().