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 ( 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 8 of file EmissionVetoHook1.h.

11  :
12  nFinalExt(nFinalIn),
13  vetoOn(vetoOnIn), vetoCount(vetoCountIn),
14  pThardMode(pThardModeIn), pTemtMode(pTemtModeIn),
15  emittedMode(emittedModeIn), pTdefMode(pTdefModeIn),
16  MPIvetoOn(MPIvetoOnIn), QEDvetoMode(QEDvetoModeIn),
17  nFinalMode(nFinalModeIn), nISRveto(0), nFSRveto(0),
18  Verbosity(VerbosityIn) {}
unsigned long int nISRveto
unsigned long int nFSRveto
EmissionVetoHook1::~EmissionVetoHook1 ( )
inlineoverride

Definition at line 19 of file EmissionVetoHook1.h.

References gather_cfg::cout, nFSRveto, and nISRveto.

19  {
20  std::cout << "Number of ISR vetoed = " << nISRveto << std::endl;
21  std::cout << "Number of FSR vetoed = " << nFSRveto << std::endl;
22  }
unsigned long int nISRveto
unsigned long int nFSRveto

Member Function Documentation

bool EmissionVetoHook1::canVetoFSREmission ( )
inlineoverride

Definition at line 33 of file EmissionVetoHook1.h.

References doVetoFSREmission(), and vetoOn.

33 { return vetoOn; }
bool EmissionVetoHook1::canVetoISREmission ( )
inlineoverride

Definition at line 30 of file EmissionVetoHook1.h.

References doVetoISREmission(), and vetoOn.

30 { return vetoOn; }
bool EmissionVetoHook1::canVetoMPIEmission ( )
inlineoverride
bool EmissionVetoHook1::canVetoMPIStep ( )
inlineoverride

Definition at line 26 of file EmissionVetoHook1.h.

26 { return true; }
bool EmissionVetoHook1::doVetoFSREmission ( int  ,
const Pythia8::Event &  e,
int  iSys,
bool  inResonance 
)
override

Definition at line 418 of file EmissionVetoHook1.cc.

References gather_cfg::cout, mps_fire::i, gen::k, min(), alignCSCRings::r, mps_update::status, AlCaHLTBitMon_QueryRunRegistry::string, and std::swap().

Referenced by canVetoFSREmission().

418  {
419 
420  // Must be radiation from the hard system
421  if (iSys != 0) return false;
422 
423  // only use for outside resonance vetos in combination with bb4l:FSREmission:veto
424  if (!inResonance && settingsPtr->flag("POWHEG:bb4l:FSREmission:veto")==1) return false;
425 
426  // If we already have accepted 'vetoCount' emissions in a row, do nothing
427  if (vetoOn && nAcceptSeq >= vetoCount) return false;
428 
429  // Pythia radiator (before and after), emitted and recoiler (after)
430  int iRecAft = e.size() - 1;
431  int iEmt = e.size() - 2;
432  int iRadAft = e.size() - 3;
433  int iRadBef = e[iEmt].mother1();
434  if ( (e[iRecAft].status() != 52 && e[iRecAft].status() != -53) ||
435  e[iEmt].status() != 51 || e[iRadAft].status() != 51) {
436  e.list();
437  fatalEmissionVeto(std::string("Couldn't find Pythia FSR emission"));
438  }
439 
440  // Behaviour based on pTemtMode:
441  // 0 - pT of emitted w.r.t. radiator before
442  // 1 - std::min(pT of emitted w.r.t. all incoming/outgoing)
443  // 2 - std::min(pT of all outgoing w.r.t. all incoming/outgoing)
444  int xSR = (pTemtMode == 0) ? 1 : -1;
445  int i = (pTemtMode == 0) ? iRadBef : -1;
446  int k = (pTemtMode == 0) ? iRadAft : -1;
447  int r = (pTemtMode == 0) ? iRecAft : -1;
448 
449  // When pTemtMode is 0 or 1, iEmt has been selected
450  double pTemt = -1.;
451  if (pTemtMode == 0 || pTemtMode == 1) {
452  // Which parton is emitted, based on emittedMode:
453  // 0 - Pythia definition of emitted
454  // 1 - Pythia definition of radiated after emission
455  // 2 - Random selection of emitted or radiated after emission
456  // 3 - Try both emitted and radiated after emission
457  int j = iRadAft;
458  if (emittedMode == 0 || (emittedMode == 2 && rndmPtr->flat() < 0.5)) j++;
459 
460  for (int jLoop = 0; jLoop < 2; jLoop++) {
461  if (jLoop == 0) pTemt = pTcalc(e, i, j, k, r, xSR);
462  else if (jLoop == 1) pTemt = std::min(pTemt, pTcalc(e, i, j, k, r, xSR));
463 
464  // For emittedMode == 3, have tried iRadAft, now try iEmt
465  if (emittedMode != 3) break;
466  if (k != -1) std::swap(j, k); else j = iEmt;
467  }
468 
469  // If pTemtMode is 2, then try all final-state partons as emitted
470  } else if (pTemtMode == 2) {
471  pTemt = pTcalc(e, i, -1, k, r, xSR);
472 
473  }
474 
475 #ifdef DBGOUTPUT
476  std::cout << "doVetoFSREmission: pTemt = " << pTemt << std::endl << std::endl;
477 #endif
478 
479  // If a Born configuration, and a colorless FS, and QEDvetoMode=2,
480  // then don't veto photons, W, or Z harder than pThard
481  bool vetoParton = (!isEmt && e[iEmt].colType()==0 && QEDvetoMode==2)
482  ? false: true;
483 
484  // Veto if pTemt > pThard
485  if (pTemt > pThard) {
486  if(!vetoParton) {
487  // Don't veto ANY emissions afterwards
488  nAcceptSeq = vetoCount-1;
489  } else {
490  nAcceptSeq = 0;
491  nFSRveto++;
492  return true;
493  }
494  }
495 
496  // Else mark that an emission has been accepted and continue
497  nAcceptSeq++;
498  accepted = true;
499  return false;
500 }
void swap(edm::DataFrameContainer &lhs, edm::DataFrameContainer &rhs)
T min(T a, T b)
Definition: MathUtil.h:58
int k[5][pyjets_maxn]
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)
bool EmissionVetoHook1::doVetoISREmission ( int  ,
const Pythia8::Event &  e,
int  iSys 
)
override

