Page Menu
Home
HEPForge
Search
Configure Global Search
Log In
Files
F7877342
No One
Temporary
Actions
View File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Flag For Later
Size
17 KB
Subscribers
None
View Options
diff --git a/Models/General/WeakCurrentDecayConstructor.cc b/Models/General/WeakCurrentDecayConstructor.cc
--- a/Models/General/WeakCurrentDecayConstructor.cc
+++ b/Models/General/WeakCurrentDecayConstructor.cc
@@ -1,296 +1,295 @@
// -*- C++ -*-
//
// WeakCurrentDecayConstructor.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2019 The Herwig Collaboration
//
// Herwig is licenced under version 3 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 WeakCurrentDecayConstructor class.
//
#include "WeakCurrentDecayConstructor.h"
#include "ThePEG/Utilities/DescribeClass.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Interface/Parameter.h"
#include "ThePEG/Interface/ParVector.h"
#include "ThePEG/Interface/RefVector.h"
#include "ThePEG/Interface/Switch.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
#include "ThePEG/Repository/EventGenerator.h"
#include "ThePEG/PDT/DecayMode.h"
#include "ThePEG/Helicity/Vertex/AbstractFFVVertex.fh"
#include "DecayConstructor.h"
using namespace Herwig;
using ThePEG::Helicity::VertexBasePtr;
IBPtr WeakCurrentDecayConstructor::clone() const {
return new_ptr(*this);
}
IBPtr WeakCurrentDecayConstructor::fullclone() const {
return new_ptr(*this);
}
void WeakCurrentDecayConstructor::doinit() {
NBodyDecayConstructorBase::doinit();
- _theModel = dynamic_ptr_cast<Ptr<Herwig::StandardModel>::pointer>
- (generator()->standardModel());
+ model_ = dynamic_ptr_cast<Ptr<Herwig::StandardModel>::pointer>(generator()->standardModel());
unsigned int isize=decayTags_.size();
if(isize!=_norm .size()||isize!=_current.size())
throw InitException() << "Invalid sizes for the decay mode vectors in "
<< " WeakCurrentDecayConstructor "
<< decayTags_.size() << " " << _norm.size() << " "
<< _current.size() << Exception::runerror;
// get the particles from the tags
for(unsigned int ix=0;ix<decayTags_.size();++ix) {
_current[ix]->init();
particles_.push_back(vector<tPDPtr>());
string tag=decayTags_[ix];
do {
string::size_type next = min(tag.find(','), tag.find(';'));
particles_.back().push_back(generator()->findParticle(tag.substr(0,next)));
if(!particles_.back().back())
throw Exception() << "Failed to find particle " << tag.substr(0,next)
<< " in DecayMode " << decayTags_[ix]
<< " in WeakCurrentDecayConstructor::doinit()"
<< Exception::runerror;
if(tag[next]==';') break;
tag = tag.substr(next+1);
}
while(true);
}
}
void WeakCurrentDecayConstructor::persistentOutput(PersistentOStream & os) const {
- os << ounit(_masscut,GeV) << decayTags_ << particles_ << _norm << _current;
+ os << ounit(massCut_,GeV) << decayTags_ << particles_ << _norm << _current;
}
void WeakCurrentDecayConstructor::persistentInput(PersistentIStream & is, int) {
- is >> iunit(_masscut,GeV) >> decayTags_ >> particles_ >> _norm >> _current;
+ is >> iunit(massCut_,GeV) >> decayTags_ >> particles_ >> _norm >> _current;
}
// The following static variable is needed for the type
// description system in ThePEG.
DescribeClass<WeakCurrentDecayConstructor,NBodyDecayConstructorBase>
describeHerwigWeakCurrentDecayConstructor("Herwig::WeakCurrentDecayConstructor", "Herwig.so");
void WeakCurrentDecayConstructor::Init() {
static ClassDocumentation<WeakCurrentDecayConstructor> documentation
("The WeakCurrentDecayConstructor class implemets the decay of BSM particles "
"to low mass hadronic states using the Weak current");
static ParVector<WeakCurrentDecayConstructor,string> interfaceDecayModes
("DecayModes",
"The decays of the weak current",
&WeakCurrentDecayConstructor::decayTags_, -1, "", "", "",
false, false, Interface::nolimits);
static ParVector<WeakCurrentDecayConstructor,double> interfaceNormalisation
("Normalisation",
"The normalisation of the different modes",
&WeakCurrentDecayConstructor::_norm, -1, 1.0, 0.0, 10.0,
false, false, Interface::limited);
static RefVector<WeakCurrentDecayConstructor,WeakCurrent> interfaceCurrent
("Current",
"The current for the decay mode",
&WeakCurrentDecayConstructor::_current, -1, false, false, true, false, false);
static Parameter<WeakCurrentDecayConstructor,Energy> interfaceMassCut
("MassCut",
"The maximum mass difference for the decay",
- &WeakCurrentDecayConstructor::_masscut, GeV, 5.0*GeV, 1.0*GeV, 10.0*GeV,
+ &WeakCurrentDecayConstructor::massCut_, GeV, 5.0*GeV, 1.0*GeV, 10.0*GeV,
false, false, Interface::limited);
}
void WeakCurrentDecayConstructor::DecayList(const set<PDPtr> & part) {
if( part.empty() ) return;
- unsigned int nv(_theModel->numberOfVertices());
+ unsigned int nv(model_->numberOfVertices());
for(set<PDPtr>::const_iterator ip=part.begin();ip!=part.end();++ip) {
for(unsigned int iv = 0; iv < nv; ++iv) {
for(unsigned int ilist = 0; ilist < 3; ++ilist) {
vector<TwoBodyDecay> decays =
- createModes(*ip, _theModel->vertex(iv),ilist);
+ createModes(*ip, model_->vertex(iv),ilist);
if(!decays.empty()) createDecayMode(decays);
}
}
}
}
vector<TwoBodyDecay> WeakCurrentDecayConstructor::createModes(const PDPtr inpart,
const VertexBasePtr vert,
unsigned int ilist) {
int id = inpart->id();
if( !vert->isIncoming(inpart) || vert->getNpoint() != 3 )
return vector<TwoBodyDecay>();
Energy m1(inpart->mass());
vector<tPDPtr> decaylist;
decaylist = vert->search(ilist,inpart);
tPDVector::size_type nd = decaylist.size();
vector<TwoBodyDecay> decays;
for( tPDVector::size_type i = 0; i < nd; i += 3 ) {
tPDPtr pa(decaylist[i]), pb(decaylist.at(i + 1)),
pc(decaylist.at(i + 2));
if( pb->id() == id ) swap(pa, pb);
if( pc->id() == id ) swap(pa, pc);
//One of the products must be a W
Energy mp(ZERO);
if( abs(pb->id()) == ParticleID::Wplus )
mp = pc->mass();
else if( abs(pc->id()) == ParticleID::Wplus )
mp = pb->mass();
else
continue;
//allowed on-shell decay and passes mass cut
if( m1 >= pb->mass() + pc->mass() ) continue;
if( m1 < mp ) continue;
- if( m1 - mp >= _masscut ) continue;
+ if( m1 - mp >= massCut_ ) continue;
//vertices are defined with all particles incoming
if( pb->CC() ) pb = pb->CC();
if( pc->CC() ) pc = pc->CC();
decays.push_back( TwoBodyDecay(inpart,pb, pc, vert) );
if(abs(decays.back().children_.second->id())!=ParticleID::Wplus)
swap(decays.back().children_.first,decays.back().children_.second);
assert(abs(decays.back().children_.second->id())==ParticleID::Wplus);
}
return decays;
}
GeneralCurrentDecayerPtr WeakCurrentDecayConstructor::createDecayer(PDPtr in, PDPtr out1,
vector<tPDPtr> outCurrent,
VertexBasePtr vertex,
WeakCurrentPtr current) {
string name;
using namespace ThePEG::Helicity::VertexType;
switch(vertex->getName()) {
case FFV :
name = "FFVCurrentDecayer";
break;
default :
ostringstream message;
message << "Invalid vertex for decays of " << in->PDGName() << " -> " << out1->PDGName()
<< " via weak current " << vertex->fullName() << "\n";
generator()->logWarning(NBodyDecayConstructorError(message.str(),
Exception::warning));
return GeneralCurrentDecayerPtr();
}
ostringstream fullname;
fullname << "/Herwig/Decays/" << name << "_" << in->PDGName() << "_"
<< out1->PDGName();
for(unsigned int ix=0;ix<outCurrent.size();++ix)
fullname << "_" << outCurrent[ix]->PDGName();
string classname = "Herwig::" + name;
GeneralCurrentDecayerPtr decayer = dynamic_ptr_cast<GeneralCurrentDecayerPtr>
(generator()->preinitCreate(classname,fullname.str()));
- decayer->setDecayInfo(in,out1,outCurrent,vertex,current,_masscut);
+ decayer->setDecayInfo(in,out1,outCurrent,vertex,current,massCut_);
// set decayer options from base class
setDecayerInterfaces(fullname.str());
// initialize the decayer
decayer->init();
// return the decayer
return decayer;
}
void WeakCurrentDecayConstructor::
createDecayMode(vector<TwoBodyDecay> & decays) {
assert(!decays.empty());
for(unsigned int ix = 0; ix < decays.size(); ++ix) {
PDVector particles(3);
particles[0] = decays[ix].parent_;
particles[1] = decays[ix].children_.first ;
bool Wplus=decays[ix].children_.second->id()==ParticleID::Wplus;
for(unsigned int iy=0;iy<_current.size();++iy) {
particles.resize(2);
vector<tPDPtr> wprod=particles_[iy];
int icharge=0;
Energy msum = particles[0]->mass()-particles[1]->mass();
for(unsigned int iz=0;iz<wprod.size();++iz) {
icharge += wprod[iz]->iCharge();
msum -=wprod[iz]->mass();
}
if(msum<=ZERO) continue;
bool cc = (Wplus&&icharge==-3)||(!Wplus&&icharge==3);
OrderedParticles outgoing;
outgoing.insert(particles[1]);
for(unsigned int iz=0;iz<wprod.size();++iz) {
if(cc&&wprod[iz]->CC()) wprod[iz]=wprod[iz]->CC();
outgoing.insert(wprod[iz]);
}
// check outgoing particles initialised
for(unsigned int iz=0;iz<wprod.size();++iz) wprod[iz]->init();
// create the tag for the decay mode
string tag = particles[0]->PDGName() + "->";
OrderedParticles::const_iterator it = outgoing.begin();
do {
tag += (**it).name();
++it;
if(it!=outgoing.end()) tag +=",";
else tag +=";";
}
while(it!=outgoing.end());
// find the decay mode
tDMPtr dm= generator()->findDecayMode(tag);
if( !dm && createDecayModes() ) {
// create the decayer
GeneralCurrentDecayerPtr decayer = createDecayer(particles[0],particles[1],
wprod,decays[ix].vertex_,
_current[iy]);
if(!decayer) continue;
// calculate the width
Energy pWidth = _norm[iy]*decayer->partialWidth(particles[0],particles[1],wprod);
if(pWidth<=ZERO) {
generator()->preinitInterface(decayer->fullName(),
"Initialize", "set","0");
continue;
}
tDMPtr ndm = generator()->preinitCreateDecayMode(tag);
if(!ndm) throw NBodyDecayConstructorError()
<< "WeakCurrentDecayConstructor::createDecayMode - Needed to create "
<< "new decaymode but one could not be created for the tag "
<< tag
<< Exception::warning;
generator()->preinitInterface(ndm, "Decayer", "set",
decayer->fullName());
generator()->preinitInterface(ndm, "Active", "set", "Yes");
setBranchingRatio(ndm, pWidth);
particles[0]->stable(false);
if(ndm->brat()<decayConstructor()->minimumBR()) {
generator()->preinitInterface(decayer->fullName(),
"Initialize", "set","0");
}
}
else if (dm) {
// create the decayer
GeneralCurrentDecayerPtr decayer = createDecayer(particles[0],particles[1],
wprod,decays[ix].vertex_,
_current[iy]);
if(!decayer) continue;
generator()->preinitInterface(dm, "Decayer", "set", decayer->fullName());
particles[0]->stable(false);
if(createDecayModes()) {
// calculate the width
Energy pWidth = _norm[iy]*decayer->partialWidth(particles[0],particles[1],wprod);
if(pWidth<=ZERO) {
generator()->preinitInterface(decayer->fullName(),
"Initialize", "set","0");
continue;
}
generator()->preinitInterface(dm, "Active", "set", "Yes");
particles[0]->width(particles[0]->width()*(1.-dm->brat()));
setBranchingRatio(dm, pWidth);
}
if(dm->brat()<decayConstructor()->minimumBR()) {
generator()->preinitInterface(decayer->fullName(),
"Initialize", "set","0");
}
}
}
}
}
diff --git a/Models/General/WeakCurrentDecayConstructor.h b/Models/General/WeakCurrentDecayConstructor.h
--- a/Models/General/WeakCurrentDecayConstructor.h
+++ b/Models/General/WeakCurrentDecayConstructor.h
@@ -1,194 +1,194 @@
// -*- C++ -*-
//
// WeakCurrentDecayConstructor.h is a part of Herwig - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2019 The Herwig Collaboration
//
// Herwig is licenced under version 3 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
#ifndef HERWIG_WeakCurrentDecayConstructor_H
#define HERWIG_WeakCurrentDecayConstructor_H
//
// This is the declaration of the WeakCurrentDecayConstructor class.
//
#include "NBodyDecayConstructorBase.h"
#include "ThePEG/Helicity/Vertex/VertexBase.h"
#include "Herwig/Decay/General/GeneralCurrentDecayer.fh"
#include "Herwig/Models/StandardModel/StandardModel.h"
#include "Herwig/Decay/WeakCurrents/WeakCurrent.h"
#include "Herwig/Decay/General/GeneralCurrentDecayer.h"
#include "TwoBodyDecay.h"
namespace Herwig {
using namespace ThePEG;
/**
* Here is the documentation of the WeakCurrentDecayConstructor class.
*
* @see \ref WeakCurrentDecayConstructorInterfaces "The interfaces"
* defined for WeakCurrentDecayConstructor.
*/
class WeakCurrentDecayConstructor: public NBodyDecayConstructorBase {
public:
/**
* The default constructor.
*/
- WeakCurrentDecayConstructor() : _masscut(5.*GeV) {}
+ WeakCurrentDecayConstructor() : massCut_(5.*GeV) {}
/**
* Function used to determine allowed decaymodes, to be implemented
* in derived class.
*@param part vector of ParticleData pointers containing particles in model
*/
virtual void DecayList(const set<PDPtr> & part);
/**
* Number of outgoing lines. Required for correct ordering (do this one last)
*/
virtual unsigned int numBodies() const { return 1000; }
/**
* Cut off
*/
- Energy massCut() const { return _masscut;}
+ Energy massCut() const { return massCut_;}
public:
/** @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);
//@}
/**
* The standard Init function used to initialize the interfaces.
* Called exactly once for each class by the class description system
* before the main function starts or
* when this class is dynamically loaded.
*/
static void Init();
protected:
/** @name Clone Methods. */
//@{
/**
* Make a simple clone of this object.
* @return a pointer to the new object.
*/
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;
//@}
protected:
/** @name Standard Interfaced functions. */
//@{
/**
* Initialize this object after the setup phase before saving an
* EventGenerator to disk.
* @throws InitException if object could not be initialized properly.
*/
virtual void doinit();
//@}
private:
/** @name Functions to create decayers and decaymodes. */
//@{
/**
* Function to create decays
* @param inpart Incoming particle
* @param vert The vertex to create decays for
* @param ilist Which list to search
* @param iv Row number in _theExistingDecayers member
* @return vector of ParticleData ptrs
*/
vector<TwoBodyDecay>
createModes(const PDPtr inpart,const VertexBasePtr vert,
unsigned int ilist);
/**
* Function to create decayer for specific vertex
* @param vert Pointer to vertex
* @param icol Integer referring to the colmun in _theExistingDecayers
* @param ivert Integer referring to the row in _theExistingDecayers
* member variable
*/
GeneralCurrentDecayerPtr createDecayer(PDPtr in, PDPtr out1,
vector<tPDPtr> outCurrent,
VertexBasePtr vertex,
WeakCurrentPtr current);
/**
* Create decay mode(s) from given part and decay modes
* @param inpart pointer to incoming particle
* @param decays list of allowed interactions
* @param decayer The decayer responsible for this decay
*/
void createDecayMode(vector<TwoBodyDecay> & decays);
//@}
private:
/**
* The assignment operator is private and must never be called.
* In fact, it should not even be implemented.
*/
WeakCurrentDecayConstructor & operator=(const WeakCurrentDecayConstructor &) = delete;
private:
/**
* Model Pointer
*/
- Ptr<Herwig::StandardModel>::pointer _theModel;
+ Ptr<Herwig::StandardModel>::pointer model_;
/**
* Cut-off on the mass difference
*/
- Energy _masscut;
+ Energy massCut_;
/**
* Tags for the modes
*/
vector<string> decayTags_;
/**
* Particles for the mode
*/
vector<vector<tPDPtr> > particles_;
/**
* Normalisation
*/
vector<double> _norm;
/**
* The current for the mode
*/
vector<WeakCurrentPtr> _current;
};
}
#endif /* HERWIG_WeakCurrentDecayConstructor_H */
File Metadata
Details
Attached
Mime Type
text/x-diff
Expires
Tue, Nov 19, 3:19 PM (1 d, 14 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
3795742
Default Alt Text
(17 KB)
Attached To
R563 testingHerwigHG
Event Timeline
Log In to Comment