Page MenuHomeHEPForge
Files F8308511

No OneTemporary
Actions

View Options

diff --git a/Decay/General/GeneralTwoBodyDecayer.cc b/Decay/General/GeneralTwoBodyDecayer.cc
--- a/Decay/General/GeneralTwoBodyDecayer.cc
+++ b/Decay/General/GeneralTwoBodyDecayer.cc
@@ -1,835 +1,837 @@
// -*- C++ -*-
//
// GeneralTwoBodyDecayer.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2019 The Herwig Collaboration
//
// Herwig is licenced under version 3 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
//
// This is the implementation of the non-inlined, non-templated member
// functions of the GeneralTwoBodyDecayer class.
//
#include "GeneralTwoBodyDecayer.h"
#include "ThePEG/Utilities/DescribeClass.h"
#include "ThePEG/Interface/ClassDocumentation.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
#include "ThePEG/PDT/DecayMode.h"
#include "ThePEG/Utilities/Exception.h"
#include "Herwig/Shower/RealEmissionProcess.h"
#include "Herwig/Utilities/Kinematics.h"
-#include "Herwig/Shower/QTilde/Dark/HiddenValleyModel.h"
+#include "Herwig/Models/HiddenValley/HiddenValleyModel.h"
using namespace Herwig;
ParticleVector GeneralTwoBodyDecayer::decay(const Particle & parent,
const tPDVector & children) const {
// return empty vector if products heavier than parent
Energy mout(ZERO);
for(tPDVector::const_iterator it=children.begin();
it!=children.end();++it) mout+=(**it).massMin();
if(mout>parent.mass()) return ParticleVector();
// generate the decay
bool cc;
int imode=modeNumber(cc,parent.dataPtr(),children);
// generate the kinematics
ParticleVector decay=generate(generateIntermediates(),cc,imode,parent);
// make the colour connections
colourConnections(parent, decay);
// return the answer
return decay;
}
void GeneralTwoBodyDecayer::doinit() {
PerturbativeDecayer::doinit();
assert( incoming_ && outgoing_.size()==2);
//create phase space mode
addMode(new_ptr(PhaseSpaceMode(incoming_,{outgoing_[0],outgoing_[1]},
maxWeight_)));
+
+ _HVmodel = dynamic_ptr_cast<tHiddenValleyPtr>(generator()->standardModel());
}
int GeneralTwoBodyDecayer::modeNumber(bool & cc, tcPDPtr parent,
const tPDVector & children) const {
long parentID = parent->id();
long id1 = children[0]->id();
long id2 = children[1]->id();
cc = false;
long out1 = outgoing_[0]->id();
long out2 = outgoing_[1]->id();
if( parentID == incoming_->id() &&
((id1 == out1 && id2 == out2) ||
(id1 == out2 && id2 == out1)) ) {
return 0;
}
else if(incoming_->CC() && parentID == incoming_->CC()->id()) {
cc = true;
if( outgoing_[0]->CC()) out1 = outgoing_[0]->CC()->id();
if( outgoing_[1]->CC()) out2 = outgoing_[1]->CC()->id();
if((id1 == out1 && id2 == out2) ||
(id1 == out2 && id2 == out1)) return 0;
}
return -1;
}
void GeneralTwoBodyDecayer::
colourConnections(const Particle & parent,
const ParticleVector & out) const {
PDT::Colour incColour(parent.data().iColour());
PDT::Colour outaColour(out[0]->data().iColour());
PDT::Colour outbColour(out[1]->data().iColour());
//incoming colour singlet
if(incColour == PDT::Colour0) {
// colour triplet-colourantitriplet
if((outaColour == PDT::Colour3 && outbColour == PDT::Colour3bar) ||
(outaColour == PDT::Colour3bar && outbColour == PDT::Colour3)) {
bool ac(out[0]->id() < 0);
out[0]->colourNeighbour(out[1],!ac);
}
//colour octet
else if(outaColour == PDT::Colour8 && outbColour == PDT::Colour8) {
out[0]->colourNeighbour(out[1]);
out[0]->antiColourNeighbour(out[1]);
}
// colour singlets
else if(outaColour == PDT::Colour0 && outbColour == PDT::Colour0) {
}
// unknown
else
throw Exception() << "Unknown outgoing colours for decaying "
<< "colour singlet in "
<< "GeneralTwoBodyDecayer::colourConnections "
<< outaColour << " " << outbColour
<< Exception::runerror;
}
//incoming colour triplet
else if(incColour == PDT::Colour3) {
// colour triplet + singlet
if(outaColour == PDT::Colour3 && outbColour == PDT::Colour0) {
out[0]->incomingColour(const_ptr_cast<tPPtr>(&parent));
}
//opposite order
else if(outaColour == PDT::Colour0 && outbColour == PDT::Colour3) {
out[1]->incomingColour(const_ptr_cast<tPPtr>(&parent));
}
// octet + triplet
else if(outaColour == PDT::Colour8 && outbColour == PDT::Colour3) {
out[0]->incomingColour(const_ptr_cast<tPPtr>(&parent));
out[1]->antiColourNeighbour(out[0]);
}
//opposite order
else if(outaColour == PDT::Colour3 && outbColour == PDT::Colour8) {
out[1]->incomingColour(const_ptr_cast<tPPtr>(&parent));
out[0]->antiColourNeighbour(out[1]);
}
else if(outaColour == PDT::Colour3bar && outaColour == PDT::Colour3bar) {
tColinePtr col[2] = {ColourLine::create(out[0],true),
ColourLine::create(out[1],true)};
parent.colourLine()->setSinkNeighbours(col[0],col[1]);
}
else
throw Exception() << "Unknown outgoing colours for decaying "
<< "colour triplet in "
<< "GeneralTwoBodyDecayer::colourConnections() "
<< outaColour << " " << outbColour
<< Exception::runerror;
}
// incoming colour anti triplet
else if(incColour == PDT::Colour3bar) {
// colour antitriplet +singlet
if(outaColour == PDT::Colour3bar && outbColour == PDT::Colour0) {
out[0]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
}
//opposite order
else if(outaColour == PDT::Colour0 && outbColour == PDT::Colour3bar) {
out[1]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
}
//octet + antitriplet
else if(outaColour == PDT::Colour3bar && outbColour == PDT::Colour8) {
out[1]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
out[0]->colourNeighbour(out[1]);
}
//opposite order
else if(outaColour == PDT::Colour8 && outbColour == PDT::Colour3bar) {
out[0]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
out[1]->colourNeighbour(out[0]);
}
else if(outaColour == PDT::Colour3 && outbColour == PDT::Colour3) {
tColinePtr col[2] = {ColourLine::create(out[0]),
ColourLine::create(out[1])};
parent.antiColourLine()->setSourceNeighbours(col[0],col[1]);
}
else
throw Exception() << "Unknown outgoing colours for decaying "
<< "colour antitriplet "
<< "in GeneralTwoBodyDecayer::colourConnections() "
<< outaColour << " " << outbColour
<< Exception::runerror;
}
//incoming colour octet
else if(incColour == PDT::Colour8) {
// triplet-antitriplet
if(outaColour == PDT::Colour3&&outbColour == PDT::Colour3bar) {
out[0]->incomingColour(const_ptr_cast<tPPtr>(&parent));
out[1]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
}
// opposite order
else if(outbColour == PDT::Colour3&&outaColour == PDT::Colour3bar) {
out[0]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
out[1]->incomingColour(const_ptr_cast<tPPtr>(&parent));
}
// neutral octet
else if(outaColour == PDT::Colour0&&outbColour == PDT::Colour8) {
out[1]->incomingColour(const_ptr_cast<tPPtr>(&parent));
out[1]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
}
else if(outbColour == PDT::Colour0&&outaColour == PDT::Colour8) {
out[0]->incomingColour(const_ptr_cast<tPPtr>(&parent));
out[0]->incomingAntiColour(const_ptr_cast<tPPtr>(&parent));
}
else
throw Exception() << "Unknown outgoing colours for decaying "
<< "colour octet "
<< "in GeneralTwoBodyDecayer::colourConnections() "
<< outaColour << " " << outbColour
<< Exception::runerror;
}
else if(incColour == PDT::Colour6) {
if(outaColour == PDT::Colour3 && outbColour == PDT::Colour3) {
tPPtr tempParent = const_ptr_cast<tPPtr>(&parent);
Ptr<MultiColour>::pointer parentColour =
dynamic_ptr_cast<Ptr<MultiColour>::pointer>
(tempParent->colourInfo());
tColinePtr line1 = const_ptr_cast<tColinePtr>(parentColour->colourLines()[0]);
line1->addColoured(dynamic_ptr_cast<tPPtr>(out[0]));
tColinePtr line2 = const_ptr_cast<tColinePtr>(parentColour->colourLines()[1]);
line2->addColoured(dynamic_ptr_cast<tPPtr>(out[1]));
}
else
throw Exception() << "Unknown outgoing colours for decaying "
<< "colour sextet "
<< "in GeneralTwoBodyDecayer::colourConnections() "
<< outaColour << " " << outbColour
<< Exception::runerror;
}
else if(incColour == PDT::Colour6bar) {
if(outaColour == PDT::Colour3bar && outbColour == PDT::Colour3bar) {
tPPtr tempParent = const_ptr_cast<tPPtr>(&parent);
Ptr<MultiColour>::pointer parentColour =
dynamic_ptr_cast<Ptr<MultiColour>::pointer>
(tempParent->colourInfo());
tColinePtr line1 = const_ptr_cast<tColinePtr>(parentColour->antiColourLines()[0]);
line1->addAntiColoured(dynamic_ptr_cast<tPPtr>(out[0]));
tColinePtr line2 = const_ptr_cast<tColinePtr>(parentColour->antiColourLines()[1]);
line2->addAntiColoured(dynamic_ptr_cast<tPPtr>(out[1]));
}
else
throw Exception() << "Unknown outgoing colours for decaying "
<< "colour anti-sextet "
<< "in GeneralTwoBodyDecayer::colourConnections() "
<< outaColour << " " << outbColour
<< Exception::runerror;
}
else
throw Exception() << "Unknown incoming colour in "
<< "GeneralTwoBodyDecayer::colourConnections() "
<< incColour
<< Exception::runerror;
}
bool GeneralTwoBodyDecayer::twoBodyMEcode(const DecayMode & dm, int & mecode,
double & coupling) const {
assert(dm.parent()->id() == incoming_->id());
ParticleMSet::const_iterator pit = dm.products().begin();
long id1 = (*pit)->id();
++pit;
long id2 = (*pit)->id();
long id1t(outgoing_[0]->id()), id2t(outgoing_[1]->id());
mecode = -1;
coupling = 1.;
if( id1 == id1t && id2 == id2t ) {
return true;
}
else if( id1 == id2t && id2 == id1t ) {
return false;
}
else
assert(false);
return false;
}
void GeneralTwoBodyDecayer::persistentOutput(PersistentOStream & os) const {
- os << incoming_ << outgoing_ << maxWeight_;
+ os << incoming_ << outgoing_ << maxWeight_ << _HVmodel;
}
void GeneralTwoBodyDecayer::persistentInput(PersistentIStream & is, int) {
- is >> incoming_ >> outgoing_ >> maxWeight_;
+ is >> incoming_ >> outgoing_ >> maxWeight_ >> _HVmodel;
}
// The following static variable is needed for the type
// description system in ThePEG.
DescribeAbstractClass<GeneralTwoBodyDecayer,PerturbativeDecayer>
describeHerwigGeneralTwoBodyDecayer("Herwig::GeneralTwoBodyDecayer", "Herwig.so");
void GeneralTwoBodyDecayer::Init() {
static ClassDocumentation<GeneralTwoBodyDecayer> documentation
("This class is designed to be a base class for all 2 body decays"
"in a general model");
}
double GeneralTwoBodyDecayer::brat(const DecayMode &, const Particle & p,
double oldbrat) const {
ParticleVector children = p.children();
if( children.size() != 2 || !p.data().widthGenerator() )
return oldbrat;
// partial width for this mode
Energy scale = p.mass();
Energy pwidth =
partialWidth( make_pair(p.dataPtr(), scale),
make_pair(children[0]->dataPtr(), children[0]->mass()),
make_pair(children[1]->dataPtr(), children[1]->mass()) );
Energy width = p.data().widthGenerator()->width(p.data(), scale);
return pwidth/width;
}
void GeneralTwoBodyDecayer::doinitrun() {
PerturbativeDecayer::doinitrun();
for(unsigned int ix=0;ix<numberModes();++ix) {
double fact = pow(1.5,int(mode(ix)->incoming().first->iSpin())-1);
mode(ix)->maxWeight(fact*mode(ix)->maxWeight());
}
}
double GeneralTwoBodyDecayer::colourFactor(tcPDPtr in, tcPDPtr out1,
tcPDPtr out2) const {
// identical particle symmetry factor
double output = out1->id()==out2->id() ? 0.5 : 1.;
// colour neutral incoming particle
if(in->iColour()==PDT::Colour0) {
// both colour neutral
if(out1->iColour()==PDT::Colour0 && out2->iColour()==PDT::Colour0)
output *= 1.;
// colour triplet/ antitriplet
else if((out1->iColour()==PDT::Colour3 && out2->iColour()==PDT::Colour3bar) ||
(out1->iColour()==PDT::Colour3bar && out2->iColour()==PDT::Colour3 ) ) {
output *= 3.;
}
// colour octet colour octet
else if(out1->iColour()==PDT::Colour8 && out2->iColour()==PDT::Colour8 ) {
output *= 8.;
}
else if((out1->iColour()==PDT::DarkColourFundamental
&& out2->iColour()==PDT::DarkColourAntiFundamental) ||
(out1->iColour()==PDT::DarkColourAntiFundamental
&& out2->iColour()==PDT::DarkColourFundamental)) {
- output *= 3.;
+ output *= _HVmodel->NC();
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour neutral particle in "
<< "GeneralTwoBodyDecayer::colourFactor() for "
<< in->PDGName() << " -> "
<< out1->PDGName() << " " << out2->PDGName()
<< Exception::runerror;
}
// triplet
else if(in->iColour()==PDT::Colour3) {
// colour triplet + neutral
if((out1->iColour()==PDT::Colour0 && out2->iColour()==PDT::Colour3) ||
(out1->iColour()==PDT::Colour3 && out2->iColour()==PDT::Colour0) ) {
output *= 1.;
}
// colour triplet + octet
else if((out1->iColour()==PDT::Colour8 && out2->iColour()==PDT::Colour3) ||
(out1->iColour()==PDT::Colour3 && out2->iColour()==PDT::Colour8) ) {
output *= 4./3.;
}
// colour anti triplet anti triplet
else if(out1->iColour()==PDT::Colour3bar &&
out2->iColour()==PDT::Colour3bar) {
output *= 2.;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour triplet particle in "
<< "GeneralTwoBodyDecayer::colourFactor() for "
<< in->PDGName() << " -> "
<< out1->PDGName() << " " << out2->PDGName()
<< Exception::runerror;
}
// anti triplet
else if(in->iColour()==PDT::Colour3bar) {
// colour anti triplet + neutral
if((out1->iColour()==PDT::Colour0 && out2->iColour()==PDT::Colour3bar ) ||
(out1->iColour()==PDT::Colour3bar && out2->iColour()==PDT::Colour0 ) ) {
output *= 1.;
}
// colour anti triplet + octet
else if((out1->iColour()==PDT::Colour8 && out2->iColour()==PDT::Colour3bar ) ||
(out1->iColour()==PDT::Colour3bar && out2->iColour()==PDT::Colour8 ) ) {
output *= 4./3.;
}
// colour triplet triplet
else if(out1->iColour()==PDT::Colour3 &&
out2->iColour()==PDT::Colour3) {
output *= 2.;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour anti triplet particle in "
<< "GeneralTwoBodyDecayer::colourFactor() for "
<< in->PDGName() << " -> "
<< out1->PDGName() << " " << out2->PDGName()
<< Exception::runerror;
}
else if(in->iColour()==PDT::Colour8) {
// colour octet + neutral
if((out1->iColour()==PDT::Colour0 && out2->iColour()==PDT::Colour8 ) ||
(out1->iColour()==PDT::Colour8 && out2->iColour()==PDT::Colour0 ) ) {
output *= 1.;
}
// colour triplet/antitriplet
else if((out1->iColour()==PDT::Colour3 && out2->iColour()==PDT::Colour3bar) ||
(out1->iColour()==PDT::Colour3bar && out2->iColour()==PDT::Colour3 ) ) {
output *= 0.5;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour octet particle in "
<< "GeneralTwoBodyDecayer::colourFactor() for "
<< in->PDGName() << " -> "
<< out1->PDGName() << " " << out2->PDGName()
<< Exception::runerror;
}
else if(in->iColour()==PDT::Colour6) {
// colour sextet -> triplet triplet
if( out1->iColour()==PDT::Colour3 && out2->iColour()==PDT::Colour3 ) {
output *= 1.;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour sextet particle in "
<< "GeneralTwoBodyDecayer::colourFactor() for "
<< in->PDGName() << " -> "
<< out1->PDGName() << " " << out2->PDGName()
<< Exception::runerror;
}
else if(in->iColour()==PDT::Colour6bar) {
// colour sextet -> triplet triplet
if( out1->iColour()==PDT::Colour3bar && out2->iColour()==PDT::Colour3bar ) {
output *= 1.;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour anti-sextet particle in "
<< "GeneralTwoBodyDecayer::colourFactor() for "
<< in->PDGName() << " -> "
<< out1->PDGName() << " " << out2->PDGName()
<< Exception::runerror;
}
else
throw Exception() << "Unknown colour "
<< in->iColour() << " for the decaying particle in "
<< "GeneralTwoBodyDecayer::colourFactor() for "
<< in->PDGName() << " -> "
<< out1->PDGName() << " " << out2->PDGName()
<< Exception::runerror;
return output;
}
Energy GeneralTwoBodyDecayer::partialWidth(PMPair inpart, PMPair outa,
PMPair outb) const {
// select the number of the mode
tPDVector children;
children.push_back(const_ptr_cast<PDPtr>(outa.first));
children.push_back(const_ptr_cast<PDPtr>(outb.first));
bool cc;
int nmode=modeNumber(cc,inpart.first,children);
tcPDPtr newchild[2] = {mode(nmode)->outgoing()[0],
mode(nmode)->outgoing()[1]};
// make the particles
Lorentz5Momentum pparent = Lorentz5Momentum(inpart.second);
PPtr parent = inpart.first->produceParticle(pparent);
vector<Lorentz5Momentum> pout(2);
double ctheta,phi;
Kinematics::generateAngles(ctheta,phi);
Kinematics::twoBodyDecay(pparent, outa.second, outb.second,
ctheta, phi,pout[0],pout[1]);
if( ( !cc && outa.first!=newchild[0]) ||
( cc && !(( outa.first->CC() && outa.first->CC() == newchild[0])||
( !outa.first->CC() && outa.first == newchild[0]) )))
swap(pout[0],pout[1]);
tPDVector out = {const_ptr_cast<tPDPtr>(outa.first),
const_ptr_cast<tPDPtr>(outb.first)};
double me = me2(-1,*parent,out,pout,Initialize);
Energy pcm = Kinematics::pstarTwoBodyDecay(inpart.second,
outa.second, outb.second);
return me/(8.*Constants::pi)*pcm;
}
void GeneralTwoBodyDecayer::decayInfo(PDPtr incoming, PDPair outgoing) {
incoming_=incoming;
outgoing_.clear();
outgoing_.push_back(outgoing.first );
outgoing_.push_back(outgoing.second);
}
double GeneralTwoBodyDecayer::matrixElementRatio(const Particle & inpart,
const ParticleVector & decay2,
const ParticleVector & decay3,
MEOption meopt,
ShowerInteraction inter) {
// calculate R/B
tPDVector outgoing = {const_ptr_cast<tPDPtr>(decay2[0]->dataPtr()),
const_ptr_cast<tPDPtr>(decay2[1]->dataPtr())};
const vector<Lorentz5Momentum> mom = {decay2[0]->momentum(),
decay2[1]->momentum()};
double B = me2 (0, inpart, outgoing,mom, meopt);
double R = threeBodyME(0, inpart, decay3, inter, meopt);
return R/B;
}
const vector<DVector> & GeneralTwoBodyDecayer::getColourFactors(const Particle & inpart,
const ParticleVector & decay,
unsigned int & nflow) {
// calculate the colour factors for the three-body decay
vector<int> sing,trip,atrip,oct,sex,asex;
for(unsigned int it=0;it<decay.size();++it) {
if (decay[it]->dataPtr()->iColour() == PDT::Colour0 ) sing. push_back(it);
else if(decay[it]->dataPtr()->iColour() == PDT::Colour3 ) trip. push_back(it);
else if(decay[it]->dataPtr()->iColour() == PDT::Colour3bar ) atrip.push_back(it);
else if(decay[it]->dataPtr()->iColour() == PDT::Colour8 ) oct. push_back(it);
else if(decay[it]->dataPtr()->iColour() == PDT::Colour6 ) sex. push_back(it);
else if(decay[it]->dataPtr()->iColour() == PDT::Colour6bar ) asex. push_back(it);
}
// identical particle symmetry factor
double symFactor=1.;
if (( sing.size()==2 && decay[ sing[0]]->id()==decay[ sing[1]]->id()) ||
( trip.size()==2 && decay[ trip[0]]->id()==decay[ trip[1]]->id()) ||
(atrip.size()==2 && decay[atrip[0]]->id()==decay[atrip[1]]->id()) ||
( oct.size()==2 && decay[ oct[0]]->id()==decay[ oct[1]]->id()) ||
( sex.size()==2 && decay[ sex[0]]->id()==decay[ sex[1]]->id()) ||
( asex.size()==2 && decay[ asex[0]]->id()==decay[ asex[1]]->id()))
symFactor /= 2.;
else if (oct.size()==3 &&
decay[oct[0]]->id()==decay[oct[1]]->id() &&
decay[oct[0]]->id()==decay[oct[2]]->id())
symFactor /= 6.;
colour_ = vector<DVector>(1,DVector(1,symFactor*1.));
// decaying colour singlet
if(inpart.dataPtr()->iColour() == PDT::Colour0) {
if(trip.size()==1 && atrip.size()==1 && oct.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*4.));
}
else if(trip.size()==1 && atrip.size()==1 && sing.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*3.));
}
else if (oct.size()==3){
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*24.));
}
else if(sing.size()==3) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor));
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour scalar particle in "
<< "GeneralTwoBodyDecayer::getColourFactors() for "
<< inpart. dataPtr()->PDGName() << " -> "
<< decay[0]->dataPtr()->PDGName() << " "
<< decay[1]->dataPtr()->PDGName() << " "
<< decay[2]->dataPtr()->PDGName()
<< Exception::runerror;
}
// decaying colour triplet
else if(inpart.dataPtr()->iColour() == PDT::Colour3) {
if(trip.size()==1 && sing.size()==1 && oct.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*4./3.));
}
else if(trip.size()==1 && sing.size()==2) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor));
}
else if(trip.size()==1 && oct.size()==2) {
nflow = 2;
colour_.clear();
colour_.resize(2,DVector(2,0.));
colour_[0][0] = symFactor*16./9.; colour_[0][1] = -symFactor*2./9.;
colour_[1][0] = -symFactor*2./9.; colour_[1][1] = symFactor*16./9.;
}
else if(atrip.size()==2 && sing.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*2.));
}
else if(atrip.size()==2 && oct.size()==1) {
nflow = 2;
colour_.clear();
colour_.resize(2,DVector(2,0.));
colour_[0][0] = 8./3.*symFactor; colour_[0][1] =-4./3.*symFactor;
colour_[1][0] = -4./3.*symFactor; colour_[1][1] = 8./3.*symFactor;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour triplet particle in "
<< "GeneralTwoBodyDecayer::getColourFactors() for "
<< inpart. dataPtr()->PDGName() << " -> "
<< decay[0]->dataPtr()->PDGName() << " "
<< decay[1]->dataPtr()->PDGName() << " "
<< decay[2]->dataPtr()->PDGName()
<< Exception::runerror;
}
// decaying colour anti-triplet
else if(inpart.dataPtr()->iColour() == PDT::Colour3bar) {
if(atrip.size()==1 && sing.size()==1 && oct.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*4./3.));
}
else if(atrip.size()==1 && sing.size()==2) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor));
}
else if(atrip.size()==1 && oct.size()==2){
nflow = 2;
colour_.clear();
colour_ .resize(2,DVector(2,0.));
colour_[0][0] = symFactor*16./9.; colour_[0][1] = -symFactor*2./9.;
colour_[1][0] = -symFactor*2./9.; colour_[1][1] = symFactor*16./9.;
}
else if(trip.size()==2 && sing.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*2.));
}
else if(trip.size()==2 && oct.size()==1) {
nflow = 2;
colour_.clear();
colour_.resize(2,DVector(2,0.));
colour_[0][0] = 8./3.*symFactor; colour_[0][1] =-4./3.*symFactor;
colour_[1][0] = -4./3.*symFactor; colour_[1][1] = 8./3.*symFactor;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for decay colour anti-triplet particle in "
<< "GeneralTwoBodyDecayer::getColourFactors() for "
<< inpart. dataPtr()->PDGName() << " -> "
<< decay[0]->dataPtr()->PDGName() << " "
<< decay[1]->dataPtr()->PDGName() << " "
<< decay[2]->dataPtr()->PDGName()
<< Exception::runerror;
}
// decaying colour octet
else if(inpart.dataPtr()->iColour() == PDT::Colour8) {
if(oct.size()==1 && trip.size()==1 && atrip.size()==1) {
nflow = 2;
colour_.clear();
colour_.resize(2,DVector(2,0.));
colour_[0][0] = symFactor*2./3. ; colour_[0][1] = -symFactor*1./12.;
colour_[1][0] = -symFactor*1./12.; colour_[1][1] = symFactor*2./3. ;
}
else if (sing.size()==1 && trip.size()==1 && atrip.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*0.5));
}
else if (oct.size()==2 && sing.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor*3.));
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for a decaying colour octet particle in "
<< "GeneralTwoBodyDecayer::getColourFactors() for "
<< inpart. dataPtr()->PDGName() << " -> "
<< decay[0]->dataPtr()->PDGName() << " "
<< decay[1]->dataPtr()->PDGName() << " "
<< decay[2]->dataPtr()->PDGName()
<< Exception::runerror;
}
// Sextet
else if(inpart.dataPtr()->iColour() == PDT::Colour6) {
if(trip.size()==2 && sing.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor));
}
else if(trip.size()==2 && oct.size()==1) {
nflow = 2;
colour_.clear();
colour_.resize(2,DVector(2,0.));
colour_[0][0] = 4./3.*symFactor; colour_[0][1] = 1./3.*symFactor;
colour_[1][0] = 1./3.*symFactor; colour_[1][1] = 4./3.*symFactor;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for a decaying colour sextet particle in "
<< "GeneralTwoBodyDecayer::getColourFactors() for "
<< inpart. dataPtr()->PDGName() << " -> "
<< decay[0]->dataPtr()->PDGName() << " "
<< decay[1]->dataPtr()->PDGName() << " "
<< decay[2]->dataPtr()->PDGName()
<< Exception::runerror;
}
// anti Sextet
else if(inpart.dataPtr()->iColour() == PDT::Colour6bar) {
if(atrip.size()==2 && sing.size()==1) {
nflow = 1;
colour_ = vector<DVector>(1,DVector(1,symFactor));
}
else if(atrip.size()==2 && oct.size()==1) {
nflow = 2;
colour_.clear();
colour_.resize(2,DVector(2,0.));
colour_[0][0] = 4./3.*symFactor; colour_[0][1] = 1./3.*symFactor;
colour_[1][0] = 1./3.*symFactor; colour_[1][1] = 4./3.*symFactor;
}
else
throw Exception() << "Unknown colour for the outgoing particles"
<< " for a decaying colour anti-sextet particle in "
<< "GeneralTwoBodyDecayer::getColourFactors() for "
<< inpart. dataPtr()->PDGName() << " -> "
<< decay[0]->dataPtr()->PDGName() << " "
<< decay[1]->dataPtr()->PDGName() << " "
<< decay[2]->dataPtr()->PDGName()
<< Exception::runerror;
}
else
throw Exception() << "Unknown colour for the decaying particle in "
<< "GeneralTwoBodyDecayer::getColourFactors() for "
<< inpart. dataPtr()->PDGName() << " -> "
<< decay[0]->dataPtr()->PDGName() << " "
<< decay[1]->dataPtr()->PDGName() << " "
<< decay[2]->dataPtr()->PDGName()
<< Exception::runerror;
return colour_;
}
const GeneralTwoBodyDecayer::CFlow &
GeneralTwoBodyDecayer::colourFlows(const Particle & inpart,
const ParticleVector & decay) {
// static initialization of commonly used colour structures
static const CFlow init = CFlow(3, CFlowPairVec(1, make_pair(0, 1.)));
static CFlow tripflow = init;
static CFlow atripflow = init;
static CFlow octflow = init;
static CFlow sexflow = init;
static CFlow epsflow = init;
static const CFlow fpflow = CFlow(4, CFlowPairVec(1, make_pair(0, 1.)));
static bool initialized = false;
if (! initialized) {
tripflow[2].resize(2, make_pair(0,1.));
tripflow[2][0] = make_pair(0, 1.);
tripflow[2][1] = make_pair(1,-1.);
tripflow[1][0] = make_pair(1, 1.);
atripflow[1].resize(2, make_pair(0,1.));
atripflow[1][0] = make_pair(0, 1.);
atripflow[1][1] = make_pair(1,-1.);
atripflow[2][0] = make_pair(1, 1.);
octflow[0].resize(2, make_pair(0,1.));
octflow[0][0] = make_pair(0,-1.);
octflow[0][1] = make_pair(1, 1.);
octflow[2][0] = make_pair(1, 1.);
sexflow[0].resize(2, make_pair(0,1.));
sexflow[0][1] = make_pair(1,1.);
sexflow[2][0] = make_pair(1,1.);
epsflow[0].resize(2, make_pair(0,-1.));
epsflow[0][0] = make_pair(0,-1.);
epsflow[0][1] = make_pair(1,-1.);
epsflow[2][0] = make_pair(1,1.);
initialized = true;
}
// main function body
int sing=0,trip=0,atrip=0,oct=0,sex=0,asex=0;
for (size_t it=0; it<decay.size(); ++it) {
switch ( decay[it]->dataPtr()->iColour() ) {
case PDT::Colour0: ++sing; break;
case PDT::Colour3: ++trip; break;
case PDT::Colour3bar: ++atrip; break;
case PDT::Colour8: ++oct; break;
case PDT::Colour6: ++sex; break;
case PDT::Colour6bar: ++asex; break;
/// @todo: handle these better
case PDT::ColourUndefined: break;
case PDT::Coloured: break;
case PDT::DarkColoured: break;
case PDT::DarkColourNeutral: break;
case PDT::DarkColourFundamental: break;
case PDT::DarkColourAntiFundamental: break;
case PDT::DarkColourAdjoint: break;
}
}
const CFlow * retval = 0;
// decaying colour triplet
if(inpart.dataPtr()->iColour() == PDT::Colour3 &&
trip==1 && oct==2) {
retval = &tripflow;
}
// decaying colour anti-triplet
else if(inpart.dataPtr()->iColour() == PDT::Colour3bar &&
atrip==1 && oct==2){
retval = &atripflow;
}
// decaying colour octet
else if(inpart.dataPtr()->iColour() == PDT::Colour8 &&
oct==1 && trip==1 && atrip==1) {
retval = &octflow;
}
// decaying colour sextet
else if(inpart.dataPtr()->iColour() ==PDT::Colour6 &&
oct==1 && trip==2) {
retval = &sexflow;
}
// decaying anti colour sextet
else if(inpart.dataPtr()->iColour() ==PDT::Colour6bar &&
oct==1 && atrip==2) {
retval = &sexflow;
}
// decaying colour triplet (eps)
else if(inpart.dataPtr()->iColour() == PDT::Colour3 &&
atrip==2 && oct==1) {
retval = &epsflow;
}
// decaying colour anti-triplet (eps)
else if(inpart.dataPtr()->iColour() == PDT::Colour3bar &&
trip==2 && oct==1){
retval = &epsflow;
}
else {
retval = &fpflow;
}
return *retval;
}
double GeneralTwoBodyDecayer::threeBodyME(const int , const Particle &,
const ParticleVector &,
ShowerInteraction, MEOption) {
throw Exception() << "Base class GeneralTwoBodyDecayer::threeBodyME() "
<< "called, should have an implementation in the inheriting class"
<< Exception::runerror;
return 0.;
}
diff --git a/Decay/General/GeneralTwoBodyDecayer.h b/Decay/General/GeneralTwoBodyDecayer.h
--- a/Decay/General/GeneralTwoBodyDecayer.h
+++ b/Decay/General/GeneralTwoBodyDecayer.h
@@ -1,295 +1,301 @@
// -*- C++ -*-
//
// GeneralTwoBodyDecayer.h is a part of Herwig - A multi-purpose Monte Carlo event generator
// Copyright (C) 2002-2019 The Herwig Collaboration
//
// Herwig is licenced under version 3 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
#ifndef HERWIG_GeneralTwoBodyDecayer_H
#define HERWIG_GeneralTwoBodyDecayer_H
//
// This is the declaration of the GeneralTwoBodyDecayer class.
//
#include "Herwig/Decay/PerturbativeDecayer.h"
#include "Herwig/Decay/PhaseSpaceMode.h"
#include "ThePEG/Helicity/Vertex/VertexBase.h"
+#include "Herwig/Models/HiddenValley/HiddenValleyModel.h"
#include "GeneralTwoBodyDecayer.fh"
namespace Herwig {
using namespace ThePEG;
using Helicity::VertexBasePtr;
/** \ingroup Decay
* The GeneralTwoBodyDecayer class is designed to be the base class
* for 2 body decays for some general model. It inherits from
* PerturbativeDecayer and implements the modeNumber() virtual function
* that is the same for all of the decays. A decayer for
* a specific spin configuration should inherit from this and implement
* the me2() and partialWidth() member functions. The colourConnections()
* member should be called from inside me2() in the inheriting decayer
* to set up the colour lines.
*
* @see \ref GeneralTwoBodyDecayerInterfaces "The interfaces"
* defined for GeneralTwoBodyDecayer.
* @see PerturbativeDecayer
*/
class GeneralTwoBodyDecayer: public PerturbativeDecayer {
public:
/** A ParticleData ptr and (possible) mass pair.*/
typedef pair<tcPDPtr, Energy> PMPair;
public:
/**
* The default constructor.
*/
GeneralTwoBodyDecayer() : maxWeight_(1.), colour_(1,DVector(1,1.))
{}
/** @name Virtual functions required by the Decayer2 class. */
//@{
/**
* For a given decay mode and a given particle instance, perform the
* decay and return the decay products. As this is the base class this
* is not implemented.
* @return The vector of particles produced in the decay.
*/
virtual ParticleVector decay(const Particle & parent,
const tPDVector & children) const;
/**
* Which of the possible decays is required
* @param cc Is this mode the charge conjugate
* @param parent The decaying particle
* @param children The decay products
*/
virtual int modeNumber(bool & cc, tcPDPtr parent,const tPDVector & children) const;
/**
* Function to return partial Width
* @param inpart The decaying particle.
* @param outa One of the decay products.
* @param outb The other decay product.
*/
virtual Energy partialWidth(PMPair inpart, PMPair outa,
PMPair outb) const;
/**
* Specify the \f$1\to2\f$ matrix element to be used in the running width
* calculation.
* @param dm The DecayMode
* @param mecode The code for the matrix element as described
* in the GenericWidthGenerator class.
* @param coupling The coupling for the matrix element.
* @return True if the the order of the particles in the
* decayer is the same as the DecayMode tag.
*/
virtual bool twoBodyMEcode(const DecayMode & dm, int & mecode,
double & coupling) const;
/**
* An overidden member to calculate a branching ratio for a certain
* particle instance.
* @param dm The DecayMode of the particle
* @param p The particle object
* @param oldbrat The branching fraction given in the DecayMode object
*/
virtual double brat(const DecayMode & dm, const Particle & p,
double oldbrat) const;
//@}
/**
* Set the information on the decay
*/
virtual void setDecayInfo(PDPtr incoming, PDPair outgoing,
vector<VertexBasePtr>,
map<ShowerInteraction,VertexBasePtr> &,
const vector<map<ShowerInteraction,VertexBasePtr> > &,
map<ShowerInteraction,VertexBasePtr>) =0;
protected:
/** @name Functions used by inheriting decayers. */
//@{
/**
* Set integration weight
* @param wgt Maximum integration weight
*/
void setWeight(double wgt) { maxWeight_ = wgt; }
/**
* Set colour connections
* @param parent Parent particle
* @param out Particle vector containing particles to
* connect colour lines
*/
void colourConnections(const Particle & parent,
const ParticleVector & out) const;
/**
* Compute the spin and colour factor
*/
double colourFactor(tcPDPtr in, tcPDPtr out1, tcPDPtr out2) const;
/**
* Calculate matrix element ratio R/B
*/
double matrixElementRatio(const Particle & inpart, const ParticleVector & decay2,
const ParticleVector & decay3, MEOption meopt,
ShowerInteraction inter);
/**
* Set the information on the decay
*/
void decayInfo(PDPtr incoming, PDPair outgoing);
//@}
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();
/**
* Initialize this object. Called in the run phase just before
* a run begins.
*/
virtual void doinitrun();
//@}
protected:
/**
* Member for the generation of additional hard radiation
*/
//@{
/**
* Return the matrix of colour factors
*/
typedef vector<pair<int,double > > CFlowPairVec;
typedef vector<CFlowPairVec> CFlow;
const vector<DVector> & getColourFactors(const Particle & inpart,
const ParticleVector & decay,
unsigned int & nflow);
const CFlow & colourFlows(const Particle & inpart,
const ParticleVector & decay);
/**
* Three-body matrix element including additional QCD radiation
*/
virtual double threeBodyME(const int , const Particle & inpart,
const ParticleVector & decay,
ShowerInteraction inter, MEOption meopt);
//@}
private:
/**
* The assignment operator is private and must never be called.
* In fact, it should not even be implemented.
*/
GeneralTwoBodyDecayer & operator=(const GeneralTwoBodyDecayer &) = delete;
private:
/**
* Store the incoming particle
*/
PDPtr incoming_;
/**
* Outgoing particles
*/
vector<PDPtr> outgoing_;
/**
* Maximum weight for integration
*/
double maxWeight_;
/**
* Store colour factors for ME calc.
*/
vector<DVector> colour_;
+ /**
+ * Pointer to the Hiddwn Valley object.
+ */
+ tHiddenValleyPtr _HVmodel;
+
};
/**
* Write a map with ShowerInteraction as the key
*/
template<typename T, typename Cmp, typename A>
inline PersistentOStream & operator<<(PersistentOStream & os,
const map<ShowerInteraction,T,Cmp,A> & m) {
os << m.size();
if(m.find(ShowerInteraction::QCD)!=m.end()) {
os << 0 << m.at(ShowerInteraction::QCD);
}
if(m.find(ShowerInteraction::QED)!=m.end()) {
os << 1 << m.at(ShowerInteraction::QED);
}
return os;
}
/**
* Read a map with ShowerInteraction as the key
*/
template <typename T, typename Cmp, typename A>
inline PersistentIStream & operator>>(PersistentIStream & is, map<ShowerInteraction,T,Cmp,A> & m) {
m.clear();
long size;
int k;
is >> size;
while ( size-- && is ) {
is >> k;
if(k==0)
is >> m[ShowerInteraction::QCD];
else if(k==1)
is >> m[ShowerInteraction::QED];
else
assert(false);
}
return is;
}
}
#endif /* HERWIG_GeneralTwoBodyDecayer_H */
diff --git a/Makefile.am b/Makefile.am
--- a/Makefile.am
+++ b/Makefile.am
@@ -1,29 +1,29 @@
SUBDIRS = include \
- Utilities PDT Decay PDF Models \
+ Utilities PDT PDF Models Decay \
Shower Hadronization MatrixElement \
UnderlyingEvent Analysis Looptools Sampling \
API lib src Doc Contrib Tests
EXTRA_DIST = GUIDELINES cat_with_cpplines
DISTCHECK_CONFIGURE_FLAGS = --enable-debug --with-thepeg=$(THEPEGPATH) --with-gsl=$(GSLPATH)
ACLOCAL_AMFLAGS = -I m4
DISTCLEANFILES = config.herwig
libclean:
find . -name '*.la' -delete
cd lib && $(MAKE) $(AM_MAKEFLAGS) clean
cd src && $(MAKE) $(AM_MAKEFLAGS) clean
tests:
cd Tests && $(MAKE) $(AM_MAKEFLAGS) tests
## ThePEG registration
unregister:
cd src && $(MAKE) $(AM_MAKEFLAGS) unregister
register:
cd src && $(MAKE) $(AM_MAKEFLAGS) register
diff --git a/Models/HiddenValley/HiddenValleyFFZPrimeVertex.cc b/Models/HiddenValley/HiddenValleyFFZPrimeVertex.cc
new file mode 100644
--- /dev/null
+++ b/Models/HiddenValley/HiddenValleyFFZPrimeVertex.cc
@@ -0,0 +1,95 @@
+// -*- C++ -*-
+//
+// This is the implementation of the non-inlined, non-templated member
+// functions of the HiddenValleyFFZPrimeVertex class.
+//
+
+#include "HiddenValleyFFZPrimeVertex.h"
+#include "ThePEG/Interface/ClassDocumentation.h"
+#include "ThePEG/Persistency/PersistentOStream.h"
+#include "ThePEG/Persistency/PersistentIStream.h"
+#include "HiddenValleyModel.h"
+#include "ThePEG/Utilities/DescribeClass.h"
+
+using namespace Herwig;
+
+void HiddenValleyFFZPrimeVertex::persistentOutput(PersistentOStream & os) const {
+ os << _gl << _gr << _gql << _gqr << _gPrime;
+}
+
+void HiddenValleyFFZPrimeVertex::persistentInput(PersistentIStream & is, int) {
+ is >> _gl >> _gr >> _gql >> _gqr >> _gPrime;
+}
+
+// Static variable needed for the type description system in ThePEG.
+DescribeClass<HiddenValleyFFZPrimeVertex,FFVVertex>
+describeHerwigHiddenValleyFFZPrimeVertex("Herwig::HiddenValleyFFZPrimeVertex", "Herwig.so");
+
+void HiddenValleyFFZPrimeVertex::Init() {
+
+ static ClassDocumentation<HiddenValleyFFZPrimeVertex> documentation
+ ("There is no documentation for the HiddenValleyFFZPrimeVertex class");
+
+}
+
+void HiddenValleyFFZPrimeVertex::setCoupling(Energy2, tcPDPtr a,tcPDPtr,tcPDPtr) {
+ norm(_gPrime);
+ // the left and right couplings
+ int iferm=abs(a->id());
+ if((iferm>=1 && iferm<=6)||(iferm>=11 &&iferm<=16)) {
+ left (_gl[iferm]);
+ right(_gr[iferm]);
+ return;
+ }
+ iferm -= ParticleID::darkg;
+ if(iferm<=int(_gql.size())) {
+ left (_gql[iferm]);
+ right(_gqr[iferm]);
+ return;
+ }
+ assert(false);
+}
+
+HiddenValleyFFZPrimeVertex::HiddenValleyFFZPrimeVertex() : _gl(17,0.0), _gr(17,0.0) {
+ orderInGem(1);
+ orderInGs(0);
+ // PDG codes for the particles
+ // the quarks
+ for(int ix=1;ix<7;++ix)
+ addToList(-ix,ix,32);
+ // the leptons
+ for(int ix=11;ix<17;++ix)
+ addToList(-ix,ix,32);
+}
+
+void HiddenValleyFFZPrimeVertex::doinit() {
+ tcHiddenValleyPtr model =
+ dynamic_ptr_cast<tcHiddenValleyPtr>(generator()->standardModel());
+ if(!model)
+ throw InitException() << "Must be using the HiddenValleyModel in "
+ << "HiddenValleyFFZPrimeVertex::doinit()"
+ << Exception::runerror;
+ _gPrime = model->gPrime();
+ // SM couplings
+ _gl.resize(17,0.);
+ _gr.resize(17,0.);
+ for(int ix=1;ix<4;++ix) {
+ _gl[2*ix-1] = model->qL();
+ _gl[2*ix ] = model->qL();
+ _gl[2*ix+9 ] = model->lL();
+ _gl[2*ix+10] = model->lL();
+ _gr[2*ix-1] = model->dR();
+ _gr[2*ix ] = model->uR();
+ _gr[2*ix+9 ] = model->lR();
+ _gr[2*ix+10] = 0. ;
+ }
+ FFVVertex::doinit();
+ _gql.resize(model->NF()+1,0.);
+ _gqr.resize(model->NF()+1,0.);
+ for(int ix=0;ix<int(model->NF());++ix) {
+ int id = ParticleID::darkg+1+ix;
+ addToList(-id,id,32);
+ _gql[ix+1] = model->qCharge()[ix];
+ _gqr[ix+1] = model->qCharge()[ix]-2.;
+ }
+}
diff --git a/Models/HiddenValley/HiddenValleyFFZPrimeVertex.h b/Models/HiddenValley/HiddenValleyFFZPrimeVertex.h
new file mode 100644
--- /dev/null
+++ b/Models/HiddenValley/HiddenValleyFFZPrimeVertex.h
@@ -0,0 +1,132 @@
+// -*- C++ -*-
+#ifndef RADIATIVEZPRIME_HiddenValleyFFZPrimeVertex_H
+#define RADIATIVEZPRIME_HiddenValleyFFZPrimeVertex_H
+//
+// This is the declaration of the HiddenValleyFFZPrimeVertex class.
+//
+
+#include "ThePEG/Helicity/Vertex/Vector/FFVVertex.h"
+
+namespace Herwig {
+
+using namespace ThePEG;
+
+/**
+ * Here is the documentation of the HiddenValleyFFZPrimeVertex class.
+ *
+ * @see \ref HiddenValleyFFZPrimeVertexInterfaces "The interfaces"
+ * defined for HiddenValleyFFZPrimeVertex.
+ */
+class HiddenValleyFFZPrimeVertex: public Helicity::FFVVertex {
+
+public:
+
+ /**
+ * The default constructor.
+ */
+ inline HiddenValleyFFZPrimeVertex();
+
+ /**
+ * Calculate the couplings.
+ * @param q2 The scale \f$q^2\f$ for the coupling at the vertex.
+ * @param part1 The ParticleData pointer for the first particle.
+ * @param part2 The ParticleData pointer for the second particle.
+ * @param part3 The ParticleData pointer for the third particle.
+ */
+ virtual void setCoupling(Energy2 q2,tcPDPtr part1,tcPDPtr part2,tcPDPtr part3);
+
+public:
+
+ /** @name Functions used by the persistent I/O system. */
+ //@{
+ /**
+ * Function used to write out object persistently.
+ * @param os the persistent output stream written to.
+ */
+ void persistentOutput(PersistentOStream & os) const;
+
+ /**
+ * Function used to read in object persistently.
+ * @param is the persistent input stream read from.
+ * @param version the version number of the object when written.
+ */
+ void persistentInput(PersistentIStream & is, int version);
+ //@}
+
+ /**
+ * The standard Init function used to initialize the interfaces.
+ * Called exactly once for each class by the class description system
+ * before the main function starts or
+ * when this class is dynamically loaded.
+ */
+ static void Init();
+
+protected:
+
+ /** @name Clone Methods. */
+ //@{
+ /**
+ * Make a simple clone of this object.
+ * @return a pointer to the new object.
+ */
+ inline 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.
+ */
+ inline virtual IBPtr fullclone() const {return new_ptr(*this);}
+ //@}
+
+protected:
+
+ /**
+ * Initialize this object after the setup phase before saving and
+ * EventGenerator to disk.
+ * @throws InitException if object could not be initialized properly.
+ */
+ virtual void doinit();
+
+private:
+
+ /**
+ * The assignment operator is private and must never be called.
+ * In fact, it should not even be implemented.
+ */
+ HiddenValleyFFZPrimeVertex & operator=(const HiddenValleyFFZPrimeVertex &);
+
+
+ /**
+ * Storage of the couplings.
+ */
+ //@{
+ /**
+ * The left couplings of the Standard Model fermions.
+ */
+ vector<double> _gl;
+
+ /**
+ * The right couplings of the Standard Model fermions.
+ */
+ vector<double> _gr;
+ /**
+ * The left couplings of the new fermions.
+ */
+ vector<double> _gql;
+
+ /**
+ * The right couplings of the new fermions.
+ */
+ vector<double> _gqr;
+
+ /**
+ * Coupling
+ */
+ double _gPrime;
+ //@}
+
+};
+
+}
+
+#endif /* RADIATIVEZPRIME_HiddenValleyFFZPrimeVertex_H */
diff --git a/Models/HiddenValley/HiddenValleyModel.cc b/Models/HiddenValley/HiddenValleyModel.cc
new file mode 100644
--- /dev/null
+++ b/Models/HiddenValley/HiddenValleyModel.cc
@@ -0,0 +1,132 @@
+// -*- C++ -*-
+//
+// This is the implementation of the non-inlined, non-templated member
+// functions of the HiddenValleyModel class.
+//
+
+#include "HiddenValleyModel.h"
+#include "ThePEG/Interface/ClassDocumentation.h"
+#include "ThePEG/Interface/Reference.h"
+#include "ThePEG/Interface/Parameter.h"
+#include "ThePEG/Interface/ParVector.h"
+#include "ThePEG/Interface/Switch.h"
+#include "ThePEG/Utilities/EnumIO.h"
+#include "ThePEG/Persistency/PersistentOStream.h"
+#include "ThePEG/Persistency/PersistentIStream.h"
+
+using namespace Herwig;
+
+HiddenValleyModel::HiddenValleyModel() : groupType_(SU), Nc_(3), Nf_(1),
+ qL_(-0.2), uR_(-0.2), dR_(0.6),
+ lL_(0.6), lR_(-0.2), gPrime_(1./7.),
+ qCharge_(1,-0.2)
+{}
+
+HiddenValleyModel::~HiddenValleyModel() {}
+
+IBPtr HiddenValleyModel::clone() const {
+ return new_ptr(*this);
+}
+
+IBPtr HiddenValleyModel::fullclone() const {
+ return new_ptr(*this);
+}
+
+void HiddenValleyModel::persistentOutput(PersistentOStream & os) const {
+ os << oenum(groupType_) << Nc_ << Nf_ << FFZPVertex_
+ << qL_ << uR_ << dR_ << lL_ << lR_ << qCharge_ << gPrime_;
+}
+
+void HiddenValleyModel::persistentInput(PersistentIStream & is, int) {
+ is >> ienum(groupType_) >> Nc_ >> Nf_ >> FFZPVertex_
+ >> qL_ >> uR_ >> dR_ >> lL_ >> lR_ >> qCharge_ >> gPrime_;
+}
+
+ClassDescription<HiddenValleyModel> HiddenValleyModel::initHiddenValleyModel;
+// Definition of the static class description member.
+
+void HiddenValleyModel::Init() {
+
+ static ClassDocumentation<HiddenValleyModel> documentation
+ ("The HiddenValleyModel class is the main class for the Hidden Valley Model.");
+
+ static Switch<HiddenValleyModel,GroupType> interfaceGroupType
+ ("GroupType",
+ "Type of the new unbroken group",
+ &HiddenValleyModel::groupType_, SU, false, false);
+ static SwitchOption interfaceGroupTypeSU
+ (interfaceGroupType,
+ "SU",
+ "Use an SU(N) group",
+ SU);
+
+ static Parameter<HiddenValleyModel,unsigned int> interfaceGroupOrder
+ ("GroupOrder",
+ "The value of the number of 'colours' for the new symmetry group",
+ &HiddenValleyModel::Nc_, 3, 1, 1000,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyModel,unsigned int> interfaceNumberOfFermions
+ ("NumberOfFermions",
+ "The number of fermions charged under the new group",
+ &HiddenValleyModel::Nf_, 1, 1, 100,
+ false, false, Interface::limited);
+
+ static Reference<HiddenValleyModel,AbstractFFVVertex> interfaceVertexFFZPrime
+ ("Vertex/FFZPrime",
+ "The vertex coupling the Zprime to the fermions",
+ &HiddenValleyModel::FFZPVertex_, false, false, true, false, false);
+
+ static ParVector<HiddenValleyModel,double> interfaceQuirkCharges
+ ("QuirkCharges",
+ "The charges under the new U(1) of the new quirks",
+ &HiddenValleyModel::qCharge_, -1, -0.2, -10., 10.0,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyModel,double> interfaceZPrimeQL
+ ("ZPrime/QL",
+ "The charge of the left-handed quarks under the new U(1)",
+ &HiddenValleyModel::qL_, -0.2, 10.0, 10.0,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyModel,double> interfaceZPrimeUR
+ ("ZPrime/UR",
+ "The charge of the right-handed up quarks under the new U(1)",
+ &HiddenValleyModel::uR_, -0.2, 10.0, 10.0,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyModel,double> interfaceZPrimeDR
+ ("ZPrime/DR",
+ "The charge of the right-handed down quarks under the new U(1)",
+ &HiddenValleyModel::uR_, 0.6, 10.0, 10.0,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyModel,double> interfaceZPrimeLL
+ ("ZPrime/LL",
+ "The charge of the left-handed leptons under the new U(1)",
+ &HiddenValleyModel::lL_, 0.6, 10.0, 10.0,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyModel,double> interfaceZPrimeLR
+ ("ZPrime/LR",
+ "The charge of the right-handed charged leptons under the new U(1)",
+ &HiddenValleyModel::lR_, -0.2, 10.0, 10.0,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyModel,double> interfaceGPrime
+ ("GPrime",
+ "The new gauge coupling",
+ &HiddenValleyModel::gPrime_, 1./7., 0.0, 10.0,
+ false, false, Interface::limited);
+
+}
+
+void HiddenValleyModel::doinit() {
+ addVertex(FFZPVertex_);
+ BSMModel::doinit();
+ FFZPVertex_->init();
+ if(qCharge_.size()!=Nf_)
+ throw InitException() << "Number of new fermions and number of new charges "
+ << "different in HiddenValleyModel::doinit()"
+ << Exception::runerror;
+}
diff --git a/Models/HiddenValley/HiddenValleyModel.h b/Models/HiddenValley/HiddenValleyModel.h
new file mode 100644
--- /dev/null
+++ b/Models/HiddenValley/HiddenValleyModel.h
@@ -0,0 +1,346 @@
+// -*- C++ -*-
+#ifndef HERWIG_HiddenValleyModel_H
+#define HERWIG_HiddenValleyModel_H
+//
+// This is the declaration of the HiddenValleyModel class.
+//
+
+#include "Herwig/Models/General/BSMModel.h"
+#include "ThePEG/Helicity/Vertex/AbstractFFVVertex.h"
+
+namespace Herwig {
+
+using namespace ThePEG;
+
+/**
+ * Here is the documentation of the HiddenValleyModel class.
+ *
+ * @see \ref HiddenValleyModelInterfaces "The interfaces"
+ * defined for HiddenValleyModel.
+ */
+class HiddenValleyModel: public BSMModel {
+
+/**
+ * Enumeration to define the type of the confining group
+ */
+ enum GroupType { SU, SO, SP };
+
+public:
+
+ /** @name Standard constructors and destructors. */
+ //@{
+ /**
+ * Default constructor
+ */
+ HiddenValleyModel();
+
+ /**
+ * Destructor
+ */
+ virtual ~HiddenValleyModel();
+ //@}
+
+public:
+
+ /**
+ * Properties of the new confining gauge group
+ */
+ //@{
+ /**
+ * The hidden colour charge of the fundamental representation
+ */
+ double CF() const {
+ switch (groupType_) {
+ case SU:
+ return 0.5*(sqr(double(Nc_))-1.)/double(Nc_);
+ case SO:
+ if(Nc_==3) {
+ return (Nc_-1);
+ }
+ else {
+ return 0.5*(Nc_-1);
+ }
+ case SP:
+ return 0.25*(Nc_+1);
+ default:
+ assert(false);
+ }
+ return 0.;
+ }
+
+ /**
+ * The hidden colour charge of the adjoint representation
+ */
+ double CA() const {
+ switch (groupType_) {
+ case SU:
+ return double(Nc_);
+ case SO:
+ if(Nc_==3) {
+ return 2.*(Nc_-2.);
+ }
+ else {
+ return (Nc_-2.);
+ }
+ case SP:
+ return 0.5*(Nc_+2.);
+ default:
+ assert(false);
+ }
+ return 0.;
+ }
+
+ /**
+ * The \f$T_R\f$ colour factor
+ */
+ double TR() const {
+ switch (groupType_) {
+ case SU: case SP:
+ return 0.5;
+ break;
+ case SO:
+ if(Nc_==3) {
+ return 2.;
+ }
+ else {
+ return 1.;
+ }
+ default:
+ assert(false);
+ }
+ return 0.;
+ }
+
+ /**
+ * The type of the group
+ */
+ GroupType groupType() const {return groupType_;}
+
+ /**
+ * The number of colours
+ */
+ unsigned int NC() const {return Nc_;}
+
+ /**
+ * The number of fermions charged under the new group
+ */
+ unsigned int NF() const {return Nf_;}
+ //@}
+
+ /**
+ * Properties of the new U(1) group
+ */
+ //@{
+ /**
+ * Charge of the left-handed quarks
+ */
+ double qL() const {return qL_;}
+
+ /**
+ * Charge of the right-handed up-type quarks
+ */
+ double uR() const {return uR_;}
+
+ /**
+ * Charge of the right-handed down-type quarks
+ */
+ double dR() const {return dR_;}
+
+ /**
+ * Charge of the left-handed leptons
+ */
+ double lL() const {return lL_;}
+
+ /**
+ * Charge of the right-handed charged leptons
+ */
+ double lR() const {return lR_;}
+
+ /**
+ * Coupling
+ */
+ double gPrime() const {return gPrime_;}
+
+ /**
+ * The couplings of the quirks
+ */
+ vector<double> qCharge() const {return qCharge_;}
+ //@}
+
+public:
+
+ /** @name Functions used by the persistent I/O system. */
+ //@{
+ /**
+ * Function used to write out object persistently.
+ * @param os the persistent output stream written to.
+ */
+ void persistentOutput(PersistentOStream & os) const;
+
+ /**
+ * Function used to read in object persistently.
+ * @param is the persistent input stream read from.
+ * @param version the version number of the object when written.
+ */
+ void persistentInput(PersistentIStream & is, int version);
+ //@}
+
+ /**
+ * The standard Init function used to initialize the interfaces.
+ * Called exactly once for each class by the class description system
+ * before the main function starts or
+ * when this class is dynamically loaded.
+ */
+ static void Init();
+
+protected:
+
+ /** @name Clone Methods. */
+ //@{
+ /**
+ * Make a simple clone of this object.
+ * @return a pointer to the new object.
+ */
+ virtual IBPtr clone() const;
+
+ /** Make a clone of this object, possibly modifying the cloned object
+ * to make it sane.
+ * @return a pointer to the new object.
+ */
+ virtual IBPtr fullclone() const;
+ //@}
+
+private:
+
+ /**
+ * The static object used to initialize the description of this class.
+ * Indicates that this is a concrete class with persistent data.
+ */
+ static ClassDescription<HiddenValleyModel> initHiddenValleyModel;
+
+ /**
+ * The assignment operator is private and must never be called.
+ * In fact, it should not even be implemented.
+ */
+ HiddenValleyModel & operator=(const HiddenValleyModel &);
+
+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:
+
+ /**
+ * Symmetry group for the new unbroken symmetry
+ */
+ //@{
+ /**
+ * Type of group
+ */
+ GroupType groupType_;
+
+ /**
+ * Order of the group
+ */
+ unsigned int Nc_;
+ //@}
+
+ /**
+ * Number of fermions
+ */
+ unsigned int Nf_;
+
+ /**
+ * Charges etc under the new \f$U(1)\f$ group
+ */
+ //@{
+ /**
+ * Charge of the left-handed quarks
+ */
+ double qL_;
+
+ /**
+ * Charge of the right-handed up-type quarks
+ */
+ double uR_;
+
+ /**
+ * Charge of the right-handed down-type quarks
+ */
+ double dR_;
+
+ /**
+ * Charge of the left-handed leptons
+ */
+ double lL_;
+
+ /**
+ * Charge of the right-handed charged leptons
+ */
+ double lR_;
+
+ /**
+ * Coupling
+ */
+ double gPrime_;
+
+ /**
+ * the vertex
+ */
+ AbstractFFVVertexPtr FFZPVertex_;
+
+ /**
+ * The couplings of the quirks
+ */
+ vector<double> qCharge_;
+ //@}
+};
+
+}
+
+#include "ThePEG/Utilities/ClassTraits.h"
+
+namespace ThePEG {
+
+/** @cond TRAITSPECIALIZATIONS */
+
+/** This template specialization informs ThePEG about the
+ * base classes of HiddenValleyModel. */
+template <>
+struct BaseClassTrait<Herwig::HiddenValleyModel,1> {
+ /** Typedef of the first base class of HiddenValleyModel. */
+ typedef Herwig::BSMModel NthBase;
+};
+
+/** This template specialization informs ThePEG about the name of
+ * the HiddenValleyModel class and the shared object where it is defined. */
+template <>
+struct ClassTraits<Herwig::HiddenValleyModel>
+ : public ClassTraitsBase<Herwig::HiddenValleyModel> {
+ /** Return a platform-independent class name */
+ static string className() { return "Herwig::HiddenValleyModel"; }
+ /**
+ * The name of a file containing the dynamic library where the class
+ * HiddenValleyModel is implemented. It may also include several, space-separated,
+ * libraries if the class HiddenValleyModel depends on other classes (base classes
+ * excepted). In this case the listed libraries will be dynamically
+ * linked in the order they are specified.
+ */
+ static string library() { return "Herwig.so"; }
+};
+
+ThePEG_DECLARE_POINTERS(Herwig::HiddenValleyModel,HiddenValleyPtr);
+
+/** @endcond */
+
+}
+
+#endif /* HERWIG_HiddenValleyModel_H */
diff --git a/Models/HiddenValley/Makefile.am b/Models/HiddenValley/Makefile.am
new file mode 100644
--- /dev/null
+++ b/Models/HiddenValley/Makefile.am
@@ -0,0 +1,17 @@
+BUILT_SOURCES = HV__all.cc
+CLEANFILES = HV__all.cc
+
+HV__all.cc : $(DIR_H_FILES) $(DIR_CC_FILES) Makefile
+ @echo "Concatenating .cc files into $@"
+ @$(top_srcdir)/cat_with_cpplines $(DIR_CC_FILES) > $@
+
+EXTRA_DIST = $(ALL_H_FILES) $(ALL_CC_FILES)
+
+DIR_H_FILES = $(addprefix $(srcdir)/,$(ALL_H_FILES))
+ALL_H_FILES = \
+HiddenValleyModel.h HiddenValleyFFZPrimeVertex.h
+
+DIR_CC_FILES = $(addprefix $(srcdir)/,$(ALL_CC_FILES))
+ALL_CC_FILES = \
+HiddenValleyModel.cc \
+HiddenValleyFFZPrimeVertex.cc
\ No newline at end of file
diff --git a/Models/Makefile.am b/Models/Makefile.am
--- a/Models/Makefile.am
+++ b/Models/Makefile.am
@@ -1,193 +1,199 @@
SUBDIRS = RSModel StandardModel General Susy UED Zprime \
Transplanckian ADD Leptoquarks Sextet TTbAsymm \
- LH LHTP Feynrules DarkMatter
+ LH LHTP Feynrules DarkMatter HiddenValley
noinst_LTLIBRARIES = libHwStandardModel.la
nodist_libHwStandardModel_la_SOURCES = \
StandardModel/SM__all.cc
libHwStandardModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/StandardModel
noinst_LTLIBRARIES += libHwModelGenerator.la
nodist_libHwModelGenerator_la_SOURCES = \
General/ModelGenerator__all.cc
libHwModelGenerator_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/General
+noinst_LTLIBRARIES += libHwHiddenValley.la
+nodist_libHwHiddenValley_la_SOURCES = \
+HiddenValley/HV__all.cc
+
+libHwHiddenValley_la_CPPFLAGS = \
+$(AM_CPPFLAGS) -I$(srcdir)/HiddenValley
if WANT_BSM
pkglib_LTLIBRARIES = \
HwRSModel.la \
HwUED.la \
HwSusy.la \
HwNMSSM.la \
HwRPV.la \
HwZprimeModel.la \
HwTransplanck.la \
HwADDModel.la \
HwLeptoquarkModel.la \
HwSextetModel.la \
HwTTbAModel.la \
HwLHModel.la \
HwDMModel.la \
HwLHTPModel.la
endif
#############
HwRSModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 12:0:0
HwRSModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/RSModel
nodist_HwRSModel_la_SOURCES = \
RSModel/RS__all.cc
#############
HwUED_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 12:0:0
HwUED_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/UED
nodist_HwUED_la_SOURCES = \
UED/UED__all.cc
#############
HwSusy_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 15:0:0
HwSusy_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/Susy
nodist_HwSusy_la_SOURCES = \
Susy/Susy__all.cc
#############
HwNMSSM_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 7:0:0
HwNMSSM_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/Susy/NMSSM
nodist_HwNMSSM_la_SOURCES = \
Susy/NMSSM/NMSSM__all.cc
#############
HwRPV_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 5:0:0
HwRPV_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/Susy/RPV
nodist_HwRPV_la_SOURCES = \
Susy/RPV/RPV__all.cc
#############
HwZprimeModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 5:0:0
HwZprimeModel_la_SOURCES = \
Zprime/ZprimeModel.cc Zprime/ZprimeModelZPQQVertex.cc
#############
HwTransplanck_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 6:0:0
HwTransplanck_la_SOURCES = \
Transplanckian/METRP2to2.cc
#############
HwADDModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 6:0:0
HwADDModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/ADD
nodist_HwADDModel_la_SOURCES = \
ADD/ADD__all.cc
#############
HwLeptoquarkModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 7:0:0
HwLeptoquarkModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/Leptoquarks
nodist_HwLeptoquarkModel_la_SOURCES = \
Leptoquarks/Leptoquark__all.cc
#############
HwSextetModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 5:0:0
HwSextetModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/Sextet
nodist_HwSextetModel_la_SOURCES = \
Sextet/Sextet__all.cc
#############
HwTTbAModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 5:0:0
HwTTbAModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/TTbAsymm
nodist_HwTTbAModel_la_SOURCES = \
TTbAsymm/TTbA__all.cc
#############
HwLHModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 7:0:0
HwLHModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/LH
nodist_HwLHModel_la_SOURCES = \
LH/LH__all.cc
#############
HwDMModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 7:0:0
HwDMModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) -I$(srcdir)/DarkMatter
nodist_HwDMModel_la_SOURCES = \
DarkMatter/DM__all.cc
#############
HwLHTPModel_la_LDFLAGS = \
$(AM_LDFLAGS) -module -version-info 7:0:0
HwLHTPModel_la_LIBADD = \
$(GSLLIBS)
HwLHTPModel_la_CPPFLAGS = \
$(AM_CPPFLAGS) $(GSLINCLUDE) -I$(srcdir)/LHTP
nodist_HwLHTPModel_la_SOURCES = \
LHTP/LHTP__all.cc
#############
diff --git a/Shower/QTilde/Couplings/HiddenValleyAlpha.cc b/Shower/QTilde/Couplings/HiddenValleyAlpha.cc
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/Couplings/HiddenValleyAlpha.cc
@@ -0,0 +1,333 @@
+// -*- C++ -*-
+//
+// HiddenValleyAlpha.cc 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.
+//
+//
+// This is the implementation of the non-inlined, non-templated member
+// functions of the HiddenValleyAlpha class.
+//
+#include "HiddenValleyAlpha.h"
+#include "Herwig/Models/HiddenValley/HiddenValleyModel.h"
+#include "ThePEG/PDT/EnumParticles.h"
+#include "ThePEG/PDT/ParticleData.h"
+#include "ThePEG/Interface/ClassDocumentation.h"
+#include "ThePEG/Interface/Switch.h"
+#include "ThePEG/Interface/Parameter.h"
+#include "ThePEG/Interface/Command.h"
+#include "ThePEG/Persistency/PersistentOStream.h"
+#include "ThePEG/Persistency/PersistentIStream.h"
+#include "ThePEG/Utilities/Throw.h"
+
+using namespace Herwig;
+
+IBPtr HiddenValleyAlpha::clone() const {
+ return new_ptr(*this);
+}
+
+IBPtr HiddenValleyAlpha::fullclone() const {
+ return new_ptr(*this);
+}
+
+void HiddenValleyAlpha::persistentOutput(PersistentOStream & os) const {
+ os << _asType << _asMaxNP << ounit(_qmin,GeV) << _nloop << _thresopt
+ << ounit(_lambdain,GeV) << _tolerance << _maxtry << _alphamin
+ << ounit(_thresholds,GeV) << ounit(_lambda,GeV) << _ca << _cf << _tr;
+}
+
+void HiddenValleyAlpha::persistentInput(PersistentIStream & is, int) {
+ is >> _asType >> _asMaxNP >> iunit(_qmin,GeV) >> _nloop >> _thresopt
+ >> iunit(_lambdain,GeV) >> _tolerance >> _maxtry >> _alphamin
+ >> iunit(_thresholds,GeV) >> iunit(_lambda,GeV) >> _ca >> _cf >> _tr;
+}
+
+ClassDescription<HiddenValleyAlpha> HiddenValleyAlpha::initHiddenValleyAlpha;
+// Definition of the static class description member.
+
+void HiddenValleyAlpha::Init() {
+
+ static ClassDocumentation<HiddenValleyAlpha> documentation
+ ("This (concrete) class describes the QCD alpha running.");
+
+ static Switch<HiddenValleyAlpha, int> intAsType
+ ("NPAlphaS",
+ "Behaviour of AlphaS in the NP region",
+ &HiddenValleyAlpha::_asType, 1, false, false);
+ static SwitchOption intAsTypeZero
+ (intAsType, "Zero","zero below Q_min", 1);
+ static SwitchOption intAsTypeConst
+ (intAsType, "Const","const as(qmin) below Q_min", 2);
+ static SwitchOption intAsTypeLin
+ (intAsType, "Linear","growing linearly below Q_min", 3);
+ static SwitchOption intAsTypeQuad
+ (intAsType, "Quadratic","growing quadratically below Q_min", 4);
+ static SwitchOption intAsTypeExx1
+ (intAsType, "Exx1", "quadratic from AlphaMaxNP down to as(Q_min)", 5);
+ static SwitchOption intAsTypeExx2
+ (intAsType, "Exx2", "const = AlphaMaxNP below Q_min", 6);
+
+ // default such that as(qmin) = 1 in the current parametrization.
+ // min = Lambda3
+ static Parameter<HiddenValleyAlpha,Energy> intQmin
+ ("Qmin", "Q < Qmin is treated with NP parametrization as of (unit [GeV])",
+ &HiddenValleyAlpha::_qmin, GeV, 0.630882*GeV, 0.330445*GeV,
+ 100.0*GeV,false,false,false);
+
+ static Parameter<HiddenValleyAlpha,double> interfaceAlphaMaxNP
+ ("AlphaMaxNP",
+ "Max value of alpha in NP region, only relevant if NPAlphaS = 5,6",
+ &HiddenValleyAlpha::_asMaxNP, 1.0, 0., 100.0,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyAlpha,unsigned int> interfaceNumberOfLoops
+ ("NumberOfLoops",
+ "The number of loops to use in the alpha_S calculation",
+ &HiddenValleyAlpha::_nloop, 2, 1, 2,
+ false, false, Interface::limited);
+
+ static Switch<HiddenValleyAlpha,bool> interfaceLambdaOption
+ ("LambdaOption",
+ "Option for the calculation of the Lambda used in the simulation from the input"
+ " Lambda_MSbar",
+ &HiddenValleyAlpha::_lambdaopt, false, false, false);
+ static SwitchOption interfaceLambdaOptionfalse
+ (interfaceLambdaOption,
+ "Same",
+ "Use the same value",
+ false);
+ static SwitchOption interfaceLambdaOptionConvert
+ (interfaceLambdaOption,
+ "Convert",
+ "Use the conversion to the Herwig scheme from NPB349, 635",
+ true);
+
+ static Parameter<HiddenValleyAlpha,Energy> interfaceLambdaDark
+ ("LambdaDark",
+ "Input value of Lambda_MSBar",
+ &HiddenValleyAlpha::_lambdain, GeV, 0.208364*GeV, 100.0*MeV, 100.0*GeV,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyAlpha,double> interfaceTolerance
+ ("Tolerance",
+ "The tolerance for discontinuities in alphaS at thresholds.",
+ &HiddenValleyAlpha::_tolerance, 1e-10, 1e-20, 1e-4,
+ false, false, Interface::limited);
+
+ static Parameter<HiddenValleyAlpha,unsigned int> interfaceMaximumIterations
+ ("MaximumIterations",
+ "The maximum number of iterations for the Newton-Raphson method to converge.",
+ &HiddenValleyAlpha::_maxtry, 100, 10, 1000,
+ false, false, Interface::limited);
+
+ static Switch<HiddenValleyAlpha,bool> interfaceThresholdOption
+ ("ThresholdOption",
+ "Whether to use the consistuent or normal masses for the thresholds",
+ &HiddenValleyAlpha::_thresopt, true, false, false);
+ static SwitchOption interfaceThresholdOptionCurrent
+ (interfaceThresholdOption,
+ "Current",
+ "Use the current masses",
+ true);
+ static SwitchOption interfaceThresholdOptionConstituent
+ (interfaceThresholdOption,
+ "Constituent",
+ "Use the constitent masses.",
+ false);
+}
+
+double HiddenValleyAlpha::value(const Energy2 scale) const {
+ pair<short,Energy> nflam;
+ Energy q = sqrt(scale);
+ double val(0.);
+ // special handling if the scale is less than Qmin
+ if (q < _qmin) {
+ nflam = getLamNfTwoLoop(_qmin);
+ double val0 = alphaS(_qmin, nflam.second, nflam.first);
+ switch (_asType) {
+ case 1:
+ // flat, zero; the default type with no NP effects.
+ val = 0.;
+ break;
+ case 2:
+ // flat, non-zero alpha_s = alpha_s(q2min).
+ val = val0;
+ break;
+ case 3:
+ // linear in q
+ val = val0*q/_qmin;
+ break;
+ case 4:
+ // quadratic in q
+ val = val0*sqr(q/_qmin);
+ break;
+ case 5:
+ // quadratic in q, starting off at asMaxNP, ending on as(qmin)
+ val = (val0 - _asMaxNP)*sqr(q/_qmin) + _asMaxNP;
+ break;
+ case 6:
+ // just asMaxNP and constant
+ val = _asMaxNP;
+ break;
+ }
+ }
+ else {
+ // the 'ordinary' case
+ nflam = getLamNfTwoLoop(q);
+ val = alphaS(q, nflam.second, nflam.first);
+ }
+ return scaleFactor() * val;
+}
+
+double HiddenValleyAlpha::overestimateValue() const {
+ return scaleFactor() * _alphamin;
+}
+
+double HiddenValleyAlpha::ratio(const Energy2 scale,double) const {
+ pair<short,Energy> nflam;
+ Energy q = sqrt(scale);
+ double val(0.);
+ // special handling if the scale is less than Qmin
+ if (q < _qmin) {
+ nflam = getLamNfTwoLoop(_qmin);
+ double val0 = alphaS(_qmin, nflam.second, nflam.first);
+ switch (_asType) {
+ case 1:
+ // flat, zero; the default type with no NP effects.
+ val = 0.;
+ break;
+ case 2:
+ // flat, non-zero alpha_s = alpha_s(q2min).
+ val = val0;
+ break;
+ case 3:
+ // linear in q
+ val = val0*q/_qmin;
+ break;
+ case 4:
+ // quadratic in q
+ val = val0*sqr(q/_qmin);
+ break;
+ case 5:
+ // quadratic in q, starting off at asMaxNP, ending on as(qmin)
+ val = (val0 - _asMaxNP)*sqr(q/_qmin) + _asMaxNP;
+ break;
+ case 6:
+ // just asMaxNP and constant
+ val = _asMaxNP;
+ break;
+ }
+ } else {
+ // the 'ordinary' case
+ nflam = getLamNfTwoLoop(q);
+ val = alphaS(q, nflam.second, nflam.first);
+ }
+ // denominator
+ return val/_alphamin;
+}
+
+Energy HiddenValleyAlpha::computeLambda(Energy match,
+ double alpha,
+ unsigned int nflav) const {
+ Energy lamtest=200.0*MeV;
+ double xtest;
+ unsigned int ntry=0;
+ do {
+ ++ntry;
+ xtest=log(sqr(match/lamtest));
+ xtest+= (alpha-alphaS(match,lamtest,nflav))/derivativealphaS(match,lamtest,nflav);
+ lamtest=match/exp(0.5*xtest);
+ }
+ while(abs(alpha-alphaS(match,lamtest,nflav)) > _tolerance && ntry < _maxtry);
+ return lamtest;
+}
+
+pair<short, Energy> HiddenValleyAlpha::getLamNfTwoLoop(Energy q) const {
+ unsigned int ix=1;
+ for(;ix<_thresholds.size();ix++) {
+ if(q<_thresholds[ix]) break;
+ }
+ --ix;
+ return pair<short,Energy>(ix+_nf_light, _lambda[ix]);
+}
+
+void HiddenValleyAlpha::doinit() {
+ ShowerAlpha::doinit();
+ // get the model for parameters
+ tcHiddenValleyPtr model = dynamic_ptr_cast<tcHiddenValleyPtr>
+ (generator()->standardModel());
+ if(!model) throw InitException() << "Must be using the HiddenValleyModel"
+ << " in HiddenValleyAlpha::doinit()"
+ << Exception::runerror;
+ // get the colour factors
+ _ca = model->CA();
+ _cf = model->CF();
+ _tr = model->TR();
+ // get the thresholds
+ _thresholds.push_back(_qmin);
+ _nf_light = 0;
+ for(unsigned int ix=1;ix<=model->NF();++ix) {
+ Energy qmass = getParticleData(ParticleID::darkg+long(ix))->mass();
+ if (qmass > _qmin) _thresholds.push_back(qmass);
+ else _nf_light++;
+ }
+ _lambda.resize(_thresholds.size());
+
+
+ unsigned int nf = _nf_light;
+
+ // Set lambda below heavy quark thresholds to input value
+ _lambda[0]=_lambdain;
+ // convert lambda to the Monte Carlo scheme if needed
+ using Constants::pi;
+ if(_lambdaopt) _lambda[0] *=exp((0.5*_ca*(67./3.-sqr(pi))-_tr*nf*10./3.)/(11*_ca-2*nf))/sqrt(2.);
+
+ // Compute the threshold matching
+ for(int ix=0;ix<int(_thresholds.size())-1;++ix) {
+ _lambda[ix+1] = computeLambda(_thresholds[ix+1],alphaS(_thresholds[ix+1],
+ _lambda[ix],ix),ix+1);
+ }
+ // compute the maximum value of as
+ if ( _asType < 5 ) _alphamin = value(sqr(_qmin)+1.0e-8*sqr(MeV)); // approx as = 1
+ else _alphamin = max(_asMaxNP, value(sqr(_qmin)+1.0e-8*sqr(MeV)));
+ // check consistency lambda_3 < qmin
+ if(_lambda[0]>_qmin)
+ Throw<InitException>() << "The value of Qmin is less than Lambda in "
+ << _qmin/GeV << " < " << _lambda[0]/GeV
+ << " HiddenValleyAlpha::doinit " << Exception::abortnow;
+}
+
+double HiddenValleyAlpha::derivativealphaS(Energy q, Energy lam, int nf) const {
+ using Constants::pi;
+ double lx = log(sqr(q/lam));
+ // N.B. b_1 is divided by 2 due Hw++ convention
+ double b0 = 11./3.*_ca - 4./3.*_tr*nf;
+ double b1 = 17./3.*sqr(_ca) - nf*_tr*(10./3.*_ca+2.*_cf);
+ if(_nloop==1)
+ return -4.*pi/(b0*sqr(lx));
+ else if(_nloop==2)
+ return -4.*pi/(b0*sqr(lx))*(1.-2.*b1/sqr(b0)/lx*(1.-2.*log(lx)));
+ else
+ assert(false);
+ return 0.;
+}
+
+double HiddenValleyAlpha::alphaS(Energy q, Energy lam, int nf) const {
+ using Constants::pi;
+ double lx(log(sqr(q/lam)));
+ // N.B. b_1 is divided by 2 due Hw++ convention
+ double b0 = 11./3.*_ca - 4./3.*_tr*nf;
+ double b1 = 17./3.*sqr(_ca) - nf*_tr*(10./3.*_ca+2.*_cf);
+ // one loop
+ if(_nloop==1)
+ return 4.*pi/(b0*lx);
+ // two loop
+ else if(_nloop==2) {
+ return 4.*pi/(b0*lx)*(1.-2.*b1/sqr(b0)*log(lx)/lx);
+ }
+ else
+ assert(false);
+ return 0.;
+}
diff --git a/Shower/QTilde/Couplings/HiddenValleyAlpha.h b/Shower/QTilde/Couplings/HiddenValleyAlpha.h
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/Couplings/HiddenValleyAlpha.h
@@ -0,0 +1,333 @@
+// -*- C++ -*-
+//
+// HiddenValleyAlpha.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_HiddenValleyAlpha_H
+#define HERWIG_HiddenValleyAlpha_H
+//
+// This is the declaration of the HiddenValleyAlpha class.
+//
+
+#include "Herwig/Shower/ShowerAlpha.h"
+#include "ThePEG/Config/Constants.h"
+
+namespace Herwig {
+
+using namespace ThePEG;
+
+/** \ingroup Shower
+ *
+ * This concrete class provides the definition of the
+ * pure virtual function value() and overestimateValue() for the
+ * strong coupling.
+ *
+ * A number of different options for the running of the coupling
+ * and its initial definition are supported.
+ *
+ * @see \ref HiddenValleyAlphaInterfaces "The interfaces"
+ * defined for HiddenValleyAlpha.
+ */
+class HiddenValleyAlpha: public ShowerAlpha {
+
+public:
+
+ /**
+ * The default constructor.
+ */
+ HiddenValleyAlpha() : ShowerAlpha(),
+ _qmin(0.630882*GeV), _asType(1), _asMaxNP(1.0),
+ _nloop(2),_thresopt(false),
+ _lambdain(0.208364*GeV),
+ _tolerance(1e-10),
+ _maxtry(100),_alphamin(0.) {}
+
+public:
+
+ /**
+ * Methods to return the coupling
+ */
+ //@{
+ /**
+ * It returns the running coupling value evaluated at the input scale
+ * multiplied by the scale factor scaleFactor().
+ * @param scale The scale
+ * @return The coupling
+ */
+ virtual double value(const Energy2 scale) const;
+
+ /**
+ * It returns the running coupling value evaluated at the input scale
+ * multiplied by the scale factor scaleFactor().
+ */
+ virtual double overestimateValue() const;
+
+ /**
+ * Return the ratio of the coupling at the scale to the overestimated value
+ */
+ virtual double ratio(const Energy2 scale,double factor=1.) const;
+
+ /**
+ * Initialize this coupling.
+ */
+ virtual void initialize() { doinit(); }
+
+ //@}
+
+ /**
+ * Get the value of \f$\Lambda_{\rm QCd}\f$
+ * @param nf number of flavours
+ */
+ Energy lambdaQCD(unsigned int nf) {
+ if (nf <= 3) return _lambda[0];
+ else if (nf==4 || nf==5) return _lambda[nf-3];
+ else return _lambda[3];
+ }
+
+public:
+
+ /** @name Functions used by the persistent I/O system. */
+ //@{
+ /**
+ * Function used to write out object persistently.
+ * @param os the persistent output stream written to.
+ */
+ void persistentOutput(PersistentOStream & os) const;
+
+ /**
+ * Function used to read in object persistently.
+ * @param is the persistent input stream read from.
+ * @param version the version number of the object when written.
+ */
+ void persistentInput(PersistentIStream & is, int version);
+ //@}
+
+ /**
+ * The standard Init function used to initialize the interfaces.
+ * Called exactly once for each class by the class description system
+ * before the main function starts or
+ * when this class is dynamically loaded.
+ */
+ static void Init();
+
+protected:
+
+ /** @name Clone Methods. */
+ //@{
+ /**
+ * Make a simple clone of this object.
+ * @return a pointer to the new object.
+ */
+ virtual IBPtr clone() const;
+
+ /** Make a clone of this object, possibly modifying the cloned object
+ * to make it sane.
+ * @return a pointer to the new object.
+ */
+ virtual IBPtr fullclone() const;
+ //@}
+
+
+protected:
+
+ /** @name Standard Interfaced functions. */
+ //@{
+ /**
+ * Initialize this object after the setup phase before saving an
+ * EventGenerator to disk.
+ * @throws InitException if object could not be initialized properly.
+ */
+ virtual void doinit();
+ //@}
+
+private:
+
+ /**
+ * Member functions which calculate the coupling
+ */
+ //@{
+ /**
+ * The 1,2,3-loop parametrization of \f$\alpha_S\f$.
+ * @param q The scale
+ * @param lam \f$\Lambda_{\rm QCD}\f$
+ * @param nf The number of flavours
+ */
+ double alphaS(Energy q, Energy lam, int nf) const;
+
+ /**
+ * The derivative of \f$\alpha_S\f$ with respect to \f$\ln(Q^2/\Lambda^2)\f$
+ * @param q The scale
+ * @param lam \f$\Lambda_{\rm QCD}\f$
+ * @param nf The number of flavours
+ */
+ double derivativealphaS(Energy q, Energy lam, int nf) const;
+
+ /**
+ * Compute the value of \f$Lambda\f$ needed to get the input value of
+ * the strong coupling at the scale given for the given number of flavours
+ * using the Newton-Raphson method
+ * @param match The scale for the coupling
+ * @param alpha The input coupling
+ * @param nflav The number of flavours
+ */
+ Energy computeLambda(Energy match, double alpha, unsigned int nflav) const;
+
+ /**
+ * Return the value of \f$\Lambda\f$ and the number of flavours at the scale.
+ * @param q The scale
+ * @return The number of flavours at the scale and \f$\Lambda\f$.
+ */
+ pair<short, Energy> getLamNfTwoLoop(Energy q) const;
+ //@}
+
+private:
+
+ /**
+ * The static object used to initialize the description of this class.
+ * Indicates that this is a concrete class with persistent data.
+ */
+ static ClassDescription<HiddenValleyAlpha> initHiddenValleyAlpha;
+
+ /**
+ * The assignment operator is private and must never be called.
+ * In fact, it should not even be implemented.
+ */
+ HiddenValleyAlpha & operator=(const HiddenValleyAlpha &);
+
+private:
+
+ /**
+ * Minimum value of the scale
+ */
+ Energy _qmin;
+
+ /**
+ * Parameter controlling the behaviour of \f$\alpha_S\f$ in the
+ * non-perturbative region.
+ */
+ int _asType;
+
+ /**
+ * Another parameter, a possible (maximum) value of alpha in the
+ * non-perturbative region.
+ */
+ double _asMaxNP;
+
+ /**
+ * Thresholds for the different number of flavours
+ */
+ vector<Energy> _thresholds;
+
+ /**
+ * \f$\Lambda\f$ for the different number of flavours
+ */
+ vector<Energy> _lambda;
+
+ /**
+ * Option for the number of loops
+ */
+ unsigned int _nloop;
+
+ /**
+ * Option for the threshold masses
+ */
+ bool _thresopt;
+
+ /**
+ * Input value of Lambda
+ */
+ Energy _lambdain;
+
+ /**
+ * Whether to convert lambda from MSbar scheme
+ */
+ bool _lambdaopt;
+
+ /**
+ * Whether to use input value of alphas(MZ) or lambda
+ */
+ bool _inopt;
+
+ /**
+ * Tolerance for discontinuities at the thresholds
+ */
+ double _tolerance;
+
+ /**
+ * Maximum number of iterations for the Newton-Raphson method to converge
+ */
+ unsigned int _maxtry;
+
+ /**
+ * Number of light flavours (below qmin)
+ */
+ unsigned int _nf_light;
+
+ /**
+ * The minimum value of the coupling
+ */
+ double _alphamin;
+
+ /**
+ * Colour factors
+ */
+ //@{
+ /**
+ * \f$C_A\f$
+ */
+ double _ca;
+ /**
+ * \f$C_A\f$
+ */
+ double _cf;
+
+ /**
+ * \f$T_R\f$
+ */
+ double _tr;
+ //@}
+
+};
+
+}
+
+#include "ThePEG/Utilities/ClassTraits.h"
+
+namespace ThePEG {
+
+/** @cond TRAITSPECIALIZATIONS */
+
+/** This template specialization informs ThePEG about the
+ * base classes of HiddenValleyAlpha. */
+template <>
+struct BaseClassTrait<Herwig::HiddenValleyAlpha,1> {
+ /** Typedef of the first base class of HiddenValleyAlpha. */
+ typedef Herwig::ShowerAlpha NthBase;
+};
+
+/** This template specialization informs ThePEG about the name of
+ * the HiddenValleyAlpha class and the shared object where it is defined. */
+template <>
+struct ClassTraits<Herwig::HiddenValleyAlpha>
+ : public ClassTraitsBase<Herwig::HiddenValleyAlpha> {
+ /** Return a platform-independent class name */
+ static string className() { return "Herwig::HiddenValleyAlpha"; }
+ /**
+ * The name of a file containing the dynamic library where the class
+ * HiddenValleyAlpha is implemented. It may also include several,
+ * space-separated, libraries if the class HiddenValleyAlpha depends on
+ * other classes (base classes excepted). In this case the listed
+ * libraries will be dynamically linked in the order they are
+ * specified.
+ */
+ static string library() { return "HwShower.so Herwig.so"; }
+};
+
+/** @endcond */
+
+}
+
+#endif /* HERWIG_HiddenValleyAlpha_H */
diff --git a/Shower/QTilde/Dark/AdditionalGaugeParticleData.fh b/Shower/QTilde/Dark/AdditionalGaugeParticleData.fh
deleted file mode 100644
--- a/Shower/QTilde/Dark/AdditionalGaugeParticleData.fh
+++ /dev/null
@@ -1,19 +0,0 @@
-// -*- C++ -*-
-//
-// This is the forward declaration of the AdditionalGaugeParticleData class.
-//
-#ifndef HERWIG_AdditionalGaugeParticleData_FH
-#define HERWIG_AdditionalGaugeParticleData_FH
-
-#include "ThePEG/Config/ThePEG.h"
-
-namespace Herwig {
-
-class AdditionalGaugeParticleData;
-
-ThePEG_DECLARE_POINTERS(Herwig::AdditionalGaugeParticleData,
- AdditionalGaugeParticleDataPtr);
-
-}
-
-#endif
diff --git a/Shower/QTilde/Dark/HalfHalfOneDarkSplitFn.cc b/Shower/QTilde/Dark/HalfHalfOneDarkSplitFn.cc
deleted file mode 100644
--- a/Shower/QTilde/Dark/HalfHalfOneDarkSplitFn.cc
+++ /dev/null
@@ -1,171 +0,0 @@
-// -*- C++ -*-
-//
-// HalfHalfOneDarkSplitFn.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
-// Copyright (C) 2002-2019 The Herwig Collaboration
-//
-// Herwig is licenced under version 3 of the GPL, see COPYING for details.
-// Please respect the MCnet academic guidelines, see GUIDELINES for details.
-//
-//
-// This is the implementation of the non-inlined, non-templated member
-// functions of the HalfHalfOneDarkSplitFn class.
-//
-
-#include "HalfHalfOneDarkSplitFn.h"
-#include "HiddenValleyModel.h"
-#include "ThePEG/PDT/ParticleData.h"
-#include "ThePEG/Interface/ClassDocumentation.h"
-#include "ThePEG/Utilities/DescribeClass.h"
-#include <cassert>
-
-using namespace Herwig;
-
-DescribeNoPIOClass<HalfHalfOneDarkSplitFn,Herwig::Sudakov1to2FormFactor>
-describeHalfHalfOneDarkSplitFn ("Herwig::HalfHalfOneDarkSplitFn","HwShower.so");
-
-void HalfHalfOneDarkSplitFn::Init() {
- // documentation
- static ClassDocumentation<HalfHalfOneDarkSplitFn> documentation
- ("The HalfHalfOneDarkSplitFn class implements the dark coloured q -> qg splitting function");
-}
-
-bool HalfHalfOneDarkSplitFn::checkColours(const IdList & ids) const {
- if(colourStructure()==TripletTripletOctet) {
- if(ids[0]!=ids[1]) return false;
- if((ids[0]->iColour()==PDT::DarkColourFundamental||
- ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
- ids[2]->iColour()==PDT::DarkColourAdjoint) return true;
- return false;
- }
- else if(colourStructure()==OctetOctetOctet) {
- for(unsigned int ix=0;ix<3;++ix) {
- if(ids[ix]->iColour()!=PDT::DarkColourAdjoint) return false;
- }
- return true;
- }
- else if(colourStructure()==OctetTripletTriplet) {
- if(ids[0]->iColour()!=PDT::DarkColourAdjoint) return false;
- if(ids[1]->iColour()==PDT::DarkColourFundamental&&
- ids[2]->iColour()==PDT::DarkColourAntiFundamental)
- return true;
- if(ids[1]->iColour()==PDT::DarkColourAntiFundamental&&
- ids[2]->iColour()==PDT::DarkColourFundamental)
- return true;
- return false;
- }
- else if(colourStructure()==TripletOctetTriplet) {
- if(ids[0]!=ids[2]) return false;
- if((ids[0]->iColour()==PDT::DarkColourFundamental||
- ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
- ids[1]->iColour()==PDT::DarkColourAdjoint) return true;
- return false;
- }
- else {
- assert(false);
- }
- return false;
-}
-
-double HalfHalfOneDarkSplitFn::P(const double z, const Energy2 t,
- const IdList &ids, const bool mass, const RhoDMatrix & ) const {
- double val = (1. + sqr(z))/(1.-z);
- if(mass) {
- Energy m = ids[0]->mass();
- val -= 2.*sqr(m)/t;
- }
- return val;
-}
-
-double HalfHalfOneDarkSplitFn::overestimateP(const double z,
- const IdList & ) const {
- return 2./(1.-z);
-}
-
-double HalfHalfOneDarkSplitFn::ratioP(const double z, const Energy2 t,
- const IdList & ids, const bool mass, const RhoDMatrix & ) const {
- double val = 1. + sqr(z);
- if(mass) {
- Energy m = ids[0]->mass();
- val -= 2.*sqr(m)*(1.-z)/t;
- }
- return 0.5*val;
-}
-
-double HalfHalfOneDarkSplitFn::integOverP(const double z,
- const IdList & ,
- unsigned int PDFfactor) const {
- switch (PDFfactor) {
- case 0:
- return -2.*Math::log1m(z);
- case 1:
- return 2.*log(z/(1.-z));
- case 2:
- return 2./(1.-z);
- case 3:
- default:
- throw Exception() << "HalfHalfOneDarkSplitFn::integOverP() invalid PDFfactor = "
- << PDFfactor << Exception::runerror;
- }
-}
-
-double HalfHalfOneDarkSplitFn::invIntegOverP(const double r, const IdList &,
- unsigned int PDFfactor) const {
- switch (PDFfactor) {
- case 0:
- return 1. - exp(- 0.5*r);
- case 1:
- return 1./(1.-exp(-0.5*r));
- case 2:
- return 1.-2./r;
- case 3:
- default:
- throw Exception() << "HalfHalfOneDarkSplitFn::invIntegOverP() invalid PDFfactor = "
- << PDFfactor << Exception::runerror;
- }
-}
-
-bool HalfHalfOneDarkSplitFn::accept(const IdList &ids) const {
- // 3 particles and in and out fermion same
- if(ids.size()!=3 || ids[0]!=ids[1]) return false;
- if(ids[0]->iSpin()!=PDT::Spin1Half ||
- ids[2]->iSpin()!=PDT::Spin1) return false;
- return checkColours(ids);
-}
-
-vector<pair<int, Complex> >
-HalfHalfOneDarkSplitFn::generatePhiForward(const double, const Energy2, const IdList & ,
- const RhoDMatrix &) {
- // no dependence on the spin density matrix, dependence on off-diagonal terms cancels
- // and rest = splitting function for Tr(rho)=1 as required by defn
- return {{ {0, 1.} }};
-}
-
-vector<pair<int, Complex> >
-HalfHalfOneDarkSplitFn::generatePhiBackward(const double, const Energy2, const IdList & ,
- const RhoDMatrix &) {
- // no dependence on the spin density matrix, dependence on off-diagonal terms cancels
- // and rest = splitting function for Tr(rho)=1 as required by defn
- return {{ {0, 1.} }};
-}
-
-DecayMEPtr HalfHalfOneDarkSplitFn::matrixElement(const double z, const Energy2 t,
- const IdList & ids, const double phi,
- bool timeLike) {
- // calculate the kernal
- DecayMEPtr kernal(new_ptr(TwoBodyDecayMatrixElement(PDT::Spin1Half,PDT::Spin1Half,PDT::Spin1)));
- Energy m = !timeLike ? ZERO : ids[0]->mass();
- double mt = m/sqrt(t);
- double root = sqrt(1.-(1.-z)*sqr(m)/z/t);
- double romz = sqrt(1.-z);
- double rz = sqrt(z);
- Complex phase = exp(Complex(0.,1.)*phi);
- (*kernal)(0,0,0) = -root/romz*phase;
- (*kernal)(1,1,2) = -conj((*kernal)(0,0,0));
- (*kernal)(0,0,2) = root/romz*z/phase;
- (*kernal)(1,1,0) = -conj((*kernal)(0,0,2));
- (*kernal)(1,0,2) = mt*(1.-z)/rz;
- (*kernal)(0,1,0) = conj((*kernal)(1,0,2));
- (*kernal)(0,1,2) = 0.;
- (*kernal)(1,0,0) = 0.;
- return kernal;
-}
diff --git a/Shower/QTilde/Dark/HalfHalfOneDarkSplitFn.h b/Shower/QTilde/Dark/HalfHalfOneDarkSplitFn.h
deleted file mode 100644
--- a/Shower/QTilde/Dark/HalfHalfOneDarkSplitFn.h
+++ /dev/null
@@ -1,189 +0,0 @@
-// -*- C++ -*-
-//
-// HalfHalfOneDarkSplitFn.h is a part of Herwig - A multi-purpose Monte Carlo event generator
-// Copyright (C) 2002-2019 The Herwig Collaboration
-//
-// Herwig is licenced under version 3 of the GPL, see COPYING for details.
-// Please respect the MCnet academic guidelines, see GUIDELINES for details.
-//
-#ifndef HERWIG_HalfHalfOneDarkSplitFn_H
-#define HERWIG_HalfHalfOneDarkSplitFn_H
-//
-// This is the declaration of the HalfHalfOneDarkSplitFn class.
-//
-
-#include "Herwig/Shower/QTilde/SplittingFunctions/Sudakov1to2FormFactor.h"
-#include "Herwig/Decay/TwoBodyDecayMatrixElement.h"
-
-namespace Herwig {
-
-using namespace ThePEG;
-
-/**\ingroup Shower
- *
- * This class provides the concrete implementation of the exact leading-order
- * splitting function for \f$\frac12\to q\frac12 1\f$.
- *
- * In this case the splitting function is given by
- * \f[P(z,t) =C\left(\frac{1+z^2}{1-z}-2\frac{m^2_q}{t}\right),\f]
- * where \f$C\f$ is the corresponding colour factor.
- * Our choice for the overestimate is
- * \f[P_{\rm over}(z) = \frac{2C}{1-z},\f]
- * therefore the integral is
- * \f[\int P_{\rm over}(z) {\rm d}z = -2C\ln(1-z),\f]
- * and its inverse is
- * \f[1-\exp\left(\frac{r}{2C}\right).\f]
- *
- * @see \ref HalfHalfOneDarkSplitFnInterfaces "The interfaces"
- * defined for HalfHalfOneDarkSplitFn.
- */
-class HalfHalfOneDarkSplitFn: public Sudakov1to2FormFactor {
-
-public:
-
- /**
- * Method to check the colours are correct
- */
- virtual bool checkColours(const IdList & ids) const;
-
- /**
- * Concrete implementation of the method to determine whether this splitting
- * function can be used for a given set of particles.
- * @param ids The PDG codes for the particles in the splitting.
- */
- virtual bool accept(const IdList & ids) const;
-
- /**
- * Methods to return the splitting function.
- */
- //@{
- /**
- * The concrete implementation of the splitting function, \f$P(z,t)\f$.
- * @param z The energy fraction.
- * @param t The scale.
- * @param ids The PDG codes for the particles in the splitting.
- * @param mass Whether or not to include the mass dependent terms
- * @param rho The spin density matrix
- */
- virtual double P(const double z, const Energy2 t, const IdList & ids,
- const bool mass, const RhoDMatrix & rho) const;
-
- /**
- * The concrete implementation of the overestimate of the splitting function,
- * \f$P_{\rm over}\f$.
- * @param z The energy fraction.
- * @param ids The PDG codes for the particles in the splitting.
- */
- virtual double overestimateP(const double z, const IdList & ids) const;
-
- /**
- * The concrete implementation of the
- * the ratio of the splitting function to the overestimate, i.e.
- * \f$P(z,t)/P_{\rm over}(z)\f$.
- * @param z The energy fraction.
- * @param t The scale.
- * @param ids The PDG codes for the particles in the splitting.
- * @param mass Whether or not to include the mass dependent terms
- * @param rho The spin density matrix
- */
- virtual double ratioP(const double z, const Energy2 t, const IdList & ids,
- const bool mass, const RhoDMatrix & rho) const;
-
- /**
- * The concrete implementation of the indefinite integral of the
- * overestimated splitting function, \f$P_{\rm over}\f$.
- * @param z The energy fraction.
- * @param ids The PDG codes for the particles in the splitting.
- * @param PDFfactor Which additional factor to include for the PDF
- * 0 is no additional factor,
- * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
- */
- virtual double integOverP(const double z, const IdList & ids,
- unsigned int PDFfactor=0) const;
-
- /**
- * The concrete implementation of the inverse of the indefinite integral.
- * @param r Value of the splitting function to be inverted
- * @param ids The PDG codes for the particles in the splitting.
- * @param PDFfactor Which additional factor to include for the PDF
- * 0 is no additional factor,
- * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
- */
- virtual double invIntegOverP(const double r, const IdList & ids,
- unsigned int PDFfactor=0) const;
- //@}
-
- /**
- * Method to calculate the azimuthal angle
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- * @return The weight
- */
- virtual vector<pair<int,Complex> >
- generatePhiForward(const double z, const Energy2 t, const IdList & ids,
- const RhoDMatrix &);
-
- /**
- * Method to calculate the azimuthal angle for backward evolution
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- * @return The weight
- */
- virtual vector<pair<int,Complex> >
- generatePhiBackward(const double z, const Energy2 t, const IdList & ids,
- const RhoDMatrix &);
-
- /**
- * Calculate the matrix element for the splitting
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- */
- virtual DecayMEPtr matrixElement(const double z, const Energy2 t,
- const IdList & ids, const double phi, bool timeLike);
-
-public:
-
- /**
- * The standard Init function used to initialize the interfaces.
- * Called exactly once for each class by the class description system
- * before the main function starts or
- * when this class is dynamically loaded.
- */
- static void Init();
-
-protected:
-
- /** @name Clone Methods. */
- //@{
- /**
- * Make a simple clone of this object.
- * @return a pointer to the new object.
- */
- virtual IBPtr clone() const {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);}
- //@}
-
-private:
-
- /**
- * The assignment operator is private and must never be called.
- * In fact, it should not even be implemented.
- */
- HalfHalfOneDarkSplitFn & operator=(const HalfHalfOneDarkSplitFn &) = delete;
-
-};
-
-}
-
-#endif /* HERWIG_HalfHalfOneDarkSplitFn_H */
diff --git a/Shower/QTilde/Dark/HiddenValleyAlpha.cc b/Shower/QTilde/Dark/HiddenValleyAlpha.cc
deleted file mode 100644
--- a/Shower/QTilde/Dark/HiddenValleyAlpha.cc
+++ /dev/null
@@ -1,333 +0,0 @@
-// -*- C++ -*-
-//
-// HiddenValleyAlpha.cc 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.
-//
-//
-// This is the implementation of the non-inlined, non-templated member
-// functions of the HiddenValleyAlpha class.
-//
-#include "HiddenValleyAlpha.h"
-#include "HiddenValleyModel.h"
-#include "ThePEG/PDT/EnumParticles.h"
-#include "ThePEG/PDT/ParticleData.h"
-#include "ThePEG/Interface/ClassDocumentation.h"
-#include "ThePEG/Interface/Switch.h"
-#include "ThePEG/Interface/Parameter.h"
-#include "ThePEG/Interface/Command.h"
-#include "ThePEG/Persistency/PersistentOStream.h"
-#include "ThePEG/Persistency/PersistentIStream.h"
-#include "ThePEG/Utilities/Throw.h"
-
-using namespace Herwig;
-
-IBPtr HiddenValleyAlpha::clone() const {
- return new_ptr(*this);
-}
-
-IBPtr HiddenValleyAlpha::fullclone() const {
- return new_ptr(*this);
-}
-
-void HiddenValleyAlpha::persistentOutput(PersistentOStream & os) const {
- os << _asType << _asMaxNP << ounit(_qmin,GeV) << _nloop << _thresopt
- << ounit(_lambdain,GeV) << _tolerance << _maxtry << _alphamin
- << ounit(_thresholds,GeV) << ounit(_lambda,GeV) << _ca << _cf << _tr;
-}
-
-void HiddenValleyAlpha::persistentInput(PersistentIStream & is, int) {
- is >> _asType >> _asMaxNP >> iunit(_qmin,GeV) >> _nloop >> _thresopt
- >> iunit(_lambdain,GeV) >> _tolerance >> _maxtry >> _alphamin
- >> iunit(_thresholds,GeV) >> iunit(_lambda,GeV) >> _ca >> _cf >> _tr;
-}
-
-ClassDescription<HiddenValleyAlpha> HiddenValleyAlpha::initHiddenValleyAlpha;
-// Definition of the static class description member.
-
-void HiddenValleyAlpha::Init() {
-
- static ClassDocumentation<HiddenValleyAlpha> documentation
- ("This (concrete) class describes the QCD alpha running.");
-
- static Switch<HiddenValleyAlpha, int> intAsType
- ("NPAlphaS",
- "Behaviour of AlphaS in the NP region",
- &HiddenValleyAlpha::_asType, 1, false, false);
- static SwitchOption intAsTypeZero
- (intAsType, "Zero","zero below Q_min", 1);
- static SwitchOption intAsTypeConst
- (intAsType, "Const","const as(qmin) below Q_min", 2);
- static SwitchOption intAsTypeLin
- (intAsType, "Linear","growing linearly below Q_min", 3);
- static SwitchOption intAsTypeQuad
- (intAsType, "Quadratic","growing quadratically below Q_min", 4);
- static SwitchOption intAsTypeExx1
- (intAsType, "Exx1", "quadratic from AlphaMaxNP down to as(Q_min)", 5);
- static SwitchOption intAsTypeExx2
- (intAsType, "Exx2", "const = AlphaMaxNP below Q_min", 6);
-
- // default such that as(qmin) = 1 in the current parametrization.
- // min = Lambda3
- static Parameter<HiddenValleyAlpha,Energy> intQmin
- ("Qmin", "Q < Qmin is treated with NP parametrization as of (unit [GeV])",
- &HiddenValleyAlpha::_qmin, GeV, 0.630882*GeV, 0.330445*GeV,
- 100.0*GeV,false,false,false);
-
- static Parameter<HiddenValleyAlpha,double> interfaceAlphaMaxNP
- ("AlphaMaxNP",
- "Max value of alpha in NP region, only relevant if NPAlphaS = 5,6",
- &HiddenValleyAlpha::_asMaxNP, 1.0, 0., 100.0,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyAlpha,unsigned int> interfaceNumberOfLoops
- ("NumberOfLoops",
- "The number of loops to use in the alpha_S calculation",
- &HiddenValleyAlpha::_nloop, 2, 1, 2,
- false, false, Interface::limited);
-
- static Switch<HiddenValleyAlpha,bool> interfaceLambdaOption
- ("LambdaOption",
- "Option for the calculation of the Lambda used in the simulation from the input"
- " Lambda_MSbar",
- &HiddenValleyAlpha::_lambdaopt, false, false, false);
- static SwitchOption interfaceLambdaOptionfalse
- (interfaceLambdaOption,
- "Same",
- "Use the same value",
- false);
- static SwitchOption interfaceLambdaOptionConvert
- (interfaceLambdaOption,
- "Convert",
- "Use the conversion to the Herwig scheme from NPB349, 635",
- true);
-
- static Parameter<HiddenValleyAlpha,Energy> interfaceLambdaDark
- ("LambdaDark",
- "Input value of Lambda_MSBar",
- &HiddenValleyAlpha::_lambdain, GeV, 0.208364*GeV, 100.0*MeV, 100.0*GeV,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyAlpha,double> interfaceTolerance
- ("Tolerance",
- "The tolerance for discontinuities in alphaS at thresholds.",
- &HiddenValleyAlpha::_tolerance, 1e-10, 1e-20, 1e-4,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyAlpha,unsigned int> interfaceMaximumIterations
- ("MaximumIterations",
- "The maximum number of iterations for the Newton-Raphson method to converge.",
- &HiddenValleyAlpha::_maxtry, 100, 10, 1000,
- false, false, Interface::limited);
-
- static Switch<HiddenValleyAlpha,bool> interfaceThresholdOption
- ("ThresholdOption",
- "Whether to use the consistuent or normal masses for the thresholds",
- &HiddenValleyAlpha::_thresopt, true, false, false);
- static SwitchOption interfaceThresholdOptionCurrent
- (interfaceThresholdOption,
- "Current",
- "Use the current masses",
- true);
- static SwitchOption interfaceThresholdOptionConstituent
- (interfaceThresholdOption,
- "Constituent",
- "Use the constitent masses.",
- false);
-}
-
-double HiddenValleyAlpha::value(const Energy2 scale) const {
- pair<short,Energy> nflam;
- Energy q = sqrt(scale);
- double val(0.);
- // special handling if the scale is less than Qmin
- if (q < _qmin) {
- nflam = getLamNfTwoLoop(_qmin);
- double val0 = alphaS(_qmin, nflam.second, nflam.first);
- switch (_asType) {
- case 1:
- // flat, zero; the default type with no NP effects.
- val = 0.;
- break;
- case 2:
- // flat, non-zero alpha_s = alpha_s(q2min).
- val = val0;
- break;
- case 3:
- // linear in q
- val = val0*q/_qmin;
- break;
- case 4:
- // quadratic in q
- val = val0*sqr(q/_qmin);
- break;
- case 5:
- // quadratic in q, starting off at asMaxNP, ending on as(qmin)
- val = (val0 - _asMaxNP)*sqr(q/_qmin) + _asMaxNP;
- break;
- case 6:
- // just asMaxNP and constant
- val = _asMaxNP;
- break;
- }
- }
- else {
- // the 'ordinary' case
- nflam = getLamNfTwoLoop(q);
- val = alphaS(q, nflam.second, nflam.first);
- }
- return scaleFactor() * val;
-}
-
-double HiddenValleyAlpha::overestimateValue() const {
- return scaleFactor() * _alphamin;
-}
-
-double HiddenValleyAlpha::ratio(const Energy2 scale,double) const {
- pair<short,Energy> nflam;
- Energy q = sqrt(scale);
- double val(0.);
- // special handling if the scale is less than Qmin
- if (q < _qmin) {
- nflam = getLamNfTwoLoop(_qmin);
- double val0 = alphaS(_qmin, nflam.second, nflam.first);
- switch (_asType) {
- case 1:
- // flat, zero; the default type with no NP effects.
- val = 0.;
- break;
- case 2:
- // flat, non-zero alpha_s = alpha_s(q2min).
- val = val0;
- break;
- case 3:
- // linear in q
- val = val0*q/_qmin;
- break;
- case 4:
- // quadratic in q
- val = val0*sqr(q/_qmin);
- break;
- case 5:
- // quadratic in q, starting off at asMaxNP, ending on as(qmin)
- val = (val0 - _asMaxNP)*sqr(q/_qmin) + _asMaxNP;
- break;
- case 6:
- // just asMaxNP and constant
- val = _asMaxNP;
- break;
- }
- } else {
- // the 'ordinary' case
- nflam = getLamNfTwoLoop(q);
- val = alphaS(q, nflam.second, nflam.first);
- }
- // denominator
- return val/_alphamin;
-}
-
-Energy HiddenValleyAlpha::computeLambda(Energy match,
- double alpha,
- unsigned int nflav) const {
- Energy lamtest=200.0*MeV;
- double xtest;
- unsigned int ntry=0;
- do {
- ++ntry;
- xtest=log(sqr(match/lamtest));
- xtest+= (alpha-alphaS(match,lamtest,nflav))/derivativealphaS(match,lamtest,nflav);
- lamtest=match/exp(0.5*xtest);
- }
- while(abs(alpha-alphaS(match,lamtest,nflav)) > _tolerance && ntry < _maxtry);
- return lamtest;
-}
-
-pair<short, Energy> HiddenValleyAlpha::getLamNfTwoLoop(Energy q) const {
- unsigned int ix=1;
- for(;ix<_thresholds.size();ix++) {
- if(q<_thresholds[ix]) break;
- }
- --ix;
- return pair<short,Energy>(ix+_nf_light, _lambda[ix]);
-}
-
-void HiddenValleyAlpha::doinit() {
- ShowerAlpha::doinit();
- // get the model for parameters
- tcHiddenValleyPtr model = dynamic_ptr_cast<tcHiddenValleyPtr>
- (generator()->standardModel());
- if(!model) throw InitException() << "Must be using the HiddenValleyModel"
- << " in HiddenValleyAlpha::doinit()"
- << Exception::runerror;
- // get the colour factors
- _ca = model->CA();
- _cf = model->CF();
- _tr = model->TR();
- // get the thresholds
- _thresholds.push_back(_qmin);
- _nf_light = 0;
- for(unsigned int ix=1;ix<=model->NF();++ix) {
- Energy qmass = getParticleData(ParticleID::darkg+long(ix))->mass();
- if (qmass > _qmin) _thresholds.push_back(qmass);
- else _nf_light++;
- }
- _lambda.resize(_thresholds.size());
-
-
- unsigned int nf = _nf_light;
-
- // Set lambda below heavy quark thresholds to input value
- _lambda[0]=_lambdain;
- // convert lambda to the Monte Carlo scheme if needed
- using Constants::pi;
- if(_lambdaopt) _lambda[0] *=exp((0.5*_ca*(67./3.-sqr(pi))-_tr*nf*10./3.)/(11*_ca-2*nf))/sqrt(2.);
-
- // Compute the threshold matching
- for(int ix=0;ix<int(_thresholds.size())-1;++ix) {
- _lambda[ix+1] = computeLambda(_thresholds[ix+1],alphaS(_thresholds[ix+1],
- _lambda[ix],ix),ix+1);
- }
- // compute the maximum value of as
- if ( _asType < 5 ) _alphamin = value(sqr(_qmin)+1.0e-8*sqr(MeV)); // approx as = 1
- else _alphamin = max(_asMaxNP, value(sqr(_qmin)+1.0e-8*sqr(MeV)));
- // check consistency lambda_3 < qmin
- if(_lambda[0]>_qmin)
- Throw<InitException>() << "The value of Qmin is less than Lambda in "
- << _qmin/GeV << " < " << _lambda[0]/GeV
- << " HiddenValleyAlpha::doinit " << Exception::abortnow;
-}
-
-double HiddenValleyAlpha::derivativealphaS(Energy q, Energy lam, int nf) const {
- using Constants::pi;
- double lx = log(sqr(q/lam));
- // N.B. b_1 is divided by 2 due Hw++ convention
- double b0 = 11./3.*_ca - 4./3.*_tr*nf;
- double b1 = 17./3.*sqr(_ca) - nf*_tr*(10./3.*_ca+2.*_cf);
- if(_nloop==1)
- return -4.*pi/(b0*sqr(lx));
- else if(_nloop==2)
- return -4.*pi/(b0*sqr(lx))*(1.-2.*b1/sqr(b0)/lx*(1.-2.*log(lx)));
- else
- assert(false);
- return 0.;
-}
-
-double HiddenValleyAlpha::alphaS(Energy q, Energy lam, int nf) const {
- using Constants::pi;
- double lx(log(sqr(q/lam)));
- // N.B. b_1 is divided by 2 due Hw++ convention
- double b0 = 11./3.*_ca - 4./3.*_tr*nf;
- double b1 = 17./3.*sqr(_ca) - nf*_tr*(10./3.*_ca+2.*_cf);
- // one loop
- if(_nloop==1)
- return 4.*pi/(b0*lx);
- // two loop
- else if(_nloop==2) {
- return 4.*pi/(b0*lx)*(1.-2.*b1/sqr(b0)*log(lx)/lx);
- }
- else
- assert(false);
- return 0.;
-}
diff --git a/Shower/QTilde/Dark/HiddenValleyAlpha.h b/Shower/QTilde/Dark/HiddenValleyAlpha.h
deleted file mode 100644
--- a/Shower/QTilde/Dark/HiddenValleyAlpha.h
+++ /dev/null
@@ -1,333 +0,0 @@
-// -*- C++ -*-
-//
-// HiddenValleyAlpha.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_HiddenValleyAlpha_H
-#define HERWIG_HiddenValleyAlpha_H
-//
-// This is the declaration of the HiddenValleyAlpha class.
-//
-
-#include "Herwig/Shower/ShowerAlpha.h"
-#include "ThePEG/Config/Constants.h"
-
-namespace Herwig {
-
-using namespace ThePEG;
-
-/** \ingroup Shower
- *
- * This concrete class provides the definition of the
- * pure virtual function value() and overestimateValue() for the
- * strong coupling.
- *
- * A number of different options for the running of the coupling
- * and its initial definition are supported.
- *
- * @see \ref HiddenValleyAlphaInterfaces "The interfaces"
- * defined for HiddenValleyAlpha.
- */
-class HiddenValleyAlpha: public ShowerAlpha {
-
-public:
-
- /**
- * The default constructor.
- */
- HiddenValleyAlpha() : ShowerAlpha(),
- _qmin(0.630882*GeV), _asType(1), _asMaxNP(1.0),
- _nloop(2),_thresopt(false),
- _lambdain(0.208364*GeV),
- _tolerance(1e-10),
- _maxtry(100),_alphamin(0.) {}
-
-public:
-
- /**
- * Methods to return the coupling
- */
- //@{
- /**
- * It returns the running coupling value evaluated at the input scale
- * multiplied by the scale factor scaleFactor().
- * @param scale The scale
- * @return The coupling
- */
- virtual double value(const Energy2 scale) const;
-
- /**
- * It returns the running coupling value evaluated at the input scale
- * multiplied by the scale factor scaleFactor().
- */
- virtual double overestimateValue() const;
-
- /**
- * Return the ratio of the coupling at the scale to the overestimated value
- */
- virtual double ratio(const Energy2 scale,double factor=1.) const;
-
- /**
- * Initialize this coupling.
- */
- virtual void initialize() { doinit(); }
-
- //@}
-
- /**
- * Get the value of \f$\Lambda_{\rm QCd}\f$
- * @param nf number of flavours
- */
- Energy lambdaQCD(unsigned int nf) {
- if (nf <= 3) return _lambda[0];
- else if (nf==4 || nf==5) return _lambda[nf-3];
- else return _lambda[3];
- }
-
-public:
-
- /** @name Functions used by the persistent I/O system. */
- //@{
- /**
- * Function used to write out object persistently.
- * @param os the persistent output stream written to.
- */
- void persistentOutput(PersistentOStream & os) const;
-
- /**
- * Function used to read in object persistently.
- * @param is the persistent input stream read from.
- * @param version the version number of the object when written.
- */
- void persistentInput(PersistentIStream & is, int version);
- //@}
-
- /**
- * The standard Init function used to initialize the interfaces.
- * Called exactly once for each class by the class description system
- * before the main function starts or
- * when this class is dynamically loaded.
- */
- static void Init();
-
-protected:
-
- /** @name Clone Methods. */
- //@{
- /**
- * Make a simple clone of this object.
- * @return a pointer to the new object.
- */
- virtual IBPtr clone() const;
-
- /** Make a clone of this object, possibly modifying the cloned object
- * to make it sane.
- * @return a pointer to the new object.
- */
- virtual IBPtr fullclone() const;
- //@}
-
-
-protected:
-
- /** @name Standard Interfaced functions. */
- //@{
- /**
- * Initialize this object after the setup phase before saving an
- * EventGenerator to disk.
- * @throws InitException if object could not be initialized properly.
- */
- virtual void doinit();
- //@}
-
-private:
-
- /**
- * Member functions which calculate the coupling
- */
- //@{
- /**
- * The 1,2,3-loop parametrization of \f$\alpha_S\f$.
- * @param q The scale
- * @param lam \f$\Lambda_{\rm QCD}\f$
- * @param nf The number of flavours
- */
- double alphaS(Energy q, Energy lam, int nf) const;
-
- /**
- * The derivative of \f$\alpha_S\f$ with respect to \f$\ln(Q^2/\Lambda^2)\f$
- * @param q The scale
- * @param lam \f$\Lambda_{\rm QCD}\f$
- * @param nf The number of flavours
- */
- double derivativealphaS(Energy q, Energy lam, int nf) const;
-
- /**
- * Compute the value of \f$Lambda\f$ needed to get the input value of
- * the strong coupling at the scale given for the given number of flavours
- * using the Newton-Raphson method
- * @param match The scale for the coupling
- * @param alpha The input coupling
- * @param nflav The number of flavours
- */
- Energy computeLambda(Energy match, double alpha, unsigned int nflav) const;
-
- /**
- * Return the value of \f$\Lambda\f$ and the number of flavours at the scale.
- * @param q The scale
- * @return The number of flavours at the scale and \f$\Lambda\f$.
- */
- pair<short, Energy> getLamNfTwoLoop(Energy q) const;
- //@}
-
-private:
-
- /**
- * The static object used to initialize the description of this class.
- * Indicates that this is a concrete class with persistent data.
- */
- static ClassDescription<HiddenValleyAlpha> initHiddenValleyAlpha;
-
- /**
- * The assignment operator is private and must never be called.
- * In fact, it should not even be implemented.
- */
- HiddenValleyAlpha & operator=(const HiddenValleyAlpha &);
-
-private:
-
- /**
- * Minimum value of the scale
- */
- Energy _qmin;
-
- /**
- * Parameter controlling the behaviour of \f$\alpha_S\f$ in the
- * non-perturbative region.
- */
- int _asType;
-
- /**
- * Another parameter, a possible (maximum) value of alpha in the
- * non-perturbative region.
- */
- double _asMaxNP;
-
- /**
- * Thresholds for the different number of flavours
- */
- vector<Energy> _thresholds;
-
- /**
- * \f$\Lambda\f$ for the different number of flavours
- */
- vector<Energy> _lambda;
-
- /**
- * Option for the number of loops
- */
- unsigned int _nloop;
-
- /**
- * Option for the threshold masses
- */
- bool _thresopt;
-
- /**
- * Input value of Lambda
- */
- Energy _lambdain;
-
- /**
- * Whether to convert lambda from MSbar scheme
- */
- bool _lambdaopt;
-
- /**
- * Whether to use input value of alphas(MZ) or lambda
- */
- bool _inopt;
-
- /**
- * Tolerance for discontinuities at the thresholds
- */
- double _tolerance;
-
- /**
- * Maximum number of iterations for the Newton-Raphson method to converge
- */
- unsigned int _maxtry;
-
- /**
- * Number of light flavours (below qmin)
- */
- unsigned int _nf_light;
-
- /**
- * The minimum value of the coupling
- */
- double _alphamin;
-
- /**
- * Colour factors
- */
- //@{
- /**
- * \f$C_A\f$
- */
- double _ca;
- /**
- * \f$C_A\f$
- */
- double _cf;
-
- /**
- * \f$T_R\f$
- */
- double _tr;
- //@}
-
-};
-
-}
-
-#include "ThePEG/Utilities/ClassTraits.h"
-
-namespace ThePEG {
-
-/** @cond TRAITSPECIALIZATIONS */
-
-/** This template specialization informs ThePEG about the
- * base classes of HiddenValleyAlpha. */
-template <>
-struct BaseClassTrait<Herwig::HiddenValleyAlpha,1> {
- /** Typedef of the first base class of HiddenValleyAlpha. */
- typedef Herwig::ShowerAlpha NthBase;
-};
-
-/** This template specialization informs ThePEG about the name of
- * the HiddenValleyAlpha class and the shared object where it is defined. */
-template <>
-struct ClassTraits<Herwig::HiddenValleyAlpha>
- : public ClassTraitsBase<Herwig::HiddenValleyAlpha> {
- /** Return a platform-independent class name */
- static string className() { return "Herwig::HiddenValleyAlpha"; }
- /**
- * The name of a file containing the dynamic library where the class
- * HiddenValleyAlpha is implemented. It may also include several,
- * space-separated, libraries if the class HiddenValleyAlpha depends on
- * other classes (base classes excepted). In this case the listed
- * libraries will be dynamically linked in the order they are
- * specified.
- */
- static string library() { return "HwShower.so HwHiddenValley.so"; }
-};
-
-/** @endcond */
-
-}
-
-#endif /* HERWIG_HiddenValleyAlpha_H */
diff --git a/Shower/QTilde/Dark/HiddenValleyFFZPrimeVertex.cc b/Shower/QTilde/Dark/HiddenValleyFFZPrimeVertex.cc
deleted file mode 100644
--- a/Shower/QTilde/Dark/HiddenValleyFFZPrimeVertex.cc
+++ /dev/null
@@ -1,93 +0,0 @@
-// -*- C++ -*-
-//
-// This is the implementation of the non-inlined, non-templated member
-// functions of the HiddenValleyFFZPrimeVertex class.
-//
-
-#include "HiddenValleyFFZPrimeVertex.h"
-#include "ThePEG/Interface/ClassDocumentation.h"
-#include "ThePEG/Persistency/PersistentOStream.h"
-#include "ThePEG/Persistency/PersistentIStream.h"
-#include "HiddenValleyModel.h"
-
-using namespace Herwig;
-
-void HiddenValleyFFZPrimeVertex::persistentOutput(PersistentOStream & os) const {
- os << _gl << _gr << _gql << _gqr << _gPrime;
-}
-
-void HiddenValleyFFZPrimeVertex::persistentInput(PersistentIStream & is, int) {
- is >> _gl >> _gr >> _gql >> _gqr >> _gPrime;
-}
-
-ClassDescription<HiddenValleyFFZPrimeVertex> HiddenValleyFFZPrimeVertex::initHiddenValleyFFZPrimeVertex;
-// Definition of the static class description member.
-
-void HiddenValleyFFZPrimeVertex::Init() {
-
- static ClassDocumentation<HiddenValleyFFZPrimeVertex> documentation
- ("There is no documentation for the HiddenValleyFFZPrimeVertex class");
-
-}
-
-void HiddenValleyFFZPrimeVertex::setCoupling(Energy2, tcPDPtr a,tcPDPtr,tcPDPtr) {
- norm(_gPrime);
- // the left and right couplings
- int iferm=abs(a->id());
- if((iferm>=1 && iferm<=6)||(iferm>=11 &&iferm<=16)) {
- left (_gl[iferm]);
- right(_gr[iferm]);
- return;
- }
- iferm -= ParticleID::darkg;
- if(iferm<=int(_gql.size())) {
- left (_gql[iferm]);
- right(_gqr[iferm]);
- return;
- }
- assert(false);
-}
-
-HiddenValleyFFZPrimeVertex::HiddenValleyFFZPrimeVertex() : _gl(17,0.0), _gr(17,0.0) {
- orderInGem(1);
- orderInGs(0);
- // PDG codes for the particles
- // the quarks
- for(int ix=1;ix<7;++ix)
- addToList(-ix,ix,32);
- // the leptons
- for(int ix=11;ix<17;++ix)
- addToList(-ix,ix,32);
-}
-
-void HiddenValleyFFZPrimeVertex::doinit() {
- tcHiddenValleyPtr model =
- dynamic_ptr_cast<tcHiddenValleyPtr>(generator()->standardModel());
- if(!model)
- throw InitException() << "Must be using the HiddenValleyModel in "
- << "HiddenValleyFFZPrimeVertex::doinit()"
- << Exception::runerror;
- _gPrime = model->gPrime();
- // SM couplings
- _gl.resize(17,0.);
- _gr.resize(17,0.);
- for(int ix=1;ix<4;++ix) {
- _gl[2*ix-1] = model->qL();
- _gl[2*ix ] = model->qL();
- _gl[2*ix+9 ] = model->lL();
- _gl[2*ix+10] = model->lL();
- _gr[2*ix-1] = model->dR();
- _gr[2*ix ] = model->uR();
- _gr[2*ix+9 ] = model->lR();
- _gr[2*ix+10] = 0. ;
- }
- FFVVertex::doinit();
- _gql.resize(model->NF()+1,0.);
- _gqr.resize(model->NF()+1,0.);
- for(int ix=0;ix<int(model->NF());++ix) {
- int id = ParticleID::darkg+1+ix;
- addToList(-id,id,32);
- _gql[ix+1] = model->qCharge()[ix];
- _gqr[ix+1] = model->qCharge()[ix]-2.;
- }
-}
diff --git a/Shower/QTilde/Dark/HiddenValleyFFZPrimeVertex.h b/Shower/QTilde/Dark/HiddenValleyFFZPrimeVertex.h
deleted file mode 100644
--- a/Shower/QTilde/Dark/HiddenValleyFFZPrimeVertex.h
+++ /dev/null
@@ -1,173 +0,0 @@
-// -*- C++ -*-
-#ifndef RADIATIVEZPRIME_HiddenValleyFFZPrimeVertex_H
-#define RADIATIVEZPRIME_HiddenValleyFFZPrimeVertex_H
-//
-// This is the declaration of the HiddenValleyFFZPrimeVertex class.
-//
-
-#include "ThePEG/Helicity/Vertex/Vector/FFVVertex.h"
-
-namespace Herwig {
-
-using namespace ThePEG;
-
-/**
- * Here is the documentation of the HiddenValleyFFZPrimeVertex class.
- *
- * @see \ref HiddenValleyFFZPrimeVertexInterfaces "The interfaces"
- * defined for HiddenValleyFFZPrimeVertex.
- */
-class HiddenValleyFFZPrimeVertex: public Helicity::FFVVertex {
-
-public:
-
- /**
- * The default constructor.
- */
- inline HiddenValleyFFZPrimeVertex();
-
- /**
- * Calculate the couplings.
- * @param q2 The scale \f$q^2\f$ for the coupling at the vertex.
- * @param part1 The ParticleData pointer for the first particle.
- * @param part2 The ParticleData pointer for the second particle.
- * @param part3 The ParticleData pointer for the third particle.
- */
- virtual void setCoupling(Energy2 q2,tcPDPtr part1,tcPDPtr part2,tcPDPtr part3);
-
-public:
-
- /** @name Functions used by the persistent I/O system. */
- //@{
- /**
- * Function used to write out object persistently.
- * @param os the persistent output stream written to.
- */
- void persistentOutput(PersistentOStream & os) const;
-
- /**
- * Function used to read in object persistently.
- * @param is the persistent input stream read from.
- * @param version the version number of the object when written.
- */
- void persistentInput(PersistentIStream & is, int version);
- //@}
-
- /**
- * The standard Init function used to initialize the interfaces.
- * Called exactly once for each class by the class description system
- * before the main function starts or
- * when this class is dynamically loaded.
- */
- static void Init();
-
-protected:
-
- /** @name Clone Methods. */
- //@{
- /**
- * Make a simple clone of this object.
- * @return a pointer to the new object.
- */
- inline 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.
- */
- inline virtual IBPtr fullclone() const {return new_ptr(*this);}
- //@}
-
-protected:
-
- /**
- * Initialize this object after the setup phase before saving and
- * EventGenerator to disk.
- * @throws InitException if object could not be initialized properly.
- */
- virtual void doinit();
-
-private:
-
- /**
- * The static object used to initialize the description of this class.
- * Indicates that this is a concrete class with persistent data.
- */
- static ClassDescription<HiddenValleyFFZPrimeVertex> initHiddenValleyFFZPrimeVertex;
-
- /**
- * The assignment operator is private and must never be called.
- * In fact, it should not even be implemented.
- */
- HiddenValleyFFZPrimeVertex & operator=(const HiddenValleyFFZPrimeVertex &);
-
-
- /**
- * Storage of the couplings.
- */
- //@{
- /**
- * The left couplings of the Standard Model fermions.
- */
- vector<double> _gl;
-
- /**
- * The right couplings of the Standard Model fermions.
- */
- vector<double> _gr;
- /**
- * The left couplings of the new fermions.
- */
- vector<double> _gql;
-
- /**
- * The right couplings of the new fermions.
- */
- vector<double> _gqr;
-
- /**
- * Coupling
- */
- double _gPrime;
- //@}
-
-};
-
-}
-
-#include "ThePEG/Utilities/ClassTraits.h"
-
-namespace ThePEG {
-
-/** @cond TRAITSPECIALIZATIONS */
-
-/** This template specialization informs ThePEG about the
- * base classes of HiddenValleyFFZPrimeVertex. */
-template <>
-struct BaseClassTrait<Herwig::HiddenValleyFFZPrimeVertex,1> {
- /** Typedef of the first base class of HiddenValleyFFZPrimeVertex. */
- typedef Helicity::FFVVertex NthBase;
-};
-
-/** This template specialization informs ThePEG about the name of
- * the HiddenValleyFFZPrimeVertex class and the shared object where it is defined. */
-template <>
-struct ClassTraits<Herwig::HiddenValleyFFZPrimeVertex>
- : public ClassTraitsBase<Herwig::HiddenValleyFFZPrimeVertex> {
- /** Return a platform-independent class name */
- static string className() { return "Herwig::HiddenValleyFFZPrimeVertex"; }
- /**
- * The name of a file containing the dynamic library where the class
- * HiddenValleyFFZPrimeVertex is implemented. It may also include several, space-separated,
- * libraries if the class HiddenValleyFFZPrimeVertex depends on other classes (base classes
- * excepted). In this case the listed libraries will be dynamically
- * linked in the order they are specified.
- */
- static string library() { return "HwShower.so HwHiddenValley.so"; }
-};
-
-/** @endcond */
-
-}
-
-#endif /* RADIATIVEZPRIME_HiddenValleyFFZPrimeVertex_H */
diff --git a/Shower/QTilde/Dark/HiddenValleyModel.cc b/Shower/QTilde/Dark/HiddenValleyModel.cc
deleted file mode 100644
--- a/Shower/QTilde/Dark/HiddenValleyModel.cc
+++ /dev/null
@@ -1,130 +0,0 @@
-// -*- C++ -*-
-//
-// This is the implementation of the non-inlined, non-templated member
-// functions of the HiddenValleyModel class.
-//
-
-#include "HiddenValleyModel.h"
-#include "ThePEG/Interface/ClassDocumentation.h"
-#include "ThePEG/Interface/Reference.h"
-#include "ThePEG/Interface/Parameter.h"
-#include "ThePEG/Interface/ParVector.h"
-#include "ThePEG/Interface/Switch.h"
-#include "ThePEG/Utilities/EnumIO.h"
-#include "ThePEG/Persistency/PersistentOStream.h"
-#include "ThePEG/Persistency/PersistentIStream.h"
-
-using namespace Herwig;
-
-HiddenValleyModel::HiddenValleyModel() : groupType_(SU), Nc_(3), Nf_(1),
- qL_(-0.2), uR_(-0.2), dR_(0.6),
- lL_(0.6), lR_(-0.2), gPrime_(1./7.),
- qCharge_(1,-0.2)
-{}
-
-IBPtr HiddenValleyModel::clone() const {
- return new_ptr(*this);
-}
-
-IBPtr HiddenValleyModel::fullclone() const {
- return new_ptr(*this);
-}
-
-void HiddenValleyModel::persistentOutput(PersistentOStream & os) const {
- os << oenum(groupType_) << Nc_ << Nf_ << FFZPVertex_
- << qL_ << uR_ << dR_ << lL_ << lR_ << qCharge_ << gPrime_;
-}
-
-void HiddenValleyModel::persistentInput(PersistentIStream & is, int) {
- is >> ienum(groupType_) >> Nc_ >> Nf_ >> FFZPVertex_
- >> qL_ >> uR_ >> dR_ >> lL_ >> lR_ >> qCharge_ >> gPrime_;
-}
-
-ClassDescription<HiddenValleyModel> HiddenValleyModel::initHiddenValleyModel;
-// Definition of the static class description member.
-
-void HiddenValleyModel::Init() {
-
- static ClassDocumentation<HiddenValleyModel> documentation
- ("The HiddenValleyModel class is the main class for the Hidden Valley Model.");
-
- static Switch<HiddenValleyModel,GroupType> interfaceGroupType
- ("GroupType",
- "Type of the new unbroken group",
- &HiddenValleyModel::groupType_, SU, false, false);
- static SwitchOption interfaceGroupTypeSU
- (interfaceGroupType,
- "SU",
- "Use an SU(N) group",
- SU);
-
- static Parameter<HiddenValleyModel,unsigned int> interfaceGroupOrder
- ("GroupOrder",
- "The value of the number of 'colours' for the new symmetry group",
- &HiddenValleyModel::Nc_, 3, 1, 1000,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyModel,unsigned int> interfaceNumberOfFermions
- ("NumberOfFermions",
- "The number of fermions charged under the new group",
- &HiddenValleyModel::Nf_, 1, 1, 100,
- false, false, Interface::limited);
-
- static Reference<HiddenValleyModel,AbstractFFVVertex> interfaceVertexFFZPrime
- ("Vertex/FFZPrime",
- "The vertex coupling the Zprime to the fermions",
- &HiddenValleyModel::FFZPVertex_, false, false, true, false, false);
-
- static ParVector<HiddenValleyModel,double> interfaceQuirkCharges
- ("QuirkCharges",
- "The charges under the new U(1) of the new quirks",
- &HiddenValleyModel::qCharge_, -1, -0.2, -10., 10.0,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyModel,double> interfaceZPrimeQL
- ("ZPrime/QL",
- "The charge of the left-handed quarks under the new U(1)",
- &HiddenValleyModel::qL_, -0.2, 10.0, 10.0,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyModel,double> interfaceZPrimeUR
- ("ZPrime/UR",
- "The charge of the right-handed up quarks under the new U(1)",
- &HiddenValleyModel::uR_, -0.2, 10.0, 10.0,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyModel,double> interfaceZPrimeDR
- ("ZPrime/DR",
- "The charge of the right-handed down quarks under the new U(1)",
- &HiddenValleyModel::uR_, 0.6, 10.0, 10.0,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyModel,double> interfaceZPrimeLL
- ("ZPrime/LL",
- "The charge of the left-handed leptons under the new U(1)",
- &HiddenValleyModel::lL_, 0.6, 10.0, 10.0,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyModel,double> interfaceZPrimeLR
- ("ZPrime/LR",
- "The charge of the right-handed charged leptons under the new U(1)",
- &HiddenValleyModel::lR_, -0.2, 10.0, 10.0,
- false, false, Interface::limited);
-
- static Parameter<HiddenValleyModel,double> interfaceGPrime
- ("GPrime",
- "The new gauge coupling",
- &HiddenValleyModel::gPrime_, 1./7., 0.0, 10.0,
- false, false, Interface::limited);
-
-}
-
-void HiddenValleyModel::doinit() {
- addVertex(FFZPVertex_);
- BSMModel::doinit();
- FFZPVertex_->init();
- if(qCharge_.size()!=Nf_)
- throw InitException() << "Number of new fermions and number of new charges "
- << "different in HiddenValleyModel::doinit()"
- << Exception::runerror;
-}
diff --git a/Shower/QTilde/Dark/HiddenValleyModel.h b/Shower/QTilde/Dark/HiddenValleyModel.h
deleted file mode 100644
--- a/Shower/QTilde/Dark/HiddenValleyModel.h
+++ /dev/null
@@ -1,338 +0,0 @@
-// -*- C++ -*-
-#ifndef HERWIG_HiddenValleyModel_H
-#define HERWIG_HiddenValleyModel_H
-//
-// This is the declaration of the HiddenValleyModel class.
-//
-
-#include "Herwig/Models/General/BSMModel.h"
-#include "ThePEG/Helicity/Vertex/AbstractFFVVertex.h"
-
-namespace Herwig {
-
-using namespace ThePEG;
-
-/**
- * Here is the documentation of the HiddenValleyModel class.
- *
- * @see \ref HiddenValleyModelInterfaces "The interfaces"
- * defined for HiddenValleyModel.
- */
-class HiddenValleyModel: public BSMModel {
-
-/**
- * Enumeration to define the type of the confining group
- */
- enum GroupType { SU, SO, SP };
-
-public:
-
- /**
- * The default constructor.
- */
- HiddenValleyModel();
-
-public:
-
- /**
- * Properties of the new confining gauge group
- */
- //@{
- /**
- * The hidden colour charge of the fundamental representation
- */
- double CF() const {
- switch (groupType_) {
- case SU:
- return 0.5*(sqr(double(Nc_))-1.)/double(Nc_);
- case SO:
- if(Nc_==3) {
- return (Nc_-1);
- }
- else {
- return 0.5*(Nc_-1);
- }
- case SP:
- return 0.25*(Nc_+1);
- default:
- assert(false);
- }
- return 0.;
- }
-
- /**
- * The hidden colour charge of the adjoint representation
- */
- double CA() const {
- switch (groupType_) {
- case SU:
- return double(Nc_);
- case SO:
- if(Nc_==3) {
- return 2.*(Nc_-2.);
- }
- else {
- return (Nc_-2.);
- }
- case SP:
- return 0.5*(Nc_+2.);
- default:
- assert(false);
- }
- return 0.;
- }
-
- /**
- * The \f$T_R\f$ colour factor
- */
- double TR() const {
- switch (groupType_) {
- case SU: case SP:
- return 0.5;
- break;
- case SO:
- if(Nc_==3) {
- return 2.;
- }
- else {
- return 1.;
- }
- default:
- assert(false);
- }
- return 0.;
- }
-
- /**
- * The type of the group
- */
- GroupType groupType() const {return groupType_;}
-
- /**
- * The number of colours
- */
- unsigned int NC() const {return Nc_;}
-
- /**
- * The number of fermions charged under the new group
- */
- unsigned int NF() const {return Nf_;}
- //@}
-
- /**
- * Properties of the new U(1) group
- */
- //@{
- /**
- * Charge of the left-handed quarks
- */
- double qL() const {return qL_;}
-
- /**
- * Charge of the right-handed up-type quarks
- */
- double uR() const {return uR_;}
-
- /**
- * Charge of the right-handed down-type quarks
- */
- double dR() const {return dR_;}
-
- /**
- * Charge of the left-handed leptons
- */
- double lL() const {return lL_;}
-
- /**
- * Charge of the right-handed charged leptons
- */
- double lR() const {return lR_;}
-
- /**
- * Coupling
- */
- double gPrime() const {return gPrime_;}
-
- /**
- * The couplings of the quirks
- */
- vector<double> qCharge() const {return qCharge_;}
- //@}
-
-public:
-
- /** @name Functions used by the persistent I/O system. */
- //@{
- /**
- * Function used to write out object persistently.
- * @param os the persistent output stream written to.
- */
- void persistentOutput(PersistentOStream & os) const;
-
- /**
- * Function used to read in object persistently.
- * @param is the persistent input stream read from.
- * @param version the version number of the object when written.
- */
- void persistentInput(PersistentIStream & is, int version);
- //@}
-
- /**
- * The standard Init function used to initialize the interfaces.
- * Called exactly once for each class by the class description system
- * before the main function starts or
- * when this class is dynamically loaded.
- */
- static void Init();
-
-protected:
-
- /** @name Clone Methods. */
- //@{
- /**
- * Make a simple clone of this object.
- * @return a pointer to the new object.
- */
- virtual IBPtr clone() const;
-
- /** Make a clone of this object, possibly modifying the cloned object
- * to make it sane.
- * @return a pointer to the new object.
- */
- virtual IBPtr fullclone() const;
- //@}
-
-private:
-
- /**
- * The static object used to initialize the description of this class.
- * Indicates that this is a concrete class with persistent data.
- */
- static ClassDescription<HiddenValleyModel> initHiddenValleyModel;
-
- /**
- * The assignment operator is private and must never be called.
- * In fact, it should not even be implemented.
- */
- HiddenValleyModel & operator=(const HiddenValleyModel &);
-
-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:
-
- /**
- * Symmetry group for the new unbroken symmetry
- */
- //@{
- /**
- * Type of group
- */
- GroupType groupType_;
-
- /**
- * Order of the group
- */
- unsigned int Nc_;
- //@}
-
- /**
- * Number of fermions
- */
- unsigned int Nf_;
-
- /**
- * Charges etc under the new \f$U(1)\f$ group
- */
- //@{
- /**
- * Charge of the left-handed quarks
- */
- double qL_;
-
- /**
- * Charge of the right-handed up-type quarks
- */
- double uR_;
-
- /**
- * Charge of the right-handed down-type quarks
- */
- double dR_;
-
- /**
- * Charge of the left-handed leptons
- */
- double lL_;
-
- /**
- * Charge of the right-handed charged leptons
- */
- double lR_;
-
- /**
- * Coupling
- */
- double gPrime_;
-
- /**
- * the vertex
- */
- AbstractFFVVertexPtr FFZPVertex_;
-
- /**
- * The couplings of the quirks
- */
- vector<double> qCharge_;
- //@}
-};
-
-}
-
-#include "ThePEG/Utilities/ClassTraits.h"
-
-namespace ThePEG {
-
-/** @cond TRAITSPECIALIZATIONS */
-
-/** This template specialization informs ThePEG about the
- * base classes of HiddenValleyModel. */
-template <>
-struct BaseClassTrait<Herwig::HiddenValleyModel,1> {
- /** Typedef of the first base class of HiddenValleyModel. */
- typedef Herwig::BSMModel NthBase;
-};
-
-/** This template specialization informs ThePEG about the name of
- * the HiddenValleyModel class and the shared object where it is defined. */
-template <>
-struct ClassTraits<Herwig::HiddenValleyModel>
- : public ClassTraitsBase<Herwig::HiddenValleyModel> {
- /** Return a platform-independent class name */
- static string className() { return "Herwig::HiddenValleyModel"; }
- /**
- * The name of a file containing the dynamic library where the class
- * HiddenValleyModel is implemented. It may also include several, space-separated,
- * libraries if the class HiddenValleyModel depends on other classes (base classes
- * excepted). In this case the listed libraries will be dynamically
- * linked in the order they are specified.
- */
- static string library() { return "HwShower.so HwHiddenValley.so"; }
-};
-
-ThePEG_DECLARE_POINTERS(Herwig::HiddenValleyModel,HiddenValleyPtr);
-
-/** @endcond */
-
-}
-
-#endif /* HERWIG_HiddenValleyModel_H */
diff --git a/Shower/QTilde/Dark/OneHalfHalfDarkSplitFn.cc b/Shower/QTilde/Dark/OneHalfHalfDarkSplitFn.cc
deleted file mode 100644
--- a/Shower/QTilde/Dark/OneHalfHalfDarkSplitFn.cc
+++ /dev/null
@@ -1,187 +0,0 @@
-// -*- C++ -*-
-//
-// OneHalfHalfDarkSplitFn.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
-// Copyright (C) 2002-2019 The Herwig Collaboration
-//
-// Herwig is licenced under version 3 of the GPL, see COPYING for details.
-// Please respect the MCnet academic guidelines, see GUIDELINES for details.
-//
-//
-// This is the implementation of the non-inlined, non-templated member
-// functions of the OneHalfHalfDarkSplitFn class.
-//
-
-#include "OneHalfHalfDarkSplitFn.h"
-#include "HiddenValleyModel.h"
-#include "ThePEG/PDT/ParticleData.h"
-#include "ThePEG/Interface/ClassDocumentation.h"
-#include "ThePEG/Utilities/DescribeClass.h"
-
-using namespace Herwig;
-
-DescribeNoPIOClass<OneHalfHalfDarkSplitFn,Herwig::Sudakov1to2FormFactor>
-describeOneHalfHalfDarkSplitFn ("Herwig::OneHalfHalfDarkSplitFn","HwShower.so");
-
-void OneHalfHalfDarkSplitFn::Init() {
-
- static ClassDocumentation<OneHalfHalfDarkSplitFn> documentation
- ("The OneHalfHalfDarkSplitFn class implements the splitting function for g->q qbar");
-
-}
-
-bool OneHalfHalfDarkSplitFn::checkColours(const IdList & ids) const {
- if(colourStructure()==TripletTripletOctet) {
- if(ids[0]!=ids[1]) return false;
- if((ids[0]->iColour()==PDT::DarkColourFundamental||
- ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
- ids[2]->iColour()==PDT::DarkColourAdjoint) return true;
- return false;
- }
- else if(colourStructure()==OctetOctetOctet) {
- for(unsigned int ix=0;ix<3;++ix) {
- if(ids[ix]->iColour()!=PDT::DarkColourAdjoint) return false;
- }
- return true;
- }
- else if(colourStructure()==OctetTripletTriplet) {
- if(ids[0]->iColour()!=PDT::DarkColourAdjoint) return false;
- if(ids[1]->iColour()==PDT::DarkColourFundamental&&
- ids[2]->iColour()==PDT::DarkColourAntiFundamental)
- return true;
- if(ids[1]->iColour()==PDT::DarkColourAntiFundamental&&
- ids[2]->iColour()==PDT::DarkColourFundamental)
- return true;
- return false;
- }
- else if(colourStructure()==TripletOctetTriplet) {
- if(ids[0]!=ids[2]) return false;
- if((ids[0]->iColour()==PDT::DarkColourFundamental||
- ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
- ids[1]->iColour()==PDT::DarkColourAdjoint) return true;
- return false;
- }
- else {
- assert(false);
- }
- return false;
-}
-
-double OneHalfHalfDarkSplitFn::P(const double z, const Energy2 t,
- const IdList &ids, const bool mass, const RhoDMatrix &) const {
- double zz = z*(1.-z);
- double val=1.-2.*zz;
- if(mass) {
- Energy m = ids[1]->mass();
- val +=2.*sqr(m)/t;
- }
- return val;
-}
-
-double OneHalfHalfDarkSplitFn::overestimateP(const double,
- const IdList &) const {
- return 1.;
-}
-
-double OneHalfHalfDarkSplitFn::ratioP(const double z, const Energy2 t,
- const IdList &ids, const bool mass, const RhoDMatrix &) const {
- double zz = z*(1.-z);
- double val = 1.-2.*zz;
- if(mass) {
- Energy m = ids[1]->mass();
- val+= 2.*sqr(m)/t;
- }
- return val;
-}
-
-double OneHalfHalfDarkSplitFn::integOverP(const double z, const IdList & ,
- unsigned int PDFfactor) const {
- switch(PDFfactor) {
- case 0:
- return z;
- case 1:
- return log(z);
- case 2:
- return -log(1.-z);
- case 3:
- return log(z/(1.-z));
- case 4:
- return 2.*sqrt(z);
- case 5:
- return (2./3.)*z*sqrt(z);
- default:
- throw Exception() << "OneHalfHalfDarkSplitFn::integOverP() invalid PDFfactor = "
- << PDFfactor << Exception::runerror;
- }
-}
-
-double OneHalfHalfDarkSplitFn::invIntegOverP(const double r,
- const IdList & ,
- unsigned int PDFfactor) const {
- switch(PDFfactor) {
- case 0:
- return r;
- case 1:
- return exp(r);
- case 2:
- return 1.-exp(-r);
- case 3:
- return 1./(1.+exp(-r));
- case 4:
- return 0.25*sqr(r);
- case 5:
- return pow(1.5*r,2./3.);
- default:
- throw Exception() << "OneHalfHalfDarkSplitFn::integOverP() invalid PDFfactor = "
- << PDFfactor << Exception::runerror;
- }
-}
-
-bool OneHalfHalfDarkSplitFn::accept(const IdList &ids) const {
- if(ids.size()!=3) return false;
- if(ids[1]!=ids[2]->CC()) return false;
- if(ids[1]->iSpin()!=PDT::Spin1Half) return false;
- if(ids[0]->iSpin()!=PDT::Spin1) return false;
- return OneHalfHalfDarkSplitFn::checkColours(ids);
-}
-
-vector<pair<int, Complex> >
-OneHalfHalfDarkSplitFn::generatePhiForward(const double z, const Energy2 t, const IdList & ids,
- const RhoDMatrix & rho) {
- assert(rho.iSpin()==PDT::Spin1);
- double modRho = abs(rho(0,2));
- Energy mq = ids[1]->mass();
- Energy2 mq2 = sqr(mq);
- double fact = z*(1.-z)-mq2/t;
- double max = 1.+2.*fact*(-1.+2.*modRho);
- vector<pair<int, Complex> > output;
- output.push_back(make_pair( 0,(rho(0,0)+rho(2,2))*(1.-2.*fact)/max));
- output.push_back(make_pair(-2,2.*fact*rho(0,2)/max));
- output.push_back(make_pair( 2,2.*fact*rho(2,0)/max));
- return output;
-}
-
-vector<pair<int, Complex> >
-OneHalfHalfDarkSplitFn::generatePhiBackward(const double, const Energy2, const IdList &,
- const RhoDMatrix & ) {
- // no dependance
- return {{ {0, 1.} }};
-}
-
-DecayMEPtr OneHalfHalfDarkSplitFn::matrixElement(const double z, const Energy2 t,
- const IdList & ids, const double phi,
- bool timeLike) {
- static const Complex ii(0.,1.);
- // calculate the kernal
- DecayMEPtr kernal(new_ptr(TwoBodyDecayMatrixElement(PDT::Spin1,PDT::Spin1Half,PDT::Spin1Half)));
- double mt = !timeLike ? ZERO : ids[1]->mass()/sqrt(t);
- double root =sqrt(1.-sqr(mt)/z/(1.-z));
- (*kernal)(0,0,0) = mt/sqrt(z*(1.-z));
- (*kernal)(2,1,1) = (*kernal)(0,0,0);
- (*kernal)(0,0,1) = -z*root*exp(-ii*phi);
- (*kernal)(2,1,0) = -conj((*kernal)(0,0,1));
- (*kernal)(0,1,0) = (1.-z)*exp(-ii*phi)*root;
- (*kernal)(2,0,1) = -conj((*kernal)(0,1,0));
- (*kernal)(0,1,1) = 0.;
- (*kernal)(2,0,0) = 0.;
- return kernal;
-}
diff --git a/Shower/QTilde/Dark/OneHalfHalfDarkSplitFn.h b/Shower/QTilde/Dark/OneHalfHalfDarkSplitFn.h
deleted file mode 100644
--- a/Shower/QTilde/Dark/OneHalfHalfDarkSplitFn.h
+++ /dev/null
@@ -1,194 +0,0 @@
-// -*- C++ -*-
-//
-// OneHalfHalfDarkSplitFn.h is a part of Herwig - A multi-purpose Monte Carlo event generator
-// Copyright (C) 2002-2019 The Herwig Collaboration
-//
-// Herwig is licenced under version 3 of the GPL, see COPYING for details.
-// Please respect the MCnet academic guidelines, see GUIDELINES for details.
-//
-#ifndef HERWIG_OneHalfHalfDarkSplitFn_H
-#define HERWIG_OneHalfHalfDarkSplitFn_H
-//
-// This is the declaration of the OneHalfHalfDarkSplitFn class.
-//
-
-#include "Herwig/Shower/QTilde/SplittingFunctions/Sudakov1to2FormFactor.h"
-#include "Herwig/Decay/TwoBodyDecayMatrixElement.h"
-
-namespace Herwig {
-
-using namespace ThePEG;
-
-/**\ingroup Shower
- *
- * This class provides the concrete implementation of the exact leading-order
- * splitting function for \f$1\to \frac12\frac12\f$.
- *
- * In this case the splitting function is given by
- * \f[P(z,t) =C\left(1-2z(1-z)+2\frac{m_q^2}{t}\right),\f]
- * where \f$C\f$ is the corresponding colour factor
- * Our choice for the overestimate is
- * \f[P_{\rm over}(z) = C,\f]
- * therefore the integral is
- * \f[\int P_{\rm over}(z) {\rm d}z =Cz,\f]
- * and its inverse is
- * \f[\frac{r}{C}\f]
- *
- * @see \ref OneHalfHalfDarkSplitFnInterfaces "The interfaces"
- * defined for OneHalfHalfDarkSplitFn.
- */
-class OneHalfHalfDarkSplitFn: public Sudakov1to2FormFactor {
-
-public:
-
- /**
- * Method to check the colours are correct
- */
- virtual bool checkColours(const IdList & ids) const;
-
- /**
- * Concrete implementation of the method to determine whether this splitting
- * function can be used for a given set of particles.
- * @param ids The PDG codes for the particles in the splitting.
- */
- virtual bool accept(const IdList & ids) const;
-
- /**
- * Methods to return the splitting function.
- */
- //@{
- /**
- * The concrete implementation of the splitting function, \f$P\f$.
- * @param z The energy fraction.
- * @param t The scale.
- * @param ids The PDG codes for the particles in the splitting.
- * @param mass Whether or not to include the mass dependent terms
- * @param rho The spin density matrix
- */
- virtual double P(const double z, const Energy2 t, const IdList & ids,
- const bool mass, const RhoDMatrix & rho) const;
-
-
- /**
- * The concrete implementation of the overestimate of the splitting function,
- * \f$P_{\rm over}\f$.
- * @param z The energy fraction.
- * @param ids The PDG codes for the particles in the splitting.
- */
- virtual double overestimateP(const double z, const IdList & ids) const;
-
- /**
- * The concrete implementation of the
- * the ratio of the splitting function to the overestimate, i.e.
- * \f$P(z,\tilde{q}^2)/P_{\rm over}(z)\f$.
- * @param z The energy fraction.
- * @param t The scale.
- * @param ids The PDG codes for the particles in the splitting.
- * @param mass Whether or not to include the mass dependent terms
- * @param rho The spin density matrix
- */
- virtual double ratioP(const double z, const Energy2 t, const IdList & ids,
- const bool mass, const RhoDMatrix & rho) const;
-
- /**
- * The concrete implementation of the indefinite integral of the
- * overestimated splitting function, \f$P_{\rm over}\f$.
- * @param z The energy fraction.
- * @param ids The PDG codes for the particles in the splitting.
- * @param PDFfactor Which additional factor to include for the PDF
- * 0 is no additional factor,
- * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
- */
- virtual double integOverP(const double z, const IdList & ids,
- unsigned int PDFfactor=0) const;
-
- /**
- * The concrete implementation of the inverse of the indefinite integral.
- * @param r Value of the splitting function to be inverted
- * @param ids The PDG codes for the particles in the splitting.
- * @param PDFfactor Which additional factor to include for the PDF
- * 0 is no additional factor,
- * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
- */
- virtual double invIntegOverP(const double r, const IdList & ids,
- unsigned int PDFfactor=0) const;
- //@}
-
- /**
- * Method to calculate the azimuthal angle
- * @param particle The particle which is branching
- * @param showerkin The ShowerKinematics object
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- * @return The weight
- */
- virtual vector<pair<int,Complex> >
- generatePhiForward(const double z, const Energy2 t, const IdList & ids,
- const RhoDMatrix &);
-
- /**
- * Method to calculate the azimuthal angle
- * @param particle The particle which is branching
- * @param showerkin The ShowerKinematics object
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- * @return The weight
- */
- virtual vector<pair<int,Complex> >
- generatePhiBackward(const double z, const Energy2 t, const IdList & ids,
- const RhoDMatrix &);
-
- /**
- * Calculate the matrix element for the splitting
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- */
- virtual DecayMEPtr matrixElement(const double z, const Energy2 t,
- const IdList & ids, const double phi, bool timeLike);
-
-public:
-
- /**
- * The standard Init function used to initialize the interfaces.
- * Called exactly once for each class by the class description system
- * before the main function starts or
- * when this class is dynamically loaded.
- */
- static void Init();
-
-protected:
-
- /** @name Clone Methods. */
- //@{
- /**
- * Make a simple clone of this object.
- * @return a pointer to the new object.
- */
- virtual IBPtr clone() const {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);}
- //@}
-
-private:
-
- /**
- * The assignment operator is private and must never be called.
- * In fact, it should not even be implemented.
- */
- OneHalfHalfDarkSplitFn & operator=(const OneHalfHalfDarkSplitFn &) = delete;
-
-};
-
-}
-
-#endif /* HERWIG_OneHalfHalfDarkSplitFn_H */
diff --git a/Shower/QTilde/Dark/OneOneOneDarkSplitFn.cc b/Shower/QTilde/Dark/OneOneOneDarkSplitFn.cc
deleted file mode 100644
--- a/Shower/QTilde/Dark/OneOneOneDarkSplitFn.cc
+++ /dev/null
@@ -1,167 +0,0 @@
-// -*- C++ -*-
-//
-// OneOneOneDarkSplitFn.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
-// Copyright (C) 2002-2019 The Herwig Collaboration
-//
-// Herwig is licenced under version 3 of the GPL, see COPYING for details.
-// Please respect the MCnet academic guidelines, see GUIDELINES for details.
-//
-//
-// This is the implementation of the non-inlined, non-templated member
-// functions of the OneOneOneDarkSplitFn class.
-//
-
-#include "OneOneOneDarkSplitFn.h"
-#include "HiddenValleyModel.h"
-#include "ThePEG/PDT/ParticleData.h"
-#include "ThePEG/Interface/ClassDocumentation.h"
-#include "ThePEG/Utilities/DescribeClass.h"
-
-using namespace Herwig;
-
-DescribeNoPIOClass<OneOneOneDarkSplitFn,Herwig::Sudakov1to2FormFactor>
-describeOneOneOneDarkSplitFn ("Herwig::OneOneOneDarkSplitFn","HwShower.so");
-
-void OneOneOneDarkSplitFn::Init() {
-
- static ClassDocumentation<OneOneOneDarkSplitFn> documentation
- ("The OneOneOneDarkSplitFn class implements the g -> gg splitting function");
-
-}
-
-bool OneOneOneDarkSplitFn::checkColours(const IdList & ids) const {
- if(colourStructure()==TripletTripletOctet) {
- if(ids[0]!=ids[1]) return false;
- if((ids[0]->iColour()==PDT::DarkColourFundamental||
- ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
- ids[2]->iColour()==PDT::DarkColourAdjoint) return true;
- return false;
- }
- else if(colourStructure()==OctetOctetOctet) {
- for(unsigned int ix=0;ix<3;++ix) {
- if(ids[ix]->iColour()!=PDT::DarkColourAdjoint) return false;
- }
- return true;
- }
- else if(colourStructure()==OctetTripletTriplet) {
- if(ids[0]->iColour()!=PDT::DarkColourAdjoint) return false;
- if(ids[1]->iColour()==PDT::DarkColourFundamental&&
- ids[2]->iColour()==PDT::DarkColourAntiFundamental)
- return true;
- if(ids[1]->iColour()==PDT::DarkColourAntiFundamental&&
- ids[2]->iColour()==PDT::DarkColourFundamental)
- return true;
- return false;
- }
- else if(colourStructure()==TripletOctetTriplet) {
- if(ids[0]!=ids[2]) return false;
- if((ids[0]->iColour()==PDT::DarkColourFundamental||
- ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
- ids[1]->iColour()==PDT::DarkColourAdjoint) return true;
- return false;
- }
- else {
- assert(false);
- }
- return false;
-}
-
-double OneOneOneDarkSplitFn::P(const double z, const Energy2,
- const IdList & , const bool, const RhoDMatrix &)const {
- return sqr(1.-z*(1.-z))/(z*(1.-z));
-}
-
-double OneOneOneDarkSplitFn::overestimateP(const double z,
- const IdList & ) const {
- return (1/z + 1/(1.-z));
-}
-
-
-double OneOneOneDarkSplitFn::ratioP(const double z, const Energy2,
- const IdList & , const bool, const RhoDMatrix &) const {
- return sqr(1.-z*(1.-z));
-}
-
-double OneOneOneDarkSplitFn::invIntegOverP(const double r,
- const IdList & ,
- unsigned int PDFfactor) const {
- switch(PDFfactor) {
- case 0:
- return 1./(1.+exp(-r));
- case 1:
- case 2:
- case 3:
- default:
- throw Exception() << "OneOneOneDarkSplitFn::integOverP() invalid PDFfactor = "
- << PDFfactor << Exception::runerror;
- }
-}
-
-double OneOneOneDarkSplitFn::integOverP(const double z, const IdList & ,
- unsigned int PDFfactor) const {
- switch(PDFfactor) {
- case 0:
- assert(z>0.&&z<1.);
- return log(z/(1.-z));
- case 1:
- case 2:
- case 3:
- default:
- throw Exception() << "OneOneOneDarkSplitFn::integOverP() invalid PDFfactor = "
- << PDFfactor << Exception::runerror;
- }
-}
-
-bool OneOneOneDarkSplitFn::accept(const IdList & ids) const {
- if(ids.size()!=3) return false;
- for(unsigned int ix=0;ix<ids.size();++ix) {
- if(ids[0]->iSpin()!=PDT::Spin1) return false;
- }
- return OneOneOneDarkSplitFn::checkColours(ids);
-}
-
-vector<pair<int, Complex> >
-OneOneOneDarkSplitFn::generatePhiForward(const double z, const Energy2, const IdList &,
- const RhoDMatrix & rho) {
- assert(rho.iSpin()==PDT::Spin1);
- double modRho = abs(rho(0,2));
- double max = 2.*z*modRho*(1.-z)+sqr(1.-(1.-z)*z)/(z*(1.-z));
- vector<pair<int, Complex> > output;
- output.push_back(make_pair( 0,(rho(0,0)+rho(2,2))*sqr(1.-(1.-z)*z)/(z*(1.-z))/max));
- output.push_back(make_pair(-2,-rho(0,2)*z*(1.-z)/max));
- output.push_back(make_pair( 2,-rho(2,0)*z*(1.-z)/max));
- return output;
-}
-
-vector<pair<int, Complex> >
-OneOneOneDarkSplitFn::generatePhiBackward(const double z, const Energy2, const IdList &,
- const RhoDMatrix & rho) {
- assert(rho.iSpin()==PDT::Spin1);
- double diag = sqr(1 - (1 - z)*z)/(1 - z)/z;
- double off = (1.-z)/z;
- double max = 2.*abs(rho(0,2))*off+diag;
- vector<pair<int, Complex> > output;
- output.push_back(make_pair( 0, (rho(0,0)+rho(2,2))*diag/max));
- output.push_back(make_pair( 2,-rho(0,2) * off/max));
- output.push_back(make_pair(-2,-rho(2,0) * off/max));
- return output;
-}
-
-DecayMEPtr OneOneOneDarkSplitFn::matrixElement(const double z, const Energy2,
- const IdList &, const double phi,
- bool) {
- // calculate the kernal
- DecayMEPtr kernal(new_ptr(TwoBodyDecayMatrixElement(PDT::Spin1,PDT::Spin1,PDT::Spin1)));
- double omz = 1.-z;
- double root = sqrt(z*omz);
- Complex phase = exp(Complex(0.,1.)*phi);
- (*kernal)(0,0,0) = phase/root;
- (*kernal)(2,2,2) = -conj((*kernal)(0,0,0));
- (*kernal)(0,0,2) = -sqr(z)/root/phase;
- (*kernal)(2,2,0) = -conj((*kernal)(0,0,2));
- (*kernal)(0,2,0) = -sqr(omz)/root/phase;
- (*kernal)(2,0,2) = -conj((*kernal)(0,2,0));
- (*kernal)(0,2,2) = 0.;
- (*kernal)(2,0,0) = 0.;
- return kernal;
-}
diff --git a/Shower/QTilde/Dark/OneOneOneDarkSplitFn.h b/Shower/QTilde/Dark/OneOneOneDarkSplitFn.h
deleted file mode 100644
--- a/Shower/QTilde/Dark/OneOneOneDarkSplitFn.h
+++ /dev/null
@@ -1,190 +0,0 @@
-// -*- C++ -*-
-//
-// OneOneOneDarkSplitFn.h is a part of Herwig - A multi-purpose Monte Carlo event generator
-// Copyright (C) 2002-2019 The Herwig Collaboration
-//
-// Herwig is licenced under version 3 of the GPL, see COPYING for details.
-// Please respect the MCnet academic guidelines, see GUIDELINES for details.
-//
-#ifndef HERWIG_OneOneOneDarkSplitFn_H
-#define HERWIG_OneOneOneDarkSplitFn_H
-//
-// This is the declaration of the OneOneOneDarkSplitFn class.
-//
-
-#include "Herwig/Shower/QTilde/SplittingFunctions/Sudakov1to2FormFactor.h"
-#include "Herwig/Decay/TwoBodyDecayMatrixElement.h"
-
-namespace Herwig {
-
-using namespace ThePEG;
-
-/** \ingroup Shower
- *
- * This class provides the concrete implementation of the exact leading-order
- * splitting function for \f$1\to 11\f$.
- *
- * In this case the splitting function is given by
- * \f[P(z) =2C*\left(\frac{z}{1-z}+\frac{1-z}{z}+z(1-z)\right),\f]
- * where \f$C=\f$ is the corresponding colour factor.
- * Our choice for the overestimate is
- * \f[P_{\rm over}(z) = 2C\left(\frac1z+\frac1{1-z}\right),\f]
- * therefore the integral is
- * \f[\int P_{\rm over}(z) {\rm d}z =2C\ln\left(\frac{z}{1-z}\right),\f]
- * and its inverse is
- * \f[\frac{\exp\left(\frac{r}{2C}\right)}{(1+\exp\left(\frac{r}{2C}\right)}\f]
- *
- *
- * @see \ref OneOneOneDarkSplitFnInterfaces "The interfaces"
- * defined for OneOneOneDarkSplitFn.
- */
-class OneOneOneDarkSplitFn: public Sudakov1to2FormFactor {
-
-public:
-
- /**
- * Method to check the colours are correct
- */
- virtual bool checkColours(const IdList & ids) const;
-
- /**
- * Concrete implementation of the method to determine whether this splitting
- * function can be used for a given set of particles.
- * @param ids The PDG codes for the particles in the splitting.
- */
- virtual bool accept(const IdList & ids) const;
-
- /**
- * Methods to return the splitting function.
- */
- //@{
- /**
- * The concrete implementation of the splitting function, \f$P(z,t)\f$.
- * @param z The energy fraction.
- * @param t The scale.
- * @param ids The PDG codes for the particles in the splitting.
- * @param mass Whether or not to include the mass dependent terms
- * @param rho The spin density matrix
- */
- virtual double P(const double z, const Energy2 t, const IdList & ids,
- const bool mass, const RhoDMatrix & rho) const;
-
- /**
- * The concrete implementation of the overestimate of the splitting function,
- * \f$P_{\rm over}\f$.
- * @param z The energy fraction.
- * @param ids The PDG codes for the particles in the splitting.
- */
- virtual double overestimateP(const double z, const IdList & ids) const;
-
- /**
- * The concrete implementation of the
- * the ratio of the splitting function to the overestimate, i.e.
- * \f$P(z,t)/P_{\rm over}(z)\f$.
- * @param z The energy fraction.
- * @param t The scale.
- * @param ids The PDG codes for the particles in the splitting.
- * @param mass Whether or not to include the mass dependent terms
- * @param rho The spin density matrix
- */
- virtual double ratioP(const double z, const Energy2 t, const IdList & ids,
- const bool mass, const RhoDMatrix & rho) const;
-
- /**
- * The concrete implementation of the indefinite integral of the
- * overestimated splitting function, \f$P_{\rm over}\f$.
- * @param z The energy fraction.
- * @param ids The PDG codes for the particles in the splitting.
- * @param PDFfactor Which additional factor to include for the PDF
- * 0 is no additional factor,
- * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
- */
- virtual double integOverP(const double z, const IdList & ids,
- unsigned int PDFfactor=0) const;
-
- /**
- * The concrete implementation of the inverse of the indefinite integral.
- * @param r Value of the splitting function to be inverted
- * @param ids The PDG codes for the particles in the splitting.
- * @param PDFfactor Which additional factor to include for the PDF
- * 0 is no additional factor,
- * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
- */
- virtual double invIntegOverP(const double r, const IdList & ids,
- unsigned int PDFfactor=0) const;
- //@}
-
- /**
- * Method to calculate the azimuthal angle for forward evolution
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- * @return The weight
- */
- virtual vector<pair<int,Complex> >
- generatePhiForward(const double z, const Energy2 t, const IdList & ids,
- const RhoDMatrix &);
-
- /**
- * Method to calculate the azimuthal angle for backward evolution
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- * @return The weight
- */
- virtual vector<pair<int,Complex> >
- generatePhiBackward(const double z, const Energy2 t, const IdList & ids,
- const RhoDMatrix &);
-
- /**
- * Calculate the matrix element for the splitting
- * @param z The energy fraction
- * @param t The scale \f$t=2p_j\cdot p_k\f$.
- * @param ids The PDG codes for the particles in the splitting.
- * @param The azimuthal angle, \f$\phi\f$.
- */
- virtual DecayMEPtr matrixElement(const double z, const Energy2 t,
- const IdList & ids, const double phi, bool timeLike);
-
-public:
-
- /**
- * The standard Init function used to initialize the interfaces.
- * Called exactly once for each class by the class description system
- * before the main function starts or
- * when this class is dynamically loaded.
- */
- static void Init();
-
-protected:
-
- /** @name Clone Methods. */
- //@{
- /**
- * Make a simple clone of this object.
- * @return a pointer to the new object.
- */
- virtual IBPtr clone() const {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);}
- //@}
-
-private:
-
- /**
- * The assignment operator is private and must never be called.
- * In fact, it should not even be implemented.
- */
- OneOneOneDarkSplitFn & operator=(const OneOneOneDarkSplitFn &) = delete;
-
-};
-
-}
-
-#endif /* HERWIG_OneOneOneDarkSplitFn_H */
diff --git a/Shower/QTilde/Makefile.am b/Shower/QTilde/Makefile.am
--- a/Shower/QTilde/Makefile.am
+++ b/Shower/QTilde/Makefile.am
@@ -1,61 +1,56 @@
-SUBDIRS = Matching SplittingFunctions/Onium
+SUBDIRS = Matching SplittingFunctions/Onium SplittingFunctions/Dark
pkglib_LTLIBRARIES = HwShower.la
HwShower_la_LDFLAGS = $(AM_LDFLAGS) -module -version-info 30:0:0
HwShower_la_SOURCES = \
Couplings/ShowerAlphaQCD.h Couplings/ShowerAlphaQCD.cc \
Couplings/ShowerAlphaQED.h Couplings/ShowerAlphaQED.cc\
Couplings/ShowerAlphaQCDNP.h Couplings/ShowerAlphaQCDNP.cc\
+Couplings/HiddenValleyAlpha.cc Couplings/HiddenValleyAlpha.h \
QTildeShowerHandler.h QTildeShowerHandler.fh QTildeShowerHandler.cc \
SplittingFunctions/HalfHalfOneSplitFn.h SplittingFunctions/HalfHalfOneSplitFn.cc\
SplittingFunctions/HalfHalfOneEWSplitFn.h SplittingFunctions/HalfHalfOneEWSplitFn.cc\
SplittingFunctions/HalfHalfZeroEWSplitFn.h SplittingFunctions/HalfHalfZeroEWSplitFn.cc\
SplittingFunctions/OneOneOneEWSplitFn.h SplittingFunctions/OneOneOneEWSplitFn.cc\
SplittingFunctions/OneOneOneQEDSplitFn.h SplittingFunctions/OneOneOneQEDSplitFn.cc\
SplittingFunctions/OneOneZeroEWSplitFn.h SplittingFunctions/OneOneZeroEWSplitFn.cc\
SplittingFunctions/ZeroZeroZeroEWSplitFn.h SplittingFunctions/ZeroZeroZeroEWSplitFn.cc\
SplittingFunctions/OneZeroZeroEWSplitFn.h SplittingFunctions/OneZeroZeroEWSplitFn.cc\
SplittingFunctions/ZeroZeroOneEWSplitFn.h SplittingFunctions/ZeroZeroOneEWSplitFn.cc\
SplittingFunctions/OneOneOneSplitFn.h SplittingFunctions/OneOneOneSplitFn.cc\
SplittingFunctions/OneOneOneMassiveSplitFn.h SplittingFunctions/OneOneOneMassiveSplitFn.cc\
SplittingFunctions/ZeroZeroOneSplitFn.h SplittingFunctions/ZeroZeroOneSplitFn.cc\
SplittingFunctions/OneHalfHalfSplitFn.h SplittingFunctions/OneHalfHalfSplitFn.cc\
SplittingFunctions/HalfOneHalfSplitFn.h SplittingFunctions/HalfOneHalfSplitFn.cc\
SplittingFunctions/CMWOneOneOneSplitFn.h SplittingFunctions/CMWOneOneOneSplitFn.cc\
SplittingFunctions/CMWHalfHalfOneSplitFn.h SplittingFunctions/CMWHalfHalfOneSplitFn.cc\
-Dark/EnumParticles.h Dark/HiddenValleyAlpha.cc Dark/HiddenValleyAlpha.h \
-Dark/HiddenValleyModel.cc Dark/HiddenValleyModel.h \
-Dark/HiddenValleyFFZPrimeVertex.cc Dark/HiddenValleyFFZPrimeVertex.h \
-Dark/OneHalfHalfDarkSplitFn.cc Dark/OneHalfHalfDarkSplitFn.h \
-Dark/OneOneOneDarkSplitFn.cc Dark/OneOneOneDarkSplitFn.h \
-Dark/HalfHalfOneDarkSplitFn.cc Dark/HalfHalfOneDarkSplitFn.h \
Kinematics/Decay_QTildeShowerKinematics1to2.cc \
Kinematics/Decay_QTildeShowerKinematics1to2.h \
Kinematics/IS_QTildeShowerKinematics1to2.cc Kinematics/IS_QTildeShowerKinematics1to2.h \
Kinematics/FS_QTildeShowerKinematics1to2.cc Kinematics/FS_QTildeShowerKinematics1to2.h \
Kinematics/KinematicHelpers.h \
Kinematics/KinematicsReconstructor.cc \
Kinematics/KinematicsReconstructor.tcc \
Kinematics/KinematicsReconstructor.h \
Kinematics/KinematicsReconstructor.fh \
Base/HardTree.cc Base/HardTree.h Base/HardTree.fh \
Base/HardBranching.h Base/HardBranching.fh Base/HardBranching.cc\
Base/PartnerFinder.h Base/PartnerFinder.fh Base/PartnerFinder.cc \
Base/ShowerVeto.h Base/ShowerVeto.fh Base/ShowerVeto.cc \
Base/FullShowerVeto.h Base/FullShowerVeto.fh Base/FullShowerVeto.cc \
SplittingFunctions/SplittingGenerator.cc SplittingFunctions/SplittingGenerator.h\
SplittingFunctions/SplittingGenerator.fh \
Base/ShowerTree.h Base/ShowerTree.fh Base/ShowerTree.cc \
ShowerConfig.h ShowerConfig.cc \
Base/Branching.h \
Base/ShowerParticle.cc Base/ShowerParticle.fh Base/ShowerParticle.h \
Kinematics/ShowerKinematics.fh Kinematics/ShowerKinematics.h Kinematics/ShowerKinematics.cc \
Kinematics/ShowerBasis.fh Kinematics/ShowerBasis.h Kinematics/ShowerBasis.cc \
Base/ShowerProgenitor.fh Base/ShowerProgenitor.h \
SplittingFunctions/SudakovFormFactor.cc SplittingFunctions/SudakovFormFactor.h SplittingFunctions/SudakovFormFactor.fh \
SplittingFunctions/SudakovCutOff.cc SplittingFunctions/SudakovCutOff.h SplittingFunctions/SudakovCutOff.fh \
SplittingFunctions/PTCutOff.cc SplittingFunctions/PTCutOff.h \
SplittingFunctions/MassCutOff.cc SplittingFunctions/MassCutOff.h \
SplittingFunctions/VariableMassCutOff.cc SplittingFunctions/VariableMassCutOff.h \
SplittingFunctions/Sudakov1to2FormFactor.h SplittingFunctions/Sudakov1to2FormFactor.fh \
SplittingFunctions/Sudakov1to2FormFactor.cc \
Base/ShowerVertex.cc Base/ShowerVertex.fh Base/ShowerVertex.h
diff --git a/Shower/QTilde/SplittingFunctions/Dark/HalfHalfOneDarkSplitFn.cc b/Shower/QTilde/SplittingFunctions/Dark/HalfHalfOneDarkSplitFn.cc
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/SplittingFunctions/Dark/HalfHalfOneDarkSplitFn.cc
@@ -0,0 +1,171 @@
+// -*- C++ -*-
+//
+// HalfHalfOneDarkSplitFn.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
+// Copyright (C) 2002-2019 The Herwig Collaboration
+//
+// Herwig is licenced under version 3 of the GPL, see COPYING for details.
+// Please respect the MCnet academic guidelines, see GUIDELINES for details.
+//
+//
+// This is the implementation of the non-inlined, non-templated member
+// functions of the HalfHalfOneDarkSplitFn class.
+//
+
+#include "HalfHalfOneDarkSplitFn.h"
+#include "Herwig/Models/HiddenValley/HiddenValleyModel.h"
+#include "ThePEG/PDT/ParticleData.h"
+#include "ThePEG/Interface/ClassDocumentation.h"
+#include "ThePEG/Utilities/DescribeClass.h"
+#include <cassert>
+
+using namespace Herwig;
+
+DescribeNoPIOClass<HalfHalfOneDarkSplitFn,Herwig::Sudakov1to2FormFactor>
+describeHalfHalfOneDarkSplitFn ("Herwig::HalfHalfOneDarkSplitFn","HwDarkShower.so");
+
+void HalfHalfOneDarkSplitFn::Init() {
+ // documentation
+ static ClassDocumentation<HalfHalfOneDarkSplitFn> documentation
+ ("The HalfHalfOneDarkSplitFn class implements the dark coloured q -> qg splitting function");
+}
+
+bool HalfHalfOneDarkSplitFn::checkColours(const IdList & ids) const {
+ if(colourStructure()==TripletTripletOctet) {
+ if(ids[0]!=ids[1]) return false;
+ if((ids[0]->iColour()==PDT::DarkColourFundamental||
+ ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
+ ids[2]->iColour()==PDT::DarkColourAdjoint) return true;
+ return false;
+ }
+ else if(colourStructure()==OctetOctetOctet) {
+ for(unsigned int ix=0;ix<3;++ix) {
+ if(ids[ix]->iColour()!=PDT::DarkColourAdjoint) return false;
+ }
+ return true;
+ }
+ else if(colourStructure()==OctetTripletTriplet) {
+ if(ids[0]->iColour()!=PDT::DarkColourAdjoint) return false;
+ if(ids[1]->iColour()==PDT::DarkColourFundamental&&
+ ids[2]->iColour()==PDT::DarkColourAntiFundamental)
+ return true;
+ if(ids[1]->iColour()==PDT::DarkColourAntiFundamental&&
+ ids[2]->iColour()==PDT::DarkColourFundamental)
+ return true;
+ return false;
+ }
+ else if(colourStructure()==TripletOctetTriplet) {
+ if(ids[0]!=ids[2]) return false;
+ if((ids[0]->iColour()==PDT::DarkColourFundamental||
+ ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
+ ids[1]->iColour()==PDT::DarkColourAdjoint) return true;
+ return false;
+ }
+ else {
+ assert(false);
+ }
+ return false;
+}
+
+double HalfHalfOneDarkSplitFn::P(const double z, const Energy2 t,
+ const IdList &ids, const bool mass, const RhoDMatrix & ) const {
+ double val = (1. + sqr(z))/(1.-z);
+ if(mass) {
+ Energy m = ids[0]->mass();
+ val -= 2.*sqr(m)/t;
+ }
+ return val;
+}
+
+double HalfHalfOneDarkSplitFn::overestimateP(const double z,
+ const IdList & ) const {
+ return 2./(1.-z);
+}
+
+double HalfHalfOneDarkSplitFn::ratioP(const double z, const Energy2 t,
+ const IdList & ids, const bool mass, const RhoDMatrix & ) const {
+ double val = 1. + sqr(z);
+ if(mass) {
+ Energy m = ids[0]->mass();
+ val -= 2.*sqr(m)*(1.-z)/t;
+ }
+ return 0.5*val;
+}
+
+double HalfHalfOneDarkSplitFn::integOverP(const double z,
+ const IdList & ,
+ unsigned int PDFfactor) const {
+ switch (PDFfactor) {
+ case 0:
+ return -2.*Math::log1m(z);
+ case 1:
+ return 2.*log(z/(1.-z));
+ case 2:
+ return 2./(1.-z);
+ case 3:
+ default:
+ throw Exception() << "HalfHalfOneDarkSplitFn::integOverP() invalid PDFfactor = "
+ << PDFfactor << Exception::runerror;
+ }
+}
+
+double HalfHalfOneDarkSplitFn::invIntegOverP(const double r, const IdList &,
+ unsigned int PDFfactor) const {
+ switch (PDFfactor) {
+ case 0:
+ return 1. - exp(- 0.5*r);
+ case 1:
+ return 1./(1.-exp(-0.5*r));
+ case 2:
+ return 1.-2./r;
+ case 3:
+ default:
+ throw Exception() << "HalfHalfOneDarkSplitFn::invIntegOverP() invalid PDFfactor = "
+ << PDFfactor << Exception::runerror;
+ }
+}
+
+bool HalfHalfOneDarkSplitFn::accept(const IdList &ids) const {
+ // 3 particles and in and out fermion same
+ if(ids.size()!=3 || ids[0]!=ids[1]) return false;
+ if(ids[0]->iSpin()!=PDT::Spin1Half ||
+ ids[2]->iSpin()!=PDT::Spin1) return false;
+ return checkColours(ids);
+}
+
+vector<pair<int, Complex> >
+HalfHalfOneDarkSplitFn::generatePhiForward(const double, const Energy2, const IdList & ,
+ const RhoDMatrix &) {
+ // no dependence on the spin density matrix, dependence on off-diagonal terms cancels
+ // and rest = splitting function for Tr(rho)=1 as required by defn
+ return {{ {0, 1.} }};
+}
+
+vector<pair<int, Complex> >
+HalfHalfOneDarkSplitFn::generatePhiBackward(const double, const Energy2, const IdList & ,
+ const RhoDMatrix &) {
+ // no dependence on the spin density matrix, dependence on off-diagonal terms cancels
+ // and rest = splitting function for Tr(rho)=1 as required by defn
+ return {{ {0, 1.} }};
+}
+
+DecayMEPtr HalfHalfOneDarkSplitFn::matrixElement(const double z, const Energy2 t,
+ const IdList & ids, const double phi,
+ bool timeLike) {
+ // calculate the kernal
+ DecayMEPtr kernal(new_ptr(TwoBodyDecayMatrixElement(PDT::Spin1Half,PDT::Spin1Half,PDT::Spin1)));
+ Energy m = !timeLike ? ZERO : ids[0]->mass();
+ double mt = m/sqrt(t);
+ double root = sqrt(1.-(1.-z)*sqr(m)/z/t);
+ double romz = sqrt(1.-z);
+ double rz = sqrt(z);
+ Complex phase = exp(Complex(0.,1.)*phi);
+ (*kernal)(0,0,0) = -root/romz*phase;
+ (*kernal)(1,1,2) = -conj((*kernal)(0,0,0));
+ (*kernal)(0,0,2) = root/romz*z/phase;
+ (*kernal)(1,1,0) = -conj((*kernal)(0,0,2));
+ (*kernal)(1,0,2) = mt*(1.-z)/rz;
+ (*kernal)(0,1,0) = conj((*kernal)(1,0,2));
+ (*kernal)(0,1,2) = 0.;
+ (*kernal)(1,0,0) = 0.;
+ return kernal;
+}
diff --git a/Shower/QTilde/SplittingFunctions/Dark/HalfHalfOneDarkSplitFn.h b/Shower/QTilde/SplittingFunctions/Dark/HalfHalfOneDarkSplitFn.h
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/SplittingFunctions/Dark/HalfHalfOneDarkSplitFn.h
@@ -0,0 +1,189 @@
+// -*- C++ -*-
+//
+// HalfHalfOneDarkSplitFn.h is a part of Herwig - A multi-purpose Monte Carlo event generator
+// Copyright (C) 2002-2019 The Herwig Collaboration
+//
+// Herwig is licenced under version 3 of the GPL, see COPYING for details.
+// Please respect the MCnet academic guidelines, see GUIDELINES for details.
+//
+#ifndef HERWIG_HalfHalfOneDarkSplitFn_H
+#define HERWIG_HalfHalfOneDarkSplitFn_H
+//
+// This is the declaration of the HalfHalfOneDarkSplitFn class.
+//
+
+#include "Herwig/Shower/QTilde/SplittingFunctions/Sudakov1to2FormFactor.h"
+#include "Herwig/Decay/TwoBodyDecayMatrixElement.h"
+
+namespace Herwig {
+
+using namespace ThePEG;
+
+/**\ingroup Shower
+ *
+ * This class provides the concrete implementation of the exact leading-order
+ * splitting function for \f$\frac12\to q\frac12 1\f$.
+ *
+ * In this case the splitting function is given by
+ * \f[P(z,t) =C\left(\frac{1+z^2}{1-z}-2\frac{m^2_q}{t}\right),\f]
+ * where \f$C\f$ is the corresponding colour factor.
+ * Our choice for the overestimate is
+ * \f[P_{\rm over}(z) = \frac{2C}{1-z},\f]
+ * therefore the integral is
+ * \f[\int P_{\rm over}(z) {\rm d}z = -2C\ln(1-z),\f]
+ * and its inverse is
+ * \f[1-\exp\left(\frac{r}{2C}\right).\f]
+ *
+ * @see \ref HalfHalfOneDarkSplitFnInterfaces "The interfaces"
+ * defined for HalfHalfOneDarkSplitFn.
+ */
+class HalfHalfOneDarkSplitFn: public Sudakov1to2FormFactor {
+
+public:
+
+ /**
+ * Method to check the colours are correct
+ */
+ virtual bool checkColours(const IdList & ids) const;
+
+ /**
+ * Concrete implementation of the method to determine whether this splitting
+ * function can be used for a given set of particles.
+ * @param ids The PDG codes for the particles in the splitting.
+ */
+ virtual bool accept(const IdList & ids) const;
+
+ /**
+ * Methods to return the splitting function.
+ */
+ //@{
+ /**
+ * The concrete implementation of the splitting function, \f$P(z,t)\f$.
+ * @param z The energy fraction.
+ * @param t The scale.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param mass Whether or not to include the mass dependent terms
+ * @param rho The spin density matrix
+ */
+ virtual double P(const double z, const Energy2 t, const IdList & ids,
+ const bool mass, const RhoDMatrix & rho) const;
+
+ /**
+ * The concrete implementation of the overestimate of the splitting function,
+ * \f$P_{\rm over}\f$.
+ * @param z The energy fraction.
+ * @param ids The PDG codes for the particles in the splitting.
+ */
+ virtual double overestimateP(const double z, const IdList & ids) const;
+
+ /**
+ * The concrete implementation of the
+ * the ratio of the splitting function to the overestimate, i.e.
+ * \f$P(z,t)/P_{\rm over}(z)\f$.
+ * @param z The energy fraction.
+ * @param t The scale.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param mass Whether or not to include the mass dependent terms
+ * @param rho The spin density matrix
+ */
+ virtual double ratioP(const double z, const Energy2 t, const IdList & ids,
+ const bool mass, const RhoDMatrix & rho) const;
+
+ /**
+ * The concrete implementation of the indefinite integral of the
+ * overestimated splitting function, \f$P_{\rm over}\f$.
+ * @param z The energy fraction.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param PDFfactor Which additional factor to include for the PDF
+ * 0 is no additional factor,
+ * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
+ */
+ virtual double integOverP(const double z, const IdList & ids,
+ unsigned int PDFfactor=0) const;
+
+ /**
+ * The concrete implementation of the inverse of the indefinite integral.
+ * @param r Value of the splitting function to be inverted
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param PDFfactor Which additional factor to include for the PDF
+ * 0 is no additional factor,
+ * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
+ */
+ virtual double invIntegOverP(const double r, const IdList & ids,
+ unsigned int PDFfactor=0) const;
+ //@}
+
+ /**
+ * Method to calculate the azimuthal angle
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ * @return The weight
+ */
+ virtual vector<pair<int,Complex> >
+ generatePhiForward(const double z, const Energy2 t, const IdList & ids,
+ const RhoDMatrix &);
+
+ /**
+ * Method to calculate the azimuthal angle for backward evolution
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ * @return The weight
+ */
+ virtual vector<pair<int,Complex> >
+ generatePhiBackward(const double z, const Energy2 t, const IdList & ids,
+ const RhoDMatrix &);
+
+ /**
+ * Calculate the matrix element for the splitting
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ */
+ virtual DecayMEPtr matrixElement(const double z, const Energy2 t,
+ const IdList & ids, const double phi, bool timeLike);
+
+public:
+
+ /**
+ * The standard Init function used to initialize the interfaces.
+ * Called exactly once for each class by the class description system
+ * before the main function starts or
+ * when this class is dynamically loaded.
+ */
+ static void Init();
+
+protected:
+
+ /** @name Clone Methods. */
+ //@{
+ /**
+ * Make a simple clone of this object.
+ * @return a pointer to the new object.
+ */
+ virtual IBPtr clone() const {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);}
+ //@}
+
+private:
+
+ /**
+ * The assignment operator is private and must never be called.
+ * In fact, it should not even be implemented.
+ */
+ HalfHalfOneDarkSplitFn & operator=(const HalfHalfOneDarkSplitFn &) = delete;
+
+};
+
+}
+
+#endif /* HERWIG_HalfHalfOneDarkSplitFn_H */
diff --git a/Shower/QTilde/SplittingFunctions/Dark/Makefile.am b/Shower/QTilde/SplittingFunctions/Dark/Makefile.am
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/SplittingFunctions/Dark/Makefile.am
@@ -0,0 +1,8 @@
+pkglib_LTLIBRARIES = HwDarkShower.la
+
+HwDarkShower_la_LDFLAGS = $(AM_LDFLAGS) -module -version-info 1:0:0
+
+HwDarkShower_la_SOURCES = \
+ OneHalfHalfDarkSplitFn.cc OneHalfHalfDarkSplitFn.h \
+ OneOneOneDarkSplitFn.cc OneOneOneDarkSplitFn.h \
+ HalfHalfOneDarkSplitFn.cc HalfHalfOneDarkSplitFn.h
diff --git a/Shower/QTilde/SplittingFunctions/Dark/OneHalfHalfDarkSplitFn.cc b/Shower/QTilde/SplittingFunctions/Dark/OneHalfHalfDarkSplitFn.cc
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/SplittingFunctions/Dark/OneHalfHalfDarkSplitFn.cc
@@ -0,0 +1,187 @@
+// -*- C++ -*-
+//
+// OneHalfHalfDarkSplitFn.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
+// Copyright (C) 2002-2019 The Herwig Collaboration
+//
+// Herwig is licenced under version 3 of the GPL, see COPYING for details.
+// Please respect the MCnet academic guidelines, see GUIDELINES for details.
+//
+//
+// This is the implementation of the non-inlined, non-templated member
+// functions of the OneHalfHalfDarkSplitFn class.
+//
+
+#include "OneHalfHalfDarkSplitFn.h"
+#include "Herwig/Models/HiddenValley/HiddenValleyModel.h"
+#include "ThePEG/PDT/ParticleData.h"
+#include "ThePEG/Interface/ClassDocumentation.h"
+#include "ThePEG/Utilities/DescribeClass.h"
+
+using namespace Herwig;
+
+DescribeNoPIOClass<OneHalfHalfDarkSplitFn,Herwig::Sudakov1to2FormFactor>
+describeOneHalfHalfDarkSplitFn ("Herwig::OneHalfHalfDarkSplitFn","HwDarkShower.so");
+
+void OneHalfHalfDarkSplitFn::Init() {
+
+ static ClassDocumentation<OneHalfHalfDarkSplitFn> documentation
+ ("The OneHalfHalfDarkSplitFn class implements the splitting function for g->q qbar");
+
+}
+
+bool OneHalfHalfDarkSplitFn::checkColours(const IdList & ids) const {
+ if(colourStructure()==TripletTripletOctet) {
+ if(ids[0]!=ids[1]) return false;
+ if((ids[0]->iColour()==PDT::DarkColourFundamental||
+ ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
+ ids[2]->iColour()==PDT::DarkColourAdjoint) return true;
+ return false;
+ }
+ else if(colourStructure()==OctetOctetOctet) {
+ for(unsigned int ix=0;ix<3;++ix) {
+ if(ids[ix]->iColour()!=PDT::DarkColourAdjoint) return false;
+ }
+ return true;
+ }
+ else if(colourStructure()==OctetTripletTriplet) {
+ if(ids[0]->iColour()!=PDT::DarkColourAdjoint) return false;
+ if(ids[1]->iColour()==PDT::DarkColourFundamental&&
+ ids[2]->iColour()==PDT::DarkColourAntiFundamental)
+ return true;
+ if(ids[1]->iColour()==PDT::DarkColourAntiFundamental&&
+ ids[2]->iColour()==PDT::DarkColourFundamental)
+ return true;
+ return false;
+ }
+ else if(colourStructure()==TripletOctetTriplet) {
+ if(ids[0]!=ids[2]) return false;
+ if((ids[0]->iColour()==PDT::DarkColourFundamental||
+ ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
+ ids[1]->iColour()==PDT::DarkColourAdjoint) return true;
+ return false;
+ }
+ else {
+ assert(false);
+ }
+ return false;
+}
+
+double OneHalfHalfDarkSplitFn::P(const double z, const Energy2 t,
+ const IdList &ids, const bool mass, const RhoDMatrix &) const {
+ double zz = z*(1.-z);
+ double val=1.-2.*zz;
+ if(mass) {
+ Energy m = ids[1]->mass();
+ val +=2.*sqr(m)/t;
+ }
+ return val;
+}
+
+double OneHalfHalfDarkSplitFn::overestimateP(const double,
+ const IdList &) const {
+ return 1.;
+}
+
+double OneHalfHalfDarkSplitFn::ratioP(const double z, const Energy2 t,
+ const IdList &ids, const bool mass, const RhoDMatrix &) const {
+ double zz = z*(1.-z);
+ double val = 1.-2.*zz;
+ if(mass) {
+ Energy m = ids[1]->mass();
+ val+= 2.*sqr(m)/t;
+ }
+ return val;
+}
+
+double OneHalfHalfDarkSplitFn::integOverP(const double z, const IdList & ,
+ unsigned int PDFfactor) const {
+ switch(PDFfactor) {
+ case 0:
+ return z;
+ case 1:
+ return log(z);
+ case 2:
+ return -log(1.-z);
+ case 3:
+ return log(z/(1.-z));
+ case 4:
+ return 2.*sqrt(z);
+ case 5:
+ return (2./3.)*z*sqrt(z);
+ default:
+ throw Exception() << "OneHalfHalfDarkSplitFn::integOverP() invalid PDFfactor = "
+ << PDFfactor << Exception::runerror;
+ }
+}
+
+double OneHalfHalfDarkSplitFn::invIntegOverP(const double r,
+ const IdList & ,
+ unsigned int PDFfactor) const {
+ switch(PDFfactor) {
+ case 0:
+ return r;
+ case 1:
+ return exp(r);
+ case 2:
+ return 1.-exp(-r);
+ case 3:
+ return 1./(1.+exp(-r));
+ case 4:
+ return 0.25*sqr(r);
+ case 5:
+ return pow(1.5*r,2./3.);
+ default:
+ throw Exception() << "OneHalfHalfDarkSplitFn::integOverP() invalid PDFfactor = "
+ << PDFfactor << Exception::runerror;
+ }
+}
+
+bool OneHalfHalfDarkSplitFn::accept(const IdList &ids) const {
+ if(ids.size()!=3) return false;
+ if(ids[1]!=ids[2]->CC()) return false;
+ if(ids[1]->iSpin()!=PDT::Spin1Half) return false;
+ if(ids[0]->iSpin()!=PDT::Spin1) return false;
+ return OneHalfHalfDarkSplitFn::checkColours(ids);
+}
+
+vector<pair<int, Complex> >
+OneHalfHalfDarkSplitFn::generatePhiForward(const double z, const Energy2 t, const IdList & ids,
+ const RhoDMatrix & rho) {
+ assert(rho.iSpin()==PDT::Spin1);
+ double modRho = abs(rho(0,2));
+ Energy mq = ids[1]->mass();
+ Energy2 mq2 = sqr(mq);
+ double fact = z*(1.-z)-mq2/t;
+ double max = 1.+2.*fact*(-1.+2.*modRho);
+ vector<pair<int, Complex> > output;
+ output.push_back(make_pair( 0,(rho(0,0)+rho(2,2))*(1.-2.*fact)/max));
+ output.push_back(make_pair(-2,2.*fact*rho(0,2)/max));
+ output.push_back(make_pair( 2,2.*fact*rho(2,0)/max));
+ return output;
+}
+
+vector<pair<int, Complex> >
+OneHalfHalfDarkSplitFn::generatePhiBackward(const double, const Energy2, const IdList &,
+ const RhoDMatrix & ) {
+ // no dependance
+ return {{ {0, 1.} }};
+}
+
+DecayMEPtr OneHalfHalfDarkSplitFn::matrixElement(const double z, const Energy2 t,
+ const IdList & ids, const double phi,
+ bool timeLike) {
+ static const Complex ii(0.,1.);
+ // calculate the kernal
+ DecayMEPtr kernal(new_ptr(TwoBodyDecayMatrixElement(PDT::Spin1,PDT::Spin1Half,PDT::Spin1Half)));
+ double mt = !timeLike ? ZERO : ids[1]->mass()/sqrt(t);
+ double root =sqrt(1.-sqr(mt)/z/(1.-z));
+ (*kernal)(0,0,0) = mt/sqrt(z*(1.-z));
+ (*kernal)(2,1,1) = (*kernal)(0,0,0);
+ (*kernal)(0,0,1) = -z*root*exp(-ii*phi);
+ (*kernal)(2,1,0) = -conj((*kernal)(0,0,1));
+ (*kernal)(0,1,0) = (1.-z)*exp(-ii*phi)*root;
+ (*kernal)(2,0,1) = -conj((*kernal)(0,1,0));
+ (*kernal)(0,1,1) = 0.;
+ (*kernal)(2,0,0) = 0.;
+ return kernal;
+}
diff --git a/Shower/QTilde/SplittingFunctions/Dark/OneHalfHalfDarkSplitFn.h b/Shower/QTilde/SplittingFunctions/Dark/OneHalfHalfDarkSplitFn.h
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/SplittingFunctions/Dark/OneHalfHalfDarkSplitFn.h
@@ -0,0 +1,194 @@
+// -*- C++ -*-
+//
+// OneHalfHalfDarkSplitFn.h is a part of Herwig - A multi-purpose Monte Carlo event generator
+// Copyright (C) 2002-2019 The Herwig Collaboration
+//
+// Herwig is licenced under version 3 of the GPL, see COPYING for details.
+// Please respect the MCnet academic guidelines, see GUIDELINES for details.
+//
+#ifndef HERWIG_OneHalfHalfDarkSplitFn_H
+#define HERWIG_OneHalfHalfDarkSplitFn_H
+//
+// This is the declaration of the OneHalfHalfDarkSplitFn class.
+//
+
+#include "Herwig/Shower/QTilde/SplittingFunctions/Sudakov1to2FormFactor.h"
+#include "Herwig/Decay/TwoBodyDecayMatrixElement.h"
+
+namespace Herwig {
+
+using namespace ThePEG;
+
+/**\ingroup Shower
+ *
+ * This class provides the concrete implementation of the exact leading-order
+ * splitting function for \f$1\to \frac12\frac12\f$.
+ *
+ * In this case the splitting function is given by
+ * \f[P(z,t) =C\left(1-2z(1-z)+2\frac{m_q^2}{t}\right),\f]
+ * where \f$C\f$ is the corresponding colour factor
+ * Our choice for the overestimate is
+ * \f[P_{\rm over}(z) = C,\f]
+ * therefore the integral is
+ * \f[\int P_{\rm over}(z) {\rm d}z =Cz,\f]
+ * and its inverse is
+ * \f[\frac{r}{C}\f]
+ *
+ * @see \ref OneHalfHalfDarkSplitFnInterfaces "The interfaces"
+ * defined for OneHalfHalfDarkSplitFn.
+ */
+class OneHalfHalfDarkSplitFn: public Sudakov1to2FormFactor {
+
+public:
+
+ /**
+ * Method to check the colours are correct
+ */
+ virtual bool checkColours(const IdList & ids) const;
+
+ /**
+ * Concrete implementation of the method to determine whether this splitting
+ * function can be used for a given set of particles.
+ * @param ids The PDG codes for the particles in the splitting.
+ */
+ virtual bool accept(const IdList & ids) const;
+
+ /**
+ * Methods to return the splitting function.
+ */
+ //@{
+ /**
+ * The concrete implementation of the splitting function, \f$P\f$.
+ * @param z The energy fraction.
+ * @param t The scale.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param mass Whether or not to include the mass dependent terms
+ * @param rho The spin density matrix
+ */
+ virtual double P(const double z, const Energy2 t, const IdList & ids,
+ const bool mass, const RhoDMatrix & rho) const;
+
+
+ /**
+ * The concrete implementation of the overestimate of the splitting function,
+ * \f$P_{\rm over}\f$.
+ * @param z The energy fraction.
+ * @param ids The PDG codes for the particles in the splitting.
+ */
+ virtual double overestimateP(const double z, const IdList & ids) const;
+
+ /**
+ * The concrete implementation of the
+ * the ratio of the splitting function to the overestimate, i.e.
+ * \f$P(z,\tilde{q}^2)/P_{\rm over}(z)\f$.
+ * @param z The energy fraction.
+ * @param t The scale.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param mass Whether or not to include the mass dependent terms
+ * @param rho The spin density matrix
+ */
+ virtual double ratioP(const double z, const Energy2 t, const IdList & ids,
+ const bool mass, const RhoDMatrix & rho) const;
+
+ /**
+ * The concrete implementation of the indefinite integral of the
+ * overestimated splitting function, \f$P_{\rm over}\f$.
+ * @param z The energy fraction.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param PDFfactor Which additional factor to include for the PDF
+ * 0 is no additional factor,
+ * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
+ */
+ virtual double integOverP(const double z, const IdList & ids,
+ unsigned int PDFfactor=0) const;
+
+ /**
+ * The concrete implementation of the inverse of the indefinite integral.
+ * @param r Value of the splitting function to be inverted
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param PDFfactor Which additional factor to include for the PDF
+ * 0 is no additional factor,
+ * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
+ */
+ virtual double invIntegOverP(const double r, const IdList & ids,
+ unsigned int PDFfactor=0) const;
+ //@}
+
+ /**
+ * Method to calculate the azimuthal angle
+ * @param particle The particle which is branching
+ * @param showerkin The ShowerKinematics object
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ * @return The weight
+ */
+ virtual vector<pair<int,Complex> >
+ generatePhiForward(const double z, const Energy2 t, const IdList & ids,
+ const RhoDMatrix &);
+
+ /**
+ * Method to calculate the azimuthal angle
+ * @param particle The particle which is branching
+ * @param showerkin The ShowerKinematics object
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ * @return The weight
+ */
+ virtual vector<pair<int,Complex> >
+ generatePhiBackward(const double z, const Energy2 t, const IdList & ids,
+ const RhoDMatrix &);
+
+ /**
+ * Calculate the matrix element for the splitting
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ */
+ virtual DecayMEPtr matrixElement(const double z, const Energy2 t,
+ const IdList & ids, const double phi, bool timeLike);
+
+public:
+
+ /**
+ * The standard Init function used to initialize the interfaces.
+ * Called exactly once for each class by the class description system
+ * before the main function starts or
+ * when this class is dynamically loaded.
+ */
+ static void Init();
+
+protected:
+
+ /** @name Clone Methods. */
+ //@{
+ /**
+ * Make a simple clone of this object.
+ * @return a pointer to the new object.
+ */
+ virtual IBPtr clone() const {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);}
+ //@}
+
+private:
+
+ /**
+ * The assignment operator is private and must never be called.
+ * In fact, it should not even be implemented.
+ */
+ OneHalfHalfDarkSplitFn & operator=(const OneHalfHalfDarkSplitFn &) = delete;
+
+};
+
+}
+
+#endif /* HERWIG_OneHalfHalfDarkSplitFn_H */
diff --git a/Shower/QTilde/SplittingFunctions/Dark/OneOneOneDarkSplitFn.cc b/Shower/QTilde/SplittingFunctions/Dark/OneOneOneDarkSplitFn.cc
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/SplittingFunctions/Dark/OneOneOneDarkSplitFn.cc
@@ -0,0 +1,167 @@
+// -*- C++ -*-
+//
+// OneOneOneDarkSplitFn.cc is a part of Herwig - A multi-purpose Monte Carlo event generator
+// Copyright (C) 2002-2019 The Herwig Collaboration
+//
+// Herwig is licenced under version 3 of the GPL, see COPYING for details.
+// Please respect the MCnet academic guidelines, see GUIDELINES for details.
+//
+//
+// This is the implementation of the non-inlined, non-templated member
+// functions of the OneOneOneDarkSplitFn class.
+//
+
+#include "OneOneOneDarkSplitFn.h"
+#include "Herwig/Models/HiddenValley/HiddenValleyModel.h"
+#include "ThePEG/PDT/ParticleData.h"
+#include "ThePEG/Interface/ClassDocumentation.h"
+#include "ThePEG/Utilities/DescribeClass.h"
+
+using namespace Herwig;
+
+DescribeNoPIOClass<OneOneOneDarkSplitFn,Herwig::Sudakov1to2FormFactor>
+describeOneOneOneDarkSplitFn ("Herwig::OneOneOneDarkSplitFn","HwDarkShower.so");
+
+void OneOneOneDarkSplitFn::Init() {
+
+ static ClassDocumentation<OneOneOneDarkSplitFn> documentation
+ ("The OneOneOneDarkSplitFn class implements the g -> gg splitting function");
+
+}
+
+bool OneOneOneDarkSplitFn::checkColours(const IdList & ids) const {
+ if(colourStructure()==TripletTripletOctet) {
+ if(ids[0]!=ids[1]) return false;
+ if((ids[0]->iColour()==PDT::DarkColourFundamental||
+ ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
+ ids[2]->iColour()==PDT::DarkColourAdjoint) return true;
+ return false;
+ }
+ else if(colourStructure()==OctetOctetOctet) {
+ for(unsigned int ix=0;ix<3;++ix) {
+ if(ids[ix]->iColour()!=PDT::DarkColourAdjoint) return false;
+ }
+ return true;
+ }
+ else if(colourStructure()==OctetTripletTriplet) {
+ if(ids[0]->iColour()!=PDT::DarkColourAdjoint) return false;
+ if(ids[1]->iColour()==PDT::DarkColourFundamental&&
+ ids[2]->iColour()==PDT::DarkColourAntiFundamental)
+ return true;
+ if(ids[1]->iColour()==PDT::DarkColourAntiFundamental&&
+ ids[2]->iColour()==PDT::DarkColourFundamental)
+ return true;
+ return false;
+ }
+ else if(colourStructure()==TripletOctetTriplet) {
+ if(ids[0]!=ids[2]) return false;
+ if((ids[0]->iColour()==PDT::DarkColourFundamental||
+ ids[0]->iColour()==PDT::DarkColourAntiFundamental) &&
+ ids[1]->iColour()==PDT::DarkColourAdjoint) return true;
+ return false;
+ }
+ else {
+ assert(false);
+ }
+ return false;
+}
+
+double OneOneOneDarkSplitFn::P(const double z, const Energy2,
+ const IdList & , const bool, const RhoDMatrix &)const {
+ return sqr(1.-z*(1.-z))/(z*(1.-z));
+}
+
+double OneOneOneDarkSplitFn::overestimateP(const double z,
+ const IdList & ) const {
+ return (1/z + 1/(1.-z));
+}
+
+
+double OneOneOneDarkSplitFn::ratioP(const double z, const Energy2,
+ const IdList & , const bool, const RhoDMatrix &) const {
+ return sqr(1.-z*(1.-z));
+}
+
+double OneOneOneDarkSplitFn::invIntegOverP(const double r,
+ const IdList & ,
+ unsigned int PDFfactor) const {
+ switch(PDFfactor) {
+ case 0:
+ return 1./(1.+exp(-r));
+ case 1:
+ case 2:
+ case 3:
+ default:
+ throw Exception() << "OneOneOneDarkSplitFn::integOverP() invalid PDFfactor = "
+ << PDFfactor << Exception::runerror;
+ }
+}
+
+double OneOneOneDarkSplitFn::integOverP(const double z, const IdList & ,
+ unsigned int PDFfactor) const {
+ switch(PDFfactor) {
+ case 0:
+ assert(z>0.&&z<1.);
+ return log(z/(1.-z));
+ case 1:
+ case 2:
+ case 3:
+ default:
+ throw Exception() << "OneOneOneDarkSplitFn::integOverP() invalid PDFfactor = "
+ << PDFfactor << Exception::runerror;
+ }
+}
+
+bool OneOneOneDarkSplitFn::accept(const IdList & ids) const {
+ if(ids.size()!=3) return false;
+ for(unsigned int ix=0;ix<ids.size();++ix) {
+ if(ids[0]->iSpin()!=PDT::Spin1) return false;
+ }
+ return OneOneOneDarkSplitFn::checkColours(ids);
+}
+
+vector<pair<int, Complex> >
+OneOneOneDarkSplitFn::generatePhiForward(const double z, const Energy2, const IdList &,
+ const RhoDMatrix & rho) {
+ assert(rho.iSpin()==PDT::Spin1);
+ double modRho = abs(rho(0,2));
+ double max = 2.*z*modRho*(1.-z)+sqr(1.-(1.-z)*z)/(z*(1.-z));
+ vector<pair<int, Complex> > output;
+ output.push_back(make_pair( 0,(rho(0,0)+rho(2,2))*sqr(1.-(1.-z)*z)/(z*(1.-z))/max));
+ output.push_back(make_pair(-2,-rho(0,2)*z*(1.-z)/max));
+ output.push_back(make_pair( 2,-rho(2,0)*z*(1.-z)/max));
+ return output;
+}
+
+vector<pair<int, Complex> >
+OneOneOneDarkSplitFn::generatePhiBackward(const double z, const Energy2, const IdList &,
+ const RhoDMatrix & rho) {
+ assert(rho.iSpin()==PDT::Spin1);
+ double diag = sqr(1 - (1 - z)*z)/(1 - z)/z;
+ double off = (1.-z)/z;
+ double max = 2.*abs(rho(0,2))*off+diag;
+ vector<pair<int, Complex> > output;
+ output.push_back(make_pair( 0, (rho(0,0)+rho(2,2))*diag/max));
+ output.push_back(make_pair( 2,-rho(0,2) * off/max));
+ output.push_back(make_pair(-2,-rho(2,0) * off/max));
+ return output;
+}
+
+DecayMEPtr OneOneOneDarkSplitFn::matrixElement(const double z, const Energy2,
+ const IdList &, const double phi,
+ bool) {
+ // calculate the kernal
+ DecayMEPtr kernal(new_ptr(TwoBodyDecayMatrixElement(PDT::Spin1,PDT::Spin1,PDT::Spin1)));
+ double omz = 1.-z;
+ double root = sqrt(z*omz);
+ Complex phase = exp(Complex(0.,1.)*phi);
+ (*kernal)(0,0,0) = phase/root;
+ (*kernal)(2,2,2) = -conj((*kernal)(0,0,0));
+ (*kernal)(0,0,2) = -sqr(z)/root/phase;
+ (*kernal)(2,2,0) = -conj((*kernal)(0,0,2));
+ (*kernal)(0,2,0) = -sqr(omz)/root/phase;
+ (*kernal)(2,0,2) = -conj((*kernal)(0,2,0));
+ (*kernal)(0,2,2) = 0.;
+ (*kernal)(2,0,0) = 0.;
+ return kernal;
+}
diff --git a/Shower/QTilde/SplittingFunctions/Dark/OneOneOneDarkSplitFn.h b/Shower/QTilde/SplittingFunctions/Dark/OneOneOneDarkSplitFn.h
new file mode 100644
--- /dev/null
+++ b/Shower/QTilde/SplittingFunctions/Dark/OneOneOneDarkSplitFn.h
@@ -0,0 +1,190 @@
+// -*- C++ -*-
+//
+// OneOneOneDarkSplitFn.h is a part of Herwig - A multi-purpose Monte Carlo event generator
+// Copyright (C) 2002-2019 The Herwig Collaboration
+//
+// Herwig is licenced under version 3 of the GPL, see COPYING for details.
+// Please respect the MCnet academic guidelines, see GUIDELINES for details.
+//
+#ifndef HERWIG_OneOneOneDarkSplitFn_H
+#define HERWIG_OneOneOneDarkSplitFn_H
+//
+// This is the declaration of the OneOneOneDarkSplitFn class.
+//
+
+#include "Herwig/Shower/QTilde/SplittingFunctions/Sudakov1to2FormFactor.h"
+#include "Herwig/Decay/TwoBodyDecayMatrixElement.h"
+
+namespace Herwig {
+
+using namespace ThePEG;
+
+/** \ingroup Shower
+ *
+ * This class provides the concrete implementation of the exact leading-order
+ * splitting function for \f$1\to 11\f$.
+ *
+ * In this case the splitting function is given by
+ * \f[P(z) =2C*\left(\frac{z}{1-z}+\frac{1-z}{z}+z(1-z)\right),\f]
+ * where \f$C=\f$ is the corresponding colour factor.
+ * Our choice for the overestimate is
+ * \f[P_{\rm over}(z) = 2C\left(\frac1z+\frac1{1-z}\right),\f]
+ * therefore the integral is
+ * \f[\int P_{\rm over}(z) {\rm d}z =2C\ln\left(\frac{z}{1-z}\right),\f]
+ * and its inverse is
+ * \f[\frac{\exp\left(\frac{r}{2C}\right)}{(1+\exp\left(\frac{r}{2C}\right)}\f]
+ *
+ *
+ * @see \ref OneOneOneDarkSplitFnInterfaces "The interfaces"
+ * defined for OneOneOneDarkSplitFn.
+ */
+class OneOneOneDarkSplitFn: public Sudakov1to2FormFactor {
+
+public:
+
+ /**
+ * Method to check the colours are correct
+ */
+ virtual bool checkColours(const IdList & ids) const;
+
+ /**
+ * Concrete implementation of the method to determine whether this splitting
+ * function can be used for a given set of particles.
+ * @param ids The PDG codes for the particles in the splitting.
+ */
+ virtual bool accept(const IdList & ids) const;
+
+ /**
+ * Methods to return the splitting function.
+ */
+ //@{
+ /**
+ * The concrete implementation of the splitting function, \f$P(z,t)\f$.
+ * @param z The energy fraction.
+ * @param t The scale.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param mass Whether or not to include the mass dependent terms
+ * @param rho The spin density matrix
+ */
+ virtual double P(const double z, const Energy2 t, const IdList & ids,
+ const bool mass, const RhoDMatrix & rho) const;
+
+ /**
+ * The concrete implementation of the overestimate of the splitting function,
+ * \f$P_{\rm over}\f$.
+ * @param z The energy fraction.
+ * @param ids The PDG codes for the particles in the splitting.
+ */
+ virtual double overestimateP(const double z, const IdList & ids) const;
+
+ /**
+ * The concrete implementation of the
+ * the ratio of the splitting function to the overestimate, i.e.
+ * \f$P(z,t)/P_{\rm over}(z)\f$.
+ * @param z The energy fraction.
+ * @param t The scale.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param mass Whether or not to include the mass dependent terms
+ * @param rho The spin density matrix
+ */
+ virtual double ratioP(const double z, const Energy2 t, const IdList & ids,
+ const bool mass, const RhoDMatrix & rho) const;
+
+ /**
+ * The concrete implementation of the indefinite integral of the
+ * overestimated splitting function, \f$P_{\rm over}\f$.
+ * @param z The energy fraction.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param PDFfactor Which additional factor to include for the PDF
+ * 0 is no additional factor,
+ * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
+ */
+ virtual double integOverP(const double z, const IdList & ids,
+ unsigned int PDFfactor=0) const;
+
+ /**
+ * The concrete implementation of the inverse of the indefinite integral.
+ * @param r Value of the splitting function to be inverted
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param PDFfactor Which additional factor to include for the PDF
+ * 0 is no additional factor,
+ * 1 is \f$1/z\f$, 2 is \f$1/(1-z)\f$ and 3 is \f$1/z/(1-z)\f$
+ */
+ virtual double invIntegOverP(const double r, const IdList & ids,
+ unsigned int PDFfactor=0) const;
+ //@}
+
+ /**
+ * Method to calculate the azimuthal angle for forward evolution
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ * @return The weight
+ */
+ virtual vector<pair<int,Complex> >
+ generatePhiForward(const double z, const Energy2 t, const IdList & ids,
+ const RhoDMatrix &);
+
+ /**
+ * Method to calculate the azimuthal angle for backward evolution
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ * @return The weight
+ */
+ virtual vector<pair<int,Complex> >
+ generatePhiBackward(const double z, const Energy2 t, const IdList & ids,
+ const RhoDMatrix &);
+
+ /**
+ * Calculate the matrix element for the splitting
+ * @param z The energy fraction
+ * @param t The scale \f$t=2p_j\cdot p_k\f$.
+ * @param ids The PDG codes for the particles in the splitting.
+ * @param The azimuthal angle, \f$\phi\f$.
+ */
+ virtual DecayMEPtr matrixElement(const double z, const Energy2 t,
+ const IdList & ids, const double phi, bool timeLike);
+
+public:
+
+ /**
+ * The standard Init function used to initialize the interfaces.
+ * Called exactly once for each class by the class description system
+ * before the main function starts or
+ * when this class is dynamically loaded.
+ */
+ static void Init();
+
+protected:
+
+ /** @name Clone Methods. */
+ //@{
+ /**
+ * Make a simple clone of this object.
+ * @return a pointer to the new object.
+ */
+ virtual IBPtr clone() const {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);}
+ //@}
+
+private:
+
+ /**
+ * The assignment operator is private and must never be called.
+ * In fact, it should not even be implemented.
+ */
+ OneOneOneDarkSplitFn & operator=(const OneOneOneDarkSplitFn &) = delete;
+
+};
+
+}
+
+#endif /* HERWIG_OneOneOneDarkSplitFn_H */
diff --git a/configure.ac b/configure.ac
--- a/configure.ac
+++ b/configure.ac
@@ -1,273 +1,275 @@
dnl Process this file with autoconf to produce a configure script.
AC_PREREQ([2.71])
AC_INIT([Herwig],[devel],[herwig@projects.hepforge.org],[Herwig])
AC_CONFIG_SRCDIR([Utilities/HerwigStrategy.cc])
AC_CONFIG_AUX_DIR([Config])
AC_CONFIG_MACRO_DIR([m4])
AC_CONFIG_HEADERS([Config/config.h])
dnl AC_PRESERVE_HELP_ORDER
AC_CANONICAL_HOST
dnl === disable debug symbols by default =====
if test "x$CXXFLAGS" = "x"; then
CXXFLAGS="-O2 -DBOOST_UBLAS_NDEBUG -Wno-deprecated-declarations -Wno-deprecated-copy"
fi
if test "x$CFLAGS" = "x"; then
CFLAGS=-O2
fi
AC_LANG([C++])
AM_INIT_AUTOMAKE([1.11 subdir-objects gnu dist-bzip2 no-dist-gzip -Wall -Wno-portability])
m4_ifdef([AM_SILENT_RULES], [AM_SILENT_RULES([yes])])
m4_ifdef([AM_PROG_AR], [AM_PROG_AR])
dnl Checks for C++ compiler. Handle C++11 flags.
AC_PROG_CXX
AX_CXX_COMPILE_STDCXX([11],[noext],[mandatory])
dnl check for POSIX
AC_CHECK_HEADER([unistd.h],[],
[AC_MSG_ERROR([Herwig needs "unistd.h". Non-POSIX systems are not supported.])])
AC_CHECK_HEADER([sys/stat.h],[],
[AC_MSG_ERROR([Herwig needs "sys/stat.h". Non-POSIX systems are not supported.])])
dnl Checks for programs.
AC_PROG_INSTALL
AC_PROG_MAKE_SET
AC_PROG_LN_S
dnl modified search order
AC_PROG_FC([gfortran g95 g77])
dnl xlf95 f95 fort ifort ifc efc pgf95 lf95 ftn xlf90 f90 pgf90 pghpf epcf90 xlf f77 frt pgf77 cf77 fort77 fl32 af77])
AC_LANG_PUSH([Fortran])
AC_MSG_CHECKING([if the Fortran compiler ($FC) works])
AC_COMPILE_IFELSE(
AC_LANG_PROGRAM([],[ print *[,]"Hello"]),
[AC_MSG_RESULT([yes])],
[AC_MSG_RESULT([no])
AC_MSG_ERROR([A Fortran compiler is required to build Herwig.])
]
)
AC_LANG_POP([Fortran])
AC_FC_WRAPPERS
LT_PREREQ([2.2.6])
LT_INIT([disable-static dlopen pic-only])
dnl ####################################
dnl ####################################
dnl for Doc/fixinterfaces.pl
AC_PATH_PROG(PERL, perl)
dnl for Models/Feynrules
AM_PATH_PYTHON([2.6],, [:])
AM_CONDITIONAL([HAVE_PYTHON], [test "x$PYTHON" != "x:"])
HERWIG_CHECK_GSL
HERWIG_CHECK_THEPEG
BOOST_REQUIRE([1.41])
BOOST_FIND_HEADER([boost/numeric/ublas/io.hpp])
dnl Boost 1.64 is missing a required header to make these work
dnl we just assume they're there if io.hpp has been found OK above
dnl BOOST_FIND_HEADER([boost/numeric/ublas/matrix.hpp])
dnl BOOST_FIND_HEADER([boost/numeric/ublas/matrix_proxy.hpp])
dnl BOOST_FIND_HEADER([boost/numeric/ublas/matrix_sparse.hpp])
dnl BOOST_FIND_HEADER([boost/numeric/ublas/symmetric.hpp])
dnl BOOST_FIND_HEADER([boost/numeric/ublas/vector.hpp])
BOOST_FIND_HEADER([boost/operators.hpp])
BOOST_TEST()
HERWIG_CHECK_VBFNLO
HERWIG_CHECK_NJET
HERWIG_CHECK_GOSAM
HERWIG_CHECK_GOSAM_CONTRIB
HERWIG_CHECK_OPENLOOPS
HERWIG_CHECK_MADGRAPH
HERWIG_CHECK_EVTGEN
HERWIG_CHECK_PYTHIA
HERWIG_COMPILERFLAGS
HERWIG_LOOPTOOLS
FASTJET_CHECK_FASTJET
HERWIG_ENABLE_MODELS
SHARED_FLAG=-shared
AM_CONDITIONAL(NEED_APPLE_FIXES,
[test "xx${host/darwin/foundit}xx" != "xx${host}xx"])
if test "xx${host/darwin/foundit}xx" != "xx${host}xx"; then
APPLE_DSO_FLAGS=-Wl,-undefined,dynamic_lookup
SHARED_FLAG=-bundle
fi
AC_SUBST([APPLE_DSO_FLAGS])
AC_SUBST([SHARED_FLAG])
AC_SUBST([PYTHON])
AC_CONFIG_FILES([UnderlyingEvent/Makefile
Models/Makefile
Models/StandardModel/Makefile
Models/RSModel/Makefile
Models/General/Makefile
Models/Susy/Makefile
Models/Susy/NMSSM/Makefile
Models/Susy/RPV/Makefile
Models/UED/Makefile
Models/LH/Makefile
Models/DarkMatter/Makefile
Models/LHTP/Makefile
Models/Transplanckian/Makefile
Models/Leptoquarks/Makefile
Models/Zprime/Makefile
Models/TTbAsymm/Makefile
Models/Feynrules/Makefile
Models/Feynrules/python/Makefile-FR
Models/ADD/Makefile
Models/Sextet/Makefile
+ Models/HiddenValley/Makefile
Decay/Makefile
Decay/FormFactors/Makefile
Decay/Tau/Makefile
Decay/Baryon/Makefile
Decay/VectorMeson/Makefile
Decay/Perturbative/Makefile
Decay/HeavyMeson/Makefile
Decay/ScalarMeson/Makefile
Decay/Dalitz/Makefile
Decay/TensorMeson/Makefile
Decay/WeakCurrents/Makefile
Decay/Partonic/Makefile
Decay/General/Makefile
Decay/Radiation/Makefile
Decay/EvtGen/Makefile
Doc/refman.conf
Doc/refman.h
PDT/Makefile
PDF/Makefile
MatrixElement/Makefile
MatrixElement/General/Makefile
MatrixElement/Lepton/Makefile
MatrixElement/Hadron/Makefile
MatrixElement/DIS/Makefile
MatrixElement/Powheg/Makefile
MatrixElement/Gamma/Makefile
MatrixElement/Reweighters/Makefile
MatrixElement/Onium/Makefile
MatrixElement/Matchbox/Makefile
MatrixElement/Matchbox/Base/Makefile
MatrixElement/Matchbox/Utility/Makefile
MatrixElement/Matchbox/Phasespace/Makefile
MatrixElement/Matchbox/Dipoles/Makefile
MatrixElement/Matchbox/InsertionOperators/Makefile
MatrixElement/Matchbox/Matching/Makefile
MatrixElement/Matchbox/Cuts/Makefile
MatrixElement/Matchbox/Scales/Makefile
MatrixElement/Matchbox/ColorFull/Makefile
MatrixElement/Matchbox/CVolver/Makefile
MatrixElement/Matchbox/Builtin/Makefile
MatrixElement/Matchbox/Builtin/Amplitudes/Makefile
MatrixElement/Matchbox/Tests/Makefile
MatrixElement/Matchbox/External/Makefile
MatrixElement/Matchbox/External/BLHAGeneric/Makefile
MatrixElement/Matchbox/External/VBFNLO/Makefile
MatrixElement/Matchbox/External/NJet/Makefile
MatrixElement/Matchbox/External/GoSam/Makefile
MatrixElement/Matchbox/External/OpenLoops/Makefile
MatrixElement/Matchbox/External/MadGraph/Makefile
MatrixElement/Matchbox/External/MadGraph/mg2herwig
MatrixElement/Matchbox/External/GoSam/gosam2herwig
MatrixElement/FxFx/Makefile
Sampling/Makefile
Sampling/CellGrids/Makefile
Shower/Makefile
Shower/QTilde/Makefile
Shower/QTilde/Matching/Makefile
Shower/QTilde/SplittingFunctions/Onium/Makefile
+ Shower/QTilde/SplittingFunctions/Dark/Makefile
Shower/Dipole/Makefile
Shower/Dipole/Base/Makefile
Shower/Dipole/Kernels/Makefile
Shower/Dipole/Kinematics/Makefile
Shower/Dipole/Utility/Makefile
Shower/Dipole/AlphaS/Makefile
Shower/Dipole/SpinCorrelations/Makefile
Utilities/Makefile
Utilities/XML/Makefile
Utilities/Statistics/Makefile
Hadronization/Makefile
lib/Makefile
include/Makefile
src/Makefile
src/defaults/Makefile
src/snippets/Makefile
src/Matchbox/Makefile
src/herwig-config
Doc/Makefile
Doc/HerwigDefaults.in
Looptools/Makefile
Analysis/Makefile
API/Makefile
src/Makefile-UserModules
src/defaults/Analysis.in
src/defaults/MatchboxDefaults.in
src/defaults/Decays.in
src/defaults/decayers.in
src/defaults/setup.gosam.in
src/Matchbox/LO-DefaultShower.in
src/Matchbox/LO-DipoleShower.in
src/Matchbox/MCatLO-DefaultShower.in
src/Matchbox/MCatLO-DipoleShower.in
src/Matchbox/LO-NoShower.in
src/Matchbox/MCatNLO-DefaultShower.in
src/Matchbox/MCatNLO-DipoleShower.in
src/Matchbox/NLO-NoShower.in
src/Matchbox/Powheg-DefaultShower.in
src/Matchbox/Powheg-DipoleShower.in
src/Merging/Makefile
Shower/Dipole/Merging/Makefile
src/defaults/MatchboxMergingDefaults.in
Contrib/Makefile
Contrib/make_makefiles.sh
Tests/Makefile
Makefile])
AC_CONFIG_FILES([Tests/python/rivet_check ],[chmod +x Tests/python/rivet_check ])
AC_CONFIG_FILES([Tests/python/LowEnergy-EE.py ],[chmod +x Tests/python/LowEnergy-EE.py ])
AC_CONFIG_FILES([Tests/python/LowEnergy-Photon.py],[chmod +x Tests/python/LowEnergy-Photon.py])
AC_CONFIG_FILES([Tests/python/OniumSigma.py ],[chmod +x Tests/python/OniumSigma.py ])
AC_CONFIG_FILES([Tests/python/OniumSplitting.py ],[chmod +x Tests/python/OniumSplitting.py ])
AC_CONFIG_FILES([Tests/python/make_input_files.py],[chmod +x Tests/python/make_input_files.py])
AC_CONFIG_FILES([Tests/python/merge-DIS ],[chmod +x Tests/python/merge-DIS ])
AC_CONFIG_FILES([Tests/python/merge-EE ],[chmod +x Tests/python/merge-EE ])
AC_CONFIG_FILES([Tests/python/merge-EE-Gamma ],[chmod +x Tests/python/merge-EE-Gamma ])
AC_CONFIG_FILES([Tests/python/merge-Fixed ],[chmod +x Tests/python/merge-Fixed ])
AC_CONFIG_FILES([Tests/python/merge-GammaGamma ],[chmod +x Tests/python/merge-GammaGamma ])
AC_CONFIG_FILES([Tests/python/merge-LHC-EW ],[chmod +x Tests/python/merge-LHC-EW ])
AC_CONFIG_FILES([Tests/python/merge-LHC-Jets ],[chmod +x Tests/python/merge-LHC-Jets ])
AC_CONFIG_FILES([Tests/python/merge-LHC-Photon ],[chmod +x Tests/python/merge-LHC-Photon ])
AC_CONFIG_FILES([Tests/python/mergeLowEnergy.py ],[chmod +x Tests/python/mergeLowEnergy.py ])
AC_CONFIG_FILES([Tests/python/merge-TVT-EW ],[chmod +x Tests/python/merge-TVT-EW ])
AC_CONFIG_FILES([Tests/python/merge-TVT-Jets ],[chmod +x Tests/python/merge-TVT-Jets ])
AC_CONFIG_FILES([Tests/python/merge-TVT-Photon ],[chmod +x Tests/python/merge-TVT-Photon ])
AC_CONFIG_FILES([Tests/python/plot-EE ],[chmod +x Tests/python/plot-EE ])
AC_CONFIG_FILES([Tests/python/R.py ],[chmod +x Tests/python/R.py ])
AC_CONFIG_FILES([Sampling/herwig-mergegrids ],[chmod +x Sampling/herwig-mergegrids ])
AC_CONFIG_LINKS([Doc/BSMlibs.in:Doc/BSMlibs.in])
AC_CONFIG_FILES([Doc/fixinterfaces.pl],[chmod +x Doc/fixinterfaces.pl])
AC_CONFIG_HEADERS([PDF/SaSPhotonPDF.cc])
HERWIG_OVERVIEW
AC_CONFIG_COMMANDS([summary],[cat config.herwig])
AC_OUTPUT
diff --git a/lib/Makefile.am b/lib/Makefile.am
--- a/lib/Makefile.am
+++ b/lib/Makefile.am
@@ -1,50 +1,51 @@
pkglib_LTLIBRARIES = Herwig.la
Herwig_la_SOURCES =
Herwig_la_LIBTOOLFLAGS = --tag=CXX
Herwig_la_LDFLAGS = $(AM_LDFLAGS) -module -version-info 27:0:0
Herwig_la_LDFLAGS += $(THEPEGLDFLAGS) $(BOOST_SYSTEM_LDFLAGS) $(BOOST_FILESYSTEM_LDFLAGS) $(FCLIBS)
Herwig_la_LIBADD = \
$(top_builddir)/Hadronization/libHwHadronization.la \
$(top_builddir)/Models/libHwStandardModel.la \
+$(top_builddir)/Models/libHwHiddenValley.la \
$(top_builddir)/Decay/libHwDecay.la \
$(top_builddir)/Decay/libHwFormFactor.la \
$(top_builddir)/Decay/libHwDecRad.la \
$(top_builddir)/Utilities/libHwUtils.la \
$(top_builddir)/Models/libHwModelGenerator.la \
$(top_builddir)/Decay/General/libHwGeneralDecay.la \
$(top_builddir)/MatrixElement/General/libHwGeneralME.la \
$(top_builddir)/MatrixElement/libHwME.la \
$(top_builddir)/MatrixElement/Reweighters/libHwReweighters.la \
$(top_builddir)/MatrixElement/Matchbox/libHwMatchbox.la \
$(top_builddir)/Decay/libHwWeakCurrent.la \
$(top_builddir)/Looptools/libHwLooptools.la \
$(top_builddir)/Shower/libHwShower.la \
$(THEPEGLIB) -ldl
dist_noinst_SCRIPTS = fix-osx-path
POSTPROCESSING = done-all-links
if NEED_APPLE_FIXES
POSTPROCESSING += apple-fixes
endif
all-local: $(POSTPROCESSING)
done-all-links: Herwig.la
find $(top_builddir) \( -name '*.so.*' -or -name '*.so' \) \
-not -name 'lib*' -not -path '$(top_builddir)/lib/*' \
-not -path '$(top_builddir)/Tests/*' \
-not -path '$(top_builddir)/.hg/*' -exec $(LN_S) -f \{\} \;
$(LN_S) -f .libs/Herwig*so* .
echo "stamp" > $@
apple-fixes: fix-osx-path done-all-links
./$<
echo "stamp" > $@
clean-local:
rm -f *.so *.so.* done-all-links apple-fixes
diff --git a/src/defaults/Shower.in b/src/defaults/Shower.in
--- a/src/defaults/Shower.in
+++ b/src/defaults/Shower.in
@@ -1,477 +1,445 @@
# -*- ThePEG-repository -*-
############################################################
# Setup of default parton shower
#
# Useful switches for users are marked near the top of
# this file.
#
# Don't edit this file directly, but reset the switches
# in your own input files!
############################################################
library HwMPI.so
library HwShower.so
library HwMatching.so
mkdir /Herwig/Shower
cd /Herwig/Shower
create Herwig::QTildeShowerHandler ShowerHandler
newdef ShowerHandler:MPIHandler /Herwig/UnderlyingEvent/MPIHandler
newdef ShowerHandler:RemDecayer /Herwig/Partons/RemnantDecayer
# use LO PDFs for Shower, can be changed later
newdef ShowerHandler:PDFA /Herwig/Partons/ShowerLOPDF
newdef ShowerHandler:PDFB /Herwig/Partons/ShowerLOPDF
newdef ShowerHandler:PDFARemnant /Herwig/Partons/RemnantPDF
newdef ShowerHandler:PDFBRemnant /Herwig/Partons/RemnantPDF
#####################################
# initial setup, don't change these!
#####################################
create Herwig::SplittingGenerator SplittingGenerator
create Herwig::ShowerAlphaQCD AlphaQCDISR
create Herwig::ShowerAlphaQCD AlphaQCDFSR
create Herwig::ShowerAlphaQED AlphaQED
set AlphaQED:CouplingSource Thompson
create Herwig::ShowerAlphaQED AlphaEW
set AlphaEW:CouplingSource MZ
-create Herwig::HiddenValleyAlpha AlphaDARK
create Herwig::PartnerFinder PartnerFinder
newdef PartnerFinder:PartnerMethod 0
newdef PartnerFinder:ScaleChoice 0
create Herwig::KinematicsReconstructor KinematicsReconstructor
newdef KinematicsReconstructor:ReconstructionOption Colour3
newdef KinematicsReconstructor:InitialStateReconOption SofterFraction
newdef KinematicsReconstructor:InitialInitialBoostOption LongTransBoost
newdef KinematicsReconstructor:FinalFinalWeight Yes
newdef /Herwig/Partons/RemnantDecayer:AlphaS AlphaQCDISR
newdef /Herwig/Partons/RemnantDecayer:AlphaEM AlphaQED
newdef ShowerHandler:PartnerFinder PartnerFinder
newdef ShowerHandler:KinematicsReconstructor KinematicsReconstructor
newdef ShowerHandler:SplittingGenerator SplittingGenerator
newdef ShowerHandler:Interactions ALL
newdef ShowerHandler:SpinCorrelations Yes
newdef ShowerHandler:SoftCorrelations Singular
##################################################################
# Intrinsic pT
#
# Recommended:
# 1.9 GeV for Tevatron W/Z production.
# 2.1 GeV for LHC W/Z production at 10 TeV
# 2.2 GeV for LHC W/Z production at 14 TeV
#
# Set all parameters to 0 to disable
##################################################################
newdef ShowerHandler:IntrinsicPtGaussian 2.2*GeV
newdef ShowerHandler:IntrinsicPtBeta 0
newdef ShowerHandler:IntrinsicPtGamma 0*GeV
newdef ShowerHandler:IntrinsicPtIptmax 0*GeV
#############################################################
# Set up truncated shower handler.
#############################################################
create Herwig::PowhegShowerHandler PowhegShowerHandler
set PowhegShowerHandler:MPIHandler /Herwig/UnderlyingEvent/MPIHandler
set PowhegShowerHandler:RemDecayer /Herwig/Partons/RemnantDecayer
newdef PowhegShowerHandler:PDFA /Herwig/Partons/ShowerLOPDF
newdef PowhegShowerHandler:PDFB /Herwig/Partons/ShowerLOPDF
newdef PowhegShowerHandler:PDFARemnant /Herwig/Partons/RemnantPDF
newdef PowhegShowerHandler:PDFBRemnant /Herwig/Partons/RemnantPDF
newdef PowhegShowerHandler:MPIHandler /Herwig/UnderlyingEvent/MPIHandler
newdef PowhegShowerHandler:RemDecayer /Herwig/Partons/RemnantDecayer
newdef PowhegShowerHandler:PDFA /Herwig/Partons/ShowerLOPDF
newdef PowhegShowerHandler:PDFB /Herwig/Partons/ShowerLOPDF
newdef PowhegShowerHandler:PDFARemnant /Herwig/Partons/RemnantPDF
newdef PowhegShowerHandler:PDFBRemnant /Herwig/Partons/RemnantPDF
newdef PowhegShowerHandler:PartnerFinder PartnerFinder
newdef PowhegShowerHandler:KinematicsReconstructor KinematicsReconstructor
newdef PowhegShowerHandler:SplittingGenerator SplittingGenerator
newdef PowhegShowerHandler:Interactions ALL
newdef PowhegShowerHandler:SpinCorrelations Yes
newdef PowhegShowerHandler:SoftCorrelations Singular
newdef PowhegShowerHandler:IntrinsicPtGaussian 2.2*GeV
newdef PowhegShowerHandler:IntrinsicPtBeta 0
newdef PowhegShowerHandler:IntrinsicPtGamma 0*GeV
newdef PowhegShowerHandler:IntrinsicPtIptmax 0*GeV
newdef PowhegShowerHandler:EvolutionScheme DotProduct
#############################################################
# End of interesting user servicable section.
#
# Anything that follows below should only be touched if you
# know what you're doing.
#
# Really.
#############################################################
#
# a few default values
newdef ShowerHandler:MECorrMode 1
newdef ShowerHandler:EvolutionScheme DotProduct
newdef AlphaQCDISR:ScaleFactor 1.0
newdef AlphaQCDISR:NPAlphaS 2
newdef AlphaQCDISR:Qmin 0.935
newdef AlphaQCDISR:NumberOfLoops 2
newdef AlphaQCDISR:AlphaIn 0.1185
newdef AlphaQCDFSR:ScaleFactor 1.0
newdef AlphaQCDFSR:NPAlphaS 2
newdef AlphaQCDFSR:Qmin 0.935
newdef AlphaQCDFSR:NumberOfLoops 2
newdef AlphaQCDFSR:AlphaIn 0.102337
-newdef AlphaDARK:ScaleFactor 1.0
-newdef AlphaDARK:NPAlphaS 2
-newdef AlphaDARK:Qmin 0.935
-newdef AlphaDARK:NumberOfLoops 2
#
#
#
# Lets set up all the splittings and Sudakovs
# cut-off
create Herwig::PTCutOff PTCutOff
newdef PTCutOff:pTmin 0.654714*GeV
# q -> q g
create Herwig::HalfHalfOneSplitFn QtoQGSudakov
newdef QtoQGSudakov:InteractionType QCD
newdef QtoQGSudakov:Alpha AlphaQCDISR
newdef QtoQGSudakov:Cutoff PTCutOff
newdef QtoQGSudakov:ColourStructure TripletTripletOctet
newdef QtoQGSudakov:AngularOrdered Yes
newdef QtoQGSudakov:StrictAO Yes
newdef QtoQGSudakov:PDFmax 1.9
-create Herwig::HalfHalfOneDarkSplitFn QtoQGDARKSudakov
-newdef QtoQGDARKSudakov:InteractionType DARK
-newdef QtoQGDARKSudakov:Alpha AlphaDARK
-newdef QtoQGDARKSudakov:Cutoff PTCutOff
-newdef QtoQGDARKSudakov:ColourStructure TripletTripletOctet
-newdef QtoQGDARKSudakov:AngularOrdered Yes
-newdef QtoQGDARKSudakov:StrictAO Yes
-newdef QtoQGDARKSudakov:PDFmax 1.9
-
# q -> q gamma
create Herwig::HalfHalfOneSplitFn QtoQGammaSudakov
newdef QtoQGammaSudakov:ColourStructure ChargedChargedNeutral
newdef QtoQGammaSudakov:AngularOrdered Yes
newdef QtoQGammaSudakov:StrictAO Yes
newdef QtoQGammaSudakov:InteractionType QED
newdef QtoQGammaSudakov:Alpha AlphaQED
newdef QtoQGammaSudakov:Cutoff PTCutOff
newdef QtoQGammaSudakov:PDFmax 1.9
# l -> l gamma
cp QtoQGammaSudakov LtoLGammaSudakov
cp PTCutOff LtoLGammaPTCutOff
# Technical parameter to stop evolution.
newdef LtoLGammaPTCutOff:pTmin 0.000001
newdef LtoLGammaSudakov:Cutoff LtoLGammaPTCutOff
# g -> g g
create Herwig::OneOneOneSplitFn GtoGGSudakov
newdef GtoGGSudakov:Alpha AlphaQCDISR
newdef GtoGGSudakov:Cutoff PTCutOff
newdef GtoGGSudakov:PDFmax 2.0
newdef GtoGGSudakov:ColourStructure OctetOctetOctet
newdef GtoGGSudakov:AngularOrdered Yes
newdef GtoGGSudakov:StrictAO Yes
newdef GtoGGSudakov:InteractionType QCD
# W -> W gamma
create Herwig::OneOneOneMassiveSplitFn WtoWGammaSudakov
newdef WtoWGammaSudakov:ColourStructure ChargedChargedNeutral
newdef WtoWGammaSudakov:Cutoff PTCutOff
newdef WtoWGammaSudakov:AngularOrdered Yes
newdef WtoWGammaSudakov:StrictAO Yes
newdef WtoWGammaSudakov:InteractionType QED
newdef WtoWGammaSudakov:Alpha AlphaQED
# g to q qbar
create Herwig::OneHalfHalfSplitFn GtoQQbarSudakov
newdef GtoQQbarSudakov:Alpha AlphaQCDISR
newdef GtoQQbarSudakov:Cutoff PTCutOff
newdef GtoQQbarSudakov:ColourStructure OctetTripletTriplet
newdef GtoQQbarSudakov:AngularOrdered Yes
newdef GtoQQbarSudakov:StrictAO Yes
newdef GtoQQbarSudakov:InteractionType QCD
newdef GtoQQbarSudakov:PDFmax 120.0
# g -> b bbar (separate for efficenicy in ISR)
cp GtoQQbarSudakov GtobbbarSudakov
newdef GtobbbarSudakov:PDFmax 40000.0
# g -> c cbar (separate for efficenicy in ISR)
cp GtoQQbarSudakov GtoccbarSudakov
newdef GtoccbarSudakov:PDFmax 2000.0
# gamma -> q qbar
create Herwig::OneHalfHalfSplitFn GammatoQQbarSudakov
newdef GammatoQQbarSudakov:Alpha AlphaQED
newdef GammatoQQbarSudakov:Cutoff PTCutOff
newdef GammatoQQbarSudakov:ColourStructure NeutralChargedCharged
newdef GammatoQQbarSudakov:AngularOrdered Yes
newdef GammatoQQbarSudakov:StrictAO Yes
newdef GammatoQQbarSudakov:InteractionType QED
newdef GammatoQQbarSudakov:PDFmax 10000.0
# q -> g q
create Herwig::HalfOneHalfSplitFn QtoGQSudakov
newdef QtoGQSudakov:Alpha AlphaQCDISR
newdef QtoGQSudakov:Cutoff PTCutOff
newdef QtoGQSudakov:ColourStructure TripletOctetTriplet
newdef QtoGQSudakov:AngularOrdered Yes
newdef QtoGQSudakov:StrictAO Yes
newdef QtoGQSudakov:InteractionType QCD
# u -> g u (for efficiency)
cp QtoGQSudakov utoGuSudakov
newdef utoGuSudakov:PDFFactor OverOneMinusZ
newdef utoGuSudakov:PDFmax 5.0
# d -> g d (for efficiency)
cp QtoGQSudakov dtoGdSudakov
newdef dtoGdSudakov:PDFFactor OverOneMinusZ
# q -> gamma q
create Herwig::HalfOneHalfSplitFn QtoGammaQSudakov
newdef QtoGammaQSudakov:Alpha AlphaQED
newdef QtoGammaQSudakov:Cutoff PTCutOff
newdef QtoGammaQSudakov:ColourStructure ChargedNeutralCharged
newdef QtoGammaQSudakov:AngularOrdered Yes
newdef QtoGammaQSudakov:StrictAO Yes
newdef QtoGammaQSudakov:InteractionType QED
create Herwig::HalfHalfOneEWSplitFn QtoQWZSudakov
newdef QtoQWZSudakov:InteractionType EW
newdef QtoQWZSudakov:ColourStructure EW
newdef QtoQWZSudakov:Alpha AlphaEW
newdef QtoQWZSudakov:PDFmax 1.9
newdef QtoQWZSudakov:Cutoff PTCutOff
create Herwig::HalfHalfOneEWSplitFn LtoLWZSudakov
newdef LtoLWZSudakov:InteractionType EW
newdef LtoLWZSudakov:ColourStructure EW
newdef LtoLWZSudakov:Alpha AlphaEW
newdef LtoLWZSudakov:PDFmax 1.9
newdef LtoLWZSudakov:Cutoff PTCutOff
create Herwig::HalfHalfZeroEWSplitFn QtoQHSudakov
newdef QtoQHSudakov:InteractionType EW
newdef QtoQHSudakov:ColourStructure EW
newdef QtoQHSudakov:Alpha AlphaEW
newdef QtoQHSudakov:PDFmax 1.9
newdef QtoQHSudakov:Cutoff PTCutOff
cp QtoQHSudakov LtoLHSudakov
create Herwig::OneOneOneEWSplitFn VtoVVSudakov
newdef VtoVVSudakov:InteractionType EW
newdef VtoVVSudakov:ColourStructure EW
newdef VtoVVSudakov:Alpha AlphaQED
newdef VtoVVSudakov:Cutoff PTCutOff
-create Herwig::OneOneOneDarkSplitFn GtoGGDARKSudakov
-newdef GtoGGDARKSudakov:InteractionType DARK
-newdef GtoGGDARKSudakov:ColourStructure OctetOctetOctet
-newdef GtoGGDARKSudakov:AngularOrdered Yes
-newdef GtoGGDARKSudakov:Alpha AlphaDARK
-newdef GtoGGDARKSudakov:StrictAO Yes
-newdef GtoGGDARKSudakov:PDFmax 2.0
-newdef GtoGGDARKSudakov:Cutoff PTCutOff
-
create Herwig::OneOneOneQEDSplitFn GammatoWWSudakov
newdef GammatoWWSudakov:InteractionType QED
newdef GammatoWWSudakov:ColourStructure NeutralChargedCharged
newdef GammatoWWSudakov:Alpha AlphaEW
newdef GammatoWWSudakov:PDFmax 1.9
newdef GammatoWWSudakov:Cutoff PTCutOff
create Herwig::OneOneZeroEWSplitFn VtoVHSudakov
newdef VtoVHSudakov:InteractionType EW
newdef VtoVHSudakov:ColourStructure EW
newdef VtoVHSudakov:Alpha AlphaEW
newdef VtoVHSudakov:PDFmax 1.9
newdef VtoVHSudakov:Cutoff PTCutOff
-create Herwig::OneHalfHalfDarkSplitFn GtoQQbarDARKSudakov
-newdef GtoQQbarDARKSudakov:InteractionType DARK
-newdef GtoQQbarDARKSudakov:ColourStructure OctetTripletTriplet
-newdef GtoQQbarDARKSudakov:AngularOrdered Yes
-newdef GtoQQbarDARKSudakov:StrictAO Yes
-newdef GtoQQbarDARKSudakov:Alpha AlphaDARK
-newdef GtoQQbarDARKSudakov:PDFmax 120.0
-newdef GtoQQbarDARKSudakov:Cutoff PTCutOff
-
create Herwig::ZeroZeroZeroEWSplitFn HtoHHSudakov
newdef HtoHHSudakov:InteractionType EW
newdef HtoHHSudakov:ColourStructure EW
newdef HtoHHSudakov:Alpha AlphaEW
newdef HtoHHSudakov:PDFmax 1.9
newdef HtoHHSudakov:Cutoff PTCutOff
#Use a different Sudakov for FS QCD splittings in order to use a different value of alphaS
cp QtoQGSudakov QtoQGSudakovFSR
cp GtoGGSudakov GtoGGSudakovFSR
cp GtoQQbarSudakov GtoQQbarSudakovFSR
cp GtobbbarSudakov GtobbbarSudakovFSR
cp GtobbbarSudakov GtoccbarSudakovFSR
set QtoQGSudakovFSR:Alpha AlphaQCDFSR
set GtoGGSudakovFSR:Alpha AlphaQCDFSR
set GtoQQbarSudakovFSR:Alpha AlphaQCDFSR
set GtobbbarSudakovFSR:Alpha AlphaQCDFSR
set GtoccbarSudakovFSR:Alpha AlphaQCDFSR
#
# Now add the final splittings
#
do SplittingGenerator:AddFinalSplitting u->u,g; QtoQGSudakovFSR
do SplittingGenerator:AddFinalSplitting d->d,g; QtoQGSudakovFSR
do SplittingGenerator:AddFinalSplitting s->s,g; QtoQGSudakovFSR
do SplittingGenerator:AddFinalSplitting c->c,g; QtoQGSudakovFSR
do SplittingGenerator:AddFinalSplitting b->b,g; QtoQGSudakovFSR
do SplittingGenerator:AddFinalSplitting t->t,g; QtoQGSudakovFSR
#
do SplittingGenerator:AddFinalSplitting g->g,g; GtoGGSudakovFSR
#
do SplittingGenerator:AddFinalSplitting g->u,ubar; GtoQQbarSudakovFSR
do SplittingGenerator:AddFinalSplitting g->d,dbar; GtoQQbarSudakovFSR
do SplittingGenerator:AddFinalSplitting g->s,sbar; GtoQQbarSudakovFSR
do SplittingGenerator:AddFinalSplitting g->c,cbar; GtoccbarSudakovFSR
do SplittingGenerator:AddFinalSplitting g->b,bbar; GtobbbarSudakovFSR
do SplittingGenerator:AddFinalSplitting g->t,tbar; GtoQQbarSudakovFSR
#
do SplittingGenerator:AddFinalSplitting gamma->u,ubar; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->d,dbar; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->s,sbar; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->c,cbar; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->b,bbar; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->t,tbar; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->e-,e+; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->mu-,mu+; GammatoQQbarSudakov
do SplittingGenerator:AddFinalSplitting gamma->tau-,tau+; GammatoQQbarSudakov
#
do SplittingGenerator:AddFinalSplitting u->u,gamma; QtoQGammaSudakov
do SplittingGenerator:AddFinalSplitting d->d,gamma; QtoQGammaSudakov
do SplittingGenerator:AddFinalSplitting s->s,gamma; QtoQGammaSudakov
do SplittingGenerator:AddFinalSplitting c->c,gamma; QtoQGammaSudakov
do SplittingGenerator:AddFinalSplitting b->b,gamma; QtoQGammaSudakov
do SplittingGenerator:AddFinalSplitting t->t,gamma; QtoQGammaSudakov
do SplittingGenerator:AddFinalSplitting e-->e-,gamma; LtoLGammaSudakov
do SplittingGenerator:AddFinalSplitting mu-->mu-,gamma; LtoLGammaSudakov
do SplittingGenerator:AddFinalSplitting tau-->tau-,gamma; LtoLGammaSudakov
do SplittingGenerator:AddFinalSplitting W+->W+,gamma; WtoWGammaSudakov
#
# Now lets add the initial splittings. Remember the form a->b,c; means
# that the current particle b is given and we backward branch to new
# particle a which is initial state and new particle c which is final state
#
do SplittingGenerator:AddInitialSplitting u->u,g; QtoQGSudakov
do SplittingGenerator:AddInitialSplitting d->d,g; QtoQGSudakov
do SplittingGenerator:AddInitialSplitting s->s,g; QtoQGSudakov
do SplittingGenerator:AddInitialSplitting c->c,g; QtoQGSudakov
do SplittingGenerator:AddInitialSplitting b->b,g; QtoQGSudakov
do SplittingGenerator:AddInitialSplitting u->u,gamma; QtoQGammaSudakov
do SplittingGenerator:AddInitialSplitting d->d,gamma; QtoQGammaSudakov
do SplittingGenerator:AddInitialSplitting s->s,gamma; QtoQGammaSudakov
do SplittingGenerator:AddInitialSplitting c->c,gamma; QtoQGammaSudakov
do SplittingGenerator:AddInitialSplitting b->b,gamma; QtoQGammaSudakov
do SplittingGenerator:AddInitialSplitting t->t,gamma; QtoQGammaSudakov
do SplittingGenerator:AddInitialSplitting g->g,g; GtoGGSudakov
#
do SplittingGenerator:AddInitialSplitting g->d,dbar; GtoQQbarSudakov
do SplittingGenerator:AddInitialSplitting g->u,ubar; GtoQQbarSudakov
do SplittingGenerator:AddInitialSplitting g->s,sbar; GtoQQbarSudakov
do SplittingGenerator:AddInitialSplitting g->c,cbar; GtoccbarSudakov
do SplittingGenerator:AddInitialSplitting g->b,bbar; GtobbbarSudakov
#
#do SplittingGenerator:AddInitialSplitting gamma->d,dbar; GammatoQQbarSudakov
#do SplittingGenerator:AddInitialSplitting gamma->u,ubar; GammatoQQbarSudakov
#do SplittingGenerator:AddInitialSplitting gamma->s,sbar; GammatoQQbarSudakov
#do SplittingGenerator:AddInitialSplitting gamma->c,cbar; GammatoQQbarSudakov
#do SplittingGenerator:AddInitialSplitting gamma->b,bbar; GammatoQQbarSudakov
do SplittingGenerator:AddInitialSplitting gamma->e-,e+; GammatoQQbarSudakov
do SplittingGenerator:AddInitialSplitting gamma->mu-,mu+; GammatoQQbarSudakov
do SplittingGenerator:AddInitialSplitting gamma->tau-,tau+; GammatoQQbarSudakov
#
do SplittingGenerator:AddInitialSplitting d->g,d; dtoGdSudakov
do SplittingGenerator:AddInitialSplitting u->g,u; utoGuSudakov
do SplittingGenerator:AddInitialSplitting s->g,s; QtoGQSudakov
do SplittingGenerator:AddInitialSplitting c->g,c; QtoGQSudakov
do SplittingGenerator:AddInitialSplitting b->g,b; QtoGQSudakov
do SplittingGenerator:AddInitialSplitting dbar->g,dbar; dtoGdSudakov
do SplittingGenerator:AddInitialSplitting ubar->g,ubar; utoGuSudakov
do SplittingGenerator:AddInitialSplitting sbar->g,sbar; QtoGQSudakov
do SplittingGenerator:AddInitialSplitting cbar->g,cbar; QtoGQSudakov
do SplittingGenerator:AddInitialSplitting bbar->g,bbar; QtoGQSudakov
#
do SplittingGenerator:AddInitialSplitting d->gamma,d; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting u->gamma,u; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting s->gamma,s; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting c->gamma,c; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting b->gamma,b; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting dbar->gamma,dbar; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting ubar->gamma,ubar; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting sbar->gamma,sbar; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting cbar->gamma,cbar; QtoGammaQSudakov
do SplittingGenerator:AddInitialSplitting bbar->gamma,bbar; QtoGammaQSudakov
#
# Electroweak
#
do SplittingGenerator:AddFinalSplitting u->u,Z0; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting d->d,Z0; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting s->s,Z0; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting c->c,Z0; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting b->b,Z0; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting t->t,Z0; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting u->u,Z0; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting d->d,Z0; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting s->s,Z0; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting c->c,Z0; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting b->b,Z0; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting t->t,Z0; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting u->d,W+; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting c->s,W+; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting d->u,W-; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting s->c,W-; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting u->d,W+; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting c->s,W+; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting d->u,W-; QtoQWZSudakov
do SplittingGenerator:AddInitialSplitting s->c,W-; QtoQWZSudakov
do SplittingGenerator:AddFinalSplitting c->c,h0; QtoQHSudakov
do SplittingGenerator:AddFinalSplitting b->b,h0; QtoQHSudakov
do SplittingGenerator:AddFinalSplitting t->t,h0; QtoQHSudakov
do SplittingGenerator:AddInitialSplitting c->c,h0; QtoQHSudakov
do SplittingGenerator:AddInitialSplitting b->b,h0; QtoQHSudakov
do SplittingGenerator:AddInitialSplitting t->t,h0; QtoQHSudakov
do SplittingGenerator:AddFinalSplitting gamma->W+,W-; GammatoWWSudakov
do SplittingGenerator:AddFinalSplitting Z0->W+,W-; VtoVVSudakov
do SplittingGenerator:AddFinalSplitting W+->W+,gamma; VtoVVSudakov
do SplittingGenerator:AddFinalSplitting W+->W+,Z0; VtoVVSudakov
do SplittingGenerator:AddFinalSplitting W+->W+,h0; VtoVHSudakov
do SplittingGenerator:AddFinalSplitting Z0->Z0,h0; VtoVHSudakov
do SplittingGenerator:AddFinalSplitting h0->h0,h0; HtoHHSudakov
#
# Electroweak l -> l V
#
#do SplittingGenerator:AddFinalSplitting e-->e-,Z0; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting mu-->mu-,Z0; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting tau-->tau-,Z0; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting nu_e->nu_e,Z0; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting nu_mu->nu_mu,Z0; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting nu_tau->nu_tau,Z0; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting e-->nu_e,W-; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting mu-->nu_mu,W-; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting tau-->nu_tau,W-; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting nu_e->e-,W+; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting nu_mu->mu-,W+; LtoLWZSudakov
#do SplittingGenerator:AddFinalSplitting nu_tau->tau-,W+; LtoLWZSudakov

File Metadata

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

Event Timeline