herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
AlpGenHandlerOL.h
1 // -*- C++ -*-
2 #ifndef HERWIG_AlpGenHandlerOL_H
3 #define HERWIG_AlpGenHandlerOL_H
4 //
5 // This is the declaration of the AlpGenHandlerOL class.
6 //
7 
8 #include "Herwig/Shower/ShowerHandler.h"
9 #include "Herwig/Shower/QTilde/QTildeShowerHandler.h"
10 #include "ThePEG/Config/Pointers.h"
11 #include "Herwig/Shower/ShowerAlpha.h"
12 #include "fastjet/PseudoJet.hh"
13 #include "fastjet/ClusterSequence.hh"
14 #include "HiggsPair.h"
15 
16 namespace Herwig {
17  class AlpGenHandlerOL;
18 }
19 
20 //declaration of thepeg ptr
21 namespace ThePEG {
22  ThePEG_DECLARE_POINTERS(Herwig::AlpGenHandlerOL,AlpGenHandlerOLPtr);
23 }
24 
25 namespace Herwig {
26 
27 using namespace ThePEG;
28 
36 
37 public:
38 
43 
44 public:
45 
52  void persistentOutput(PersistentOStream & os) const;
53 
59  void persistentInput(PersistentIStream & is, int version);
61 
68  static void Init();
69 
70 protected:
71 
77  virtual void dofinish();
78 
84  virtual void doinit();
85 
90  virtual void doinitrun();
92 
93 public:
98  virtual bool showerHardProcessVeto() const;
99 
100 
101 protected:
102 
109  virtual IBPtr clone() const;
110 
115  virtual IBPtr fullclone() const;
117 
118 private:
119 
120  /*
121  * Run MLM jet-parton matching on the 'extra' jets.
122  */
123 
124  bool lightJetPartonVeto();
125 
126  /*
127  * Function that calculates deltaR between a parton and a jet
128  */
129 
130  double partonJetDeltaR(ThePEG::tPPtr partonptr, LorentzMomentum jetmom) const;
131 
138  void calini_m() const;
139 
140 
145  void calsim_m() const;
146 
150  void getFastJets(double rjet, Energy ejcut, double etajcut) const;
151 
152 
161  void getjet_m(double rjet, Energy ejcut, double etajcut) const;
162 
172  void caldel_m() const;
173 
183  void caldel_hvq() const;
184 
189  void getPreshowerParticles() const;
190 
196  void getShoweredParticles() const;
197 
204  void doSanityChecks(int debugLevel) const;
205 
210  void getDescendents(PPtr theParticle) const;
211 
216  void getTopRadiation(PPtr theParticle) const;
217 
222  ParticleVector pTsort(ParticleVector unsortedVec);
223  pair< vector<Energy>, vector<Lorentz5Momentum> > ETsort(vector<Energy> unsortedetjet, vector<Lorentz5Momentum> unsortedVec);
224 
225  /*
226  * A function that prints a vector of Lorentz5Momenta in a fancy way
227  */
228  void printMomVec(vector<Lorentz5Momentum> momVec);
229 
230 
231  /*
232  * A probability function for varying etclus_ about the mean value
233  */
234  Energy etclusran_(double petc) const;
235 
236 private:
237 
243 
248  AlpGenHandlerOL & operator=(const AlpGenHandlerOL &) = delete;
249 
250 private:
251 
257 
263 
271 
277 
283 
289 
298 
304 
305  /*
306  * The shower progenitors
307  */
308 
309  mutable PPtr theProgenitor;
310  mutable PPtr theLastProgenitor;
311 
317 
323 
327  mutable ShowerAlphaPtr alphaS_;
328 
336  mutable Energy pdfScale_;
337 
341  mutable Energy2 sHat_;
342 
347  mutable double alphaSME_;
349 
350  /*
351  * Number of rapidity segments of the calorimeter.
352  */
353  mutable unsigned int ncy_;
354 
355  /*
356  * Number of phi segments of the calorimeter.
357  */
358  mutable unsigned int ncphi_;
359 
360  /*
361  * Heavy flavour in WQQ,ZQQ,2Q etc (4=c, 5=b, 6=t).
362  */
363  mutable int ihvy_;
364 
365  /*
366  * Number of photons in the AlpGen process.
367  */
368  mutable int nph_;
369 
370  /*
371  * Number of higgses in the AlpGen process.
372  */
373  mutable int nh_;
374 
375  /*
376  * Jet ET cut to apply in jet clustering (in merging).
377  */
378  mutable Energy etclus_;
379 
380  /*
381  * Mean Jet ET cut to apply in jet clustering (in merging).
382  */
383  mutable Energy etclusmean_;
384 
385  /*
386  * maximum deviation from mean Jet ET cut to apply in jet clustering (in merging).
387  */
388  mutable Energy epsetclus_;
389 
390  /*
391  * type of smoothing function to use
392  */
393  mutable unsigned int smoothingtype_;
394 
395 
396 
397 
398 
399  /*
400  * Cone size used in jet clustering (in merging).
401  */
402  mutable double rclus_;
403 
404  /*
405  * Max |eta| for jets in clustering (in merging).
406  */
407  double etaclmax_;
408 
409  /*
410  * Default 1.5 factor used to decide if a jet matches a parton
411  * in merging: if DR(parton,jet)<rclusfactor*rclus the parton
412  * and jet are said to have been matched.
413  */
414  double rclusfactor_;
415 
416  /*
417  * The AlpGen hard process code. Relation to the AlpGen process names:
418  * 1: wqq, 2: zqq, 3: wjet, 4: zjet, 5: vbjet, 6: 2Q, 8: QQh, 9: Njet,
419  * 10: wcjet, 11: phjet, 12: hjet, 13: top, 14: wphjet, 15: wphqq,
420  * 16: 2Qph.
421  */
422  mutable int ihrd_;
423 
424  /*
425  * The number of light jets in the AlpGen process (i.e. the 'extra' ones).
426  */
427  mutable int njets_;
428 
429  /*
430  * Mimimum parton-parton R-sep used for generation (used for hvq merging).
431  */
432  mutable double drjmin_;
433 
434  /*
435  * This flags that the highest multiplicity ME-level process is
436  * being processed.
437  */
438  mutable bool highestMultiplicity_;
439 
440  /*
441  * This is the highest number of jets to be included in the matching.
442  * This implies that exclusive rates will be calculated up to highestNjets_-1 and
443  * inclusive for highestNjets_.
444  */
445  mutable int highestNjets_;
446 
447  /*
448  * This flags that the highest NLO multiplicity ME-level process is
449  * being processed.
450  */
451  mutable bool highestNLOMultiplicity_;
452 
453  /*
454  * This flags whether the etclus_ (merging scale) should be fixed or variable according to a prob. distribution around the mean
455  */
456  mutable bool etclusfixed_;
457 
458  /*
459  * The forwards rapidity span of the calorimeter.
460  */
461  mutable double ycmax_;
462 
463  /*
464  * The backwards rapidity span of the calorimeter.
465  */
466  mutable double ycmin_;
467 
468  /*
469  * The jet algorithm used for parton-jet matching in the MLM procedure.
470  */
471  mutable int jetAlgorithm_;
472 
473  /*
474  * Allows the vetoing to be turned off completely - just for convenience.
475  */
476  mutable bool vetoIsTurnedOff_;
477 
478  /*
479  * Switch between original and OpenLoops implementations of the handler
480  *
481  */
482  mutable int vetoType_;
483 
484  /*
485  * Signals that the LH file being read-in is a NLO (Powheg one).
486  */
487  mutable bool inputIsNLO_;
488 
489  /*
490  * Cosine of phi values of calorimeter cell centres.
491  * Goes phi~=0 to phi~=2*pi (index = 0 ---> ncphi).
492  * ==> Cosine goes from +1 ---> +1 (index = 0 ---> ncphi).
493  */
494  mutable vector<double> cphcal_;
495 
496  /*
497  * Sine of phi values of calorimeter cell centres.
498  * Goes phi~=0 to phi~=2*pi (index = 0 ---> ncphi).
499  * ==> Sine goes 0 -> 1 -> 0 -> -1 -> 0 (index = 0 ---> ncphi).
500  */
501  mutable vector<double> sphcal_;
502 
503  /*
504  * Cosine of theta values of calorimeter cell centres in Y.
505  * Goes bwds th~=pi to fwds th~=0 (index = 0 ---> ncy).
506  * ==> Cosine goes from -1 ---> +1 (index = 0 ---> ncy).
507  */
508  mutable vector<double> cthcal_;
509 
510  /*
511  * Sine of theta values of calorimeter cell centres in Y.
512  * Goes bwds th~=pi to fwds th~=0 (index = 0 ---> ncy).
513  * ==> Sine goes from 0 ---> +1 ---> 0 (index = 0 ---> ncy).
514  */
515  mutable vector<double> sthcal_;
516 
517  /*
518  * Transverse energy deposit in a given calorimeter cell.
519  * First array index corresponds to rapidity index of cell,
520  * second array index corresponds to phi cell index.
521  */
522  mutable vector<vector<Energy> > et_;
523 
524  /*
525  * For a given calorimeter cell this holds the index of the jet
526  * that the cell was clustered into.
527  */
528  mutable vector<vector<int> > jetIdx_;
529 
530  /*
531  * Vector holding the Lorentz 5 momenta of each jet.
532  */
533  mutable vector<Lorentz5Momentum> pjet_;
534 
535  /*
536  * Vector holding the list of FS particles resulting from
537  * the particle input to getDescendents.
538  */
539  mutable ParticleVector tmpList_;
540 
541  /*
542  * Variables for the C++ translation of the calini_m(), calsim_m(),
543  * getjet_m(...) and caldel_m() functions
544  */
545  mutable vector<Energy> etjet_;
546  mutable double dely_, delphi_;
547 
548 };
549 
550 }
551 
552 #include "ThePEG/Utilities/ClassTraits.h"
553 
554 namespace ThePEG {
555 
560 template <>
561 struct BaseClassTrait<Herwig::AlpGenHandlerOL,1> {
563  typedef Herwig::QTildeShowerHandler NthBase;
564 };
565 
568 template <>
569 struct ClassTraits<Herwig::AlpGenHandlerOL>
570  : public ClassTraitsBase<Herwig::AlpGenHandlerOL> {
572  static string className() { return "Herwig::AlpGenHandlerOL"; }
580  static string library() { return "AlpGenHandlerOL.so"; }
581 };
582 
585 }
586 
587 #endif /* HERWIG_AlpGenHandlerOL_H */
ThePEG::Ptr< Particle >::pointer PPtr
ParticleVector showeredFSPsToDelete_
Final-state outgoing partICLEs after shower of hard process to_be_removed from showeredFSPs_ prior to...
The QTildeShowerHandler class.
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
PPair showeredRems_
Final-state remnants after shower of hard process (look for remnants initially in showeredFSPs_)...
tPVector particlesToCluster_
ONLY the final-state particles from showeredFSPs_ (and maybe also showeredRems_) that are supposed to...
Energy pdfScale_
Information extracted from the XComb object.
Here is the documentation of the AlpGenHandlerOL class.
tPVector showeredFSPs_
Final-state outgoing partICLEs after shower of hard process (eventHandler()->currentEvent()->getFinal...
PPair showeredISHs_
Initial-state incoming hadrons after shower of hard process (eventHandler()->currentEvent()->incoming...
ShowerAlphaPtr alphaS_
Pointer to the object calculating the strong coupling.
ParticleVector preshowerFSPsToDelete_
Final-state outgoing partICLEs prior to showering to_be_removed from preShowerFSPs_ prior to the ligh...
PPair preshowerISPs_
Initial-state incoming partons prior to showering (i.e.
double alphaSME_
Constant alphaS used to generate LH events - if not already using CKKW scale (ickkw = 1 in AlpGen for...
pair< PPtr, PPtr > PPair
ThePEG::Ptr< Particle >::transient_pointer tPPtr
-*- C++ -*-
ParticleVector partonsToMatch_
ONLY the final-state partons from preshowerFSPs_ that are supposed to enter the jet-parton matching...
Energy2 sHat_
Centre of mass energy.
ParticleVector preshowerFSPs_
Final-state outgoing partICLEs prior to showering (i.e.
vector< tPPtr > tPVector
PPair showeredISPs_
Initial-state incoming partons after shower of hard process (look for partonic children of showeredIS...
vector< PPtr > ParticleVector
static ClassDescription< AlpGenHandlerOL > initAlpGenHandlerOL
The static object used to initialize the description of this class.