herwig is hosted by Hepforge, IPPP Durham
Herwig  7.1.5
HwRemDecayer.h
1 // -*- C++ -*-
2 //
3 // HwRemDecayer.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2017 The Herwig Collaboration
5 //
6 // Herwig is licenced under version 3 of the GPL, see COPYING for details.
7 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
8 //
9 #ifndef HERWIG_HwRemDecayer_H
10 #define HERWIG_HwRemDecayer_H
11 //
12 // This is the declaration of the HwRemDecayer class.
13 //
14 
15 #include "ThePEG/PDT/RemnantDecayer.h"
16 #include "ThePEG/Handlers/EventHandler.h"
17 #include "ThePEG/Repository/EventGenerator.h"
18 #include "ThePEG/EventRecord/SubProcess.h"
19 #include "ThePEG/PDF/BeamParticleData.h"
20 #include "Herwig/Shower/Core/Couplings/ShowerAlpha.h"
21 #include "Herwig/PDT/StandardMatchers.h"
23 #include "HwRemDecayer.fh"
24 
25 namespace Herwig {
26 using namespace ThePEG;
57 
58 public:
59 
61  typedef vector<pair<tPPtr, tPPtr> > PartnerMap;
62 
63 public:
64 
68  HwRemDecayer() : allowTop_(false), multiPeriph_(false), quarkPair_(false),
69  ptmin_(-1.*GeV), beta_(ZERO),
70  maxtrySoft_(10),
71  colourDisrupt_(1.0),
72  ladderbFactor_(0.0),
73  ladderPower_(-0.08),
74  ladderNorm_(1.0),
75  gaussWidth_(0.1),
76  valOfN_(0),
77  initTotRap_(0),
78  _kinCutoff(0.75*GeV),
79  _forcedSplitScale(2.5*GeV),
80  _range(1.1), _zbin(0.05),_ybin(0.),
81  _nbinmax(100), DISRemnantOpt_(0),
82  pomeronStructure_(0), mg_(ZERO) {}
83 
91  virtual bool accept(const DecayMode &) const {
92  return true;
93  }
94 
99  virtual bool canHandle(tcPDPtr particle, tcPDPtr parton) const;
100 
105  virtual bool multiCapable() const {
106  return true;
107  }
108 
116  virtual ParticleVector decay(const DecayMode & dm, const Particle & p, Step & step) const;
118 
119 public:
120 
126 
133  void persistentOutput(PersistentOStream & os) const;
134 
140  void persistentInput(PersistentIStream & is, int version);
142 
149  static void Init();
150 
154  void initialize(pair<tRemPPtr, tRemPPtr> rems, tPPair beam, Step & step,
155  Energy forcedSplitScale);
156 
162  void initSoftInteractions(Energy ptmin, InvEnergy2 beta);
163 
171  void doSplit(pair<tPPtr, tPPtr> partons, pair<tcPDFPtr, tcPDFPtr> pdfs, bool first);
172 
180  void finalize(double colourDisrupt=0.0, unsigned int softInt=0);
181 
185  void findChildren(tPPtr,vector<PPtr> &) const;
186 
187 protected:
188 
195  virtual IBPtr clone() const {return new_ptr(*this);}
196 
201  virtual IBPtr fullclone() const {return new_ptr(*this);}
203 
204 protected:
205 
213  virtual void doinit() {
215  _ybin=0.25/_zbin;
216  mg_ = getParticleData(ParticleID::g)->constituentMass();
217  }
219 
220 private:
221 
227 
232  HwRemDecayer & operator=(const HwRemDecayer &) = delete;
233 
234 public:
235 
240  struct HadronContent {
241 
245  inline void extract(int id) {
246  for(unsigned int i=0; i<flav.size(); i++) {
247  if(id == sign*flav[i]){
248  if(hadron->id() == ParticleID::gamma ||
249  (hadron->id() == ParticleID::pomeron && pomeronStructure==1) ||
250  hadron->id() == ParticleID::reggeon) {
251  flav[0] = id;
252  flav[1] = -id;
253  extracted = 0;
254  flav.resize(2);
255  }
256  else if (hadron->id() == ParticleID::pomeron && pomeronStructure==0) {
257  extracted = 0;
258  }
259  else {
260  extracted = i;
261  }
262  break;
263  }
264  }
265  }
266 
271  long RemID() const;
272 
277  bool isSeaQuark(tcPPtr parton) const {
278  return ((parton->id() != ParticleID::g) && ( !isValenceQuark(parton) ) );
279  }
280 
284  bool isValenceQuark(tcPPtr parton) const {
285  return isValenceQuark(parton->id());
286  }
287 
292  bool isSeaQuarkData(tcPDPtr partonData) const {
293  return ((partonData->id() != ParticleID::g) && ( !isValenceQuarkData(partonData) ) );
294  }
295 
299  bool isValenceQuarkData(tcPDPtr partonData) const {
300  int id(sign*partonData->id());
301  return find(flav.begin(),flav.end(),id) != flav.end();
302  }
303 
307  bool isValenceQuark(int id) const {
308  return find(flav.begin(),flav.end(),sign*id) != flav.end();
309  }
310 
312  vector<int> flav;
313 
316 
318  int sign;
319 
322 
324  unsigned int pomeronStructure;
325  };
326 
330  const pair<HadronContent, HadronContent>& content() const {
331  return theContent;
332  }
333 
337  HadronContent getHadronContent(tcPPtr hadron) const;
338 
343  theContent.first = getHadronContent(beam.first);
344  theContent.second = getHadronContent(beam.second);
345  }
346 
347 private:
348 
361  void split(tPPtr parton, HadronContent & content, tRemPPtr rem,
362  Lorentz5Momentum & used, PartnerMap & partners, tcPDFPtr pdf, bool first);
363 
370  void mergeColour(tPPtr p1, tPPtr p2, bool anti) const;
371 
378  void fixColours(PartnerMap partners, bool anti, double disrupt) const;
379 
383  void setRemMasses() const;
384 
390  PPtr finalSplit(const tRemPPtr rem, long remID,
391  Lorentz5Momentum usedMomentum) const {
392  // Create the remnant and set its momentum, also reset all of the decay
393  // products from the hadron
394  PPtr remnant = new_ptr(Particle(getParticleData(remID)));
395  Lorentz5Momentum prem(rem->momentum()-usedMomentum);
396  prem.setMass(getParticleData(remID)->constituentMass());
397  prem.rescaleEnergy();
398  remnant->set5Momentum(prem);
399  // Add the remnant to the step, but don't do colour connections
400  thestep->addDecayProduct(rem,remnant,false);
401  return remnant;
402  }
403 
404 
417  PPtr forceSplit(const tRemPPtr rem, long child, Energy &oldQ, double &oldx,
418  Lorentz5Momentum &pf, Lorentz5Momentum &p,
419  HadronContent & content) const;
420 
425  bool isPartonic(tPPtr parton) const;
426 
429 
433  Energy softPt() const;
434 
439  void softKinematics(Lorentz5Momentum &r1, Lorentz5Momentum &r2,
440  Lorentz5Momentum &g1, Lorentz5Momentum &g2) const;
441 
445  void doSoftInteractions(unsigned int N){
446  if(!multiPeriph_){
447  doSoftInteractions_old(N);}
448  else{
449  doSoftInteractions_multiPeriph(N);
450  }
451  }
452 
456  void doSoftInteractions_old(unsigned int N);
457 
461  void doSoftInteractions_multiPeriph(unsigned int N);
462 
469  tPPtr addParticle(tcPPtr parent, long id, Lorentz5Momentum p) const;
471 
476  pair<bool, bool> theanti;
477 
481  pair<double, double> theX;
482 
484  pair<HadronContent, HadronContent> theContent;
485 
487  pair<Lorentz5Momentum, Lorentz5Momentum> theUsed;
488 
493  pair<PartnerMap, PartnerMap> theMaps;
494 
501 
506  pair<RemPPtr, RemPPtr> theRems;
507 
511  mutable tcPPtr theBeam;
512 
516  mutable Ptr<BeamParticleData>::const_pointer theBeamData;
517 
521  mutable tcPDFPtr _pdf;
522 
523 private:
524 
528  bool allowTop_;
529 
534 
539 
542 
547 
551  Energy ptmin_;
552 
556  InvEnergy2 beta_;
557 
562  unsigned int maxtrySoft_;
563 
569 
575 
579  double ladderPower_;
580 
584  double ladderNorm_;
585 
591  double gaussWidth_;
592 
597  double valOfN_;
598 
603  double initTotRap_;
604 
606 
609 
613  Energy _kinCutoff;
614 
619 
623  double _range;
624 
628  double _zbin;
629 
633  double _ybin;
634 
638  int _nbinmax;
639 
643  ShowerAlphaPtr _alphaS;
644 
648  ShowerAlphaPtr _alphaEM;
649 
653  unsigned int DISRemnantOpt_;
654 
658  unsigned int pomeronStructure_;
660 
664  Energy mg_;
665 
666 };
667 
668 
669 }
670 
671 #include "ThePEG/Utilities/ClassTraits.h"
672 
673 namespace ThePEG {
674 
679 template <>
680 struct BaseClassTrait<Herwig::HwRemDecayer,1> {
682  typedef RemnantDecayer NthBase;
683 };
684 
687 template <>
688 struct ClassTraits<Herwig::HwRemDecayer>
689  : public ClassTraitsBase<Herwig::HwRemDecayer> {
691  static string className() { return "Herwig::HwRemDecayer"; }
699  static string library() { return "HwShower.so"; }
700 };
701 
704 }
705 
706 #endif /* HERWIG_HwRemDecayer_H */
virtual void doinit()
Initialize this object after the setup phase before saving an EventGenerator to disk.
Definition: HwRemDecayer.h:213
ThePEG::Ptr< Particle >::pointer PPtr
Energy _kinCutoff
The kinematic cut-off.
Definition: HwRemDecayer.h:613
vector< pair< tPPtr, tPPtr > > PartnerMap
Typedef to store information about colour partners.
Definition: HwRemDecayer.h:61
const pair< HadronContent, HadronContent > & content() const
Return the hadron content objects for the incoming particles.
Definition: HwRemDecayer.h:330
bool multiPeriph_
Switch to control using multiperipheral kinemaics.
Definition: HwRemDecayer.h:533
tcPDFPtr _pdf
The PDF for the current initial-state shower.
Definition: HwRemDecayer.h:521
int extracted
The array index of the extracted particle.
Definition: HwRemDecayer.h:315
pair< double, double > theX
variable to sum up the x values of the extracted particles
Definition: HwRemDecayer.h:481
unsigned int DISRemnantOpt_
Option for the DIS remnant.
Definition: HwRemDecayer.h:653
int sign
-1 if the particle is an anti-particle.
Definition: HwRemDecayer.h:318
StepPtr thestep
Variable to hold a pointer to the current step.
Definition: HwRemDecayer.h:500
static ClassDescription< HwRemDecayer > initHwRemDecayer
The static object used to initialize the description of this class.
Definition: HwRemDecayer.h:226
HwRemDecayer()
The default constructor.
Definition: HwRemDecayer.h:68
ThePEG::Ptr< PDFBase >::transient_const_pointer tcPDFPtr
void extract(int id)
manually extract the valence flavour id.
Definition: HwRemDecayer.h:245
vector< int > flav
The valence flavours of the corresponding baryon.
Definition: HwRemDecayer.h:312
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
Energy _forcedSplitScale
The PDF freezing scale as set in ShowerHandler.
Definition: HwRemDecayer.h:618
pair< PartnerMap, PartnerMap > theMaps
Pair of PartnerMap&#39;s to store the particles, which will be colour connected in the end...
Definition: HwRemDecayer.h:493
bool isValenceQuark(int id) const
Method to determine whether parton is a valence quark.
Definition: HwRemDecayer.h:307
unsigned int pomeronStructure
Pomeron treatment.
Definition: HwRemDecayer.h:324
pair< HadronContent, HadronContent > theContent
Pair of HadronContent structs to know about the quark content of the beams.
Definition: HwRemDecayer.h:484
InvEnergy2 beta_
slope of the soft pt-spectrum: dN/dp_T = N p_T exp(-beta*p_T^2)
Definition: HwRemDecayer.h:556
void doSoftInteractions(unsigned int N)
Create N soft gluon interactions.
Definition: HwRemDecayer.h:445
virtual bool accept(const DecayMode &) const
Check if this decayer can perfom the decay specified by the given decay mode.
Definition: HwRemDecayer.h:91
unsigned int pomeronStructure_
Option for the treatment of the pomeron structure.
Definition: HwRemDecayer.h:658
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
Definition: HwRemDecayer.h:201
pair< tPPtr, tPPtr > tPPair
pair< bool, bool > theanti
A flag which indicates, whether the extracted valence quark was a anti particle.
Definition: HwRemDecayer.h:476
Iterator find(IteratorRange< Iterator > r, const T &t)
virtual bool multiCapable() const
Return true if this decayed can extract more than one parton from a particle.
Definition: HwRemDecayer.h:105
virtual IBPtr clone() const
Make a simple clone of this object.
Definition: HwRemDecayer.h:195
bool isValenceQuarkData(tcPDPtr partonData) const
Method to determine whether parton is a valence quark.
Definition: HwRemDecayer.h:299
void setMass(Energy a)
bool quarkPair_
True if kinematics is to be calculated for quarks.
Definition: HwRemDecayer.h:538
PPtr finalSplit(const tRemPPtr rem, long remID, Lorentz5Momentum usedMomentum) const
This creates a parton from the remaining flavours of the hadron.
Definition: HwRemDecayer.h:390
Simple struct to store info about baryon quark and di-quark constituents.
Definition: HwRemDecayer.h:240
tcPDPtr hadron
The ParticleData objects of the hadron.
Definition: HwRemDecayer.h:321
bool isSeaQuark(tcPPtr parton) const
Method to determine whether parton is a quark from the sea.
Definition: HwRemDecayer.h:277
double ladderbFactor_
Variable to store the additive factor of the multiperipheral ladder multiplicity. ...
Definition: HwRemDecayer.h:574
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
double gaussWidth_
Variable to store the gaussian width of the fluctuation of the longitudinal momentum fraction...
Definition: HwRemDecayer.h:591
bool allowTop_
Switch to control handling of top quarks in proton.
Definition: HwRemDecayer.h:528
double ladderNorm_
Variable of the parameterization of the ladder multiplicity.
Definition: HwRemDecayer.h:584
ThePEG::Ptr< Particle >::transient_pointer tPPtr
double initTotRap_
Variable to store the initial total rapidity between of the remnants.
Definition: HwRemDecayer.h:603
ThePEG::Ptr< Particle >::transient_const_pointer tcPPtr
Ptr< BeamParticleData >::const_pointer theBeamData
the beam data
Definition: HwRemDecayer.h:516
tPPair softRems_
Pair of soft Remnant pointers, i.e.
Definition: HwRemDecayer.h:546
struct that is used to catch exceptions which are thrown due to energy conservation issues of additio...
Definition: HwRemDecayer.h:125
double colourDisrupt_
Variable to store the relative number of colour disrupted connections to additional soft subprocesses...
Definition: HwRemDecayer.h:568
Energy ptmin_
ptcut of the UE model
Definition: HwRemDecayer.h:551
tcPPtr theBeam
The beam particle data for the current incoming hadron.
Definition: HwRemDecayer.h:511
double _zbin
Size of the bins in z for the interpolation.
Definition: HwRemDecayer.h:628
bool isSeaQuarkData(tcPDPtr partonData) const
Method to determine whether parton is a quark from the sea.
Definition: HwRemDecayer.h:292
-*- C++ -*-
virtual void doinit()
int _nbinmax
Maximum number of bins for the z interpolation.
Definition: HwRemDecayer.h:638
pair< Lorentz5Momentum, Lorentz5Momentum > theUsed
Pair of Lorentz5Momentum to keep track of the forced splitting product momenta.
Definition: HwRemDecayer.h:487
double sign(double x)
Small helper.
Definition: RandomHelpers.h:30
double ladderPower_
Variable of the parameterization of the ladder multiplicity.
Definition: HwRemDecayer.h:579
unsigned int maxtrySoft_
Maximum number of attempts for the regeneration of an additional soft scattering, before the number o...
Definition: HwRemDecayer.h:562
Energy mg_
The gluon constituent mass.
Definition: HwRemDecayer.h:664
void setHadronContent(tPPair beam)
Set the hadron contents.
Definition: HwRemDecayer.h:342
The HwRemDecayer class is responsible for the decay of the remnants.
Definition: HwRemDecayer.h:56
double _ybin
Size of the bins in y for the interpolation.
Definition: HwRemDecayer.h:633
double valOfN_
Variable to store the current total multiplicity of a ladder.
Definition: HwRemDecayer.h:597
vector< PPtr > ParticleVector
bool isValenceQuark(tcPPtr parton) const
Method to determine whether parton is a valence quark.
Definition: HwRemDecayer.h:284
ThePEG::Ptr< Step >::pointer StepPtr
constexpr ZeroUnit ZERO
pair< RemPPtr, RemPPtr > theRems
Pair of Remnant pointers.
Definition: HwRemDecayer.h:506
ShowerAlphaPtr _alphaEM
Pointer to the object calculating the QED coupling.
Definition: HwRemDecayer.h:648
double _range
Range for emission.
Definition: HwRemDecayer.h:623
ShowerAlphaPtr _alphaS
Pointer to the object calculating the QCD coupling.
Definition: HwRemDecayer.h:643