Definition at line 357 of file EmissionVetoHook1.cc.

References gather_cfg::cout, mps_fire::i, gen::k, alignCSCRings::r, mps_update::status, and AlCaHLTBitMon_QueryRunRegistry::string.

Referenced by canVetoISREmission().

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

Definition at line 506 of file EmissionVetoHook1.cc.

References gather_cfg::cout.

Referenced by canVetoMPIEmission().

506  {
507  if (MPIvetoOn) {
508  if (e[e.size() - 1].pT() > pTMPI) {
509 #ifdef DBGOUTPUT
510  std::cout << "doVetoMPIEmission: pTnow = " << e[e.size() - 1].pT()
511  << ", pTMPI = " << pTMPI << std::endl << std::endl;
512 #endif
513  return true;
514  }
515  }
516  return false;
517 }
bool EmissionVetoHook1::doVetoMPIStep ( int  nMPI,
const Pythia8::Event &  e 
)
override

Definition at line 255 of file EmissionVetoHook1.cc.

References funct::abs(), KineDebug3::count(), gather_cfg::cout, plotBeamSpotDB::first, mps_fire::i, plotBeamSpotDB::last, and AlCaHLTBitMon_QueryRunRegistry::string.

Referenced by numberVetoMPIStep().

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

Definition at line 8 of file EmissionVetoHook1.cc.

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

Referenced by canVetoMPIEmission().

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

Definition at line 27 of file EmissionVetoHook1.h.

References doVetoMPIStep(), and MillePedeFileConverter_cfg::e.

27 { return 1; }
double EmissionVetoHook1::pTcalc ( const Pythia8::Event &  e,
int  i,
int  j,
int  k,
int  r,
int  xSRin 
)

Definition at line 121 of file EmissionVetoHook1.cc.

References gather_cfg::cout, and min().

Referenced by canVetoMPIEmission().

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

Definition at line 78 of file EmissionVetoHook1.cc.

References gather_cfg::cout, MillePedeFileConverter_cfg::e, AlCaHLTBitMon_ParallelJobs::p, and mathSSE::sqrt().

Referenced by canVetoMPIEmission().

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

Definition at line 21 of file EmissionVetoHook1.cc.

References funct::abs(), gather_cfg::cout, class-composition::Q, Validation_hcalonly_cfi::sign, mathSSE::sqrt(), and globals_cff::x1.

Referenced by canVetoMPIEmission().

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

Member Data Documentation

bool EmissionVetoHook1::accepted
private

Definition at line 55 of file EmissionVetoHook1.h.

int EmissionVetoHook1::emittedMode
private

Definition at line 51 of file EmissionVetoHook1.h.

bool EmissionVetoHook1::isEmt
private

Definition at line 55 of file EmissionVetoHook1.h.

int EmissionVetoHook1::MPIvetoOn
private

Definition at line 51 of file EmissionVetoHook1.h.

Referenced by canVetoMPIEmission().

int EmissionVetoHook1::nAcceptSeq
private

Definition at line 57 of file EmissionVetoHook1.h.

int EmissionVetoHook1::nFinal
private

Definition at line 53 of file EmissionVetoHook1.h.

int EmissionVetoHook1::nFinalExt
private

Definition at line 51 of file EmissionVetoHook1.h.

int EmissionVetoHook1::nFinalMode
private

Definition at line 51 of file EmissionVetoHook1.h.

unsigned long int EmissionVetoHook1::nFSRveto
private

Definition at line 59 of file EmissionVetoHook1.h.

Referenced by ~EmissionVetoHook1().

unsigned long int EmissionVetoHook1::nISRveto
private

Definition at line 59 of file EmissionVetoHook1.h.

Referenced by ~EmissionVetoHook1().

int EmissionVetoHook1::pTdefMode
private

Definition at line 51 of file EmissionVetoHook1.h.

int EmissionVetoHook1::pTemtMode
private

Definition at line 51 of file EmissionVetoHook1.h.

double EmissionVetoHook1::pThard
private

Definition at line 54 of file EmissionVetoHook1.h.

int EmissionVetoHook1::pThardMode
private

Definition at line 51 of file EmissionVetoHook1.h.

double EmissionVetoHook1::pTMPI
private

Definition at line 54 of file EmissionVetoHook1.h.

int EmissionVetoHook1::QEDvetoMode
private

Definition at line 51 of file EmissionVetoHook1.h.

int EmissionVetoHook1::Verbosity
private

Definition at line 60 of file EmissionVetoHook1.h.

int EmissionVetoHook1::vetoCount
private

Definition at line 51 of file EmissionVetoHook1.h.

int EmissionVetoHook1::vetoOn
private

Definition at line 51 of file EmissionVetoHook1.h.

Referenced by canVetoFSREmission(), and canVetoISREmission().