No OneTemporary
Actions

Size

51 KB

Subscribers

None

View Options

	diff --git a/Decay/DecayIntegrator.cc b/Decay/DecayIntegrator.cc
	--- a/Decay/DecayIntegrator.cc
	+++ b/Decay/DecayIntegrator.cc
	@@ -1,302 +1,302 @@
	// -- C++ --
	//
	// DecayIntegrator.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 DecayIntegrator class.
	//

	#include "DecayIntegrator.h"
	#include "PhaseSpaceMode.h"
	#include "ThePEG/Interface/Reference.h"
	#include "ThePEG/Interface/Switch.h"
	#include "ThePEG/Interface/ClassDocumentation.h"
	#include "ThePEG/EventRecord/Particle.h"
	#include "ThePEG/Repository/UseRandom.h"
	#include "ThePEG/Repository/EventGenerator.h"
	#include "ThePEG/Utilities/DescribeClass.h"
	#include "ThePEG/Persistency/PersistentOStream.h"
	#include "ThePEG/Persistency/PersistentIStream.h"
	#include "Herwig/PDT/WidthCalculatorBase.h"

	using namespace Herwig;

	void DecayIntegrator::persistentOutput(PersistentOStream & os) const {
	os << modes_ << nIter_ << nPoint_ << nTry_
	<< photonGen_ << generateInter_ << ounit(eps_,GeV) << warnings_;
	}

	void DecayIntegrator::persistentInput(PersistentIStream & is, int) {
	is >> modes_ >> nIter_ >> nPoint_ >> nTry_
	>> photonGen_ >> generateInter_ >> iunit(eps_,GeV) >> warnings_;
	}

	// The following static variable is needed for the type
	// description system in ThePEG.
	DescribeAbstractClass<DecayIntegrator,HwDecayerBase>
	describeHerwigDecayIntegrator("Herwig::DecayIntegrator", "Herwig.so");

	void DecayIntegrator::Init() {
	static ClassDocumentation<DecayIntegrator> documentation
	("The DecayIntegrator class is a base decayer class "
	"including a multi-channel integrator.");

	static Parameter<DecayIntegrator,unsigned int> interfaceIteration
	("Iteration",
	"Number of iterations for the initialization of the phase space",
	&DecayIntegrator::nIter_, 10, 0, 100,
	false, false, true);

	static Parameter<DecayIntegrator,unsigned int> interfacePoints
	("Points",
	"number of phase space points to generate in the initialisation.",
	&DecayIntegrator::nPoint_, 10000, 1, 1000000000,
	false, false, true);

	static Parameter<DecayIntegrator,unsigned int> interfaceNtry
	("Ntry",
	"Number of attempts to generate the decay",
	&DecayIntegrator::nTry_, 500, 0, 100000,
	false, false, true);

	static Reference<DecayIntegrator,DecayRadiationGenerator> interfacePhotonGenerator
	("PhotonGenerator",
	"Object responsible for generating photons in the decay.",
	&DecayIntegrator::photonGen_, false, false, true, true, false);

	static Switch<DecayIntegrator,bool> interfaceGenerateIntermediates
	("GenerateIntermediates",
	"Whether or not to include intermediate particles in the output",
	&DecayIntegrator::generateInter_, false, false, false);
	static SwitchOption interfaceGenerateIntermediatesNoIntermediates
	(interfaceGenerateIntermediates,
	"No",
	"Don't include the intermediates",
	false);
	static SwitchOption interfaceGenerateIntermediatesIncludeIntermediates
	(interfaceGenerateIntermediates,
	"Yes",
	"include the intermediates",
	true);

	static Switch<DecayIntegrator, bool> InterfacePhaseSpaceWarning
	("PhaseSpaceWarning",
	"Switch on/off text warnings in PhaseSpaceMode class",
	- &DecayIntegrator::warnings_, true, false, false);
	+ &DecayIntegrator::warnings_, false, false, false);
	static SwitchOption on
	(InterfacePhaseSpaceWarning,"on","turn on the warnings", true);
	static SwitchOption off
	(InterfacePhaseSpaceWarning,"off","turn off the warnings", false);

	}

	double DecayIntegrator::oneLoopVirtualME(unsigned int ,
	const Particle &,
	const ParticleVector &) {
	throw Exception()
	<< "DecayIntegrator::oneLoopVirtualME() called. This should"
	<< " have been overidden in an inheriting class if it is used"
	<< Exception::runerror;
	}

	InvEnergy2 DecayIntegrator::realEmissionME(unsigned int,
	const Particle &,
	ParticleVector &,
	unsigned int,
	double, double,
	const LorentzRotation &,
	const LorentzRotation &) {
	throw Exception()
	<< "DecayIntegrator::realEmmisionME() called. This should"
	<< " have been overidden in an inheriting class if it is used"
	<< Exception::runerror;
	}

	ParticleVector DecayIntegrator::decay(const Particle & parent,
	const tPDVector & children) const {
	// return empty vector if products heavier than parent
	Energy mout(ZERO);
	for(tPDPtr pd : children) mout += pd->massMin();
	if(mout>parent.mass()) return ParticleVector();
	// generate the decay
	bool cc;
	iMode_ = modeNumber(cc,parent.dataPtr(),children);
	if(numberModes()==0) return ParticleVector();
	return modes_[iMode_]->generateDecay(parent,this,generateInter_,cc);
	}

	void DecayIntegrator::doinitrun() {
	HwDecayerBase::doinitrun();
	if ( initialize() && Debug::level > 1 )
	CurrentGenerator::current().log() << "Start of the initialisation for "
	<< name() << "\n";
	for(unsigned int ix=0;ix<modes_.size();++ix) {
	if(!modes_[ix]) continue;
	modes_[ix]->initrun();
	iMode_=ix;
	modes_[ix]->initializePhaseSpace(initialize(),this);
	}
	}

	// output the information for the database
	void DecayIntegrator::dataBaseOutput(ofstream & output,bool header) const {
	// header for MySQL
	if(header) output << "update decayers set parameters=\"";
	output << "newdef " << name() << ":Iteration " << nIter_ << "\n";
	output << "newdef " << name() << ":Ntry " << nTry_ << "\n";
	output << "newdef " << name() << ":Points " << nPoint_ << "\n";
	//if(_photongen){;}
	output << "newdef " << name() << ":GenerateIntermediates " << generateInter_ << " \n";
	// footer for MySQL
	if(header) {
	output << "\n\" where BINARY ThePEGName=\"" << fullName() << "\";\n";
	}
	}

	// set the code for the partial width
	void DecayIntegrator::setPartialWidth(const DecayMode & dm, int imode) {
	int ifound = findMode(dm);
	if(ifound>=0) modes_[ifound]->setPartialWidth(imode);
	}

	WidthCalculatorBasePtr
	DecayIntegrator::threeBodyMEIntegrator(const DecayMode &) const {
	return WidthCalculatorBasePtr();
	}

	int DecayIntegrator::findMode(const DecayMode & dm) {
	int imode(-1);
	vector<int> extid;
	bool found(false);
	int id;
	unsigned int ix(0),iy,N,iz,tmax,nmatched;
	if(modes_.size()==0) return -1;
	do {
	if(!modes_[ix]) {
	++ix;
	continue;
	}
	tcPDPtr in = modes_[ix]->incoming().first;
	tcPDPtr cc = modes_[ix]->incoming().first->CC();
	tmax=1;if(!cc){++tmax;}
	for(iz=0;iz<tmax;++iz) {
	extid.clear();
	// check the parent
	if(dm.parent()!=in && dm.parent()!=cc) continue;
	if(dm.parent()->id()==in->id()&&iz==0) {
	for(iy=0,N=modes_[ix]->numberOfParticles();iy<N;++iy) {
	extid.push_back(modes_[ix]->outgoing()[iy]->id());
	}
	}
	else if(dm.parent()->id()==in->id()&&iz==1) {
	for(iy=0,N=modes_[ix]->numberOfParticles();iy<N;++iy) {
	tcPDPtr cc2=modes_[ix]->outgoing()[iy]->CC();
	extid.push_back( cc2 ? cc2->id() : modes_[ix]->outgoing()[iy]->id());
	}
	}
	else if(cc&&dm.parent()->id()==cc->id()) {
	for(iy=0,N=modes_[ix]->numberOfParticles();iy<N;++iy) {
	tcPDPtr cc2 = modes_[ix]->outgoing()[iy]->CC();
	extid.push_back( cc2 ? cc2->id() : modes_[ix]->outgoing()[iy]->id());
	}
	}
	// if the parents match
	if(!extid.empty()) {
	vector<bool> matched(extid.size(),false);
	bool done;
	nmatched=0;
	ParticleMSet::const_iterator pit = dm.products().begin();
	do {
	id=(**pit).id();
	done=false;
	iy=0;
	do {
	if(id==extid[iy]&&!matched[iy]) {
	matched[iy]=true;
	++nmatched;
	done=true;
	}
	++iy;
	}
	while(iy<extid.size()&&!done);
	++pit;
	}
	while(pit!=dm.products().end());
	if(nmatched==extid.size()) {
	imode=ix;
	found=true;
	}
	}
	}
	++ix;
	}
	while(!found&&ix<modes_.size());
	return imode;
	}

	// the matrix element to be integrated for the me
	double DecayIntegrator::threeBodyMatrixElement(const int,const Energy2,
	const Energy2,
	const Energy2,const Energy2,
	const Energy, const Energy,
	const Energy) const {
	throw Exception()
	<< "Calling the virtual DecayIntegrator::threeBodyMatrixElement"
	<< "method. This must be overwritten in the classes "
	<< "inheriting from DecayIntegrator where it is needed"
	<< Exception::runerror;
	}

	// the differential three body decay rate with one integral performed
	InvEnergy DecayIntegrator::threeBodydGammads(const int, const Energy2,
	const Energy2,
	const Energy, const Energy,
	const Energy) const {
	throw Exception()
	<< "Calling the virtual DecayIntegrator::threeBodydGammads()"
	<<"method. This must be overwritten in the classes "
	<< "inheriting from DecayIntegrator where it is needed"
	<< Exception::runerror;
	}

	// generate the momenta for the decay
	ParticleVector DecayIntegrator::generate(bool inter,bool cc,
	const unsigned int & imode,
	const Particle & inpart) const {
	iMode_=imode;
	return modes_[imode]->generateDecay(inpart,this,inter,cc);
	}

	void DecayIntegrator::addMode(PhaseSpaceModePtr mode) const {
	modes_.push_back(mode);
	if(mode) mode->init();
	}

	ostream & Herwig::operator<<(ostream & os, const DecayIntegrator & decay) {
	os << "The integrator has " << decay.modes_.size() << " modes" << endl;
	for(unsigned int ix=0;ix<decay.modes_.size();++ix) {
	os << "Information on mode " << ix << endl;
	os << *(decay.modes_[ix]);
	}
	return os;
	}

	// reset the properities of all intermediates
	void DecayIntegrator::resetIntermediate(tcPDPtr part, Energy mass, Energy width) {
	if(!part) return;
	for(unsigned int ix=0,N=modes_.size();ix<N;++ix) {
	modes_[ix]->resetIntermediate(part,mass,width);
	}
	}

	Energy DecayIntegrator::initializePhaseSpaceMode(unsigned int imode,bool init, bool onShell) const{
	tcPhaseSpaceModePtr cmodeptr=mode(imode);
	tPhaseSpaceModePtr modeptr = const_ptr_cast<tPhaseSpaceModePtr>(cmodeptr);
	modeptr->init();
	return modeptr->initializePhaseSpace(init,this,onShell);
	}
	diff --git a/Decay/DecayIntegrator.h b/Decay/DecayIntegrator.h
	--- a/Decay/DecayIntegrator.h
	+++ b/Decay/DecayIntegrator.h
	@@ -1,526 +1,526 @@
	// -- C++ --
	//
	// DecayIntegrator.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_DecayIntegrator_H
	#define Herwig_DecayIntegrator_H
	//
	// This is the declaration of the DecayIntegrator class.
	//

	#include "DecayIntegrator.fh"
	#include "HwDecayerBase.h"
	#include "PhaseSpaceMode.fh"
	#include "Herwig/PDT/WidthCalculatorBase.fh"
	#include "Radiation/DecayRadiationGenerator.h"
	#include <Herwig/Decay/DecayVertex.h>

	namespace Herwig {

	using namespace ThePEG;

	/** \ingroup Decay
	* \class DecayIntegrator
	* \brief Main class for Decayers implementing multi-channel phase space integration.
	* \author Peter Richardson
	*
	* This class is designed to be the base class for Herwig decays including
	* the implementation of a multichannel decayer or n-body phase space decays.
	*
	* The <code>DecayIntegrator</code> class inherits from ThePEG's Decayer class
	* and makes use of the <code>PhaseSpaceMode</code> class to specify a number
	* of decay modes.
	*
	* Additional modes can be added using the addMode method. In practice the
	* phase space channels for a particular mode are usually constructed in the
	* doinit member of a Decayer and then the modes added to the Decayer.
	*
	* For the majority of the decays currently implemented the
	* phase-space integration has been optimised and the maximum weight set.
	* If the parameters of the decay model are changed the Initialize interface
	* can be used to optimise the integration and calculate the maximum weight.
	*
	* In classes inheriting from this the me2() member which gives the matrix element
	* squared must be implemented. This should be combined with the setting of the
	* phase space channels, and the setting of which channels to use and their
	* initial weights in the doinit() member. The different decay modes should then
	* be initialized in the initrun() member if needed. The generate member can then
	* be called from the decay() member to generate a phase-space configuration for a
	* decay.
	*
	* @see DecayPhaseSpaceMode
	* @see DecayPhaseSpaceChannel
	* @see \ref DecayIntegratorInterfaces "The interfaces"
	* defined for DecayIntegrator.
	*/
	class DecayIntegrator: public HwDecayerBase {

	public:

	/**
	* and DecayPhaseMode
	*/
	friend class PhaseSpaceMode;

	/**
	* Enum for the matrix element option
	*/
	enum MEOption {Initialize,Calculate,Terminate};

	public:

	/**
	* The default constructor.
	*/
	DecayIntegrator() : nIter_(10), nPoint_(10000), nTry_(500),
	generateInter_(false), iMode_(-1),
	- realME_(false), virtualME_(false), eps_(ZERO), warnings_(true)
	+ realME_(false), virtualME_(false), eps_(ZERO), warnings_(false)
	{}

	public:

	/**
	* Check if this decayer can perfom the decay for a particular mode.
	* Uses the modeNumber member but can be overridden
	* @param parent The decaying particle
	* @param children The decay products
	*/
	virtual bool accept(tcPDPtr parent, const tPDVector & children) const {
	bool cc;
	return modeNumber(cc,parent,children)>=0;
	}

	/**
	* For a given decay mode and a given particle instance, perform the
	* decay and return the decay products. As this is the base class this
	* is not implemented.
	* @return The vector of particles produced in the decay.
	*/
	virtual ParticleVector decay(const Particle & parent,
	const tPDVector & children) const;

	/**
	* Which of the possible decays is required
	* @param cc Is this mode the charge conjugate
	* @param parent The decaying particle
	* @param children The decay products
	*/
	virtual int modeNumber(bool & cc, tcPDPtr parent,
	const tPDVector & children) const = 0;

	/**
	* The mode being used for this decay
	*/
	int imode() const {return iMode_;}

	/**
	* Add a phase-space mode to the list
	* @param mode The mode being added.
	*/
	void addMode(PhaseSpaceModePtr mode) const;

	/**
	* Return the matrix element squared for a given mode and phase-space channel.
	* @param ichan The channel we are calculating the matrix element for.
	* @param part The decaying Particle.
	* @param outgoing The particles produced in the decay
	* @param momenta The momenta of the particles produced in the decay
	* @param meopt Option for the calculation of the matrix element
	* @return The matrix element squared for the phase-space configuration.
	*/
	virtual double me2(const int ichan, const Particle & part,
	const tPDVector & outgoing,
	const vector<Lorentz5Momentum> & momenta,
	MEOption meopt) const = 0;

	/**
	* Construct the SpinInfos for the particles produced in the decay
	*/
	virtual void constructSpinInfo(const Particle & part,
	ParticleVector outgoing) const = 0;

	/**
	* Output the setup information for the particle database
	* @param os The stream to output the information to
	* @param header Whether or not to output the information for MySQL
	*/
	virtual void dataBaseOutput(ofstream & os,bool header) const;

	/**
	* Set the code for the partial width. Finds the partial width in the
	* GenericWidthGenerator class which corresponds to the decay mode.
	* @param dm The DecayMode
	* @param imode The mode.
	*/
	void setPartialWidth(const DecayMode & dm, int imode);
	/**
	* Specify the \f$1\to2\f$ matrix element to be used in the running width calculation.
	* @param mecode The code for the matrix element as described
	* in the GenericWidthGenerator class.
	* @param coupling The coupling for the matrix element.
	* @return True or False if this mode can be handled.
	*/
	virtual bool twoBodyMEcode(const DecayMode &, int & mecode,
	double & coupling) const {
	coupling = 1.;
	mecode = -1;
	return false;
	}

	/**
	* Method to return an object to calculate the 3 (or higher body) partial width
	* @param dm The DecayMode
	* @return A pointer to a WidthCalculatorBase object capable of calculating the width
	*/
	virtual WidthCalculatorBasePtr threeBodyMEIntegrator(const DecayMode & dm) const;

	/**
	* The matrix element to be integrated for the three-body decays as a function
	* of the invariant masses of pairs of the outgoing particles.
	* @param imode The mode for which the matrix element is needed.
	* @param q2 The scale, \e i.e. the mass squared of the decaying particle.
	* @param s3 The invariant mass squared of particles 1 and 2, \f$s_3=m^2_{12}\f$.
	* @param s2 The invariant mass squared of particles 1 and 3, \f$s_2=m^2_{13}\f$.
	* @param s1 The invariant mass squared of particles 2 and 3, \f$s_1=m^2_{23}\f$.
	* @param m1 The mass of the first outgoing particle.
	* @param m2 The mass of the second outgoing particle.
	* @param m3 The mass of the third outgoing particle.
	* @return The matrix element
	*/
	virtual double threeBodyMatrixElement(const int imode, const Energy2 q2,
	const Energy2 s3, const Energy2 s2,
	const Energy2 s1, const Energy m1,
	const Energy m2, const Energy m3) const;

	/**
	* The differential three body decay rate with one integral performed.
	* @param imode The mode for which the matrix element is needed.
	* @param q2 The scale, \e i.e. the mass squared of the decaying particle.
	* @param s The invariant mass which still needs to be integrate over.
	* @param m1 The mass of the first outgoing particle.
	* @param m2 The mass of the second outgoing particle.
	* @param m3 The mass of the third outgoing particle.
	* @return The differential rate \f$\frac{d\Gamma}{ds}\f$
	*/
	virtual InvEnergy threeBodydGammads(const int imode, const Energy2 q2,
	const Energy2 s,
	const Energy m1, const Energy m2,
	const Energy m3) const;

	/**
	* Finds the phase-space mode corresponding to a given decay mode
	* @param dm The DecayMode
	*/
	int findMode(const DecayMode & dm);

	public:

	/**
	* Members for the generation of QED radiation in the decays
	*/
	//@{
	/**
	* Use the DecayRadiationGenerator to generate photons in the decay.
	* @param p The Particle instance being decayed
	* @param children The decay products
	* @return A particle vector containing the decay products after the generation
	* of photons.
	*/
	ParticleVector generatePhotons(const Particle & p,ParticleVector children) {
	return photonGen_->generatePhotons(p,children,this);
	}

	/**
	* check if photons can be generated in the decay
	*/
	bool canGeneratePhotons() {return photonGen_;}

	/**
	* The one-loop virtual correction.
	* @param imode The mode required.
	* @param part The decaying particle.
	* @param products The decay products including the radiated photon.
	* @return Whether the correction is implemented
	*/
	virtual double oneLoopVirtualME(unsigned int imode,
	const Particle & part,
	const ParticleVector & products);

	/**
	* Whether or not the one loop matrix element is implemented
	*/
	bool hasOneLoopME() {return virtualME_;}

	/**
	* The real emission matrix element
	* @param imode The mode required
	* @param part The decaying particle
	* @param products The decay products including the radiated photon
	* @param iemitter The particle which emitted the photon
	* @param ctheta The cosine of the polar angle between the photon and the
	* emitter
	* @param stheta The sine of the polar angle between the photon and the
	* emitter
	* @param rot1 Rotation from rest frame to frame for real emission
	* @param rot2 Rotation to place emitting particle along z
	*/
	virtual InvEnergy2 realEmissionME(unsigned int imode,
	const Particle & part,
	ParticleVector & products,
	unsigned int iemitter,
	double ctheta, double stheta,
	const LorentzRotation & rot1,
	const LorentzRotation & rot2);

	/**
	* Whether or not the real emission matrix element is implemented
	*/
	bool hasRealEmissionME() {return realME_;}
	//@}

	public:

	/**
	* The output operator is a friend, this is mainly for debugging
	*/
	friend ostream & operator<<(ostream & os, const DecayIntegrator & decay);


	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 Standard Interfaced functions. */
	//@{
	/**
	* Initialize this object. Called in the run phase just before
	* a run begins.
	*/
	virtual void doinitrun();
	//@}

	protected:

	/**
	* Generate the momenta for the decay
	* @param inter Generate the intermediates produced in the decay as well as the
	* final particles.
	* @param cc Is this the mode defined or its charge conjugate.
	* @param imode The mode being generated.
	* @param inpart The decaying particle.
	* @return The particles produced inthe decay.
	*/
	ParticleVector generate(bool inter,bool cc, const unsigned int & imode,
	const Particle & inpart) const;

	/**
	* Set the mode being use for this decay.
	*/
	void imode(int in) { iMode_ = in;}

	/**
	* Set the helicity matrix element for the decay.
	*/
	void ME(DecayMEPtr in) const { matrixElement_ = in;}

	/**
	* The helicity amplitude matrix element for spin correlations.
	*/
	DecayMEPtr ME() const {return matrixElement_;}

	/**
	* Reset the properities of all intermediates.
	* @param part The intermediate particle being reset.
	* @param mass The mass of the particle.
	* @param width The width of the particle.
	*/
	void resetIntermediate(tcPDPtr part, Energy mass, Energy width);

	/**
	* Initialize the phase-space mode
	* @param imode The mode
	* @param init Whether or not to perform the initialization
	*/
	Energy initializePhaseSpaceMode(unsigned int imode,bool init, bool onShell=false) const;

	protected:

	/**
	* Methods to set variables in inheriting classes
	*/
	//@{
	/**
	* Set whether or not the intermediates are included
	*/
	void generateIntermediates(bool in) {generateInter_=in;}

	/**
	* Set whether or not the intermediates are included
	*/
	bool generateIntermediates() const {return generateInter_;}

	/**
	* Whether or not the one loop matrix element is implemented
	*/
	void hasOneLoopME(bool in) {virtualME_=in;}

	/**
	* Whether or not the real emission matrix element is implemented
	*/
	void hasRealEmissionME(bool in) {realME_=in;}

	/**
	* Set the epsilon parameter
	*/
	void epsilonPS(Energy in) {eps_=in;}

	/**
	* Clear the models
	*/
	void clearModes() {modes_.clear();}

	protected:

	/**
	* Number of decay modes
	*/
	unsigned int numberModes() const {return modes_.size();}

	/**
	* Pointer to a mode
	*/
	tPhaseSpaceModePtr mode(unsigned int ix) {
	return modes_[ix];
	}
	/**
	* Pointer to a mode
	*/
	tcPhaseSpaceModePtr mode(unsigned int ix) const {
	return modes_[ix];
	}

	public:

	bool warnings() const {
	return warnings_;
	}

	private:

	/**
	* Private and non-existent assignment operator.
	*/
	DecayIntegrator & operator=(const DecayIntegrator &) = delete;

	/**
	* Parameters for the integration
	*/
	//@{
	/**
	* Number of iterations for th initialization.
	*/
	unsigned int nIter_;

	/**
	* Number of points for initialisation
	*/
	unsigned int nPoint_;

	/**
	* number of attempts to generate the decay
	*/
	unsigned int nTry_;

	/**
	* List of the decay modes
	*/
	mutable vector<PhaseSpaceModePtr> modes_;

	//@}

	/**
	* Whether to include the intermediates whne outputing the results.
	*/
	bool generateInter_;

	/**
	* Pointer to the object generating the QED radiation in the decay
	*/
	DecayRadiationGeneratorPtr photonGen_;

	/**
	* mode currently being generated
	*/
	mutable int iMode_;

	/**
	* The helicity matrix element for the current decay
	*/
	mutable DecayMEPtr matrixElement_;

	/**
	* Whether or not the real photon emission matrix element exists
	*/
	bool realME_;

	/**
	* Whether or not the one-loop matrix element exists
	*/
	bool virtualME_;

	/**
	* Epsilon parameter for phase-space integration
	*/
	Energy eps_;

	/**
	* option for turinh on/off log warnings in Phase class
	*/
	bool warnings_;

	protected:

	/**
	* Exception for this class and those inheriting from it
	*/
	class DecayIntegratorError: public Exception {};

	};

	/**
	* Output information on the DecayIntegrator for debugging purposes
	*/
	ostream & operator<<(ostream &, const DecayIntegrator &);

	}

	#endif /* Herwig_DecayIntegrator_H */
	diff --git a/src/defaults/Hadronization.in b/src/defaults/Hadronization.in
	--- a/src/defaults/Hadronization.in
	+++ b/src/defaults/Hadronization.in
	@@ -1,132 +1,133 @@
	# -- ThePEG-repository --

	############################################################
	# Setup of default hadronization
	#
	# There are no user servicable parts inside.
	#
	# Anything that follows below should only be touched if you
	# know what you're doing.
	#############################################################

	cd /Herwig/Particles
	create ThePEG::ParticleData Cluster
	setup Cluster 81 Cluster 0.00990 0.0 0.0 0.0 0 0 0 1
	create ThePEG::ParticleData Remnant
	setup Remnant 82 Remnant 0.00990 0.0 0.0 0.0 0 0 0 1

	mkdir /Herwig/Hadronization
	cd /Herwig/Hadronization

	create Herwig::ClusterHadronizationHandler ClusterHadHandler
	create Herwig::PartonSplitter PartonSplitter
	create Herwig::ClusterFinder ClusterFinder
	create Herwig::ColourReconnector ColourReconnector
	create Herwig::ClusterFissioner ClusterFissioner
	create Herwig::LightClusterDecayer LightClusterDecayer
	create Herwig::ClusterDecayer ClusterDecayer
	create Herwig::HwppSelector HadronSelector

	newdef ClusterHadHandler:PartonSplitter PartonSplitter
	newdef ClusterHadHandler:ClusterFinder ClusterFinder
	newdef ClusterHadHandler:ColourReconnector ColourReconnector
	newdef ClusterHadHandler:ClusterFissioner ClusterFissioner
	newdef ClusterHadHandler:LightClusterDecayer LightClusterDecayer
	newdef ClusterHadHandler:ClusterDecayer ClusterDecayer
	newdef ClusterHadHandler:MinVirtuality2 0.1*GeV2
	newdef ClusterHadHandler:MaxDisplacement 1.0e-10*millimeter
	newdef ClusterHadHandler:UnderlyingEventHandler NULL

	newdef ClusterFissioner:HadronSelector HadronSelector
	newdef LightClusterDecayer:HadronSelector HadronSelector
	newdef ClusterDecayer:HadronSelector HadronSelector
	newdef ClusterFinder:HadronSelector HadronSelector

	# ColourReconnector Default Parameters
	newdef ColourReconnector:ColourReconnection Yes
	newdef ColourReconnector:Algorithm Baryonic

	# Statistical CR Parameters:
	newdef ColourReconnector:AnnealingFactor 0.9
	newdef ColourReconnector:AnnealingSteps 50
	newdef ColourReconnector:TriesPerStepFactor 5.0
	newdef ColourReconnector:InitialTemperature 0.1

	# Plain and Baryonic CR Paramters
	newdef ColourReconnector:ReconnectionProbability 0.95
	newdef ColourReconnector:ReconnectionProbabilityBaryonic 0.7

	# BaryonicMesonic and BaryonicMesonic CR Paramters
	newdef ColourReconnector:ReconnectionProbability3Mto3M 0.5
	newdef ColourReconnector:ReconnectionProbability3MtoBBbar 0.5
	newdef ColourReconnector:ReconnectionProbabilityBbarBto3M 0.5
	newdef ColourReconnector:ReconnectionProbability2Bto2B 0.05
	newdef ColourReconnector:ReconnectionProbabilityMBtoMB 0.5
	newdef ColourReconnector:StepFactor 1.0
	newdef ColourReconnector:MesonToBaryonFactor 1.333

	# General Parameters and switches
	newdef ColourReconnector:MaxDistance 1.0e50
	newdef ColourReconnector:OctetTreatment All
	newdef ColourReconnector:CR2BeamClusters No
	newdef ColourReconnector:Junction Yes
	newdef ColourReconnector:PrePlainCR No
	newdef ColourReconnector:LocalCR No
	newdef ColourReconnector:CausalCR No
	# Debugging
	newdef ColourReconnector:Debug No

	# set ClusterFissioner parameters
	-set /Herwig/Hadronization/ClusterFissioner:ProbablityPowerFactor 1.65
	-set /Herwig/Hadronization/ClusterFissioner:ProbablityShift 0.00
	-set /Herwig/Hadronization/ClusterFissioner:KineticThresholdShift 2.75
	+set /Herwig/Hadronization/ClusterFissioner:KinematicThreshold Dynamic
	+set /Herwig/Hadronization/ClusterFissioner:ProbablityPowerFactor 3.3462
	+set /Herwig/Hadronization/ClusterFissioner:ProbablityShift 6.0442
	+set /Herwig/Hadronization/ClusterFissioner:KineticThresholdShift -1.5321

	# Clustering parameters for light quarks
	-newdef ClusterFissioner:ClMaxLight 3.649
	-newdef ClusterFissioner:ClPowLight 2.780
	-newdef ClusterFissioner:PSplitLight 0.899
	+newdef ClusterFissioner:ClMaxLight 3.2344
	+newdef ClusterFissioner:ClPowLight 2.6464
	+newdef ClusterFissioner:PSplitLight 0.7235
	newdef ClusterDecayer:ClDirLight 1
	newdef ClusterDecayer:ClSmrLight 0.78

	# Clustering parameters for b-quarks
	newdef ClusterFissioner:ClMaxBottom 3.757
	newdef ClusterFissioner:ClPowBottom 0.547
	newdef ClusterFissioner:PSplitBottom 0.625
	newdef ClusterDecayer:ClDirBottom 1
	newdef ClusterDecayer:ClSmrBottom 0.078
	newdef HadronSelector:SingleHadronLimitBottom 0.000

	# Clustering parameters for c-quarks
	newdef ClusterFissioner:ClMaxCharm 3.950
	newdef ClusterFissioner:ClPowCharm 2.559
	newdef ClusterFissioner:PSplitCharm 0.994
	newdef ClusterDecayer:ClDirCharm 1
	newdef ClusterDecayer:ClSmrCharm 0.163
	newdef HadronSelector:SingleHadronLimitCharm 0.000

	# Clustering parameters for exotic quarks
	# (e.g. hadronizing Susy particles)
	newdef ClusterFissioner:ClMaxExotic 2.7*GeV
	newdef ClusterFissioner:ClPowExotic 1.46
	newdef ClusterFissioner:PSplitExotic 1.00
	newdef ClusterDecayer:ClDirExotic 1
	newdef ClusterDecayer:ClSmrExotic 0.
	newdef HadronSelector:SingleHadronLimitExotic 0.

	#
	newdef PartonSplitter:SplitPwtSquark 0.824135
	newdef PartonSplitter:Split uds

	#
	newdef HadronSelector:PwtDquark 1.0
	newdef HadronSelector:PwtUquark 1.0
	-newdef HadronSelector:PwtSquark 0.291717
	+newdef HadronSelector:PwtSquark 0.35743
	newdef HadronSelector:PwtCquark 0.0
	newdef HadronSelector:PwtBquark 0.0
	-newdef HadronSelector:PwtDIquark 0.298
	-newdef HadronSelector:SngWt 0.74
	-newdef HadronSelector:DecWt 0.62
	+newdef HadronSelector:PwtDIquark 0.36452
	+newdef HadronSelector:SngWt 0.88022
	+newdef HadronSelector:DecWt 0.34622
	newdef HadronSelector:Mode 1
	newdef HadronSelector:BelowThreshold All

	create Herwig::SpinHadronizer SpinHadronizer
	diff --git a/src/defaults/Shower.in b/src/defaults/Shower.in
	--- a/src/defaults/Shower.in
	+++ b/src/defaults/Shower.in
	@@ -1,464 +1,464 @@
	# -- ThePEG-repository --

	############################################################
	# Setup of default parton shower
	#
	# Useful switches for users are marked near the top of
	# this file.
	#
	# Don't edit this file directly, but reset the switches
	# in your own input files!
	############################################################
	library HwMPI.so
	library HwShower.so
	library HwMatching.so

	mkdir /Herwig/Shower
	cd /Herwig/Shower

	create Herwig::QTildeShowerHandler ShowerHandler
	newdef ShowerHandler:MPIHandler /Herwig/UnderlyingEvent/MPIHandler
	newdef ShowerHandler:RemDecayer /Herwig/Partons/RemnantDecayer
	# use LO PDFs for Shower, can be changed later
	newdef ShowerHandler:PDFA /Herwig/Partons/ShowerLOPDF
	newdef ShowerHandler:PDFB /Herwig/Partons/ShowerLOPDF
	newdef ShowerHandler:PDFARemnant /Herwig/Partons/RemnantPDF
	newdef ShowerHandler:PDFBRemnant /Herwig/Partons/RemnantPDF

	#####################################
	# initial setup, don't change these!
	#####################################
	create Herwig::SplittingGenerator SplittingGenerator
	create Herwig::ShowerAlphaQCD AlphaQCD
	create Herwig::ShowerAlphaQCD AlphaQCDFSR
	create Herwig::ShowerAlphaQED AlphaQED
	set AlphaQED:CouplingSource Thompson
	create Herwig::ShowerAlphaQED AlphaEW
	set AlphaEW:CouplingSource MZ
	create Herwig::PartnerFinder PartnerFinder
	newdef PartnerFinder:PartnerMethod 0
	newdef PartnerFinder:ScaleChoice 0
	create Herwig::KinematicsReconstructor KinematicsReconstructor
	newdef KinematicsReconstructor:ReconstructionOption Colour3
	newdef KinematicsReconstructor:InitialStateReconOption SofterFraction
	newdef KinematicsReconstructor:InitialInitialBoostOption LongTransBoost
	newdef KinematicsReconstructor:FinalFinalWeight Yes

	newdef /Herwig/Partons/RemnantDecayer:AlphaS AlphaQCD
	newdef /Herwig/Partons/RemnantDecayer:AlphaEM AlphaQED

	newdef ShowerHandler:PartnerFinder PartnerFinder
	newdef ShowerHandler:KinematicsReconstructor KinematicsReconstructor
	newdef ShowerHandler:SplittingGenerator SplittingGenerator
	newdef ShowerHandler:Interactions QEDQCD
	newdef ShowerHandler:SpinCorrelations Yes
	newdef ShowerHandler:SoftCorrelations Singular

	##################################################################
	# Intrinsic pT
	#
	# Recommended:
	# 1.9 GeV for Tevatron W/Z production.
	# 2.1 GeV for LHC W/Z production at 10 TeV
	# 2.2 GeV for LHC W/Z production at 14 TeV
	#
	# Set all parameters to 0 to disable
	##################################################################
	newdef ShowerHandler:IntrinsicPtGaussian 2.2*GeV
	newdef ShowerHandler:IntrinsicPtBeta 0
	newdef ShowerHandler:IntrinsicPtGamma 0*GeV
	newdef ShowerHandler:IntrinsicPtIptmax 0*GeV
	#############################################################
	# Set up truncated shower handler.
	#############################################################

	create Herwig::PowhegShowerHandler PowhegShowerHandler
	set PowhegShowerHandler:MPIHandler /Herwig/UnderlyingEvent/MPIHandler
	set PowhegShowerHandler:RemDecayer /Herwig/Partons/RemnantDecayer
	newdef PowhegShowerHandler:PDFA /Herwig/Partons/ShowerLOPDF
	newdef PowhegShowerHandler:PDFB /Herwig/Partons/ShowerLOPDF
	newdef PowhegShowerHandler:PDFARemnant /Herwig/Partons/RemnantPDF
	newdef PowhegShowerHandler:PDFBRemnant /Herwig/Partons/RemnantPDF
	newdef PowhegShowerHandler:MPIHandler /Herwig/UnderlyingEvent/MPIHandler
	newdef PowhegShowerHandler:RemDecayer /Herwig/Partons/RemnantDecayer
	newdef PowhegShowerHandler:PDFA /Herwig/Partons/ShowerLOPDF
	newdef PowhegShowerHandler:PDFB /Herwig/Partons/ShowerLOPDF
	newdef PowhegShowerHandler:PDFARemnant /Herwig/Partons/RemnantPDF
	newdef PowhegShowerHandler:PDFBRemnant /Herwig/Partons/RemnantPDF
	newdef PowhegShowerHandler:PartnerFinder PartnerFinder
	newdef PowhegShowerHandler:KinematicsReconstructor KinematicsReconstructor
	newdef PowhegShowerHandler:SplittingGenerator SplittingGenerator
	newdef PowhegShowerHandler:Interactions QEDQCD
	newdef PowhegShowerHandler:SpinCorrelations Yes
	newdef PowhegShowerHandler:SoftCorrelations Singular
	newdef PowhegShowerHandler:IntrinsicPtGaussian 2.2*GeV
	newdef PowhegShowerHandler:IntrinsicPtBeta 0
	newdef PowhegShowerHandler:IntrinsicPtGamma 0*GeV
	newdef PowhegShowerHandler:IntrinsicPtIptmax 0*GeV
	newdef PowhegShowerHandler:EvolutionScheme DotProduct

	#############################################################
	# End of interesting user servicable section.
	#
	# Anything that follows below should only be touched if you
	# know what you're doing.
	#
	# Really.
	#############################################################
	#
	# a few default values
	newdef ShowerHandler:MECorrMode 1
	newdef ShowerHandler:EvolutionScheme DotProduct
	newdef AlphaQCD:ScaleFactor 1.0
	newdef AlphaQCD:NPAlphaS 2
	newdef AlphaQCD:Qmin 0.935
	newdef AlphaQCD:NumberOfLoops 2
	-newdef AlphaQCD:AlphaIn 0.126
	+newdef AlphaQCD:AlphaIn 0.1185
	newdef AlphaQCDFSR:ScaleFactor 1.0
	newdef AlphaQCDFSR:NPAlphaS 2
	newdef AlphaQCDFSR:Qmin 0.935
	newdef AlphaQCDFSR:NumberOfLoops 2
	newdef AlphaQCDFSR:AlphaIn 0.1186
	#
	#
	# Lets set up all the splittings
	create Herwig::HalfHalfOneSplitFn QtoQGammaSplitFn
	set QtoQGammaSplitFn:InteractionType QED
	set QtoQGammaSplitFn:ColourStructure ChargedChargedNeutral
	set QtoQGammaSplitFn:AngularOrdered Yes
	set QtoQGammaSplitFn:StrictAO Yes

	create Herwig::HalfHalfOneSplitFn QtoQGSplitFn
	newdef QtoQGSplitFn:InteractionType QCD
	newdef QtoQGSplitFn:ColourStructure TripletTripletOctet
	set QtoQGSplitFn:AngularOrdered Yes
	set QtoQGSplitFn:StrictAO Yes

	create Herwig::OneOneOneSplitFn GtoGGSplitFn
	newdef GtoGGSplitFn:InteractionType QCD
	newdef GtoGGSplitFn:ColourStructure OctetOctetOctet
	set GtoGGSplitFn:AngularOrdered Yes
	set GtoGGSplitFn:StrictAO Yes

	create Herwig::OneOneOneMassiveSplitFn WtoWGammaSplitFn
	newdef WtoWGammaSplitFn:InteractionType QED
	newdef WtoWGammaSplitFn:ColourStructure ChargedChargedNeutral
	set WtoWGammaSplitFn:AngularOrdered Yes
	set WtoWGammaSplitFn:StrictAO Yes

	create Herwig::OneHalfHalfSplitFn GtoQQbarSplitFn
	newdef GtoQQbarSplitFn:InteractionType QCD
	newdef GtoQQbarSplitFn:ColourStructure OctetTripletTriplet
	set GtoQQbarSplitFn:AngularOrdered Yes
	set GtoQQbarSplitFn:StrictAO Yes

	create Herwig::OneHalfHalfSplitFn GammatoQQbarSplitFn
	newdef GammatoQQbarSplitFn:InteractionType QED
	newdef GammatoQQbarSplitFn:ColourStructure NeutralChargedCharged
	set GammatoQQbarSplitFn:AngularOrdered Yes
	set GammatoQQbarSplitFn:StrictAO Yes

	create Herwig::HalfOneHalfSplitFn QtoGQSplitFn
	newdef QtoGQSplitFn:InteractionType QCD
	newdef QtoGQSplitFn:ColourStructure TripletOctetTriplet
	set QtoGQSplitFn:AngularOrdered Yes
	set QtoGQSplitFn:StrictAO Yes

	create Herwig::HalfOneHalfSplitFn QtoGammaQSplitFn
	newdef QtoGammaQSplitFn:InteractionType QED
	newdef QtoGammaQSplitFn:ColourStructure ChargedNeutralCharged
	set QtoGammaQSplitFn:AngularOrdered Yes
	set QtoGammaQSplitFn:StrictAO Yes

	create Herwig::HalfHalfOneEWSplitFn QtoQWZSplitFn
	newdef QtoQWZSplitFn:InteractionType EW
	newdef QtoQWZSplitFn:ColourStructure EW

	create Herwig::HalfHalfOneEWSplitFn LtoLWZSplitFn
	newdef LtoLWZSplitFn:InteractionType EW
	newdef LtoLWZSplitFn:ColourStructure EW

	create Herwig::HalfHalfZeroEWSplitFn QtoQHSplitFn
	newdef QtoQHSplitFn:InteractionType EW
	newdef QtoQHSplitFn:ColourStructure EW

	create Herwig::OneOneOneEWSplitFn VtoVVSplitFn
	newdef VtoVVSplitFn:InteractionType EW
	newdef VtoVVSplitFn:ColourStructure EW

	create Herwig::OneOneOneQEDSplitFn GammatoWWSplitFn
	newdef GammatoWWSplitFn:InteractionType QED
	newdef GammatoWWSplitFn:ColourStructure ChargedNeutralCharged

	create Herwig::OneOneZeroEWSplitFn VtoVHSplitFn
	newdef VtoVHSplitFn:InteractionType EW
	newdef VtoVHSplitFn:ColourStructure EW

	#
	# Now the Sudakovs
	create Herwig::PTCutOff PTCutOff
	-newdef PTCutOff:pTmin 0.958*GeV
	+newdef PTCutOff:pTmin 0.7473*GeV

	create Herwig::SudakovFormFactor SudakovCommon
	newdef SudakovCommon:Alpha AlphaQCD
	newdef SudakovCommon:Cutoff PTCutOff
	newdef SudakovCommon:PDFmax 1.0

	cp SudakovCommon QtoQGSudakov
	newdef QtoQGSudakov:SplittingFunction QtoQGSplitFn
	newdef QtoQGSudakov:PDFmax 1.9

	cp SudakovCommon QtoQGammaSudakov
	set QtoQGammaSudakov:SplittingFunction QtoQGammaSplitFn
	set QtoQGammaSudakov:Alpha AlphaQED
	set QtoQGammaSudakov:PDFmax 1.9

	cp QtoQGammaSudakov LtoLGammaSudakov
	cp PTCutOff LtoLGammaPTCutOff
	# Technical parameter to stop evolution.
	set LtoLGammaPTCutOff:pTmin 0.000001
	set LtoLGammaSudakov:Cutoff LtoLGammaPTCutOff

	cp SudakovCommon QtoQWZSudakov
	set QtoQWZSudakov:SplittingFunction QtoQWZSplitFn
	set QtoQWZSudakov:Alpha AlphaEW
	set QtoQWZSudakov:PDFmax 1.9

	cp SudakovCommon LtoLWZSudakov
	set LtoLWZSudakov:SplittingFunction LtoLWZSplitFn
	set LtoLWZSudakov:Alpha AlphaEW
	set LtoLWZSudakov:PDFmax 1.9

	cp SudakovCommon QtoQHSudakov
	set QtoQHSudakov:SplittingFunction QtoQHSplitFn
	set QtoQHSudakov:Alpha AlphaEW
	set QtoQHSudakov:PDFmax 1.9

	cp QtoQHSudakov LtoLHSudakov

	cp SudakovCommon VtoVVSudakov
	set VtoVVSudakov:SplittingFunction VtoVVSplitFn
	set VtoVVSudakov:Alpha AlphaQED

	cp SudakovCommon GammatoWWSudakov
	set GammatoWWSudakov:SplittingFunction GammatoWWSplitFn
	set GammatoWWSudakov:Alpha AlphaEW
	set GammatoWWSudakov:PDFmax 1.9



	cp SudakovCommon VtoVHSudakov
	set VtoVHSudakov:SplittingFunction VtoVHSplitFn
	set VtoVHSudakov:Alpha AlphaEW
	set VtoVHSudakov:PDFmax 1.9

	cp SudakovCommon GtoGGSudakov
	newdef GtoGGSudakov:SplittingFunction GtoGGSplitFn
	newdef GtoGGSudakov:PDFmax 2.0

	cp SudakovCommon WtoWGammaSudakov
	newdef WtoWGammaSudakov:SplittingFunction WtoWGammaSplitFn
	set WtoWGammaSudakov:Alpha AlphaQED

	cp SudakovCommon GtoQQbarSudakov
	newdef GtoQQbarSudakov:SplittingFunction GtoQQbarSplitFn
	newdef GtoQQbarSudakov:PDFmax 120.0

	cp SudakovCommon GammatoQQbarSudakov
	newdef GammatoQQbarSudakov:SplittingFunction GammatoQQbarSplitFn
	set GammatoQQbarSudakov:Alpha AlphaQED
	newdef GammatoQQbarSudakov:PDFmax 10000.0

	cp SudakovCommon GtobbbarSudakov
	newdef GtobbbarSudakov:SplittingFunction GtoQQbarSplitFn
	newdef GtobbbarSudakov:PDFmax 40000.0

	cp SudakovCommon GtoccbarSudakov
	newdef GtoccbarSudakov:SplittingFunction GtoQQbarSplitFn
	newdef GtoccbarSudakov:PDFmax 2000.0

	cp SudakovCommon QtoGQSudakov
	newdef QtoGQSudakov:SplittingFunction QtoGQSplitFn

	cp SudakovCommon QtoGammaQSudakov
	newdef QtoGammaQSudakov:SplittingFunction QtoGammaQSplitFn
	set QtoGammaQSudakov:Alpha AlphaQED
	newdef QtoGammaQSudakov:PDFmax 2000.0

	cp SudakovCommon utoGuSudakov
	newdef utoGuSudakov:SplittingFunction QtoGQSplitFn
	newdef utoGuSudakov:PDFFactor OverOneMinusZ
	newdef utoGuSudakov:PDFmax 5.0

	cp SudakovCommon dtoGdSudakov
	newdef dtoGdSudakov:SplittingFunction QtoGQSplitFn
	newdef dtoGdSudakov:PDFFactor OverOneMinusZ


	#Use a different Sudakov for FS QCD splittings in order to use a different value of alphaS
	cp QtoQGSudakov QtoQGSudakovFSR
	cp GtoGGSudakov GtoGGSudakovFSR
	cp GtoQQbarSudakov GtoQQbarSudakovFSR
	cp GtobbbarSudakov GtobbbarSudakovFSR
	cp GtobbbarSudakov GtoccbarSudakovFSR

	set QtoQGSudakovFSR:Alpha AlphaQCDFSR
	set GtoGGSudakovFSR:Alpha AlphaQCDFSR
	set GtoQQbarSudakovFSR:Alpha AlphaQCDFSR
	set GtobbbarSudakovFSR:Alpha AlphaQCDFSR
	set GtoccbarSudakovFSR:Alpha AlphaQCDFSR

	#
	# Now add the final splittings
	#
	do SplittingGenerator:AddFinalSplitting u->u,g; QtoQGSudakovFSR
	do SplittingGenerator:AddFinalSplitting d->d,g; QtoQGSudakovFSR
	do SplittingGenerator:AddFinalSplitting s->s,g; QtoQGSudakovFSR
	do SplittingGenerator:AddFinalSplitting c->c,g; QtoQGSudakovFSR
	do SplittingGenerator:AddFinalSplitting b->b,g; QtoQGSudakovFSR
	do SplittingGenerator:AddFinalSplitting t->t,g; QtoQGSudakovFSR
	#
	do SplittingGenerator:AddFinalSplitting g->g,g; GtoGGSudakovFSR
	#
	do SplittingGenerator:AddFinalSplitting g->u,ubar; GtoQQbarSudakovFSR
	do SplittingGenerator:AddFinalSplitting g->d,dbar; GtoQQbarSudakovFSR
	do SplittingGenerator:AddFinalSplitting g->s,sbar; GtoQQbarSudakovFSR
	do SplittingGenerator:AddFinalSplitting g->c,cbar; GtoccbarSudakovFSR
	do SplittingGenerator:AddFinalSplitting g->b,bbar; GtobbbarSudakovFSR
	do SplittingGenerator:AddFinalSplitting g->t,tbar; GtoQQbarSudakovFSR
	#
	do SplittingGenerator:AddFinalSplitting gamma->u,ubar; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->d,dbar; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->s,sbar; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->c,cbar; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->b,bbar; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->t,tbar; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->e-,e+; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->mu-,mu+; GammatoQQbarSudakov
	do SplittingGenerator:AddFinalSplitting gamma->tau-,tau+; GammatoQQbarSudakov
	#
	do SplittingGenerator:AddFinalSplitting u->u,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddFinalSplitting d->d,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddFinalSplitting s->s,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddFinalSplitting c->c,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddFinalSplitting b->b,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddFinalSplitting t->t,gamma; QtoQGammaSudakov

	do SplittingGenerator:AddFinalSplitting e-->e-,gamma; LtoLGammaSudakov
	do SplittingGenerator:AddFinalSplitting mu-->mu-,gamma; LtoLGammaSudakov
	do SplittingGenerator:AddFinalSplitting tau-->tau-,gamma; LtoLGammaSudakov

	do SplittingGenerator:AddFinalSplitting W+->W+,gamma; WtoWGammaSudakov
	#
	# Now lets add the initial splittings. Remember the form a->b,c; means
	# that the current particle b is given and we backward branch to new
	# particle a which is initial state and new particle c which is final state
	#
	do SplittingGenerator:AddInitialSplitting u->u,g; QtoQGSudakov
	do SplittingGenerator:AddInitialSplitting d->d,g; QtoQGSudakov
	do SplittingGenerator:AddInitialSplitting s->s,g; QtoQGSudakov
	do SplittingGenerator:AddInitialSplitting c->c,g; QtoQGSudakov
	do SplittingGenerator:AddInitialSplitting b->b,g; QtoQGSudakov
	do SplittingGenerator:AddInitialSplitting u->u,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddInitialSplitting d->d,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddInitialSplitting s->s,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddInitialSplitting c->c,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddInitialSplitting b->b,gamma; QtoQGammaSudakov
	do SplittingGenerator:AddInitialSplitting t->t,gamma; QtoQGammaSudakov

	do SplittingGenerator:AddInitialSplitting g->g,g; GtoGGSudakov
	#
	do SplittingGenerator:AddInitialSplitting g->d,dbar; GtoQQbarSudakov
	do SplittingGenerator:AddInitialSplitting g->u,ubar; GtoQQbarSudakov
	do SplittingGenerator:AddInitialSplitting g->s,sbar; GtoQQbarSudakov
	do SplittingGenerator:AddInitialSplitting g->c,cbar; GtoccbarSudakov
	do SplittingGenerator:AddInitialSplitting g->b,bbar; GtobbbarSudakov
	#
	#do SplittingGenerator:AddInitialSplitting gamma->d,dbar; GammatoQQbarSudakov
	#do SplittingGenerator:AddInitialSplitting gamma->u,ubar; GammatoQQbarSudakov
	#do SplittingGenerator:AddInitialSplitting gamma->s,sbar; GammatoQQbarSudakov
	#do SplittingGenerator:AddInitialSplitting gamma->c,cbar; GammatoQQbarSudakov
	#do SplittingGenerator:AddInitialSplitting gamma->b,bbar; GammatoQQbarSudakov
	do SplittingGenerator:AddInitialSplitting gamma->e-,e+; GammatoQQbarSudakov
	do SplittingGenerator:AddInitialSplitting gamma->mu-,mu+; GammatoQQbarSudakov
	do SplittingGenerator:AddInitialSplitting gamma->tau-,tau+; GammatoQQbarSudakov
	#
	do SplittingGenerator:AddInitialSplitting d->g,d; dtoGdSudakov
	do SplittingGenerator:AddInitialSplitting u->g,u; utoGuSudakov
	do SplittingGenerator:AddInitialSplitting s->g,s; QtoGQSudakov
	do SplittingGenerator:AddInitialSplitting c->g,c; QtoGQSudakov
	do SplittingGenerator:AddInitialSplitting b->g,b; QtoGQSudakov
	do SplittingGenerator:AddInitialSplitting dbar->g,dbar; dtoGdSudakov
	do SplittingGenerator:AddInitialSplitting ubar->g,ubar; utoGuSudakov
	do SplittingGenerator:AddInitialSplitting sbar->g,sbar; QtoGQSudakov
	do SplittingGenerator:AddInitialSplitting cbar->g,cbar; QtoGQSudakov
	do SplittingGenerator:AddInitialSplitting bbar->g,bbar; QtoGQSudakov
	#
	do SplittingGenerator:AddInitialSplitting d->gamma,d; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting u->gamma,u; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting s->gamma,s; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting c->gamma,c; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting b->gamma,b; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting dbar->gamma,dbar; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting ubar->gamma,ubar; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting sbar->gamma,sbar; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting cbar->gamma,cbar; QtoGammaQSudakov
	do SplittingGenerator:AddInitialSplitting bbar->gamma,bbar; QtoGammaQSudakov

	#
	# Electroweak
	#
	do SplittingGenerator:AddFinalSplitting u->u,Z0; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting d->d,Z0; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting s->s,Z0; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting c->c,Z0; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting b->b,Z0; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting t->t,Z0; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting u->u,Z0; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting d->d,Z0; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting s->s,Z0; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting c->c,Z0; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting b->b,Z0; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting t->t,Z0; QtoQWZSudakov

	do SplittingGenerator:AddFinalSplitting u->d,W+; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting c->s,W+; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting d->u,W-; QtoQWZSudakov
	do SplittingGenerator:AddFinalSplitting s->c,W-; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting u->d,W+; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting c->s,W+; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting d->u,W-; QtoQWZSudakov
	do SplittingGenerator:AddInitialSplitting s->c,W-; QtoQWZSudakov

	do SplittingGenerator:AddFinalSplitting c->c,h0; QtoQHSudakov
	do SplittingGenerator:AddFinalSplitting b->b,h0; QtoQHSudakov
	do SplittingGenerator:AddFinalSplitting t->t,h0; QtoQHSudakov
	do SplittingGenerator:AddInitialSplitting c->c,h0; QtoQHSudakov
	do SplittingGenerator:AddInitialSplitting b->b,h0; QtoQHSudakov
	do SplittingGenerator:AddInitialSplitting t->t,h0; QtoQHSudakov

	do SplittingGenerator:AddFinalSplitting gamma->W+,W-; GammatoWWSudakov

	do SplittingGenerator:AddFinalSplitting Z0->W+,W-; VtoVVSudakov
	do SplittingGenerator:AddFinalSplitting W+->W+,gamma; VtoVVSudakov
	do SplittingGenerator:AddFinalSplitting W+->W+,Z0; VtoVVSudakov

	do SplittingGenerator:AddFinalSplitting W+->W+,h0; VtoVHSudakov
	do SplittingGenerator:AddFinalSplitting Z0->Z0,h0; VtoVHSudakov

	#
	# Electroweak l -> l V
	#
	#do SplittingGenerator:AddFinalSplitting e-->e-,Z0; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting mu-->mu-,Z0; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting tau-->tau-,Z0; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting nu_e->nu_e,Z0; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting nu_mu->nu_mu,Z0; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting nu_tau->nu_tau,Z0; LtoLWZSudakov

	#do SplittingGenerator:AddFinalSplitting e-->nu_e,W-; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting mu-->nu_mu,W-; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting tau-->nu_tau,W-; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting nu_e->e-,W+; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting nu_mu->mu-,W+; LtoLWZSudakov
	#do SplittingGenerator:AddFinalSplitting nu_tau->tau-,W+; LtoLWZSudakov

File Metadata

Mime Type: text/x-diff
Expires: Wed, May 14, 10:18 AM (1 d, 17 h)
Storage Engine: blob
Storage Format: Raw Data
Storage Handle: 5091205
Default Alt Text: (51 KB)

No OneTemporaryActions

View Options

File Metadata

Event Timeline

No OneTemporary
Actions