No OneTemporary
Actions

Size

45 KB

Subscribers

None

View Options

	diff --git a/Models/General/NBodyDecayConstructorBase.cc b/Models/General/NBodyDecayConstructorBase.cc
	--- a/Models/General/NBodyDecayConstructorBase.cc
	+++ b/Models/General/NBodyDecayConstructorBase.cc
	@@ -1,416 +1,428 @@
	// -- C++ --
	//
	// NBodyDecayConstructorBase.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 NBodyDecayConstructorBase class.
	//

	#include "NBodyDecayConstructorBase.h"
	#include "ThePEG/Repository/EventGenerator.h"
	#include "ThePEG/PDT/DecayMode.h"
	#include "ThePEG/Interface/Parameter.h"
	#include "ThePEG/Interface/Switch.h"
	#include "ThePEG/Interface/RefVector.h"
	#include "ThePEG/Interface/ClassDocumentation.h"
	#include "ThePEG/Persistency/PersistentOStream.h"
	#include "ThePEG/Persistency/PersistentIStream.h"
	#include "DecayConstructor.h"
	#include "Herwig++/Models/StandardModel/StandardModel.h"

	using namespace Herwig;
	using namespace ThePEG;

	void NBodyDecayConstructorBase::persistentOutput(PersistentOStream & os ) const {
	os << init_ << iteration_ << points_ << info_ << decayConstructor_
	- << removeOnShell_ << excludedVector_ << excludedSet_
	- << excludeEffective_ << includeTopOnShell_;
	+ << removeOnShell_ << excludeEffective_ << includeTopOnShell_
	+ << excludedVerticesVector_ << excludedVerticesSet_
	+ << excludedParticlesVector_ << excludedParticlesSet_;
	}

	void NBodyDecayConstructorBase::persistentInput(PersistentIStream & is , int) {
	is >> init_ >> iteration_ >> points_ >> info_ >> decayConstructor_
	- >> removeOnShell_ >> excludedVector_ >> excludedSet_
	- >> excludeEffective_ >> includeTopOnShell_;
	+ >> removeOnShell_ >> excludeEffective_ >> includeTopOnShell_
	+ >> excludedVerticesVector_ >> excludedVerticesSet_
	+ >> excludedParticlesVector_ >> excludedParticlesSet_;
	}

	AbstractClassDescription<NBodyDecayConstructorBase>
	NBodyDecayConstructorBase::initNBodyDecayConstructorBase;
	// Definition of the static class description member.

	void NBodyDecayConstructorBase::Init() {

	static ClassDocumentation<NBodyDecayConstructorBase> documentation
	("The NBodyDecayConstructorBase class is the base class for the automatic"
	"construction of the decay modes");

	static Switch<NBodyDecayConstructorBase,bool> interfaceInitializeDecayers
	("InitializeDecayers",
	"Initialize new decayers",
	&NBodyDecayConstructorBase::init_, true, false, false);
	static SwitchOption interfaceInitializeDecayersInitializeDecayersOn
	(interfaceInitializeDecayers,
	"Yes",
	"Initialize new decayers to find max weights",
	true);
	static SwitchOption interfaceInitializeDecayersoff
	(interfaceInitializeDecayers,
	"No",
	"Use supplied weights for integration",
	false);

	static Parameter<NBodyDecayConstructorBase,int> interfaceInitIteration
	("InitIteration",
	"Number of iterations to optimise integration weights",
	&NBodyDecayConstructorBase::iteration_, 1, 0, 10,
	false, false, true);

	static Parameter<NBodyDecayConstructorBase,int> interfaceInitPoints
	("InitPoints",
	"Number of points to generate when optimising integration",
	&NBodyDecayConstructorBase::points_, 1000, 100, 100000000,
	false, false, true);

	static Switch<NBodyDecayConstructorBase,bool> interfaceOutputInfo
	("OutputInfo",
	"Whether to output information about the decayers",
	&NBodyDecayConstructorBase::info_, false, false, false);
	static SwitchOption interfaceOutputInfoOff
	(interfaceOutputInfo,
	"No",
	"Do not output information regarding the created decayers",
	false);
	static SwitchOption interfaceOutputInfoOn
	(interfaceOutputInfo,
	"Yes",
	"Output information regarding the decayers",
	true);

	static Switch<NBodyDecayConstructorBase,bool> interfaceCreateDecayModes
	("CreateDecayModes",
	"Whether to create the ThePEG::DecayMode objects as well as the decayers",
	&NBodyDecayConstructorBase::createModes_, true, false, false);
	static SwitchOption interfaceCreateDecayModesOn
	(interfaceCreateDecayModes,
	"Yes",
	"Create the ThePEG::DecayMode objects",
	true);
	static SwitchOption interfaceCreateDecayModesOff
	(interfaceCreateDecayModes,
	"No",
	"Only create the Decayer objects",
	false);

	static Switch<NBodyDecayConstructorBase,unsigned int> interfaceRemoveOnShell
	("RemoveOnShell",
	"Remove on-shell diagrams as should be treated as a sequence of 1->2 decays",
	&NBodyDecayConstructorBase::removeOnShell_, 1, false, false);
	static SwitchOption interfaceRemoveOnShellYes
	(interfaceRemoveOnShell,
	"Yes",
	"Remove the diagrams if neither the production of decay or the intermediate"
	" can happen",
	1);
	static SwitchOption interfaceRemoveOnShellNo
	(interfaceRemoveOnShell,
	"No",
	"Never remove the intermediate",
	0);
	static SwitchOption interfaceRemoveOnShellProduction
	(interfaceRemoveOnShell,
	"Production",
	"Remove the diagram if the on-shell production of the intermediate is allowed",
	2);

	static RefVector<NBodyDecayConstructorBase,VertexBase> interfaceExcludedVertices
	("ExcludedVertices",
	"Vertices which are not included in the three-body decayers",
	- &NBodyDecayConstructorBase::excludedVector_, -1, false, false, true, true, false);
	+ &NBodyDecayConstructorBase::excludedVerticesVector_,
	+ -1, false, false, true, true, false);
	+
	+ static RefVector<NBodyDecayConstructorBase,ParticleData> interfaceExcludedIntermediates
	+ ("ExcludedIntermediates",
	+ "Excluded intermediate particles",
	+ &NBodyDecayConstructorBase::excludedParticlesVector_,
	+ -1, false, false, true, true, false);

	static Switch<NBodyDecayConstructorBase,bool> interfaceExcludeEffectiveVertices
	("ExcludeEffectiveVertices",
	"Exclude effectice vertices",
	&NBodyDecayConstructorBase::excludeEffective_, true, false, false);
	static SwitchOption interfaceExcludeEffectiveVerticesYes
	(interfaceExcludeEffectiveVertices,
	"Yes",
	"Exclude the effective vertices",
	true);
	static SwitchOption interfaceExcludeEffectiveVerticesNo
	(interfaceExcludeEffectiveVertices,
	"No",
	"Don't exclude the effective vertices",
	false);

	static Parameter<NBodyDecayConstructorBase,double> interfaceMinReleaseFraction
	("MinReleaseFraction",
	"The minimum energy release for a three-body decay, as a "
	"fraction of the parent mass.",
	&NBodyDecayConstructorBase::minReleaseFraction_, 1e-3, 0.0, 1.0,
	false, false, Interface::limited);

	static Switch<NBodyDecayConstructorBase,unsigned int> interfaceMaximumGaugeBosons
	("MaximumGaugeBosons",
	"Maximum number of electroweak gauge bosons"
	" to be produced as decay products",
	&NBodyDecayConstructorBase::maxBoson_, 1, false, false);
	static SwitchOption interfaceMaximumGaugeBosonsNone
	(interfaceMaximumGaugeBosons,
	"None",
	"Produce no W/Zs",
	0);
	static SwitchOption interfaceMaximumGaugeBosonsSingle
	(interfaceMaximumGaugeBosons,
	"Single",
	"Produce at most one W/Zs",
	1);
	static SwitchOption interfaceMaximumGaugeBosonsDouble
	(interfaceMaximumGaugeBosons,
	"Double",
	"Produce at most two W/Zs",
	2);
	static SwitchOption interfaceMaximumGaugeBosonsTriple
	(interfaceMaximumGaugeBosons,
	"Triple",
	"Produce at most three W/Zs",
	3);

	static Switch<NBodyDecayConstructorBase,unsigned int> interfaceMaximumNewParticles
	("MaximumNewParticles",
	"Maximum number of particles from the list of "
	"decaying particles to be allowed as decay products",
	&NBodyDecayConstructorBase::maxList_, 0, false, false);
	static SwitchOption interfaceMaximumNewParticlesNone
	(interfaceMaximumNewParticles,
	"None",
	"No particles from the list",
	0);
	static SwitchOption interfaceMaximumNewParticlesSingle
	(interfaceMaximumNewParticles,
	"Single",
	"A single particle from the list",
	1);
	static SwitchOption interfaceMaximumNewParticlesDouble
	(interfaceMaximumNewParticles,
	"Double",
	"Two particles from the list",
	2);
	static SwitchOption interfaceMaximumNewParticlesTriple
	(interfaceMaximumNewParticles,
	"Triple",
	"Four particles from the list",
	3);

	static Switch<NBodyDecayConstructorBase,bool> interfaceIncludeOnShellTop
	("IncludeOnShellTop",
	"Include the on-shell diagrams involving t -> bW",
	&NBodyDecayConstructorBase::includeTopOnShell_, false, false, false);
	static SwitchOption interfaceIncludeOnShellTopYes
	(interfaceIncludeOnShellTop,
	"Yes",
	"Inlude them",
	true);
	static SwitchOption interfaceIncludeOnShellTopNo
	(interfaceIncludeOnShellTop,
	"No",
	"Don't include them",
	true);
	}

	void NBodyDecayConstructorBase::setBranchingRatio(tDMPtr dm, Energy pwidth) {
	//Need width and branching ratios for all currently created decay modes
	PDPtr parent = const_ptr_cast<PDPtr>(dm->parent());
	DecaySet modes = parent->decayModes();
	if( modes.empty() ) return;
	double dmbrat(0.);
	if( modes.size() == 1 ) {
	parent->width(pwidth);
	if( pwidth > ZERO ) parent->cTau(hbarc/pwidth);
	dmbrat = 1.;
	}
	else {
	Energy currentwidth(parent->width());
	Energy newWidth(currentwidth + pwidth);
	parent->width(newWidth);
	if( newWidth > ZERO ) parent->cTau(hbarc/newWidth);
	//need to reweight current branching fractions if there are any
	double factor = newWidth > ZERO ? double(currentwidth/newWidth) : 0.;
	for(DecaySet::const_iterator dit = modes.begin();
	dit != modes.end(); ++dit) {
	if( *dit == dm \|\| !(**dit).on() ) continue;
	double newbrat = ((*dit).brat())factor;
	ostringstream brf;
	brf << setprecision(13)<< newbrat;
	generator()->preinitInterface(*dit, "BranchingRatio",
	"set", brf.str());
	}
	//set brat for current mode
	dmbrat = newWidth > ZERO ? double(pwidth/newWidth) : 0.;
	}
	ostringstream br;
	br << setprecision(13) << dmbrat;
	generator()->preinitInterface(dm, "BranchingRatio",
	"set", br.str());
	}

	void NBodyDecayConstructorBase::setDecayerInterfaces(string fullname) const {
	if( initialize() ) {
	ostringstream value;
	value << initialize();
	generator()->preinitInterface(fullname, "Initialize", "set",
	value.str());
	value.str("");
	value << iteration();
	generator()->preinitInterface(fullname, "Iteration", "set",
	value.str());
	value.str("");
	value << points();
	generator()->preinitInterface(fullname, "Points", "set",
	value.str());
	}
	// QED stuff if needed
	if(decayConstructor()->QEDGenerator())
	generator()->preinitInterface(fullname, "PhotonGenerator", "set",
	decayConstructor()->QEDGenerator()->fullName());
	string outputmodes;
	if( info() ) outputmodes = string("Output");
	else outputmodes = string("NoOutput");
	generator()->preinitInterface(fullname, "OutputModes", "set",
	outputmodes);
	}

	void NBodyDecayConstructorBase::doinit() {
	Interfaced::doinit();
	- excludedSet_ = set<VertexBasePtr>(excludedVector_.begin(),
	- excludedVector_.end());
	+ excludedVerticesSet_ = set<VertexBasePtr>(excludedVerticesVector_.begin(),
	+ excludedVerticesVector_.end());
	+ excludedParticlesSet_ = set<PDPtr>(excludedParticlesVector_.begin(),
	+ excludedParticlesVector_.end());
	if(removeOnShell_==0&&numBodies()>2)
	generator()->log() << "Warning: Including diagrams with on-shell "
	<< "intermediates in " << numBodies() << "-body BSM decays, this"
	<< " can lead to double counting and is not"
	<< " recommended unless you really know what you are doing\n"
	<< "This can be switched off using\n set "
	<< fullName() << ":RemoveOnShell Yes\n";
	}

	void NBodyDecayConstructorBase::DecayList(const set<PDPtr> & particles) {
	if( particles.empty() ) return;
	// cast the StandardModel to the Hw++ one to get the vertices
	tHwSMPtr model = dynamic_ptr_cast<tHwSMPtr>(generator()->standardModel());
	unsigned int nv(model->numberOfVertices());
	// loop over the particles and create the modes
	for(set<PDPtr>::const_iterator ip=particles.begin();
	ip!=particles.end();++ip) {
	// get the decaying particle
	tPDPtr parent = *ip;
	// first create prototype 1->2 decays
	std::queue<PrototypeVertexPtr> prototypes;
	for(unsigned int iv = 0; iv < nv; ++iv) {
	VertexBasePtr vertex = model->vertex(iv);
	if(excluded(vertex)) continue;
	PrototypeVertex::createPrototypes(parent, vertex, prototypes);
	}
	// now expand them into potential decay modes
	list<vector<PrototypeVertexPtr> > modes;
	while(!prototypes.empty()) {
	// get the first prototype from the queue
	PrototypeVertexPtr proto = prototypes.front();
	prototypes.pop();
	// multiplcity too low
	if(proto->npart<numBodies()) {
	// loop over all vertices and expand
	for(unsigned int iv = 0; iv < nv; ++iv) {
	VertexBasePtr vertex = model->vertex(iv);
	if(excluded(vertex)) continue;
	- PrototypeVertex::expandPrototypes(proto,vertex,prototypes);
	+ PrototypeVertex::expandPrototypes(proto,vertex,prototypes,
	+ excludedParticlesSet_);
	}
	}
	// multiplcity too high disgard
	else if(proto->npart>numBodies()) {
	continue;
	}
	// right multiplicity
	else {
	// check it's kinematical allowed for physical masses
	if( proto->incomingMass() < proto->outgoingMass() ) continue;
	// and for constituent masses
	Energy outgoingMass = proto->outgoingConstituentMass();
	if( proto->incomingMass() < proto->outgoingConstituentMass() ) continue;
	// remove processes which aren't kinematically allowed within
	// the release fraction
	if( proto->incomingMass() - outgoingMass <=
	minReleaseFraction_ * proto->incomingMass() ) continue;
	// check the external particles
	if(!proto->checkExternal()) continue;
	// check if first piece on-shell
	bool onShell = proto->canBeOnShell(removeOnShell(),proto->incomingMass(),true);
	// special treatment for three-body top decays
	if(onShell) {
	if(includeTopOnShell_ &&
	abs(proto->incoming->id())==ParticleID::t && proto->npart==3) {
	unsigned int nprimary=0;
	bool foundW=false, foundQ=false;
	for(OrderedVertices::const_iterator it = proto->outgoing.begin();
	it!=proto->outgoing.end();++it) {
	if(abs(it->first->id())==ParticleID::Wplus)
	foundW = true;
	if(abs(it->first->id())==ParticleID::b \|\|
	abs(it->first->id())==ParticleID::s \|\|
	abs(it->first->id())==ParticleID::d)
	foundQ = true;
	++nprimary;
	}
	if(!(foundW&&foundQ&&nprimary==2)) continue;
	}
	else continue;
	}
	// check if should be added to an existing decaymode
	bool added = false;
	for(list<vector<PrototypeVertexPtr> >::iterator it=modes.begin();
	it!=modes.end();++it) {
	// is the same ?
	if(!(it)[0]->sameDecay(proto)) continue;
	// it is the same
	added = true;
	// check if it is a duplicate
	bool already = false;
	for(unsigned int iz = 0; iz < (*it).size(); ++iz) {
	if( (it)[iz] == *proto) {
	already = true;
	break;
	}
	}
	if(!already) (*it).push_back(proto);
	break;
	}
	if(!added) modes.push_back(vector<PrototypeVertexPtr>(1,proto));
	}
	}
	// now look at the decay modes
	for(list<vector<PrototypeVertexPtr> >::iterator mit = modes.begin();
	mit!=modes.end();++mit) {
	// count the number of gauge bosons and particles from the list
	unsigned int nlist(0),nbos(0);
	for(OrderedParticles::const_iterator it=(*mit)[0]->outPart.begin();
	it!=(*mit)[0]->outPart.end();++it) {
	if(abs((**it).id()) == ParticleID::Wplus \|\|
	abs((**it).id()) == ParticleID::Z0) ++nbos;
	if(particles.find(*it)!=particles.end()) ++nlist;
	if((it).CC() && particles.find((it).CC())!=particles.end()) ++nlist;
	}
	// if too many ignore the mode
	if(nbos > maxBoson_ \|\| nlist > maxList_) continue;
	// create the decay
	createDecayMode(*mit);
	}
	}
	}

	void NBodyDecayConstructorBase::createDecayMode(vector<PrototypeVertexPtr> &) {
	throw Exception() << "In NBodyDecayConstructorBase::createDecayMode() which"
	<< " should have be overridden in an inheriting class"
	<< Exception::abortnow;
	}
	diff --git a/Models/General/NBodyDecayConstructorBase.h b/Models/General/NBodyDecayConstructorBase.h
	--- a/Models/General/NBodyDecayConstructorBase.h
	+++ b/Models/General/NBodyDecayConstructorBase.h
	@@ -1,324 +1,334 @@
	// -- C++ --
	//
	// NBodyDecayConstructorBase.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_NBodyDecayConstructorBase_H
	#define HERWIG_NBodyDecayConstructorBase_H
	//
	// This is the declaration of the NBodyDecayConstructorBase class.
	//

	#include "ThePEG/Interface/Interfaced.h"
	#include "ThePEG/Utilities/Exception.h"
	#include "ThePEG/PDT/ParticleData.h"
	#include "NBodyDecayConstructorBase.fh"
	#include "PrototypeVertex.h"
	#include "DecayConstructor.fh"

	namespace Herwig {

	using namespace ThePEG;

	/**
	* This is the base class for NBodyDecayConstructors. An N-body
	* decay constructor should inherit from this and implement the
	* DecayList virtual funtcion to create the decays and decayers.
	*
	* @see \ref NBodyDecayConstructorBaseInterfaces "The interfaces"
	* defined for NBodyDecayConstructor.
	*/
	class NBodyDecayConstructorBase: public Interfaced {

	public:

	/**
	* The default constructor.
	*/
	NBodyDecayConstructorBase() :
	init_(true),iteration_(1), points_(1000), info_(false),
	createModes_(true), removeOnShell_(1), excludeEffective_(true),
	minReleaseFraction_(1e-3), maxBoson_(1), maxList_(1),
	includeTopOnShell_(false) {}

	/**
	* Function used to determine allowed decaymodes, to be implemented
	* in derived class.
	* @param particles vector of ParticleData pointers containing
	* particles in model
	*/
	virtual void DecayList(const set<PDPtr> & particles);

	/**
	* Number of outgoing lines. Required for correct ordering.
	*/
	virtual unsigned int numBodies() const = 0;

	/**
	* Set the pointer to the DecayConstrcutor
	*/
	void decayConstructor(tDecayConstructorPtr d) {
	decayConstructor_ = d;
	}

	protected:

	/**
	* Method to set up the decay mode, should be overidden in inheriting class
	*/
	virtual void createDecayMode(vector<PrototypeVertexPtr> & mode);


	/**
	* Set the branching ratio of this mode. This requires
	* calculating a new width for the decaying particle and reweighting
	* the current branching fractions.
	* @param dm The decaymode for which to set the branching ratio
	* @param pwidth The calculated width of the mode
	*/
	void setBranchingRatio(tDMPtr dm, Energy pwidth);

	/**
	* Set the interfaces of the decayers depending on the flags stored.
	* @param name Fullname of the decayer in the EventGenerator
	* including the path
	*/
	void setDecayerInterfaces(string name) const;

	/**
	* Whether to initialize decayers or not
	*/
	bool initialize() const { return init_; }

	/**
	* Number of iterations if initializing (default 1)
	*/
	int iteration() const { return iteration_; }

	/**
	* Number of points to do in initialization
	*/
	int points() const { return points_; }

	/**
	* Whether to output information on the decayers
	*/
	bool info() const { return info_; }

	/**
	* Whether to create the DecayModes as well as the Decayer objects
	*/
	bool createDecayModes() const { return createModes_; }

	/**
	* Maximum number of electroweak gauge bosons
	*/
	unsigned int maximumGaugeBosons() const { return maxBoson_;}

	/**
	* Maximum number of particles from the list whose decays we are calculating
	*/
	unsigned int maximumList() const { return maxList_;}

	/**
	* Minimum energy release fraction
	*/
	double minimumReleaseFraction() const {return minReleaseFraction_;}

	/**
	* Get the pointer to the DecayConstructor object
	*/
	tDecayConstructorPtr decayConstructor() const {
	return decayConstructor_;
	}

	/**
	* Option for on-shell particles
	*/
	unsigned int removeOnShell() const { return removeOnShell_; }

	/**
	* Check if a vertex is excluded
	*/
	bool excluded(VertexBasePtr vertex) const {
	// skip an effective vertex
	if( excludeEffective_ &&
	int(vertex->orderInGs() + vertex->orderInGem()) != int(vertex->getNpoint())-2)
	return true;
	// check if explicitly forbidden
	- return excludedSet_.find(vertex)!=excludedSet_.end();
	+ return excludedVerticesSet_.find(vertex)!=excludedVerticesSet_.end();
	}

	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 after the setup phase before saving an
	* EventGenerator to disk.
	* @throws InitException if object could not be initialized properly.
	*/
	virtual void doinit();
	//@}

	private:

	/**
	* The static object used to initialize the description of this class.
	* Indicates that this is an abstract class with persistent data.
	*/
	static AbstractClassDescription<NBodyDecayConstructorBase>
	initNBodyDecayConstructorBase;

	/**
	* The assignment operator is private and must never be called.
	* In fact, it should not even be implemented.
	*/
	NBodyDecayConstructorBase & operator=(const NBodyDecayConstructorBase &);

	private:

	/**
	* Whether to initialize decayers or not
	*/
	bool init_;

	/**
	* Number of iterations if initializing (default 1)
	*/
	int iteration_;

	/**
	* Number of points to do in initialization
	*/
	int points_;

	/**
	* Whether to output information on the decayers
	*/
	bool info_;

	/**
	* Whether to create the DecayModes as well as the Decayer objects
	*/
	bool createModes_;

	/**
	* Whether or not to remove on-shell diagrams
	*/
	unsigned int removeOnShell_;

	/**
	* Excluded Vertices
	*/
	- vector<VertexBasePtr> excludedVector_;
	+ vector<VertexBasePtr> excludedVerticesVector_;

	/**
	* Excluded Vertices
	*/
	- set<VertexBasePtr> excludedSet_;
	+ set<VertexBasePtr> excludedVerticesSet_;
	+
	+ /**
	+ * Excluded Particles
	+ */
	+ vector<PDPtr> excludedParticlesVector_;
	+
	+ /**
	+ * Excluded Particles
	+ */
	+ set<PDPtr> excludedParticlesSet_;

	/**
	* Whether or not to exclude effective vertices
	*/
	bool excludeEffective_;

	/**
	* A pointer to the DecayConstructor object
	*/
	tDecayConstructorPtr decayConstructor_;

	/**
	* The minimum energy release for a three-body decay as a
	* fraction of the parent mass
	*/
	double minReleaseFraction_;

	/**
	* Maximum number of EW gauge bosons
	*/
	unsigned int maxBoson_;

	/**
	* Maximum number of particles from the decaying particle list
	*/
	unsigned int maxList_;

	/**
	* Include on-shell for \f$t\to b W\f$
	*/
	bool includeTopOnShell_;
	};

	/** An Exception class that can be used by all inheriting classes to
	* indicate a setup problem. */
	class NBodyDecayConstructorError : public Exception {

	public:

	NBodyDecayConstructorError() : Exception() {}

	NBodyDecayConstructorError(const string & str,
	Severity sev) : Exception(str,sev)
	{}
	};

	}

	#include "ThePEG/Utilities/ClassTraits.h"

	namespace ThePEG {

	/** @cond TRAITSPECIALIZATIONS */

	/** This template specialization informs ThePEG about the
	* base classes of NBodyDecayConstructorBase. */
	template <>
	struct BaseClassTrait<Herwig::NBodyDecayConstructorBase,1> {
	/** Typedef of the first base class of NBodyDecayConstructorBase. */
	typedef Interfaced NthBase;
	};

	/** This template specialization informs ThePEG about the name of
	* the NBodyDecayConstructorBase class and the shared object where it is defined. */
	template <>
	struct ClassTraits<Herwig::NBodyDecayConstructorBase>
	: public ClassTraitsBase<Herwig::NBodyDecayConstructorBase> {
	/** Return a platform-independent class name */
	static string className() { return "Herwig::NBodyDecayConstructorBase"; }
	};

	/** @endcond */

	}

	#endif /* HERWIG_NBodyDecayConstructorBase_H */
	diff --git a/Models/General/PrototypeVertex.cc b/Models/General/PrototypeVertex.cc
	--- a/Models/General/PrototypeVertex.cc
	+++ b/Models/General/PrototypeVertex.cc
	@@ -1,247 +1,251 @@
	// -- C++ --
	//
	// This is the implementation of the non-inlined, non-templated member
	// functions of the FourBodyDecayConstructor class.
	//

	#include "ThePEG/Helicity/Vertex/VertexBase.h"
	#include "ThePEG/PDT/StandardMatchers.h"
	#include "ThePEG/Utilities/Throw.h"
	#include "ThePEG/Utilities/Exception.h"
	#include "PrototypeVertex.h"

	using namespace Herwig;

	void PrototypeVertex::createPrototypes(tPDPtr inpart, VertexBasePtr vertex,
	std::queue<PrototypeVertexPtr> & prototypes) {
	int id = inpart->id();
	if(!vertex->isIncoming(inpart)) return;
	for(unsigned int list=0;list<vertex->getNpoint();++list) {
	tPDVector decaylist = vertex->search(list, inpart);
	tPDVector::size_type nd = decaylist.size();
	for( tPDVector::size_type i = 0; i < nd; i += vertex->getNpoint() ) {
	tPDVector pout(decaylist.begin()+i,
	decaylist.begin()+i+vertex->getNpoint());
	OrderedVertices out;
	for(unsigned int ix=1;ix<pout.size();++ix) {
	if(pout[ix]->id() == id ) swap(pout[0], pout[ix]);
	if(pout[ix]->CC()) pout[ix] = pout[ix]->CC();
	out.insert(make_pair(pout[ix],PrototypeVertexPtr()));
	}
	if(vertex->getNpoint()==3) {
	// remove radiation
	if((pout[0]==pout[1] && (pout[2]->id()==ParticleID::gamma\|\|
	pout[2]->id()==ParticleID::g\|\|
	pout[2]->id()==ParticleID::Z0)) \|\|
	(pout[0]==pout[2] && (pout[1]->id()==ParticleID::gamma\|\|
	pout[1]->id()==ParticleID::g\|\|
	pout[1]->id()==ParticleID::Z0)))
	continue;
	}
	prototypes.push(new_ptr(PrototypeVertex(inpart,out,vertex,
	int(vertex->getNpoint())-1)));
	}
	}
	}

	PrototypeVertexPtr PrototypeVertex::replicateTree(PrototypeVertexPtr parent,
	PrototypeVertexPtr oldChild,
	PrototypeVertexPtr & newChild) {
	PrototypeVertexPtr newParent =
	new_ptr(PrototypeVertex(parent->incoming,OrderedVertices(),
	parent->vertex,parent->npart));
	for(OrderedVertices::const_iterator it = parent->outgoing.begin();
	it!=parent->outgoing.end();++it) {
	PrototypeVertexPtr child = it->second ?
	replicateTree(it->second,oldChild,newChild) :
	PrototypeVertexPtr();
	newParent->outgoing.insert(make_pair(it->first,child));
	if(child) child->parent = newParent;
	if(it->second==oldChild) newChild=child;
	}
	if(oldChild==parent) newChild=newParent;
	return newParent;
	}

	void PrototypeVertex::expandPrototypes(PrototypeVertexPtr proto, VertexBasePtr vertex,
	- std::queue<PrototypeVertexPtr> & prototypes) {
	+ std::queue<PrototypeVertexPtr> & prototypes,
	+ const set<PDPtr> & excluded) {
	for(OrderedVertices::const_iterator it = proto->outgoing.begin();
	it!=proto->outgoing.end();++it) {
	if(it->second) {
	- expandPrototypes(it->second,vertex,prototypes);
	+ expandPrototypes(it->second,vertex,prototypes,excluded);
	}
	else {
	if(!vertex->isIncoming(it->first)) continue;
	+ if(excluded.find(it->first)!=excluded.end()) continue;
	+ if(it->first->CC() &&
	+ excluded.find(it->first->CC())!=excluded.end()) continue;
	int id = it->first->id();
	PrototypeVertexPtr parent=proto;
	while(parent->parent) parent=parent->parent;
	for(unsigned int il = 0; il < vertex->getNpoint(); ++il) {
	tPDVector decaylist = vertex->search(il,it->first );
	tPDVector::size_type nd = decaylist.size();
	for( tPDVector::size_type i = 0; i < nd; i += vertex->getNpoint() ) {
	tPDVector pout(decaylist.begin()+i,
	decaylist.begin()+i+vertex->getNpoint());
	OrderedVertices outgoing;
	for(unsigned int iy=1;iy<pout.size();++iy) {
	if(pout[iy]->id() == id ) swap(pout[0], pout[iy]);
	if(pout[iy]->CC()) pout[iy] = pout[iy]->CC();
	outgoing.insert(make_pair(pout[iy],PrototypeVertexPtr()));
	}
	if(vertex->getNpoint()==3) {
	if((pout[0]==pout[1] && (pout[2]->id()==ParticleID::gamma\|\|
	pout[2]->id()==ParticleID::g\|\|
	pout[2]->id()==ParticleID::Z0)) \|\|
	(pout[0]==pout[2] && (pout[1]->id()==ParticleID::gamma\|\|
	pout[1]->id()==ParticleID::g\|\|
	pout[1]->id()==ParticleID::Z0)))
	continue;
	// remove weak decays of quarks other than top
	if(StandardQCDPartonMatcher::Check(pout[0]->id()) &&
	((StandardQCDPartonMatcher::Check(pout[1]->id()) &&
	abs(pout[2]->id())==ParticleID::Wplus)\|\|
	(StandardQCDPartonMatcher::Check(pout[2]->id())&&
	abs(pout[1]->id())==ParticleID::Wplus))) continue;
	// remove weak decays of leptons
	if((abs(pout[0]->id())>=11&&abs(pout[0]->id())<=16) &&
	(((abs(pout[1]->id())>=11&&abs(pout[1]->id())<=16) &&
	abs(pout[2]->id())==ParticleID::Wplus)\|\|
	((abs(pout[2]->id())>=11&&abs(pout[2]->id())<=16)&&
	abs(pout[1]->id())==ParticleID::Wplus))) continue;
	}
	PrototypeVertexPtr newBranch =
	new_ptr(PrototypeVertex(it->first,
	outgoing,vertex,int(vertex->getNpoint())-1));
	PrototypeVertexPtr newChild;
	PrototypeVertexPtr newVertex = replicateTree(parent,proto,newChild);
	newBranch->parent = newChild;
	OrderedVertices::iterator kt = newChild->outgoing.begin();
	for(OrderedVertices::const_iterator jt = proto->outgoing.begin();
	jt!=it;++jt,++kt) {;}
	pair< tPDPtr, PrototypeVertexPtr > newPair = make_pair(kt->first,newBranch);
	newChild->outgoing.erase(kt);
	newChild->outgoing.insert(newPair);
	newChild->incrementN(newBranch->npart-1);
	prototypes.push(newVertex);
	}
	}
	}
	}
	}

	bool PrototypeVertex::canBeOnShell(unsigned int opt,Energy maxMass,bool first) {
	Energy outMass = outgoingMass();
	if(!first) {
	bool in = maxMass>incomingMass();
	bool out = incomingMass()>outMass;
	if(opt!=0) {
	if( in && ( out \|\| opt==2 ) ) return true;
	}
	else if (incoming->width() == ZERO) {
	tPrototypeVertexPtr original = this;
	while(original->parent) {
	original=original->parent;
	};
	ostringstream name;
	name << original->incoming->PDGName() << " -> ";
	for(OrderedParticles::const_iterator it = original->outPart.begin();
	it!= original->outPart.end();++it)
	name << (**it).PDGName() << " ";
	Throw<InitException>()
	<< "Trying to include on-shell diagram in decay"
	<< name.str() << "including on-shell "
	<< incoming->PDGName() << " with zero width.\n"
	<< "You should make sure that the width for the intermediate is either"
	<< " read from an SLHA file or the intermediate is included in the "
	<< "DecayParticles list of the ModelGenerator.\n"
	<< Exception::runerror;
	}
	}
	else maxMass = incomingMass();
	// check the decay products
	for(OrderedVertices::const_iterator it = outgoing.begin();
	it!=outgoing.end();++it) {
	if(!it->second) continue;
	Energy newMax = maxMass-outMass+it->second->outgoingMass();
	if(it->second->canBeOnShell(opt,newMax,false)) {
	if(first) possibleOnShell = true;
	return true;
	}
	}
	return false;
	}

	NBDiagram::NBDiagram(PrototypeVertexPtr proto)
	: incoming(proto->incoming), outgoing(proto->outPart),
	vertex(proto->vertex) {
	// create the vertices
	for(OrderedVertices::const_iterator it=proto->outgoing.begin();
	it!=proto->outgoing.end();++it) {
	vertices.push_back(make_pair(it->first,NBVertex(it->second)));
	}
	// now let's re-order so that branchings are at the end
	for(list<pair<PDPtr,NBVertex> >::iterator it=vertices.begin();
	it!=vertices.end();++it) {
	if(!it->second.incoming) continue;
	list<pair<PDPtr,NBVertex> >::iterator jt=it;
	for( ; jt!=vertices.end();++jt) {
	if(jt==it) continue;
	if(!jt->second.incoming) {
	break;
	}
	}
	if(jt!=vertices.end()) {
	list<pair<PDPtr,NBVertex> >::iterator kt = it;
	while(kt!=jt) {
	list<pair<PDPtr,NBVertex> >::iterator lt = kt;
	++lt;
	swap(kt,lt);
	kt=lt;
	}
	}
	}
	// finally work out the channel and the ordering
	unsigned int loc(outgoing.size());
	for(list<pair<PDPtr,NBVertex> >::iterator it=vertices.begin();
	it!=vertices.end();++it) {
	unsigned int ipart(0);
	for( ; ipart<outgoing.size(); ++ipart) {
	}


	--loc;
	}

	exit(0);
	}

	NBVertex::NBVertex(PrototypeVertexPtr proto) {
	if(!proto) return;
	incoming = proto->incoming;
	outgoing = proto->outPart;
	vertex = proto->vertex;
	for(OrderedVertices::const_iterator it=proto->outgoing.begin();
	it!=proto->outgoing.end();++it) {
	vertices.push_back(make_pair(it->first,NBVertex(it->second)));
	}
	// now let's re-order so that branchings are at the end
	for(list<pair<PDPtr,NBVertex> >::iterator it=vertices.begin();
	it!=vertices.end();++it) {
	if(!it->second.incoming) continue;
	list<pair<PDPtr,NBVertex> >::iterator jt=it;
	for( ; jt!=vertices.end();++jt) {
	if(jt==it) continue;
	if(!jt->second.incoming) {
	break;
	}
	}
	if(jt!=vertices.end()) {
	list<pair<PDPtr,NBVertex> >::iterator kt = it;
	while(kt!=jt) {
	list<pair<PDPtr,NBVertex> >::iterator lt = kt;
	++lt;
	swap(kt,lt);
	kt=lt;
	}
	}
	}
	}
	diff --git a/Models/General/PrototypeVertex.h b/Models/General/PrototypeVertex.h
	--- a/Models/General/PrototypeVertex.h
	+++ b/Models/General/PrototypeVertex.h
	@@ -1,439 +1,440 @@
	// -- C++ --
	//
	// PrototypeVertex.h is a part of Herwig++ - A multi-purpose Monte Carlo event generator
	// Copyright (C) 2002-2007 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_PrototypeVertex_H
	#define HERWIG_PrototypeVertex_H
	#include <queue>
	#include "ThePEG/Helicity/Vertex/VertexBase.h"
	#include "ThePEG/Persistency/PersistentOStream.h"
	#include "ThePEG/Persistency/PersistentIStream.h"
	#include "ThePEG/Utilities/EnumIO.h"
	//
	// This is the declaration of the PrototypeVertex class.
	//

	namespace Herwig {
	using namespace ThePEG;
	using Helicity::VertexBasePtr;

	class PrototypeVertex;
	ThePEG_DECLARE_POINTERS(Herwig::PrototypeVertex,PrototypeVertexPtr);

	/** Pair of int,double */
	typedef pair<unsigned int, double> CFPair;


	/**
	* A struct to order the particles in the same way as in the DecayMode's
	*/
	struct ParticleOrdering {
	/**
	* Operator for the ordering
	* @param p1 The first ParticleData object
	* @param p2 The second ParticleData object
	*/
	bool operator()(PDPtr p1, PDPtr p2) {
	return abs(p1->id()) > abs(p2->id()) \|\|
	( abs(p1->id()) == abs(p2->id()) && p1->id() > p2->id() ) \|\|
	( p1->id() == p2->id() && p1->fullName() > p2->fullName() );
	}
	};

	/**
	* A struct to order the particles in the same way as in the DecayMode's
	*/
	struct VertexOrdering {
	/**
	* Operator for the ordering
	* @param p1 The first ParticleData object
	* @param p2 The second ParticleData object
	*/
	bool operator()(pair< tPDPtr, PrototypeVertexPtr > p1,
	pair< tPDPtr, PrototypeVertexPtr > p2) {
	return abs(p1.first->id()) > abs(p2.first->id()) \|\|
	( abs(p1.first->id()) == abs(p2.first->id()) && p1.first->id() > p2.first->id() ) \|\|
	( p1.first->id() == p2.first->id() && p1.first->fullName() > p2.first->fullName() );
	}
	};

	typedef multiset<pair< tPDPtr, PrototypeVertexPtr >,VertexOrdering > OrderedVertices;

	/**
	* A set of ParticleData objects ordered as for the DecayMode's
	*/
	typedef multiset<PDPtr,ParticleOrdering> OrderedParticles;

	/**
	* Storage of a potenital n-body decay
	*/
	class PrototypeVertex : public Base {

	public:

	/**
	* Default Constructor
	*/
	PrototypeVertex() : npart(0), possibleOnShell(false) {}

	/**
	* Constructor
	*/
	PrototypeVertex(tPDPtr in, OrderedVertices out,
	VertexBasePtr v, int n) :
	incoming(in), outgoing(out), vertex(v), npart(n),
	possibleOnShell(false) {}

	/**
	* Incoming particle
	*/
	tPDPtr incoming;

	/**
	* Outgoing particles
	*/
	OrderedVertices outgoing;

	/**
	* The vertex for the interaction
	*/
	VertexBasePtr vertex;

	/**
	* The parent of the vertex
	*/
	tPrototypeVertexPtr parent;

	/**
	* Number of particles
	*/
	unsigned int npart;

	/**
	* Outgoing particles
	*/
	mutable OrderedParticles outPart;

	/**
	* Can have on-shell intermediates
	*/
	bool possibleOnShell;

	/**
	* Increment the number of particles
	*/
	void incrementN(int in) {
	npart += in;
	if(parent) parent->incrementN(in);
	}

	/**
	* Mass of the incoming particle
	*/
	Energy incomingMass() {
	return incoming->mass();
	}

	/**
	* Total mass of all the outgoing particles
	*/
	Energy outgoingMass() {
	Energy mass(ZERO);
	for(OrderedVertices::const_iterator it = outgoing.begin();
	it!=outgoing.end();++it) {
	mass += it->second ?
	it->second->outgoingMass() : it->first->mass();
	}
	return mass;
	}

	/**
	* Total constituent mass of all the outgoing particles
	*/
	Energy outgoingConstituentMass() {
	Energy mass(ZERO);
	for(OrderedVertices::const_iterator it = outgoing.begin();
	it!=outgoing.end();++it) {
	mass += it->second ?
	it->second->outgoingConstituentMass() : it->first->constituentMass();
	}
	return mass;
	}

	/**
	* Check the external particles
	*/
	bool checkExternal(bool first=true) {
	if(outPart.empty()) setOutgoing();
	if(first&&outPart.find(incoming)!=outPart.end()) return false;
	bool output = true;
	for(OrderedVertices::const_iterator it = outgoing.begin();
	it!=outgoing.end();++it) {
	if(it->second&& !it->second->checkExternal(false)) output = false;
	}
	return output;
	}

	/**
	* Set the outgoing particles
	*/
	void setOutgoing() const {
	assert(outPart.empty());
	for(OrderedVertices::const_iterator it = outgoing.begin();
	it!=outgoing.end();++it) {
	if(it->second) {
	it->second->setOutgoing();
	outPart.insert(it->second->outPart.begin(),
	it->second->outPart.end());
	}
	else
	outPart.insert(it->first);
	}
	}

	/**
	* Are there potential on-shell intermediates?
	*/
	bool canBeOnShell(unsigned int opt,Energy maxMass,bool first);

	/**
	* Check if same external particles
	*/
	bool sameDecay(const PrototypeVertex & x) const {
	if(incoming != x.incoming) return false;
	if(outPart.empty()) setOutgoing();
	if(x.outPart.empty()) x.setOutgoing();
	OrderedParticles::const_iterator cit = outPart.begin();
	OrderedParticles::const_iterator cjt = x.outPart.begin();
	if(x.npart!=npart) return false;
	for(;cit!=outPart.end();++cit,++cjt) {
	if(cit!=cjt) return false;
	}
	return true;
	}

	/**
	* Create a \f$1\to2\f$ prototype
	*/
	static void createPrototypes(tPDPtr inpart, VertexBasePtr vertex,
	std::queue<PrototypeVertexPtr> & prototypes);

	/**
	* Expand the prototypes by adding more legs
	*/
	static void expandPrototypes(PrototypeVertexPtr proto, VertexBasePtr vertex,
	- std::queue<PrototypeVertexPtr> & prototypes);
	+ std::queue<PrototypeVertexPtr> & prototypes,
	+ const set<PDPtr> & excluded);

	/**
	* Copy the whole structure with a new branching
	*/
	static PrototypeVertexPtr replicateTree(PrototypeVertexPtr parent,
	PrototypeVertexPtr oldChild,
	PrototypeVertexPtr & newChild);

	};

	/**
	* Output to a stream
	*/
	inline ostream & operator<<(ostream & os, const PrototypeVertex & diag) {
	os << diag.incoming->PDGName() << " -> ";
	bool seq=false;
	for(OrderedVertices::const_iterator it = diag.outgoing.begin();
	it!=diag.outgoing.end();++it) {
	os << it->first->PDGName() << " ";
	if(it->second) seq = true;
	}
	os << " decays via "
	<< diag.vertex->fullName() << " in a "
	<< diag.npart << "-body decay\n";
	if(!seq) return os;
	os << "Followed by\n";
	for(OrderedVertices::const_iterator it = diag.outgoing.begin();
	it!=diag.outgoing.end();++it) {
	if(it->second) os << *it->second;
	}
	return os;
	}

	/**
	* Test whether two diagrams are identical.
	*/
	inline bool operator==(const PrototypeVertex & x, const PrototypeVertex & y) {
	if(x.incoming != y.incoming) return false;
	if(x.vertex != y.vertex) return false;
	if(x.npart != y.npart) return false;
	if(x.outgoing.size() != y.outgoing.size()) return false;
	OrderedVertices::const_iterator xt = x.outgoing.begin();
	OrderedVertices::const_iterator yt = y.outgoing.begin();
	for(;xt!=x.outgoing.end();++xt,++yt) {
	if(xt->first != yt->first) return false;
	if(xt->second && yt->second) {
	if((xt->second)==(yt->second)) continue;
	else return false;
	}
	else if(xt->second \|\| yt->second)
	return false;
	}
	return true;
	}

	/**
	* A simple vertex for the N-body diagram
	*/
	struct NBVertex {

	/**
	* Constructor taking a prototype vertex as the arguments
	*/
	NBVertex(PrototypeVertexPtr proto = PrototypeVertexPtr() );

	/**
	* Incoming particle
	*/
	tPDPtr incoming;

	/**
	* Outgoing particles
	*/
	mutable OrderedParticles outgoing;

	/**
	* The vertices
	*/
	list<pair<PDPtr,NBVertex> > vertices;

	/**
	* The vertex
	*/
	VertexBasePtr vertex;
	};

	/**
	* The NBDiagram struct contains information about a \f$1\ton\f$ decay
	* that has been automatically generated.
	*/
	struct NBDiagram {

	/**
	* Enumeration for the channel type
	*/
	enum Channel {UNDEFINED = -1,
	TBchannel23=0, TBchannel13=1,
	TBchannel12=2, TBfourPoint=3};

	/**
	* Standard Constructor
	*/
	NBDiagram() : channelType(UNDEFINED) {}

	/**
	* Constructor taking a prototype vertex as the arguments*/
	NBDiagram(PrototypeVertexPtr proto);

	/**
	* Incoming particle
	*/
	tPDPtr incoming;

	/**
	* The type of channel
	*/
	Channel channelType;

	/**
	* Outgoing particles
	*/
	mutable OrderedParticles outgoing;

	/**
	* The vertices
	*/
	list<pair<PDPtr,NBVertex> > vertices;

	/**
	* The vertex for the parent branching
	*/
	VertexBasePtr vertex;

	/** Store colour flow at \f$N_c=3\f$ information */
	mutable vector<CFPair> colourFlow;

	/** Store colour flow at \f$N_c=\infty\f$ information */
	mutable vector<CFPair> largeNcColourFlow;
	};
	// /**
	// * Test whether this and x are the same decay
	// * @param x The other process to check
	// */
	// bool sameDecay(const TBDiagram & x) const {
	// if(ids[0] != x.ids[0]) return false;
	// bool used[4]={false,false,false,false};
	// for(unsigned int ix=1;ix<4;++ix) {
	// bool found=false;
	// for(unsigned int iy=1;iy<4;++iy) {
	// if(used[iy]) continue;
	// if(ids[ix]==x.ids[iy]) {
	// used[iy]=true;
	// found=true;
	// break;
	// }
	// }
	// if(!found) return false;
	// }
	// return true;
	// }

	/**
	* Output operator to allow the structure to be persistently written
	* @param os The output stream
	* @param x The TBVertex
	*/
	inline PersistentOStream & operator<<(PersistentOStream & os,
	const NBVertex & x) {
	os << x.incoming << x.outgoing << x.vertices << x.vertex;
	return os;
	}

	/**
	* Input operator to allow persistently written data to be read in
	* @param is The input stream
	* @param x The TBVertex
	*/
	inline PersistentIStream & operator>>(PersistentIStream & is,
	NBVertex & x) {
	is >> x.incoming >> x.outgoing >> x.vertices >> x.vertex;
	return is;
	}

	/**
	* Output operator to allow the structure to be persistently written
	* @param os The output stream
	* @param x The TBDiagram
	*/
	inline PersistentOStream & operator<<(PersistentOStream & os,
	const NBDiagram & x) {
	os << x.incoming << oenum(x.channelType) << x.outgoing << x.vertices << x.vertex
	<< x.colourFlow << x.largeNcColourFlow;
	return os;
	}

	/**
	* Input operator to allow persistently written data to be read in
	* @param is The input stream
	* @param x The TBDiagram
	*/
	inline PersistentIStream & operator>>(PersistentIStream & is,
	NBDiagram & x) {
	is >> x.incoming >> ienum(x.channelType) >> x.outgoing >> x.vertices >> x.vertex
	>> x.colourFlow >> x.largeNcColourFlow;
	return is;
	}

	}

	#endif /* HERWIG_PrototypeVertex_H */

File Metadata

Mime Type: text/x-diff
Expires: Wed, May 14, 10:04 AM (1 d, 13 h)
Storage Engine: blob
Storage Format: Raw Data
Storage Handle: 5104679
Default Alt Text: (45 KB)

No OneTemporaryActions

View Options

File Metadata

Event Timeline

No OneTemporary
Actions