herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
GenericMassGenerator.h
1 // -*- C++ -*-
2 //
3 // GenericMassGenerator.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2019 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_GenericMassGenerator_H
10 #define HERWIG_GenericMassGenerator_H
11 //
12 // This is the declaration of the GenericMassGenerator class.
13 //
14 #include "ThePEG/PDT/MassGenerator.h"
15 #include "ThePEG/Repository/EventGenerator.h"
16 #include "ThePEG/PDT/ParticleData.h"
17 #include "GenericMassGenerator.fh"
18 #include "ThePEG/PDT/WidthGenerator.h"
19 #include "GenericWidthGenerator.fh"
20 #include "ThePEG/Repository/CurrentGenerator.h"
21 
22 namespace Herwig {
23 using namespace ThePEG;
24 
28  class ModelGenerator;
29 
57 
61  friend class ModelGenerator;
62 
63 public:
64 
69 
73  virtual ~GenericMassGenerator();
74 
81  void persistentOutput(PersistentOStream & os) const;
82 
88  void persistentInput(PersistentIStream & is, int version);
90 
94  static void Init();
95 
96 public:
97 
103  bool accept(const ParticleData & part) const;
104 
112  Energy mass(const ParticleData & part) const {
113  return mass(part,lowerMass_,upperMass_);
114  }
115 
123  Energy mass(const ParticleData & part,
124  const Energy low,const Energy upp) const {
125  if(upp<low) return low;
126  Energy output;
127  int ntry=0; double wgt=0.;
128  do {
129  ++ntry;
130  output=mass(wgt,part,low,upp,3);
131  if(wgt>maxWgt_) maxWgt_=wgt;
132  }
133  while(maxWgt_*(UseRandom::rnd())>wgt&&ntry<nGenerate_);
134  return (ntry>=nGenerate_) ? mass_ : output;
135  }
136 
144  Energy mass(double & wgt, const ParticleData & part,
145  double r=UseRandom::rnd()) const {
146  return mass(wgt,part,lowerMass_,upperMass_,r);
147  }
148 
158  Energy mass(double & wgt, const ParticleData & part,
159  const Energy low,const Energy upp,
160  double r=UseRandom::rnd()) const {
161  return mass(wgt,part,low,upp,BWShape_,r);
162  }
163 
169  virtual double weight(Energy q) const {
170  return weight(q,BWShape_);
171  }
172 
176  virtual InvEnergy2 BreitWignerWeight(Energy q) {
177  return BreitWignerWeight(q,BWShape_);
178  }
180 
184  virtual void dataBaseOutput(ofstream &,bool);
185 
186 public:
187 
195  pair<Energy,Energy> width(Energy q,int shape) const;
196 
200  Energy lowerLimit() const {return lowerMass_;}
201 
205  Energy upperLimit() const {return upperMass_;}
206 
210  Energy nominalMass() const {return mass_;}
211 
215  Energy nominalWidth() const {return width_;}
216 
217 protected:
218 
228  virtual Energy mass(double & wgt, const ParticleData & ,
229  const Energy low,const Energy upp, int shape,
230  double r=UseRandom::rnd()) const {
231  // calculate the mass square using fixed width BW
232  Energy lo=max(low,lowerMass_),up=min(upp,upperMass_);
233  double rhomin=atan2((lo*lo-mass2_),mWidth_);
234  double rhomax=atan2((up*up-mass2_),mWidth_)-rhomin;
235  double rho=rhomin+rhomax*r;
236  Energy2 q2 = mass2_+mWidth_*tan(rho);
237  Energy q = sqrt(q2);
238  wgt = rhomax*weight(q,shape);
239  // return the mass
240  return q;
241  }
242 
251  Energy mass(double & wgt, const ParticleData & part, int shape,
252  double r=UseRandom::rnd()) const {
253  return mass(wgt,part,lowerMass_,upperMass_,shape,r);
254  }
255 
262  inline virtual double weight(Energy q, int shape) const {
263  Energy2 q2 = q*q;
264  Energy4 sq=sqr(q2-mass2_);
265  pair<Energy,Energy> gam=width(q,shape);
266  // finish the calculation of the width
267  Energy2 num;
268  if(shape==2) num = mass_*gam.first;
269  else if(shape==3) num = mass_*gam.first;
270  else num = q *gam.first;
271  Energy4 den = (shape==2) ? sq+mass2_*gam.second*gam.second : sq+q2*gam.second*gam.second;
272  return num/den*(sq+mWidth_*mWidth_)/Constants::pi/mWidth_;
273  }
274 
278  virtual InvEnergy2 BreitWignerWeight(Energy q, int shape) const {
279  Energy2 q2 = q*q;
280  Energy4 sq=sqr(q2-mass2_);
281  pair<Energy,Energy> gam=width(q,shape);
282  // finish the calculation of the width
283  Energy2 num;
284  if(shape==2) num = mass_*gam.first;
285  else if(shape==3) num = mass_*gam.first;
286  else num = q *gam.first;
287  Energy4 den = (shape==2) ? sq+mass2_*gam.second*gam.second : sq+q2*gam.second*gam.second;
288  return num/den/Constants::pi;
289  }
290 
294  tcPDPtr particle() const {return particle_;}
295 
299  void particle(tPDPtr in) {particle_ = in;}
300 
301 protected:
302 
309  virtual IBPtr clone() const;
310 
315  virtual IBPtr fullclone() const;
317 
318 protected:
319 
327  virtual void doinit();
328 
333  virtual void dofinish();
334 
344  virtual void rebind(const TranslationMap & trans)
345  ;
346 
352  virtual IVector getReferences();
354 
355 private:
356 
360  GenericMassGenerator & operator=(const GenericMassGenerator &) = delete;
361 
362 private:
363 
367  void setParticle(string);
368 
372  string getParticle() const;
373 
374 private:
375 
379  mutable double maxWgt_;
380 
384  int BWShape_;
385 
390 
391 private:
392 
397 
401  Energy lowerMass_;
402 
406  Energy upperMass_;
407 
411  Energy mass_;
412 
416  Energy width_;
417 
421  Energy2 mass2_;
422 
426  Energy2 mWidth_;
427 
432 
437 
441  bool output_;
442 
446  WidthGeneratorPtr widthGen_;
447 
451  GenericWidthGeneratorPtr widthGenB_;
452 
456  bool widthOpt_;
457 
458 };
459 
460 }
461 
462 
463 #endif /* HERWIG_GenericMassGenerator_H */
Energy mass(const ParticleData &part) const
Generate a mass using the default limits.
Energy mass_
Mass of the particle.
double sqrt(int x)
int BWShape_
parameter controlling the shape of the Breit-Wigner
Energy lowerLimit() const
Lower limit on the mass.
ThePEG::Ptr< ParticleData >::transient_pointer tPDPtr
bool output_
Whether or not to output the data to a file.
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
Energy mass(double &wgt, const ParticleData &part, double r=UseRandom::rnd()) const
Return a mass with the weight using the default limits.
constexpr auto sqr(const T &x) -> decltype(x *x)
virtual InvEnergy2 BreitWignerWeight(Energy q, int shape) const
Return the full weight.
virtual InvEnergy2 BreitWignerWeight(Energy q)
Return the full weight.
int nGenerate_
Number of attempts to generate the mass.
static double rnd()
Energy nominalMass() const
Default mass.
void particle(tPDPtr in)
Set the particle.
GenericWidthGeneratorPtr widthGenB_
Pointer to the width generator.
Energy2 mass2_
Mass of the particle squared.
Energy upperLimit() const
Upper limit on the mass.
virtual double weight(Energy q) const
Weight for the factor.
tcPDPtr particle() const
Accesss to the particle.
Energy mass(double &wgt, const ParticleData &part, int shape, double r=UseRandom::rnd()) const
Return a mass with the weight using the default limits.
Energy2 mWidth_
Mass of the particle times the width.
Energy lowerMass_
Lower limit on the particle&#39;s mass.
Energy width_
Width of the particle.
int nInitial_
Number of weights to generate when initializing.
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
virtual double weight(Energy q, int shape) const
Weight for the factor.
tPDPtr particle_
Pointer to the particle.
Energy nominalWidth() const
Default Width.
double maxWgt_
The maximum weight for unweighting when generating the mass.
-*- C++ -*-
bool widthOpt_
Option for the treatment of the width.
This class is designed to store the particles in some model and then call the appropriate function to...
bool initialize_
Whether or not to initialize the GenericMassGenerator.
Energy mass(double &wgt, const ParticleData &part, const Energy low, const Energy upp, double r=UseRandom::rnd()) const
Return a mass with the weight using the specified limits.
vector< IBPtr > IVector
Energy mass(const ParticleData &part, const Energy low, const Energy upp) const
Generate a mass using specified limits.
virtual Energy mass(double &wgt, const ParticleData &, const Energy low, const Energy upp, int shape, double r=UseRandom::rnd()) const
Return a mass with the weight using the specified limits.
Energy upperMass_
Upper limit on the particle&#39;s mass.
WidthGeneratorPtr widthGen_
Pointer to the width generator.
The GenericMassGenerator class is a simple class for the generation of particle masses in Herwig...