PYTHIA  8.312
UserHooks.h
1 // UserHooks.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2024 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 
6 // Header file to allow user access to program at different stages.
7 // UserHooks: almost empty base class, with user to write the rela code.
8 // MyUserHooks: derived class, only intended as an example.
9 
10 #ifndef Pythia8_UserHooks_H
11 #define Pythia8_UserHooks_H
12 
13 #include "Pythia8/Event.h"
14 #include "Pythia8/PartonSystems.h"
15 #include "Pythia8/PhysicsBase.h"
16 #include "Pythia8/PythiaStdlib.h"
17 #include "Pythia8/SigmaProcess.h"
18 
19 namespace Pythia8 {
20 
21 //==========================================================================
22 
23 // Forward references.
24 class PhaseSpace;
25 class StringEnd;
26 class HadronLevel;
27 
28 //==========================================================================
29 
30 // UserHooks is base class for user access to program execution.
31 
32 class UserHooks : public PhysicsBase {
33 
34 public:
35 
36  // Destructor.
37  virtual ~UserHooks() {}
38 
39  // Initialisation after beams have been set by Pythia::init().
40  virtual bool initAfterBeams() { return true; }
41 
42  // Possibility to modify cross section of process.
43  virtual bool canModifySigma() {return false;}
44 
45  // Multiplicative factor modifying the cross section of a hard process.
46  virtual double multiplySigmaBy(const SigmaProcess* sigmaProcessPtr,
47  const PhaseSpace* phaseSpacePtr, bool inEvent);
48 
49  // Possibility to bias selection of events, compensated by a weight.
50  virtual bool canBiasSelection() {return false;}
51 
52  // Multiplicative factor in the phase space selection of a hard process.
53  virtual double biasSelectionBy(const SigmaProcess* sigmaProcessPtr,
54  const PhaseSpace* phaseSpacePtr, bool inEvent);
55 
56  // Event weight to compensate for selection weight above.
57  virtual double biasedSelectionWeight() {return 1./selBias;}
58 
59  // Possibility to veto event after process-level selection.
60  virtual bool canVetoProcessLevel() {return false;}
61 
62  // Decide whether to veto current process or not, based on process record.
63  // Usage: doVetoProcessLevel( process).
64  virtual bool doVetoProcessLevel(Event& ) {return false;}
65 
66  // Possibility to set low energy cross sections.
67  // Usage: canSetLowEnergySigma(idA, idB).
68  virtual bool canSetLowEnergySigma(int, int) const {return false;}
69 
70  // Set low energy cross sections for ids where canSetLowEnergySigma is true.
71  // Usage: doSetLowEnergySigma(idA, idB, eCM, mA, mB).
72  virtual double doSetLowEnergySigma(int, int, double, double, double) const {
73  return 0.;}
74 
75  // Possibility to veto resonance decay chain.
76  virtual bool canVetoResonanceDecays() {return false;}
77 
78  // Decide whether to veto current resonance decay chain or not, based on
79  // process record. Usage: doVetoProcessLevel( process).
80  virtual bool doVetoResonanceDecays(Event& ) {return false;}
81 
82  // Possibility to veto MPI + ISR + FSR evolution and kill event,
83  // making decision at a fixed pT scale. Useful for MLM-style matching.
84  virtual bool canVetoPT() {return false;}
85 
86  // Transverse-momentum scale for veto test.
87  virtual double scaleVetoPT() {return 0.;}
88 
89  // Decide whether to veto current event or not, based on event record.
90  // Usage: doVetoPT( iPos, event), where iPos = 0: no emissions so far;
91  // iPos = 1/2/3 joint evolution, latest step was MPI/ISR/FSR;
92  // iPos = 4: FSR only afterwards; iPos = 5: FSR in resonance decay.
93  virtual bool doVetoPT( int , const Event& ) {return false;}
94 
95  // Possibility to veto MPI + ISR + FSR evolution and kill event,
96  // making decision after fixed number of ISR or FSR steps.
97  virtual bool canVetoStep() {return false;}
98 
99  // Up to how many ISR + FSR steps of hardest interaction should be checked.
100  virtual int numberVetoStep() {return 1;}
101 
102  // Decide whether to veto current event or not, based on event record.
103  // Usage: doVetoStep( iPos, nISR, nFSR, event), where iPos as above,
104  // nISR and nFSR number of emissions so far for hard interaction only.
105  virtual bool doVetoStep( int , int , int , const Event& ) {return false;}
106 
107  // Possibility to veto MPI + ISR + FSR evolution and kill event,
108  // making decision after fixed number of MPI steps.
109  virtual bool canVetoMPIStep() {return false;}
110 
111  // Up to how many MPI steps should be checked.
112  virtual int numberVetoMPIStep() {return 1;}
113 
114  // Decide whether to veto current event or not, based on event record.
115  // Usage: doVetoMPIStep( nMPI, event), where nMPI is number of MPI's so far.
116  virtual bool doVetoMPIStep( int , const Event& ) {return false;}
117 
118  // Possibility to veto event after ISR + FSR + MPI in parton level,
119  // but before beam remnants and resonance decays.
120  virtual bool canVetoPartonLevelEarly() {return false;}
121 
122  // Decide whether to veto current partons or not, based on event record.
123  // Usage: doVetoPartonLevelEarly( event).
124  virtual bool doVetoPartonLevelEarly( const Event& ) {return false;}
125 
126  // Retry same ProcessLevel with a new PartonLevel after a veto in
127  // doVetoPT, doVetoStep, doVetoMPIStep or doVetoPartonLevelEarly
128  // if you overload this method to return true.
129  virtual bool retryPartonLevel() {return false;}
130 
131  // Possibility to veto event after parton-level selection.
132  virtual bool canVetoPartonLevel() {return false;}
133 
134  // Decide whether to veto current partons or not, based on event record.
135  // Usage: doVetoPartonLevel( event).
136  virtual bool doVetoPartonLevel( const Event& ) {return false;}
137 
138  // Possibility to set initial scale in TimeShower for resonance decay.
139  virtual bool canSetResonanceScale() {return false;}
140 
141  // Initial scale for TimeShower evolution.
142  // Usage: scaleResonance( iRes, event), where iRes is location
143  // of decaying resonance in the event record.
144  virtual double scaleResonance( int, const Event& ) {return 0.;}
145 
146  // Possibility to veto an emission in the ISR machinery.
147  virtual bool canVetoISREmission() {return false;}
148 
149  // Decide whether to veto current emission or not, based on event record.
150  // Usage: doVetoISREmission( sizeOld, event, iSys) where sizeOld is size
151  // of event record before current emission-to-be-scrutinized was added,
152  // and iSys is the system of the radiation (according to PartonSystems).
153  virtual bool doVetoISREmission( int, const Event&, int ) {return false;}
154 
155  // Possibility to veto an emission in the FSR machinery.
156  virtual bool canVetoFSREmission() {return false;}
157 
158  // Decide whether to veto current emission or not, based on event record.
159  // Usage: doVetoFSREmission( sizeOld, event, iSys, inResonance) where
160  // sizeOld is size of event record before current emission-to-be-scrutinized
161  // was added, iSys is the system of the radiation (according to
162  // PartonSystems), and inResonance is true if the emission takes place in a
163  // resonance decay.
164  virtual bool doVetoFSREmission( int, const Event&, int, bool = false )
165  {return false;}
166 
167  // Possibility to veto an MPI.
168  virtual bool canVetoMPIEmission() { return false; }
169 
170  // Decide whether to veto an MPI based on event record.
171  // Usage: doVetoMPIEmission( sizeOld, event) where sizeOld
172  // is size of event record before the current MPI.
173  virtual bool doVetoMPIEmission( int, const Event &) { return false; }
174 
175  // Possibility to reconnect colours from resonance decay systems.
176  virtual bool canReconnectResonanceSystems() { return false; }
177 
178  // Do reconnect colours from resonance decay systems.
179  // Usage: doVetoFSREmission( oldSizeEvt, event)
180  // where oldSizeEvent is the event size before resonance decays.
181  // Should normally return true, while false means serious failure.
182  // Value of PartonLevel:earlyResDec determines where method is called.
183  virtual bool doReconnectResonanceSystems( int, Event &) {return true;}
184 
185  // Can change fragmentation parameters.
186  virtual bool canChangeFragPar() { return false;}
187 
188  // Set initial ends of a string to be fragmented. This is done once
189  // for each string. Note that the second string end may be zero in case
190  // we are hadronising a string piece leading to a junction.
191  virtual void setStringEnds( const StringEnd*, const StringEnd*,
192  vector<int>) {}
193 
194  // Do change fragmentation parameters.
195  // Input: flavPtr, zPtr, pTPtr, idEnd, m2Had, iParton and posEnd (or
196  // negEnd).
197  virtual bool doChangeFragPar( StringFlav*, StringZ*, StringPT*, int,
198  double, vector<int>, const StringEnd* ) { return false;}
199 
200  // Do a veto on a hadron just before it is added to the final state.
201  // The StringEnd from which the the hadron was produced is included
202  // for information.
203  virtual bool doVetoFragmentation( Particle, const StringEnd *)
204  { return false;}
205 
206  // Do a veto on a hadron just before it is added to the final state
207  // (final two hadron case).
209  const StringEnd*, const StringEnd* ) { return false;}
210 
211  // Possibility to veto an event after hadronization based
212  // on event contents. Works as an early trigger to avoid
213  // running the time consuming rescattering process on
214  // uninteresting events.
215  virtual bool canVetoAfterHadronization() {return false;}
216 
217  // Do the actual veto after hadronization.
218  virtual bool doVetoAfterHadronization(const Event& ) {return false;}
219 
220  // Can set the overall impact parameter for the MPI treatment.
221  virtual bool canSetImpactParameter() const { return false; }
222 
223  // Set the overall impact parameter for the MPI treatment.
224  virtual double doSetImpactParameter() { return 0.0; }
225 
226  // Custom processing at the end of HadronLevel::next.
227  virtual bool onEndHadronLevel(HadronLevel&, Event&) { return true; }
228 
229 protected:
230 
231  // Constructor.
233 
234  // After initInfoPtr, initialize workEvent
235  virtual void onInitInfoPtr() override {
236  // Set smart pointer to null, in order to avoid circular dependency
237  userHooksPtr = nullptr;
238  workEvent.init("(work event)", particleDataPtr);
239  }
240 
241  // omitResonanceDecays omits resonance decay chains from process record.
242  void omitResonanceDecays(const Event& process, bool finalOnly = false);
243 
244  // subEvent extracts currently resolved partons in the hard process.
245  void subEvent(const Event& event, bool isHardest = true);
246 
247  // Have one event object around as work area.
249 
250  // User-imposed selection bias.
251  double selBias = 1.;
252 
253  // Bookkept quantities for boosted event weights.
254  double enhancedEventWeight = {}, pTEnhanced = {}, wtEnhanced = {};
255 
256 };
257 
258 //==========================================================================
259 
260 // SuppressSmallPT is a derived class for user access to program execution.
261 // It is a simple example, illustrating how to suppress the cross section
262 // of 2 -> 2 processes by a factor pT^4 / (pT0^2 + pT^2)^2, with pT0 input,
263 // and also modify alpha_strong scale similarly.
264 
265 class SuppressSmallPT : public UserHooks {
266 
267 public:
268 
269  // Constructor.
270  SuppressSmallPT( double pT0timesMPIIn = 1., int numberAlphaSIn = 0,
271  bool useSameAlphaSasMPIIn = true) : pT20(0.) {isInit = false;
272  pT0timesMPI = pT0timesMPIIn; numberAlphaS = numberAlphaSIn;
273  useSameAlphaSasMPI = useSameAlphaSasMPIIn;}
274 
275  // Possibility to modify cross section of process.
276  virtual bool canModifySigma() {return true;}
277 
278  // Multiplicative factor modifying the cross section of a hard process.
279  // Usage: inEvent is true for event generation, false for initialization.
280  virtual double multiplySigmaBy(const SigmaProcess* sigmaProcessPtr,
281  const PhaseSpace* phaseSpacePtr, bool );
282 
283 private:
284 
285  // Save input properties and the squared pT0 scale.
286  bool isInit, useSameAlphaSasMPI;
287  int numberAlphaS;
288  double pT0timesMPI, pT20;
289 
290  // Alpha_strong calculation.
291  AlphaStrong alphaS;
292 
293 };
294 
295 //==========================================================================
296 
297 // UserHooksVector implements a vector of UserHooks and is itself a UserHooks.
298 
299 class UserHooksVector: public UserHooks {
300 
301 public:
302 
303  // The default constructor is private, and should only be used
304  // internally in Pythia.
306  friend class Pythia;
307 
308  // Destructor.
309  virtual ~UserHooksVector() {}
310 
311  // Initialisation after beams have been set by Pythia::init().
312  // Check that there are no (obvious) clashes.
313  virtual bool initAfterBeams() {
314  int nCanSetResonanceScale = 0;
315  int nCanChangeFragPar = 0;
316  int nCanSetImpactParameter = 0;
317  for ( int i = 0, N = hooks.size(); i < N; ++i ) {
318  registerSubObject(*hooks[i]);
319  if ( !hooks[i]->initAfterBeams() ) return false;
320  if (hooks[i]->canSetResonanceScale()) ++nCanSetResonanceScale;
321  if (hooks[i]->canChangeFragPar()) ++nCanChangeFragPar;
322  if (hooks[i]->canSetImpactParameter()) ++nCanSetImpactParameter;
323  }
324  if (nCanSetResonanceScale > 1) {
325  loggerPtr->ERROR_MSG(
326  "multiple UserHooks with canSetResonanceScale() not allowed");
327  return false;
328  }
329  if (nCanChangeFragPar > 1) {
330  loggerPtr->ERROR_MSG(
331  "multiple UserHooks with canChangeFragPar() not allowed");
332  return false;
333  }
334  if (nCanSetImpactParameter > 1) {
335  loggerPtr->ERROR_MSG(
336  "multiple UserHooks with canSetImpactParameter() not allowed");
337  return false;
338  }
339  return true;
340  }
341 
342  // Possibility to modify cross section of process.
343  virtual bool canModifySigma() {
344  for ( int i = 0, N = hooks.size(); i < N; ++i )
345  if ( hooks[i]->canModifySigma() ) return true;
346  return false;
347  }
348 
349  // Multiplicative factor modifying the cross section of a hard process.
350  virtual double multiplySigmaBy(const SigmaProcess* sigmaProcessPtr,
351  const PhaseSpace* phaseSpacePtr, bool inEvent) {
352  double f = 1.0;
353  for ( int i = 0, N = hooks.size(); i < N; ++i )
354  if ( hooks[i]->canModifySigma() )
355  f *= hooks[i]->multiplySigmaBy(sigmaProcessPtr, phaseSpacePtr,inEvent);
356  return f;
357  }
358 
359  // Possibility to bias selection of events, compensated by a weight.
360  virtual bool canBiasSelection() {
361  for ( int i = 0, N = hooks.size(); i < N; ++i )
362  if ( hooks[i]->canBiasSelection() ) return true;
363  return false;
364  }
365 
366  // Multiplicative factor in the phase space selection of a hard process.
367  virtual double biasSelectionBy(const SigmaProcess* sigmaProcessPtr,
368  const PhaseSpace* phaseSpacePtr, bool inEvent) {
369  double f = 1.0;
370  for ( int i = 0, N = hooks.size(); i < N; ++i )
371  if ( hooks[i]->canBiasSelection() )
372  f *= hooks[i]->biasSelectionBy(sigmaProcessPtr, phaseSpacePtr,
373  inEvent);
374  return f;
375  }
376 
377  // Event weight to compensate for selection weight above.
378  virtual double biasedSelectionWeight() {
379  double f = 1.0;
380  for ( int i = 0, N = hooks.size(); i < N; ++i )
381  if ( hooks[i]->canBiasSelection() )
382  f *= hooks[i]->biasedSelectionWeight();
383  return f;
384  }
385 
386  // Possibility to veto event after process-level selection.
387  virtual bool canVetoProcessLevel() {
388  for ( int i = 0, N = hooks.size(); i < N; ++i )
389  if ( hooks[i]->canVetoProcessLevel() ) return true;
390  return false;
391  }
392 
393  // Decide whether to veto current process or not, based on process record.
394  // Usage: doVetoProcessLevel( process).
395  virtual bool doVetoProcessLevel(Event& e) {
396  for ( int i = 0, N = hooks.size(); i < N; ++i )
397  if ( hooks[i]->canVetoProcessLevel() &&
398  hooks[i]->doVetoProcessLevel(e) ) return true;
399  return false;
400  }
401 
402  // Possibility to veto resonance decay chain.
403  virtual bool canVetoResonanceDecays() {
404  for ( int i = 0, N = hooks.size(); i < N; ++i )
405  if ( hooks[i]->canVetoResonanceDecays() ) return true;
406  return false;
407  }
408 
409  // Decide whether to veto current resonance decay chain or not, based on
410  // process record. Usage: doVetoProcessLevel( process).
411  virtual bool doVetoResonanceDecays(Event& e) {
412  for ( int i = 0, N = hooks.size(); i < N; ++i )
413  if ( hooks[i]->canVetoResonanceDecays() &&
414  hooks[i]->doVetoResonanceDecays(e) ) return true;
415  return false;
416  }
417 
418  // Possibility to veto MPI + ISR + FSR evolution and kill event,
419  // making decision at a fixed pT scale. Useful for MLM-style matching.
420  virtual bool canVetoPT() {
421  for ( int i = 0, N = hooks.size(); i < N; ++i )
422  if ( hooks[i]->canVetoPT() ) return true;
423  return false;
424  }
425 
426  // Transverse-momentum scale for veto test.
427  virtual double scaleVetoPT() {
428  double s = 0.0;
429  for ( int i = 0, N = hooks.size(); i < N; ++i )
430  if ( hooks[i]->canVetoPT() ) s = max(s, hooks[i]->scaleVetoPT());
431  return s;
432  }
433 
434  // Decide whether to veto current event or not, based on event record.
435  // Usage: doVetoPT( iPos, event), where iPos = 0: no emissions so far;
436  // iPos = 1/2/3 joint evolution, latest step was MPI/ISR/FSR;
437  // iPos = 4: FSR only afterwards; iPos = 5: FSR in resonance decay.
438  virtual bool doVetoPT( int iPos, const Event& e) {
439  for ( int i = 0, N = hooks.size(); i < N; ++i )
440  if ( hooks[i]->canVetoPT() && hooks[i]->doVetoPT(iPos, e) ) return true;
441  return false;
442  }
443 
444  // Possibility to veto MPI + ISR + FSR evolution and kill event,
445  // making decision after fixed number of ISR or FSR steps.
446  virtual bool canVetoStep() {
447  for ( int i = 0, N = hooks.size(); i < N; ++i )
448  if ( hooks[i]->canVetoStep() ) return true;
449  return false;
450  }
451 
452  // Up to how many ISR + FSR steps of hardest interaction should be checked.
453  virtual int numberVetoStep() {
454  int n = 1;
455  for ( int i = 0, N = hooks.size(); i < N; ++i )
456  if ( hooks[i]->canVetoStep() ) n = max(n, hooks[i]->numberVetoStep());
457  return n;
458  }
459 
460  // Decide whether to veto current event or not, based on event record.
461  // Usage: doVetoStep( iPos, nISR, nFSR, event), where iPos as above,
462  // nISR and nFSR number of emissions so far for hard interaction only.
463  virtual bool doVetoStep( int iPos, int nISR, int nFSR, const Event& e) {
464  for ( int i = 0, N = hooks.size(); i < N; ++i )
465  if ( hooks[i]->canVetoStep()
466  && hooks[i]->doVetoStep(iPos, nISR, nFSR, e) ) return true;
467  return false;
468  }
469 
470  // Possibility to veto MPI + ISR + FSR evolution and kill event,
471  // making decision after fixed number of MPI steps.
472  virtual bool canVetoMPIStep() {
473  for ( int i = 0, N = hooks.size(); i < N; ++i )
474  if ( hooks[i]->canVetoMPIStep() ) return true;
475  return false;
476  }
477 
478  // Up to how many MPI steps should be checked.
479  virtual int numberVetoMPIStep() {
480  int n = 1;
481  for ( int i = 0, N = hooks.size(); i < N; ++i )
482  if ( hooks[i]->canVetoMPIStep() )
483  n = max(n, hooks[i]->numberVetoMPIStep());
484  return n;
485  }
486 
487  // Decide whether to veto current event or not, based on event record.
488  // Usage: doVetoMPIStep( nMPI, event), where nMPI is number of MPI's so far.
489  virtual bool doVetoMPIStep( int nMPI, const Event& e) {
490  for ( int i = 0, N = hooks.size(); i < N; ++i )
491  if ( hooks[i]->canVetoMPIStep() && hooks[i]->doVetoMPIStep(nMPI, e) )
492  return true;
493  return false;
494  }
495 
496  // Possibility to veto event after ISR + FSR + MPI in parton level,
497  // but before beam remnants and resonance decays.
498  virtual bool canVetoPartonLevelEarly() {
499  for ( int i = 0, N = hooks.size(); i < N; ++i )
500  if ( hooks[i]->canVetoPartonLevelEarly() ) return true;
501  return false;
502  }
503 
504  // Decide whether to veto current partons or not, based on event record.
505  // Usage: doVetoPartonLevelEarly( event).
506  virtual bool doVetoPartonLevelEarly( const Event& e) {
507  for ( int i = 0, N = hooks.size(); i < N; ++i )
508  if ( hooks[i]->canVetoPartonLevelEarly()
509  && hooks[i]->doVetoPartonLevelEarly(e) ) return true;
510  return false;
511  }
512 
513  // Retry same ProcessLevel with a new PartonLevel after a veto in
514  // doVetoPT, doVetoStep, doVetoMPIStep or doVetoPartonLevelEarly
515  // if you overload this method to return true.
516  virtual bool retryPartonLevel() {
517  for ( int i = 0, N = hooks.size(); i < N; ++i )
518  if ( hooks[i]->retryPartonLevel() ) return true;
519  return false;
520  }
521 
522  // Possibility to veto event after parton-level selection.
523  virtual bool canVetoPartonLevel() {
524  for ( int i = 0, N = hooks.size(); i < N; ++i )
525  if ( hooks[i]->canVetoPartonLevel() ) return true;
526  return false;
527  }
528 
529  // Decide whether to veto current partons or not, based on event record.
530  // Usage: doVetoPartonLevel( event).
531  virtual bool doVetoPartonLevel( const Event& e) {
532  for ( int i = 0, N = hooks.size(); i < N; ++i )
533  if ( hooks[i]->canVetoPartonLevel()
534  && hooks[i]->doVetoPartonLevel(e) ) return true;
535  return false;
536  }
537 
538  // Possibility to set initial scale in TimeShower for resonance decay.
539  virtual bool canSetResonanceScale() {
540  for ( int i = 0, N = hooks.size(); i < N; ++i )
541  if ( hooks[i]->canSetResonanceScale() ) return true;
542  return false;
543  }
544 
545  // Initial scale for TimeShower evolution.
546  // Usage: scaleResonance( iRes, event), where iRes is location
547  // of decaying resonance in the event record.
548  virtual double scaleResonance( int iRes, const Event& e) {
549  double s = 0.0;
550  for ( int i = 0, N = hooks.size(); i < N; ++i )
551  if ( hooks[i]->canSetResonanceScale() )
552  s = max(s, hooks[i]->scaleResonance(iRes, e));
553  return s;
554  }
555 
556  // Possibility to veto an emission in the ISR machinery.
557  virtual bool canVetoISREmission() {
558  for ( int i = 0, N = hooks.size(); i < N; ++i )
559  if ( hooks[i]->canVetoISREmission() ) return true;
560  return false;
561  }
562 
563  // Decide whether to veto current emission or not, based on event record.
564  // Usage: doVetoISREmission( sizeOld, event, iSys) where sizeOld is size
565  // of event record before current emission-to-be-scrutinized was added,
566  // and iSys is the system of the radiation (according to PartonSystems).
567  virtual bool doVetoISREmission( int sizeOld, const Event& e, int iSys) {
568  for ( int i = 0, N = hooks.size(); i < N; ++i )
569  if ( hooks[i]->canVetoISREmission()
570  && hooks[i]->doVetoISREmission(sizeOld, e, iSys) ) return true;
571  return false;
572  }
573 
574  // Possibility to veto an emission in the FSR machinery.
575  virtual bool canVetoFSREmission() {
576  for ( int i = 0, N = hooks.size(); i < N; ++i )
577  if ( hooks[i]->canVetoFSREmission() ) return true;
578  return false;
579  }
580 
581  // Decide whether to veto current emission or not, based on event record.
582  // Usage: doVetoFSREmission( sizeOld, event, iSys, inResonance) where
583  // sizeOld is size of event record before current emission-to-be-scrutinized
584  // was added, iSys is the system of the radiation (according to
585  // PartonSystems), and inResonance is true if the emission takes place in a
586  // resonance decay.
587  virtual bool doVetoFSREmission(int sizeOld, const Event& e,
588  int iSys, bool inResonance = false ) {
589  for ( int i = 0, N = hooks.size(); i < N; ++i )
590  if ( hooks[i]->canVetoFSREmission()
591  && hooks[i]->doVetoFSREmission(sizeOld, e, iSys, inResonance) )
592  return true;
593  return false;
594  }
595 
596  // Possibility to veto an MPI.
597  virtual bool canVetoMPIEmission() {
598  for ( int i = 0, N = hooks.size(); i < N; ++i )
599  if ( hooks[i]->canVetoMPIEmission() ) return true;
600  return false;
601  }
602 
603  // Decide whether to veto an MPI based on event record.
604  // Usage: doVetoMPIEmission( sizeOld, event) where sizeOld
605  // is size of event record before the current MPI.
606  virtual bool doVetoMPIEmission( int sizeOld, const Event & e) {
607  for ( int i = 0, N = hooks.size(); i < N; ++i )
608  if ( hooks[i]->canVetoMPIEmission()
609  && hooks[i]->doVetoMPIEmission(sizeOld, e) )
610  return true;
611  return false;
612  }
613 
614  // Possibility to reconnect colours from resonance decay systems.
616  for ( int i = 0, N = hooks.size(); i < N; ++i )
617  if ( hooks[i]->canReconnectResonanceSystems() ) return true;
618  return false;
619  }
620 
621  // Do reconnect colours from resonance decay systems.
622  // Usage: doVetoFSREmission( oldSizeEvt, event)
623  // where oldSizeEvent is the event size before resonance decays.
624  // Should normally return true, while false means serious failure.
625  // Value of PartonLevel:earlyResDec determines where method is called.
626  virtual bool doReconnectResonanceSystems( int j, Event & e) {
627  for ( int i = 0, N = hooks.size(); i < N; ++i )
628  if ( hooks[i]->canReconnectResonanceSystems()
629  && hooks[i]->doReconnectResonanceSystems(j, e) ) return true;
630  return false;
631  }
632 
633  // Can change fragmentation parameters.
634  virtual bool canChangeFragPar() {
635  for ( int i = 0, N = hooks.size(); i < N; ++i )
636  if ( hooks[i]->canChangeFragPar() ) return true;
637  return false;
638  }
639 
640  // Set initial ends of a string to be fragmented. This is done once
641  // for each string. Note that the second string end may be zero in
642  // case we are hadronising a string piece leading to a junction.
643  virtual void setStringEnds( const StringEnd* pos, const StringEnd* neg,
644  vector<int> iPart) {
645  for ( int i = 0, N = hooks.size(); i < N; ++i )
646  hooks[i]->setStringEnds( pos, neg, iPart);
647  }
648 
649  // Do change fragmentation parameters.
650  // Input: flavPtr, zPtr, pTPtr, idEnd, m2Had, iParton and posEnd (or
651  // negEnd).
652  virtual bool doChangeFragPar( StringFlav* sfIn, StringZ* zIn,
653  StringPT* ptIn, int idIn, double mIn, vector<int> parIn,
654  const StringEnd* endIn) {
655  for ( int i = 0, N = hooks.size(); i < N; ++i )
656  if ( hooks[i]->canChangeFragPar()
657  && hooks[i]->doChangeFragPar(sfIn, zIn, ptIn, idIn,
658  mIn, parIn, endIn) ) return true;
659  return false;}
660 
661  // Do a veto on a hadron just before it is added to the final state.
662  virtual bool doVetoFragmentation(Particle p, const StringEnd* nowEnd) {
663  for ( int i = 0, N = hooks.size(); i < N; ++i )
664  if ( hooks[i]->canChangeFragPar()
665  && hooks[i]->doVetoFragmentation(p, nowEnd) ) return true;
666  return false;
667  }
668 
669  virtual bool doVetoFragmentation(Particle p1, Particle p2,
670  const StringEnd* e1, const StringEnd* e2) {
671  for ( int i = 0, N = hooks.size(); i < N; ++i )
672  if ( hooks[i]->canChangeFragPar()
673  && hooks[i]->doVetoFragmentation(p1, p2, e1, e2) ) return true;
674  return false;
675  }
676 
677  // Possibility to veto an event after hadronization based on event
678  // contents. Works as an early trigger to avoid running the time
679  // consuming rescattering process on uninteresting events.
680  virtual bool canVetoAfterHadronization() {
681  for ( int i = 0, N = hooks.size(); i < N; ++i )
682  if ( hooks[i]->canVetoAfterHadronization() ) return true;
683  return false;
684  }
685 
686  // Do the actual veto after hadronization.
687  virtual bool doVetoAfterHadronization(const Event& e) {
688  for ( int i = 0, N = hooks.size(); i < N; ++i )
689  if ( hooks[i]->canVetoAfterHadronization()
690  && hooks[i]->doVetoAfterHadronization(e) ) return true;
691  return false;
692  }
693 
694  // Can set the overall impact parameter for the MPI treatment.
695  virtual bool canSetImpactParameter() const {
696  for ( int i = 0, N = hooks.size(); i < N; ++i )
697  if ( hooks[i]->canSetImpactParameter() ) return true;
698  return false;
699  }
700 
701  // Set the overall impact parameter for the MPI treatment.
702  virtual double doSetImpactParameter() {
703  for ( int i = 0, N = hooks.size(); i < N; ++i )
704  if ( hooks[i]->canSetImpactParameter() )
705  return hooks[i]->doSetImpactParameter();
706  return 0.0;
707  }
708 
709  // The vector of user hooks.
710  vector< shared_ptr<UserHooks> > hooks = {};
711 
712 };
713 
714 //==========================================================================
715 
716 } // end namespace Pythia8
717 
718 #endif // Pythia8_UserHooks_H
virtual bool canBiasSelection()
Possibility to bias selection of events, compensated by a weight.
Definition: UserHooks.h:50
virtual int numberVetoStep()
Up to how many ISR + FSR steps of hardest interaction should be checked.
Definition: UserHooks.h:453
virtual double scaleVetoPT()
Transverse-momentum scale for veto test.
Definition: UserHooks.h:87
virtual int numberVetoStep()
Up to how many ISR + FSR steps of hardest interaction should be checked.
Definition: UserHooks.h:100
virtual bool canVetoPT()
Definition: UserHooks.h:420
virtual double scaleVetoPT()
Transverse-momentum scale for veto test.
Definition: UserHooks.h:427
virtual bool canVetoPartonLevelEarly()
Definition: UserHooks.h:120
virtual double multiplySigmaBy(const SigmaProcess *sigmaProcessPtr, const PhaseSpace *phaseSpacePtr, bool inEvent)
Multiplicative factor modifying the cross section of a hard process.
Definition: UserHooks.cc:27
Definition: StandardModel.h:23
UserHooksVector()
Definition: UserHooks.h:305
Definition: PhysicsBase.h:27
virtual bool canModifySigma()
Possibility to modify cross section of process.
Definition: UserHooks.h:276
void registerSubObject(PhysicsBase &pb)
Register a sub object that should have its information in sync with this.
Definition: PhysicsBase.cc:56
virtual double doSetImpactParameter()
Set the overall impact parameter for the MPI treatment.
Definition: UserHooks.h:702
The Event class holds all info on the generated event.
Definition: Event.h:453
UserHooks()
Constructor.
Definition: UserHooks.h:232
virtual bool doVetoISREmission(int sizeOld, const Event &e, int iSys)
Definition: UserHooks.h:567
virtual bool retryPartonLevel()
Definition: UserHooks.h:129
virtual bool onEndHadronLevel(HadronLevel &, Event &)
Custom processing at the end of HadronLevel::next.
Definition: UserHooks.h:227
SigmaProcess is the base class for cross section calculations.
Definition: SigmaProcess.h:86
virtual bool doVetoPT(int iPos, const Event &e)
Definition: UserHooks.h:438
virtual bool doVetoPartonLevel(const Event &e)
Definition: UserHooks.h:531
virtual bool doVetoMPIStep(int nMPI, const Event &e)
Definition: UserHooks.h:489
virtual double biasedSelectionWeight()
Event weight to compensate for selection weight above.
Definition: UserHooks.h:378
virtual bool doChangeFragPar(StringFlav *, StringZ *, StringPT *, int, double, vector< int >, const StringEnd *)
Definition: UserHooks.h:197
virtual bool canVetoProcessLevel()
Possibility to veto event after process-level selection.
Definition: UserHooks.h:387
virtual void setStringEnds(const StringEnd *, const StringEnd *, vector< int >)
Definition: UserHooks.h:191
The StringPT class is used to select select transverse momenta.
Definition: FragmentationFlavZpT.h:322
virtual bool canVetoMPIEmission()
Possibility to veto an MPI.
Definition: UserHooks.h:597
virtual bool doVetoStep(int iPos, int nISR, int nFSR, const Event &e)
Definition: UserHooks.h:463
virtual bool canVetoStep()
Definition: UserHooks.h:446
virtual bool canVetoAfterHadronization()
Definition: UserHooks.h:680
virtual bool canSetResonanceScale()
Possibility to set initial scale in TimeShower for resonance decay.
Definition: UserHooks.h:139
virtual bool doVetoMPIStep(int, const Event &)
Definition: UserHooks.h:116
virtual bool doVetoFSREmission(int, const Event &, int, bool=false)
Definition: UserHooks.h:164
virtual bool doVetoPartonLevelEarly(const Event &)
Definition: UserHooks.h:124
virtual ~UserHooks()
Destructor.
Definition: UserHooks.h:37
virtual int numberVetoMPIStep()
Up to how many MPI steps should be checked.
Definition: UserHooks.h:112
virtual void setStringEnds(const StringEnd *pos, const StringEnd *neg, vector< int > iPart)
Definition: UserHooks.h:643
virtual bool doVetoFSREmission(int sizeOld, const Event &e, int iSys, bool inResonance=false)
Definition: UserHooks.h:587
virtual bool doVetoResonanceDecays(Event &)
Definition: UserHooks.h:80
virtual bool doChangeFragPar(StringFlav *sfIn, StringZ *zIn, StringPT *ptIn, int idIn, double mIn, vector< int > parIn, const StringEnd *endIn)
Definition: UserHooks.h:652
virtual bool doVetoResonanceDecays(Event &e)
Definition: UserHooks.h:411
Definition: UserHooks.h:265
virtual double biasedSelectionWeight()
Event weight to compensate for selection weight above.
Definition: UserHooks.h:57
virtual bool doVetoPT(int, const Event &)
Definition: UserHooks.h:93
ParticleData * particleDataPtr
Pointer to the particle data table.
Definition: PhysicsBase.h:84
virtual bool doVetoAfterHadronization(const Event &e)
Do the actual veto after hadronization.
Definition: UserHooks.h:687
virtual bool initAfterBeams()
Definition: UserHooks.h:313
virtual bool canVetoISREmission()
Possibility to veto an emission in the ISR machinery.
Definition: UserHooks.h:147
The StringZ class is used to sample the fragmentation function f(z).
Definition: FragmentationFlavZpT.h:264
virtual bool canVetoResonanceDecays()
Possibility to veto resonance decay chain.
Definition: UserHooks.h:403
virtual double scaleResonance(int, const Event &)
Definition: UserHooks.h:144
virtual bool canVetoPartonLevel()
Possibility to veto event after parton-level selection.
Definition: UserHooks.h:523
Definition: HadronLevel.h:45
virtual bool canVetoMPIEmission()
Possibility to veto an MPI.
Definition: UserHooks.h:168
UserHooksPtr userHooksPtr
Definition: PhysicsBase.h:124
virtual void onInitInfoPtr() override
After initInfoPtr, initialize workEvent.
Definition: UserHooks.h:235
virtual bool doVetoFragmentation(Particle, const StringEnd *)
Definition: UserHooks.h:203
virtual bool canVetoResonanceDecays()
Possibility to veto resonance decay chain.
Definition: UserHooks.h:76
virtual bool doVetoProcessLevel(Event &e)
Definition: UserHooks.h:395
void omitResonanceDecays(const Event &process, bool finalOnly=false)
omitResonanceDecays omits resonance decay chains from process record.
Definition: UserHooks.cc:132
virtual bool canSetResonanceScale()
Possibility to set initial scale in TimeShower for resonance decay.
Definition: UserHooks.h:539
virtual bool doVetoPartonLevelEarly(const Event &e)
Definition: UserHooks.h:506
virtual bool canChangeFragPar()
Can change fragmentation parameters.
Definition: UserHooks.h:186
Definition: StringFragmentation.h:33
double selBias
User-imposed selection bias.
Definition: UserHooks.h:251
virtual bool canModifySigma()
Possibility to modify cross section of process.
Definition: UserHooks.h:343
virtual bool canModifySigma()
Possibility to modify cross section of process.
Definition: UserHooks.h:43
virtual bool doVetoFragmentation(Particle p, const StringEnd *nowEnd)
Do a veto on a hadron just before it is added to the final state.
Definition: UserHooks.h:662
virtual bool canReconnectResonanceSystems()
Possibility to reconnect colours from resonance decay systems.
Definition: UserHooks.h:176
virtual int numberVetoMPIStep()
Up to how many MPI steps should be checked.
Definition: UserHooks.h:479
SuppressSmallPT(double pT0timesMPIIn=1., int numberAlphaSIn=0, bool useSameAlphaSasMPIIn=true)
Constructor.
Definition: UserHooks.h:270
virtual bool doVetoFragmentation(Particle, Particle, const StringEnd *, const StringEnd *)
Definition: UserHooks.h:208
virtual bool canBiasSelection()
Possibility to bias selection of events, compensated by a weight.
Definition: UserHooks.h:360
virtual bool doReconnectResonanceSystems(int j, Event &e)
Definition: UserHooks.h:626
virtual double biasSelectionBy(const SigmaProcess *sigmaProcessPtr, const PhaseSpace *phaseSpacePtr, bool inEvent)
Multiplicative factor in the phase space selection of a hard process.
Definition: UserHooks.cc:79
virtual bool canSetLowEnergySigma(int, int) const
Definition: UserHooks.h:68
virtual ~UserHooksVector()
Destructor.
Definition: UserHooks.h:309
virtual bool canVetoProcessLevel()
Possibility to veto event after process-level selection.
Definition: UserHooks.h:60
virtual bool doVetoMPIEmission(int sizeOld, const Event &e)
Definition: UserHooks.h:606
virtual bool canSetImpactParameter() const
Can set the overall impact parameter for the MPI treatment.
Definition: UserHooks.h:221
Definition: Event.h:32
virtual bool canVetoStep()
Definition: UserHooks.h:97
virtual bool canVetoPartonLevel()
Possibility to veto event after parton-level selection.
Definition: UserHooks.h:132
virtual bool doReconnectResonanceSystems(int, Event &)
Definition: UserHooks.h:183
virtual bool canChangeFragPar()
Can change fragmentation parameters.
Definition: UserHooks.h:634
double enhancedEventWeight
Bookkept quantities for boosted event weights.
Definition: UserHooks.h:254
virtual bool canVetoMPIStep()
Definition: UserHooks.h:472
virtual bool canVetoFSREmission()
Possibility to veto an emission in the FSR machinery.
Definition: UserHooks.h:575
Logger * loggerPtr
Pointer to logger.
Definition: PhysicsBase.h:87
virtual bool doVetoISREmission(int, const Event &, int)
Definition: UserHooks.h:153
virtual bool retryPartonLevel()
Definition: UserHooks.h:516
virtual bool canSetImpactParameter() const
Can set the overall impact parameter for the MPI treatment.
Definition: UserHooks.h:695
virtual bool canVetoPT()
Definition: UserHooks.h:84
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
The Pythia class contains the top-level routines to generate an event.
Definition: Pythia.h:71
virtual bool doVetoMPIEmission(int, const Event &)
Definition: UserHooks.h:173
virtual bool canVetoAfterHadronization()
Definition: UserHooks.h:215
Event workEvent
Have one event object around as work area.
Definition: UserHooks.h:248
UserHooksVector implements a vector of UserHooks and is itself a UserHooks.
Definition: UserHooks.h:299
The StringFlav class is used to select quark and hadron flavours.
Definition: FragmentationFlavZpT.h:84
UserHooks is base class for user access to program execution.
Definition: UserHooks.h:32
void subEvent(const Event &event, bool isHardest=true)
subEvent extracts currently resolved partons in the hard process.
Definition: UserHooks.cc:182
Definition: PhaseSpace.h:42
virtual bool canVetoFSREmission()
Possibility to veto an emission in the FSR machinery.
Definition: UserHooks.h:156
void init(string headerIn="", ParticleData *particleDataPtrIn=0, int startColTagIn=100)
Initialize header for event listing, particle data table, and colour.
Definition: Event.h:467
virtual bool initAfterBeams()
Initialisation after beams have been set by Pythia::init().
Definition: UserHooks.h:40
virtual double biasSelectionBy(const SigmaProcess *sigmaProcessPtr, const PhaseSpace *phaseSpacePtr, bool inEvent)
Multiplicative factor in the phase space selection of a hard process.
Definition: UserHooks.h:367
virtual double multiplySigmaBy(const SigmaProcess *sigmaProcessPtr, const PhaseSpace *phaseSpacePtr, bool inEvent)
Multiplicative factor modifying the cross section of a hard process.
Definition: UserHooks.h:350
virtual bool canReconnectResonanceSystems()
Possibility to reconnect colours from resonance decay systems.
Definition: UserHooks.h:615
virtual bool doVetoProcessLevel(Event &)
Definition: UserHooks.h:64
virtual bool canVetoMPIStep()
Definition: UserHooks.h:109
virtual bool doVetoStep(int, int, int, const Event &)
Definition: UserHooks.h:105
virtual bool doVetoPartonLevel(const Event &)
Definition: UserHooks.h:136
virtual bool canVetoPartonLevelEarly()
Definition: UserHooks.h:498
virtual bool doVetoFragmentation(Particle p1, Particle p2, const StringEnd *e1, const StringEnd *e2)
Definition: UserHooks.h:669
virtual bool doVetoAfterHadronization(const Event &)
Do the actual veto after hadronization.
Definition: UserHooks.h:218
virtual bool canVetoISREmission()
Possibility to veto an emission in the ISR machinery.
Definition: UserHooks.h:557
virtual double scaleResonance(int iRes, const Event &e)
Definition: UserHooks.h:548
virtual double doSetLowEnergySigma(int, int, double, double, double) const
Definition: UserHooks.h:72
virtual double doSetImpactParameter()
Set the overall impact parameter for the MPI treatment.
Definition: UserHooks.h:224