No OneTemporary
Actions

Size

58 KB

Subscribers

None

View Options

	diff --git a/Decay/WeakCurrents/FivePionCurrent.cc b/Decay/WeakCurrents/FivePionCurrent.cc
	--- a/Decay/WeakCurrents/FivePionCurrent.cc
	+++ b/Decay/WeakCurrents/FivePionCurrent.cc
	@@ -1,770 +1,728 @@
	// -- C++ --
	//
	// FivePionCurrent.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
	// Copyright (C) 2002-2017 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 FivePionCurrent class.
	//

	#include "FivePionCurrent.h"
	#include "ThePEG/Utilities/DescribeClass.h"
	#include "ThePEG/Interface/Parameter.h"
	#include "ThePEG/Interface/Switch.h"
	#include "ThePEG/Interface/ClassDocumentation.h"
	#include "ThePEG/Persistency/PersistentOStream.h"
	#include "ThePEG/Persistency/PersistentIStream.h"
	#include "ThePEG/Helicity/WaveFunction/ScalarWaveFunction.h"

	using namespace Herwig;
	using namespace ThePEG::Helicity;

	FivePionCurrent::FivePionCurrent() {
	// set the number of modes
	addDecayMode(2,-1);
	addDecayMode(2,-1);
	addDecayMode(2,-1);
	setInitialModes(3);
	// masses of the intermediates
	_rhomass = 776*MeV;
	_a1mass = 1260*MeV;
	_omegamass = 782*MeV;
	_sigmamass = 800*MeV;
	// widths of the intermediates
	_rhowidth = 150*MeV;
	_a1width = 400*MeV;
	_omegawidth = 8.5*MeV;
	_sigmawidth = 600*MeV;
	// use local values of the resonance parameters
	_localparameters=true;
	// include the rho Breit-Wigners in omega decay
	_rhoomega = true;
	// Normalisation parameters for the different currents
	_c =4.*GeV2;
	_c0=3.;
	// various meson coupling constants
	_fomegarhopi=0.07/MeV;
	_grhopipi=6.0;
	_garhopi=6.*GeV;
	_faaf=4.*GeV;
	_ffpipi=5.*GeV;
	_presigma = ZERO;
	_preomega = ZERO;
	}

	inline void FivePionCurrent::doinit() {
	- WeakDecayCurrent::doinit();
	+ WeakCurrent::doinit();
	if(!_localparameters) {
	_rhomass = getParticleData(ParticleID::rhominus)->mass();
	_rhowidth = getParticleData(ParticleID::rhominus)->width();
	_omegamass = getParticleData(ParticleID::omega)->mass();
	_omegawidth = getParticleData(ParticleID::omega)->width();
	_sigmamass = getParticleData(9000221)->mass();
	_sigmawidth = getParticleData(9000221)->width();
	_a1mass = getParticleData(ParticleID::a_1minus)->mass();
	_a1width = getParticleData(ParticleID::a_1minus)->width();
	}
	// prefactors
	_presigma = _c/sqr(sqr(_a1mass)_sigmamass_rhomass)_faaf_ffpipi*
	_garhopi*_grhopipi;
	_preomega = _c0_fomegarhopisqr(_grhopipi/(sqr(_rhomass)*_omegamass));
	}

	void FivePionCurrent::persistentOutput(PersistentOStream & os) const {
	static const InvEnergy7 InvGeV7 = pow<-7,1>(GeV);
	static const InvEnergy3 InvGeV3 = pow<-3,1>(GeV);
	os << ounit(_rhomass,GeV) << ounit(_a1mass,GeV) << ounit(_omegamass,GeV)
	<< ounit(_sigmamass,GeV) << ounit(_rhowidth,GeV)
	<< ounit(_a1width,GeV) << ounit(_omegawidth,GeV) << ounit(_sigmawidth,GeV)
	<< _localparameters << ounit(_c,GeV2) << _c0
	<< ounit(_fomegarhopi,1/GeV) << _grhopipi << ounit(_garhopi,GeV)
	<< ounit(_faaf,GeV) << ounit(_ffpipi,GeV)
	<< ounit(_preomega,InvGeV7) << ounit(_presigma,InvGeV3) << _rhoomega;
	}

	void FivePionCurrent::persistentInput(PersistentIStream & is, int) {
	static const InvEnergy7 InvGeV7 = pow<-7,1>(GeV);
	static const InvEnergy3 InvGeV3 = pow<-3,1>(GeV);
	is >> iunit(_rhomass,GeV) >> iunit(_a1mass,GeV) >> iunit(_omegamass,GeV)
	>> iunit(_sigmamass,GeV) >> iunit(_rhowidth,GeV)
	>> iunit(_a1width,GeV) >> iunit(_omegawidth,GeV) >> iunit(_sigmawidth,GeV)
	>> _localparameters >> iunit(_c,GeV2) >> _c0
	>> iunit(_fomegarhopi,1/GeV) >> _grhopipi >> iunit(_garhopi,GeV)
	>> iunit(_faaf,GeV) >> iunit(_ffpipi,GeV)
	>> iunit(_preomega,InvGeV7) >> iunit(_presigma,InvGeV3) >> _rhoomega;
	}

	// The following static variable is needed for the type
	// description system in ThePEG.
	-DescribeClass<FivePionCurrent,WeakDecayCurrent>
	+DescribeClass<FivePionCurrent,WeakCurrent>
	describeHerwigFivePionCurrent("Herwig::FivePionCurrent", "HwWeakCurrents.so");

	void FivePionCurrent::Init() {

	static ClassDocumentation<FivePionCurrent> documentation
	("The FivePionCurrent class implements the model of hep-ph/0602162",
	"The model of \\cite{Kuhn:2006nw} was used for the hadronic five pion current.",
	"\\bibitem{Kuhn:2006nw} J.~H.~Kuhn and Z.~Was, hep-ph/0602162, (2006).");

	static Parameter<FivePionCurrent,Energy> interfaceRhoMass
	("RhoMass",
	"The mass of the rho meson",
	&FivePionCurrent::_rhomass, MeV, 776MeV, 500MeV, 1000*MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceA1Mass
	("A1Mass",
	"The mass of the a_1 meson",
	&FivePionCurrent::_a1mass, MeV, 1260MeV, 1000MeV, 1500*MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceOmegaMass
	("OmegaMass",
	"The mass of the omega meson",
	&FivePionCurrent::_omegamass, MeV, 782MeV, 600MeV, 900*MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceSigmaMass
	("SigmaMass",
	"The mass of the sigma meson",
	&FivePionCurrent::_sigmamass, MeV, 800MeV, 400MeV, 1200*MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceRhoWidth
	("RhoWidth",
	"The width of the rho meson",
	&FivePionCurrent::_rhowidth, MeV, 150MeV, 100MeV, 300*MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceA1Width
	("A1Width",
	"The width of the a_1 meson",
	&FivePionCurrent::_a1width, MeV, 400MeV, 100MeV, 800*MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceOmegaWidth
	("OmegaWidth",
	"The width of the omega meson",
	&FivePionCurrent::_omegawidth, MeV, 8.5MeV, 1.0MeV, 20.0*MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceSigmaWidth
	("SigmaWidth",
	"The width of the sigma meson",
	&FivePionCurrent::_sigmawidth, MeV, 600MeV, 100MeV, 1200*MeV,
	false, false, Interface::limited);

	static Switch<FivePionCurrent,bool> interfaceLocalParameters
	("LocalParameters",
	"Use local values of the meson masses and widths or those from the"
	" ParticleData objects.",
	&FivePionCurrent::_localparameters, true, false, false);
	static SwitchOption interfaceLocalParametersLocal
	(interfaceLocalParameters,
	"Local",
	"Use local values",
	true);
	static SwitchOption interfaceLocalParametersParticleData
	(interfaceLocalParameters,
	"ParticleData",
	"Use values from the particle data objects",
	false);

	static Switch<FivePionCurrent,bool> interfaceRhoOmega
	("RhoOmega",
	"Option for the treatment of the rho Breit-Wigners in the omega decay",
	&FivePionCurrent::_rhoomega, true, false, false);
	static SwitchOption interfaceRhoOmegaInclude
	(interfaceRhoOmega,
	"Yes",
	"Include the rho Breit-Wigners",
	true);
	static SwitchOption interfaceRhoOmegaOmit
	(interfaceRhoOmega,
	"No",
	"Don't include the rho Breit-Wigners",
	false);

	static Parameter<FivePionCurrent,Energy2> interfaceC
	("C",
	"The normalisation parameter for the a_1 sigma current",
	&FivePionCurrent::_c, GeV2, 4.0GeV2, 0.1GeV2, 20.0*GeV2,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,double> interfaceC0
	("C0",
	"The normalisation constant for the omega-rho current",
	&FivePionCurrent::_c0, 3., 0.1, 10.0,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,InvEnergy> interfacegomegarhopi
	("fomegarhopi",
	"The coupling of omega-rho-pi",
	&FivePionCurrent::_fomegarhopi, 1./MeV, 0.07/MeV, 0.01/MeV, 0.2/MeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,double> interfacegrhopipi
	("grhopipi",
	"The coupling for rho-pi-pi",
	&FivePionCurrent::_grhopipi, 6.0, 1.0, 20.0,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfacegarhopi
	("garhopi",
	"The coupling of a-rho-pi",
	&FivePionCurrent::_garhopi, GeV, 6.0GeV, 0.1GeV, 20.0*GeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfacefaaf
	("faaf",
	"The coupling of a-a-f",
	&FivePionCurrent::_faaf, GeV, 4.0GeV, 0.1GeV, 20.0*GeV,
	false, false, Interface::limited);

	static Parameter<FivePionCurrent,Energy> interfaceffpipi
	("ffpipi",
	"The coupling of f-pi-pi",
	&FivePionCurrent::_ffpipi, GeV, 5.0GeV, 0.1GeV, 20.0*GeV,
	false, false, Interface::limited);
	}

	bool FivePionCurrent::accept(vector<int> id) {
	bool allowed(false);
	// check five products
	if(id.size()!=5){return false;}
	int npiminus=0,npiplus=0,npi0=0;
	for(unsigned int ix=0;ix<id.size();++ix) {
	if(id[ix]==ParticleID:: piplus){++npiplus;}
	else if(id[ix]==ParticleID::piminus){++npiminus;}
	else if(id[ix]==ParticleID::pi0){++npi0;}
	}
	if(npiplus>npiminus) swap(npiplus,npiminus);
	if( npiminus==3&&npiplus==2&&npi0==0) allowed=true;
	else if(npiminus==2&&npiplus==1&&npi0==2) allowed=true;
	else if(npiminus==1&&npiplus==0&&npi0==4) allowed=true;
	return allowed;
	}

	-bool FivePionCurrent::createMode(int icharge, unsigned int imode,
	- DecayPhaseSpaceModePtr mode,
	+bool FivePionCurrent::createMode(int icharge, tcPDPtr resonance,
	+ IsoSpin::IsoSpin Itotal, IsoSpin::I3 i3,
	+ unsigned int imode,PhaseSpaceModePtr mode,
	unsigned int iloc,unsigned int ires,
	- DecayPhaseSpaceChannelPtr phase,Energy upp) {
	+ PhaseSpaceChannel phase, Energy upp ) {
	// check the charge
	if(abs(icharge)!=3) return false;
	+ // check the total isospin
	+ if(Itotal!=IsoSpin::IUnknown) {
	+ if(Itotal!=IsoSpin::IOne) return false;
	+ }
	+ // check I_3
	+ if(i3!=IsoSpin::I3Unknown) {
	+ switch(i3) {
	+ case IsoSpin::I3Zero:
	+ return false;
	+ break;
	+ case IsoSpin::I3One:
	+ if(icharge ==-3) return false;
	+ break;
	+ case IsoSpin::I3MinusOne:
	+ if(icharge ==3) return false;
	+ default:
	+ return false;
	+ }
	+ }
	// check that the modes are kinematical allowed
	Energy min(ZERO);
	// 3 pi- 2pi+
	if(imode==0) {
	min=5.*getParticleData(ParticleID::piplus)->mass();
	}
	// 2pi- pi+ 2pi0
	else if(imode==1) {
	min=3.*getParticleData(ParticleID::piplus)->mass()
	+2.*getParticleData(ParticleID::pi0)->mass();
	}
	// pi- 4pi0
	else {
	min= getParticleData(ParticleID::piplus)->mass()
	+4.*getParticleData(ParticleID::pi0)->mass();
	}
	if(min>upp) return false;
	// intermediates for the channels
	tPDPtr omega(getParticleData(ParticleID::omega)),rhop,rhom,
	rho0(getParticleData(ParticleID::rho0)),a1m,a10(getParticleData(ParticleID::a_10)),
	sigma(getParticleData(9000221));
	if(icharge==3) {
	rhop = getParticleData(ParticleID::rhominus);
	rhom = getParticleData(ParticleID::rhoplus);
	a1m = getParticleData(ParticleID::a_1plus);
	}
	else {
	rhop = getParticleData(ParticleID::rhoplus);
	rhom = getParticleData(ParticleID::rhominus);
	a1m = getParticleData(ParticleID::a_1minus);
	}
	- DecayPhaseSpaceChannelPtr newchannel;
	+ if(resonance && resonance !=a1m) return false;
	// all charged mode
	if(imode==0) {
	if(sigma) {
	// 1st two a_1 sigma channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+1, iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+4,ires+2,iloc+5,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+1,
	+ ires+4,iloc+2,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+4,ires+2,iloc+5,
	+ ires+1,a1m ,ires+3,rho0 ,ires+3,iloc+2,
	+ ires+4,iloc+1,ires+4,iloc+3));
	// 2nd two a_1 sigma channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc );
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+1, iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc );
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+4, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+4,ires+2,iloc+1,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+5,
	+ ires+4,iloc+2,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+4,ires+2,iloc+1,
	+ ires+1,a1m ,ires+3,rho0 ,ires+3,iloc+2,
	+ ires+4,iloc+5,ires+4,iloc+3));
	// 3rd two a_1 sigma channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+4, iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+4,ires+2,iloc+2,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+1,
	+ ires+4,iloc+5,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+4,ires+2,iloc+2,
	+ ires+1,a1m ,ires+3,rho0 ,ires+3,iloc+5,
	+ ires+4,iloc+1,ires+4,iloc+3));
	// 4th two a_1 sigma channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+1, iloc+3);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+3);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+3,ires+2,iloc+5,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+1,
	+ ires+4,iloc+2,ires+4,iloc+4));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+3,ires+2,iloc+5,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+2,
	+ ires+4,iloc+1,ires+4,iloc+4));
	// 5th two a_1 sigma channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc );
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+1, iloc+3);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc );
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+4, iloc+3);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+3,ires+2,iloc+1,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+5,
	+ ires+4,iloc+2,ires+4,iloc+4));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+3,ires+2,iloc+1,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+2,
	+ ires+4,iloc+5,ires+4,iloc+4));
	// 6th two a_1 sigma channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc+4, iloc+3);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+3);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+3,ires+2,iloc+2,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+1,
	+ ires+4,iloc+5,ires+4,iloc+4));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma ,ires+2,iloc+3,ires+2,iloc+2,
	+ ires+1,a1m ,ires+3, rho0 ,ires+3,iloc+5,
	+ ires+4,iloc+1,ires+4,iloc+4));
	}
	}
	// 2 neutral mode
	else if(imode==1) {
	// first three omega channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+4);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+2);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+1);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+4);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rhop ,0,0.0, iloc , iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+4);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+4,ires+2,iloc+5,
	+ ires+1,omega,ires+3, rho0,ires+3,iloc+3,
	+ ires+4,iloc+1,ires+4,iloc+2));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+4,ires+2,iloc+5,
	+ ires+1,omega,ires+3, rhop,ires+3,iloc+2,
	+ ires+4,iloc+1,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+4,ires+2,iloc+5,
	+ ires+1,omega,ires+3, rhom,ires+3,iloc+1,
	+ ires+4,iloc+2,ires+4,iloc+3));
	// second three omega channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+4);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+2);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+3);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+4);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+3);
	- newchannel->addIntermediate(rhop ,0,0.0, iloc , iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+4);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+2,ires+2,iloc+5,
	+ ires+1,omega,ires+3, rho0,ires+3,iloc+3,
	+ ires+4,iloc+1,ires+4,iloc+4));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+2,ires+2,iloc+5,
	+ ires+1,omega,ires+3, rhop,ires+3,iloc+4,
	+ ires+4,iloc+1,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+2,ires+2,iloc+5,
	+ ires+1,omega,ires+3, rhom,ires+3,iloc+1,
	+ ires+4,iloc+4,ires+4,iloc+3));
	// third three omega channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+1);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rhop ,0,0.0, iloc , iloc+4);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+4);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+4,ires+2,iloc+3,
	+ ires+1,omega,ires+3, rho0,ires+3,iloc+5,
	+ ires+4,iloc+1,ires+4,iloc+2));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+4,ires+2,iloc+3,
	+ ires+1,omega,ires+3, rhop,ires+3,iloc+2,
	+ ires+4,iloc+1,ires+4,iloc+5));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+4,ires+2,iloc+3,
	+ ires+1,omega,ires+3, rhom,ires+3,iloc+1,
	+ ires+4,iloc+2,ires+4,iloc+5));
	// fourth three omega channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+3);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc+3);
	- newchannel->addIntermediate(rhop ,0,0.0, iloc , iloc+4);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- newchannel->addIntermediate(omega ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+4);
	- mode->addChannel(newchannel);
	- // first two sigma a_1 channels
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+2,ires+2,iloc+3,
	+ ires+1,omega,ires+3, rho0,ires+3,iloc+5,
	+ ires+4,iloc+1,ires+4,iloc+4));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+2,ires+2,iloc+3,
	+ ires+1,omega,ires+3, rhop,ires+3,iloc+4,
	+ ires+4,iloc+1,ires+4,iloc+5));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, rhom,ires+2,iloc+2,ires+2,iloc+3,
	+ ires+1,omega,ires+3, rhom,ires+3,iloc+1,
	+ ires+4,iloc+4,ires+4,iloc+5));
	if(sigma) {
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+4);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- mode->addChannel(newchannel);
	+ // first two sigma a_1 channels
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma,ires+2,iloc+1,ires+2,iloc+2,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+3,
	+ ires+4,iloc+4,ires+4,iloc+5));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma,ires+2,iloc+1,ires+2,iloc+2,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+5,
	+ ires+4,iloc+4,ires+4,iloc+3));
	// second two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+3);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+2);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+4);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+3);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma,ires+2,iloc+1,ires+2,iloc+4,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+3,
	+ ires+4,iloc+2,ires+4,iloc+5));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma,ires+2,iloc+1,ires+2,iloc+4,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+5,
	+ ires+4,iloc+2,ires+4,iloc+3));
	// third two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+3);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+2, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+3);
	- newchannel->addIntermediate(rho0 ,0,0.0, iloc , iloc+1);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma,ires+2,iloc+3,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rho0,ires+3,iloc+2,
	+ ires+4,iloc+1,ires+4,iloc+4));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1,sigma,ires+2,iloc+3,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rho0,ires+3,iloc+4,
	+ ires+4,iloc+1,ires+4,iloc+2));
	}
	}
	// 4 neutral mode
	else {
	if(sigma) {
	// 1st two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc , iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- mode->addChannel(newchannel);
	- // 2nd two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc );
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+4, iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc );
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+4,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+2,
	+ ires+4,iloc+1,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+4,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+1,
	+ ires+4,iloc+2,ires+4,iloc+3));
	+ // // 2nd two sigma a_1 channels
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+4,ires+2,iloc+1,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+2,
	+ ires+4,iloc+5,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+4,ires+2,iloc+1,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+5,
	+ ires+4,iloc+2,ires+4,iloc+3));
	// 3rd two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+1);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+3);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+1, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+1,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+2,
	+ ires+4,iloc+4,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+1,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+4,
	+ ires+4,iloc+2,ires+4,iloc+3));
	// 4th two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+1, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+3);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc , iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+1, iloc+4);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+2,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+4,
	+ ires+4,iloc+1,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+2,ires+2,iloc+5,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+1,
	+ ires+4,iloc+4,ires+4,iloc+3));
	// 5th two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc , iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc+3, iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc );
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+4, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+4,ires+2,iloc+2,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+5,
	+ ires+4,iloc+1,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+4,ires+2,iloc+2,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+1,
	+ ires+4,iloc+5,ires+4,iloc+3));
	// 6th two sigma a_1 channels
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+3);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+4, iloc+2);
	- mode->addChannel(newchannel);
	- newchannel=new_ptr(DecayPhaseSpaceChannel(*phase));
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-1,-ires-2);
	- newchannel->addIntermediate(sigma ,0,0.0, iloc , iloc+1);
	- newchannel->addIntermediate(a1m ,0,0.0,-ires-3, iloc+4);
	- newchannel->addIntermediate(rhom ,0,0.0, iloc+3, iloc+2);
	- mode->addChannel(newchannel);
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+1,ires+2,iloc+2,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+4,
	+ ires+4,iloc+5,ires+4,iloc+3));
	+ mode->addChannel((PhaseSpaceChannel(phase),ires,a1m,
	+ ires+1, sigma,ires+2,iloc+1,ires+2,iloc+2,
	+ ires+1, a1m,ires+3, rhom,ires+3,iloc+5,
	+ ires+4,iloc+4,ires+4,iloc+3));
	}
	}
	// reset the parameters of the resonances if using local values
	if(_localparameters) {
	mode->resetIntermediate(rhom,_rhomass,_rhowidth);
	mode->resetIntermediate(rhop,_rhomass,_rhowidth);
	mode->resetIntermediate(rho0,_rhomass,_rhowidth);
	mode->resetIntermediate(omega,_omegamass,_omegawidth);
	mode->resetIntermediate(a1m,_a1mass,_a1width);
	mode->resetIntermediate(a10,_a1mass,_a1width);
	if(sigma) mode->resetIntermediate(sigma,_sigmamass,_sigmawidth);
	}
	// return if successful
	return true;
	}

	// the particles produced by the current
	tPDVector FivePionCurrent::particles(int icharge, unsigned int imode,int,int) {
	// particle data objects for the pions
	tPDPtr piplus (getParticleData(ParticleID::piplus ));
	tPDPtr pi0 (getParticleData(ParticleID::pi0 ));
	tPDPtr piminus(getParticleData(ParticleID::piminus));
	if(icharge==3) swap(piplus,piminus);
	tPDVector output(5);
	// all charged
	if(imode==0) {
	output[0]=piminus;
	output[1]=piminus;
	output[2]=piplus;;
	output[3]=piplus;
	output[4]=piminus;
	}
	// two neutral
	else if(imode==1) {
	output[0]=piplus;
	output[1]=piminus;
	output[2]=pi0;;
	output[3]=piminus;
	output[4]=pi0;
	}
	// four neutral
	else {
	output[0]=pi0;
	output[1]=pi0;
	output[2]=piminus;;
	output[3]=pi0;
	output[4]=pi0;
	}
	return output;
	}

	// the decay mode
	unsigned int FivePionCurrent::decayMode(vector<int> idout) {
	unsigned int npi(0);
	for(unsigned int ix=0;ix<idout.size();++ix) {
	if(abs(idout[ix])==ParticleID::pi0) ++npi;
	}
	return npi/2;
	}

	// output the information for the database
	void FivePionCurrent::dataBaseOutput(ofstream & output,bool header,bool create) const {
	if(header) output << "update decayers set parameters=\"";
	if(create) output << "create Herwig::FivePionCurrent " << name()
	<< " HwWeakCurrents.so\n";
	output << "newdef " << name() << ":RhoMass " << _rhomass/MeV << "\n";
	output << "newdef " << name() << ":A1Mass " << _a1mass/MeV << "\n";
	output << "newdef " << name() << ":SigmaMass " << _sigmamass/MeV << "\n";
	output << "newdef " << name() << ":OmegaMass " << _omegamass/MeV << "\n";
	output << "newdef " << name() << ":RhoWidth " << _rhowidth/MeV << "\n";
	output << "newdef " << name() << ":A1Width " << _a1width/MeV << "\n";
	output << "newdef " << name() << ":SigmaWidth " << _sigmawidth/MeV << "\n";
	output << "newdef " << name() << ":OmegaWidth " << _omegawidth/MeV << "\n";
	output << "newdef " << name() << ":LocalParameters " << _localparameters << "\n";
	output << "newdef " << name() << ":RhoOmega " << _rhoomega << "\n";
	output << "newdef " << name() << ":C " << _c/GeV2 << "\n";
	output << "newdef " << name() << ":C0 " << _c0 << "\n";
	output << "newdef " << name() << ":fomegarhopi " <<_fomegarhopi*MeV << "\n";
	output << "newdef " << name() << ":grhopipi " <<_grhopipi << "\n";
	output << "newdef " << name() << ":garhopi " << _garhopi/GeV << "\n";
	output << "newdef " << name() << ":faaf " <<_faaf/GeV << "\n";
	output << "newdef " << name() << ":ffpipi " << _ffpipi/GeV << "\n";
	- WeakDecayCurrent::dataBaseOutput(output,false,false);
	+ WeakCurrent::dataBaseOutput(output,false,false);
	if(header) output << "\n\" where BINARY ThePEGName=\"" << fullName() << "\";\n";
	}

	+void FivePionCurrent::constructSpinInfo(ParticleVector decay) const {
	+ for(unsigned int ix=0;ix<decay.size();++ix)
	+ ScalarWaveFunction::constructSpinInfo(decay[ix],outgoing,true);
	+}
	+
	vector<LorentzPolarizationVectorE>
	-FivePionCurrent::current(const int imode,const int ichan,
	- Energy & scale,const ParticleVector & decay,
	- DecayIntegrator::MEOption meopt) const {
	+FivePionCurrent::current(tcPDPtr,
	+ IsoSpin::IsoSpin Itotal, IsoSpin::I3 i3,
	+ const int imode, const int ichan,Energy & scale,
	+ const tPDVector & outgoing,
	+ const vector<Lorentz5Momentum> & momenta,
	+ DecayIntegrator2::MEOption) const {
	+ // check the isospin
	+ if(Itotal!=IsoSpin::IUnknown && Itotal!=IsoSpin::IOne)
	+ return vector<LorentzPolarizationVectorE>();
	+ int icharge = outgoing[0]->iCharge()+outgoing[1]->iCharge();
	+ // check I_3
	+ if(i3!=IsoSpin::I3Unknown) {
	+ switch(i3) {
	+ case IsoSpin::I3Zero:
	+ return vector<LorentzPolarizationVectorE>();
	+ break;
	+ case IsoSpin::I3One:
	+ if(icharge ==-3) return vector<LorentzPolarizationVectorE>();
	+ break;
	+ case IsoSpin::I3MinusOne:
	+ if(icharge ==3) return vector<LorentzPolarizationVectorE>();
	+ default:
	+ return vector<LorentzPolarizationVectorE>();
	+ }
	+ }
	useMe();
	- if(meopt==DecayIntegrator::Terminate) {
	- for(unsigned int ix=0;ix<5;++ix)
	- ScalarWaveFunction::constructSpinInfo(decay[ix],outgoing,true);
	- return vector<LorentzPolarizationVectorE>(1,LorentzPolarizationVectorE());
	- }
	LorentzVector<complex<InvEnergy2> > output;
	- Lorentz5Momentum q1(decay[0]->momentum());
	- Lorentz5Momentum q2(decay[1]->momentum());
	- Lorentz5Momentum q3(decay[2]->momentum());
	- Lorentz5Momentum q4(decay[3]->momentum());
	- Lorentz5Momentum q5(decay[4]->momentum());
	+ Lorentz5Momentum q1(momenta[0]);
	+ Lorentz5Momentum q2(momenta[1]);
	+ Lorentz5Momentum q3(momenta[2]);
	+ Lorentz5Momentum q4(momenta[3]);
	+ Lorentz5Momentum q5(momenta[4]);
	// total momentum
	Lorentz5Momentum Q(q1+q2+q3+q4+q5);
	Q.rescaleMass();
	scale=Q.mass();
	// decide which decay mode
	if(imode==0) {
	if(ichan<0) {
	output=
	a1SigmaCurrent(0,Q,q1,q2,q3,q4,q5)+
	a1SigmaCurrent(0,Q,q5,q2,q3,q4,q1)+
	a1SigmaCurrent(0,Q,q1,q5,q3,q4,q2)+
	a1SigmaCurrent(0,Q,q1,q2,q4,q3,q5)+
	a1SigmaCurrent(0,Q,q5,q2,q4,q3,q1)+
	a1SigmaCurrent(0,Q,q1,q5,q4,q3,q2);
	}
	else if(ichan==0 ) output=a1SigmaCurrent(2,Q,q1,q2,q3,q4,q5);
	else if(ichan==1 ) output=a1SigmaCurrent(1,Q,q1,q2,q3,q4,q5);
	else if(ichan==2 ) output=a1SigmaCurrent(2,Q,q5,q2,q3,q4,q1);
	else if(ichan==3 ) output=a1SigmaCurrent(1,Q,q5,q2,q3,q4,q1);
	else if(ichan==4 ) output=a1SigmaCurrent(2,Q,q1,q5,q3,q4,q2);
	else if(ichan==5 ) output=a1SigmaCurrent(1,Q,q1,q5,q3,q4,q2);
	else if(ichan==6 ) output=a1SigmaCurrent(2,Q,q1,q2,q4,q3,q5);
	else if(ichan==7 ) output=a1SigmaCurrent(1,Q,q1,q2,q4,q3,q5);
	else if(ichan==8 ) output=a1SigmaCurrent(2,Q,q5,q2,q4,q3,q1);
	else if(ichan==9 ) output=a1SigmaCurrent(1,Q,q5,q2,q4,q3,q1);
	else if(ichan==10) output=a1SigmaCurrent(2,Q,q1,q5,q4,q3,q2);
	else if(ichan==11) output=a1SigmaCurrent(1,Q,q1,q5,q4,q3,q2);
	// identical particle symmetry factor
	output/=sqrt(12.);
	}
	else if(imode==1) {
	if(ichan<0) {
	output=
	rhoOmegaCurrent(0,Q,q1,q2,q3,q4,q5)
	+rhoOmegaCurrent(0,Q,q1,q4,q3,q2,q5)
	+rhoOmegaCurrent(0,Q,q1,q2,q5,q4,q3)
	+rhoOmegaCurrent(0,Q,q1,q4,q5,q2,q3)
	+a1SigmaCurrent(0,Q,q2,q4,q1,q3,q5)
	+a1SigmaCurrent(0,Q,q3,q5,q2,q1,q4)
	+a1SigmaCurrent(0,Q,q3,q5,q4,q1,q2);
	}
	else if(ichan==0 ) output=rhoOmegaCurrent(3,Q,q1,q2,q3,q4,q5);
	else if(ichan==1 ) output=rhoOmegaCurrent(2,Q,q1,q2,q3,q4,q5);
	else if(ichan==2 ) output=rhoOmegaCurrent(1,Q,q1,q2,q3,q4,q5);
	else if(ichan==3 ) output=rhoOmegaCurrent(3,Q,q1,q4,q3,q2,q5);
	else if(ichan==4 ) output=rhoOmegaCurrent(2,Q,q1,q4,q3,q2,q5);
	else if(ichan==5 ) output=rhoOmegaCurrent(1,Q,q1,q4,q3,q2,q5);
	else if(ichan==6 ) output=rhoOmegaCurrent(3,Q,q1,q2,q5,q4,q3);
	else if(ichan==7 ) output=rhoOmegaCurrent(2,Q,q1,q2,q5,q4,q3);
	else if(ichan==8 ) output=rhoOmegaCurrent(1,Q,q1,q2,q5,q4,q3);
	else if(ichan==9 ) output=rhoOmegaCurrent(3,Q,q1,q4,q5,q2,q3);
	else if(ichan==10) output=rhoOmegaCurrent(2,Q,q1,q4,q5,q2,q3);
	else if(ichan==11) output=rhoOmegaCurrent(1,Q,q1,q4,q5,q2,q3);
	else if(ichan==12) output=a1SigmaCurrent(2,Q,q3,q5,q4,q1,q2);
	else if(ichan==13) output=a1SigmaCurrent(1,Q,q3,q5,q4,q1,q2);
	else if(ichan==14) output=a1SigmaCurrent(2,Q,q3,q5,q2,q1,q4);
	else if(ichan==15) output=a1SigmaCurrent(1,Q,q3,q5,q2,q1,q4);
	else if(ichan==16) output=a1SigmaCurrent(2,Q,q2,q4,q1,q3,q5);
	else if(ichan==17) output=a1SigmaCurrent(1,Q,q2,q4,q1,q3,q5);
	// identical particle symmetry factor
	output/=2.;
	}
	else if(imode==2) {
	if(ichan<0) {
	output=
	a1SigmaCurrent(0,Q,q1,q2,q3,q4,q5)+
	a1SigmaCurrent(0,Q,q5,q2,q3,q4,q1)+
	a1SigmaCurrent(0,Q,q2,q4,q3,q1,q5)+
	a1SigmaCurrent(0,Q,q1,q4,q3,q2,q5)+
	a1SigmaCurrent(0,Q,q1,q5,q3,q4,q2)+
	a1SigmaCurrent(0,Q,q4,q5,q3,q1,q2);
	}
	else if(ichan==0 ) output=a1SigmaCurrent(1,Q,q1,q2,q3,q4,q5);
	else if(ichan==1 ) output=a1SigmaCurrent(2,Q,q1,q2,q3,q4,q5);
	else if(ichan==2 ) output=a1SigmaCurrent(1,Q,q5,q2,q3,q4,q1);
	else if(ichan==3 ) output=a1SigmaCurrent(2,Q,q5,q2,q3,q4,q1);
	else if(ichan==4 ) output=a1SigmaCurrent(2,Q,q2,q4,q3,q1,q5);
	else if(ichan==5 ) output=a1SigmaCurrent(1,Q,q2,q4,q3,q1,q5);
	else if(ichan==6 ) output=a1SigmaCurrent(1,Q,q1,q4,q3,q2,q5);
	else if(ichan==7 ) output=a1SigmaCurrent(2,Q,q1,q4,q3,q2,q5);
	else if(ichan==8 ) output=a1SigmaCurrent(1,Q,q1,q5,q3,q4,q2);
	else if(ichan==9 ) output=a1SigmaCurrent(2,Q,q1,q5,q3,q4,q2);
	else if(ichan==10) output=a1SigmaCurrent(2,Q,q4,q5,q3,q1,q2);
	else if(ichan==11) output=a1SigmaCurrent(1,Q,q4,q5,q3,q1,q2);
	// identical particle symmetry factor
	output/=sqrt(24.);
	}
	else {
	throw DecayIntegratorError() << "Unknown decay mode in the "
	<< "FivePionCurrent::"
	<< "hadronCurrent()" << Exception::abortnow;
	}
	// normalise and return the current
	return vector<LorentzPolarizationVectorE>(1, output * pow<3,1>(scale));
	}
	diff --git a/Decay/WeakCurrents/FivePionCurrent.h b/Decay/WeakCurrents/FivePionCurrent.h
	--- a/Decay/WeakCurrents/FivePionCurrent.h
	+++ b/Decay/WeakCurrents/FivePionCurrent.h
	@@ -1,408 +1,428 @@
	// -- C++ --
	//
	// FivePionCurrent.h is a part of Herwig - A multi-purpose Monte Carlo event generator
	// Copyright (C) 2002-2017 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_FivePionCurrent_H
	#define HERWIG_FivePionCurrent_H
	//
	// This is the declaration of the FivePionCurrent class.
	//

	-#include "WeakDecayCurrent.h"
	+#include "WeakCurrent.h"
	#include "ThePEG/Helicity/epsilon.h"

	namespace Herwig {

	using namespace ThePEG;

	/**
	* Here is the documentation of the FivePionCurrent class.
	*
	* @see \ref FivePionCurrentInterfaces "The interfaces"
	* defined for FivePionCurrent.
	*/
	-class FivePionCurrent: public WeakDecayCurrent {
	+class FivePionCurrent: public WeakCurrent {

	public:

	/**
	* The default constructor.
	*/
	FivePionCurrent();

	/** @name Methods for the construction of the phase space integrator. */
	//@{
	/**
	* Complete the construction of the decay mode for integration.classes inheriting
	* from this one.
	- * @param icharge The total charge of the outgoing particles in the current.
	- * @param imode The mode in the current being asked for.
	- * @param mode The phase space mode for the integration
	- * @param iloc The location of the of the first particle from the current in
	- * the list of outgoing particles.
	- * @param ires The location of the first intermediate for the current.
	- * @param phase The prototype phase space channel for the integration.
	- * @param upp The maximum possible mass the particles in the current are
	- * allowed to have.
	+ * This method is purely virtual and must be implemented in the classes inheriting
	+ * from WeakCurrent.
	+ * @param icharge The total charge of the outgoing particles in the current.
	+ * @param resonance If specified only include terms with this particle
	+ * @param Itotal If specified the total isospin of the current
	+ * @param I3 If specified the thrid component of isospin
	+ * @param imode The mode in the current being asked for.
	+ * @param mode The phase space mode for the integration
	+ * @param iloc The location of the of the first particle from the current in
	+ * the list of outgoing particles.
	+ * @param ires The location of the first intermediate for the current.
	+ * @param phase The prototype phase space channel for the integration.
	+ * @param upp The maximum possible mass the particles in the current are
	+ * allowed to have.
	* @return Whether the current was sucessfully constructed.
	*/
	- virtual bool createMode(int icharge,unsigned int imode,DecayPhaseSpaceModePtr mode,
	+ virtual bool createMode(int icharge, tcPDPtr resonance,
	+ IsoSpin::IsoSpin Itotal, IsoSpin::I3 i3,
	+ unsigned int imode,PhaseSpaceModePtr mode,
	unsigned int iloc,unsigned int ires,
	- DecayPhaseSpaceChannelPtr phase,Energy upp);
	+ PhaseSpaceChannel phase, Energy upp );

	/**
	* The particles produced by the current.
	* @param icharge The total charge of the particles in the current.
	* @param imode The mode for which the particles are being requested
	* @param iq The PDG code for the quark
	* @param ia The PDG code for the antiquark
	* @return The external particles for the current.
	*/
	virtual tPDVector particles(int icharge, unsigned int imode, int iq, int ia);
	//@}

	/**
	* Hadronic current. This method is purely virtual and must be implemented in
	* all classes inheriting from this one.
	+ * @param resonance If specified only include terms with this particle
	+ * @param Itotal If specified the total isospin of the current
	+ * @param I3 If specified the thrid component of isospin
	* @param imode The mode
	* @param ichan The phase-space channel the current is needed for.
	* @param scale The invariant mass of the particles in the current.
	- * @param decay The decay products
	+ * @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 current.
	*/
	- virtual vector<LorentzPolarizationVectorE>
	- current(const int imode,const int ichan,Energy & scale,
	- const ParticleVector & decay, DecayIntegrator::MEOption meopt) const;
	+ virtual vector<LorentzPolarizationVectorE>
	+ current(tcPDPtr resonance,
	+ IsoSpin::IsoSpin Itotal, IsoSpin::I3 i3,
	+ const int imode, const int ichan,Energy & scale,
	+ const tPDVector & outgoing,
	+ const vector<Lorentz5Momentum> & momenta,
	+ DecayIntegrator2::MEOption meopt) const;
	+
	+ /**
	+ * Construct the SpinInfo for the decay products
	+ */
	+ void constructSpinInfo(ParticleVector decay) const;

	/**
	* Accept the decay.
	* @param id The id's of the particles in the current.
	* @return Can this current have the external particles specified.
	*/
	virtual bool accept(vector<int> id);

	/**
	* Return the decay mode number for a given set of particles in the current.
	* @param id The id's of the particles in the current.
	* @return The number of the mode
	*/
	virtual unsigned int decayMode(vector<int> id);

	/**
	* 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
	* @param create Whether or not to add a statement creating the object
	*/
	virtual void dataBaseOutput(ofstream & os,bool header,bool create) const;

	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:

	/**
	* Methods to calculate the Breit-Wigner distributions for the various
	* mesons.
	*/
	//@{
	/**
	* Breit-Wigner for the \f$\rho\f$.
	* @param scale The virtual mass
	*/
	Complex rhoBreitWigner(Energy2 scale) const {
	Energy2 m2=sqr(_rhomass);
	return m2/(m2-scale-Complex(0.,1.)_rhomass_rhowidth);
	}

	/**
	* Breit-Wigner for the \f$a_1\f$.
	* @param scale The virtual mass
	*/
	Complex a1BreitWigner(Energy2 scale) const {
	Energy2 m2=sqr(_a1mass);
	return m2/(m2-scale-Complex(0.,1.)_a1mass_a1width);
	}

	/**
	* Breit-Wigner for the \f$\omega\f$.
	* @param scale The virtual mass
	*/
	Complex omegaBreitWigner(Energy2 scale) const {
	Energy2 m2=sqr(_omegamass);
	return m2/(m2-scale-Complex(0.,1.)_omegamass_omegawidth);
	}

	/**
	* Breit-Wigner for the \f$\sigma\f$.
	* @param scale The virtual mass
	*/
	Complex sigmaBreitWigner(Energy2 scale) const {
	Energy2 m2=sqr(_sigmamass);
	return m2/(m2-scale-Complex(0.,1.)_sigmamass_sigmawidth);
	}
	//@}

	/**
	* Currents for the different channels
	*/
	//@{
	/**
	* The \f$\rho\omega\f$ current
	* @param iopt Option for the inclusion of \f$\rho\f$ Breit-Wigner terms in the
	* \f$\omega\f$ decay piece
	* @param Q The total momentum for the current
	* @param q1 The first momentum
	* @param q2 The first momentum
	* @param q3 The first momentum
	* @param q4 The first momentum
	* @param q5 The first momentum
	*/
	LorentzVector<complex<InvEnergy2> >
	rhoOmegaCurrent(unsigned int iopt,
	const Lorentz5Momentum & Q,
	const Lorentz5Momentum & q1,
	const Lorentz5Momentum & q2,
	const Lorentz5Momentum & q3,
	const Lorentz5Momentum & q4,
	const Lorentz5Momentum & q5) const {
	// prefactor
	complex<InvEnergy7> pre(_preomegaa1BreitWigner(Q.m2())
	omegaBreitWigner((q1+q2+q3).m2())*
	rhoBreitWigner((q4+q5).m2()));
	// omega piece
	Complex omega(-1.);
	if(_rhoomega) {
	if(iopt==1) omega=rhoBreitWigner((q2+q3).m2());
	else if(iopt==2) omega=rhoBreitWigner((q1+q3).m2());
	else if(iopt==3) omega=rhoBreitWigner((q1+q2).m2());
	else
	omega=rhoBreitWigner((q2+q3).m2())+rhoBreitWigner((q1+q3).m2())+
	rhoBreitWigner((q1+q2).m2());
	}
	LorentzVector<complex<Energy3> > omegacurrent(Helicity::epsilon(q1,q2,q3));
	LorentzVector<complex<InvEnergy2> > output =
	pre * omega * Helicity::epsilon(q4-q5,omegacurrent,Q);
	return output;
	}

	/**
	* The \f$a_1\sigma\f$ current
	* @param iopt Option for the inclusion of \f$\rho\f$ Breit-Wigner terms in the
	* \f$a_1\f$ decay piece
	* @param Q The total momentum for the current
	* @param q1 The first momentum
	* @param q2 The first momentum
	* @param q3 The first momentum
	* @param q4 The first momentum
	* @param q5 The first momentum
	*/
	LorentzVector<complex<InvEnergy2> >
	a1SigmaCurrent(unsigned int iopt,
	const Lorentz5Momentum & Q,
	const Lorentz5Momentum & q1,
	const Lorentz5Momentum & q2,
	const Lorentz5Momentum & q3,
	const Lorentz5Momentum & q4,
	const Lorentz5Momentum & q5) const {
	Lorentz5Momentum pa1(q1+q2+q3);pa1.rescaleMass();
	Energy2 ma12(pa1.m2());
	complex<InvEnergy3> pre(_presigmaa1BreitWigner(Q.m2())a1BreitWigner(ma12)*
	sigmaBreitWigner((q4+q5).m2()));
	Energy2 pdot[2]={q2(q1-q3),q1(q2-q3)};
	LorentzPolarizationVectorE rho[2] =
	{(pdot[0]/ma12pa1-q1+q3)rhoBreitWigner((q1+q3).m2()),
	(pdot[1]/ma12pa1-q2+q3)rhoBreitWigner((q2+q3).m2())};
	LorentzPolarizationVectorE total;
	if(iopt==1) total = rho[0];
	else if(iopt==2) total = rho[1];
	else total = rho[0]+rho[1];
	Complex qdot = total * Q / Q.m2();
	LorentzPolarizationVectorE cq(Q);
	cq = cq * qdot;
	cq -= total;
	return pre * cq;
	}
	//@}

	protected:

	/** @name Clone Methods. */
	//@{
	/**
	* Make a simple clone of this object.
	* @return a pointer to the new object.
	*/
	virtual IBPtr clone() const {return new_ptr(*this);}

	/** Make a clone of this object, possibly modifying the cloned object
	* to make it sane.
	* @return a pointer to the new object.
	*/
	virtual IBPtr fullclone() const {return new_ptr(*this);}
	//@}

	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 assignment operator is private and must never be called.
	* In fact, it should not even be implemented.
	*/
	FivePionCurrent & operator=(const FivePionCurrent &);

	private:

	/**
	* The masses and widths of the intermediate particles
	*/
	//@{
	/**
	* The mass of the \f$\rho\f$ for the current.
	*/
	Energy _rhomass;

	/**
	* The mass of the \f$a_1\f$ for the current.
	*/
	Energy _a1mass;

	/**
	* The mass of the \f$\omega\f$ for the current.
	*/
	Energy _omegamass;

	/**
	* The mass of the \f$\sigma\f$ for the current.
	*/
	Energy _sigmamass;

	/**
	* The width for the \f$\rho\f$.
	*/
	Energy _rhowidth;

	/**
	* The \f$a_1\f$ width
	*/
	Energy _a1width;

	/**
	* The \f$\omega\f$ width.
	*/
	Energy _omegawidth;

	/**
	* The \f$\sigma\f$ width.
	*/
	Energy _sigmawidth;
	//@}

	/**
	* use local values of the particle masses
	*/
	bool _localparameters;

	/**
	* Option for the treatment of \f$\rho\f$ Breit-Wigners in \f$\omega\f$ decay
	*/
	bool _rhoomega;

	/**
	* Normalisation parameters for the different currents
	*/
	//@{
	/**
	* The \f$c\f$ parameter
	*/
	Energy2 _c;

	/**
	* The \f$c_0\f$ parameter
	*/
	double _c0;

	/**
	* The \f$f_{\omega\rho\pi}\f$ parameter
	*/
	InvEnergy _fomegarhopi;

	/**
	* The \f$g_{\rho\pi\pi}\f$ parameter
	*/
	double _grhopipi;

	/**
	* The \f$G_{a\rho\pi}\f$ parameter
	*/
	Energy _garhopi;

	/**
	* The \f$f_{aaf}\f$ parameter
	*/
	Energy _faaf;

	/**
	* The \f$f_{f\pi\pi}\f$ parameter
	*/
	Energy _ffpipi;
	//@}

	/**
	* Values cached to avoid unnessacary calculations
	*/
	//@{
	/**
	* Prefactor for the \f$\rho\omega\f$ current
	*/
	InvEnergy7 _preomega;

	/**
	* Prefactor for the \f$a_1\sigma\f$ current
	*/
	InvEnergy3 _presigma;
	//@}
	};

	}

	#endif /* HERWIG_FivePionCurrent_H */

File Metadata

Mime Type: text/x-diff
Expires: Sat, Dec 21, 12:45 PM (1 d, 21 h)
Storage Engine: blob
Storage Format: Raw Data
Storage Handle: 4022814
Default Alt Text: (58 KB)

No OneTemporaryActions

View Options

File Metadata

Event Timeline

No OneTemporary
Actions