Page MenuHomeHEPForge
Files F7879302

No OneTemporary
Actions

View Options

diff --git a/Models/ADD/GravitonMassGenerator.cc b/Models/ADD/GravitonMassGenerator.cc
--- a/Models/ADD/GravitonMassGenerator.cc
+++ b/Models/ADD/GravitonMassGenerator.cc
@@ -1,91 +1,92 @@
// -*- C++ -*-
//
// This is the implementation of the non-inlined, non-templated member
// functions of the GravitonMassGenerator class.
//
#include "GravitonMassGenerator.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Interface/Parameter.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
+#include "Herwig/PDT/GenericWidthGenerator.h"
#include "ADDModel.h"
using namespace Herwig;
GravitonMassGenerator::GravitonMassGenerator()
: prefactor_(0.), delta_(2), md_(1000.*GeV), mMin_(MeV)
{}
IBPtr GravitonMassGenerator::clone() const {
return new_ptr(*this);
}
IBPtr GravitonMassGenerator::fullclone() const {
return new_ptr(*this);
}
void GravitonMassGenerator::persistentOutput(PersistentOStream & os) const {
os << prefactor_ << delta_ << ounit(md_,GeV) << ounit(mMin_,GeV);
}
void GravitonMassGenerator::persistentInput(PersistentIStream & is, int) {
is >> prefactor_ >> delta_ >> iunit(md_,GeV) >> iunit(mMin_,GeV);
}
ClassDescription<GravitonMassGenerator>
GravitonMassGenerator::initGravitonMassGenerator;
// Definition of the static class description member.
void GravitonMassGenerator::Init() {
static ClassDocumentation<GravitonMassGenerator> documentation
("The GravitonMassGenerator class generates the mass for external gravitions "
"in the ADD model.");
static Parameter<GravitonMassGenerator,Energy> interfaceMinimumMass
("MinimumMass",
"Minimum gravition mass to avoid numerical problems",
&GravitonMassGenerator::mMin_, GeV, MeV, MeV, GeV,
false, false, Interface::limited);
}
void GravitonMassGenerator::doinit() {
GenericMassGenerator::doinit();
tcHwADDPtr hwADD = dynamic_ptr_cast<tcHwADDPtr>(generator()->standardModel());
if(!hwADD)
throw Exception() << "Must have ADDModel in GravitonMassGenerator::doinit()"
<< Exception::runerror;
delta_ = hwADD->delta();
md_ = hwADD->MD();
// calculate the prefactor
prefactor_ = sqr(hwADD->MPlanckBar()/md_);
// even no of dimensions
if(delta_%2==0) {
unsigned int n = delta_/2;
prefactor_ *= 2.*pow(Constants::pi,int(n));
for(unsigned int ix=1;ix<n;++ix) {
prefactor_ /= double(ix);
}
}
// odd number of dimensions
else {
unsigned int n = (delta_-1)/2;
prefactor_ *= 2.*pow(Constants::pi,int(n));
for(unsigned int ix=0;ix<n;++ix) {
prefactor_ /= double(ix)+0.5;
}
}
}
Energy GravitonMassGenerator::mass(double & wgt, const ParticleData & ,
const Energy low,const Energy upp, int,
double r) const {
Energy low2 = max(mMin_,low);
double rlow = pow(double(low2/md_),int(delta_))/double(delta_);
double rupp = pow(double(upp /md_),int(delta_))/double(delta_);
double rho = rlow + (rupp-rlow)*r;
wgt = (rupp-rlow)*prefactor_;
return pow(double(delta_)*rho,1./double(delta_))*md_;
}
diff --git a/PDT/GenericMassGenerator.cc b/PDT/GenericMassGenerator.cc
--- a/PDT/GenericMassGenerator.cc
+++ b/PDT/GenericMassGenerator.cc
@@ -1,232 +1,274 @@
// -*- C++ -*-
//
// GenericMassGenerator.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2011 The Herwig Collaboration
//
// Herwig is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
//
// This is the implementation of the non-inlined, non-templated member
// functions of the GenericMassGenerator class.
//
#include "GenericMassGenerator.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Interface/Parameter.h"
#include "ThePEG/Interface/Switch.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
#include "ThePEG/Utilities/Throw.h"
#include "ThePEG/Utilities/Rebinder.h"
+#include "GenericWidthGenerator.h"
using namespace Herwig;
using namespace ThePEG;
+GenericMassGenerator::GenericMassGenerator()
+ : maxWgt_(),BWShape_(0),nGenerate_(100),
+ lowerMass_(),upperMass_(),
+ mass_(),width_(),mass2_(),mWidth_(),
+ nInitial_(1000),
+ initialize_(false), output_(false), widthOpt_(false) {
+}
+
+GenericMassGenerator::~GenericMassGenerator() {
+}
+
+IBPtr GenericMassGenerator::clone() const {return new_ptr(*this);}
+IBPtr GenericMassGenerator::fullclone() const {return new_ptr(*this);}
+
void GenericMassGenerator::persistentOutput(PersistentOStream & os) const {
os << particle_
<< ounit(lowerMass_,GeV) << ounit(upperMass_,GeV) << maxWgt_
<< BWShape_ << nGenerate_
<< ounit(mass_,GeV) << ounit(width_,GeV)
<< ounit(mass2_,GeV2) << ounit(mWidth_,GeV2)
- << nInitial_ << initialize_ << output_ << widthGen_ << widthOpt_;
+ << nInitial_ << initialize_ << output_ << widthGen_ << widthGenB_
+ << widthOpt_;
}
void GenericMassGenerator::persistentInput(PersistentIStream & is, int) {
is >> particle_
>> iunit(lowerMass_,GeV) >> iunit(upperMass_,GeV) >> maxWgt_
>> BWShape_ >> nGenerate_
>> iunit(mass_,GeV) >> iunit(width_ ,GeV)
>> iunit(mass2_,GeV2) >> iunit(mWidth_ ,GeV2)
- >> nInitial_ >> initialize_ >> output_ >> widthGen_ >> widthOpt_;
+ >> nInitial_ >> initialize_ >> output_ >> widthGen_ >> widthGenB_
+ >> widthOpt_;
}
ClassDescription<GenericMassGenerator> GenericMassGenerator::initGenericMassGenerator;
// Definition of the static class description member.
void GenericMassGenerator::Init() {
static ClassDocumentation<GenericMassGenerator> documentation
("The GenericMassGenerator class is the main class for"
" mass generation in Herwig.");
static Parameter<GenericMassGenerator,string> interfaceParticle
("Particle",
"The particle data object for this class",
0, "", true, false,
&GenericMassGenerator::setParticle,
&GenericMassGenerator::getParticle);
static Switch<GenericMassGenerator,bool> interfaceInitialize
("Initialize",
"Initialize the calculation of the maximum weight etc",
&GenericMassGenerator::initialize_, false, false, false);
static SwitchOption interfaceInitializeInitialization
(interfaceInitialize,
"Yes",
"Do the initialization",
true);
static SwitchOption interfaceInitializeNoInitialization
(interfaceInitialize,
"No",
"Don't do the initalization",
false);
static Switch<GenericMassGenerator,bool> interfaceOutput
("Output",
"Output the setup",
&GenericMassGenerator::output_, false, false, false);
static SwitchOption interfaceOutputYes
(interfaceOutput,
"Yes",
"Output the data",
true);
static SwitchOption interfaceOutputNo
(interfaceOutput,
"No",
"Don't output the data",
false);
static Switch<GenericMassGenerator,int> interfaceBreitWignerShape
("BreitWignerShape",
"Controls the shape of the mass distribution generated",
&GenericMassGenerator::BWShape_, 0, false, false);
static SwitchOption interfaceBreitWignerShapeDefault
(interfaceBreitWignerShape,
"Default",
"Running width with q in numerator and denominator width factor",
0);
static SwitchOption interfaceBreitWignerShapeFixedWidth
(interfaceBreitWignerShape,
"FixedWidth",
"Use a fixed width",
1);
static SwitchOption interfaceBreitWignerShapeNoq
(interfaceBreitWignerShape,
"Noq",
"Use M rather than q in the numerator and denominator width factor",
2);
static SwitchOption interfaceBreitWignerShapeNoNumerator
(interfaceBreitWignerShape,
"NoNumerator",
"Neglect the numerator factors",
3);
static Parameter<GenericMassGenerator,double> interfaceMaximumWeight
("MaximumWeight",
"The maximum weight for the unweighting",
&GenericMassGenerator::maxWgt_, 1.0, 0.0, 1000.0,
false, false, true);
static Parameter<GenericMassGenerator,int> interfaceNGenerate
("NGenerate",
"The number of tries to generate the mass",
&GenericMassGenerator::nGenerate_, 100, 0, 10000,
false, false, true);
static Parameter<GenericMassGenerator,int> interfaceNInitial
("NInitial",
"Number of tries for the initialisation",
&GenericMassGenerator::nInitial_, 1000, 0, 100000,
false, false, true);
static Switch<GenericMassGenerator,bool> interfaceWidthOption
("WidthOption",
"Which width to use",
&GenericMassGenerator::widthOpt_, false, false, false);
static SwitchOption interfaceWidthOptionNominalWidth
(interfaceWidthOption,
"NominalWidth",
"Use the normal width from the particle data object",
false);
static SwitchOption interfaceWidthOptionPhysicalWidth
(interfaceWidthOption,
"PhysicalWidth",
"Use the width calculated at the on-shell mass",
true);
}
bool GenericMassGenerator::accept(const ParticleData & in) const {
if(!particle_) return false;
return in.id() == particle_->id() ||
( particle_->CC() && particle_->CC()->id() == in.id() );
}
void GenericMassGenerator::doinit() {
MassGenerator::doinit();
if(particle_->massGenerator()!=this) return;
// the width generator
particle_->init();
widthGen_=particle_->widthGenerator();
- if(widthGen_){widthGen_->init();}
+ if(widthGen_) widthGen_->init();
+ widthGenB_ = dynamic_ptr_cast<GenericWidthGeneratorPtr>(widthGen_);
// local storage of particle properties for speed
mass_=particle_->mass();
width_= widthOpt_ ? particle_->generateWidth(mass_) : particle_->width();
mass2_=mass_*mass_;
mWidth_=mass_*width_;
lowerMass_ = mass_-particle_->widthLoCut();
upperMass_ = mass_+particle_->widthUpCut();
// print out messagw if doing the initialisation
if(initialize_) {
// zero the maximum weight
maxWgt_=0.;
// storage of variable for the loop
double wgt=0.;
// perform the initialisation
// (mass() contains a ThePEG::Random call)
for(int ix=0; ix<nInitial_; ++ix) {
mass(wgt, *particle_, 3);
if ( wgt > maxWgt_ )
maxWgt_ = wgt;
}
}
}
void GenericMassGenerator::dataBaseOutput(ofstream & output, bool header) {
if(header) output << "update Mass_Generators set parameters=\"";
output << "newdef " << fullName() << ":BreitWignerShape " << BWShape_ << "\n";
output << "newdef " << fullName() << ":MaximumWeight " << maxWgt_ << "\n";
output << "newdef " << fullName() << ":NGenerate " << nGenerate_ << "\n";
output << "newdef " << fullName() << ":WidthOption " << widthOpt_ << "\n";
if(header) output << "\n\" where BINARY ThePEGFullName=\""
<< fullName() << "\";" << endl;
}
void GenericMassGenerator::dofinish() {
if(output_) {
string fname = CurrentGenerator::current().filename() +
string("-") + name() + string(".output");
ofstream output(fname.c_str());
dataBaseOutput(output,true);
}
MassGenerator::dofinish();
}
void GenericMassGenerator::setParticle(string p) {
if ( (particle_ = Repository::GetPtr<tPDPtr>(p)) ) return;
particle_ = Repository::findParticle(StringUtils::basename(p));
if ( ! particle_ )
Throw<InterfaceException>()
<< "Could not set Particle interface "
<< "for the object \"" << name()
<< "\". Particle \"" << StringUtils::basename(p) << "\" not found."
<< Exception::runerror;
}
string GenericMassGenerator::getParticle() const {
return particle_ ? particle_->fullName() : "";
}
void GenericMassGenerator::rebind(const TranslationMap & trans)
{
particle_ = trans.translate(particle_);
MassGenerator::rebind(trans);
}
IVector GenericMassGenerator::getReferences() {
IVector ret = MassGenerator::getReferences();
ret.push_back(particle_);
return ret;
}
+
+pair<Energy,Energy> GenericMassGenerator::width(Energy q,int shape) const {
+ Energy gam;
+ if(shape==3 || BWShape_==1 || !widthGen_ ) {
+ gam = width_;
+ }
+ else if(!widthGenB_) {
+ gam = widthGen_->width(*particle_,q);
+ }
+ else {
+ return widthGenB_->width(q,*particle_);
+ }
+ pair<Energy,Energy> output;
+ output.first = width_;
+ output.second = width_;
+ // take BR's into account
+ double brSum=0.;
+ for(DecaySet::const_iterator it = particle_->decayModes().begin();
+ it!=particle_->decayModes().end();++it) {
+ if((**it).on()) brSum += (**it).brat();
+ }
+ output.first *= brSum;
+ return output;
+}
diff --git a/PDT/GenericMassGenerator.h b/PDT/GenericMassGenerator.h
--- a/PDT/GenericMassGenerator.h
+++ b/PDT/GenericMassGenerator.h
@@ -1,496 +1,499 @@
// -*- C++ -*-
//
// GenericMassGenerator.h is a part of Herwig - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2011 The Herwig Collaboration
//
// Herwig is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
#ifndef HERWIG_GenericMassGenerator_H
#define HERWIG_GenericMassGenerator_H
//
// This is the declaration of the GenericMassGenerator class.
//
#include "ThePEG/PDT/MassGenerator.h"
#include "ThePEG/Repository/EventGenerator.h"
#include "ThePEG/PDT/ParticleData.h"
+#include "GenericMassGenerator.fh"
#include "ThePEG/PDT/WidthGenerator.h"
-#include "GenericMassGenerator.fh"
+#include "GenericWidthGenerator.fh"
#include "ThePEG/Repository/CurrentGenerator.h"
namespace Herwig {
using namespace ThePEG;
/**
* Declare ModelGenerator class as must be friend to set the particle
*/
class ModelGenerator;
/** \ingroup PDT
*
* The <code>GenericMassGenerator</code> class is a simple class for the
* generation of particle masses in Herwig. It inherits from the
* <code>MassGenerator</code> class of ThePEG and implements a Breit-Wigner
* using the width generator to give the running width.
*
* In general the width generator will be an instance of the
* <code>GenericWidthGenerator</code> class which uses the Herwig decayers
* based on the <code>DecayIntegrator</code> class to define the shape of the
* running width.
*
* This class is designed so that the weight
*
* \f[\int dm^2 \frac{m\Gamma(m)}{(m^2-M^2)^2+m^2\Gamma^2(m)}\f]
*
* can be included in the production of the particle to take off-shell effects into
* account. This is the default form of the weight.
* The numerator and running of the width can
* be changed using the BreitWignerShape interface.
*
* @see MassGenerator
* @see DecayIntegrator
* @see GenericWidthGenerator
*
*/
class GenericMassGenerator: public MassGenerator {
/**
* ModelGenerator class as must be friend to set the particle
*/
friend class ModelGenerator;
public:
/**
* Default constructor
*/
- GenericMassGenerator()
- : maxWgt_(),BWShape_(0),nGenerate_(100),
- lowerMass_(),upperMass_(),
- mass_(),width_(),mass2_(),mWidth_(),
- nInitial_(1000),
- initialize_(false), output_(false), widthOpt_(false) {}
+ GenericMassGenerator();
+
+ /**
+ * Destructor
+ */
+ virtual ~GenericMassGenerator();
/** @name Functions used by the persistent I/O system. */
//@{
/**
* Function used to write out object persistently.
* @param os the persistent output stream written to.
*/
void persistentOutput(PersistentOStream & os) const;
/**
* Function used to read in object persistently.
* @param is the persistent input stream read from.
* @param version the version number of the object when written.
*/
void persistentInput(PersistentIStream & is, int version);
//@}
/**
* Standard Init function used to initialize the interfaces.
*/
static void Init();
public:
/**
* Return true if this mass generator can handle the given particle type.
* @param part The particle data pointer of the particle.
* @return True ig this class can handle the particle and false otherwise
*/
bool accept(const ParticleData & part) const;
/** @name Members to generate the mass of a particle instance */
//@{
/**
* Generate a mass using the default limits.
* @param part The particle data pointer of the particle.
* @return The mass of the particle instance.
*/
Energy mass(const ParticleData & part) const {
return mass(part,lowerMass_,upperMass_);
}
/**
* Generate a mass using specified limits
* @param part The particle data pointer of the particle.
* @param low The lower limit on the particle's mass.
* @param upp The upper limit on the particle's mass.
* @return The mass of the particle instance.
*/
Energy mass(const ParticleData & part,
const Energy low,const Energy upp) const {
if(upp<low) return low;
Energy output;
int ntry=0; double wgt=0.;
do {
++ntry;
output=mass(wgt,part,low,upp,3);
if(wgt>maxWgt_) maxWgt_=wgt;
}
while(maxWgt_*(UseRandom::rnd())>wgt&&ntry<nGenerate_);
return (ntry>=nGenerate_) ? mass_ : output;
}
/**
* Return a mass with the weight using the default limits.
* @param part The particle data pointer of the particle.
* @param wgt The weight for this mass.
* @param r The random number used for the weight
* @return The mass of the particle instance.
*/
Energy mass(double & wgt, const ParticleData & part,
double r=UseRandom::rnd()) const {
return mass(wgt,part,lowerMass_,upperMass_,r);
}
/**
* Return a mass with the weight using the specified limits.
* @param part The particle data pointer of the particle.
* @param low The lower limit on the particle's mass.
* @param upp The upper limit on the particle's mass.
* @param wgt The weight for this mass.
* @param r The random number used for the weight
* @return The mass of the particle instance.
*/
Energy mass(double & wgt, const ParticleData & part,
const Energy low,const Energy upp,
double r=UseRandom::rnd()) const {
return mass(wgt,part,low,upp,BWShape_,r);
}
/**
* Weight for the factor.
* @param q The mass of the instance
* @return The weight.
*/
virtual double weight(Energy q) const {
return weight(q,BWShape_);
}
/**
* Return the full weight
*/
virtual InvEnergy2 BreitWignerWeight(Energy q) {
return BreitWignerWeight(q,BWShape_);
}
//@}
/**
* Output the initialisation info for the database
*/
virtual void dataBaseOutput(ofstream &,bool);
public:
/** @name Access to particle properties */
//@{
/**
* The running width.
* @param q The mass for the calculation of the running width
* @return The running width.
*/
- Energy width(Energy q) const {
- return (BWShape_==1||!widthGen_) ?
- width_ : widthGen_->width(*particle_,q);
- }
+ pair<Energy,Energy> width(Energy q,int shape) const;
/**
* Lower limit on the mass
*/
Energy lowerLimit() const {return lowerMass_;}
/**
* Upper limit on the mass
*/
Energy upperLimit() const {return upperMass_;}
/**
* Default mass
*/
Energy nominalMass() const {return mass_;}
/**
* Default Width
*/
Energy nominalWidth() const {return width_;}
protected:
/**
* Return a mass with the weight using the specified limits.
* @param low The lower limit on the particle's mass.
* @param upp The upper limit on the particle's mass.
* @param wgt The weight for this mass.
* @param shape The type of shape to use
* @param r The random number used for the weight
* @return The mass of the particle instance.
*/
virtual Energy mass(double & wgt, const ParticleData & ,
const Energy low,const Energy upp, int shape,
double r=UseRandom::rnd()) const {
// calculate the mass square using fixed width BW
Energy lo=max(low,lowerMass_),up=min(upp,upperMass_);
double rhomin=atan2((lo*lo-mass2_),mWidth_);
double rhomax=atan2((up*up-mass2_),mWidth_)-rhomin;
double rho=rhomin+rhomax*r;
Energy2 q2 = mass2_+mWidth_*tan(rho);
Energy q = sqrt(q2);
wgt = rhomax*weight(q,shape);
// return the mass
return q;
}
/**
* Return a mass with the weight using the default limits.
* @param part The particle data pointer of the particle.
* @param wgt The weight for this mass.
* @param shape The type of shape to use
* @param r The random number used for the weight
* @return The mass of the particle instance.
*/
Energy mass(double & wgt, const ParticleData & part, int shape,
double r=UseRandom::rnd()) const {
return mass(wgt,part,lowerMass_,upperMass_,shape,r);
}
/**
* Weight for the factor.
* @param q The mass of the instance
* @param shape The type of shape to use as for the BreitWignerShape interface
* @return The weight.
*/
inline virtual double weight(Energy q, int shape) const {
Energy2 q2 = q*q;
Energy4 sq=sqr(q2-mass2_);
- Energy gam=width(q);
+ pair<Energy,Energy> gam=width(q,shape);
// finish the calculation of the width
Energy2 num;
- if(shape==2) num = mass_*gam;
- else if(shape==3) num = mass_*width_;
- else num = q*gam;
- Energy4 den = (shape==2) ? sq+mass2_*gam*gam : sq+q2*gam*gam;
+ if(shape==2) num = mass_*gam.first;
+ else if(shape==3) num = mass_*gam.first;
+ else num = q *gam.first;
+ Energy4 den = (shape==2) ? sq+mass2_*gam.second*gam.second : sq+q2*gam.second*gam.second;
return num/den*(sq+mWidth_*mWidth_)/Constants::pi/mWidth_;
}
/**
* Return the full weight
*/
virtual InvEnergy2 BreitWignerWeight(Energy q, int shape) const {
Energy2 q2 = q*q;
Energy4 sq=sqr(q2-mass2_);
- Energy gam=width(q);
+ pair<Energy,Energy> gam=width(q,shape);
// finish the calculation of the width
Energy2 num;
- if(shape==2) num = mass_*gam;
- else if(shape==3) num = mass_*width_;
- else num = q*gam;
- Energy4 den = (shape==2) ? sq+mass2_*gam*gam : sq+q2*gam*gam;
+ if(shape==2) num = mass_*gam.first;
+ else if(shape==3) num = mass_*gam.first;
+ else num = q *gam.first;
+ Energy4 den = (shape==2) ? sq+mass2_*gam.second*gam.second : sq+q2*gam.second*gam.second;
return num/den/Constants::pi;
}
/**
* Accesss to the particle
*/
tcPDPtr particle() const {return particle_;}
/**
* Set the particle
*/
void particle(tPDPtr in) {particle_ = in;}
protected:
/** @name Clone Methods. */
//@{
/**
* Make a simple clone of this object.
* @return a pointer to the new object.
*/
- virtual IBPtr clone() const {return new_ptr(*this);}
+ virtual IBPtr clone() const;
/** Make a clone of this object, possibly modifying the cloned object
* to make it sane.
* @return a pointer to the new object.
*/
- virtual IBPtr fullclone() const {return new_ptr(*this);}
+ virtual IBPtr fullclone() const;
//@}
protected:
/** @name Standard Interfaced functions. */
//@{
/**
* Initialize this object after the setup phase before saving and
* EventGenerator to disk.
* @throws InitException if object could not be initialized properly.
*/
virtual void doinit();
/**
* Finalize this object. Called in the run phase just after a
* run has ended. Used eg. to write out statistics.
*/
virtual void dofinish();
/**
* Rebind pointer to other Interfaced objects. Called in the setup phase
* after all objects used in an EventGenerator has been cloned so that
* the pointers will refer to the cloned objects afterwards.
* @param trans a TranslationMap relating the original objects to
* their respective clones.
* @throws RebindException if no cloned object was found for a given
* pointer.
*/
virtual void rebind(const TranslationMap & trans)
;
/**
* Return a vector of all pointers to Interfaced objects used in this
* object.
* @return a vector of pointers.
*/
virtual IVector getReferences();
//@}
private:
/**
* Describe a concrete class with persistent data.
*/
static ClassDescription<GenericMassGenerator> initGenericMassGenerator;
/**
* Private and non-existent assignment operator.
*/
GenericMassGenerator & operator=(const GenericMassGenerator &);
private:
/**
* Helper function for the interface
*/
void setParticle(string);
/**
* Helper function for the interface
*/
string getParticle() const;
private:
/**
* The maximum weight for unweighting when generating the mass.
*/
mutable double maxWgt_;
/**
* parameter controlling the shape of the Breit-Wigner
*/
int BWShape_;
/**
* Number of attempts to generate the mass.
*/
int nGenerate_;
private:
/**
* Pointer to the particle
*/
tPDPtr particle_;
/**
* Lower limit on the particle's mass
*/
Energy lowerMass_;
/**
* Upper limit on the particle's mass
*/
Energy upperMass_;
/**
* Mass of the particle
*/
Energy mass_;
/**
* Width of the particle
*/
Energy width_;
/**
* Mass of the particle squared.
*/
Energy2 mass2_;
/**
* Mass of the particle times the width.
*/
Energy2 mWidth_;
/**
* Number of weights to generate when initializing
*/
int nInitial_;
/**
* Whether or not to initialize the GenericMassGenerator
*/
bool initialize_;
/**
* Whether or not to output the data to a file
*/
bool output_;
/**
* Pointer to the width generator
*/
WidthGeneratorPtr widthGen_;
/**
+ * Pointer to the width generator
+ */
+ GenericWidthGeneratorPtr widthGenB_;
+
+ /**
* Option for the treatment of the width
*/
bool widthOpt_;
};
}
#include "ThePEG/Utilities/ClassTraits.h"
namespace ThePEG {
/** @cond TRAITSPECIALIZATIONS */
/**
* The following template specialization informs ThePEG about the
* base class of GenericMassGenerator.
*/
template <>
struct BaseClassTrait<Herwig::GenericMassGenerator,1> {
/** Typedef of the base class of GenericMassGenerator. */
typedef MassGenerator NthBase;
};
/**
* The following template specialization informs ThePEG about the
* name of this class and the shared object where it is defined.
*/
template <>
struct ClassTraits<Herwig::GenericMassGenerator>
: public ClassTraitsBase<Herwig::GenericMassGenerator> {
/** Return the class name. */
static string className() { return "Herwig::GenericMassGenerator"; }
};
/** @endcond */
}
#endif /* HERWIG_GenericMassGenerator_H */
diff --git a/PDT/ScalarMassGenerator.cc b/PDT/ScalarMassGenerator.cc
--- a/PDT/ScalarMassGenerator.cc
+++ b/PDT/ScalarMassGenerator.cc
@@ -1,163 +1,164 @@
// -*- C++ -*-
//
// ScalarMassGenerator.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2011 The Herwig Collaboration
//
// Herwig is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
//
// This is the implementation of the non-inlined, non-templated member
// functions of the ScalarMassGenerator class.
//
#include "ScalarMassGenerator.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Interface/ParVector.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
+#include "GenericWidthGenerator.h"
using namespace Herwig;
void ScalarMassGenerator::persistentOutput(PersistentOStream & os) const {
os << ounit(_coupling,GeV) << ounit(_mass1,GeV) << ounit(_mass2,GeV)
<< ounit(_m2plus,GeV2) << ounit(_m2minus,GeV2);
}
void ScalarMassGenerator::persistentInput(PersistentIStream & is, int) {
is >> iunit(_coupling,GeV) >> iunit(_mass1,GeV) >> iunit(_mass2,GeV)
>> iunit(_m2plus,GeV2) >> iunit(_m2minus,GeV2);
}
ClassDescription<ScalarMassGenerator> ScalarMassGenerator::initScalarMassGenerator;
// Definition of the static class description member.
void ScalarMassGenerator::Init() {
static ClassDocumentation<ScalarMassGenerator> documentation
("The ScalarMassGenerator class inherits from the "
"GenericMassGenerator class and includes finite width effects "
"for the scalar f_0 and a_0 mesons.",
"Finite width effects for the scalar $f_0$ and $a_0$ mesons follow \\cite{Flatte:1976xu}.",
"%\\cite{Flatte:1976xu}\n"
"\\bibitem{Flatte:1976xu}\n"
" S.~M.~Flatte,\n"
" ``Coupled - Channel Analysis Of The Pi Eta And K Anti-K Systems Near K Anti-K\n"
" Threshold,''\n"
" Phys.\\ Lett.\\ B {\\bf 63}, 224 (1976).\n"
" %%CITATION = PHLTA,B63,224;%%\n"
);
static ParVector<ScalarMassGenerator,Energy> interfacecoupling
("Coupling",
"The coupling",
&ScalarMassGenerator::_coupling,
GeV, 0, ZERO, ZERO, Constants::MaxEnergy, false, false, true);
static ParVector<ScalarMassGenerator,Energy> interfacemass1
("Mass1",
"The mass for first particle",
&ScalarMassGenerator::_mass1,
GeV, 0, ZERO, ZERO, Constants::MaxEnergy, false, false, true);
static ParVector<ScalarMassGenerator,Energy> interfacemass2
("Mass2",
"The mass for second particle",
&ScalarMassGenerator::_mass2,
GeV, 0, ZERO, ZERO, Constants::MaxEnergy, false, false, true);
}
void ScalarMassGenerator::dataBaseOutput(ofstream & output,bool header) {
if(header) output << "update Mass_Generators set parameters=\"";
GenericMassGenerator::dataBaseOutput(output,false);
for(unsigned int ix=0;ix<_coupling.size();++ix)
output << "insert " << fullName() << ":Coupling "
<< ix << " " << _coupling[ix]/GeV << "\n";
for(unsigned int ix=0;ix<_mass1.size();++ix)
output << "insert " << fullName() << ":Mass1 "
<< ix << " " << _mass1[ix]/GeV << "\n";
for(unsigned int ix=0;ix<_mass2.size();++ix)
output << "insert " << fullName()
<< ":Mass2 " << ix << " " << _mass2[ix]/GeV << "\n";
if(header) output << "\n\" where BINARY ThePEGName=\""
<< fullName() << "\";" << endl;
}
void ScalarMassGenerator::doinit() {
if(_coupling.size()!=_mass1.size()||
_coupling.size()!=_mass2.size())
throw InitException() << "Parameter vectors have inconsistent sizes in "
<< "ScalarMassGenerator::doinit()"
<< Exception::runerror;
// initialise the local variables
for(unsigned int ix=0;ix<_mass1.size();++ix) {
_m2plus .push_back(sqr(_mass1[ix]+_mass2[ix]));
_m2minus.push_back(sqr(_mass1[ix]-_mass2[ix]));
}
// rest of the initialisation is handled in the base class
GenericMassGenerator::doinit();
}
double ScalarMassGenerator::weight(Energy q,int shape) const {
useMe();
Energy2 q2 = sqr(q);
Energy2 mass2 = sqr(nominalMass());
Energy2 mwidth= nominalMass()*nominalWidth();
Energy2 gamma[2]={0.*MeV2,ZERO};
if(shape==1) {
for(unsigned int ix=0;ix<_coupling.size();++ix) {
double lambda = (mass2-_m2plus[ix])*(mass2-_m2minus[ix])/sqr(mass2);
if(lambda>=0.) gamma[0]+=sqr(_coupling[ix])*sqrt( lambda)*q/nominalMass();
else gamma[1]+=sqr(_coupling[ix])*sqrt(-lambda)*q/nominalMass();
}
}
else {
for(unsigned int ix=0;ix<_coupling.size();++ix) {
double lambda = (q2-_m2plus[ix])*(q2-_m2minus[ix])/sqr(q2);
if(lambda>=0.) gamma[0]+=sqr(_coupling[ix])*sqrt( lambda);
else gamma[1]+=sqr(_coupling[ix])*sqrt(-lambda);
}
}
Energy2 numer(ZERO);
if(shape==2) numer = nominalMass()*gamma[0]/q;
else if(shape==3) numer = nominalMass()*nominalWidth();
else numer = gamma[0];
complex<Energy2> denom = (shape==2) ?
mass2 - q2 +mass2/q2*complex<Energy2>(gamma[1],-gamma[0]) :
mass2 - q2 +complex<Energy2>(gamma[1],-gamma[0]);
// complex denominantor
Energy4 den = real(denom*conj(denom));
return numer/den*(sqr(mass2-q2)+sqr(mwidth))/Constants::pi/mwidth;
}
InvEnergy2 ScalarMassGenerator::BreitWignerWeight(Energy q, int shape) const {
useMe();
Energy2 q2 = sqr(q);
Energy2 mass2 = sqr(nominalMass());
Energy2 gamma[2]={0.*MeV2,ZERO};
if(shape==1) {
for(unsigned int ix=0;ix<_coupling.size();++ix) {
double lambda = (mass2-_m2plus[ix])*(mass2-_m2minus[ix])/sqr(mass2);
if(lambda>=0.) gamma[0]+=sqr(_coupling[ix])*sqrt( lambda)*q/nominalMass();
else gamma[1]+=sqr(_coupling[ix])*sqrt(-lambda)*q/nominalMass();
}
}
else {
for(unsigned int ix=0;ix<_coupling.size();++ix) {
double lambda = (q2-_m2plus[ix])*(q2-_m2minus[ix])/sqr(q2);
if(lambda>=0.) gamma[0]+=sqr(_coupling[ix])*sqrt( lambda);
else gamma[1]+=sqr(_coupling[ix])*sqrt(-lambda);
}
}
Energy2 numer(ZERO);
if(shape==2) numer = nominalMass()*gamma[0]/q;
else if(shape==3) numer = nominalMass()*nominalWidth();
else numer = gamma[0];
complex<Energy2> denom = (shape==2) ?
mass2 - q2 +mass2/q2*complex<Energy2>(gamma[1],-gamma[0]) :
mass2 - q2 +complex<Energy2>(gamma[1],-gamma[0]);
// complex denominantor
Energy4 den = real(denom*conj(denom));
return numer/den/Constants::pi;
}

File Metadata

Mime Type
text/x-diff
Expires
Tue, Nov 19, 7:51 PM (1 d, 7 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
3792764
Default Alt Text
(33 KB)

Event Timeline